hdiff output

r30147/climber.f90 2016-03-16 18:34:19.683535112 +0000 r30146/climber.f90 2016-03-16 18:34:21.251551181 +0000
  1:   1: 
  2:   2: 
  3: SUBROUTINE CLIMBERINT(QSTART,QEND)  3: SUBROUTINE CLIMBERINT(QSTART,QEND)
  4: USE PORFUNCS  4: USE PORFUNCS
  5: USE CLIMBERCOMMONS  5: USE CLIMBERCOMMONS
  6: USE COMMONS, ONLY: NATOMS,DEBUG,PARAM1, PARAM2, PARAM3, NOPT   6: USE COMMONS, ONLY: NATOMS,DEBUG,PARAM1, PARAM2, PARAM3, NOPT 
  7: USE KEY, ONLY: BULKT, TWOD, RIGIDBODY, MUPDATE, BFGSSTEPS, CLIMBERSTEPS, CLIMBERSPRING, CLIMBERCONV  7: USE KEY, ONLY: BULKT, TWOD, RIGIDBODY, MUPDATE, BFGSSTEPS, CLIMBERSTEPS, CLIMBERSPRING, CLIMBERCONV
  8: IMPLICIT NONE  8: IMPLICIT NONE
  9: DOUBLE PRECISION :: QSTART(3*NATOMS),QEND(3*NATOMS), NEWXYZ(3*NATOMS)  9: DOUBLE PRECISION :: QSTART(3*NATOMS),QEND(3*NATOMS), NEWXYZ(3*NATOMS)
 10: DOUBLE PRECISION :: STARTDIST , DIST2, RMAT(3,3) , DISTANCE,  EDUMMY, RMS2, RMS, EREAL, VNEW(NOPT) 10: DOUBLE PRECISION :: STARTDIST , DIST2, RMAT(3,3) , DISTANCE,  EDUMMY, RMS2, RMS, EREAL, VNEW(NOPT)
 11: INTEGER :: NSTRUCTURES,ITDONE,NCYCLES 11: INTEGER :: NSTRUCTURES,ITDONE
 12: LOGICAL :: PTEST, MFLAG 12: LOGICAL :: PTEST, MFLAG
 13: DOUBLE PRECISION, PARAMETER :: RED_SPRING=0.95,INC_SPRING=1.05 13: DOUBLE PRECISION, PARAMETER :: RED_SPRING=0.5,INC_SPRING=1.5
 14: LOGICAL KNOWE, KNOWG, KNOWH 14: 
 15: COMMON /KNOWN/ KNOWE, KNOWG, KNOWH 
 16:  
 17: KNOWE=.FALSE. 
 18: KNOWG=.FALSE. 
 19: KNOWH=.FALSE. 
 20: PTEST=.FALSE. 15: PTEST=.FALSE.
 21: ALLOCATE(SAVEDXYZ(2*CLIMBERSTEPS*(3*NATOMS))) 16: ALLOCATE(SAVEDXYZ(2*CLIMBERSTEPS*(3*NATOMS)))
 22: ALLOCATE(MOVINGXYZ(3*NATOMS)) 17: ALLOCATE(MOVINGXYZ(3*NATOMS))
 23: ALLOCATE(CONSTRATOMS(3*NATOMS)) 18: ALLOCATE(CONSTRATOMS(3*NATOMS))
 24: WRITE(*,'(A)') ' climber> Set up climber run.' 
 25: CALL CLIMBERINPUT() !read in all possible contraints 19: CALL CLIMBERINPUT() !read in all possible contraints
 26: !now get all constraints that are of the right distances after minpermdist 20: !now get all constraints that are of the right distances after minpermdist
 27: CALL MINPERMDIST(QSTART,QEND,NATOMS,DEBUG,PARAM1,PARAM2,PARAM3,BULKT,TWOD,STARTDIST,DIST2,RIGIDBODY,RMAT) 21: CALL MINPERMDIST(QSTART,QEND,NATOMS,DEBUG,PARAM1,PARAM2,PARAM3,BULKT,TWOD,STARTDIST,DIST2,RIGIDBODY,RMAT)
 28: CALL CLIMBERPAIRCONSTR(QSTART,QEND) 22: CALL CLIMBERPAIRCONSTR(QSTART,QEND)
 29: WRITE(*,'(A,F15.5)') ' climber> Initial distance between start and target structure: ', STARTDIST 
 30: NSTRUCTURES=1 23: NSTRUCTURES=1
 31: SAVEDXYZ((NSTRUCTURES-1)*3*NATOMS+1:NSTRUCTURES*3*NATOMS)=QSTART(1:3*NATOMS) 24: SAVEDXYZ((NSTRUCTURES-1)*3*NATOMS+1:NSTRUCTURES*3*NATOMS)=QSTART(1:3*NATOMS)
 32: MOVINGXYZ(1:3*NATOMS)=QSTART(1:3*NATOMS) 
 33: NCYCLES = 0 
 34: DO 25: DO
 35: 56 IF (CLIMBERSPRING.GT.1.0D10) RETURN 26: 56 NEWXYZ=MOVINGXYZ
 36:    NCYCLES = NCYCLES + 1 
 37:    IF (NCYCLES.GT.500) RETURN 
 38:    NEWXYZ(1:3*NATOMS)=MOVINGXYZ(1:3*NATOMS) 
 39:   ! CALL POTENTIAL(NEWXYZ,EREAL,MOVINGXYZ,.TRUE.,.FALSE.,RMS,PTEST,.FALSE.) 
 40:   ! WRITE(*,'(A,F15.5)') ' climber> potential: ',EREAL 
 41:    CALL MYLBFGS(NOPT,MUPDATE,NEWXYZ,.FALSE., & 27:    CALL MYLBFGS(NOPT,MUPDATE,NEWXYZ,.FALSE., &
 42:   &         MFLAG,EDUMMY,RMS2,EREAL,RMS,BFGSSTEPS,.TRUE.,ITDONE,PTEST,VNEW,.TRUE.,.FALSE.) 28:   &         MFLAG,EDUMMY,RMS2,EREAL,RMS,BFGSSTEPS,.TRUE.,ITDONE,PTEST,VNEW,.TRUE.,.FALSE.)
 43:    IF (MFLAG) WRITE(*,'(A,G20.10,A,I8)') ' climber> converged minimum with energy: ',EREAL,' iterations: ',ITDONE 29:    IF (MFLAG) WRITE(*,'(A,I6,A,G20.10)') 'climber> converged minimum with energy: ',EDUMMY,' iterations: ',ITDONE
 44:    IF (.NOT.MFLAG) THEN 30:    IF (.NOT.MFLAG) WRITE(*,'(A)') 'climber> optimisation did not converge'
 45:       WRITE(*,'(A)') ' climber> optimisation did not converge, increase Kspring' 31: 
 46:       CLIMBERSPRING=CLIMBERSPRING*INC_SPRING 
 47:       WRITE(*,'(A,F15.5)') ' climber> Kspring = ',CLIMBERSPRING 
 48:       GOTO 56 
 49:    ENDIF 
 50:    CALL MINPERMDIST(NEWXYZ,QEND,NATOMS,DEBUG,PARAM1,PARAM2,PARAM3,BULKT,TWOD,DISTANCE,DIST2,RIGIDBODY,RMAT) 32:    CALL MINPERMDIST(NEWXYZ,QEND,NATOMS,DEBUG,PARAM1,PARAM2,PARAM3,BULKT,TWOD,DISTANCE,DIST2,RIGIDBODY,RMAT)
 51:    WRITE(*,'(A,F15.5)') ' climber> Distance of new minimum and target: ',DISTANCE 
 52:    IF (DISTANCE.LE.CLIMBERCONV) THEN 33:    IF (DISTANCE.LE.CLIMBERCONV) THEN
 53:       NSTRUCTURES=NSTRUCTURES+1 34:       NSTRUCTURES=NSTRUCTURES+1
 54:       SAVEDXYZ((NSTRUCTURES-1)*3*NATOMS+1:NSTRUCTURES*3*NATOMS)=NEWXYZ(1:3*NATOMS) 35:       SAVEDXYZ((NSTRUCTURES-1)*3*NATOMS+1:NSTRUCTURES*3*NATOMS)=NEWXYZ(1:3*NATOMS)
 55:       GOTO 57 36:       GOTO 57
 56:    ENDIF 37:    ENDIF
 57:    IF (ABS(DISTANCE-STARTDIST).GE.1.5*STARTDIST/CLIMBERSTEPS) THEN  38:    IF (ABS(DISTANCE-STARTDIST).GE.1.5*STARTDIST/CLIMBERSTEPS) THEN 
 58:       WRITE(*,'(A,F15.5,A)') ' climber> structure accepted, Kspring=',CLIMBERSPRING,', decrease Kspring for next move' 39:       WRITE(*,'(A,F15.5,A)') 'climber> structure accepted, Kspring=',CLIMBERSPRING,', decrease Kspring for next move'
 59:       !save structure 40:       !save structure
 60:       NSTRUCTURES=NSTRUCTURES+1 41:       NSTRUCTURES=NSTRUCTURES+1
 61:       SAVEDXYZ((NSTRUCTURES-1)*3*NATOMS+1:NSTRUCTURES*3*NATOMS)=NEWXYZ(1:3*NATOMS) 42:       SAVEDXYZ((NSTRUCTURES-1)*3*NATOMS+1:NSTRUCTURES*3*NATOMS)=NEWXYZ(1:3*NATOMS)
 62:       MOVINGXYZ(1:3*NATOMS)=NEWXYZ(1:3*NATOMS) 43:       MOVINGXYZ=NEWXYZ
 63:       CLIMBERSPRING=CLIMBERSPRING*RED_SPRING 44:       CLIMBERSPRING=CLIMBERSPRING*RED_SPRING
 64:       WRITE(*,'(A,F15.5)') ' climber> Kspring = ',CLIMBERSPRING 
 65:       GOTO 56 45:       GOTO 56
 66:    ELSE 46:    ELSE
 67:       WRITE(*,'(A)') ' climber> structure not accepted, increase Kspring' 47:       WRITE(*,'(A)') 'climber> structure not accepted, increase Kspring'
 68:       CLIMBERSPRING=CLIMBERSPRING*INC_SPRING 48:       CLIMBERSPRING=CLIMBERSPRING*INC_SPRING
 69:       WRITE(*,'(A,F15.5)') ' climber> Kspring = ',CLIMBERSPRING 
 70:       GOTO 56 49:       GOTO 56
 71:    ENDIF 50:    ENDIF
  51: 57 WRITE(*,'(A)') 'climber> structures converged'
  52:    NSTRUCTURES=NSTRUCTURES+1
  53:    SAVEDXYZ((NSTRUCTURES-1)*3*NATOMS+1:NSTRUCTURES*3*NATOMS)=QEND(1:3*NATOMS)
  54:    CALL CLIMBEROUT(NSTRUCTURES)
  55:    RETURN
 72: ENDDO 56: ENDDO
 73: 57 WRITE(*,'(A)') ' climber> structures converged' 
 74: NSTRUCTURES=NSTRUCTURES+1 
 75: SAVEDXYZ((NSTRUCTURES-1)*3*NATOMS+1:NSTRUCTURES*3*NATOMS)=QEND(1:3*NATOMS) 
 76: CALL CLIMBEROUT(NSTRUCTURES) 
 77: RETURN 
 78: END SUBROUTINE CLIMBERINT 57: END SUBROUTINE CLIMBERINT
 79:  58: 
 80: SUBROUTINE CLIMBEROUT(NSTRUCTURES) 59: SUBROUTINE CLIMBEROUT(NSTRUCTURES)
 81: USE COMMONS, ONLY: NATOMS 60: USE COMMONS, ONLY: NATOMS
 82: USE PORFUNCS 61: USE PORFUNCS
 83: USE CLIMBERCOMMONS 62: USE CLIMBERCOMMONS
 84: IMPLICIT NONE 63: IMPLICIT NONE
 85: INTEGER :: J1,J2, GETUNIT, CLUNIT, NSTRUCTURES 64: INTEGER :: J1,J2, GETUNIT, CLUNIT, NSTRUCTURES
 86:  65: 
 87: CLUNIT=GETUNIT() 66: CLUNIT=GETUNIT()
102: USE PORFUNCS 81: USE PORFUNCS
103: IMPLICIT NONE 82: IMPLICIT NONE
104: INTEGER INPUNIT,GETUNIT,J1,NDUMMY,J2,J3 83: INTEGER INPUNIT,GETUNIT,J1,NDUMMY,J2,J3
105: INTEGER , PARAMETER :: NWORDS=20 84: INTEGER , PARAMETER :: NWORDS=20
106: CHARACTER(100) :: ENTRY 85: CHARACTER(100) :: ENTRY
107: CHARACTER(25) :: ENTRIES(NWORDS)='' 86: CHARACTER(25) :: ENTRIES(NWORDS)=''
108:  87: 
109: INPUNIT=GETUNIT() 88: INPUNIT=GETUNIT()
110: OPEN(INPUNIT,FILE='climber.constraints',STATUS='OLD',ACTION='READ') 89: OPEN(INPUNIT,FILE='climber.constraints',STATUS='OLD',ACTION='READ')
111: READ(INPUNIT,'(I8)') NBBCONSTR 90: READ(INPUNIT,'(I8)') NBBCONSTR
112: WRITE(*,'(A,I8)') ' climber> Number of C alpha for constraints: ',NBBCONSTR 91: WRITE(*,'(A,I8)') 'climber> Number of C alpha for constraints: ',NBBCONSTR
113: ALLOCATE(BBCONSTR(NBBCONSTR)) 92: ALLOCATE(BBCONSTR(NBBCONSTR))
114: J1=NBBCONSTR/10 93: J1=NBBCONSTR/10
115: J1=J1+(NBBCONSTR-J1*10) 94: J1=J1+(NBBCONSTR-J1*10)
116: NDUMMY=1 95: NDUMMY=1
117: DO J2=1,J1 96: DO J2=1,J1
118:    READ(INPUNIT,'(A)') ENTRY 97:    READ(INPUNIT,'(A)') ENTRY
119:    CALL READ_LINE(ENTRY,NWORDS,ENTRIES) 98:    CALL READ_LINE(ENTRY,NWORDS,ENTRIES)
120:    J3=1 99:    J3=1
121:    DO WHILE(J3.LE.10)100:    DO WHILE(J3.LE.10)
122:       IF (NDUMMY.LE.NBBCONSTR) THEN101:       IF (NDUMMY.LE.NBBCONSTR) THEN
123:          READ(ENTRIES(J3),'(I8)') BBCONSTR(NDUMMY)102:          READ(ENTRIES(J3),'(I8)') BBCONSTR(NDUMMY)
124:          NDUMMY=NDUMMY+1103:          NDUMMY=NDUMMY+1
125:       END IF104:       END IF
126:       J3=J3+1105:       J3=J3+1
127:    END DO106:    END DO
128: END DO107: END DO
129: READ(INPUNIT,'(I8)') NSCCONSTR108: READ(INPUNIT,'(I8)') NSCCONSTR
130: WRITE(*,'(A,I8)') ' climber> Number of heavy atoms in sidechains for constraints: ', NSCCONSTR109: WRITE(*,'(A,I8)') 'climber> Number of heavy atoms in sidechains for constraints: ', NSCCONSTR
131: ALLOCATE(SCCONSTR(NSCCONSTR))110: ALLOCATE(SCCONSTR(NSCCONSTR))
132: J1=NSCCONSTR/10111: J1=NSCCONSTR/10
133: J1=J1+(NSCCONSTR-J1*10)112: J1=J1+(NSCCONSTR-J1*10)
134: NDUMMY=1113: NDUMMY=1
135: DO J2=1,J1114: DO J2=1,J1
136:    READ(INPUNIT,'(A)',END=99) ENTRY115:    READ(INPUNIT,'(A)') ENTRY
137:    CALL READ_LINE(ENTRY,NWORDS,ENTRIES)116:    CALL READ_LINE(ENTRY,NWORDS,ENTRIES)
138:    J3=1117:    J3=1
139:    DO WHILE(J3.LE.10)118:    DO WHILE(J3.LE.10)
140:       IF (NDUMMY.LE.NSCCONSTR) THEN119:       IF (NDUMMY.LE.NSCCONSTR) THEN
141:          READ(ENTRIES(J3),'(I8)') SCCONSTR(NDUMMY)120:          READ(ENTRIES(J3),'(I8)') SCCONSTR(NDUMMY)
142:          NDUMMY=NDUMMY+1121:          NDUMMY=NDUMMY+1
143:       END IF122:       END IF
144:       J3=J3+1123:       J3=J3+1
145:    END DO124:    END DO
146: END DO125: END DO
147: 99 CLOSE(INPUNIT)126: CLOSE(INPUNIT)
148: RETURN127: RETURN
149: END SUBROUTINE CLIMBERINPUT128: END SUBROUTINE CLIMBERINPUT
150: 129: 
151: SUBROUTINE CLIMBERPAIRCONSTR(QSTART,QEND)130: SUBROUTINE CLIMBERPAIRCONSTR(QSTART,QEND)
152: USE COMMONS, ONLY: NATOMS,DEBUG131: USE COMMONS, ONLY: NATOMS,DEBUG
153: USE CLIMBERCOMMONS132: USE CLIMBERCOMMONS
154: IMPLICIT NONE133: IMPLICIT NONE
155: INTEGER :: J1,J2,J3,J4134: INTEGER :: J1,J2,J3,J4
156: DOUBLE PRECISION :: DS,DF,QSTART(3*NATOMS),QEND(3*NATOMS)135: DOUBLE PRECISION :: DS,DF,QSTART(3*NATOMS),QEND(3*NATOMS)
157: 136: 
158: CONSTRATOMS(1:3*NATOMS)=.FALSE.137: CONSTRATOMS(1:3*NATOMS)=.FALSE.
159: ALLOCATE(PAIRCONSTR(100000,2))138: ALLOCATE(PAIRCONSTR(10000,2))
160: ALLOCATE(TARGETDIST(100000))139: ALLOCATE(TARGETDIST(10000))
161: NCONSTR=0140: NCONSTR=0
162: DO J1=1,(NBBCONSTR-1)141: DO J1=1,(NBBCONSTR-1)
163:    J3=BBCONSTR(J1)142:    J3=BBCONSTR(J1)
164:    DO J2=(J1+1),NBBCONSTR143:    DO J2=(J1+1),NBBCONSTR
165:       J4=BBCONSTR(J2)144:       J4=BBCONSTR(J2)
166:       DS=SQRT((QSTART(3*(J3-1)+1)-QSTART(3*(J4-1)+1))**2 &145:       DS=SQRT((QSTART(3*(J3-1)+1)-QSTART(3*(J4-1)+1))**2 &
167:   &          +(QSTART(3*(J3-1)+2)-QSTART(3*(J4-1)+2))**2 &146:   &          +(QSTART(3*(J3-1)+2)-QSTART(3*(J4-1)+2))**2 &
168:   &          +(QSTART(3*(J3-1)+3)-QSTART(3*(J4-1)+3))**2) 147:   &          +(QSTART(3*(J3-1)+3)-QSTART(3*(J4-1)+3))**2) 
169:       DF=SQRT((QEND(3*(J3-1)+1)-QEND(3*(J4-1)+1))**2 &148:       DF=SQRT((QEND(3*(J3-1)+1)-QEND(3*(J4-1)+1))**2 &
170:   &          +(QEND(3*(J3-1)+2)-QEND(3*(J4-1)+2))**2 &149:   &          +(QEND(3*(J3-1)+2)-QEND(3*(J4-1)+2))**2 &
171:   &          +(QEND(3*(J3-1)+3)-QEND(3*(J4-1)+3))**2) 150:   &          +(QEND(3*(J3-1)+3)-QEND(3*(J4-1)+3))**2) 
172:       !IF (DEBUG) WRITE(*,'(A,2I8,A,2F15.5)') ' climber> Atoms: ',J3,J4,' distances in endpoints: ',DS,DF151:       IF (DEBUG) WRITE(*,'(A,2I8,A,2F15.5)') 'climber> Atoms: ',J3,J4,' distances in endpoints: ',DS,DF
173:       IF ((DS.GT.10.0).AND.(DF.GT.10.0)) THEN152:       IF ((DS.GT.10.0).AND.(DF.GT.10.0)) THEN
174:          NCONSTR=NCONSTR+1153:          NCONSTR=NCONSTR+1
175:          PAIRCONSTR(NCONSTR,1)=J3154:          PAIRCONSTR(NCONSTR,1)=J3
176:          PAIRCONSTR(NCONSTR,2)=J4155:          PAIRCONSTR(NCONSTR,2)=J4
177:          TARGETDIST(NCONSTR)=DF156:          TARGETDIST(NCONSTR)=DF
178:          WRITE(*,'(A,2I8,A,I8)')  ' climber> New constraint for atoms: ',J3,J4,' Number of constraints: ',NCONSTR157:          WRITE(*,'(A,2I8,A,I8)')  'climber> New constraint for atoms: ',J3,J4,' Number of constraints: ',NCONSTR
179:          CONSTRATOMS(J3)=.TRUE.158:          CONSTRATOMS(J3)=.TRUE.
180:          CONSTRATOMS(J4)=.TRUE.159:          CONSTRATOMS(J4)=.TRUE.
181:       END IF160:       END IF
182:    END DO161:    END DO
183: END DO162: END DO
184: DO J1=1,(NSCCONSTR-1)163: DO J1=1,(NSCCONSTR-1)
185:    J3=SCCONSTR(J1)164:    J3=SCCONSTR(J1)
186:    DO J2=(J1+1),NSCCONSTR165:    DO J2=(J1+1),NSCCONSTR
187:       J4=SCCONSTR(J2)166:       J4=SCCONSTR(J2)
188:       DS=SQRT((QSTART(3*(J3-1)+1)-QSTART(3*(J4-1)+1))**2 &167:       DS=SQRT((QSTART(3*(J3-1)+1)-QSTART(3*(J4-1)+1))**2 &
189:   &          +(QSTART(3*(J3-1)+2)-QSTART(3*(J4-1)+2))**2 &168:   &          +(QSTART(3*(J3-1)+2)-QSTART(3*(J4-1)+2))**2 &
190:   &          +(QSTART(3*(J3-1)+3)-QSTART(3*(J4-1)+3))**2)169:   &          +(QSTART(3*(J3-1)+3)-QSTART(3*(J4-1)+3))**2)
191:       DF=SQRT((QEND(3*(J3-1)+1)-QEND(3*(J4-1)+1))**2 &170:       DF=SQRT((QEND(3*(J3-1)+1)-QEND(3*(J4-1)+1))**2 &
192:   &          +(QEND(3*(J3-1)+2)-QEND(3*(J4-1)+2))**2 &171:   &          +(QEND(3*(J3-1)+2)-QEND(3*(J4-1)+2))**2 &
193:   &          +(QEND(3*(J3-1)+3)-QEND(3*(J4-1)+3))**2)    172:   &          +(QEND(3*(J3-1)+3)-QEND(3*(J4-1)+3))**2)    
194:       !IF (DEBUG) WRITE(*,'(A,2I8,A,2F15.5)') ' climber> Atoms: ',J3,J4,' distances in endpoints: ',DS,DF173:       IF (DEBUG) WRITE(*,'(A,2I8,A,2F15.5)') 'climber> Atoms: ',J3,J4,' distances in endpoints: ',DS,DF
195:       IF ((DS.LT.10.0).AND.(DF.LT.10.0)) THEN174:       IF ((DS.LT.10.0).AND.(DF.LT.10.0)) THEN
196:          NCONSTR=NCONSTR+1175:          NCONSTR=NCONSTR+1
197:          PAIRCONSTR(NCONSTR,1)=J3176:          PAIRCONSTR(NCONSTR,1)=J3
198:          PAIRCONSTR(NCONSTR,2)=J4177:          PAIRCONSTR(NCONSTR,2)=J4
199:          TARGETDIST(NCONSTR)=DF178:          TARGETDIST(NCONSTR)=DF
200:          WRITE(*,'(A,2I8,A,I8)')  ' climber> New constraint for atoms: ',J3,J4,' Number of constraints: ',NCONSTR179:          WRITE(*,'(A,2I8,A,I8)')  'climber> New constraint for atoms: ',J3,J4,' Number of constraints: ',NCONSTR
201:          CONSTRATOMS(J3)=.TRUE.180:          CONSTRATOMS(J3)=.TRUE.
202:          CONSTRATOMS(J4)=.TRUE.181:          CONSTRATOMS(J4)=.TRUE.
203:       END IF182:       END IF
204:    END DO183:    END DO
205: END DO184: END DO
206: WRITE(*,'(A,I8,A)') ' climber> There are ',NCONSTR,' constrained atom pairs.' 
207: RETURN 
208: END SUBROUTINE CLIMBERPAIRCONSTR185: END SUBROUTINE CLIMBERPAIRCONSTR
209: 186: 
210: SUBROUTINE CLIMBERENERGY(XYZ,ERESTRAINT)187: SUBROUTINE CLIMBERENERGY(XYZ,ERESTRAINT)
211: USE COMMONS, ONLY: NATOMS,DEBUG188: USE COMMONS, ONLY: NATOMS,DEBUG
212: USE CLIMBERCOMMONS189: USE CLIMBERCOMMONS
213: USE KEY,ONLY:CLIMBERSPRING190: USE KEY,ONLY:CLIMBERSPRING
214: IMPLICIT NONE191: IMPLICIT NONE
215: DOUBLE PRECISION, INTENT(IN) :: XYZ(3*NATOMS)192: DOUBLE PRECISION, INTENT(IN) :: XYZ(3*NATOMS)
216: DOUBLE PRECISION :: DIST, ERESTRAINT, DISTSUM193: DOUBLE PRECISION :: DIST, ERESTRAINT, DISTSUM
217: INTEGER :: J1194: INTEGER :: J1
218: 195: 
219: DISTSUM=0196: DISTSUM=0
220: DO J1=1,NCONSTR197: DO J1=1,NCONSTR
221:    DIST=SQRT((XYZ(3*(PAIRCONSTR(J1,1)-1)+1)-XYZ(3*(PAIRCONSTR(J1,2)-1)+1))**2 + &198:    DIST=SQRT((XYZ(3*(PAIRCONSTR(J1,1)-1)+1)-XYZ(3*(PAIRCONSTR(J1,2)-1)+1))**2 + &
222:  &           (XYZ(3*(PAIRCONSTR(J1,1)-1)+2)-XYZ(3*(PAIRCONSTR(J1,2)-1)+2))**2 + &199:  &           (XYZ(3*(PAIRCONSTR(J1,1)-1)+2)-XYZ(3*(PAIRCONSTR(J1,2)-1)+2))**2 + &
223:  &           (XYZ(3*(PAIRCONSTR(J1,1)-1)+3)-XYZ(3*(PAIRCONSTR(J1,2)-1)+3))**2)200:  &           (XYZ(3*(PAIRCONSTR(J1,1)-1)+3)-XYZ(3*(PAIRCONSTR(J1,2)-1)+3))**2)
224:    IF (DEBUG) THEN201:    IF (DEBUG) THEN
225:       WRITE(*,'(A,2I8,A,F15.5)') ' climberenergy> Atoms: ',(PAIRCONSTR(J1,1)),(PAIRCONSTR(J1,2)),', distance: ', DIST202:       WRITE(*,'(A,2I8,A,F15.5)') 'climberenergy> Atoms: ',(PAIRCONSTR(J1,1)),(PAIRCONSTR(J1,2)),', distance: ', DIST
226:    ENDIF203:    ENDIF
227:    DISTSUM=DISTSUM+(DIST-TARGETDIST(J1))**2204:    DISTSUM=DISTSUM+(DIST-TARGETDIST(J1))**2
228: ENDDO205: ENDDO
229: ERESTRAINT=CLIMBERSPRING*SQRT(DISTSUM)206: ERESTRAINT=CLIMBERSPRING*SQRT(DISTSUM)
230: RETURN207: RETURN
231: END SUBROUTINE CLIMBERENERGY208: END SUBROUTINE CLIMBERENERGY
232: 209: 
233: 210: 
234: SUBROUTINE CLIMBERGRADIENT(XYZ,CLIMBER_GRADATOMS)211: SUBROUTINE CLIMBERGRADIENT(XYZ,CLIMBER_GRADATOMS)
235: USE CLIMBERCOMMONS212: USE CLIMBERCOMMONS
268: 245: 
269: 246: 
270: SUBROUTINE CLIMBER_ENERGY_GRADIENT(ATOMS,COORDS,CLIMBER_ENERGY,CLIMBER_GRADATOMS)247: SUBROUTINE CLIMBER_ENERGY_GRADIENT(ATOMS,COORDS,CLIMBER_ENERGY,CLIMBER_GRADATOMS)
271: 248: 
272: USE COMMONS, ONLY: NATOMS,DEBUG249: USE COMMONS, ONLY: NATOMS,DEBUG
273: IMPLICIT NONE250: IMPLICIT NONE
274: DOUBLE PRECISION :: COORDS(3*NATOMS), CLIMBER_ENERGY, CLIMBER_GRADATOMS(3*NATOMS)251: DOUBLE PRECISION :: COORDS(3*NATOMS), CLIMBER_ENERGY, CLIMBER_GRADATOMS(3*NATOMS)
275: INTEGER :: ATOMS252: INTEGER :: ATOMS
276: 253: 
277: CALL CLIMBERENERGY(COORDS,CLIMBER_ENERGY)254: CALL CLIMBERENERGY(COORDS,CLIMBER_ENERGY)
278: IF (DEBUG) WRITE(*,'(A,G20.10)') ' climber> Restraint energy: ',CLIMBER_ENERGY255: IF (DEBUG) WRITE(*,'(A,G20.10)') 'climber> Restraint energy: ',CLIMBER_ENERGY
279: CALL CLIMBERGRADIENT(COORDS,CLIMBER_GRADATOMS)256: CALL CLIMBERGRADIENT(COORDS,CLIMBER_GRADATOMS)
280: RETURN257: RETURN
281: END SUBROUTINE CLIMBER_ENERGY_GRADIENT258: END SUBROUTINE CLIMBER_ENERGY_GRADIENT
282: 259: 
283: 260: 
284: SUBROUTINE READ_LINE(LINE,NWORDS,WORDSOUT)261: SUBROUTINE READ_LINE(LINE,NWORDS,WORDSOUT)
285:       CHARACTER(*), INTENT(IN) :: LINE262:       CHARACTER(*), INTENT(IN) :: LINE
286:       INTEGER, INTENT(IN) :: NWORDS263:       INTEGER, INTENT(IN) :: NWORDS
287:       CHARACTER(*), DIMENSION(NWORDS), INTENT(OUT) :: WORDSOUT264:       CHARACTER(*), DIMENSION(NWORDS), INTENT(OUT) :: WORDSOUT
288:       INTEGER:: J1,START_IND,END_IND,J2265:       INTEGER:: J1,START_IND,END_IND,J2


r30147/fetchz.f 2016-03-16 18:34:19.879537110 +0000 r30146/fetchz.f 2016-03-16 18:34:21.487553665 +0000
216:      &                  ' fetchz> Maximum number of steps for pathway minimisation=',PATHSDSTEPS216:      &                  ' fetchz> Maximum number of steps for pathway minimisation=',PATHSDSTEPS
217:       ENDIF217:       ENDIF
218: 218: 
219:       IF (GRADSQ) WRITE(*,'(A)') ' fetchz> Using the modulus gradient as the objective function'219:       IF (GRADSQ) WRITE(*,'(A)') ' fetchz> Using the modulus gradient as the objective function'
220:       PRINT*220:       PRINT*
221:       CALL FLUSH(6)221:       CALL FLUSH(6)
222: 222: 
223:       IF (AMBERT.OR.NABT.OR.AMBER12T) THEN223:       IF (AMBERT.OR.NABT.OR.AMBER12T) THEN
224:          NOPT=3*NATOMS224:          NOPT=3*NATOMS
225:          WRITE(*,'(A,I4,A,I4,A)') ' fetchz> ',NOPT,' Cartesian coordinates will be optimised for ',NATOMS,' AMBER atoms'225:          WRITE(*,'(A,I4,A,I4,A)') ' fetchz> ',NOPT,' Cartesian coordinates will be optimised for ',NATOMS,' AMBER atoms'
226:          IF (TWOENDS.OR.CONNECTT.OR.NEBT.OR.NEWNEBT.OR.DRAGT.OR.GUESSPATHT.OR.CLIMBERT226:          IF (TWOENDS.OR.CONNECTT.OR.NEBT.OR.NEWNEBT.OR.DRAGT.OR.GUESSPATHT
227:      $     .OR.MECCANOT.OR.MORPHT.OR.GREATCIRCLET.OR.GROWSTRINGT.OR.BHINTERPT.OR.BISECTT.OR.CLOSESTALIGNMENT) THEN227:      $     .OR.MECCANOT.OR.MORPHT.OR.GREATCIRCLET.OR.GROWSTRINGT.OR.BHINTERPT.OR.BISECTT.OR.CLOSESTALIGNMENT) THEN
228:             OPEN (UNIT=7,FILE=FINSTRING,STATUS='OLD')228:             OPEN (UNIT=7,FILE=FINSTRING,STATUS='OLD')
229:           IF(AMBERT.OR.NABT.OR.AMBER12T) THEN      ! read coordinates from file finish (containing only coordinates)229:           IF(AMBERT.OR.NABT.OR.AMBER12T) THEN      ! read coordinates from file finish (containing only coordinates)
230:             DO J1=1,NATOMS230:             DO J1=1,NATOMS
231:                READ(7,*)  FIN(3*(J1-1)+1),FIN(3*(J1-1)+2),FIN(3*(J1-1)+3)231:                READ(7,*)  FIN(3*(J1-1)+1),FIN(3*(J1-1)+2),FIN(3*(J1-1)+3)
232:             ENDDO232:             ENDDO
233:             CLOSE(7)233:             CLOSE(7)
234:           ELSE234:           ELSE
235:             DO J1=1,NATOMS235:             DO J1=1,NATOMS
236:                READ(7,*) CD1, CD2, IDUM1, IDUM2, FIN(3*(J1-1)+1),FIN(3*(J1-1)+2),FIN(3*(J1-1)+3)236:                READ(7,*) CD1, CD2, IDUM1, IDUM2, FIN(3*(J1-1)+1),FIN(3*(J1-1)+2),FIN(3*(J1-1)+3)


r30147/key.f90 2016-03-16 18:34:20.071539089 +0000 r30146/key.f90 2016-03-16 18:34:21.679555640 +0000
295:       !ds656> Mie field(s) for modelling substrate effects295:       !ds656> Mie field(s) for modelling substrate effects
296:       LOGICAL :: MIEFT, MIEF_CUTT, MIEF_PBCT296:       LOGICAL :: MIEFT, MIEF_CUTT, MIEF_PBCT
297:       CHARACTER(LEN=130) :: MIEF_FILENAME297:       CHARACTER(LEN=130) :: MIEF_FILENAME
298:       INTEGER :: MIEF_NSITES,MIEF_N,MIEF_M298:       INTEGER :: MIEF_NSITES,MIEF_N,MIEF_M
299:       DOUBLE PRECISION :: MIEF_BOX(3), MIEF_RCUT299:       DOUBLE PRECISION :: MIEF_BOX(3), MIEF_RCUT
300:       DOUBLE PRECISION, ALLOCATABLE :: MIEF_EPS(:), MIEF_SIG(:), &300:       DOUBLE PRECISION, ALLOCATABLE :: MIEF_EPS(:), MIEF_SIG(:), &
301:            MIEF_SITES(:,:), MIEF_U_RCUT(:), MIEF_DUDR_RCUT(:)301:            MIEF_SITES(:,:), MIEF_U_RCUT(:), MIEF_DUDR_RCUT(:)
302:       302:       
303: ! AMBER 12 variables303: ! AMBER 12 variables
304:       LOGICAL :: AMBER12T304:       LOGICAL :: AMBER12T
305:       LOGICAL :: CHIRALENDPOINTS305:       LOGICAL :: CLIMBERT
306:       LOGICAL :: CLIMBERT,CLIMBERINIT 
307:       INTEGER :: CLIMBERSTEPS306:       INTEGER :: CLIMBERSTEPS
308:       DOUBLE PRECISION :: CLIMBERCONV,CLIMBERSPRING307:       DOUBLE PRECISION :: CLIMBERCONV,CLIMBERSPRING
309:       INTEGER, DIMENSION(:,:), ALLOCATABLE :: BONDS308:       INTEGER, DIMENSION(:,:), ALLOCATABLE :: BONDS
310: 309: 
311:       DOUBLE PRECISION, ALLOCATABLE ::  MLPDAT(:,:)310:       DOUBLE PRECISION, ALLOCATABLE ::  MLPDAT(:,:)
312:       INTEGER, ALLOCATABLE ::  MLPOUTCOME(:)311:       INTEGER, ALLOCATABLE ::  MLPOUTCOME(:)
313: 312: 
314: 313: 
315: END MODULE KEY314: END MODULE KEY


r30147/keywords.f 2016-03-16 18:34:20.267541103 +0000 r30146/keywords.f 2016-03-16 18:34:21.875557637 +0000
881:          MCPATHDOBLOCK=0881:          MCPATHDOBLOCK=0
882:          MCPATHNEGLECT=0.0D0882:          MCPATHNEGLECT=0.0D0
883:          SSHT=.FALSE.883:          SSHT=.FALSE.
884:          NCPU=0884:          NCPU=0
885:          MCPATHGWS=0.3D0885:          MCPATHGWS=0.3D0
886:          MCPATHGWQ=0.00025D0886:          MCPATHGWQ=0.00025D0
887:          COLLAGENOP=.FALSE.887:          COLLAGENOP=.FALSE.
888: 888: 
889:          DUMPFRQST=.FALSE.889:          DUMPFRQST=.FALSE.
890: 890: 
891:          CHIRALENDPOINTS=.TRUE. 
892:  
893:          CUDAT=.FALSE.891:          CUDAT=.FALSE.
894:          CUDAPOT=' '892:          CUDAPOT=' '
895:          CUDATIMET=.FALSE.893:          CUDATIMET=.FALSE.
896: 894: 
897:          CLIMBERT=.FALSE.895:          CLIMBERT=.FALSE.
898:          CLIMBERINIT=.FALSE. 
899:          CLIMBERSTEPS=20896:          CLIMBERSTEPS=20
900:          CLIMBERCONV=0.2897:          CLIMBERCONV=0.2
901:          CLIMBERSPRING=5.0898:          CLIMBERSPRING=5.0
902: 899: 
903: !900: !
904: ! Neural network potential901: ! Neural network potential
905: !902: !
906:          MLP3T=.FALSE.903:          MLP3T=.FALSE.
907:          MLPB3T=.FALSE.904:          MLPB3T=.FALSE.
908:          MLPPROB=.FALSE.905:          MLPPROB=.FALSE.
1656:          ELSE IF (WORD.EQ.'INE_NEW') THEN1653:          ELSE IF (WORD.EQ.'INE_NEW') THEN
1657:             CLSTRINGTST=.TRUE.1654:             CLSTRINGTST=.TRUE.
1658:             CLSTRINGT=.TRUE.1655:             CLSTRINGT=.TRUE.
1659:             EVOLVESTRINGT = .TRUE.1656:             EVOLVESTRINGT = .TRUE.
1660:             CALL READF(STTSRMSCONV)1657:             CALL READF(STTSRMSCONV)
1661:             CALL READI(ST_TSSTEP)1658:             CALL READI(ST_TSSTEP)
1662:             CALL READF(LAN_DIST)1659:             CALL READF(LAN_DIST)
1663:             CALL READI(LANSTEP)1660:             CALL READI(LANSTEP)
1664:             CALL READF(LANCONV)1661:             CALL READF(LANCONV)
1665:             CALL READF(LANFACTOR)1662:             CALL READF(LANFACTOR)
 1663: ! kr366 climber input
 1664:          ELSE IF (WORD.EQ.'CLIMBER') THEN
 1665:             CLIMBERT=.TRUE.
 1666:             CALL READI(CLIMBERSTEPS)
 1667:             CALL READF(CLIMBERCONV)
 1668:             CALL READF(CLIMBERSPRING)
1666: 1669: 
1667:          ELSE IF (WORD.EQ.'CLSTRING') THEN1670:          ELSE IF (WORD.EQ.'CLSTRING') THEN
1668:             CLSTRINGT=.TRUE.1671:             CLSTRINGT=.TRUE.
1669:             EVOLVESTRINGT = .TRUE.1672:             EVOLVESTRINGT = .TRUE.
1670:          ELSE IF (WORD.EQ.'COLLAGENOP') THEN1673:          ELSE IF (WORD.EQ.'COLLAGENOP') THEN
1671:             COLLAGENOP=.TRUE.1674:             COLLAGENOP=.TRUE.
1672:             CALL READI(COLLINDICES(1))1675:             CALL READI(COLLINDICES(1))
1673:             CALL READI(COLLINDICES(2))1676:             CALL READI(COLLINDICES(2))
1674:             CALL READI(COLLINDICES(3))1677:             CALL READI(COLLINDICES(3))
1675:             CALL READI(COLLINDICES(4))1678:             CALL READI(COLLINDICES(4))
2042:             CHRIGIDT=.TRUE.2045:             CHRIGIDT=.TRUE.
2043:             CALL READF(PTRANS)2046:             CALL READF(PTRANS)
2044:             CALL READF(TRANSMAX)2047:             CALL READF(TRANSMAX)
2045:             CALL READF(PROT)2048:             CALL READF(PROT)
2046:             CALL READF(ROTMAX)2049:             CALL READF(ROTMAX)
2047: ! 2050: ! 
2048: ! CISTRANS is a CHARMM/AMBER9/NAB related keyword, which allows cis-trans isomerisation of the peptide bond .2051: ! CISTRANS is a CHARMM/AMBER9/NAB related keyword, which allows cis-trans isomerisation of the peptide bond .
2049: ! 2052: ! 
2050:          ELSE IF (WORD.EQ.'CISTRANS') THEN2053:          ELSE IF (WORD.EQ.'CISTRANS') THEN
2051:             CISTRANS=.TRUE.2054:             CISTRANS=.TRUE.
2052:  
2053: ! kr366 climber input 
2054:          ELSE IF (WORD.EQ.'CLIMBER') THEN 
2055:             CLIMBERT=.TRUE. 
2056:             CLIMBERINIT=.TRUE. 
2057:             CALL READI(CLIMBERSTEPS) 
2058:             CALL READF(CLIMBERCONV) 
2059:             CALL READF(CLIMBERSPRING) 
2060:             WRITE(*,'(A,I8)') ' keywords> Climber interpolation steps: ',CLIMBERSTEPS 
2061:             WRITE(*,'(2(A,F10.5))') ' keywords> Convergence: ',CLIMBERCONV, ' initial Kspring: ', CLIMBERSPRING 
2062:  
2063:  
2064: ! 2055: ! 
2065: ! Sometimes have to modify the cold fusion limit when using high electric fields2056: ! Sometimes have to modify the cold fusion limit when using high electric fields
2066: ! 2057: ! 
2067:          ELSE IF (WORD.EQ.'COLDFUSION') THEN2058:          ELSE IF (WORD.EQ.'COLDFUSION') THEN
2068:             IF (NITEMS.GT.1) call READF(COLDFUSIONLIMIT)2059:             IF (NITEMS.GT.1) call READF(COLDFUSIONLIMIT)
2069: ! 2060: ! 
2070: ! Connect initial minimum in odata to final minimum in file finish - maximum2061: ! Connect initial minimum in odata to final minimum in file finish - maximum
2071: ! number of transiiton states=NCONNECT. Obsolete - use NEWCONNECT instead.2062: ! number of transiiton states=NCONNECT. Obsolete - use NEWCONNECT instead.
2072: ! 2063: ! 
2073:          ELSE IF (WORD.EQ.'CONNECT') THEN2064:          ELSE IF (WORD.EQ.'CONNECT') THEN
4140: ! CHARMM related keyword to reject transition states4131: ! CHARMM related keyword to reject transition states
4141: ! that connect two minima with different omega angles, i.e. to prevent cis-trans peptide4132: ! that connect two minima with different omega angles, i.e. to prevent cis-trans peptide
4142: ! isomerisation.4133: ! isomerisation.
4143: ! 4134: ! 
4144:          ELSE IF (WORD.EQ.'UACHIRAL') THEN4135:          ELSE IF (WORD.EQ.'UACHIRAL') THEN
4145:             UACHIRAL=.TRUE.4136:             UACHIRAL=.TRUE.
4146: 4137: 
4147:          ELSE IF (WORD.EQ.'NOCHIRALCHECKS') THEN4138:          ELSE IF (WORD.EQ.'NOCHIRALCHECKS') THEN
4148:             TURNOFFCHECKCHIRALITY=.TRUE.4139:             TURNOFFCHECKCHIRALITY=.TRUE.
4149: 4140: 
4150:          ELSE IF (WORD.EQ.'NOCHIRALENDPOINTS') THEN 
4151:             CHIRALENDPOINTS=.FALSE. 
4152:             WRITE(*,*) ' keywords> No chirality checks for endpoints' 
4153:  
4154:          ELSE IF (WORD.EQ.'NOCISTRANS') THEN4141:          ELSE IF (WORD.EQ.'NOCISTRANS') THEN
4155:             NOCISTRANS=.TRUE.   ! is used in connect.f4142:             NOCISTRANS=.TRUE.   ! is used in connect.f
4156:             CHECKOMEGAT=.TRUE.  ! is used in NEWNEB4143:             CHECKOMEGAT=.TRUE.  ! is used in NEWNEB
4157:             IF (NITEMS.GT.1) CALL READF(MINOMEGA)4144:             IF (NITEMS.GT.1) CALL READF(MINOMEGA)
4158: 4145: 
4159:             IF (NITEMS.GT.2) THEN4146:             IF (NITEMS.GT.2) THEN
4160:                CALL READU(WW)4147:                CALL READU(WW)
4161:                IF (TRIM(ADJUSTL(WW))=='RNA') THEN4148:                IF (TRIM(ADJUSTL(WW))=='RNA') THEN
4162:                   NOCISTRANSRNA = .TRUE.4149:                   NOCISTRANSRNA = .TRUE.
4163:                   write(*,*) ' keywords> NOCISTRANSRNA set to .TRUE.'4150:                   write(*,*) ' keywords> NOCISTRANSRNA set to .TRUE.'


r30147/mylbfgs.f 2016-03-16 18:34:20.463543116 +0000 r30146/mylbfgs.f 2016-03-16 18:34:22.071559650 +0000
117: ! hk286 - local rigid body                                                      117: ! hk286 - local rigid body                                                      
118:       DOUBLE PRECISION :: XCOORDS(3*NATOMS), XRIGIDCOORDS(DEGFREEDOMS)118:       DOUBLE PRECISION :: XCOORDS(3*NATOMS), XRIGIDCOORDS(DEGFREEDOMS)
119: ! hk286 - Thomson problem                                                      119: ! hk286 - Thomson problem                                                      
120:       DOUBLE PRECISION :: ROT(3,3), ROTINV(3,3)120:       DOUBLE PRECISION :: ROT(3,3), ROTINV(3,3)
121:       INTEGER NFAILTHOMSON121:       INTEGER NFAILTHOMSON
122:       LOGICAL GIMBALT122:       LOGICAL GIMBALT
123: ! cs778 reset lbfgs when charmm updating nonbond list123: ! cs778 reset lbfgs when charmm updating nonbond list
124:       LOGICAL RESETLBFGS124:       LOGICAL RESETLBFGS
125:       COMMON /CRESETLBFGS/ RESETLBFGS125:       COMMON /CRESETLBFGS/ RESETLBFGS
126: 126: 
127:  
128:       DOINTERNALSTRANSFORM = intwrap_useinternals()127:       DOINTERNALSTRANSFORM = intwrap_useinternals()
129:       if(DOINTERNALSTRANSFORM) print *, "Doing optimization in internals"128:       if(DOINTERNALSTRANSFORM) print *, "Doing optimization in internals"
130: 129: 
131:       IF (.NOT.ALLOCATED(DIAG)) ALLOCATE(DIAG(N))       ! SAVE doesn't work otherwise for pgi and ifort130:       IF (.NOT.ALLOCATED(DIAG)) ALLOCATE(DIAG(N))       ! SAVE doesn't work otherwise for pgi and ifort
132:       IF (.NOT.ALLOCATED(W)) ALLOCATE(W(N*(2*M+1)+2*M)) ! SAVE doesn't work otherwise for pgi and ifort131:       IF (.NOT.ALLOCATED(W)) ALLOCATE(W(N*(2*M+1)+2*M)) ! SAVE doesn't work otherwise for pgi and ifort
133:       IF (SIZE(W,1).NE.N*(2*M+1)+2*M) THEN ! mustn't call mylbfgs with changing number of variables!!!132:       IF (SIZE(W,1).NE.N*(2*M+1)+2*M) THEN ! mustn't call mylbfgs with changing number of variables!!!
134:          PRINT '(A,I10,A,I10,A)', ' mylbfgs> ERROR, dimension of W=',SIZE(W,1),' but N*(2*M+1)+2*M=',N*(2*M+1)+2*M133:          PRINT '(A,I10,A,I10,A)', ' mylbfgs> ERROR, dimension of W=',SIZE(W,1),' but N*(2*M+1)+2*M=',N*(2*M+1)+2*M
135:          STOP134:          STOP
136:       ENDIF135:       ENDIF
137: 136: 
235:             c(3,J2+nres)=X(6*(J2-1)+6)234:             c(3,J2+nres)=X(6*(J2-1)+6)
236:          END DO235:          END DO
237:          CALL UPDATEDC236:          CALL UPDATEDC
238:          CALL int_from_cart(.true.,.false.)237:          CALL int_from_cart(.true.,.false.)
239: C238: C
240: C jmc 9/1/03 Need this call to chainbuild to calculate xrot,xloc matrices which are used239: C jmc 9/1/03 Need this call to chainbuild to calculate xrot,xloc matrices which are used
241: C in the calculation of d(bond vectors) by d(internals).240: C in the calculation of d(bond vectors) by d(internals).
242: C241: C
243:          CALL chainbuild242:          CALL chainbuild
244:       ENDIF243:       ENDIF
 244: 
245:       IF ((.NOT.KNOWE).OR.(.NOT.KNOWG)) THEN245:       IF ((.NOT.KNOWE).OR.(.NOT.KNOWG)) THEN
246:          IF(AMBERT.OR.NABT) irespa2=ITDONE+1246:          IF(AMBERT.OR.NABT) irespa2=ITDONE+1
247:          RMS=0.0D0247:          RMS=0.0D0
248:          CALL POTENTIAL(X,ENERGY,GSAVE,.TRUE.,GRADSQ,RMS,.FALSE.,.FALSE.)248:          CALL POTENTIAL(X,ENERGY,GSAVE,.TRUE.,GRADSQ,RMS,.FALSE.,.FALSE.)
249:          REALRMS=RMS249:          REALRMS=RMS
250:          IF (ENERGY.LT.COLDFUSIONLIMIT) then250:          IF (ENERGY.LT.COLDFUSIONLIMIT) then
251:             WRITE(*,'(A,2G20.10)') ' mylbfgs> Cold fusion diagnosed - step discarded, energy, limit=',ENERGY,coldFusionLimit251:             WRITE(*,'(A,2G20.10)') ' mylbfgs> Cold fusion diagnosed - step discarded, energy, limit=',ENERGY,coldFusionLimit
252:             ENERGY=1.0d60252:             ENERGY=1.0d60
253:             EREAL=1.0d60253:             EREAL=1.0d60
254:             RMS=1.0d1254:             RMS=1.0d1
464:             ENDIF464:             ENDIF
465:          ENDIF465:          ENDIF
466:          DO J1=1,3*NATOMS466:          DO J1=1,3*NATOMS
467:             G(J1)=VEC2(J1)467:             G(J1)=VEC2(J1)
468: C           G(J1)=VEC2(J1)/ENERGY468: C           G(J1)=VEC2(J1)/ENERGY
469:          ENDDO469:          ENDDO
470:       ELSE470:       ELSE
471:          EREAL=ENERGY471:          EREAL=ENERGY
472:          REALRMS=RMS472:          REALRMS=RMS
473:       ENDIF473:       ENDIF
 474: 
474:       IF (PTEST) WRITE(*,'(A,2G20.10,A,I6,A)') 475:       IF (PTEST) WRITE(*,'(A,2G20.10,A,I6,A)') 
475:      1             ' mylbfgs> Energy and RMS force=',ENERGY,RMS,' after ',ITDONE,' steps'476:      1             ' mylbfgs> Energy and RMS force=',ENERGY,RMS,' after ',ITDONE,' steps'
476:       WRITE(ESTRING,16) ' mylbfgs> Energy for last cycle=',ENERGY,' '477:       WRITE(ESTRING,16) ' mylbfgs> Energy for last cycle=',ENERGY,' '
477: 16       FORMAT(A,27X,F20.10,A)478: 16       FORMAT(A,27X,F20.10,A)
478: 479: 
479: 10    CALL FLUSH(6)480: 10    CALL FLUSH(6)
480: 481: 
481: C     add additional check for gradient in internals. Otherwise rigid bonds might never converge.482: C     add additional check for gradient in internals. Otherwise rigid bonds might never converge.
482: C     GINT with bonds projected out is never transferred back to carthesians -> check in internals483: C     GINT with bonds projected out is never transferred back to carthesians -> check in internals
483:       MFLAG=.FALSE.484:       MFLAG=.FALSE.


r30147/OPTIM.F 2016-03-16 18:34:19.487533078 +0000 r30146/OPTIM.F 2016-03-16 18:34:20.859547186 +0000
 28:       USE ZWK 28:       USE ZWK
 29:       USE MODCHARMM 29:       USE MODCHARMM
 30:       USE MODUNRES 30:       USE MODUNRES
 31:       USE NEWNEBMODULE 31:       USE NEWNEBMODULE
 32:       USE NEWCONNECTMODULE 32:       USE NEWCONNECTMODULE
 33:       USE KEYNEB, NNNIMAGE=>NIMAGE 33:       USE KEYNEB, NNNIMAGE=>NIMAGE
 34:       USE MODGUESS 34:       USE MODGUESS
 35:       USE MODMEC 35:       USE MODMEC
 36:       USE PORFUNCS 36:       USE PORFUNCS
 37:       USE CHARUTILS 37:       USE CHARUTILS
 38:       USE CHIRALITY 
 39:       USE INTCOMMONS, ONLY : INTINTERPT, DESMINT, NATINT, INTMINPERMT !msb50 remove last 38:       USE INTCOMMONS, ONLY : INTINTERPT, DESMINT, NATINT, INTMINPERMT !msb50 remove last
 40:       USE INTCUTILS, ONLY : INTSETUP, INTCLEANUP ! msb50 last 39:       USE INTCUTILS, ONLY : INTSETUP, INTCLEANUP ! msb50 last
 41:       USE MODAMBER9, ONLY : DUMPMODEN,CHIARRAY1,GOODSTRUCTURE1 40:       USE MODAMBER9, ONLY : DUMPMODEN,CHIARRAY1,GOODSTRUCTURE1
 42: !     USE BENCHMARKS, ONLY : MINBMT, MINBM 41: !     USE BENCHMARKS, ONLY : MINBMT, MINBM
 43: ! hk286 42: ! hk286
 44:       USE GENRIGID 43:       USE GENRIGID
 45:  44: 
 46:       IMPLICIT NONE 45:       IMPLICIT NONE
 47: ! subroutine parameters   46: ! subroutine parameters  
 48:       INTEGER F1,F2 47:       INTEGER F1,F2
 55:      2  DIST, OVEC(3), H1VEC(3), H2VEC(3), Q(NOPT), EINITIAL, EFINAL,  54:      2  DIST, OVEC(3), H1VEC(3), H2VEC(3), Q(NOPT), EINITIAL, EFINAL, 
 56:      3  ETIME, FTIME, STIME, DPRAND, DCOORDS(NOPT), INTFREEZETOLSAVE, 55:      3  ETIME, FTIME, STIME, DPRAND, DCOORDS(NOPT), INTFREEZETOLSAVE,
 57:      4  ETS, EPLUS, EMINUS, SLENGTH, DISP, GAMMA, NTILDE,  56:      4  ETS, EPLUS, EMINUS, SLENGTH, DISP, GAMMA, NTILDE, 
 58:      5  FRQSTS(NOPT), FRQSPLUS(NOPT), FRQSMINUS(NOPT), QMINUS(NOPT), DISTSF 57:      5  FRQSTS(NOPT), FRQSPLUS(NOPT), FRQSMINUS(NOPT), QMINUS(NOPT), DISTSF
 59:       DOUBLE PRECISION THTEMP(NOPT) 58:       DOUBLE PRECISION THTEMP(NOPT)
 60:       CHARACTER ESTRING*87, GPSTRING*80, NSTRING*80, FSTRING*80, FNAME*13, FNAMEV*18,  59:       CHARACTER ESTRING*87, GPSTRING*80, NSTRING*80, FSTRING*80, FNAME*13, FNAMEV*18, 
 61:      1          ITSTRING*22, EOFSSTRING*15 60:      1          ITSTRING*22, EOFSSTRING*15
 62:       CHARACTER(LEN=80) FNAMEF 61:       CHARACTER(LEN=80) FNAMEF
 63:       CHARACTER(LEN=20) EFNAME 62:       CHARACTER(LEN=20) EFNAME
 64:       DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: QW 63:       DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: QW
 65:       LOGICAL LDUMMY, BFGSTSSAVE, MFLAG, POTCALL, GOODENDPOINT 64:       LOGICAL LDUMMY, BFGSTSSAVE, MFLAG, POTCALL
 66:       INTEGER FRAME 65:       INTEGER FRAME
 67:       LOGICAL PVFLAG 66:       LOGICAL PVFLAG
 68:       COMMON /PVF/ PVFLAG 67:       COMMON /PVF/ PVFLAG
 69: C     COMMON /VN/ VNEW   !  common SV was also deleted 68: C     COMMON /VN/ VNEW   !  common SV was also deleted
 70:       COMMON /STRINGS/ ESTRING, GPSTRING, NSTRING, FSTRING 69:       COMMON /STRINGS/ ESTRING, GPSTRING, NSTRING, FSTRING
 71:       COMMON /PCALL/ NPCALL, ECALL, FCALL, SCALL, ETIME, FTIME, STIME 70:       COMMON /PCALL/ NPCALL, ECALL, FCALL, SCALL, ETIME, FTIME, STIME
 72:       LOGICAL PATHT, DRAGT 71:       LOGICAL PATHT, DRAGT
 73:       INTEGER NPATHFRAME, NCDONE 72:       INTEGER NPATHFRAME, NCDONE
 74:       COMMON /RUNTYPE/ DRAGT, PATHT, NPATHFRAME 73:       COMMON /RUNTYPE/ DRAGT, PATHT, NPATHFRAME
 75:       LOGICAL CONNECTT, DUMPPATH, READPATH, CALCRATES, STOPFIRST 74:       LOGICAL CONNECTT, DUMPPATH, READPATH, CALCRATES, STOPFIRST
125: C     IF (CONNECTT.AND.NEWNEBT) THEN124: C     IF (CONNECTT.AND.NEWNEBT) THEN
126: C        PRINT*,'WARNING - cannot use old connect with new neb, changing to old neb'125: C        PRINT*,'WARNING - cannot use old connect with new neb, changing to old neb'
127: C        NEWNEBT=.FALSE.126: C        NEWNEBT=.FALSE.
128: C        NEBT=.TRUE.127: C        NEBT=.TRUE.
129: C     ENDIF128: C     ENDIF
130:       IF ((FILTH2.EQ.0).AND.(FILTH.NE.0)) WRITE(FILTHSTR,'(I10)') FILTH ! otherwise FILTHSTR isn;t set correctly.129:       IF ((FILTH2.EQ.0).AND.(FILTH.NE.0)) WRITE(FILTHSTR,'(I10)') FILTH ! otherwise FILTHSTR isn;t set correctly.
131:       IF (REPELTST) ALLOCATE(REPELTS(NOPT,100)) ! PREVIOUS TS GEOMETRIES TO AVOID130:       IF (REPELTST) ALLOCATE(REPELTS(NOPT,100)) ! PREVIOUS TS GEOMETRIES TO AVOID
132:       IF (CHECKINDEX) ALLOCATE(VECCHK(NOPT,MAX(NUSEEV,HINDEX,1))) ! vectors to orthogonise to131:       IF (CHECKINDEX) ALLOCATE(VECCHK(NOPT,MAX(NUSEEV,HINDEX,1))) ! vectors to orthogonise to
133:       ALLOCATE(ZWORK(NOPT,MAX(NUSEEV,HINDEX,1)))                  ! partial eigenvectors storage132:       ALLOCATE(ZWORK(NOPT,MAX(NUSEEV,HINDEX,1)))                  ! partial eigenvectors storage
134:       IF (TWOENDS.OR.CONNECTT.OR.NEWNEBT.OR.DRAGT.OR.GUESSPATHT.OR.MECCANOT.OR.MORPHT.OR.GREATCIRCLET.OR.BHINTERPT.OR.BISECTT 133:       IF (TWOENDS.OR.CONNECTT.OR.NEWNEBT.OR.DRAGT.OR.GUESSPATHT.OR.MECCANOT.OR.MORPHT.OR.GREATCIRCLET.OR.BHINTERPT.OR.BISECTT 
135:      & .OR.CLOSESTALIGNMENT.OR.GROWSTRINGT.OR.CLIMBERT) THEN134:      & .OR.CLOSESTALIGNMENT.OR.GROWSTRINGT) THEN
136:          ALLOCATE(FIN(NOPT),START(NOPT))135:          ALLOCATE(FIN(NOPT),START(NOPT))
137:       ENDIF136:       ENDIF
138: 137: 
139:       NPCALL=0138:       NPCALL=0
140:       ECALL=0139:       ECALL=0
141:       FCALL=0140:       FCALL=0
142:       SCALL=0141:       SCALL=0
143:       ETIME=0142:       ETIME=0
144:       FTIME=0143:       FTIME=0
145:       STIME=0144:       STIME=0
517:                C(3,J1+NRES)=FIN(6*(J1-1)+6)516:                C(3,J1+NRES)=FIN(6*(J1-1)+6)
518:             ENDDO517:             ENDDO
519:             CALL UPDATEDC518:             CALL UPDATEDC
520:             CALL INT_FROM_CART(.TRUE.,.FALSE.)519:             CALL INT_FROM_CART(.TRUE.,.FALSE.)
521:             CALL GEOM_TO_VAR(NINTS,FIN(1:NINTS))520:             CALL GEOM_TO_VAR(NINTS,FIN(1:NINTS))
522: C           CALL UNGUESSPATH(Q,FIN,NINTS,EDIFFTOL,NATOMS)521: C           CALL UNGUESSPATH(Q,FIN,NINTS,EDIFFTOL,NATOMS)
523:             CALL GUESSPATH(Q,FIN,NINTS,EDIFFTOL,NATOMS)522:             CALL GUESSPATH(Q,FIN,NINTS,EDIFFTOL,NATOMS)
524:          ELSE523:          ELSE
525:             CALL GUESSPATH(Q,FIN,NOPT,EDIFFTOL,NATOMS)524:             CALL GUESSPATH(Q,FIN,NOPT,EDIFFTOL,NATOMS)
526:          ENDIF525:          ENDIF
527:       ELSE IF (CONNECTT.OR.BHINTERPT.OR.BISECTT.OR.CLIMBERT) THEN526:       ELSE IF (CONNECTT.OR.BHINTERPT.OR.BISECTT) THEN
528: C        IF (TWOENDS) THEN527: C        IF (TWOENDS) THEN
529: C           CALL CONNECTTWO528: C           CALL CONNECTTWO
530: C        ELSE529: C        ELSE
531:          IF (NEWCONNECTT.OR.BHINTERPT.OR.BISECTT.OR.CLIMBERT) THEN530:          IF (NEWCONNECTT.OR.BHINTERPT.OR.BISECTT) THEN
532:             IF (UNRST) THEN531:             IF (UNRST) THEN
533:                 DO J1=1,NRES532:                 DO J1=1,NRES
534:                    C(1,J1)=Q(6*(J1-1)+1)533:                    C(1,J1)=Q(6*(J1-1)+1)
535:                    C(2,J1)=Q(6*(J1-1)+2)534:                    C(2,J1)=Q(6*(J1-1)+2)
536:                    C(3,J1)=Q(6*(J1-1)+3)535:                    C(3,J1)=Q(6*(J1-1)+3)
537:                    C(1,J1+NRES)=Q(6*(J1-1)+4)536:                    C(1,J1+NRES)=Q(6*(J1-1)+4)
538:                    C(2,J1+NRES)=Q(6*(J1-1)+5)537:                    C(2,J1+NRES)=Q(6*(J1-1)+5)
539:                    C(3,J1+NRES)=Q(6*(J1-1)+6)538:                    C(3,J1+NRES)=Q(6*(J1-1)+6)
540:                 ENDDO539:                 ENDDO
541:                 CALL UPDATEDC540:                 CALL UPDATEDC
638:                      CALL GEOPT(FNAMEF,EFNAME,FIN,VECS,MFLAG,ENERGY,EVALMIN,VNEW)637:                      CALL GEOPT(FNAMEF,EFNAME,FIN,VECS,MFLAG,ENERGY,EVALMIN,VNEW)
639:                   ENDIF638:                   ENDIF
640:                   REOPTIMISEENDPOINTS=.FALSE.639:                   REOPTIMISEENDPOINTS=.FALSE.
641:                   GOTO 555640:                   GOTO 555
642:                ELSE641:                ELSE
643:                   CALL TSUMMARY642:                   CALL TSUMMARY
644:                   CALL FLUSH(6)643:                   CALL FLUSH(6)
645:                   GOTO 765644:                   GOTO 765
646:                ENDIF645:                ENDIF
647:             ENDIF646:             ENDIF
648:             IF (AMBER12T.AND.CHIRALENDPOINTS) THEN 
649:                CALL CHIRALITY_CHECK(FIN,GOODENDPOINT) 
650:                IF (.NOT.GOODENDPOINT) THEN 
651:                   PRINT '(A)', ' OPTIM> Chirality change in endpoints' 
652:                   STOP 
653:                ENDIF 
654:             ENDIF 
655:             IF (CLIMBERT) THEN647:             IF (CLIMBERT) THEN
656:                CLIMBERINIT=.FALSE. 
657:                CALL CLIMBERINT(Q,FIN)648:                CALL CLIMBERINT(Q,FIN)
658:                RETURN649:                RETURN
659:             ENDIF650:             ENDIF
660:             CALL MINPERMDIST(Q,FIN,NATOMS,DEBUG,PARAM1,PARAM2,PARAM3,BULKT,TWOD,DIST,DIST2,RIGIDBODY,RMAT)651:             CALL MINPERMDIST(Q,FIN,NATOMS,DEBUG,PARAM1,PARAM2,PARAM3,BULKT,TWOD,DIST,DIST2,RIGIDBODY,RMAT)
661:             IF (PERMDISTINIT) PERMDIST=.FALSE.652:             IF (PERMDISTINIT) PERMDIST=.FALSE.
662:             IF (ATOMMATCHINIT) ATOMMATCHDIST=.FALSE.653:             IF (ATOMMATCHINIT) ATOMMATCHDIST=.FALSE.
663:             IF (BISECTT) THEN654:             IF (BISECTT) THEN
664:                CALL BISECT_OPT(NATOMS,EINITIAL,Q,EFINAL,FIN,DIST)655:                CALL BISECT_OPT(NATOMS,EINITIAL,Q,EFINAL,FIN,DIST)
665:             ELSE656:             ELSE
666:                IF (ALLOCATED(SAVES)) DEALLOCATE(SAVES)657:                IF (ALLOCATED(SAVES)) DEALLOCATE(SAVES)


r30147/potential.f 2016-03-16 18:34:20.663545171 +0000 r30146/potential.f 2016-03-16 18:34:22.271561707 +0000
3045:             END IF3045:             END IF
3046:             IF (CUDAT) THEN3046:             IF (CUDAT) THEN
3047:                ! This call copies CPU coordinates to GPU, calculates energy/gradient and copies energy/gradient back to CPU3047:                ! This call copies CPU coordinates to GPU, calculates energy/gradient and copies energy/gradient back to CPU
3048:                CALL CUDA_ENEGRAD_WRAPPER(NATOMS, COORDS, ENERGY, GRADATOMS)3048:                CALL CUDA_ENEGRAD_WRAPPER(NATOMS, COORDS, ENERGY, GRADATOMS)
3049:             ELSE3049:             ELSE
3050:                CALL AMBER12_ENERGY_AND_GRADIENT(NATOMS,3050:                CALL AMBER12_ENERGY_AND_GRADIENT(NATOMS,
3051:      &                                          COORDS,3051:      &                                          COORDS,
3052:      &                                          ENERGY,3052:      &                                          ENERGY,
3053:      &                                          GRADATOMS,3053:      &                                          GRADATOMS,
3054:      &                                          ENERGY_DECOMP)3054:      &                                          ENERGY_DECOMP)
3055:                IF (CLIMBERT.AND.(.NOT.(CLIMBERINIT))) THEN3055:                IF (CLIMBERT) THEN
3056:                   CALL CLIMBER_ENERGY_GRADIENT(NATOMS,COORDS,CLIMBER_ENERGY,CLIMBER_GRADATOMS)3056:                   CALL CLIMBER_ENERGY_GRADIENT(NATOMS,COORDS,CLIMBER_ENERGY,CLIMBER_GRADATOMS)
3057:                   IF(DEBUG) WRITE(*,'(A,F15.5,A,F15.5)') ' potential> AMBER energy: ',ENERGY,', climber energy: ', CLIMBER_ENERGY 
3058:                   GRADATOMS=GRADATOMS+CLIMBER_GRADATOMS3057:                   GRADATOMS=GRADATOMS+CLIMBER_GRADATOMS
3059:                   ENERGY=ENERGY+CLIMBER_ENERGY3058:                   ENERGY=ENERGY+CLIMBER_ENERGY
3060:                ENDIF3059:                ENDIF
3061:             END IF3060:             END IF
3062:             VNEW(1:3*NATOMS) = GRADATOMS(:)3061:             VNEW(1:3*NATOMS) = GRADATOMS(:)
3063: ! Calculate the numerical hessian3062: ! Calculate the numerical hessian
3064:             IF (STEST) THEN3063:             IF (STEST) THEN
3065:                IF (.NOT. ALLOCATED(HESS)) ALLOCATE(HESS(3*NATOMS, 3*NATOMS))3064:                IF (.NOT. ALLOCATED(HESS)) ALLOCATE(HESS(3*NATOMS, 3*NATOMS))
3066:                CALL AMBER12_NUM_HESS(NATOMS, COORDS, DELTA=1.0D-5, HESSIAN=HESS(:, :))3065:                CALL AMBER12_NUM_HESS(NATOMS, COORDS, DELTA=1.0D-5, HESSIAN=HESS(:, :))
3067:             END IF3066:             END IF


legend
Lines Added 
Lines changed
 Lines Removed

hdiff - version: 2.1.0