X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;f=source%2Funres%2FMD.F90;h=3d540f5167388d6b06a7c4c17beb171b40d4589c;hb=3ff3c7823d59315586fae4842324638986cfa569;hp=a3451ec450441c10e8be8afa7fc202f4e24f567a;hpb=692e3a1a84c499761a035b29d793a6552dbcd3eb;p=unres4.git diff --git a/source/unres/MD.F90 b/source/unres/MD.F90 index a3451ec..3d540f5 100644 --- a/source/unres/MD.F90 +++ b/source/unres/MD.F90 @@ -1087,6 +1087,7 @@ use comm_gucio use control, only:tcpu use control_data + use minimm, only:minim_dc #ifdef MPI include 'mpif.h' integer :: ierror,ierrcode @@ -1114,17 +1115,17 @@ integer :: count,rstcount !ilen, !el external ilen character(len=50) :: tytul - integer :: maxcount_scale = 20 + integer :: maxcount_scale = 30 !el common /gucio/ cm !el real(kind=8),dimension(6*nres) :: stochforcvec !(MAXRES6) maxres6=6*maxres !el common /stochcalc/ stochforcvec - integer :: itime,icount_scale,itime_scal,i,j,ifac_time - logical :: scale + integer :: itime,icount_scale,itime_scal,i,j,ifac_time,iretcode + logical :: scalel real(kind=8) :: epdrift,tt0,fac_time ! if (.not.allocated(stochforcvec)) allocate(stochforcvec(6*nres)) !(MAXRES6) maxres6=6*maxres - scale=.true. + scalel=.true. icount_scale=0 if (lang.eq.1) then call sddir_precalc @@ -1141,8 +1142,9 @@ #endif endif itime_scal=0 - do while (scale) + do while (scalel) icount_scale=icount_scale+1 +! write(iout,*) "icount_scale",icount_scale,scalel if (icount_scale.gt.maxcount_scale) then write (iout,*) & "ERROR: too many attempts at scaling down the time step. ",& @@ -1204,8 +1206,19 @@ call etotal(potEcomp) ! AL 4/17/17: Reduce the steps if NaNs occurred. if (potEcomp(0).gt.0.99e18 .or. isnan(potEcomp(0)).gt.0) then - d_time=d_time/2 + d_time=d_time/10.0 + if (icount_scale.gt.15) then + write (iout,*) "Tu jest problem",potEcomp(0),d_time +! call gen_rand_conf(1,*335) +! call minim_dc(potEcomp(0),iretcode,100) + +! call zerograd +! call etotal(potEcomp) +! write(iout,*) "needed to repara,",potEcomp + endif cycle +! 335 write(iout,*) "Failed genrand" +! cycle endif ! end change if (large.and. mod(itime,ntwe).eq.0) & @@ -1230,9 +1243,13 @@ call max_accel amax=amax/(itime_scal+1)**2 call predict_edrift(epdrift) +! write(iout,*) "amax=",amax,damax,epdrift,edriftmax,amax/(itime_scal+1) + scalel=.false. +! write (iout,*) "before enter if",scalel,icount_scale if (amax/(itime_scal+1).gt.damax .or. epdrift.gt.edriftmax) then +! write(iout,*) "I enter if" ! Maximum acceleration or maximum predicted energy drift exceeded, rescale the time step - scale=.true. + scalel=.true. ifac_time=dmax1(dlog(amax/damax),dlog(epdrift/edriftmax)) & /dlog(2.0d0)+1 itime_scal=itime_scal+ifac_time @@ -1348,7 +1365,6 @@ endif #endif endif - scale=.false. endif enddo ! Calculate the kinetic and the total energy and the kinetic temperature @@ -2320,6 +2336,7 @@ character(len=16) :: form integer :: IERROR,ERRCODE #endif + integer :: iranmin,itrial,itmp ! include 'COMMON.SETUP' ! include 'COMMON.CONTROL' ! include 'COMMON.VAR' @@ -2535,6 +2552,85 @@ if(me.eq.king.or..not.out1file) & write(iout,*) 'SUMSL return code is',iretcode,' eval ',nfun endif + if(iranconf.ne.0) then +!c 8/22/17 AL Loop to produce a low-energy random conformation + DO iranmin=1,40 + if (overlapsc) then + if(me.eq.king.or..not.out1file) & + write (iout,*) 'Calling OVERLAP_SC' + call overlap_sc(fail) + endif !endif overlap + + if (searchsc) then + call sc_move(2,nres-1,10,1d10,nft_sc,etot) + print *,'SC_move',nft_sc,etot + if(me.eq.king.or..not.out1file) & + write(iout,*) 'SC_move',nft_sc,etot + endif + + if(dccart)then + print *, 'Calling MINIM_DC' + call minim_dc(etot,iretcode,nfun) + call int_from_cart1(.false.) + else + call geom_to_var(nvar,varia) + print *,'Calling MINIMIZE.' + call minimize(etot,varia,iretcode,nfun) + call var_to_geom(nvar,varia) + endif + if(me.eq.king.or..not.out1file) & + write(iout,*) 'SUMSL return code is',iretcode,' eval ',nfun + + if (isnan(etot) .or. etot.gt.4.0d6) then + write (iout,*) "Energy too large",etot, & + " trying another random conformation" + do itrial=1,100 + itmp=1 + call gen_rand_conf(itmp,*30) + goto 40 + 30 write (iout,*) 'Failed to generate random conformation', & + ', itrial=',itrial + write (*,*) 'Processor:',me, & + ' Failed to generate random conformation',& + ' itrial=',itrial + call intout +#ifdef AIX + call flush_(iout) +#else + call flush(iout) +#endif + enddo + write (iout,'(a,i3,a)') 'Processor:',me, & + ' error in generating random conformation.' + write (*,'(a,i3,a)') 'Processor:',me, & + ' error in generating random conformation.' + call flush(iout) +#ifdef MPI +! call MPI_Abort(mpi_comm_world,error_msg,ierrcode) + call MPI_Abort(MPI_COMM_WORLD,IERROR,ERRCODE) +#else + stop +#endif + 40 continue + else + goto 44 + endif + ENDDO + + write (iout,'(a,i3,a)') 'Processor:',me, & + ' failed to generate a low-energy random conformation.' + write (*,'(a,i3,a,f10.3)') 'Processor:',me, & + ' failed to generate a low-energy random conformation.',etot + call flush(iout) + call intout +#ifdef MPI +! call MPI_Abort(mpi_comm_world,error_msg,ierrcode) + call MPI_Abort(MPI_COMM_WORLD,IERROR,ERRCODE) +#else + stop +#endif + 44 continue + endif endif call chainbuild_cart call kinetic(EK)