X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;f=source%2Funres%2Fsrc-HCD-5D%2Funres.F;h=978fa59f764dea9eb45d3e5537505280d2199cb6;hb=f9783c27de185ff3e1d6e874a37a1a8c4d3c93a3;hp=76bd2808cdc71755c87ba8246389b1fb335af0c0;hpb=48ae9e01d2dd6571fa2cca6c704dc04f86e5fd7b;p=unres.git diff --git a/source/unres/src-HCD-5D/unres.F b/source/unres/src-HCD-5D/unres.F index 76bd280..978fa59 100644 --- a/source/unres/src-HCD-5D/unres.F +++ b/source/unres/src-HCD-5D/unres.F @@ -20,7 +20,7 @@ CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC include 'COMMON.GEO' include 'COMMON.HEADER' include 'COMMON.CONTROL' - include 'COMMON.CONTACTS' +c include 'COMMON.CONTACTS' include 'COMMON.CHAIN' include 'COMMON.VAR' include 'COMMON.IOUNITS' @@ -59,7 +59,8 @@ c call memmon_print_usage() C Read force field parameters and job setup data call readrtns C -c write (iout,*) "After readrtns" + write (iout,*) "After readrtns" + call flush(iout) call cartprint call intout if (me.eq.king .or. .not. out1file) then @@ -73,8 +74,8 @@ c write (iout,*) "After readrtns" call flush(iout) C if (modecalc.eq.-2) then - call test - stop +c call test +c stop else if (modecalc.eq.-1) then write(iout,*) "call check_sc_map next" call check_bond @@ -86,9 +87,14 @@ C Fine-grain slaves just do energy and gradient components. call ergastulum ! slave workhouse in Latin else #endif + if (indpdb.eq.0 .and. .not.read_cart) call chainbuild + if (indpdb.ne.0 .or. read_cart) then + dc(1,0)=c(1,1) + dc(2,0)=c(2,1) + dc(3,0)=c(3,1) + endif if (modecalc.eq.0) then write (iout,*) "Calling exec_eeval_or_minim" - call cartprint call exec_eeval_or_minim else if (modecalc.eq.1) then call exec_regularize @@ -199,7 +205,7 @@ c--------------------------------------------------------------------------- include 'COMMON.GEO' include 'COMMON.HEADER' include 'COMMON.CONTROL' - include 'COMMON.CONTACTS' +c include 'COMMON.CONTACTS' include 'COMMON.CHAIN' include 'COMMON.VAR' include 'COMMON.IOUNITS' @@ -217,9 +223,14 @@ c--------------------------------------------------------------------------- double precision rms,frac,frac_nn,co double precision varia(maxvar) double precision time00,time1,time_ene,evals +#ifdef LBFGS + character*9 status + integer niter + common /lbfgstat/ status,niter,nfun +#endif integer ilen - if (indpdb.eq.0) call chainbuild - if (indpdb.ne.0) then + if (indpdb.eq.0 .and. .not.read_cart) call chainbuild + if (indpdb.ne.0 .or. read_cart) then dc(1,0)=c(1,1) dc(2,0)=c(2,1) dc(3,0)=c(3,1) @@ -269,10 +280,18 @@ c print *,'after hairpin' c print *,'after secondary' if (minim) then crc overlap test + if (indpdb.ne.0 .and. .not.dccart) then + call bond_regular + call chainbuild_extconf + call etotal(energy(0)) + write (iout,*) "After bond regularization" + call enerprint(energy(0)) + endif + if (overlapsc) then - print *, 'Calling OVERLAP_SC' + write (iout,*) 'Calling OVERLAP_SC' call overlap_sc(fail) - print *,"After overlap_sc" + write (iout,*) "After overlap_sc" endif if (searchsc) then @@ -290,12 +309,8 @@ crc overlap test #endif call minim_dc(etot,iretcode,nfun) else - if (indpdb.ne.0) then - call bond_regular - call chainbuild_extconf - endif call geom_to_var(nvar,varia) - print *,'Calling MINIMIZE.' +c print *,'Calling MINIMIZE.' #ifdef MPI time1=MPI_WTIME() #else @@ -303,7 +318,11 @@ crc overlap test #endif call minimize(etot,varia,iretcode,nfun) endif +#ifdef LBFGS + print *,'LBFGS return code is',status,' eval ',nfun +#else print *,'SUMSL return code is',iretcode,' eval ',nfun +#endif #ifdef MPI evals=nfun/(MPI_WTIME()-time1) #else @@ -318,22 +337,23 @@ crc overlap test call etotal(energy(0)) etot = energy(0) call enerprint(energy(0)) - call intout - if (out_int) call briefout(0,etot) - if (out_cart) then - cartname=prefix(:ilen(prefix))//'.x' - potE=etot - call cartoutx(0.0d0) - endif - if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.) +#ifdef LBFGS + write (iout,'(a,a9)') 'LBFGS return code:',status + write (iout,'(a,i20)') '# of energy evaluations:',nfun+1 + write (iout,'(a,f16.3)')'# of energy evaluations/sec:',evals +#else write (iout,'(a,i3)') 'SUMSL return code:',iretcode write (iout,'(a,i20)') '# of energy evaluations:',nfun+1 write (iout,'(a,f16.3)')'# of energy evaluations/sec:',evals - else - print *,'refstr=',refstr - if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.) - call briefout(0,etot) +#endif + endif + if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.) + if (out_int) call briefout(0,etot) + if (out_cart) then + cartname=prefix(:ilen(prefix))//'.x' + potE=etot + call cartoutx(0.0d0) endif if (outpdb) call pdbout(etot,titel(:50),ipdb) if (outmol2) call mol2out(etot,titel(:32)) @@ -353,7 +373,7 @@ c--------------------------------------------------------------------------- include 'COMMON.GEO' include 'COMMON.HEADER' include 'COMMON.CONTROL' - include 'COMMON.CONTACTS' +c include 'COMMON.CONTACTS' include 'COMMON.CHAIN' include 'COMMON.VAR' include 'COMMON.IOUNITS' @@ -364,6 +384,11 @@ c--------------------------------------------------------------------------- double precision energy(0:n_ene) double precision etot,rms,frac,frac_nn,co integer iretcode +#ifdef LBFGS + character*9 status + integer niter,nfun + common /lbfgstat/ status,niter,nfun +#endif call gen_dist_constr call sc_conf @@ -378,7 +403,11 @@ c--------------------------------------------------------------------------- if (outpdb) call pdbout(etot,titel(:50),ipdb) if (outmol2) call mol2out(etot,titel(:32)) if (refstr) call rms_nac_nnc(rms,frac,frac_nn,co,.true.) +#ifdef LBFGS + write (iout,'(a,a9)') 'LBFGS return code:',status +#else write (iout,'(a,i3)') 'SUMSL return code:',iretcode +#endif return end c--------------------------------------------------------------------------- @@ -440,7 +469,7 @@ c--------------------------------------------------------------------------- include 'COMMON.GEO' include 'COMMON.HEADER' include 'COMMON.CONTROL' - include 'COMMON.CONTACTS' +c include 'COMMON.CONTACTS' include 'COMMON.CHAIN' include 'COMMON.VAR' include 'COMMON.IOUNITS' @@ -739,7 +768,7 @@ c--------------------------------------------------------------------------- include 'COMMON.GEO' include 'COMMON.HEADER' include 'COMMON.CONTROL' - include 'COMMON.CONTACTS' +c include 'COMMON.CONTACTS' include 'COMMON.CHAIN' include 'COMMON.VAR' include 'COMMON.IOUNITS' @@ -825,7 +854,8 @@ c enddo goto (10,20,30) icheckgrad 10 call check_ecartint return - 20 call check_cartgrad + 20 write (iout,*) + & "Checking the gradient of Cartesian coordinates disabled." return 30 call check_eint return