hdiff output

r32919/geopt.f 2017-07-04 12:30:19.004672169 +0100 r32918/geopt.f 2017-07-04 12:30:19.812683011 +0100
194:                  CALL CUDA_BFGSTS_WRAPPER(NSTEPS,Q,ENERGY,MFLAG,RMS,EVALMIN,VECS,ITDONE)194:                  CALL CUDA_BFGSTS_WRAPPER(NSTEPS,Q,ENERGY,MFLAG,RMS,EVALMIN,VECS,ITDONE)
195:              ELSE195:              ELSE
196:                  CALL BFGSTS(NSTEPS,Q,ENERGY,VNEW,MFLAG,RMS,EVALMIN,EVALMAX,VECS,ITDONE,.TRUE.,.TRUE.)196:                  CALL BFGSTS(NSTEPS,Q,ENERGY,VNEW,MFLAG,RMS,EVALMIN,EVALMAX,VECS,ITDONE,.TRUE.,.TRUE.)
197:              END IF197:              END IF
198:              IF (MFLAG) THEN198:              IF (MFLAG) THEN
199:                  WRITE(*,'(1X,A,I6)') 'Converged to TS (number of iterations):     ',ITDONE199:                  WRITE(*,'(1X,A,I6)') 'Converged to TS (number of iterations):     ',ITDONE
200:                  write(*,*) "Energy:", ENERGY200:                  write(*,*) "Energy:", ENERGY
201:              ELSE201:              ELSE
202:                  WRITE(*,'(1X,A,I6)') 'DID NOT converge to TS (number of iterations):     ',ITDONE202:                  WRITE(*,'(1X,A,I6)') 'DID NOT converge to TS (number of iterations):     ',ITDONE
203:              END IF203:              END IF
204:              WRITE(*,'(A,F20.10)') ' geopt> Energy of TS is:                    ',ENERGY204:              WRITE(*,'(A,F20.10)') 'Energy of TS is:                    ',ENERGY
 205: 
205:           ENDIF206:           ENDIF
206: 207: 
207: C208: C
208: C  We may now have the energy and gradient, so the first energy and gradient call in209: C  We may now have the energy and gradient, so the first energy and gradient call in
209: C  mylbfgs for BFGSSTEP may be unnecessary210: C  mylbfgs for BFGSSTEP may be unnecessary
210: C211: C
211:           IF (BFGSSTEP) THEN212:           IF (BFGSSTEP) THEN
212: C213: C
213: C  Switch to appropriate minimization after pushoff.214: C  Switch to appropriate minimization after pushoff.
214: C215: C


r32919/golden.f90 2017-07-04 12:30:19.268675732 +0100 r32918/golden.f90 2017-07-04 12:30:20.072686486 +0100
 24: USE KEY, ONLY : PUSHOPTMAX, PUSHOPTCONV 24: USE KEY, ONLY : PUSHOPTMAX, PUSHOPTCONV
 25: IMPLICIT NONE 25: IMPLICIT NONE
 26:  26: 
 27: DOUBLE PRECISION, PARAMETER :: GRATIO=1.618033989D0 27: DOUBLE PRECISION, PARAMETER :: GRATIO=1.618033989D0
 28: DOUBLE PRECISION ASTEP, BSTEP, CSTEP, DSTEP, RMS, VECS(NOPT), QINIT(NOPT), ENERGYC, ENERGYD 28: DOUBLE PRECISION ASTEP, BSTEP, CSTEP, DSTEP, RMS, VECS(NOPT), QINIT(NOPT), ENERGYC, ENERGYD
 29: DOUBLE PRECISION QC(NOPT), QD(NOPT), VNEW(NOPT), ETS, POPT, PENERGY, PINIT 29: DOUBLE PRECISION QC(NOPT), QD(NOPT), VNEW(NOPT), ETS, POPT, PENERGY, PINIT
 30: INTEGER J1 30: INTEGER J1
 31: ! 31: !
 32: ! Initially astep=0 the ts at QINIT(J1), bstep=PUSHOFF, corresponding to current Q(J1) 32: ! Initially astep=0 the ts at QINIT(J1), bstep=PUSHOFF, corresponding to current Q(J1)
 33: ! 33: !
  34: PRINT '(A,6G20.10)','VECS(1:6)=',VECS(1:6)
 34: ASTEP=0.0D0 35: ASTEP=0.0D0
 35: BSTEP=PINIT 36: BSTEP=PINIT
 36: CSTEP = BSTEP - (BSTEP - ASTEP) / GRATIO 37: CSTEP = BSTEP - (BSTEP - ASTEP) / GRATIO
 37: DSTEP = ASTEP + (BSTEP - ASTEP) / GRATIO 38: DSTEP = ASTEP + (BSTEP - ASTEP) / GRATIO
 38: IF (DEBUG) WRITE(6,'(A,4G20.10,A)') 'path> golden a, b, c, d=',ASTEP,BSTEP,CSTEP,DSTEP,' initially' 39: IF (DEBUG) WRITE(6,'(A,4G20.10,A)') 'path> golden a, b, c, d=',ASTEP,BSTEP,CSTEP,DSTEP,' initially'
 39: DO J1=1,PUSHOPTMAX 40: DO J1=1,PUSHOPTMAX
 40:    QC(1:NOPT)=QINIT(1:NOPT)+CSTEP*VECS(1:NOPT) 41:    QC(1:NOPT)=QINIT(1:NOPT)+CSTEP*VECS(1:NOPT)
 41:    QD(1:NOPT)=QINIT(1:NOPT)+DSTEP*VECS(1:NOPT) 42:    QD(1:NOPT)=QINIT(1:NOPT)+DSTEP*VECS(1:NOPT)
 42:    CALL POTENTIAL(QC,ENERGYC,VNEW,.FALSE.,.FALSE.,RMS,.FALSE.,.FALSE.) 43:    CALL POTENTIAL(QC,ENERGYC,VNEW,.FALSE.,.FALSE.,RMS,.FALSE.,.FALSE.)
 43:    CALL POTENTIAL(QD,ENERGYD,VNEW,.FALSE.,.FALSE.,RMS,.FALSE.,.FALSE.) 44:    CALL POTENTIAL(QD,ENERGYD,VNEW,.FALSE.,.FALSE.,RMS,.FALSE.,.FALSE.)
 44:    IF (DEBUG) WRITE(6,'(A,4G25.15)') 'path> energy at ts, C and D and diff=', & 45:    IF (DEBUG) WRITE(6,'(A,4G25.15)') 'path> energy at ts, C and D and diff=', &
 45:      &                    ETS,ENERGYC,ENERGYD,MIN(ENERGYC,ENERGYD)-ETS 46:      &                    ETS,ENERGYC,ENERGYD,MIN(ENERGYC,ENERGYD)-ETS
 46:    IF (ENERGYC.LT.ENERGYD) THEN 47:    IF (ENERGYC.LT.ENERGYD) THEN
 47:       BSTEP=DSTEP 48:       BSTEP=DSTEP
 48:    ELSE 49:    ELSE
 49:       ASTEP=CSTEP 50:       ASTEP=CSTEP
 50:    ENDIF 51:    ENDIF
 51:    CSTEP = BSTEP - (BSTEP - ASTEP) / GRATIO 52:    CSTEP = BSTEP - (BSTEP - ASTEP) / GRATIO
 52:    DSTEP = ASTEP + (BSTEP - ASTEP) / GRATIO 53:    DSTEP = ASTEP + (BSTEP - ASTEP) / GRATIO
 53: !  IF (DEBUG) WRITE(6,'(A,4G20.10)') 'path> golden a, b, c, d=',ASTEP,BSTEP,CSTEP,DSTEP 54:    IF (DEBUG) WRITE(6,'(A,4G20.10)') 'path> golden a, b, c, d=',ASTEP,BSTEP,CSTEP,DSTEP
 54:    IF (ABS(ASTEP-BSTEP).LT.PUSHOPTCONV) THEN 55:    IF (ABS(ASTEP-BSTEP).LT.PUSHOPTCONV) THEN
 55:       ASTEP=(ASTEP+BSTEP)/2.0D0 56:       ASTEP=(ASTEP+BSTEP)/2.0D0
  57: !     Q(1:NOPT)=QINIT(1:NOPT)+ASTEP*VECS(1:NOPT)
  58: !     STEP(1:NOPT)=ASTEP*VECS(1:NOPT)
 56:       EXIT 59:       EXIT
 57:    ENDIF 60:    ENDIF
 58: ENDDO 61: ENDDO
 59:  62: 
 60: PENERGY=MIN(ENERGYC,ENERGYD) 63: PENERGY=MIN(ENERGYC,ENERGYD)
 61: POPT=ASTEP 64: POPT=ASTEP
 62:   65:  
 63: END SUBROUTINE GOLDEN 66: END SUBROUTINE GOLDEN
 64:  67: 
 65:  68: 


r32919/path.f 2017-07-04 12:30:19.544679423 +0100 r32918/path.f 2017-07-04 12:30:20.344690127 +0100
163: !           sn402: to compute the pushoff step from the transition state?163: !           sn402: to compute the pushoff step from the transition state?
164:             CALL BFGSTS(NSTEPS,Q,ENERGY,VNEW,MFLAG,RMS,EVTS,EVALMAX,VECS,ITDONE,POTCALL,PTEST)164:             CALL BFGSTS(NSTEPS,Q,ENERGY,VNEW,MFLAG,RMS,EVTS,EVALMAX,VECS,ITDONE,POTCALL,PTEST)
165:          ENDIF165:          ENDIF
166:          ETS=ENERGY166:          ETS=ENERGY
167: !        IF (DEBUG) PRINT*,'ts step off plus points in path:'167: !        IF (DEBUG) PRINT*,'ts step off plus points in path:'
168: !        IF (DEBUG) WRITE(*,'(3G20.10)') (Q(J1),J1=1,NOPT)168: !        IF (DEBUG) WRITE(*,'(3G20.10)') (Q(J1),J1=1,NOPT)
169:          DO J1=1,NOPT169:          DO J1=1,NOPT
170:             STEP(J1)=Q(J1)-QINIT(J1)170:             STEP(J1)=Q(J1)-QINIT(J1)
171:          ENDDO171:          ENDDO
172:          IF (PUSHOPTT) THEN  ! set REDOTS to zero? maybe not172:          IF (PUSHOPTT) THEN  ! set REDOTS to zero? maybe not
 173: ! !
 174: ! ! Initially astep=0 the ts at QINIT(J1), bstep=PUSHOFF, corresponding to current Q(J1)
 175: ! !
 176: !             GOLDEN=(SQRT(5.0D0) + 1.0D0) / 2.0D0
 177: !             ASTEP=0.0D0
 178: !             BSTEP=PUSHOFF
 179: !             CSTEP = BSTEP - (BSTEP - ASTEP) / GOLDEN
 180: !             DSTEP = ASTEP + (BSTEP - ASTEP) / GOLDEN
 181: !             IF (DEBUG) WRITE(6,'(A,4G20.10,A)') 'path> golden a, b, c, d=',ASTEP,BSTEP,CSTEP,DSTEP,' initially'
 182: !             DO J1=1,100
 183: !                QC(1:NOPT)=QINIT(1:NOPT)+CSTEP*VECS(1:NOPT)
 184: !                QD(1:NOPT)=QINIT(1:NOPT)+DSTEP*VECS(1:NOPT)
 185: !                CALL POTENTIAL(QC,ENERGYC,VNEW,.FALSE.,.FALSE.,RMS,.FALSE.,.FALSE.)
 186: !                CALL POTENTIAL(QD,ENERGYD,VNEW,.FALSE.,.FALSE.,RMS,.FALSE.,.FALSE.)
 187: !                IF (DEBUG) WRITE(6,'(A,4G25.15)') 'path> energy at ts, C and D and diff=',
 188: !      &                    ETS,ENERGYC,ENERGYD,MIN(ENERGYC,ENERGYD)-ETS  
 189: !                IF (ENERGYC.LT.ENERGYD) THEN
 190: !                   BSTEP=DSTEP
 191: !                ELSE
 192: !                   ASTEP=CSTEP
 193: !                ENDIF
 194: !                CSTEP = BSTEP - (BSTEP - ASTEP) / GOLDEN
 195: !                DSTEP = ASTEP + (BSTEP - ASTEP) / GOLDEN
 196: !                IF (DEBUG) WRITE(6,'(A,4G20.10)') 'path> golden a, b, c, d=',ASTEP,BSTEP,CSTEP,DSTEP
 197: !                IF (ABS(ASTEP-BSTEP).LT.1.0D-6) THEN
 198: !                   ASTEP=(ASTEP+BSTEP)/2.0D0
 199: !                   Q(1:NOPT)=QINIT(1:NOPT)+ASTEP*VECS(1:NOPT)
 200: !                   STEP(1:NOPT)=ASTEP*VECS(1:NOPT)
 201: !                   EXIT
 202: !                ENDIF
 203: !             ENDDO
173:             PINIT=PUSHOFF204:             PINIT=PUSHOFF
174:             CALL GOLDEN(VECS,QINIT,ETS,PINIT,POPT,PENERGY)205:             CALL GOLDEN(VECS,QINIT,ETS,PINIT,POPT,PENERGY)
175:             Q(1:NOPT)=QINIT(1:NOPT)+POPT*VECS(1:NOPT)206:             Q(1:NOPT)=QINIT(1:NOPT)+POPT*VECS(1:NOPT)
176:             STEP(1:NOPT)=POPT*VECS(1:NOPT)207:             STEP(1:NOPT)=POPT*VECS(1:NOPT)
177:             WRITE(6,'(A,3G25.15)') ' path> golden section + pushoff, energy and delta energy: ',208:             WRITE(6,'(A,3G25.15)') 'path> golden section + pushoff, energy and delta energy: ',
178:      &                              POPT,PENERGY,PENERGY-ETS209:      &                              POPT,PENERGY,PENERGY-ETS
179:          ENDIF210:          ENDIF
180:          211:          
181:          IF (RKMIN) RMS=1.0D0212:          IF (RKMIN) RMS=1.0D0
182:          MFLAG=.FALSE.213:          MFLAG=.FALSE.
183:          BFGSSTEP=.FALSE.214:          BFGSSTEP=.FALSE.
184:          BFGSTST=.FALSE.215:          BFGSTST=.FALSE.
185:          IF (BFGSMINT) THEN216:          IF (BFGSMINT) THEN
186:             IF (UNRST.OR.(CHRMMT.AND.INTMINT)) THEN217:             IF (UNRST.OR.(CHRMMT.AND.INTMINT)) THEN
187:                 CALL MYLBFGS(NINTS,MUPDATE,Q,.FALSE.,MFLAG,ENERGY,RMS2,EREAL,RMS,BFGSSTEPS,218:                 CALL MYLBFGS(NINTS,MUPDATE,Q,.FALSE.,MFLAG,ENERGY,RMS2,EREAL,RMS,BFGSSTEPS,
273:                ELSE304:                ELSE
274:                   CALL MYLBFGS(NOPT,MUPDATE,Q,.FALSE.,MFLAG,ENERGY,RMS2,EREAL,RMS,BFGSSTEPS,305:                   CALL MYLBFGS(NOPT,MUPDATE,Q,.FALSE.,MFLAG,ENERGY,RMS2,EREAL,RMS,BFGSSTEPS,
275:      1                         .TRUE.,ITDONE,PTEST,VNEW,.FALSE.,.FALSE.)306:      1                         .TRUE.,ITDONE,PTEST,VNEW,.FALSE.,.FALSE.)
276:                ENDIF307:                ENDIF
277:                NSTEPPLUS=NSTEPPLUS+ITDONE308:                NSTEPPLUS=NSTEPPLUS+ITDONE
278:             ENDIF309:             ENDIF
279:          ENDIF310:          ENDIF
280:       ELSE  ! Another method for finding the pushoff step311:       ELSE  ! Another method for finding the pushoff step
281:          CALL EFOL(Q,MFLAG,1,ENERGY,ITDONE,EVTS,PTEST,FRQSTS,0)312:          CALL EFOL(Q,MFLAG,1,ENERGY,ITDONE,EVTS,PTEST,FRQSTS,0)
282:          ETS=ENERGY313:          ETS=ENERGY
283:          VECS(1:NOPT)=0.0D0 
284:          DUMMY=0.0D0 
285:          DO J1=1,NOPT314:          DO J1=1,NOPT
286:             STEP(J1)=Q(J1)-QINIT(J1)315:             STEP(J1)=Q(J1)-QINIT(J1)
287:             VECS(J1)=STEP(J1) 
288:             DUMMY=DUMMY+VECS(J1)**2 
289:          ENDDO 
290:          DUMMY=SQRT(DUMMY) 
291:          DO J1=1,NOPT 
292:             VECS(J1)=VECS(J1)/DUMMY 
293:          ENDDO316:          ENDDO
294:          IF (PUSHOPTT) THEN317:          IF (PUSHOPTT) THEN
295:             PINIT=PUSHOFF318:             PINIT=PUSHOFF
296:             CALL GOLDEN(VECS,QINIT,ETS,PINIT,POPT,PENERGY)  319:             CALL GOLDEN(VECS,QINIT,ETS,PINIT,POPT,PENERGY)
297:             Q(1:NOPT)=QINIT(1:NOPT)+POPT*VECS(1:NOPT)320:             Q(1:NOPT)=QINIT(1:NOPT)+POPT*VECS(1:NOPT)
298:             STEP(1:NOPT)=POPT*VECS(1:NOPT)321:             STEP(1:NOPT)=POPT*VECS(1:NOPT)
299:             WRITE(6,'(A,3G25.15)') ' path> golden section + pushoff, energy and delta energy: ',322:             WRITE(6,'(A,3G25.15)') 'path> golden section + pushoff, energy and delta energy: ',
300:      &                              POPT,PENERGY,PENERGY-ETS323:      &                              POPT,PENERGY,PENERGY-ETS
301:          ENDIF324:          ENDIF
302:          KNOWE=.FALSE.325:          KNOWE=.FALSE.
303:          KNOWG=.FALSE. ! we don`t know the gradient at the point in Q because VNEW isn`t passed from efol?326:          KNOWG=.FALSE. ! we don`t know the gradient at the point in Q because VNEW isn`t passed from efol?
304:          IF (BFGSMINT) THEN327:          IF (BFGSMINT) THEN
305: C           NOSHIFT=.FALSE.328: C           NOSHIFT=.FALSE.
306:             IF (CHRMMT.AND.INTMINT) THEN329:             IF (CHRMMT.AND.INTMINT) THEN
307:                CALL MYLBFGS(NINTS,MUPDATE,Q,.FALSE.,MFLAG,ENERGY,RMS2,EREAL,RMS,BFGSSTEPS,330:                CALL MYLBFGS(NINTS,MUPDATE,Q,.FALSE.,MFLAG,ENERGY,RMS2,EREAL,RMS,BFGSSTEPS,
308:      1                      .TRUE.,ITDONE,PTEST,VNEW,.FALSE.,.FALSE.)331:      1                      .TRUE.,ITDONE,PTEST,VNEW,.FALSE.,.FALSE.)
309:             ELSE332:             ELSE
466:       DO J1=1,NOPT489:       DO J1=1,NOPT
467:          QPLUS(J1)=Q(J1)490:          QPLUS(J1)=Q(J1)
468:          IF (.NOT.UNRST) Q(J1)=QINIT(J1)-STEP(J1)491:          IF (.NOT.UNRST) Q(J1)=QINIT(J1)-STEP(J1)
469:          IF (HYBRIDMINT) Q(J1)=QINIT(J1)492:          IF (HYBRIDMINT) Q(J1)=QINIT(J1)
470:       ENDDO493:       ENDDO
471:       EPLUS=ENERGY494:       EPLUS=ENERGY
472: !     IF (DEBUG) PRINT*,'ts step off minus points in path:'495: !     IF (DEBUG) PRINT*,'ts step off minus points in path:'
473: !     IF (DEBUG) WRITE(*,'(3G20.10)') (Q(J1),J1=1,NOPT)496: !     IF (DEBUG) WRITE(*,'(3G20.10)') (Q(J1),J1=1,NOPT)
474: 497: 
475:       IF (PUSHOPTT) THEN 498:       IF (PUSHOPTT) THEN 
 499: ! !
 500: ! ! Initially astep=0 the ts at QINIT(J1), bstep=PUSHOFF, corresponding to current Q(J1)
 501: ! !
 502: !          ASTEP=0.0D0
 503: !          BSTEP=-PUSHOFF
 504: !          CSTEP = BSTEP - (BSTEP - ASTEP) / GOLDEN
 505: !          DSTEP = ASTEP + (BSTEP - ASTEP) / GOLDEN
 506: !          IF (DEBUG) WRITE(6,'(A,4G20.10,A)') 'path> golden a, b, c, d=',ASTEP,BSTEP,CSTEP,DSTEP,' initially'
 507: !          DO J1=1,100
 508: !             QC(1:NOPT)=QINIT(1:NOPT)+CSTEP*VECS(1:NOPT)
 509: !             QD(1:NOPT)=QINIT(1:NOPT)+DSTEP*VECS(1:NOPT)
 510: !             CALL POTENTIAL(QC,ENERGYC,VNEW,.FALSE.,.FALSE.,RMS,.FALSE.,.FALSE.)
 511: !             CALL POTENTIAL(QD,ENERGYD,VNEW,.FALSE.,.FALSE.,RMS,.FALSE.,.FALSE.)
 512: !             IF (DEBUG) WRITE(6,'(A,4G25.15)') 'path> energy at ts, C and D and diff=',ETS,ENERGYC,ENERGYD,MIN(ENERGYC,ENERGYD)-ETS  
 513: !             IF (ENERGYC.LT.ENERGYD) THEN
 514: !                BSTEP=DSTEP
 515: !             ELSE
 516: !                ASTEP=CSTEP
 517: !             ENDIF
 518: !             CSTEP = BSTEP - (BSTEP - ASTEP) / GOLDEN
 519: !             DSTEP = ASTEP + (BSTEP - ASTEP) / GOLDEN
 520: !             IF (DEBUG) WRITE(6,'(A,4G20.10)') 'path> golden a, b, c, d=',ASTEP,BSTEP,CSTEP,DSTEP
 521: !             IF (ABS(ASTEP-BSTEP).LT.1.0D-6) THEN
 522: !                ASTEP=(ASTEP+BSTEP)/2.0D0
 523: !                Q(1:NOPT)=QINIT(1:NOPT)+ASTEP*VECS(1:NOPT)
 524: !                STEP(1:NOPT)=ASTEP*VECS(1:NOPT)
 525: !                EXIT
 526: !             ENDIF
 527: !          ENDDO
 528: 
476:          PINIT=-PUSHOFF529:          PINIT=-PUSHOFF
477:          CALL GOLDEN(VECS,QINIT,ETS,PINIT,POPT,PENERGY)530:          CALL GOLDEN(VECS,QINIT,ETS,PINIT,POPT,PENERGY)
478:          Q(1:NOPT)=QINIT(1:NOPT)+POPT*VECS(1:NOPT)531:          Q(1:NOPT)=QINIT(1:NOPT)+POPT*VECS(1:NOPT)
479:          STEP(1:NOPT)=POPT*VECS(1:NOPT)532:          STEP(1:NOPT)=POPT*VECS(1:NOPT)
480:          WRITE(6,'(A,3G25.15)') ' path> golden section - pushoff, energy and delta energy: ',533:          WRITE(6,'(A,3G25.15)') 'path> golden section - pushoff, energy and delta energy: ',
481:      &                              POPT,PENERGY,PENERGY-ETS534:      &                              POPT,PENERGY,PENERGY-ETS
482:       ENDIF535:       ENDIF
483: 536: 
484:       IF (UNRST) THEN ! jmc new intstep stuff 537:       IF (UNRST) THEN ! jmc new intstep stuff 
485:          DO J1=1,nres538:          DO J1=1,nres
486:             c(1,J1)=QINIT(6*(J1-1)+1)539:             c(1,J1)=QINIT(6*(J1-1)+1)
487:             c(2,J1)=QINIT(6*(J1-1)+2)540:             c(2,J1)=QINIT(6*(J1-1)+2)
488:             c(3,J1)=QINIT(6*(J1-1)+3)541:             c(3,J1)=QINIT(6*(J1-1)+3)
489:             c(1,J1+nres)=QINIT(6*(J1-1)+4)542:             c(1,J1+nres)=QINIT(6*(J1-1)+4)
490:             c(2,J1+nres)=QINIT(6*(J1-1)+5)543:             c(2,J1+nres)=QINIT(6*(J1-1)+5)


legend
Lines Added 
Lines changed
 Lines Removed

hdiff - version: 2.1.0