hdiff output

r29244/intlbfgs.f90 2015-11-17 23:32:39.184225579 +0000 r29243/intlbfgs.f90 2015-11-17 23:32:39.376228154 +0000
118:   &      DIAG(3*NATOMS*INTIMAGE), STP(3*NATOMS*INTIMAGE), SEARCHSTEP(0:INTMUPDATE,NOPT*INTIMAGE), &118:   &      DIAG(3*NATOMS*INTIMAGE), STP(3*NATOMS*INTIMAGE), SEARCHSTEP(0:INTMUPDATE,NOPT*INTIMAGE), &
119:   &      GDIF(0:INTMUPDATE,NOPT*INTIMAGE),GLAST(NOPT*INTIMAGE), XSAVE(NOPT*INTIMAGE), &119:   &      GDIF(0:INTMUPDATE,NOPT*INTIMAGE),GLAST(NOPT*INTIMAGE), XSAVE(NOPT*INTIMAGE), &
120:   &      XYZ(NOPT*(INTIMAGE+2)), GGG(NOPT*(INTIMAGE+2)), CHECKG(NOPT*INTIMAGE), IMGFREEZE(INTIMAGE), &120:   &      XYZ(NOPT*(INTIMAGE+2)), GGG(NOPT*(INTIMAGE+2)), CHECKG(NOPT*INTIMAGE), IMGFREEZE(INTIMAGE), &
121:   &      EEE(INTIMAGE+2), STEPIMAGE(INTIMAGE))121:   &      EEE(INTIMAGE+2), STEPIMAGE(INTIMAGE))
122: 122: 
123: SWITCHED=.FALSE.123: SWITCHED=.FALSE.
124: INTIMAGESAVE=INTIMAGE124: INTIMAGESAVE=INTIMAGE
125: NBACKTRACK=1125: NBACKTRACK=1
126: CALL MYCPU_TIME(STIME,.FALSE.)126: CALL MYCPU_TIME(STIME,.FALSE.)
127: PRINT '(A,I6)',' intlbfgs> Maximum number of steps for constraint potential phase is ',INTSTEPS1127: PRINT '(A,I6)',' intlbfgs> Maximum number of steps for constraint potential phase is ',INTSTEPS1
128: WRITE(*,'(A,I6)') ' intlbfgs> Updates: ',INTMUPDATE 
129: PREVGRAD=1.0D100128: PREVGRAD=1.0D100
130: ADDATOM=.FALSE.129: ADDATOM=.FALSE.
131: INTTST=.TRUE.  ! must set this before any possible exit130: INTTST=.TRUE.  ! must set this before any possible exit
132: NFAIL=0131: NFAIL=0
133: IMGFREEZE(1:INTIMAGE)=.FALSE.132: IF (FREEZENODEST) IMGFREEZE(1:INTIMAGE)=.FALSE.
134: D=NOPT*INTIMAGE133: D=NOPT*INTIMAGE
135: U=INTMUPDATE134: U=INTMUPDATE
136: NITERDONE=1135: NITERDONE=1
137: 136: 
138: IF ( D<=0 ) THEN137: IF ( D<=0 ) THEN
139:    PRINT *, 'd is not positive, d=',d138:    PRINT *, 'd is not positive, d=',d
140:    CALL TSUMMARY139:    CALL TSUMMARY
141:    STOP140:    STOP
142: ENDIF141: ENDIF
143: IF ( U<=0 ) THEN142: IF ( U<=0 ) THEN
161: X=>XYZ(NOPT+1:NOPT*(INTIMAGE+1))160: X=>XYZ(NOPT+1:NOPT*(INTIMAGE+1))
162: G=>GGG(NOPT+1:NOPT*(INTIMAGE+1))161: G=>GGG(NOPT+1:NOPT*(INTIMAGE+1))
163: !162: !
164: ! Initialise XYZ163: ! Initialise XYZ
165: !164: !
166: XYZ(1:NOPT)=QSTART(1:NOPT)165: XYZ(1:NOPT)=QSTART(1:NOPT)
167: XYZ(NOPT*(INTIMAGE+1)+1:NOPT*(INTIMAGE+2))=QFINISH(1:NOPT)166: XYZ(NOPT*(INTIMAGE+1)+1:NOPT*(INTIMAGE+2))=QFINISH(1:NOPT)
168: DO J1=1,INTIMAGE+2167: DO J1=1,INTIMAGE+2
169:    XYZ((J1-1)*NOPT+1:J1*NOPT)=((INTIMAGE+2-J1)*QSTART(1:NOPT)+(J1-1)*QFINISH(1:NOPT))/(INTIMAGE+1)168:    XYZ((J1-1)*NOPT+1:J1*NOPT)=((INTIMAGE+2-J1)*QSTART(1:NOPT)+(J1-1)*QFINISH(1:NOPT))/(INTIMAGE+1)
170: ENDDO169: ENDDO
171:       WRITE(*,'(A)') 'intlbfgs> here Z' 
172:       WRITE(*,'(6G20.10)') XYZ(3*(398-1)+1:3*(398-1)+3), & 
173:   &                             XYZ((INTIMAGE+1)*3*NATOMS+3*(398-1)+1:(INTIMAGE+1)*3*NATOMS+3*(398-1)+3) 
174:       WRITE(*,'(6G20.10)') XYZ(3*(400-1)+1:3*(400-1)+3), & 
175:   &                             XYZ((INTIMAGE+1)*3*NATOMS+3*(400-1)+1:(INTIMAGE+1)*3*NATOMS+3*(400-1)+3) 
176:       WRITE(*,'(6G20.10)') QSTART(1:6) 
177:       WRITE(*,'(6G20.10)') QFINISH(1:6) 
178: 170: 
179: NQCIFREEZE=0171: NQCIFREEZE=0
180: IF (FREEZE) THEN172: IF (FREEZE) THEN
181:    PRINT '(A)',' intlbfgs> ERROR *** QCI has not been coded for frozen atoms yet'173:    PRINT '(A)',' intlbfgs> ERROR *** QCI has not been coded for frozen atoms yet'
182:    STOP     174:    STOP     
183: ENDIF175: ENDIF
184: IF (ALLOCATED(INTFROZEN)) DEALLOCATE(INTFROZEN)176: IF (ALLOCATED(INTFROZEN)) DEALLOCATE(INTFROZEN)
185: ALLOCATE(INTFROZEN(NATOMS))177: ALLOCATE(INTFROZEN(NATOMS))
186: INTFROZEN(1:NATOMS)=.FALSE.178: INTFROZEN(1:NATOMS)=.FALSE.
187: DLIST(1:NATOMS)=-1179: DLIST(1:NATOMS)=-1
188: DMOVED(1:NATOMS)=1.0D100180: DMOVED(1:NATOMS)=1.0D100
189: IF (INTFREEZET) THEN181: IF (INTFREEZET) THEN
190:    DUMMY=INTFREEZETOL**2182:    DUMMY=INTFREEZETOL**2
191:    PRINT '(A,6G20.10)',' intlbfgs> INTFREEZETOL,DUMMY=',INTFREEZETOL,DUMMY 
192:    DO J1=1,NATOMS183:    DO J1=1,NATOMS
193:       DF=(XYZ(3*(J1-1)+1)-XYZ((INTIMAGE+1)*3*NATOMS+3*(J1-1)+1))**2 &184:       DF=(XYZ(3*(J1-1)+1)-XYZ((INTIMAGE+1)*3*NATOMS+3*(J1-1)+1))**2 &
194:   &     +(XYZ(3*(J1-1)+2)-XYZ((INTIMAGE+1)*3*NATOMS+3*(J1-1)+2))**2 &185:   &     +(XYZ(3*(J1-1)+2)-XYZ((INTIMAGE+1)*3*NATOMS+3*(J1-1)+2))**2 &
195:   &     +(XYZ(3*(J1-1)+3)-XYZ((INTIMAGE+1)*3*NATOMS+3*(J1-1)+3))**2186:   &     +(XYZ(3*(J1-1)+3)-XYZ((INTIMAGE+1)*3*NATOMS+3*(J1-1)+3))**2
196:       IF (DF.LT.DUMMY) THEN187:       IF (DF.LT.DUMMY) THEN
197:          NQCIFREEZE=NQCIFREEZE+1188:          NQCIFREEZE=NQCIFREEZE+1
198:          INTFROZEN(J1)=.TRUE.189:          INTFROZEN(J1)=.TRUE.
199:          IF (DEBUG) PRINT '(A,I6,A,F12.6,A,I6)',' intlbfgs> atom ',J1,' moves less than threshold: dist^2=',DF,' total=',NQCIFREEZE190: !        IF (DEBUG) PRINT '(A,I6,A,F12.6,A,I6)',' intlbfgs> atom ',J1,' moves less than threshold: dist^2=',DF,' total=',NQCIFREEZE
200:       ENDIF191:       ENDIF
201:       sortd: DO J2=1,J1192:       sortd: DO J2=1,J1
202:          IF (DF.LT.DMOVED(J2)) THEN193:          IF (DF.LT.DMOVED(J2)) THEN
203:             DO J3=J1,J2+1,-1194:             DO J3=J1,J2+1,-1
204:                DMOVED(J3)=DMOVED(J3-1)195:                DMOVED(J3)=DMOVED(J3-1)
205:                DLIST(J3)=DLIST(J3-1)196:                DLIST(J3)=DLIST(J3-1)
206:             ENDDO197:             ENDDO
207:             DMOVED(J2)=DF198:             DMOVED(J2)=DF
208:             DLIST(J2)=J1199:             DLIST(J2)=J1
209:             EXIT sortd200:             EXIT sortd
210:          ENDIF201:          ENDIF
211:       ENDDO sortd202:       ENDDO sortd
212:    ENDDO203:    ENDDO
213:    PRINT '(A,I6,A,F12.6,A,I6)',' intlbfgs> Total number of atoms moving less than threshold=',NQCIFREEZE204:    PRINT '(A,I6,A,F12.6,A,I6)',' intlbfgs> Total number of atoms moving less than threshold=',NQCIFREEZE
214: ENDIF205: ENDIF
215: 206: 
216:       WRITE(*,'(6G20.10)') XYZ(3*(398-1)+1:3*(398-1)+3), & 
217:   &                             XYZ((INTIMAGE+1)*3*NATOMS+3*(398-1)+1:(INTIMAGE+1)*3*NATOMS+3*(398-1)+3) 
218:       WRITE(*,'(6G20.10)') XYZ(3*(400-1)+1:3*(400-1)+3), & 
219:   &                             XYZ((INTIMAGE+1)*3*NATOMS+3*(400-1)+1:(INTIMAGE+1)*3*NATOMS+3*(400-1)+3) 
220:  
221: IF (NATOMS-NQCIFREEZE.LT.INTFREEZEMIN) THEN207: IF (NATOMS-NQCIFREEZE.LT.INTFREEZEMIN) THEN
222:    DO J1=NATOMS,NATOMS-INTFREEZEMIN+1,-1208:    DO J1=NATOMS,NATOMS-INTFREEZEMIN+1,-1
223:       INTFROZEN(DLIST(J1))=.FALSE.209:       INTFROZEN(DLIST(J1))=.FALSE.
224:    ENDDO210:    ENDDO
225:    NQCIFREEZE=NATOMS-INTFREEZEMIN211:    NQCIFREEZE=NATOMS-INTFREEZEMIN
226:    PRINT '(A,I6,A)',' intlbfgs> Freezing ',NQCIFREEZE,' atoms'212:    PRINT '(A,I6,A)',' intlbfgs> Freezing ',NQCIFREEZE,' atoms'
227: ENDIF213: ENDIF
228: 214: 
229: NLASTGOODE=0215: NLASTGOODE=0
230: LASTGOODE=1.0D100216: LASTGOODE=1.0D100
256:    INTIMAGE=0242:    INTIMAGE=0
257:    XYZ(1:NOPT)=QSTART(1:NOPT)243:    XYZ(1:NOPT)=QSTART(1:NOPT)
258:    XYZ(NOPT*(INTIMAGE+1)+1:NOPT*(INTIMAGE+2))=QFINISH(1:NOPT)244:    XYZ(NOPT*(INTIMAGE+1)+1:NOPT*(INTIMAGE+2))=QFINISH(1:NOPT)
259:    DO J1=1,INTIMAGE+2245:    DO J1=1,INTIMAGE+2
260:       XYZ((J1-1)*NOPT+1:J1*NOPT)=((INTIMAGE+2-J1)*QSTART(1:NOPT)+(J1-1)*QFINISH(1:NOPT))/(INTIMAGE+1)246:       XYZ((J1-1)*NOPT+1:J1*NOPT)=((INTIMAGE+2-J1)*QSTART(1:NOPT)+(J1-1)*QFINISH(1:NOPT))/(INTIMAGE+1)
261:    ENDDO247:    ENDDO
262:    GOTO 678248:    GOTO 678
263: ENDIF249: ENDIF
264: 250: 
265: NACTIVE=0251: NACTIVE=0
266: TURNONORDER(1:NATOMS)=0 
267: ATOMACTIVE(1:NATOMS)=.FALSE.252: ATOMACTIVE(1:NATOMS)=.FALSE.
268: IF (INTFREEZET) THEN253: IF (INTFREEZET) THEN
269:    DO J1=1,NATOMS254:    DO J1=1,NATOMS
270:       IF (INTFROZEN(J1)) THEN255:       IF (INTFROZEN(J1)) THEN
271: ! 256: ! 
272: ! linear interpolation 257: ! linear interpolation 
273: ! 258: ! 
274:          DO J2=2,INTIMAGE+1259:          DO J2=2,INTIMAGE+1
275:             XYZ((J2-1)*3*NATOMS+3*(J1-1)+1:(J2-1)*3*NATOMS+3*(J1-1)+3)= &260:             XYZ((J2-1)*3*NATOMS+3*(J1-1)+1:(J2-1)*3*NATOMS+3*(J1-1)+3)= &
276:   &            (INTIMAGE-J2+2)*XYZ(3*(J1-1)+1:3*(J1-1)+3)/(INTIMAGE+1) &261:   &            (INTIMAGE-J2+2)*XYZ(3*(J1-1)+1:3*(J1-1)+3)/(INTIMAGE+1) &
308:    XSAVE(1:D)=X(1:D)293:    XSAVE(1:D)=X(1:D)
309:    GOTO 567294:    GOTO 567
310: ENDIF295: ENDIF
311: DO J1=1,NCONSTRAINT296: DO J1=1,NCONSTRAINT
312:    DF=SQRT((XYZ(3*(CONI(J1)-1)+1)-XYZ((INTIMAGE+1)*3*NATOMS+3*(CONI(J1)-1)+1))**2 &297:    DF=SQRT((XYZ(3*(CONI(J1)-1)+1)-XYZ((INTIMAGE+1)*3*NATOMS+3*(CONI(J1)-1)+1))**2 &
313:   &       +(XYZ(3*(CONI(J1)-1)+2)-XYZ((INTIMAGE+1)*3*NATOMS+3*(CONI(J1)-1)+2))**2 &298:   &       +(XYZ(3*(CONI(J1)-1)+2)-XYZ((INTIMAGE+1)*3*NATOMS+3*(CONI(J1)-1)+2))**2 &
314:   &       +(XYZ(3*(CONI(J1)-1)+3)-XYZ((INTIMAGE+1)*3*NATOMS+3*(CONI(J1)-1)+3))**2)&299:   &       +(XYZ(3*(CONI(J1)-1)+3)-XYZ((INTIMAGE+1)*3*NATOMS+3*(CONI(J1)-1)+3))**2)&
315:   &  +SQRT((XYZ(3*(CONJ(J1)-1)+1)-XYZ((INTIMAGE+1)*3*NATOMS+3*(CONJ(J1)-1)+1))**2 &300:   &  +SQRT((XYZ(3*(CONJ(J1)-1)+1)-XYZ((INTIMAGE+1)*3*NATOMS+3*(CONJ(J1)-1)+1))**2 &
316:   &       +(XYZ(3*(CONJ(J1)-1)+2)-XYZ((INTIMAGE+1)*3*NATOMS+3*(CONJ(J1)-1)+2))**2 &301:   &       +(XYZ(3*(CONJ(J1)-1)+2)-XYZ((INTIMAGE+1)*3*NATOMS+3*(CONJ(J1)-1)+2))**2 &
317:   &       +(XYZ(3*(CONJ(J1)-1)+3)-XYZ((INTIMAGE+1)*3*NATOMS+3*(CONJ(J1)-1)+3))**2)302:   &       +(XYZ(3*(CONJ(J1)-1)+3)-XYZ((INTIMAGE+1)*3*NATOMS+3*(CONJ(J1)-1)+3))**2)
318:    IF (J1.EQ.3505) THEN 
319:       WRITE(*,'(A,3I10)') 'intlbfgs> J1,CONI(J1),CONJ(J1)=',J1,CONI(J1),CONJ(J1) 
320:       WRITE(*,'(6G20.10)') XYZ(3*(CONI(J1)-1)+1:3*(CONI(J1)-1)+3), & 
321:   &                             XYZ((INTIMAGE+1)*3*NATOMS+3*(CONI(J1)-1)+1:(INTIMAGE+1)*3*NATOMS+3*(CONI(J1)-1)+3) 
322:       WRITE(*,'(6G20.10)') XYZ(3*(CONJ(J1)-1)+1:3*(CONJ(J1)-1)+3), & 
323:   &                             XYZ((INTIMAGE+1)*3*NATOMS+3*(CONJ(J1)-1)+1:(INTIMAGE+1)*3*NATOMS+3*(CONJ(J1)-1)+3) 
324:    ENDIF 
325:    IF (DF.LT.DUMMY) THEN303:    IF (DF.LT.DUMMY) THEN
326:       NBEST=J1304:       NBEST=J1
327:       DUMMY=DF305:       DUMMY=DF
328:    ENDIF306:    ENDIF
329: ENDDO307: ENDDO
330: IF (DEBUG) PRINT '(A,I6,A,2I6,A,F15.5)',' intlbfgs> Smallest overall motion for constraint ',NBEST,' atoms ', &308: IF (DEBUG) PRINT '(A,I6,A,2I6,A,F15.5)',' intlbfgs> Smallest overall motion for constraint ',NBEST,' atoms ', &
331:   &                           CONI(NBEST),CONJ(NBEST),' distance=',DUMMY309:   &                           CONI(NBEST),CONJ(NBEST),' distance=',DUMMY
332: 310: 
 311: TURNONORDER(1:NATOMS)=0
333: NTRIES(1:NATOMS)=1312: NTRIES(1:NATOMS)=1
334: IF (ALLOCATED(CONACTIVE)) DEALLOCATE(CONACTIVE)313: IF (ALLOCATED(CONACTIVE)) DEALLOCATE(CONACTIVE)
335: IF (ALLOCATED(NITSTART)) DEALLOCATE(NITSTART)314: IF (ALLOCATED(NITSTART)) DEALLOCATE(NITSTART)
336: ALLOCATE(CONACTIVE(NCONSTRAINT),NITSTART(NCONSTRAINT))315: ALLOCATE(CONACTIVE(NCONSTRAINT),NITSTART(NCONSTRAINT))
337: CONACTIVE(1:NCONSTRAINT)=.FALSE.316: CONACTIVE(1:NCONSTRAINT)=.FALSE.
338: CONACTIVE(NBEST)=.TRUE.317: CONACTIVE(NBEST)=.TRUE.
339: NITSTART(NBEST)=1318: NITSTART(NBEST)=1
340: ATOMACTIVE(CONI(NBEST))=.TRUE.319: ATOMACTIVE(CONI(NBEST))=.TRUE.
341: ATOMACTIVE(CONJ(NBEST))=.TRUE.320: ATOMACTIVE(CONJ(NBEST))=.TRUE.
342: IF (.NOT.INTFROZEN(CONI(NBEST))) THEN321: IF (.NOT.INTFROZEN(CONI(NBEST))) THEN
570:                DMAX=DUMMY549:                DMAX=DUMMY
571:                JMAX=J1550:                JMAX=J1
572:             ENDIF551:             ENDIF
573:             IF (DUMMY.LT.DMIN) THEN552:             IF (DUMMY.LT.DMIN) THEN
574:                DMIN=DUMMY553:                DMIN=DUMMY
575:                JMIN=J1554:                JMIN=J1
576:             ENDIF555:             ENDIF
577: !           IF (DEBUG) PRINT '(A,I6,A,I6,A,G20.10)',' intlbfgs> distance between images ', &556: !           IF (DEBUG) PRINT '(A,I6,A,I6,A,G20.10)',' intlbfgs> distance between images ', &
578: ! &                                                  J1,' and ',J1+1,' is ',DUMMY557: ! &                                                  J1,' and ',J1+1,' is ',DUMMY
579:          ENDDO558:          ENDDO
580:          IF ((DMAX.GT.IMSEPMAX).AND.(INTIMAGE.LT.MAXINTIMAGE).AND.(.NOT.SWITCHED)) THEN559:          IF ((DMAX.GT.IMSEPMAX).AND.(INTIMAGE.LT.MAXINTIMAGE)) THEN
581:             PRINT '(A,I6,A,I6)',' intlbfgs> Add an image between ',JMAX,' and ',JMAX+1560: !           PRINT '(A,I6,A,I6)',' intlbfgs> Add an image between ',JMAX,' and ',JMAX+1
582:             ALLOCATE(DPTMP(3*NATOMS*(INTIMAGE+2)))561:             ALLOCATE(DPTMP(3*NATOMS*(INTIMAGE+2)))
583:             DPTMP(1:3*NATOMS*(INTIMAGE+2))=XYZ(1:3*NATOMS*(INTIMAGE+2))562:             DPTMP(1:3*NATOMS*(INTIMAGE+2))=XYZ(1:3*NATOMS*(INTIMAGE+2))
584:             DEALLOCATE(XYZ)563:             DEALLOCATE(XYZ)
585:             ALLOCATE(XYZ(3*NATOMS*(INTIMAGE+3)))564:             ALLOCATE(XYZ(3*NATOMS*(INTIMAGE+3)))
586:             XYZ(1:3*NATOMS*JMAX)=DPTMP(1:3*NATOMS*JMAX)565:             XYZ(1:3*NATOMS*JMAX)=DPTMP(1:3*NATOMS*JMAX)
587:             XYZ(3*NATOMS*JMAX+1:3*NATOMS*(JMAX+1))=(DPTMP(3*NATOMS*(JMAX-1)+1:3*NATOMS*JMAX) &566:             XYZ(3*NATOMS*JMAX+1:3*NATOMS*(JMAX+1))=(DPTMP(3*NATOMS*(JMAX-1)+1:3*NATOMS*JMAX) &
588:   &                                               + DPTMP(3*NATOMS*JMAX+1:3*NATOMS*(JMAX+1)))/2.0D0567:   &                                               + DPTMP(3*NATOMS*JMAX+1:3*NATOMS*(JMAX+1)))/2.0D0
589:             XYZ(3*NATOMS*(JMAX+1)+1:3*NATOMS*(INTIMAGE+3))=DPTMP(3*NATOMS*JMAX+1:3*NATOMS*(INTIMAGE+2))568:             XYZ(3*NATOMS*(JMAX+1)+1:3*NATOMS*(INTIMAGE+3))=DPTMP(3*NATOMS*JMAX+1:3*NATOMS*(INTIMAGE+2))
590: !569: !
591: ! Save step-taking memories in SEARCHSTEP and GDIF.570: ! Save step-taking memories in SEARCHSTEP and GDIF.
597:             D2TMP(0:INTMUPDATE,1:NOPT*INTIMAGE)=SEARCHSTEP(0:INTMUPDATE,1:NOPT*INTIMAGE)576:             D2TMP(0:INTMUPDATE,1:NOPT*INTIMAGE)=SEARCHSTEP(0:INTMUPDATE,1:NOPT*INTIMAGE)
598:             DEALLOCATE(SEARCHSTEP)577:             DEALLOCATE(SEARCHSTEP)
599:             ALLOCATE(SEARCHSTEP(0:INTMUPDATE,1:NOPT*(INTIMAGE+1)))578:             ALLOCATE(SEARCHSTEP(0:INTMUPDATE,1:NOPT*(INTIMAGE+1)))
600:             DO J1=0,INTMUPDATE579:             DO J1=0,INTMUPDATE
601:                IF (JMAX.GT.1) SEARCHSTEP(J1,1:3*NATOMS*(JMAX-1))=D2TMP(J1,1:3*NATOMS*(JMAX-1))580:                IF (JMAX.GT.1) SEARCHSTEP(J1,1:3*NATOMS*(JMAX-1))=D2TMP(J1,1:3*NATOMS*(JMAX-1))
602:                IF (JMAX.LT.INTIMAGE+1) SEARCHSTEP(J1,3*NATOMS*JMAX+1:3*NATOMS*(INTIMAGE+1))= &581:                IF (JMAX.LT.INTIMAGE+1) SEARCHSTEP(J1,3*NATOMS*JMAX+1:3*NATOMS*(INTIMAGE+1))= &
603:   &                 D2TMP(J1,3*NATOMS*(JMAX-1)+1:3*NATOMS*INTIMAGE)582:   &                 D2TMP(J1,3*NATOMS*(JMAX-1)+1:3*NATOMS*INTIMAGE)
604:                SEARCHSTEP(J1,3*NATOMS*(JMAX-1)+1:3*NATOMS*JMAX)= &583:                SEARCHSTEP(J1,3*NATOMS*(JMAX-1)+1:3*NATOMS*JMAX)= &
605:   &                             D2TMP(J1,3*NATOMS*(MIN(JMAX,INTIMAGE)-1)+1:3*NATOMS*MIN(JMAX,INTIMAGE))584:   &                             D2TMP(J1,3*NATOMS*(MIN(JMAX,INTIMAGE)-1)+1:3*NATOMS*MIN(JMAX,INTIMAGE))
606:             ENDDO585:             ENDDO
607: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 
608:             SEARCHSTEP(0:INTMUPDATE,1:(3*NATOMS)*(INTIMAGE+1))=0.0D0 
609: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 
610:             D2TMP(0:INTMUPDATE,1:NOPT*INTIMAGE)=GDIF(0:INTMUPDATE,1:NOPT*INTIMAGE)586:             D2TMP(0:INTMUPDATE,1:NOPT*INTIMAGE)=GDIF(0:INTMUPDATE,1:NOPT*INTIMAGE)
611:             DEALLOCATE(GDIF)587:             DEALLOCATE(GDIF)
612:             ALLOCATE(GDIF(0:INTMUPDATE,1:NOPT*(INTIMAGE+1)))588:             ALLOCATE(GDIF(0:INTMUPDATE,1:NOPT*(INTIMAGE+1)))
613:             DO J1=0,INTMUPDATE589:             DO J1=0,INTMUPDATE
614:                IF (JMAX.GT.1) GDIF(J1,1:3*NATOMS*(JMAX-1))=D2TMP(J1,1:3*NATOMS*(JMAX-1))590:                IF (JMAX.GT.1) GDIF(J1,1:3*NATOMS*(JMAX-1))=D2TMP(J1,1:3*NATOMS*(JMAX-1))
615:                IF (JMAX.LT.INTIMAGE+1) GDIF(J1,3*NATOMS*JMAX+1:3*NATOMS*(INTIMAGE+1))= &591:                IF (JMAX.LT.INTIMAGE+1) GDIF(J1,3*NATOMS*JMAX+1:3*NATOMS*(INTIMAGE+1))= &
616:   &                 D2TMP(J1,3*NATOMS*(JMAX-1)+1:3*NATOMS*INTIMAGE)592:   &                 D2TMP(J1,3*NATOMS*(JMAX-1)+1:3*NATOMS*INTIMAGE)
617:                GDIF(J1,3*NATOMS*(JMAX-1)+1:3*NATOMS*JMAX)= &593:                GDIF(J1,3*NATOMS*(JMAX-1)+1:3*NATOMS*JMAX)= &
618:   &                       D2TMP(J1,3*NATOMS*(MIN(JMAX,INTIMAGE)-1)+1:3*NATOMS*MIN(JMAX,INTIMAGE))594:   &                       D2TMP(J1,3*NATOMS*(MIN(JMAX,INTIMAGE)-1)+1:3*NATOMS*MIN(JMAX,INTIMAGE))
619:             ENDDO595:             ENDDO
620: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 
621:             GDIF(0:INTMUPDATE,1:(3*NATOMS)*(INTIMAGE+1))=0.0D0 
622: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 
623:             DEALLOCATE(D2TMP)596:             DEALLOCATE(D2TMP)
624: 597: 
625:             DEALLOCATE(TRUEEE,EEETMP,MYGTMP,GTMP,GGG, &598:             DEALLOCATE(TRUEEE,EEETMP,MYGTMP,GTMP,GGG, &
626:   &                    DIAG,STP,GLAST,XSAVE,EEE,STEPIMAGE,CHECKG,IMGFREEZE)599:   &                    DIAG,STP,GLAST,XSAVE,EEE,STEPIMAGE,CHECKG,IMGFREEZE)
627:             ALLOCATE(TRUEEE(INTIMAGE+3), &600:             ALLOCATE(TRUEEE(INTIMAGE+3), &
628:   &                  EEETMP(INTIMAGE+3), MYGTMP(3*NATOMS*(INTIMAGE+1)), &601:   &                  EEETMP(INTIMAGE+3), MYGTMP(3*NATOMS*(INTIMAGE+1)), &
629:   &                  GTMP(3*NATOMS*(INTIMAGE+1)), &602:   &                  GTMP(3*NATOMS*(INTIMAGE+1)), &
630:   &                  DIAG(3*NATOMS*(INTIMAGE+1)), STP(3*NATOMS*(INTIMAGE+1)), &603:   &                  DIAG(3*NATOMS*(INTIMAGE+1)), STP(3*NATOMS*(INTIMAGE+1)), &
631:   &                  GLAST(NOPT*(INTIMAGE+1)), &604:   &                  GLAST(NOPT*(INTIMAGE+1)), &
632:   &                  XSAVE(NOPT*(INTIMAGE+1)), CHECKG(NOPT*(INTIMAGE+1)), IMGFREEZE(INTIMAGE+1), &605:   &                  XSAVE(NOPT*(INTIMAGE+1)), CHECKG(NOPT*(INTIMAGE+1)), IMGFREEZE(INTIMAGE+1), &
651:             D=NOPT*INTIMAGE624:             D=NOPT*INTIMAGE
652:             CALL CHECKREP(INTIMAGE,XYZ,NOPT,0,1)625:             CALL CHECKREP(INTIMAGE,XYZ,NOPT,0,1)
653:             IF (CHECKCONINT) THEN626:             IF (CHECKCONINT) THEN
654:                CALL CONGRAD2(NMAXINT,NMININT,ETOTAL,XYZ,GGG,EEE,IMGFREEZE,RMS)627:                CALL CONGRAD2(NMAXINT,NMININT,ETOTAL,XYZ,GGG,EEE,IMGFREEZE,RMS)
655:             ELSE628:             ELSE
656:                CALL CONGRAD(NMAXINT,NMININT,ETOTAL,XYZ,GGG,EEE,IMGFREEZE,RMS)629:                CALL CONGRAD(NMAXINT,NMININT,ETOTAL,XYZ,GGG,EEE,IMGFREEZE,RMS)
657:             ENDIF630:             ENDIF
658: !           GOTO 864631: !           GOTO 864
659:          ELSEIF ((DMIN.LT.IMSEPMIN).AND.(INTIMAGE.GT.1)) THEN632:          ELSEIF ((DMIN.LT.IMSEPMIN).AND.(INTIMAGE.GT.1)) THEN
660:             IF (JMIN.EQ.1) JMIN=2633:             IF (JMIN.EQ.1) JMIN=2
661:             PRINT '(A,I6,A,I6)',' intlbfgs> Remove image ',JMIN634: !           PRINT '(A,I6,A,I6)',' intlbfgs> Remove image ',JMIN
662:             ALLOCATE(DPTMP(3*NATOMS*(INTIMAGE+2)))635:             ALLOCATE(DPTMP(3*NATOMS*(INTIMAGE+2)))
663:             DPTMP(1:3*NATOMS*(INTIMAGE+2))=XYZ(1:3*NATOMS*(INTIMAGE+2))636:             DPTMP(1:3*NATOMS*(INTIMAGE+2))=XYZ(1:3*NATOMS*(INTIMAGE+2))
664:             DEALLOCATE(XYZ)637:             DEALLOCATE(XYZ)
665:             ALLOCATE(XYZ(3*NATOMS*(INTIMAGE+1)))638:             ALLOCATE(XYZ(3*NATOMS*(INTIMAGE+1)))
666:             XYZ(1:3*NATOMS*(JMIN-1))=DPTMP(1:3*NATOMS*(JMIN-1))639:             XYZ(1:3*NATOMS*(JMIN-1))=DPTMP(1:3*NATOMS*(JMIN-1))
667:             XYZ(3*NATOMS*(JMIN-1)+1:3*NATOMS*(INTIMAGE+1))=DPTMP(3*NATOMS*JMIN+1:3*NATOMS*(INTIMAGE+2))640:             XYZ(3*NATOMS*(JMIN-1)+1:3*NATOMS*(INTIMAGE+1))=DPTMP(3*NATOMS*JMIN+1:3*NATOMS*(INTIMAGE+2))
668: 641: 
669:             DEALLOCATE(DPTMP)642:             DEALLOCATE(DPTMP)
670: !643: !
671: ! Save step-taking memories in SEARCHSTEP and GDIF.644: ! Save step-taking memories in SEARCHSTEP and GDIF.
674: !647: !
675:             ALLOCATE(D2TMP(0:INTMUPDATE,1:NOPT*INTIMAGE))648:             ALLOCATE(D2TMP(0:INTMUPDATE,1:NOPT*INTIMAGE))
676:             D2TMP(0:INTMUPDATE,1:NOPT*INTIMAGE)=SEARCHSTEP(0:INTMUPDATE,1:NOPT*INTIMAGE)649:             D2TMP(0:INTMUPDATE,1:NOPT*INTIMAGE)=SEARCHSTEP(0:INTMUPDATE,1:NOPT*INTIMAGE)
677:             DEALLOCATE(SEARCHSTEP)650:             DEALLOCATE(SEARCHSTEP)
678:             ALLOCATE(SEARCHSTEP(0:INTMUPDATE,1:NOPT*(INTIMAGE-1)))651:             ALLOCATE(SEARCHSTEP(0:INTMUPDATE,1:NOPT*(INTIMAGE-1)))
679:             DO J1=0,INTMUPDATE652:             DO J1=0,INTMUPDATE
680:                SEARCHSTEP(J1,1:3*NATOMS*(JMIN-2))=D2TMP(J1,1:3*NATOMS*(JMIN-2))653:                SEARCHSTEP(J1,1:3*NATOMS*(JMIN-2))=D2TMP(J1,1:3*NATOMS*(JMIN-2))
681:                SEARCHSTEP(J1,3*NATOMS*(JMIN-2)+1:3*NATOMS*(INTIMAGE-1))= &654:                SEARCHSTEP(J1,3*NATOMS*(JMIN-2)+1:3*NATOMS*(INTIMAGE-1))= &
682:   &                     D2TMP(J1,3*NATOMS*(JMIN-1)+1:3*NATOMS*INTIMAGE)655:   &                     D2TMP(J1,3*NATOMS*(JMIN-1)+1:3*NATOMS*INTIMAGE)
683:             ENDDO656:             ENDDO
684: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 
685:             SEARCHSTEP(0:INTMUPDATE,1:(3*NATOMS)*(INTIMAGE-1))=0.0D0 
686: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 
687:             D2TMP(0:INTMUPDATE,1:NOPT*INTIMAGE)=GDIF(0:INTMUPDATE,1:NOPT*INTIMAGE)657:             D2TMP(0:INTMUPDATE,1:NOPT*INTIMAGE)=GDIF(0:INTMUPDATE,1:NOPT*INTIMAGE)
688:             DEALLOCATE(GDIF)658:             DEALLOCATE(GDIF)
689:             ALLOCATE(GDIF(0:INTMUPDATE,1:NOPT*(INTIMAGE-1)))659:             ALLOCATE(GDIF(0:INTMUPDATE,1:NOPT*(INTIMAGE-1)))
690:             DO J1=0,INTMUPDATE660:             DO J1=0,INTMUPDATE
691:                GDIF(J1,1:3*NATOMS*(JMIN-2))=D2TMP(J1,1:3*NATOMS*(JMIN-2))661:                GDIF(J1,1:3*NATOMS*(JMIN-2))=D2TMP(J1,1:3*NATOMS*(JMIN-2))
692:                GDIF(J1,3*NATOMS*(JMIN-2)+1:3*NATOMS*(INTIMAGE-1))= &662:                GDIF(J1,3*NATOMS*(JMIN-2)+1:3*NATOMS*(INTIMAGE-1))= &
693:   &                     D2TMP(J1,3*NATOMS*(JMIN-1)+1:3*NATOMS*INTIMAGE)663:   &                     D2TMP(J1,3*NATOMS*(JMIN-1)+1:3*NATOMS*INTIMAGE)
694:             ENDDO664:             ENDDO
695: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 
696:             GDIF(0:INTMUPDATE,1:(3*NATOMS)*(INTIMAGE-1))=0.0D0 
697: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 
698:             DEALLOCATE(D2TMP)665:             DEALLOCATE(D2TMP)
699: 666: 
700:             DEALLOCATE(TRUEEE,EEETMP,MYGTMP,GTMP,GGG, &667:             DEALLOCATE(TRUEEE,EEETMP,MYGTMP,GTMP,GGG, &
701:   &                    DIAG,STP,GLAST,XSAVE,EEE,STEPIMAGE,CHECKG,IMGFREEZE)668:   &                    DIAG,STP,GLAST,XSAVE,EEE,STEPIMAGE,CHECKG,IMGFREEZE)
702:             ALLOCATE(TRUEEE(INTIMAGE+1),&669:             ALLOCATE(TRUEEE(INTIMAGE+1),&
703:   &                  EEETMP(INTIMAGE+1), MYGTMP(3*NATOMS*(INTIMAGE-1)), &670:   &                  EEETMP(INTIMAGE+1), MYGTMP(3*NATOMS*(INTIMAGE-1)), &
704:   &                  GTMP(3*NATOMS*(INTIMAGE-1)), &671:   &                  GTMP(3*NATOMS*(INTIMAGE-1)), &
705:   &                  DIAG(3*NATOMS*(INTIMAGE-1)), STP(3*NATOMS*(INTIMAGE-1)), &672:   &                  DIAG(3*NATOMS*(INTIMAGE-1)), STP(3*NATOMS*(INTIMAGE-1)), &
706:   &                  GLAST(NOPT*(INTIMAGE-1)), &673:   &                  GLAST(NOPT*(INTIMAGE-1)), &
707:   &                  XSAVE(NOPT*(INTIMAGE-1)), CHECKG(NOPT*(INTIMAGE-1)), IMGFREEZE(INTIMAGE-1), &674:   &                  XSAVE(NOPT*(INTIMAGE-1)), CHECKG(NOPT*(INTIMAGE-1)), IMGFREEZE(INTIMAGE-1), &
949:                  IF (NREPULSIVE.GT.NREPMAX) CALL REPDOUBLE916:                  IF (NREPULSIVE.GT.NREPMAX) CALL REPDOUBLE
950:                  REPI(NREPULSIVE)=J1917:                  REPI(NREPULSIVE)=J1
951:                  REPJ(NREPULSIVE)=J2918:                  REPJ(NREPULSIVE)=J2
952: ! 919: ! 
953: ! Use the minimum of the end point distances and INTCONSTRAINREPCUT for each contact.920: ! Use the minimum of the end point distances and INTCONSTRAINREPCUT for each contact.
954: !921: !
955:                  REPCUT(NREPULSIVE)=MIN(DMIN-1.0D-3,INTCONSTRAINREPCUT)922:                  REPCUT(NREPULSIVE)=MIN(DMIN-1.0D-3,INTCONSTRAINREPCUT)
956:    ENDDO rep2923:    ENDDO rep2
957: ENDDO924: ENDDO
958: 925: 
 926: 
959: !        NBACKTRACK=MAX(MIN(MIN(1.0D0*(NBACKTRACK+1),1.0D0*50),0.1D0*(NACTIVE-2-NQCIFREEZE)),1)927: !        NBACKTRACK=MAX(MIN(MIN(1.0D0*(NBACKTRACK+1),1.0D0*50),0.1D0*(NACTIVE-2-NQCIFREEZE)),1)
960:          NBACKTRACK=MAX(MIN(MIN(1.0D0*(NBACKTRACK+1),1.0D0*50),1.0D0*(NACTIVE-2-NQCIFREEZE)),1.0D0)928:          NBACKTRACK=MAX(MIN(MIN(1.0D0*(NBACKTRACK+1),1.0D0*50),1.0D0*(NACTIVE-2-NQCIFREEZE)),1.0D0)
961: !        IF (DEBUG) PRINT '(A,I6)',' intlbfgs> Number of atoms to backtrack is now ',NBACKTRACK929: !        IF (DEBUG) PRINT '(A,I6)',' intlbfgs> Number of atoms to backtrack is now ',NBACKTRACK
962:          NDUMMY=0930:          NDUMMY=0
963:          DO J1=1,NATOMS931:          DO J1=1,NATOMS
964:             IF (ATOMACTIVE(J1)) NDUMMY=NDUMMY+1932:             IF (ATOMACTIVE(J1)) NDUMMY=NDUMMY+1
965:          ENDDO933:          ENDDO
966:          IF (NDUMMY.NE.NACTIVE) THEN934:          IF (NDUMMY.NE.NACTIVE) THEN
967:             PRINT '(A,I6)',' intlbfgs> ERROR *** inconsistency in number of active atoms. ',NDUMMY,' should be ',NACTIVE935:             PRINT '(A,I6)',' intlbfgs> ERROR *** inconsistency in number of active atoms. ',NDUMMY,' should be ',NACTIVE
968:             DO J1=1,NATOMS936:             DO J1=1,NATOMS
1032: !   &                                 GLAST(J2),(EPLUS-EMINUS)/(2.0D0*DIFF), &1000: !   &                                 GLAST(J2),(EPLUS-EMINUS)/(2.0D0*DIFF), &
1033: !   &                                 (EPLUS-EMINUS)/(2.0D0*DIFF*GLAST(J2))1001: !   &                                 (EPLUS-EMINUS)/(2.0D0*DIFF*GLAST(J2))
1034: !             ELSE1002: !             ELSE
1035: !                WRITE(*,'(A,3I8,3G20.10)') 'OK    ',(J2-1)/NOPT+1,(J2-NOPT*((J2-1)/NOPT)-1)/3+1,J2, &1003: !                WRITE(*,'(A,3I8,3G20.10)') 'OK    ',(J2-1)/NOPT+1,(J2-NOPT*((J2-1)/NOPT)-1)/3+1,J2, &
1036: !   &                                       GLAST(J2),(EPLUS-EMINUS)/(2.0D0*DIFF), &1004: !   &                                       GLAST(J2),(EPLUS-EMINUS)/(2.0D0*DIFF), &
1037: !   &                                       (EPLUS-EMINUS)/(2.0D0*DIFF*GLAST(J2))1005: !   &                                       (EPLUS-EMINUS)/(2.0D0*DIFF*GLAST(J2))
1038: !             ENDIF1006: !             ENDIF
1039: !          ENDIF1007: !          ENDIF
1040: !       ENDDO1008: !       ENDDO
1041: !    ENDIF1009: !    ENDIF
 1010: 
1042:    IF (EXITSTATUS > 0) THEN  1011:    IF (EXITSTATUS > 0) THEN  
1043:       IF ((.NOT.SWITCHED).AND.(EXITSTATUS.EQ.1)) THEN ! add active atom or restart with true potential on1012:       IF ((.NOT.SWITCHED).AND.(EXITSTATUS.EQ.1)) THEN ! add active atom or restart with true potential on
1044:          IF (ETOTAL/INTIMAGE.GT.MAXCONE) GOTO 7771013:          IF (ETOTAL/INTIMAGE.GT.MAXCONE) GOTO 777
1045:          IF (NACTIVE.LT.NATOMS) THEN 1014:          IF (NACTIVE.LT.NATOMS) THEN 
1046:             ADDATOM=.TRUE.1015:             ADDATOM=.TRUE.
1047:             GOTO 7771016:             GOTO 777
1048:          ENDIF1017:          ENDIF
1049:          CALL MYCPU_TIME(FTIME,.FALSE.)1018:          CALL MYCPU_TIME(FTIME,.FALSE.)
1050:          PRINT '(A,I6,A,F12.6,A,I6,A,F10.1)',' intlbfgs> switch on true potential at step ',NITERDONE, &1019:          PRINT '(A,I6,A,F12.6,A,I6,A,F10.1)',' intlbfgs> switch on true potential at step ',NITERDONE, &
1051:   &                                     ' fraction=',INTCONFRAC,' images=',INTIMAGE,' time=',FTIME-STIME1020:   &                                     ' fraction=',INTCONFRAC,' images=',INTIMAGE,' time=',FTIME-STIME
1076:    ENDIF1045:    ENDIF
1077:    777 CONTINUE1046:    777 CONTINUE
1078: !1047: !
1079: ! Compute the new step and gradient change1048: ! Compute the new step and gradient change
1080: !1049: !
1081:    NPT=POINT*D1050:    NPT=POINT*D
1082:    SEARCHSTEP(POINT,:) = STP*SEARCHSTEP(POINT,:)1051:    SEARCHSTEP(POINT,:) = STP*SEARCHSTEP(POINT,:)
1083:    GDIF(POINT,:)=G-GTMP1052:    GDIF(POINT,:)=G-GTMP
1084:    POINT=POINT+1; IF (POINT==INTMUPDATE) POINT=01053:    POINT=POINT+1; IF (POINT==INTMUPDATE) POINT=0
1085: 1054: 
1086:    IF (DUMPINTXYZ.AND.MOD(NITERDONE,DUMPINTXYZFREQ)==0) CALL INTRWG(NACTIVE,NITERDONE,INTIMAGE,XYZ)1055:    IF (DUMPINTXYZ.AND.MOD(NITERDONE,DUMPINTXYZFREQ)==0) CALL RWG(NITERDONE,INTIMAGE,XYZ)
1087:    IF (DUMPINTEOS.AND.MOD(NITERDONE,DUMPINTEOSFREQ)==0) CALL WRITEPROFILE(NITERDONE,EEE,INTIMAGE)1056:    IF (DUMPINTEOS.AND.MOD(NITERDONE,DUMPINTEOSFREQ)==0) CALL WRITEPROFILE(NITERDONE,EEE,INTIMAGE)
1088:    PREVGRAD=RMS1057:    PREVGRAD=RMS
1089: 1058: 
1090:    NITERDONE=NITERDONE+11059:    NITERDONE=NITERDONE+1
1091:    IF (NITERDONE.GT.NSTEPSMAX) EXIT1060:    IF (NITERDONE.GT.NSTEPSMAX) EXIT
1092:    IF ((NIMAGEFREEZE.EQ.INTIMAGE).AND.(NACTIVE.EQ.NATOMS)) THEN1061:    IF ((NIMAGEFREEZE.EQ.INTIMAGE).AND.(NACTIVE.EQ.NATOMS)) THEN
1093:       IF (SWITCHED) THEN1062:       IF (SWITCHED) THEN
1094:          EXIT1063:          EXIT
1095:       ELSE1064:       ELSE
1096:          CALL MYCPU_TIME(FTIME,.FALSE.)1065:          CALL MYCPU_TIME(FTIME,.FALSE.)
1670:          NREPCUT(NNREPULSIVE)=REPCUT(JJ)1639:          NREPCUT(NNREPULSIVE)=REPCUT(JJ)
1671:          GOTO 2461640:          GOTO 246
1672:       ENDIF1641:       ENDIF
1673:    ENDDO 1642:    ENDDO 
1674: 246 CONTINUE1643: 246 CONTINUE
1675: ENDDO1644: ENDDO
1676: IF (DEBUG) PRINT '(A,2I8)',' checkrep> number of active repulsions and total=',NNREPULSIVE,NREPULSIVE1645: IF (DEBUG) PRINT '(A,2I8)',' checkrep> number of active repulsions and total=',NNREPULSIVE,NREPULSIVE
1677: 1646: 
1678: END SUBROUTINE CHECKREP1647: END SUBROUTINE CHECKREP
1679: 1648: 
1680: SUBROUTINE INTRWG(NACTIVE,NITER,INTIMAGE,XYZ)1649: SUBROUTINE RWG(NITER,INTIMAGE,XYZ)
1681: USE PORFUNCS1650: USE PORFUNCS
1682: USE KEY,ONLY: FILTH,FILTHSTR,STOCKT,AMHT,SEQ,NUMGLY,STOCKAAT, RBAAT,ATOMACTIVE1651: USE KEY,ONLY: FILTH,FILTHSTR,STOCKT,AMHT,SEQ,NUMGLY,STOCKAAT, RBAAT
1683: USE COMMONS, ONLY: ZSYM, NRBSITES 1652: USE COMMONS, ONLY: ZSYM, NRBSITES 
1684: USE AMHGLOBALS, ONLY : NMRES1653: USE AMHGLOBALS, ONLY : NMRES
1685: USE COMMONS, ONLY: NATOMS, NOPT1654: USE COMMONS, ONLY: NATOMS, NOPT
1686: IMPLICIT NONE1655: IMPLICIT NONE
1687: CHARACTER(LEN=10) :: XYZFILE   = 'int.xyz   '1656: CHARACTER(LEN=10) :: XYZFILE   = 'int.xyz   '
1688: CHARACTER(LEN=12) :: RBXYZFILE = 'rbint.xyz   '1657: CHARACTER(LEN=12) :: RBXYZFILE = 'rbint.xyz   '
1689: INTEGER,INTENT(IN) :: NITER1658: INTEGER,INTENT(IN) :: NITER
1690: INTEGER :: J1,J2,J3,GLY_COUNT,INTIMAGE,NACTIVE1659: INTEGER :: J1,J2,GLY_COUNT,INTIMAGE
1691: CHARACTER(LEN=80) :: FILENAME,FILENAME2,DUMMYS,DUMMYS21660: CHARACTER(LEN=80) :: FILENAME,FILENAME2,DUMMYS,DUMMYS2
1692: DOUBLE PRECISION XYZ(NOPT*(INTIMAGE+2))1661: DOUBLE PRECISION XYZ(NOPT*(INTIMAGE+2))
1693: 1662: 
1694: IF (FILTH.EQ.0) THEN1663: IF (FILTH.EQ.0) THEN
1695:    FILENAME=XYZFILE1664:    FILENAME=XYZFILE
1696:    IF (RBAAT) FILENAME2=RBXYZFILE1665:    IF (RBAAT) FILENAME2=RBXYZFILE
1697: ELSE1666: ELSE
1698:    FILENAME=TRIM(XYZFILE)//'.'//TRIM(ADJUSTL(FILTHSTR))1667:    FILENAME=TRIM(XYZFILE)//'.'//TRIM(ADJUSTL(FILTHSTR))
1699:    IF (RBAAT) FILENAME2=TRIM(RBXYZFILE)//'.'//TRIM(ADJUSTL(FILTHSTR))1668:    IF (RBAAT) FILENAME2=TRIM(RBXYZFILE)//'.'//TRIM(ADJUSTL(FILTHSTR))
1700: ENDIF 1669: ENDIF 
1764:      &                                  XYZ((J2-1)*NOPT+9*(J1-1)+3-GLY_COUNT*3)1733:      &                                  XYZ((J2-1)*NOPT+9*(J1-1)+3-GLY_COUNT*3)
1765:             WRITE(993,'(a5,1x,3f20.10)') 'C2   ',XYZ((J2-1)*NOPT+9*(J1-1)+4-GLY_COUNT*3),XYZ((J2-1)*NOPT+9*(J1-1)+5-GLY_COUNT*3), &1734:             WRITE(993,'(a5,1x,3f20.10)') 'C2   ',XYZ((J2-1)*NOPT+9*(J1-1)+4-GLY_COUNT*3),XYZ((J2-1)*NOPT+9*(J1-1)+5-GLY_COUNT*3), &
1766:      &                                  XYZ((J2-1)*NOPT+9*(J1-1)+6-GLY_COUNT*3)1735:      &                                  XYZ((J2-1)*NOPT+9*(J1-1)+6-GLY_COUNT*3)
1767:             WRITE(993,'(a5,1x,3f20.10)') 'O    ',XYZ((J2-1)*NOPT+9*(J1-1)+7-GLY_COUNT*3),XYZ((J2-1)*NOPT+9*(J1-1)+8-GLY_COUNT*3), &1736:             WRITE(993,'(a5,1x,3f20.10)') 'O    ',XYZ((J2-1)*NOPT+9*(J1-1)+7-GLY_COUNT*3),XYZ((J2-1)*NOPT+9*(J1-1)+8-GLY_COUNT*3), &
1768:      &                                  XYZ((J2-1)*NOPT+9*(J1-1)+9-GLY_COUNT*3)1737:      &                                  XYZ((J2-1)*NOPT+9*(J1-1)+9-GLY_COUNT*3)
1769:          ENDIF1738:          ENDIF
1770:       ENDDO1739:       ENDDO
1771:    ENDDO1740:    ENDDO
1772: ELSE1741: ELSE
1773:    DO J2=1,INTIMAGE+21742:    DO J2=1,INTIMAGE+2
1774: !     WRITE(993,'(i4/)') NACTIVE1743:       WRITE(993,'(i4/)') natoms
1775:       WRITE(993,'(i4/)') NATOMS1744:       WRITE(993,'(a5,1x,3f20.10)') (ZSYM((j1+2)/3),xyz( (j2-1)*Nopt+j1),&
1776:       DO J3=1,NATOMS1745:     & XYZ((J2-1)*NOPT+J1+1), XYZ((J2-1)*NOPT+J1+2),J1=1,NOPT,3)
1777:          IF (ATOMACTIVE(J3)) THEN 
1778:             WRITE(993,'(A5,1X,3F20.10)') 'LA   ',XYZ((J2-1)*3*NATOMS+3*(J3-1)+1),XYZ((J2-1)*3*NATOMS+3*(J3-1)+2), & 
1779:   &                                                                      XYZ((J2-1)*3*NATOMS+3*(J3-1)+3) 
1780:          ELSE 
1781:             WRITE(993,'(A5,1X,3F20.10)') 'DU   ',XYZ((J2-1)*3*NATOMS+3*(J3-1)+1),XYZ((J2-1)*3*NATOMS+3*(J3-1)+2), & 
1782:   &                                                                      XYZ((J2-1)*3*NATOMS+3*(J3-1)+3) 
1783:          ENDIF 
1784:       ENDDO 
1785:    ENDDO1746:    ENDDO
1786: ENDIF1747: ENDIF
1787: 1748: 
1788: PRINT *, 'rwg> Interpolated image coordinates were saved to xyz file "'//TRIM(FILENAME)//'"'1749: PRINT *, 'rwg> Interpolated image coordinates were saved to xyz file "'//TRIM(FILENAME)//'"'
1789: 1750: 
1790: CLOSE(UNIT=993)1751: CLOSE(UNIT=993)
1791: END SUBROUTINE INTRWG1752: END SUBROUTINE RWG
1792: 1753: 
1793: SUBROUTINE WRITEPROFILE(NITER,EEE,INTIMAGE)1754: SUBROUTINE WRITEPROFILE(NITER,EEE,INTIMAGE)
1794: USE KEY,ONLY: FILTH,FILTHSTR1755: USE KEY,ONLY: FILTH,FILTHSTR
1795: IMPLICIT NONE 1756: IMPLICIT NONE 
1796: INTEGER,INTENT(IN) :: NITER, INTIMAGE1757: INTEGER,INTENT(IN) :: NITER, INTIMAGE
1797: INTEGER :: I,UNIT1758: INTEGER :: I,UNIT
1798: DOUBLE PRECISION :: EEE(INTIMAGE+2)1759: DOUBLE PRECISION :: EEE(INTIMAGE+2)
1799: CHARACTER(LEN=20) :: FILENAME1760: CHARACTER(LEN=20) :: FILENAME
1800: 1761: 
1801: UNIT=9921762: UNIT=992
2475:             ENDIF2436:             ENDIF
2476:          ENDIF2437:          ENDIF
2477:          IF (INTFROZEN(CONIFIX(J1)).AND.INTFROZEN(CONJFIX(J1))) CYCLE2438:          IF (INTFROZEN(CONIFIX(J1)).AND.INTFROZEN(CONJFIX(J1))) CYCLE
2478:          J2=J2+12439:          J2=J2+1
2479:          CONI(J2)=CONIFIX(J1)2440:          CONI(J2)=CONIFIX(J1)
2480:          CONJ(J2)=CONJFIX(J1)2441:          CONJ(J2)=CONJFIX(J1)
2481:          CONDISTREF(J2)=CONDISTREFFIX(J1)2442:          CONDISTREF(J2)=CONDISTREFFIX(J1)
2482:          CONCUT(J2)=CONCUTFIX(J1)2443:          CONCUT(J2)=CONCUTFIX(J1)
2483:       ENDDO2444:       ENDDO
2484:       NCONSTRAINT=J22445:       NCONSTRAINT=J2
2485:       PRINT '(A,I6,A)',' checkperc> After allowing for frozen atoms there are ',NCONSTRAINT,' constraints'2446: !     PRINT '(A,I6,A)',' checkperc> After allowing for frozen atoms there are ',NCONSTRAINT,' constraints'
2486:       RETURN 2447:       RETURN 
2487:    ELSE2448:    ELSE
2488: !2449: !
2489: ! Put reference minima in optimal permutational alignment with reference minimum one.2450: ! Put reference minima in optimal permutational alignment with reference minimum one.
2490: !2451: !
2491:       DO J2=2,NCONGEOM2452:       DO J2=2,NCONGEOM
2492:          CALL MINPERMDIST(CONGEOM(1,1:3*NATOMS),CONGEOM(J2,1:3*NATOMS),NATOMS,DEBUG, &2453:          CALL MINPERMDIST(CONGEOM(1,1:3*NATOMS),CONGEOM(J2,1:3*NATOMS),NATOMS,DEBUG, &
2493:   &                       PARAM1,PARAM2,PARAM3,BULKT,TWOD,D,DIST2,RIGIDBODY,RMAT)2454:   &                       PARAM1,PARAM2,PARAM3,BULKT,TWOD,D,DIST2,RIGIDBODY,RMAT)
2494:       ENDDO2455:       ENDDO
2495:    ENDIF2456:    ENDIF
2496:    ALLOCATE(CONIFIX(INTCONMAX),CONJFIX(INTCONMAX),CONCUTFIX(INTCONMAX),CONDISTREFFIX(INTCONMAX))2457:    ALLOCATE(CONIFIX(INTCONMAX),CONJFIX(INTCONMAX),CONCUTFIX(INTCONMAX),CONDISTREFFIX(INTCONMAX))
2497: ENDIF2458: ENDIF
2498: 2459: 
2499: 51   NCONSTRAINT=0 2460: 51   NCONSTRAINT=0 
2500: MAXCONDIST=-1.0D02461: MAXCONDIST=-1.0D0
2501: MINCONDIST=1.0D1002462: MINCONDIST=1.0D100
2502: IF (NCONGEOM.LT.2) THEN 2463: IF (NCONGEOM.LT.2) THEN 
2503:    DO J2=1,NATOMS2464:    DO J2=1,NATOMS
2504:       DO J3=J2+1,NATOMS2465:       DO J3=J2+1,NATOMS
2505: 2466: 
 2467:          IF ((J3.EQ.701).AND.(J2.EQ.700)) PRINT '(A)','checkperc> doing atoms 700 and 701'
2506:          IF (J3-J2.GT.INTCONSEP) CYCLE ! forbid constraints corresponding to atoms distant in sequence2468:          IF (J3-J2.GT.INTCONSEP) CYCLE ! forbid constraints corresponding to atoms distant in sequence
2507:          IF (INTFROZEN(J2).AND.INTFROZEN(J3)) CYCLE ! no constraints between intfrozen atoms2469:          IF (INTFROZEN(J2).AND.INTFROZEN(J3)) CYCLE ! no constraints between intfrozen atoms
2508:          DS=SQRT((LXYZ(3*(J2-1)+1)-LXYZ(3*(J3-1)+1))**2 &2470:          DS=SQRT((LXYZ(3*(J2-1)+1)-LXYZ(3*(J3-1)+1))**2 &
2509:   &             +(LXYZ(3*(J2-1)+2)-LXYZ(3*(J3-1)+2))**2 &2471:   &             +(LXYZ(3*(J2-1)+2)-LXYZ(3*(J3-1)+2))**2 &
2510:   &             +(LXYZ(3*(J2-1)+3)-LXYZ(3*(J3-1)+3))**2) 2472:   &             +(LXYZ(3*(J2-1)+3)-LXYZ(3*(J3-1)+3))**2) 
 2473:          IF ((J3.EQ.701).AND.(J2.EQ.700)) PRINT '(A,G20.10)','checkperc> DS=',DS
2511:          IF (DS.GT.5.0D0) CYCLE ! don't allow constraints if either endpoint separation is too large DJW2474:          IF (DS.GT.5.0D0) CYCLE ! don't allow constraints if either endpoint separation is too large DJW
2512: !        IF (DS.GT.15.0D0) CYCLE ! don't allow constraints if either endpoint separation is too large DJW2475: !        IF (DS.GT.15.0D0) CYCLE ! don't allow constraints if either endpoint separation is too large DJW
2513:          DF=SQRT((LXYZ(3*NATOMS+3*(J2-1)+1)-LXYZ(3*NATOMS+3*(J3-1)+1))**2 &2476:          DF=SQRT((LXYZ(3*NATOMS+3*(J2-1)+1)-LXYZ(3*NATOMS+3*(J3-1)+1))**2 &
2514:   &             +(LXYZ(3*NATOMS+3*(J2-1)+2)-LXYZ(3*NATOMS+3*(J3-1)+2))**2 &2477:   &             +(LXYZ(3*NATOMS+3*(J2-1)+2)-LXYZ(3*NATOMS+3*(J3-1)+2))**2 &
2515:   &             +(LXYZ(3*NATOMS+3*(J2-1)+3)-LXYZ(3*NATOMS+3*(J3-1)+3))**2) 2478:   &             +(LXYZ(3*NATOMS+3*(J2-1)+3)-LXYZ(3*NATOMS+3*(J3-1)+3))**2) 
2516:          IF (DF.GT.5.0D0) CYCLE ! don't allow constraints if either endpoint separation is too large DJW2479:          IF (DF.GT.5.0D0) CYCLE ! don't allow constraints if either endpoint separation is too large DJW
2517: !        IF (DF.GT.15.0D0) CYCLE ! don't allow constraints if either endpoint separation is too large DJW2480: !        IF (DF.GT.15.0D0) CYCLE ! don't allow constraints if either endpoint separation is too large DJW
2518: !        IF (2.0D0*ABS(DS-DF)/(DS+DF).LT.LINTCONSTRAINTTOL) THEN2481: !        IF (2.0D0*ABS(DS-DF)/(DS+DF).LT.LINTCONSTRAINTTOL) THEN
 2482:          IF ((J3.EQ.701).AND.(J2.EQ.700)) PRINT '(A,G20.10)','checkperc> DF=',DF
 2483:          IF ((J3.EQ.701).AND.(J2.EQ.700)) PRINT '(A,2G20.10)','checkperc> ABS(DS-DF), &
 2484:   &                                      LINTCONSTRAINTTOL=',ABS(DS-DF),LINTCONSTRAINTTOL
2519:          IF (ABS(DS-DF).LT.LINTCONSTRAINTTOL) THEN2485:          IF (ABS(DS-DF).LT.LINTCONSTRAINTTOL) THEN
2520: !2486: !
2521: !  Add constraint for this distance to the list.2487: !  Add constraint for this distance to the list.
2522: !2488: !
2523:             NCONSTRAINT=NCONSTRAINT+12489:             NCONSTRAINT=NCONSTRAINT+1
2524: !           PRINT '(A,2I6,A,I6)','checkperc> Adding constraint for atoms ',J2,J3,'  total=',NCONSTRAINT2490: !           PRINT '(A,2I6,A,I6)','checkperc> Adding constraint for atoms ',J2,J3,'  total=',NCONSTRAINT
2525:          IF ((J3.EQ.701).AND.(J2.EQ.700)) PRINT '(A,2I6,A,I6)','checkperc> Adding constraint for atoms ', &2491:          IF ((J3.EQ.701).AND.(J2.EQ.700)) PRINT '(A,2I6,A,I6)','checkperc> Adding constraint for atoms ', &
2526:   &                                                        J2,J3,'  total=',NCONSTRAINT2492:   &                                                        J2,J3,'  total=',NCONSTRAINT
2527:             IF (NCONSTRAINT.GT.INTCONMAX) CALL CONDOUBLE2493:             IF (NCONSTRAINT.GT.INTCONMAX) CALL CONDOUBLE
2528:             CONI(NCONSTRAINT)=J22494:             CONI(NCONSTRAINT)=J2
2690:           ALPHA(CP+1) = RHO1(CP+1) * SQ2656:           ALPHA(CP+1) = RHO1(CP+1) * SQ
2691:           GTMP(1:D)        = -ALPHA(CP+1)*GDIF(CP,1:D) + GTMP(1:D)2657:           GTMP(1:D)        = -ALPHA(CP+1)*GDIF(CP,1:D) + GTMP(1:D)
2692:      ENDDO2658:      ENDDO
2693:               2659:               
2694:      GTMP(1:D)=DIAG(1)*GTMP(1:D)2660:      GTMP(1:D)=DIAG(1)*GTMP(1:D)
2695: 2661: 
2696:      DO I=1,BOUND2662:      DO I=1,BOUND
2697:           YR= DOT_PRODUCT( GDIF(CP,1:D) , GTMP )2663:           YR= DOT_PRODUCT( GDIF(CP,1:D) , GTMP )
2698:           BETA= RHO1(CP+1)*YR2664:           BETA= RHO1(CP+1)*YR
2699:           BETA= ALPHA(CP+1)-BETA2665:           BETA= ALPHA(CP+1)-BETA
2700: !         WRITE(*,'(A,I8,4G20.10)') 'makestep> I,YR,BETA,RHO1,ALPHA=',I,YR,BETA,RHO1(CP+1),ALPHA(CP+1) 
2701:           GTMP(1:D) = BETA*SEARCHSTEP(CP,1:D) + GTMP(1:D)2666:           GTMP(1:D) = BETA*SEARCHSTEP(CP,1:D) + GTMP(1:D)
2702:           CP=CP+12667:           CP=CP+1
2703: !         IF (CP==M) CP=02668: !         IF (CP==M) CP=0
2704:           IF (CP==INTMUPDATE) CP=02669:           IF (CP==INTMUPDATE) CP=0
2705:      ENDDO2670:      ENDDO
2706:               2671:               
2707:      STP(1:D) = 1.0D02672:      STP(1:D) = 1.0D0
2708: ENDIF MAIN2673: ENDIF MAIN
2709: 2674: 
2710: !  Store the new search direction2675: !  Store the new search direction


legend
Lines Added 
Lines changed
 Lines Removed

hdiff - version: 2.1.0