hdiff output

r33495/congrad.f90 2017-11-17 17:30:12.337105345 +0000 r33494/congrad.f90 2017-11-17 17:30:13.009114255 +0000
 10: ! 10: !
 11: !   You should have received a copy of the GNU General Public License 11: !   You should have received a copy of the GNU General Public License
 12: !   along with this program; if not, write to the Free Software 12: !   along with this program; if not, write to the Free Software
 13: !   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA 13: !   Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
 14: ! 14: !
 15: SUBROUTINE CONGRAD(NMAXINT,NMININT,ETOTAL,XYZ,GGG,EEE,IMGFREEZE,RMS) 15: SUBROUTINE CONGRAD(NMAXINT,NMININT,ETOTAL,XYZ,GGG,EEE,IMGFREEZE,RMS)
 16: USE KEY, ONLY: FROZEN, FREEZE, NREPI, NREPJ, NNREPULSIVE, & 16: USE KEY, ONLY: FROZEN, FREEZE, NREPI, NREPJ, NNREPULSIVE, &
 17:   &            NCONSTRAINT, CONI, CONJ, INTCONSTRAINTDEL, CONDISTREF, INTCONSTRAINTREP, CONDISTREFLOCAL, & 17:   &            NCONSTRAINT, CONI, CONJ, INTCONSTRAINTDEL, CONDISTREF, INTCONSTRAINTREP, CONDISTREFLOCAL, &
 18:   &            CONACTIVE, INTCONSTRAINREPCUT, NREPCUT,INTIMAGE, KINT, IMSEPMAX, ATOMACTIVE, JMAXCON, & 18:   &            CONACTIVE, INTCONSTRAINREPCUT, NREPCUT,INTIMAGE, KINT, IMSEPMAX, ATOMACTIVE, JMAXCON, &
 19:   &            INTFREEZET, INTFROZEN, CONCUTLOCAL, CONCUT, CONCUTABST, CONCUTABS, CONCUTFRACT, CONCUTFRAC, & 19:   &            INTFREEZET, INTFROZEN, CONCUTLOCAL, CONCUT, CONCUTABST, CONCUTABS, CONCUTFRACT, CONCUTFRAC, &
 20:   &  FREEZENODEST, INTSPRINGACTIVET, INTMINFAC, NCONOFF, CONOFFLIST, CONOFFTRIED, INTCONMAX, ECON, EREP, ESPRING, & 20:   &  FREEZENODEST, INTSPRINGACTIVET, INTMINFAC, NCONOFF, CONOFFLIST, CONOFFTRIED, INTCONMAX, ECON, EREP, ESPRING
 21:   &  CONVERGECONTEST, CONVERGEREPTEST 
 22: USE COMMONS, ONLY: NATOMS, NOPT, DEBUG 21: USE COMMONS, ONLY: NATOMS, NOPT, DEBUG
 23: USE PORFUNCS 22: USE PORFUNCS
 24: IMPLICIT NONE 23: IMPLICIT NONE
 25:             24:            
 26: INTEGER :: J1,J2,NI2,NI1,NJ2,NJ1,NMAXINT,NMININT,NREPINT(INTIMAGE+2),ISTAT,NINTMIN,NINTMIN2,MYUNIT,JMAX,IMAX,JMAXNOFF,IMAXNOFF 25: INTEGER :: J1,J2,NI2,NI1,NJ2,NJ1,NMAXINT,NMININT,NREPINT(INTIMAGE+2),ISTAT,NINTMIN,NINTMIN2,MYUNIT,JMAX,IMAX,JMAXNOFF,IMAXNOFF
 27: DOUBLE PRECISION :: ETOTAL, RMS, EMAX, EMAXNOFF 26: DOUBLE PRECISION :: ETOTAL, RMS, EMAX, EMAXNOFF
 28: INTEGER JJMAX(INTIMAGE+2) 27: INTEGER JJMAX(INTIMAGE+2)
 29: DOUBLE PRECISION  EEMAX(INTIMAGE+2) 28: DOUBLE PRECISION  EEMAX(INTIMAGE+2)
 30: DOUBLE PRECISION R1AX,R1AY,R1AZ,R2AX,R2AY,R2AZ,R1BX,R1BY,R1BZ,R2BX,R2BY,R2BZ,D2,D1 29: DOUBLE PRECISION R1AX,R1AY,R1AZ,R2AX,R2AY,R2AZ,R1BX,R1BY,R1BZ,R2BX,R2BY,R2BZ,D2,D1
 31: DOUBLE PRECISION G1(3),G2(3),DINT,G1INT(3),G2INT(3) 30: DOUBLE PRECISION G1(3),G2(3),DINT,G1INT(3),G2INT(3)
 32: DOUBLE PRECISION DUMMY, REPGRAD(3), D12, DSQ2, DSQ1, DSQI 31: DOUBLE PRECISION DUMMY, REPGRAD(3), INTCONST, D12, DSQ2, DSQ1, DSQI
 33: DOUBLE PRECISION CONE(INTIMAGE+2), REPE(INTIMAGE+2),MAXINT,MININT,REPEINT(INTIMAGE+2),RMSIM(INTIMAGE+2) 32: DOUBLE PRECISION CONE(INTIMAGE+2), REPE(INTIMAGE+2),MAXINT,MININT,REPEINT(INTIMAGE+2),RMSIM(INTIMAGE+2)
 34: LOGICAL NOINT 33: LOGICAL NOINT
 35: DOUBLE PRECISION XYZ((3*NATOMS)*(INTIMAGE+2)), GGG((3*NATOMS)*(INTIMAGE+2)), EEE(INTIMAGE+2), CCLOCAL 34: DOUBLE PRECISION XYZ((3*NATOMS)*(INTIMAGE+2)), GGG((3*NATOMS)*(INTIMAGE+2)), EEE(INTIMAGE+2), CCLOCAL
 36: LOGICAL IMGFREEZE(INTIMAGE) 35: LOGICAL IMGFREEZE(INTIMAGE)
 37: DOUBLE PRECISION DPLUS, SPGRAD(3), DCUT, r1amr1bdr2amr2b,r1apr2bmr2amr1bsq,CUTMAX,DISTMAX 36: DOUBLE PRECISION DPLUS, SPGRAD(3), DCUT, r1amr1bdr2amr2b,r1apr2bmr2amr1bsq,CUTMAX,DISTMAX,DIST2MAX
 38: DOUBLE PRECISION CONDMAX, CONREFMAX, CONCUTMAX, DUMMY2, CCLOCAL2 37: DOUBLE PRECISION CONDMAX, CONREFMAX, CONCUTMAX
 39:  38: 
 40: EEE(1:INTIMAGE+2)=0.0D0 39: EEE(1:INTIMAGE+2)=0.0D0
 41: CONE(1:INTIMAGE+2)=0.0D0 40: CONE(1:INTIMAGE+2)=0.0D0
 42: REPE(1:INTIMAGE+2)=0.0D0 41: REPE(1:INTIMAGE+2)=0.0D0
 43: REPEINT(1:INTIMAGE+2)=0.0D0 42: REPEINT(1:INTIMAGE+2)=0.0D0
 44: NREPINT(1:INTIMAGE+2)=0 43: NREPINT(1:INTIMAGE+2)=0
 45: GGG(1:(3*NATOMS)*(INTIMAGE+2))=0.0D0 44: GGG(1:(3*NATOMS)*(INTIMAGE+2))=0.0D0
 46: ECON=0.0D0; EREP=0.0D0 45: ECON=0.0D0; EREP=0.0D0
 47: NINTMIN=0 46: NINTMIN=0
 48: NINTMIN2=0 47: NINTMIN2=0
 88: !              IF (J2.EQ.1) PRINT '(A,I6,A)','J1=',J1,' IMGFREEZE(J1-2)=T and IMGFREEZE(J1-1)=T cycle' 87: !              IF (J2.EQ.1) PRINT '(A,I6,A)','J1=',J1,' IMGFREEZE(J1-2)=T and IMGFREEZE(J1-1)=T cycle'
 89:                CYCLE 88:                CYCLE
 90:             ENDIF 89:             ENDIF
 91:          ENDIF 90:          ENDIF
 92:       ENDIF 91:       ENDIF
 93:       NI1=(3*NATOMS)*(J1-1)+3*(CONI(J2)-1) 92:       NI1=(3*NATOMS)*(J1-1)+3*(CONI(J2)-1)
 94:       NJ1=(3*NATOMS)*(J1-1)+3*(CONJ(J2)-1) 93:       NJ1=(3*NATOMS)*(J1-1)+3*(CONJ(J2)-1)
 95:       R2AX=XYZ(NI1+1); R2AY=XYZ(NI1+2); R2AZ=XYZ(NI1+3) 94:       R2AX=XYZ(NI1+1); R2AY=XYZ(NI1+2); R2AZ=XYZ(NI1+3)
 96:       R2BX=XYZ(NJ1+1); R2BY=XYZ(NJ1+2); R2BZ=XYZ(NJ1+3) 95:       R2BX=XYZ(NJ1+1); R2BY=XYZ(NJ1+2); R2BZ=XYZ(NJ1+3)
 97:       D2=SQRT((R2AX-R2BX)**2+(R2AY-R2BY)**2+(R2AZ-R2BZ)**2) 96:       D2=SQRT((R2AX-R2BX)**2+(R2AY-R2BY)**2+(R2AZ-R2BZ)**2)
 98:       DUMMY=D2-CONDISTREFLOCAL(J2)   97:       IF (ABS(D2-CONDISTREFLOCAL(J2)).GT.CCLOCAL) THEN 
 99:       DUMMY2=DUMMY**2 98:          DUMMY=D2-CONDISTREFLOCAL(J2)  
100:       CCLOCAL2=CCLOCAL**2 
101: ! 
102: ! Reduced form for penalty function and gradient - multiply through by CCLOCAL**2/INTCONSTRAINTDEL 
103: ! Should save CCLOCAL**2 
104: ! Could save DUMMY**2-CCLOCAL**2 outside IF block and base the test on difference of squares 
105: ! 
106:       IF (DUMMY2.GT.CCLOCAL2) THEN  
107:          G2(1)=(R2AX-R2BX)/D2 99:          G2(1)=(R2AX-R2BX)/D2
108:          G2(2)=(R2AY-R2BY)/D2100:          G2(2)=(R2AY-R2BY)/D2
109:          G2(3)=(R2AZ-R2BZ)/D2101:          G2(3)=(R2AZ-R2BZ)/D2
110: !102:          REPGRAD(1:3)=2*INTCONSTRAINTDEL*((DUMMY/CCLOCAL)**2-1.0D0)*DUMMY*G2(1:3)
111: !        REPGRAD(1:3)=2*INTCONSTRAINTDEL*((DUMMY/CCLOCAL)**2-1.0D0)*DUMMY*G2(1:3)103:          DUMMY=INTCONSTRAINTDEL*(DUMMY**2-CCLOCAL**2)**2/(2.0D0*CCLOCAL**2)
112: !        DUMMY=INTCONSTRAINTDEL*(DUMMY**2-CCLOCAL**2)**2/(2.0D0*CCLOCAL**2) 
113: ! 
114:          REPGRAD(1:3)=2*(DUMMY2-CCLOCAL2)*DUMMY*G2(1:3) 
115:          DUMMY=(DUMMY2-CCLOCAL2)**2/2.0D0 
116:          EEE(J1)=EEE(J1)  +DUMMY104:          EEE(J1)=EEE(J1)  +DUMMY
117:          ECON=ECON        +DUMMY105:          ECON=ECON        +DUMMY
118:          IF (DUMMY.GT.EMAX) THEN106:          IF (DUMMY.GT.EMAX) THEN
119:             IMAX=J1107:             IMAX=J1
120:             JMAX=J2108:             JMAX=J2
121:             EMAX=DUMMY109:             EMAX=DUMMY
122:             CONDMAX=D2 
123:             CONREFMAX=CONDISTREFLOCAL(J2) 
124:             CONCUTMAX=CCLOCAL 
125:          ENDIF110:          ENDIF
126:          IF (DUMMY.GT.EMAXNOFF) THEN111:          IF (DUMMY.GT.EMAXNOFF) THEN
127:             IF (.NOT.CONOFFTRIED(J2)) THEN112:             IF (.NOT.CONOFFTRIED(J2)) THEN
128:                IMAXNOFF=J1113:                IMAXNOFF=J1
129:                JMAXNOFF=J2114:                JMAXNOFF=J2
130:                EMAXNOFF=DUMMY115:                EMAXNOFF=DUMMY
131:             ENDIF116:             ENDIF
132:          ENDIF117:          ENDIF
133:          IF (DUMMY.GT.EEMAX(J1)) THEN118:          IF (DUMMY.GT.EEMAX(J1)) THEN
134:             JJMAX(J1)=J2119:             JJMAX(J1)=J2
135:             EEMAX(J1)=DUMMY120:             EEMAX(J1)=DUMMY
136:          ENDIF121:          ENDIF
137:          CONE(J1)=CONE(J1)+DUMMY122:          CONE(J1)=CONE(J1)+DUMMY
138:          GGG(NI1+1:NI1+3)=GGG(NI1+1:NI1+3)+REPGRAD(1:3)123:          GGG(NI1+1:NI1+3)=GGG(NI1+1:NI1+3)+REPGRAD(1:3)
139:          GGG(NJ1+1:NJ1+3)=GGG(NJ1+1:NJ1+3)-REPGRAD(1:3)124:          GGG(NJ1+1:NJ1+3)=GGG(NJ1+1:NJ1+3)-REPGRAD(1:3)
140:       ENDIF125:       ENDIF
141: !     WRITE(MYUNIT,'(A,2I6,5G20.10)') 'J1,J2,D2,CONDISTREFLOCAL,CCLOCAL,EEE,CONE=',J1,J2,D2,CONDISTREFLOCAL(J2),CCLOCAL,EEE(J1),CONE(J1)126: !     WRITE(MYUNIT,'(A,2I6,5G20.10)') 'J1,J2,D2,CONDISTREFLOCAL,CCLOCAL,EEE,CONE=',J1,J2,D2,CONDISTREFLOCAL(J2),CCLOCAL,EEE(J1),CONE(J1)
142:    ENDDO127:    ENDDO
143: ENDDO128: ENDDO
144: IF (JMAX.GT.0) THEN129: IF (JMAX.GT.0) THEN
145:    WRITE(*,'(A,I6,A,I6,A,2I6,A,G20.10,A,3G20.10)') ' congrad> Highest constraint contribution for any image in image ',IMAX, &130:    WRITE(*,'(A,I6,A,I6,A,2I6,A,3G20.10)') ' congrad> Highest constraint contribution for any image in image ',IMAX, &
146:  & ' constraint ',JMAX, &131:  & ' constraint ',JMAX, &
147:  &                              ' atoms ',CONI(JMAX),CONJ(JMAX),' value=',EMAX,' d,ref,cutoff=',CONDMAX,CONREFMAX,CONCUTMAX132:  &                              ' atoms ',CONI(JMAX),CONJ(JMAX),'d,ref,cutoff=',CONDMAX,CONREFMAX,CONCUTMAX
148: 133: 
149: ENDIF134: ENDIF
150: CONVERGECONTEST=EMAX/INTCONSTRAINTDEL 
151: IF (JMAXNOFF.GT.0) THEN135: IF (JMAXNOFF.GT.0) THEN
152:    WRITE(*,'(A,I6,A,I6,A,2I8,A,G20.10,A,I8)') ' congrad> Highest constraint contribution never turned off for any image in image ',IMAXNOFF, &136:    WRITE(*,'(A,I6,A,I6,A,2I8,A,G20.10,A,I8)') ' congrad> Highest constraint contribution never turned off for any image in image ',IMAXNOFF, &
153:  & ' constraint ',JMAXNOFF, &137:  & ' constraint ',JMAXNOFF, &
154:  &                              ' atoms ',CONI(JMAX),CONJ(JMAX),' value=',EMAX,' off=',NCONOFF138:  &                              ' atoms ',CONI(JMAX),CONJ(JMAX),' value=',EMAX,' off=',NCONOFF
155: ELSEIF (JMAX.GT.0) THEN139: ELSEIF (JMAX.GT.0) THEN
156:    JMAXNOFF=JMAX140:    JMAXNOFF=JMAX
157:    WRITE(*,'(A,I6,A,I6,A,2I8,A,I8)') ' congrad> *** Using highest constraint contribution for any image, setting NCONOFF to 0'141:    WRITE(*,'(A,I6,A,I6,A,2I8,A,I8)') ' congrad> *** Using highest constraint contribution for any image, setting NCONOFF to 0'
158:    CONOFFTRIED(1:INTCONMAX)=.FALSE.142:    CONOFFTRIED(1:INTCONMAX)=.FALSE.
159:    NCONOFF=0143:    NCONOFF=0
160: ENDIF144: ENDIF
165: 149: 
166: ! INTCONST=INTCONSTRAINREPCUT**13150: ! INTCONST=INTCONSTRAINREPCUT**13
167: 151: 
168: EMAX=-1.0D200152: EMAX=-1.0D200
169: EEMAX(1:INTIMAGE+2)=-1.0D200153: EEMAX(1:INTIMAGE+2)=-1.0D200
170: JJMAX(1:INTIMAGE+2)=-1154: JJMAX(1:INTIMAGE+2)=-1
171: JMAX=-1155: JMAX=-1
172: IMAX=-1156: IMAX=-1
173: 157: 
174: DO J2=1,NNREPULSIVE158: DO J2=1,NNREPULSIVE
175: ! !  INTCONST=NREPCUT(J2)**13159: !  INTCONST=NREPCUT(J2)**13
176: !  INTCONST=NREPCUT(J2)**3160:    INTCONST=NREPCUT(J2)**3
177:    DO J1=2,INTIMAGE+2161:    DO J1=2,INTIMAGE+2
178: !  DO J1=1,INTIMAGE+2 ! can change when zero energies are confirmed for end images162: !  DO J1=1,INTIMAGE+2 ! can change when zero energies are confirmed for end images
179:       IF (FREEZENODEST) THEN163:       IF (FREEZENODEST) THEN
180:          IF (J1.EQ.2) THEN164:          IF (J1.EQ.2) THEN
181:             IF (IMGFREEZE(1)) CYCLE165:             IF (IMGFREEZE(1)) CYCLE
182:          ELSE IF (J1.EQ.INTIMAGE+2) THEN166:          ELSE IF (J1.EQ.INTIMAGE+2) THEN
183:             IF (IMGFREEZE(INTIMAGE)) CYCLE167:             IF (IMGFREEZE(INTIMAGE)) CYCLE
184:          ELSE168:          ELSE
185:             IF (IMGFREEZE(J1-2).AND.IMGFREEZE(J1-1)) CYCLE169:             IF (IMGFREEZE(J1-2).AND.IMGFREEZE(J1-1)) CYCLE
186:          ENDIF170:          ENDIF
219:          D2=SQRT(DSQ2)203:          D2=SQRT(DSQ2)
220:          G2(1:3)=G2(1:3)/D2204:          G2(1:3)=G2(1:3)/D2
221:          G1(1:3)=G1(1:3)/D1205:          G1(1:3)=G1(1:3)/D1
222:       ELSE206:       ELSE
223:          CALL MINMAXD2R(D2,D1,DINT,DSQ2,DSQ1,DSQI,G1,G2,G1INT,G2INT,NOINT,.FALSE.,NREPCUT(J2),r1amr1bdr2amr2b,r1apr2bmr2amr1bsq)207:          CALL MINMAXD2R(D2,D1,DINT,DSQ2,DSQ1,DSQI,G1,G2,G1INT,G2INT,NOINT,.FALSE.,NREPCUT(J2),r1amr1bdr2amr2b,r1apr2bmr2amr1bsq)
224:       ENDIF208:       ENDIF
225: !209: !
226: ! Skip image INTIMAGE+2 - no non-zero gradients on other images and no energy contributions.210: ! Skip image INTIMAGE+2 - no non-zero gradients on other images and no energy contributions.
227: !211: !
228:       IF ((D2.LT.NREPCUT(J2)).AND.(J1.LT.INTIMAGE+2)) THEN ! terms for image J1 - non-zero derivatives only for J1212:       IF ((D2.LT.NREPCUT(J2)).AND.(J1.LT.INTIMAGE+2)) THEN ! terms for image J1 - non-zero derivatives only for J1
229: ! !        D12=DSQ2**6213: !        D12=DSQ2**6
230:          D12=DSQ2214:          D12=DSQ2
231: ! !        DUMMY=INTCONSTRAINTREP*(1.0D0/D12+(12.0D0*D2-13.0D0*INTCONSTRAINREPCUT)/INTCONST)215: !        DUMMY=INTCONSTRAINTREP*(1.0D0/D12+(12.0D0*D2-13.0D0*INTCONSTRAINREPCUT)/INTCONST)
232: ! !        DUMMY=INTCONSTRAINTREP*(1.0D0/D12+(12.0D0*D2-13.0D0*NREPCUT(J2))/INTCONST)216: !        DUMMY=INTCONSTRAINTREP*(1.0D0/D12+(12.0D0*D2-13.0D0*NREPCUT(J2))/INTCONST)
233: !        DUMMY=INTCONSTRAINTREP*(1.0D0/D12+(2.0D0*D2-3.0D0*NREPCUT(J2))/INTCONST)217:          DUMMY=INTCONSTRAINTREP*(1.0D0/D12+(2.0D0*D2-3.0D0*NREPCUT(J2))/INTCONST)
234: !        DUMMY=(1.0D0/D12+(2.0D0*D2-3.0D0*NREPCUT(J2))/INTCONST)*NREPCUT(J2)**2 
235:          DUMMY=DCUT/D12+2.0D0*D2/NREPCUT(J2)-3.0D0 
236:          EEE(J1)=EEE(J1)+DUMMY218:          EEE(J1)=EEE(J1)+DUMMY
237: !        IF (PRINTE.AND.(DUMMY.GT.1.0D-4)) THEN219: !        IF (PRINTE.AND.(DUMMY.GT.1.0D-4)) THEN
238: !           WRITE(*, '(A,2I6,2L5,G20.10)') 'R1 NREPI,NREPJ,INTFROZEN(NREPI),INTFROZEN(NREPJ),DUMMY=', &220: !           WRITE(*, '(A,2I6,2L5,G20.10)') 'R1 NREPI,NREPJ,INTFROZEN(NREPI),INTFROZEN(NREPJ),DUMMY=', &
239: ! &                                     NREPI(J2),NREPJ(J2),INTFROZEN(NREPI(J2)),INTFROZEN(NREPJ(J2)),DUMMY221: ! &                                     NREPI(J2),NREPJ(J2),INTFROZEN(NREPI(J2)),INTFROZEN(NREPJ(J2)),DUMMY
240: !        ENDIF222: !        ENDIF
241:          IF (DUMMY.GT.EMAX) THEN223:          IF (DUMMY.GT.EMAX) THEN
242:             IMAX=J1224:             IMAX=J1
243:             JMAX=J2225:             JMAX=J2
 226:             DIST2MAX=D12
244:             CUTMAX=NREPCUT(J2)227:             CUTMAX=NREPCUT(J2)
245:             DISTMAX=D2228:             DISTMAX=D2
246:             EMAX=DUMMY229:             EMAX=DUMMY
247:          ENDIF230:          ENDIF
248:          IF (DUMMY.GT.EEMAX(J1)) THEN231:          IF (DUMMY.GT.EEMAX(J1)) THEN
249:             JJMAX(J1)=J2232:             JJMAX(J1)=J2
250:             EEMAX(J1)=DUMMY233:             EEMAX(J1)=DUMMY
251:          ENDIF234:          ENDIF
252:          REPE(J1)=REPE(J1)+DUMMY235:          REPE(J1)=REPE(J1)+DUMMY
253:          EREP=EREP+DUMMY236:          EREP=EREP+DUMMY
254: ! !        DUMMY=-12.0D0*INTCONSTRAINTREP*(1.0D0/(D2*D12)-1.0D0/INTCONST)237: !        DUMMY=-12.0D0*INTCONSTRAINTREP*(1.0D0/(D2*D12)-1.0D0/INTCONST)
255: !        DUMMY=-2.0D0*INTCONSTRAINTREP*(1.0D0/(D2*D12)-1.0D0/INTCONST)238:          DUMMY=-2.0D0*INTCONSTRAINTREP*(1.0D0/(D2*D12)-1.0D0/INTCONST)
256: ! 
257: !        DUMMY=-2.0D0*(1.0D0/(D2*D12)-1.0D0/INTCONST)*NREPCUT(J2)**2 
258:          DUMMY=-2.0D0*(DCUT/(D2*D12)-1.0D0/NREPCUT(J2)) 
259:          REPGRAD(1:3)=DUMMY*G2(1:3)239:          REPGRAD(1:3)=DUMMY*G2(1:3)
260: !        WRITE(*, '(A,4I6,G15.5)') 'min J1,J2,REPI,REPJ,REPGRAD=',J1,J2,NREPI(J2),NREPJ(J2), &240: !        WRITE(*, '(A,4I6,G15.5)') 'min J1,J2,REPI,REPJ,REPGRAD=',J1,J2,NREPI(J2),NREPJ(J2), &
261: ! &                              SQRT(REPGRAD(1)**2+REPGRAD(2)**2+REPGRAD(3)**2)241: ! &                              SQRT(REPGRAD(1)**2+REPGRAD(2)**2+REPGRAD(3)**2)
262:          GGG(NI2+1:NI2+3)=GGG(NI2+1:NI2+3)+REPGRAD(1:3)242:          GGG(NI2+1:NI2+3)=GGG(NI2+1:NI2+3)+REPGRAD(1:3)
263:          GGG(NJ2+1:NJ2+3)=GGG(NJ2+1:NJ2+3)-REPGRAD(1:3)243:          GGG(NJ2+1:NJ2+3)=GGG(NJ2+1:NJ2+3)-REPGRAD(1:3)
264:       ENDIF244:       ENDIF
265: !245: !
266: ! For internal minima we are counting edges. 246: ! For internal minima we are counting edges. 
267: ! Edge J1 is between images J1-1 and J1, starting from J1=2.247: ! Edge J1 is between images J1-1 and J1, starting from J1=2.
268: ! Energy contributions are shared evenly, except for248: ! Energy contributions are shared evenly, except for
269: ! edge 1, which is assigned to image 2, and edge INTIMAGE+1, which249: ! edge 1, which is assigned to image 2, and edge INTIMAGE+1, which
270: ! is assigned to image INTIMAGE+1. Gradients are set to zero for250: ! is assigned to image INTIMAGE+1. Gradients are set to zero for
271: ! the end images.251: ! the end images.
272: !252: !
273:       DUMMY=0.0D0253:       DUMMY=0.0D0
274:       IF ((.NOT.NOINT).AND.(DINT.LT.NREPCUT(J2))) THEN254:       IF ((.NOT.NOINT).AND.(DINT.LT.NREPCUT(J2))) THEN
275:          NINTMIN2=NINTMIN2+1255:          NINTMIN2=NINTMIN2+1
276: !        D12=DSQI**6256: !        D12=DSQI**6
277:          D12=DSQI257:          D12=DSQI
278: ! !        DUMMY=INTCONSTRAINTREP*(1.0D0/D12+(12.0D0*DINT-13.0D0*NREPCUT(J2))/INTCONST)258: !        DUMMY=INTCONSTRAINTREP*(1.0D0/D12+(12.0D0*DINT-13.0D0*NREPCUT(J2))/INTCONST)
279: !        DUMMY=INTMINFAC*INTCONSTRAINTREP*(1.0D0/D12+(2.0D0*DINT-3.0D0*NREPCUT(J2))/INTCONST)259:          DUMMY=INTMINFAC*INTCONSTRAINTREP*(1.0D0/D12+(2.0D0*DINT-3.0D0*NREPCUT(J2))/INTCONST)
280: !        DUMMY=INTMINFAC*(1.0D0/D12+(2.0D0*DINT-3.0D0*NREPCUT(J2))/INTCONST)*NREPCUT(J2)**2 
281:          DUMMY=INTMINFAC*(DCUT/D12+2.0D0*DINT/NREPCUT(J2)-3.0D0) 
282:          IF (J1.EQ.2) THEN260:          IF (J1.EQ.2) THEN
283:             EEE(J1)=EEE(J1)+DUMMY261:             EEE(J1)=EEE(J1)+DUMMY
284:             REPEINT(J1)=REPEINT(J1)+DUMMY262:             REPEINT(J1)=REPEINT(J1)+DUMMY
285:             NREPINT(J1)=NREPINT(J1)+1263:             NREPINT(J1)=NREPINT(J1)+1
286:          ELSE IF (J1.LT.INTIMAGE+2) THEN264:          ELSE IF (J1.LT.INTIMAGE+2) THEN
287:             EEE(J1)=EEE(J1)+DUMMY/2.0D0265:             EEE(J1)=EEE(J1)+DUMMY/2.0D0
288:             EEE(J1-1)=EEE(J1-1)+DUMMY/2.0D0266:             EEE(J1-1)=EEE(J1-1)+DUMMY/2.0D0
289:             REPEINT(J1)=REPEINT(J1)+DUMMY/2.0D0267:             REPEINT(J1)=REPEINT(J1)+DUMMY/2.0D0
290:             REPEINT(J1-1)=REPEINT(J1-1)+DUMMY/2.0D0268:             REPEINT(J1-1)=REPEINT(J1-1)+DUMMY/2.0D0
291:             NREPINT(J1)=NREPINT(J1)+1269:             NREPINT(J1)=NREPINT(J1)+1
293:          ELSE IF (J1.EQ.INTIMAGE+2) THEN271:          ELSE IF (J1.EQ.INTIMAGE+2) THEN
294:             EEE(J1-1)=EEE(J1-1)+DUMMY272:             EEE(J1-1)=EEE(J1-1)+DUMMY
295:             REPEINT(J1-1)=REPEINT(J1-1)+DUMMY273:             REPEINT(J1-1)=REPEINT(J1-1)+DUMMY
296:             NREPINT(J1-1)=NREPINT(J1-1)+1274:             NREPINT(J1-1)=NREPINT(J1-1)+1
297:          ENDIF275:          ENDIF
298:          EREP=EREP+DUMMY276:          EREP=EREP+DUMMY
299:          IF (DUMMY.GT.EMAX) THEN277:          IF (DUMMY.GT.EMAX) THEN
300:             IMAX=J1278:             IMAX=J1
301:             JMAX=J2279:             JMAX=J2
302:             EMAX=DUMMY280:             EMAX=DUMMY
303:             EMAX=DUMMY 
304:             DISTMAX=DINT281:             DISTMAX=DINT
 282:             DIST2MAX=D12
305:             CUTMAX=NREPCUT(J2)283:             CUTMAX=NREPCUT(J2)
306:          ENDIF284:          ENDIF
307:          IF (DUMMY.GT.EEMAX(J1)) THEN285:          IF (DUMMY.GT.EEMAX(J1)) THEN
308:             JJMAX(J1)=J2286:             JJMAX(J1)=J2
309:             EEMAX(J1)=DUMMY287:             EEMAX(J1)=DUMMY
310:          ENDIF288:          ENDIF
311: ! !        DUMMY=-12.0D0*INTCONSTRAINTREP*(1.0D0/(DINT*D12)-1.0D0/INTCONST)289: !        DUMMY=-12.0D0*INTCONSTRAINTREP*(1.0D0/(DINT*D12)-1.0D0/INTCONST)
312: !        DUMMY=-2.0D0*INTCONSTRAINTREP*(1.0D0/(DINT*D12)-1.0D0/INTCONST)290:          DUMMY=-2.0D0*INTCONSTRAINTREP*(1.0D0/(DINT*D12)-1.0D0/INTCONST)
313: !        DUMMY=-2.0D0*(1.0D0/(DINT*D12)-1.0D0/INTCONST)*NREPCUT(J2)**2 
314:          DUMMY=-2.0D0*(DCUT/(DINT*D12)-1.0D0/NREPCUT(J2)) 
315:          REPGRAD(1:3)=INTMINFAC*DUMMY*G1INT(1:3)291:          REPGRAD(1:3)=INTMINFAC*DUMMY*G1INT(1:3)
316: !        PRINT '(A,4I6,2G15.5)','in1 J1,J2,REPI,REPJ,REPGRAD,NREPCUT=',J1,J2,NREPI(J2),NREPJ(J2), &292: !        PRINT '(A,4I6,2G15.5)','in1 J1,J2,REPI,REPJ,REPGRAD,NREPCUT=',J1,J2,NREPI(J2),NREPJ(J2), &
317: ! &                              SQRT(REPGRAD(1)**2+REPGRAD(2)**2+REPGRAD(3)**2),NREPCUT(J2)293: ! &                              SQRT(REPGRAD(1)**2+REPGRAD(2)**2+REPGRAD(3)**2),NREPCUT(J2)
318: !294: !
319: ! Gradient contributions for image J1-1295: ! Gradient contributions for image J1-1
320: !296: !
321:          GGG(NI1+1:NI1+3)=GGG(NI1+1:NI1+3)+REPGRAD(1:3)297:          GGG(NI1+1:NI1+3)=GGG(NI1+1:NI1+3)+REPGRAD(1:3)
322:          GGG(NJ1+1:NJ1+3)=GGG(NJ1+1:NJ1+3)-REPGRAD(1:3)298:          GGG(NJ1+1:NJ1+3)=GGG(NJ1+1:NJ1+3)-REPGRAD(1:3)
323:          REPGRAD(1:3)=INTMINFAC*DUMMY*G2INT(1:3)299:          REPGRAD(1:3)=INTMINFAC*DUMMY*G2INT(1:3)
324: !        PRINT '(A,4I6,2G15.5)','in1 J1,J2,REPI,REPJ,REPGRAD,NREPCUT=',J1,J2,NREPI(J2),NREPJ(J2), &300: !        PRINT '(A,4I6,2G15.5)','in1 J1,J2,REPI,REPJ,REPGRAD,NREPCUT=',J1,J2,NREPI(J2),NREPJ(J2), &
327: ! Gradient contributions for image J1303: ! Gradient contributions for image J1
328: !304: !
329:          GGG(NI2+1:NI2+3)=GGG(NI2+1:NI2+3)+REPGRAD(1:3)305:          GGG(NI2+1:NI2+3)=GGG(NI2+1:NI2+3)+REPGRAD(1:3)
330:          GGG(NJ2+1:NJ2+3)=GGG(NJ2+1:NJ2+3)-REPGRAD(1:3)306:          GGG(NJ2+1:NJ2+3)=GGG(NJ2+1:NJ2+3)-REPGRAD(1:3)
331:       ENDIF307:       ENDIF
332:    ENDDO308:    ENDDO
333: ENDDO309: ENDDO
334: IF (JMAX.GT.0) THEN310: IF (JMAX.GT.0) THEN
335:    WRITE(*,'(A,I6,A,I6,A,2I6,A,G20.10,A,3G20.10)') ' congrad> Highest repulsive  contribution for any image in image ',IMAX, &  311:    WRITE(*,'(A,I6,A,I6,A,2I6,A,G20.10,A,3G20.10)') ' congrad> Highest repulsive  contribution for any image in image ',IMAX, &  
336:  &  ' pair index ', &312:  &  ' pair index ', &
337:  &                                JMAX,' atoms ',NREPI(JMAX),NREPJ(JMAX),' value=',EMAX,' d,cutoff=',DISTMAX,CUTMAX 313:  &                                JMAX,' atoms ',NREPI(JMAX),NREPJ(JMAX),' value=',EMAX,' d,d^2,cutoff=',DISTMAX,DIST2MAX,CUTMAX 
338: ENDIF314: ENDIF
339: CONVERGEREPTEST=EMAX/INTCONSTRAINTDEL 
340: 315: 
341: !316: !
342: ! Spring energy. Set EEE(J1) and ESPRING dividing up the pairwise317: ! Spring energy. Set EEE(J1) and ESPRING dividing up the pairwise
343: ! energy terms between images except for the end points.318: ! energy terms between images except for the end points.
344: !319: !
345: ESPRING=0.0D0320: ESPRING=0.0D0
346: IF (KINT.NE.0.0D0) THEN321: IF (KINT.NE.0.0D0) THEN
347:    DO J1=1,INTIMAGE+1 ! sum over edges from J1 to J1+1322:    DO J1=1,INTIMAGE+1 ! sum over edges from J1 to J1+1
348:       NI1=(3*NATOMS)*(J1-1)323:       NI1=(3*NATOMS)*(J1-1)
349:       NI2=(3*NATOMS)*J1324:       NI2=(3*NATOMS)*J1
373:             IF ((.NOT.INTSPRINGACTIVET).OR.ATOMACTIVE(J2)) THEN 348:             IF ((.NOT.INTSPRINGACTIVET).OR.ATOMACTIVE(J2)) THEN 
374:                SPGRAD(1:3)=DUMMY*(XYZ(NI1+3*(J2-1)+1:NI1+3*(J2-1)+3)-XYZ(NI2+3*(J2-1)+1:NI2+3*(J2-1)+3))349:                SPGRAD(1:3)=DUMMY*(XYZ(NI1+3*(J2-1)+1:NI1+3*(J2-1)+3)-XYZ(NI2+3*(J2-1)+1:NI2+3*(J2-1)+3))
375:                GGG(NI1+3*(J2-1)+1:NI1+3*(J2-1)+3)=GGG(NI1+3*(J2-1)+1:NI1+3*(J2-1)+3)+SPGRAD(1:3)350:                GGG(NI1+3*(J2-1)+1:NI1+3*(J2-1)+3)=GGG(NI1+3*(J2-1)+1:NI1+3*(J2-1)+3)+SPGRAD(1:3)
376:                GGG(NI2+3*(J2-1)+1:NI2+3*(J2-1)+3)=GGG(NI2+3*(J2-1)+1:NI2+3*(J2-1)+3)-SPGRAD(1:3)351:                GGG(NI2+3*(J2-1)+1:NI2+3*(J2-1)+3)=GGG(NI2+3*(J2-1)+1:NI2+3*(J2-1)+3)-SPGRAD(1:3)
377:             ENDIF352:             ENDIF
378:          ENDDO353:          ENDDO
379: !     ENDIF354: !     ENDIF
380:    ENDDO355:    ENDDO
381: ENDIF356: ENDIF
382: WRITE(*,'(A,I6,A,I6,A,2I6)') ' congrad> Highest spring  contribution for any image in image ',IMAX357: WRITE(*,'(A,I6,A,I6,A,2I6)') ' congrad> Highest spring  contribution for any image in image ',IMAX
383: IF (DEBUG) WRITE(*, '(A,3G20.10)') 'congrad> ECON,EREP,ESPRING=',ECON,EREP,ESPRING358: IF (DEBUG) WRITE(*, '(A,3G20.10)') 'congrad2> ECON,EREP,ESPRING=',ECON,EREP,ESPRING
384: !359: !
385: ! Set gradients on frozen atoms to zero.360: ! Set gradients on frozen atoms to zero.
386: !361: !
387: IF (FREEZE) THEN362: IF (FREEZE) THEN
388:    DO J1=2,INTIMAGE+1  363:    DO J1=2,INTIMAGE+1  
389:       DO J2=1,NATOMS364:       DO J2=1,NATOMS
390:          IF (FROZEN(J2)) THEN365:          IF (FROZEN(J2)) THEN
391:             GGG((3*NATOMS)*(J1-1)+3*(J2-1)+1)=0.0D0366:             GGG((3*NATOMS)*(J1-1)+3*(J2-1)+1)=0.0D0
392:             GGG((3*NATOMS)*(J1-1)+3*(J2-1)+2)=0.0D0367:             GGG((3*NATOMS)*(J1-1)+3*(J2-1)+2)=0.0D0
393:             GGG((3*NATOMS)*(J1-1)+3*(J2-1)+3)=0.0D0368:             GGG((3*NATOMS)*(J1-1)+3*(J2-1)+3)=0.0D0
547: !522: !
548: ! This version of congrad tests for an internal minimum in the523: ! This version of congrad tests for an internal minimum in the
549: ! constraint distances as well as the repulsions.524: ! constraint distances as well as the repulsions.
550: !525: !
551: SUBROUTINE CONGRAD2(NMAXINT,NMININT,ETOTAL,XYZ,GGG,EEE,IMGFREEZE,RMS)526: SUBROUTINE CONGRAD2(NMAXINT,NMININT,ETOTAL,XYZ,GGG,EEE,IMGFREEZE,RMS)
552: USE KEY, ONLY: FROZEN, FREEZE, NREPI, NREPJ, NNREPULSIVE, &527: USE KEY, ONLY: FROZEN, FREEZE, NREPI, NREPJ, NNREPULSIVE, &
553:   &            NCONSTRAINT, CONI, CONJ, INTCONSTRAINTDEL, CONDISTREF, INTCONSTRAINTREP, CONDISTREFLOCAL, &528:   &            NCONSTRAINT, CONI, CONJ, INTCONSTRAINTDEL, CONDISTREF, INTCONSTRAINTREP, CONDISTREFLOCAL, &
554:   &            CONACTIVE, INTCONSTRAINREPCUT, NREPCUT,FREEZENODEST, INTIMAGE, ATOMACTIVE, KINT, IMSEPMAX, &529:   &            CONACTIVE, INTCONSTRAINREPCUT, NREPCUT,FREEZENODEST, INTIMAGE, ATOMACTIVE, KINT, IMSEPMAX, &
555:   &            INTFREEZET, INTFROZEN, REPI, REPJ, CONCUT, CONCUTLOCAL, &530:   &            INTFREEZET, INTFROZEN, REPI, REPJ, CONCUT, CONCUTLOCAL, &
556:   &            CONCUTABS, CONCUTABST, CONCUTFRAC, CONCUTFRACT, INTMINFAC, INTSPRINGACTIVET, CHECKCONINT, JMAXCON, &531:   &            CONCUTABS, CONCUTABST, CONCUTFRAC, CONCUTFRACT, INTMINFAC, INTSPRINGACTIVET, CHECKCONINT, JMAXCON, &
557:   &            NCONOFF, CONOFFTRIED, INTCONMAX, KINTENDS, ECON, EREP, ESPRING, &532:   &            NCONOFF, CONOFFTRIED, INTCONMAX, KINTENDS, ECON, EREP, ESPRING
558:   &  CONVERGECONTEST, CONVERGEREPTEST 
559: USE COMMONS, ONLY: NATOMS, NOPT, DEBUG533: USE COMMONS, ONLY: NATOMS, NOPT, DEBUG
560: IMPLICIT NONE534: IMPLICIT NONE
561:            535:            
562: INTEGER :: J1,J2,NI2,NI1,NJ2,NJ1,NMAXINT,NMININT,NCONINT(INTIMAGE+2),NREPINT(INTIMAGE+2),JMAX,IMAX,JMAXNOFF,IMAXNOFF536: INTEGER :: J1,J2,NI2,NI1,NJ2,NJ1,NMAXINT,NMININT,NCONINT(INTIMAGE+2),NREPINT(INTIMAGE+2),JMAX,IMAX,JMAXNOFF,IMAXNOFF
563: DOUBLE PRECISION :: ETOTAL, RMS, EMAX, EMAXNOFF537: DOUBLE PRECISION :: ETOTAL, RMS, EMAX, EMAXNOFF
564: INTEGER JJMAX(INTIMAGE+2)538: INTEGER JJMAX(INTIMAGE+2)
565: DOUBLE PRECISION  EEMAX(INTIMAGE+2)539: DOUBLE PRECISION  EEMAX(INTIMAGE+2)
566: DOUBLE PRECISION R1AX,R1AY,R1AZ,R2AX,R2AY,R2AZ,R1BX,R1BY,R1BZ,R2BX,R2BY,R2BZ,D2,D1540: DOUBLE PRECISION R1AX,R1AY,R1AZ,R2AX,R2AY,R2AZ,R1BX,R1BY,R1BZ,R2BX,R2BY,R2BZ,D2,D1
567: DOUBLE PRECISION G1(3),G2(3),DINT,G1INT(3),G2INT(3)541: DOUBLE PRECISION G1(3),G2(3),DINT,G1INT(3),G2INT(3)
568: DOUBLE PRECISION DUMMY, REPGRAD(3), INTCONST, D12, DSQ2, DSQ1, DSQI542: DOUBLE PRECISION DUMMY, REPGRAD(3), INTCONST, D12, DSQ2, DSQ1, DSQI
803: ! &                              SQRT(REPGRAD(1)**2+REPGRAD(2)**2+REPGRAD(3)**2)777: ! &                              SQRT(REPGRAD(1)**2+REPGRAD(2)**2+REPGRAD(3)**2)
804:          GGG(NI2+1:NI2+3)=GGG(NI2+1:NI2+3)+REPGRAD(1:3)778:          GGG(NI2+1:NI2+3)=GGG(NI2+1:NI2+3)+REPGRAD(1:3)
805:          GGG(NJ2+1:NJ2+3)=GGG(NJ2+1:NJ2+3)-REPGRAD(1:3)779:          GGG(NJ2+1:NJ2+3)=GGG(NJ2+1:NJ2+3)-REPGRAD(1:3)
806: !        WRITE(*, '(A,3I7,9G13.5)') 'J1,NI2,NJ2,INTMINFAC,INTCONSTRAINTDEL,DUMMY,GGG=',J1,NI2,NJ2,INTMINFAC,INTCONSTRAINTDEL, & 780: !        WRITE(*, '(A,3I7,9G13.5)') 'J1,NI2,NJ2,INTMINFAC,INTCONSTRAINTDEL,DUMMY,GGG=',J1,NI2,NJ2,INTMINFAC,INTCONSTRAINTDEL, & 
807: ! &            DUMMY, GGG(NI2+1:NI2+3),GGG(NJ2+1:NJ2+3)   781: ! &            DUMMY, GGG(NI2+1:NI2+3),GGG(NJ2+1:NJ2+3)   
808: !        WRITE(*,'(A,2G20.10)') 'in intmin block EEE(J1),EEE(J1-1)=',EEE(J1),EEE(J1-1)782: !        WRITE(*,'(A,2G20.10)') 'in intmin block EEE(J1),EEE(J1-1)=',EEE(J1),EEE(J1-1)
809:       ENDIF783:       ENDIF
810:    ENDDO784:    ENDDO
811: ENDDO785: ENDDO
812: IF (JMAX.GT.0) THEN786: IF (JMAX.GT.0) THEN
813:    WRITE(*,'(A,I6,A,I6,A,2I6,A,G20.10,A,3G20.10)') ' congrad> Highest constraint contribution for any image in image ',IMAX, &787:    WRITE(*,'(A,I6,A,I6,A,2I6,A,3G20.10)') ' congrad2> Highest constraint contribution for any image in image ',IMAX, &
814:  & ' constraint ',JMAX, &788:  & ' constraint ',JMAX, &
815:  &                              ' atoms ',CONI(JMAX),CONJ(JMAX),' value=',EMAX,' d,ref,cutoff=',CONDMAX,CONREFMAX,CONCUTMAX789:  &                              ' atoms ',CONI(JMAX),CONJ(JMAX),' d,ref,cutoff=',CONDMAX,CONREFMAX,CONCUTMAX
816: 790: 
817: ENDIF791: ENDIF
818: CONVERGECONTEST=EMAX/INTCONSTRAINTDEL 
819: IF (JMAXNOFF.GT.0) THEN792: IF (JMAXNOFF.GT.0) THEN
820:    WRITE(*,'(A,I6,A,I6,A,2I8,A,I8)') ' congrad2> Highest constraint contribution never turned off for any image in image ',IMAXNOFF, &793:    WRITE(*,'(A,I6,A,I6,A,2I8,A,I8)') ' congrad2> Highest constraint contribution never turned off for any image in image ',IMAXNOFF, &
821:  & ' constraint ',JMAXNOFF, &794:  & ' constraint ',JMAXNOFF, &
822:  &                              ' atoms ',CONI(JMAX),CONJ(JMAX),' off=',NCONOFF795:  &                              ' atoms ',CONI(JMAX),CONJ(JMAX),' off=',NCONOFF
823: ELSEIF (JMAX.GT.0) THEN796: ELSEIF (JMAX.GT.0) THEN
824:    JMAXNOFF=JMAX797:    JMAXNOFF=JMAX
825:    WRITE(*,'(A,I6,A,I6,A,2I8,A,I8)') ' congrad2> *** Using highest constraint contribution for any image, setting NCONOFF to 0'798:    WRITE(*,'(A,I6,A,I6,A,2I8,A,I8)') ' congrad2> *** Using highest constraint contribution for any image, setting NCONOFF to 0'
826:    CONOFFTRIED(1:INTCONMAX)=.FALSE.799:    CONOFFTRIED(1:INTCONMAX)=.FALSE.
827:    NCONOFF=0800:    NCONOFF=0
828: ENDIF801: ENDIF
1011:          GGG(NI2+1:NI2+3)=GGG(NI2+1:NI2+3)+REPGRAD(1:3)984:          GGG(NI2+1:NI2+3)=GGG(NI2+1:NI2+3)+REPGRAD(1:3)
1012:          GGG(NJ2+1:NJ2+3)=GGG(NJ2+1:NJ2+3)-REPGRAD(1:3)985:          GGG(NJ2+1:NJ2+3)=GGG(NJ2+1:NJ2+3)-REPGRAD(1:3)
1013:       ENDIF986:       ENDIF
1014:    ENDDO987:    ENDDO
1015: ENDDO988: ENDDO
1016: IF (JMAX.GT.0) THEN989: IF (JMAX.GT.0) THEN
1017:    WRITE(*,'(A,I6,A,I6,A,2I6)') ' congrad2> Highest repulsive  contribution for any image in image ',IMAX, &990:    WRITE(*,'(A,I6,A,I6,A,2I6)') ' congrad2> Highest repulsive  contribution for any image in image ',IMAX, &
1018:  &  ' pair index ', &991:  &  ' pair index ', &
1019:  &                                JMAX,' atoms ',NREPI(JMAX),NREPJ(JMAX)992:  &                                JMAX,' atoms ',NREPI(JMAX),NREPJ(JMAX)
1020: ENDIF993: ENDIF
1021: CONVERGEREPTEST=EMAX/INTCONSTRAINTREP 
1022: !994: !
1023: ! Spring energy. Set EEE(J1) and ESPRING dividing up the pairwise995: ! Spring energy. Set EEE(J1) and ESPRING dividing up the pairwise
1024: ! energy terms between images except for the end points.996: ! energy terms between images except for the end points.
1025: !997: !
1026: ESPRING=0.0D0998: ESPRING=0.0D0
1027: EMAX=0.0D0999: EMAX=0.0D0
1028: IMAX=01000: IMAX=0
1029: IF (KINT.NE.0.0D0) THEN1001: IF (KINT.NE.0.0D0) THEN
1030:    DO J1=1,INTIMAGE+1 ! sum over edges from J1 to J1+11002:    DO J1=1,INTIMAGE+1 ! sum over edges from J1 to J1+1
1031:       NI1=(3*NATOMS)*(J1-1)1003:       NI1=(3*NATOMS)*(J1-1)


r33495/intlbfgs.f90 2017-11-17 17:30:12.557108261 +0000 r33494/intlbfgs.f90 2017-11-17 17:30:13.233117222 +0000
 24:      & CONDISTREFLOCAL, INTCONFRAC, CONACTIVE, REPI, & 24:      & CONDISTREFLOCAL, INTCONFRAC, CONACTIVE, REPI, &
 25:      & REPJ, NREPMAX, ATOMACTIVE, NCONSTRAINTON, CONION, CONJON, CONDISTREFLOCALON, CONDISTREFON, & 25:      & REPJ, NREPMAX, ATOMACTIVE, NCONSTRAINTON, CONION, CONJON, CONDISTREFLOCALON, CONDISTREFON, &
 26:      & NREPCUT, REPCUT, CHECKCONINT, INTCONSTEPS, INTRELSTEPS, MAXCONE, COLDFUSIONLIMIT, & 26:      & NREPCUT, REPCUT, CHECKCONINT, INTCONSTEPS, INTRELSTEPS, MAXCONE, COLDFUSIONLIMIT, &
 27:      & INTSTEPS1, DUMPINTXYZ, DUMPINTXYZFREQ, DUMPINTEOS, DUMPINTEOSFREQ, & 27:      & INTSTEPS1, DUMPINTXYZ, DUMPINTXYZFREQ, DUMPINTEOS, DUMPINTEOSFREQ, &
 28:      & IMSEPMIN, IMSEPMAX, MAXINTIMAGE, INTFREEZET, INTFREEZETOL, FREEZE, & 28:      & IMSEPMIN, IMSEPMAX, MAXINTIMAGE, INTFREEZET, INTFREEZETOL, FREEZE, &
 29:      & INTFROZEN, CHECKREPINTERVAL, NNREPULSIVE, INTFREEZEMIN, INTIMAGECHECK, & 29:      & INTFROZEN, CHECKREPINTERVAL, NNREPULSIVE, INTFREEZEMIN, INTIMAGECHECK, &
 30:      & CONCUT, CONCUTLOCAL, KINT, REPIFIX, REPJFIX, NREPULSIVEFIX, & 30:      & CONCUT, CONCUTLOCAL, KINT, REPIFIX, REPJFIX, NREPULSIVEFIX, &
 31:      & NCONSTRAINTFIX, CONIFIX, CONJFIX, QCIPERMCHECK, QCIPERMCHECKINT, BULKT, TWOD, RIGIDBODY, & 31:      & NCONSTRAINTFIX, CONIFIX, CONJFIX, QCIPERMCHECK, QCIPERMCHECKINT, BULKT, TWOD, RIGIDBODY, &
 32:      & QCIADDREP, QCIXYZ, WHOLEDNEB, QCIIMAGE, FROZEN, QCIRESTART, NPERMGROUP, NPERMSIZE, PERMGROUP, NSETS, SETS, & 32:      & QCIADDREP, QCIXYZ, WHOLEDNEB, QCIIMAGE, FROZEN, QCIRESTART, NPERMGROUP, NPERMSIZE, PERMGROUP, NSETS, SETS, &
 33:      & PERMDIST, LOCALPERMCUT, QCILPERMDIST, QCIPDINT, QCIPERMCUT, QCIAMBERT, BONDS, DOBACK, & 33:      & PERMDIST, LOCALPERMCUT, QCILPERMDIST, QCIPDINT, QCIPERMCUT, QCIAMBERT, BONDS, DOBACK, &
 34:      & QCIRESET, QCIRESETINT1, QCIRESETINT2, JMAXCON, NCONOFF, EREP, ECON, ESPRING, CONVERGECONTEST, CONVERGEREPTEST 34:      & QCIRESET, QCIRESETINT1, QCIRESETINT2, JMAXCON, NCONOFF, EREP, ECON, ESPRING
 35: USE COMMONS, ONLY: NATOMS, DEBUG, PARAM1, PARAM2, PARAM3 35: USE COMMONS, ONLY: NATOMS, DEBUG, PARAM1, PARAM2, PARAM3
 36: USE MODCHARMM, ONLY : CHRMMT 36: USE MODCHARMM, ONLY : CHRMMT
 37: USE CHIRALITY 37: USE CHIRALITY
 38:  38: 
 39: IMPLICIT NONE  39: IMPLICIT NONE 
 40:  40: 
 41: DOUBLE PRECISION, INTENT(IN) :: QSTART(3*NATOMS), QFINISH(3*NATOMS)  ! The two end points 41: DOUBLE PRECISION, INTENT(IN) :: QSTART(3*NATOMS), QFINISH(3*NATOMS)  ! The two end points
 42: INTEGER D, U 42: INTEGER D, U
 43: DOUBLE PRECISION DIST, DIST2, RMAT(3,3), SUMEEE, SUMEEE2, SIGMAEEE, NEIGHBOUR_COORDS(12), CENTRE_COORDS(3) 43: DOUBLE PRECISION DIST, DIST2, RMAT(3,3), SUMEEE, SUMEEE2, SIGMAEEE, NEIGHBOUR_COORDS(12), CENTRE_COORDS(3)
 44: DOUBLE PRECISION DMAX, DF, DMIN, LOCALSTEP, ADMAX, DUMMYX, DUMMYY, DUMMYZ 44: DOUBLE PRECISION DMAX, DF, DMIN, LOCALSTEP, ADMAX, DUMMYX, DUMMYY, DUMMYZ
659: ! Don;t want to redistribute images before even taking a step, so don;t call CHECKSEP.659: ! Don;t want to redistribute images before even taking a step, so don;t call CHECKSEP.
660: ! Must call CHECKREP to initialise NNREULSIVE, NREPI, NREPJ, etc. SEGV otherwise on second cycle!660: ! Must call CHECKREP to initialise NNREULSIVE, NREPI, NREPJ, etc. SEGV otherwise on second cycle!
661: !661: !
662: ! To take BH-type steps in the QCI space, jump back here. Leave SWITCHED true.662: ! To take BH-type steps in the QCI space, jump back here. Leave SWITCHED true.
663: !663: !
664: BESTWORST=1.0D100664: BESTWORST=1.0D100
665: 9876 CONTINUE665: 9876 CONTINUE
666: CALL CHECKREP(INTIMAGE,XYZ,(3*NATOMS),0,1)666: CALL CHECKREP(INTIMAGE,XYZ,(3*NATOMS),0,1)
667: 667: 
668: 668: 
669: !               DIFF=1.0D-7669: !               DIFF=1.0D-4
670: !               PRINT*,'analytic and numerical gradients:'670: !               PRINT*,'analytic and numerical gradients:'
671: !               CALL CONGRAD2(NMAXINT,NMININT,ETOTAL,XYZ,GGG,EEE,IMGFREEZE,RMS)671: !               CALL CONGRAD2(NMAXINT,NMININT,ETOTAL,XYZ,GGG,EEE,IMGFREEZE,RMS)
672: !               DO J1=2,INTIMAGE+1672: !               DO J1=2,INTIMAGE+1
673: !                  DO J2=1,3*NATOMS673: !                  DO J2=1,3*NATOMS
674: !                     J3=3*NATOMS*(J1-1)+J2674: !                     J3=3*NATOMS*(J1-1)+J2
675: !                     XYZ(J3)=XYZ(J3)+DIFF675: !                     XYZ(J3)=XYZ(J3)+DIFF
676: !                     CALL CONGRAD2(NMAXINT,NMININT,EPLUS,XYZ,VPLUS,EEE,IMGFREEZE,RMS)676: !                     CALL CONGRAD2(NMAXINT,NMININT,EPLUS,XYZ,VPLUS,EEE,IMGFREEZE,RMS)
677: !                     XYZ(J3)=XYZ(J3)-2.0D0*DIFF677: !                     XYZ(J3)=XYZ(J3)-2.0D0*DIFF
678: !                     CALL CONGRAD2(NMAXINT,NMININT,EMINUS,XYZ,VMINUS,EEE,IMGFREEZE,RMS)678: !                     CALL CONGRAD2(NMAXINT,NMININT,EMINUS,XYZ,VMINUS,EEE,IMGFREEZE,RMS)
679: !                     XYZ(J3)=XYZ(J3)+DIFF679: !                     XYZ(J3)=XYZ(J3)+DIFF
745: !745: !
746: ! Are we stuck? If so, try resetting problem atoms to previous image.746: ! Are we stuck? If so, try resetting problem atoms to previous image.
747: !747: !
748: IF (QCIRESET) THEN748: IF (QCIRESET) THEN
749: !  IF ((SWITCHED.AND.(MOD(NITERDONE-1,QCIRESETINT2).EQ.0)).OR.((.NOT.SWITCHED).AND.(MOD(NITERDONE-1,QCIRESETINT1).EQ.0))) THEN749: !  IF ((SWITCHED.AND.(MOD(NITERDONE-1,QCIRESETINT2).EQ.0)).OR.((.NOT.SWITCHED).AND.(MOD(NITERDONE-1,QCIRESETINT1).EQ.0))) THEN
750: !  PRINT *,'intlbfgs> NITERDONE,NLASTGOODE,QCIRESETINT1=',NITERDONE,NLASTGOODE,QCIRESETINT1750: !  PRINT *,'intlbfgs> NITERDONE,NLASTGOODE,QCIRESETINT1=',NITERDONE,NLASTGOODE,QCIRESETINT1
751:    IF ((.NOT.SWITCHED).AND.(NITERDONE-NLASTGOODE.GT.QCIRESETINT1)) THEN ! .AND.(NITERDONE-NLASTCHANGE.GT.QCIRESETINT1)) THEN751:    IF ((.NOT.SWITCHED).AND.(NITERDONE-NLASTGOODE.GT.QCIRESETINT1)) THEN ! .AND.(NITERDONE-NLASTCHANGE.GT.QCIRESETINT1)) THEN
752:       IF (DEBUG) CALL INTRWG(NACTIVE,NITERDONE,INTIMAGE,XYZ,TURNONORDER,NCONOFF)752:       IF (DEBUG) CALL INTRWG(NACTIVE,NITERDONE,INTIMAGE,XYZ,TURNONORDER,NCONOFF)
753:       IF (DEBUG) CALL WRITEPROFILE(NITERDONE,EEE,INTIMAGE)753:       IF (DEBUG) CALL WRITEPROFILE(NITERDONE,EEE,INTIMAGE)
754:       TURNOFF=.TRUE.754:       TURNOFF=.TRUE.
755: !     IF (EREP.GT.2.0D0*ECON) THEN755:       IF (EREP.GT.2.0D0*ECON) THEN
756: !        INTCONSTRAINTREP=INTCONSTRAINTREP/2.0D0756:          INTCONSTRAINTREP=INTCONSTRAINTREP/2.0D0
757: !        WRITE(*,'(A,G20.10)') 'intlbfgs> Interpolation seems to be stuck. Reducing repulsive prefactor to ',INTCONSTRAINTREP757:          WRITE(*,'(A,G20.10)') 'intlbfgs> Interpolation seems to be stuck. Reducing repulsive prefactor to ',INTCONSTRAINTREP
758: !        TURNOFF=.FALSE.758:          TURNOFF=.FALSE.
759: !     ENDIF759:       ENDIF
760:       IF (ESPRING.GT.2.0D0*ECON) THEN760:       IF (ESPRING.GT.2.0D0*ECON) THEN
761:          KINT=KINT/2.0D0761:          KINT=KINT/2.0D0
762:          WRITE(*,'(A,G20.10)') 'intlbfgs> Interpolation seems to be stuck. Reducing spring constant to ',KINT762:          WRITE(*,'(A,G20.10)') 'intlbfgs> Interpolation seems to be stuck. Reducing spring constant to ',KINT
763:          TURNOFF=.FALSE.763:          TURNOFF=.FALSE.
764:       ENDIF764:       ENDIF
765:       IF (TURNOFF) THEN765:       IF (TURNOFF) THEN
766:          WRITE(*,'(2(A,I6))') 'intlbfgs> Interpolation seems to be stuck. Dump images. Turn off worst constraint ',JMAXCON, &766:          WRITE(*,'(2(A,I6))') 'intlbfgs> Interpolation seems to be stuck. Dump images. Turn off worst constraint ',JMAXCON, &
767:   &                        ' total off=',NCONOFF+1767:   &                        ' total off=',NCONOFF+1
768:          IF ((JMAXCON.LT.1).OR.(JMAXCON.GT.NCONSTRAINT)) THEN768:          IF ((JMAXCON.LT.1).OR.(JMAXCON.GT.NCONSTRAINT)) THEN
769:             WRITE(*,'(A)') 'intlbfgs> *** ERROR *** constraint index out of allowed range'769:             WRITE(*,'(A)') 'intlbfgs> *** ERROR *** constraint index out of allowed range'
1068: !!                  ENDDO myrep21068: !!                  ENDDO myrep2
1069: !!               ENDDO1069: !!               ENDDO
1070: !!               WRITE(*,'(A,I6,A)') ' intlbfgs> Now it looks like there are ',NREPULSIVE,' possible repulsions before adding new atom'1070: !!               WRITE(*,'(A,I6,A)') ' intlbfgs> Now it looks like there are ',NREPULSIVE,' possible repulsions before adding new atom'
1071: !!!!!!!!!!!!!!!DEBUG DJW !!!!!!!!!!!1071: !!!!!!!!!!!!!!!DEBUG DJW !!!!!!!!!!!
1072: 1072: 
1073:       IF (NCONOFF.GT.0) THEN1073:       IF (NCONOFF.GT.0) THEN
1074:          CONACTIVE(CONOFFLIST(NCONOFF))=.TRUE.1074:          CONACTIVE(CONOFFLIST(NCONOFF))=.TRUE.
1075:          WRITE(*,'(2(A,I6))') 'intlbfgs> Turn back on constraint ',CONOFFLIST(NCONOFF),' total off=',NCONOFF-11075:          WRITE(*,'(2(A,I6))') 'intlbfgs> Turn back on constraint ',CONOFFLIST(NCONOFF),' total off=',NCONOFF-1
1076:          NCONOFF=NCONOFF-11076:          NCONOFF=NCONOFF-1
1077:       ELSE1077:       ELSE
1078: !        IF (DEBUG) WRITE(*,'(A)') 'intlbfgs> dump state before doaddatom index -2'1078:          IF (DEBUG) WRITE(*,'(A)') 'intlbfgs> dump state before doaddatom index -2'
1079: !        IF (DEBUG) CALL INTRWG2(NACTIVE,-2,INTIMAGE,XYZ,TURNONORDER,NCONOFF)1079:          IF (DEBUG) CALL INTRWG2(NACTIVE,-2,INTIMAGE,XYZ,TURNONORDER,NCONOFF)
1080:          CALL DOADDATOM(NCONSTRAINT,NTRIES,NEWATOM,IMGFREEZE,INTIMAGE,XYZ,EEE,GGG,TURNONORDER,NITERDONE,NACTIVE,AABACK,BACKDONE)  1080:          CALL DOADDATOM(NCONSTRAINT,NTRIES,NEWATOM,IMGFREEZE,INTIMAGE,XYZ,EEE,GGG,TURNONORDER,NITERDONE,NACTIVE,AABACK,BACKDONE)  
1081: !        IF (DEBUG) WRITE(*,'(A)') 'intlbfgs> dump state after doaddatom index -1'1081:          IF (DEBUG) WRITE(*,'(A)') 'intlbfgs> dump state after doaddatom index -1'
1082: !        IF (DEBUG) CALL INTRWG2(NACTIVE,-1,INTIMAGE,XYZ,TURNONORDER,NCONOFF)1082:          IF (DEBUG) CALL INTRWG2(NACTIVE,-1,INTIMAGE,XYZ,TURNONORDER,NCONOFF)
1083:       ENDIF1083:       ENDIF
1084:       NLASTGOODE=NITERDONE1084:       NLASTGOODE=NITERDONE
1085:       LASTGOODE=ETOTAL1085:       LASTGOODE=ETOTAL
1086:    ENDIF1086:    ENDIF
1087:    GTMP(1:D)=0.0D01087:    GTMP(1:D)=0.0D0
1088:    CALL MAKESTEP(NITERUSE,POINT,DIAG,INTIMAGE,SEARCHSTEP,G,GTMP,STP,GDIF,NPT,D,RHO1,ALPHA)1088:    CALL MAKESTEP(NITERUSE,POINT,DIAG,INTIMAGE,SEARCHSTEP,G,GTMP,STP,GDIF,NPT,D,RHO1,ALPHA)
1089: !1089: !
1090: ! If the number of images has changed since G was declared then G is not the same1090: ! If the number of images has changed since G was declared then G is not the same
1091: ! size as Gtmp and Dot_Product cannot be used.1091: ! size as Gtmp and Dot_Product cannot be used.
1092: !1092: !
1540: !    IF (ABS(MAXEEE-SUMEEE).GT.3.0D0*SIGMAEEE) THEN1540: !    IF (ABS(MAXEEE-SUMEEE).GT.3.0D0*SIGMAEEE) THEN
1541: ! !  IF (MAXEEE.GT.1.0D2*MINEEE) THEN1541: ! !  IF (MAXEEE.GT.1.0D2*MINEEE) THEN
1542:       WRITE(*,'(A,I8,A,G20.10,A,G20.10,A)') ' intlbfgs> Highest image ',JMAXEEE,' energy ',MAXEEE,' is ',ABS(MAXEEE-SUMEEE)/SIGMAEEE, &1542:       WRITE(*,'(A,I8,A,G20.10,A,G20.10,A)') ' intlbfgs> Highest image ',JMAXEEE,' energy ',MAXEEE,' is ',ABS(MAXEEE-SUMEEE)/SIGMAEEE, &
1543:   &                        ' sigma from the mean'1543:   &                        ' sigma from the mean'
1544: !       REMOVEIMAGE=.TRUE.1544: !       REMOVEIMAGE=.TRUE.
1545: !       IF (NITERDONE-NLASTGOODE.LT.20) REMOVEIMAGE=.FALSE.1545: !       IF (NITERDONE-NLASTGOODE.LT.20) REMOVEIMAGE=.FALSE.
1546: !    ENDIF1546: !    ENDIF
1547: 1547: 
1548:    IF (DEBUG) THEN1548:    IF (DEBUG) THEN
1549: !     WRITE(*,'(A,2G20.10)') ' intlbfgs> mean and sigma of image energies=',SUMEEE,SIGMAEEE1549: !     WRITE(*,'(A,2G20.10)') ' intlbfgs> mean and sigma of image energies=',SUMEEE,SIGMAEEE
1550: !     WRITE(*,'(A,I6,2G20.10,3(G20.10,I8))') ' intlbfgs> steps: ',NITERDONE,ETOTAL/INTIMAGE,RMS,STEPTOT,NACTIVE, &1550:       WRITE(*,'(A,I6,2G20.10,3(G20.10,I8))') ' intlbfgs> steps: ',NITERDONE,ETOTAL/INTIMAGE,RMS,STEPTOT,NACTIVE, &
1551: ! &                                                        MAXEEE,JMAXEEE,MAXRMS,JMAXRMS1551:   &                                                        MAXEEE,JMAXEEE,MAXRMS,JMAXRMS
1552:       WRITE(*,'(A,I6,4G20.10,I10,I6)') ' intlbfgs> steps: ',NITERDONE,CONVERGECONTEST,CONVERGEREPTEST,RMS,STEPTOT,NACTIVE,INTIMAGE+2 
1553:       CALL FLUSH(6)1552:       CALL FLUSH(6)
1554:    ENDIF1553:    ENDIF
1555: 1554: 
1556:    IF (.NOT.SWITCHED) THEN1555:    IF (.NOT.SWITCHED) THEN
1557: !     IF ((NITERDONE-NLASTGOODE.GT.INTRELSTEPS).AND.((ETOTAL.GT.LASTGOODE).OR.(ETOTAL/INTIMAGE.GT.MAXCONE*1.0D8))) THEN1556: !     IF ((NITERDONE-NLASTGOODE.GT.INTRELSTEPS).AND.((ETOTAL.GT.LASTGOODE).OR.(ETOTAL/INTIMAGE.GT.MAXCONE*1.0D8))) THEN
1558:       IF (.FALSE.) THEN ! no backtracking1557:       IF (.FALSE.) THEN ! no backtracking
1559:          WRITE(*,'(2(A,I6))') ' intlbfgs> Backtracking ',NBACKTRACK,' steps, current active atoms=',NACTIVE1558:          WRITE(*,'(2(A,I6))') ' intlbfgs> Backtracking ',NBACKTRACK,' steps, current active atoms=',NACTIVE
1560:          NTRIES(NEWATOM)=NTRIES(NEWATOM)+11559:          NTRIES(NEWATOM)=NTRIES(NEWATOM)+1
1561:          IF (FREEZENODEST) IMGFREEZE(1:INTIMAGE)=.FALSE.1560:          IF (FREEZENODEST) IMGFREEZE(1:INTIMAGE)=.FALSE.
1562: !1561: !
1656:             CALL CONGRAD2(NMAXINT,NMININT,ETOTAL,XYZ,GGG,EEE,IMGFREEZE,RMS)1655:             CALL CONGRAD2(NMAXINT,NMININT,ETOTAL,XYZ,GGG,EEE,IMGFREEZE,RMS)
1657:          ELSE1656:          ELSE
1658:             CALL CONGRAD(NMAXINT,NMININT,ETOTAL,XYZ,GGG,EEE,IMGFREEZE,RMS)1657:             CALL CONGRAD(NMAXINT,NMININT,ETOTAL,XYZ,GGG,EEE,IMGFREEZE,RMS)
1659:          ENDIF1658:          ENDIF
1660:       ENDIF1659:       ENDIF
1661:       LASTGOODE=ETOTAL1660:       LASTGOODE=ETOTAL
1662:    ENDIF1661:    ENDIF
1663: 1662: 
1664:    EXITSTATUS=01663:    EXITSTATUS=0
1665:    INTDGUESS=DIAG(1) ! should be ok for subsequent runs of the same system DJW1664:    INTDGUESS=DIAG(1) ! should be ok for subsequent runs of the same system DJW
1666:    IF ((.NOT.SWITCHED).AND.(MAXRMS<=INTRMSTOL).AND.(NITERDONE>1).AND.(CONVERGECONTEST.LT.MAXCONE).AND.(CONVERGEREPTEST.LT.MAXCONE)) EXITSTATUS=1 1665:    IF ((.NOT.SWITCHED).AND.(MAXRMS<=INTRMSTOL).AND.NITERDONE>1) EXITSTATUS=1 
1667:    IF (SWITCHED.AND.(MAXRMS<=CONVR).AND.NITERDONE>1) EXITSTATUS=1 1666:    IF (SWITCHED.AND.(MAXRMS<=CONVR).AND.NITERDONE>1) EXITSTATUS=1 
1668:    IF (NITERDONE==NSTEPSMAX) EXITSTATUS=21667:    IF (NITERDONE==NSTEPSMAX) EXITSTATUS=2
1669:    IF ((.NOT.SWITCHED).AND.(MOD(NITERDONE,INTRELSTEPS).EQ.0)) EXITSTATUS=1 ! Add an atom every INTRELSTEPS !!! DJW1668:    IF ((.NOT.SWITCHED).AND.(MOD(NITERDONE,INTRELSTEPS).EQ.0)) EXITSTATUS=1 ! Add an atom every INTRELSTEPS !!! DJW
1670: !  PRINT '(A,2G20.10,3I8)','MAXRMS,INTRMSTOL,NITERDONE,NITERDONE,NSTEPSMAX=',MAXRMS,INTRMSTOL,NITERDONE,NITERDONE,NSTEPSMAX1669: !  PRINT '(A,2G20.10,3I8)','MAXRMS,INTRMSTOL,NITERDONE,NITERDONE,NSTEPSMAX=',MAXRMS,INTRMSTOL,NITERDONE,NITERDONE,NSTEPSMAX
1671: 1670: 
1672:    IF (EXITSTATUS > 0) THEN  1671:    IF (EXITSTATUS > 0) THEN  
1673:       IF ((.NOT.SWITCHED).AND.(EXITSTATUS.EQ.1)) THEN ! add active atom or restart with true potential on1672:       IF ((.NOT.SWITCHED).AND.(EXITSTATUS.EQ.1)) THEN ! add active atom or restart with true potential on
1674: !        IF (ETOTAL/INTIMAGE.GT.MAXCONE*MAX(0.3D0,NACTIVE*1.0D0/(NATOMS*1.0D0))) GOTO 7771673: !        IF (ETOTAL/INTIMAGE.GT.MAXCONE*MAX(0.3D0,NACTIVE*1.0D0/(NATOMS*1.0D0))) GOTO 777
1675: !        PRINT '(A,3G20.10)','MAXEEE,MAXCONE,scaled=',MAXEEE,MAXCONE,MAXCONE*MAX(0.3D0,NACTIVE*1.0D0/(NATOMS*1.0D0))1674:          PRINT '(A,3G20.10)','MAXEEE,MAXCONE,scaled=',MAXEEE,MAXCONE,MAXCONE*MAX(0.3D0,NACTIVE*1.0D0/(NATOMS*1.0D0))
1676: !        IF (MAXEEE.GT.MAXCONE*MAX(0.3D0,NACTIVE*1.0D0/(NATOMS*1.0D0))) GOTO 7771675:          IF (MAXEEE.GT.MAXCONE*MAX(0.3D0,NACTIVE*1.0D0/(NATOMS*1.0D0))) GOTO 777
1677:          IF (NACTIVE.LT.NATOMS) THEN 1676:          IF (NACTIVE.LT.NATOMS) THEN 
1678:             ADDATOM=.TRUE.1677:             ADDATOM=.TRUE.
1679:             GOTO 7771678:             GOTO 777
1680:          ENDIF1679:          ENDIF
1681:          CALL MYCPU_TIME(FTIME,.FALSE.)1680:          CALL MYCPU_TIME(FTIME,.FALSE.)
1682:          WRITE(*,'(A,I6,A,F12.6,A,I6,A,G20.10)') ' intlbfgs> switch on true potential at step ',NITERDONE, &1681:          WRITE(*,'(A,I6,A,F12.6,A,I6,A,G20.10)') ' intlbfgs> switch on true potential at step ',NITERDONE, &
1683:   &                                     ' fraction=',INTCONFRAC,' images=',INTIMAGE,' time=',FTIME-STIME1682:   &                                     ' fraction=',INTCONFRAC,' images=',INTIMAGE,' time=',FTIME-STIME
1684:          CALL INTRWG(NACTIVE,NITERDONE,INTIMAGE,XYZ,TURNONORDER,NCONOFF)1683:          CALL INTRWG(NACTIVE,NITERDONE,INTIMAGE,XYZ,TURNONORDER,NCONOFF)
1685:          CALL WRITEPROFILE(NITERDONE,EEE,INTIMAGE)1684:          CALL WRITEPROFILE(NITERDONE,EEE,INTIMAGE)
1686:          WRITE(*,'(A,I6,A,F15.6)') ' intlbfgs> Allowing ',INTCONSTEPS,' further optimization steps'1685:          WRITE(*,'(A,I6,A,F15.6)') ' intlbfgs> Allowing ',INTCONSTEPS,' further optimization steps'
1711: ! Compute the new step and gradient change1710: ! Compute the new step and gradient change
1712: !1711: !
1713:    NPT=POINT*D1712:    NPT=POINT*D
1714:    SEARCHSTEP(POINT,:) = STP*SEARCHSTEP(POINT,:)1713:    SEARCHSTEP(POINT,:) = STP*SEARCHSTEP(POINT,:)
1715:    GDIF(POINT,:)=G-GTMP1714:    GDIF(POINT,:)=G-GTMP
1716:    1715:    
1717:    POINT=POINT+1; IF (POINT==INTMUPDATE) POINT=01716:    POINT=POINT+1; IF (POINT==INTMUPDATE) POINT=0
1718: 1717: 
1719:    IF (DUMPINTXYZ.AND.MOD(NITERDONE,DUMPINTXYZFREQ)==0) CALL INTRWG(NACTIVE,NITERDONE,INTIMAGE,XYZ,TURNONORDER,NCONOFF)1718:    IF (DUMPINTXYZ.AND.MOD(NITERDONE,DUMPINTXYZFREQ)==0) CALL INTRWG(NACTIVE,NITERDONE,INTIMAGE,XYZ,TURNONORDER,NCONOFF)
1720:    IF (DUMPINTEOS.AND.MOD(NITERDONE,DUMPINTEOSFREQ)==0) CALL WRITEPROFILE(NITERDONE,EEE,INTIMAGE)1719:    IF (DUMPINTEOS.AND.MOD(NITERDONE,DUMPINTEOSFREQ)==0) CALL WRITEPROFILE(NITERDONE,EEE,INTIMAGE)
1721: !  IF (NITERDONE.GT.0) CALL INTRWG2(NACTIVE,NITERDONE,INTIMAGE,XYZ,TURNONORDER,NCONOFF) !!! DEBUG DJW1720:    IF (NITERDONE.GT.0) CALL INTRWG2(NACTIVE,NITERDONE,INTIMAGE,XYZ,TURNONORDER,NCONOFF) !!! DEBUG DJW
1722: 1721: 
1723:    NITERDONE=NITERDONE+11722:    NITERDONE=NITERDONE+1
1724:    NITERUSE=NITERUSE+11723:    NITERUSE=NITERUSE+1
1725: 1724: 
1726:    IF (NITERDONE.GT.NSTEPSMAX) EXIT1725:    IF (NITERDONE.GT.NSTEPSMAX) EXIT
1727:    IF (NACTIVE.EQ.NATOMS) THEN1726:    IF (NACTIVE.EQ.NATOMS) THEN
1728:       IF (.NOT.SWITCHED) THEN1727:       IF (.NOT.SWITCHED) THEN
1729:          CALL MYCPU_TIME(FTIME,.FALSE.)1728:          CALL MYCPU_TIME(FTIME,.FALSE.)
1730:          WRITE(*,'(A,I6,A,F12.6,A,I6,A,F10.1)') ' intlbfgs> switch on true potential at step ',NITERDONE, &1729:          WRITE(*,'(A,I6,A,F12.6,A,I6,A,F10.1)') ' intlbfgs> switch on true potential at step ',NITERDONE, &
1731:   &                                     ' fraction=',INTCONFRAC,' images=',INTIMAGE,' time=',FTIME-STIME1730:   &                                     ' fraction=',INTCONFRAC,' images=',INTIMAGE,' time=',FTIME-STIME


r33495/key.f90 2017-11-17 17:30:12.777111179 +0000 r33494/key.f90 2017-11-17 17:30:13.453120140 +0000
108:      &        CPCONSTRAINTTOL, CPCONSTRAINTDEL, CPCONSTRAINTREP, CPCONSTRAINREPCUT, CPCONFRAC, &108:      &        CPCONSTRAINTTOL, CPCONSTRAINTDEL, CPCONSTRAINTREP, CPCONSTRAINREPCUT, CPCONFRAC, &
109:      &        INTLJTOL, INTLJDEL, INTLJEPS, IMSEPMIN, IMSEPMAX, TRAPK, MINOVERLAP, &109:      &        INTLJTOL, INTLJDEL, INTLJEPS, IMSEPMIN, IMSEPMAX, TRAPK, MINOVERLAP, &
110:      &        INTFREEZETOL, LOCALPERMCUT, LOCALPERMCUT2, LOCALPERMCUTINC, CHECKREPCUTOFF, CONCUTABS, &110:      &        INTFREEZETOL, LOCALPERMCUT, LOCALPERMCUT2, LOCALPERMCUTINC, CHECKREPCUTOFF, CONCUTABS, &
111:      &        CONCUTFRAC, ENDNUMHESSDELTA, DNEBEFRAC, QCHEMSCALE, KAA, SIGMAAA, QUIPATOMMASS, TEMPERATURE1, &111:      &        CONCUTFRAC, ENDNUMHESSDELTA, DNEBEFRAC, QCHEMSCALE, KAA, SIGMAAA, QUIPATOMMASS, TEMPERATURE1, &
112:      &        DISTORTINST,DELTAINST,MOLPROSCALE,COVER,STTSRMSCONV,LAN_DIST,LANCONV,LANFACTOR, &112:      &        DISTORTINST,DELTAINST,MOLPROSCALE,COVER,STTSRMSCONV,LAN_DIST,LANCONV,LANFACTOR, &
113:      &        STOCKEXP, JPARAM, MCPATHTEMP, MCPATHDMAX, MCPATHSTEP, MCPATHACCRATIO, BIASFAC, &113:      &        STOCKEXP, JPARAM, MCPATHTEMP, MCPATHDMAX, MCPATHSTEP, MCPATHACCRATIO, BIASFAC, &
114:      &        MCADDDEV, MCPATHQMIN, MCPATHQMAX, RPHQMIN, RPHQMAX, RPHTEMP, TWISTF, TWISTREF, MCPATHADDREF, &114:      &        MCADDDEV, MCPATHQMIN, MCPATHQMAX, RPHQMIN, RPHQMAX, RPHTEMP, TWISTF, TWISTREF, MCPATHADDREF, &
115:      &        MCPATHGWS, MCPATHGWQ, MCPATHNEGLECT, MCPATHTOL, FRAMESDIFF,TMRATIO, INTMINFAC, MLPLAMBDA, COLL_TOL, KLIM, SCA, &115:      &        MCPATHGWS, MCPATHGWQ, MCPATHNEGLECT, MCPATHTOL, FRAMESDIFF,TMRATIO, INTMINFAC, MLPLAMBDA, COLL_TOL, KLIM, SCA, &
116:      &        NEBMAXERISE, GDSQ, FLATEDIFF, QCIADDREPCUT, QCIADDREPEPS, QCIRADSHIFT, INTCONCUT, MLQLAMBDA, FRQCONV, FRQCONV2, &116:      &        NEBMAXERISE, GDSQ, FLATEDIFF, QCIADDREPCUT, QCIADDREPEPS, QCIRADSHIFT, INTCONCUT, MLQLAMBDA, FRQCONV, FRQCONV2, &
117:      &        TANHFAC, LJADDCUTOFF,LJADDREFNORM,MAXIMFACTOR, PUSHOPTCONV, QCICYCDIST, QCIPERMCUT, KERNELWIDTH, QUIPLATT(3,3), &117:      &        TANHFAC, LJADDCUTOFF,LJADDREFNORM,MAXIMFACTOR, PUSHOPTCONV, QCICYCDIST, QCIPERMCUT, KERNELWIDTH, QUIPLATT(3,3), &
118:      &        ECON, EREP, ESPRING, CONVERGECONTEST, CONVERGEREPTEST118:      &        ECON, EREP, ESPRING
119: 119: 
120: !     sf344120: !     sf344
121:       DOUBLE PRECISION :: PCUTOFF,PYA11(3),PYA21(3),PYA12(3),PYA22(3),PEPSILON1(3),PSCALEFAC1(2),PSCALEFAC2(2), &121:       DOUBLE PRECISION :: PCUTOFF,PYA11(3),PYA21(3),PYA12(3),PYA22(3),PEPSILON1(3),PSCALEFAC1(2),PSCALEFAC2(2), &
122:      &                     PEPSILONATTR(2),PSIGMAATTR(2), PYOVERLAPTHRESH, LJSITECOORDS(3), LJGSITESIGMA, LJGSITEEPS, &122:      &                     PEPSILONATTR(2),PSIGMAATTR(2), PYOVERLAPTHRESH, LJSITECOORDS(3), LJGSITESIGMA, LJGSITEEPS, &
123:      &                     PYLOCALSTEP(2),PYCFTHRESH,PYGRAVITYC1,PYGRAVITYC2,PERCCUT123:      &                     PYLOCALSTEP(2),PYCFTHRESH,PYGRAVITYC1,PYGRAVITYC2,PERCCUT
124:  124:  
125:       DOUBLE PRECISION, ALLOCATABLE :: PYADDEPS(:,:)125:       DOUBLE PRECISION, ALLOCATABLE :: PYADDEPS(:,:)
126:       DOUBLE PRECISION, ALLOCATABLE :: SHIFTL(:)126:       DOUBLE PRECISION, ALLOCATABLE :: SHIFTL(:)
127:       LOGICAL, ALLOCATABLE :: uniaxarray(:)127:       LOGICAL, ALLOCATABLE :: uniaxarray(:)
128:       DOUBLE PRECISION, ALLOCATABLE :: PYADDREP(:,:), PYADDATT(:,:)128:       DOUBLE PRECISION, ALLOCATABLE :: PYADDREP(:,:), PYADDATT(:,:)


legend
Lines Added 
Lines changed
 Lines Removed

hdiff - version: 2.1.0