X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;f=source%2Funres%2Fsrc_MD-M-newcorr%2Fminim_jlee.F;fp=source%2Funres%2Fsrc_MD-M-newcorr%2Fminim_jlee.F;h=d83b15bda91f9bb9b3e85ee2ff7125b2ef902b4e;hb=7308760ff07636ef6b1ee28d8c3a67a23c14b34b;hp=0000000000000000000000000000000000000000;hpb=9a54ab407f6d0d9d564d52763b3e2136450b9ffc;p=unres.git diff --git a/source/unres/src_MD-M-newcorr/minim_jlee.F b/source/unres/src_MD-M-newcorr/minim_jlee.F new file mode 100644 index 0000000..d83b15b --- /dev/null +++ b/source/unres/src_MD-M-newcorr/minim_jlee.F @@ -0,0 +1,435 @@ + subroutine minim_jlee +c controls minimization and sorting routines + implicit real*8 (a-h,o-z) + include 'DIMENSIONS' + parameter (liv=60,lv=(77+maxvar*(maxvar+17)/2)) + include 'COMMON.VAR' + include 'COMMON.IOUNITS' + include 'COMMON.MINIM' + include 'COMMON.CONTROL' + include 'mpif.h' + external func,gradient,fdum + real ran1,ran2,ran3 + include 'COMMON.SETUP' + include 'COMMON.GEO' + include 'COMMON.FFIELD' + include 'COMMON.SBRIDGE' + include 'COMMON.DISTFIT' + include 'COMMON.CHAIN' + dimension muster(mpi_status_size) + dimension var(maxvar),erg(mxch*(mxch+1)/2+1) + dimension var2(maxvar) + integer iffr(maxres),ihpbt(maxdim),jhpbt(maxdim) + double precision d(maxvar),v(1:lv+1),garbage(maxvar) + double precision energia(0:n_ene),time0s,time1s + dimension indx(9),info(12) + dimension iv(liv) + dimension idum(1),rdum(1) + dimension icont(2,maxcont) + logical check_var,fail + integer iloop(2) + common /przechowalnia/ v + data rad /1.745329252d-2/ +c receive # of start +! print *,'Processor',me,' calling MINIM_JLEE maxfun',maxfun, +! & ' maxmin',maxmin,' tolf',tolf,' rtolf',rtolf + nhpb0=nhpb + 10 continue + time0s=MPI_WTIME() +c print *, 'MINIM_JLEE: ',me,' is waiting' + call mpi_recv(info,12,mpi_integer,king,idint,CG_COMM, + * muster,ierr) + time1s=MPI_WTIME() + write (iout,'(a12,f10.4,a4)')'Waiting for ',time1s-time0s,' sec' + call flush(iout) + n=info(1) +c print *, 'MINIM_JLEE: ',me,' received: ',n + +crc if (ierr.ne.0) go to 100 +c if # = 0, return + if (n.eq.0) then + write (iout,*) 'Finishing minim_jlee - signal',n,' from master' + call flush(iout) + return + endif + + nfun=0 + IF (n.lt.0) THEN + call mpi_recv(var,nvar,mpi_double_precision, + * king,idreal,CG_COMM,muster,ierr) + call mpi_recv(iffr,nres,mpi_integer, + * king,idint,CG_COMM,muster,ierr) + call mpi_recv(var2,nvar,mpi_double_precision, + * king,idreal,CG_COMM,muster,ierr) + ELSE +c receive initial values of variables + call mpi_recv(var,nvar,mpi_double_precision, + * king,idreal,CG_COMM,muster,ierr) +crc if (ierr.ne.0) go to 100 + ENDIF + + if(vdisulf.and.info(2).ne.-1) then + if(info(4).ne.0)then + call mpi_recv(ihpbt,info(4),mpi_integer, + * king,idint,CG_COMM,muster,ierr) + call mpi_recv(jhpbt,info(4),mpi_integer, + * king,idint,CG_COMM,muster,ierr) + endif + endif + + IF (n.lt.0) THEN + n=-n + nhpb=nhpb0 + link_start=1 + link_end=nhpb + call init_int_table + call contact_cp(var,var2,iffr,nfun,n) + ENDIF + + if(vdisulf.and.info(2).ne.-1) then + nss=0 + if(info(4).ne.0)then +cd write(iout,*) 'SS=',info(4),'N=',info(1),'IT=',info(2) + call var_to_geom(nvar,var) + call chainbuild + do i=1,info(4) + if (dist(ihpbt(i),jhpbt(i)).lt.7.0) then + nss=nss+1 + ihpb(nss)=ihpbt(i) + jhpb(nss)=jhpbt(i) +cd write(iout,*) 'SS mv=',info(3), +cd & ihpb(nss)-nres,jhpb(nss)-nres, +cd & dist(ihpb(nss),jhpb(nss)) + dhpb(nss)=dbr + forcon(nss)=fbr + else +cd write(iout,*) 'rm SS mv=',info(3), +cd & ihpbt(i)-nres,jhpbt(i)-nres,dist(ihpbt(i),jhpbt(i)) + endif + enddo + endif + nhpb=nss + link_start=1 + link_end=nhpb + call init_int_table + endif + + if (info(3).eq.14) then + write(iout,*) 'calling local_move',info(7),info(8) + call local_move_init(.false.) + call var_to_geom(nvar,var) + call local_move(info(7),info(8),20d0,50d0) + call geom_to_var(nvar,var) + endif + + + if (info(3).eq.16) then + write(iout,*) 'calling beta_slide',info(7),info(8), + & info(10), info(11), info(12) + call var_to_geom(nvar,var) + call beta_slide(info(7),info(8),info(10),info(11),info(12) + & ,nfun,n) + call geom_to_var(nvar,var) + endif + + + if (info(3).eq.17) then + write(iout,*) 'calling beta_zip',info(7),info(8) + call var_to_geom(nvar,var) + call beta_zip(info(7),info(8),nfun,n) + call geom_to_var(nvar,var) + endif + + +crc overlap test + + if (overlapsc) then + + call var_to_geom(nvar,var) + call chainbuild + call etotal(energia(0)) + nfun=nfun+1 + if (energia(1).eq.1.0d20) then + info(3)=-info(3) + write (iout,'(a,1pe14.5)')'#OVERLAP evdw=1d20',energia(1) + call overlap_sc(fail) + if(.not.fail) then + call geom_to_var(nvar,var) + call etotal(energia(0)) + nfun=nfun+1 + write (iout,'(a,1pe14.5)')'#OVERLAP evdw after',energia(1) + else + v(10)=1.0d20 + iv(1)=-1 + goto 201 + endif + endif + endif + + if (searchsc) then + call var_to_geom(nvar,var) + call sc_move(2,nres-1,1,10d0,nft_sc,etot) + call geom_to_var(nvar,var) +cd write(iout,*) 'sc_move',nft_sc,etot + endif + + if (check_var(var,info)) then + v(10)=1.0d21 + iv(1)=6 + goto 201 + endif + + +crc + +! write (iout,*) 'MINIM_JLEE: Processor',me,' nvar',nvar +! write (iout,'(8f10.4)') (var(i),i=1,nvar) +! write (*,*) 'MINIM_JLEE: Processor',me,' received nvar',nvar +! write (*,'(8f10.4)') (var(i),i=1,nvar) + + do i=1,nvar + garbage(i)=var(i) + enddo + + call deflt(2,iv,liv,lv,v) +* 12 means fresh start, dont call deflt + iv(1)=12 +* max num of fun calls + if (maxfun.eq.0) maxfun=500 + iv(17)=maxfun +* max num of iterations + if (maxmin.eq.0) maxmin=1000 + iv(18)=maxmin +* controls output + iv(19)=2 +* selects output unit +cd iv(21)=iout + iv(21)=0 +* 1 means to print out result + iv(22)=0 +cd iv(22)=1 +* 1 means to print out summary stats + iv(23)=0 +* 1 means to print initial x and d + iv(24)=0 + +c if(me.eq.3.and.n.eq.255) then +c print *,' CHUJ: stoi' +c iv(21)=6 +c iv(22)=1 +c iv(23)=1 +c iv(24)=1 +c endif + +* min val for v(radfac) default is 0.1 + v(24)=0.1D0 +* max val for v(radfac) default is 4.0 + v(25)=2.0D0 +c v(25)=4.0D0 +* check false conv if (act fnctn decrease) .lt. v(26)*(exp decrease) +* the sumsl default is 0.1 + v(26)=0.1D0 +* false conv if (act fnctn decrease) .lt. v(34) +* the sumsl default is 100*machep + v(34)=v(34)/100.0D0 +* absolute convergence + if (tolf.eq.0.0D0) tolf=1.0D-4 + v(31)=tolf +* relative convergence + if (rtolf.eq.0.0D0) rtolf=1.0D-4 + v(32)=rtolf +* controls initial step size + v(35)=1.0D-1 +* large vals of d correspond to small components of step + do i=1,nphi + d(i)=1.0D-1 + enddo + do i=nphi+1,nvar + d(i)=1.0D-1 + enddo +c minimize energy +! write (iout,*) 'Processor',me,' nvar',nvar +! write (iout,*) 'Variables BEFORE minimization:' +! write (iout,'(8f10.4)') (rad2deg*var(i),i=1,nvar) + +c print *, 'MINIM_JLEE: ',me,' before SUMSL ' + + call func(nvar,var,nf,eee,idum,rdum,fdum) + nfun=nfun+1 + if(eee.ge.1.0d20) then +c print *,'MINIM_JLEE: ',me,' CHUJ NASTAPIL' +c print *,' energy before SUMSL =',eee +c print *,' aborting local minimization' + iv(1)=-1 + v(10)=eee + go to 201 + endif + +ct time0s=MPI_WTIME() + call sumsl(nvar,d,var,func,gradient,iv,liv,lv,v,idum,rdum,fdum) +ct write(iout,*) 'sumsl time=',MPI_WTIME()-time0s,iv(7),v(10) +c print *, 'MINIM_JLEE: ',me,' after SUMSL ' + +c find which conformation was returned from sumsl + nfun=nfun+iv(7) +! print *,'Processor',me,' iv(17)',iv(17),' iv(18)',iv(18),' nf',nf, +! & ' retcode',iv(1),' energy',v(10),' tolf',v(31),' rtolf',v(32) +c if (iv(1).ne.4 .or. nf.le.1) then +c write (*,*) 'Processor',me,' something bad in SUMSL',iv(1),nf +c write (*,*) 'Initial Variables' +c write (*,'(8f10.4)') (rad2deg*garbage(i),i=1,nvar) +c write (*,*) 'Variables' +c write (*,'(8f10.4)') (rad2deg*var(i),i=1,nvar) +c write (*,*) 'Vector d' +c write (*,'(8f10.4)') (d(i),i=1,nvar) +c write (iout,*) 'Processor',me,' something bad in SUMSL', +c & iv(1),nf +c write (iout,*) 'Initial Variables' +c write (iout,'(8f10.4)') (rad2deg*garbage(i),i=1,nvar) +c write (iout,*) 'Variables' +c write (iout,'(8f10.4)') (rad2deg*var(i),i=1,nvar) +c write (iout,*) 'Vector d' +c write (iout,'(8f10.4)') (d(i),i=1,nvar) +c endif +c if (nf.lt.iv(6)-1) then +c recalculate intra- and interchain energies +c call func(nvar,var,nf,v(10),iv,v,fdum) +c else if (nf.eq.iv(6)-1) then +c regenerate conformation +c call var_to_geom(nvar,var) +c call chainbuild +c endif +c change origin and axes to standard ECEPP format +c call var_to_geom(nvar,var) +! write (iout,*) 'MINIM_JLEE after minim: Processor',me,' nvar',nvar +! write (iout,'(8f10.4)') (var(i),i=1,nvar) +! write (iout,*) 'Energy:',v(10) +c send back output +c print *, 'MINIM_JLEE: ',me,' minimized: ',n + 201 continue + indx(1)=n +c return code: 6-gradient 9-number of ftn evaluation, etc + indx(2)=iv(1) +c total # of ftn evaluations (for iwf=0, it includes all minimizations). + indx(3)=nfun + indx(4)=info(2) + indx(5)=info(3) + indx(6)=nss + indx(7)=info(5) + indx(8)=info(6) + indx(9)=info(9) + call mpi_send(indx,9,mpi_integer,king,idint,CG_COMM, + * ierr) +c send back energies +c al & cc +c calculate contact order +#ifdef CO_BIAS + call contact(.false.,ncont,icont,co) + erg(1)=v(10)-1.0d2*co +#else + erg(1)=v(10) +#endif + j=1 + call mpi_send(erg,j,mpi_double_precision,king,idreal, + * CG_COMM,ierr) +#ifdef CO_BIAS + call mpi_send(co,j,mpi_double_precision,king,idreal, + * CG_COMM,ierr) +#endif +c send back values of variables + call mpi_send(var,nvar,mpi_double_precision, + * king,idreal,CG_COMM,ierr) +! print * , 'MINIM_JLEE: Processor',me,' send erg and var ' + + if(vdisulf.and.info(2).ne.-1.and.nss.ne.0) then +cd call intout +cd call chainbuild +cd call etotal(energia(0)) +cd etot=energia(0) +cd call enerprint(energia(0)) + call mpi_send(ihpb,nss,mpi_integer, + * king,idint,CG_COMM,ierr) + call mpi_send(jhpb,nss,mpi_integer, + * king,idint,CG_COMM,ierr) + endif + + go to 10 + 100 print *, ' error in receiving message from emperor', me + call mpi_abort(mpi_comm_world,ierror,ierrcode) + return + 200 print *, ' error in sending message to emperor' + call mpi_abort(mpi_comm_world,ierror,ierrcode) + return + 300 print *, ' error in communicating with emperor' + call mpi_abort(mpi_comm_world,ierror,ierrcode) + return + 956 format (' initial energy could not be calculated',41x) + 957 format (80x) + 965 format (' convergence code ',i2,' # of function calls ', + * i4,' # of gradient calls ',i4,10x) + 975 format (' energy ',1p,e12.4,' scaled gradient ',e11.3,32x) + end + + logical function check_var(var,info) + implicit real*8 (a-h,o-z) + include 'DIMENSIONS' + include 'COMMON.VAR' + include 'COMMON.IOUNITS' + include 'COMMON.GEO' + include 'COMMON.SETUP' + dimension var(maxvar) + dimension info(3) +C AL ------- + check_var=.false. + do i=nphi+ntheta+1,nphi+ntheta+nside +! Check the side chain "valence" angles alpha + if (var(i).lt.1.0d-7) then + write (iout,*) 'CHUJ NASTAPIL ABSOLUTNY!!!!!!!!!!!!' + write (iout,*) 'Processor',me,'received bad variables!!!!' + write (iout,*) 'Variables' + write (iout,'(8f10.4)') (rad2deg*var(j),j=1,nvar) + write (iout,*) 'Continuing calculations at this point', + & ' could destroy the results obtained so far... ABORTING!!!!!!' + write (iout,'(a19,i5,f10.4,a4,2i4,a3,i3)') + & 'valence angle alpha',i-nphi-ntheta,var(i), + & 'n it',info(1),info(2),'mv ',info(3) + write (*,*) 'CHUJ NASTAPIL ABSOLUTNY!!!!!!!!!!!!' + write (*,*) 'Processor',me,'received bad variables!!!!' + write (*,*) 'Variables' + write (*,'(8f10.4)') (rad2deg*var(j),j=1,nvar) + write (*,*) 'Continuing calculations at this point', + & ' could destroy the results obtained so far... ABORTING!!!!!!' + write (*,'(a19,i5,f10.4,a4,2i4,a3,i3)') + & 'valence angle alpha',i-nphi-ntheta,var(i), + & 'n it',info(1),info(2),'mv ',info(3) + check_var=.true. + return + endif + enddo +! Check the backbone "valence" angles theta + do i=nphi+1,nphi+ntheta + if (var(i).lt.1.0d-7) then + write (iout,*) 'CHUJ NASTAPIL ABSOLUTNY!!!!!!!!!!!!' + write (iout,*) 'Processor',me,'received bad variables!!!!' + write (iout,*) 'Variables' + write (iout,'(8f10.4)') (rad2deg*var(j),j=1,nvar) + write (iout,*) 'Continuing calculations at this point', + & ' could destroy the results obtained so far... ABORTING!!!!!!' + write (iout,'(a19,i5,f10.4,a4,2i4,a3,i3)') + & 'valence angle theta',i-nphi,var(i), + & 'n it',info(1),info(2),'mv ',info(3) + write (*,*) 'CHUJ NASTAPIL ABSOLUTNY!!!!!!!!!!!!' + write (*,*) 'Processor',me,'received bad variables!!!!' + write (*,*) 'Variables' + write (*,'(8f10.4)') (rad2deg*var(j),j=1,nvar) + write (*,*) 'Continuing calculations at this point', + & ' could destroy the results obtained so far... ABORTING!!!!!!' + write (*,'(a19,i5,f10.4,a4,2i4,a3,i3)') + & 'valence angle theta',i-nphi,var(i), + & 'n it',info(1),info(2),'mv ',info(3) + check_var=.true. + return + endif + enddo + return + end