hdiff output

r31704/key.f90 2017-01-21 10:35:11.818250000 +0000 r31703/key.f90 2017-01-21 10:35:13.470250000 +0000
 47:      &        QSPCFWT, QTIP4PFT, CFUSIONT, DUMPINTXYZ, DUMPINTEOS, INTLJT, INTTST, EYTRAPT, OHCELLT, MKTRAPT, & 47:      &        QSPCFWT, QTIP4PFT, CFUSIONT, DUMPINTXYZ, DUMPINTEOS, INTLJT, INTTST, EYTRAPT, OHCELLT, MKTRAPT, &
 48:      &        INTFREEZET, LPERMDIST, CHECKNEGATIVET, CHECKOVERLAPT, ACK1, ACK2, CONDATT, USERPOTT, & 48:      &        INTFREEZET, LPERMDIST, CHECKNEGATIVET, CHECKOVERLAPT, ACK1, ACK2, CONDATT, USERPOTT, &
 49:      &        CONCUTFRACT, CONCUTABST, ENDNUMHESS2, CHARMMDFTBT, PAIRCOLOURT, REVERSEUPHILLT, WHOLEDNEB, & 49:      &        CONCUTFRACT, CONCUTABST, ENDNUMHESS2, CHARMMDFTBT, PAIRCOLOURT, REVERSEUPHILLT, WHOLEDNEB, &
 50:      &        NONEBMAX, READMASST, ONEDAPBCT, ONEDPBCT, INVTONEDPBCT, INVTTWODPBCT, TWODAPBCT, TWODPBCT, THREEDAPBCT, & 50:      &        NONEBMAX, READMASST, ONEDAPBCT, ONEDPBCT, INVTONEDPBCT, INVTTWODPBCT, TWODAPBCT, TWODPBCT, THREEDAPBCT, &
 51:      &        THREEDPBCT, FOURDAPBCT, FOURDPBCT, MODEDOWNT, CHEMSHIFT, TTM3T, & 51:      &        THREEDPBCT, FOURDAPBCT, FOURDPBCT, MODEDOWNT, CHEMSHIFT, TTM3T, &
 52:      &        NOINVERSION, INVERTPT, KNOWVECS, PMPATHT, AAORIENTT, MULTIJOBT, QUIPARGSTRT, QUIPPARAMST, HESSDUMPT, & 52:      &        NOINVERSION, INVERTPT, KNOWVECS, PMPATHT, AAORIENTT, MULTIJOBT, QUIPARGSTRT, QUIPPARAMST, HESSDUMPT, &
 53:      &        CLASSICALRATEST, TSPLITTINGT, HESSREADT, INSTANTONOPTT,INSTANTONSTARTDUMPT,VARSTEPOPTT, MOLPRO, REAXFFT, & 53:      &        CLASSICALRATEST, TSPLITTINGT, HESSREADT, INSTANTONOPTT,INSTANTONSTARTDUMPT,VARSTEPOPTT, MOLPRO, REAXFFT, &
 54:      &        EIGENONLY,OVERCONV, GLJT,CLSTRINGT,CLSTRINGTST, PHI4MODT, EX1DT, MCPATHT, MCBIAST, RPHT, TWISTT, MCPATH2T, & 54:      &        EIGENONLY,OVERCONV, GLJT,CLSTRINGT,CLSTRINGTST, PHI4MODT, EX1DT, MCPATHT, MCBIAST, RPHT, TWISTT, MCPATH2T, &
 55:      &        PBST, SSHT, GAUSSIAN03, CPPNEBT, CUDAT, CUDATIMET, TRUSTMODET,MODELOST, METRICTENSOR, INTSPRINGACTIVET, & 55:      &        PBST, SSHT, GAUSSIAN03, CPPNEBT, CUDAT, CUDATIMET, TRUSTMODET,MODELOST, METRICTENSOR, INTSPRINGACTIVET, &
 56:      &        PERMGUESS, QCIPERMCHECK, DUMPFRQST, MULTIPOTT, MLP3T, MLPB3T, DUMPBESTPATH, ALIGNRBST, AVOID_COLLISIONS, MLPPROB, & 56:      &        PERMGUESS, QCIPERMCHECK, DUMPFRQST, MULTIPOTT, MLP3T, MLPB3T, DUMPBESTPATH, ALIGNRBST, AVOID_COLLISIONS, MLPPROB, &
 57:      &        MALONALDEHYDE, MLPNEWREG, DJWRBT, STEALTHYT, STEALTV, LJADDT, MLPB3NEWT, MLPVB3T, & 57:      &        MALONALDEHYDE, MLPNEWREG, DJWRBT, STEALTHYT, STEALTV, LJADDT, MLPB3NEWT, &
 58:      &        QCIPOTT, QCIPOT2T, QCIRADSHIFTT, QCINOREPINT, QCIAMBERT, SLERPT, NOTRANSROTT, MAXGAPT, BULKBOXT, GDSQT, FLATTESTT, & 58:      &        QCIPOTT, QCIPOT2T, QCIRADSHIFTT, QCINOREPINT, QCIAMBERT, SLERPT, NOTRANSROTT, MAXGAPT, BULKBOXT, GDSQT, FLATTESTT, &
 59:      &        MLQT, MLQPROB, LJADD2T, MACROIONT, NOREGBIAS 59:      &        MLQT, MLQPROB, LJADD2T, MACROIONT
 60:  60: 
 61: ! sy349 > for testing the flatpath after dneb 61: ! sy349 > for testing the flatpath after dneb
 62:       !LOGICAL, ALLOCATABLE :: FLATPATHT(:) 62:       !LOGICAL, ALLOCATABLE :: FLATPATHT(:)
 63:       LOGICAL FLATPATHT 63:       LOGICAL FLATPATHT
 64:  64: 
 65: ! bf269 > polymer in a pore (non-bonding (LJ) energy from neighbours is not subtracted) 65: ! bf269 > polymer in a pore (non-bonding (LJ) energy from neighbours is not subtracted)
 66:       LOGICAL :: PORE8T = .FALSE. ! add 8th power cylindrical pore to the potential? 66:       LOGICAL :: PORE8T = .FALSE. ! add 8th power cylindrical pore to the potential?
 67:       INTEGER :: PORE8_AXIS = 3 ! principal axis of the cylindric pore (1:x, 2:y, 3:z) 67:       INTEGER :: PORE8_AXIS = 3 ! principal axis of the cylindric pore (1:x, 2:y, 3:z)
 68:       DOUBLE PRECISION :: PORE8_ENERGY = 1.0d1 ! energy of the pore when radius = 1 68:       DOUBLE PRECISION :: PORE8_ENERGY = 1.0d1 ! energy of the pore when radius = 1
 69:       LOGICAL :: HARMPOLYT = .FALSE. ! add harmonic bonds between the beads 69:       LOGICAL :: HARMPOLYT = .FALSE. ! add harmonic bonds between the beads


r31704/keywords.f 2017-01-21 10:35:12.098250000 +0000 r31703/keywords.f 2017-01-21 10:35:13.706250000 +0000
943: !943: !
944:          STEALTHYT=.FALSE.944:          STEALTHYT=.FALSE.
945:          STEALTV=.FALSE.945:          STEALTV=.FALSE.
946: 946: 
947: !947: !
948: ! Neural network potential948: ! Neural network potential
949: !949: !
950:          MLP3T=.FALSE.950:          MLP3T=.FALSE.
951:          MLPB3T=.FALSE.951:          MLPB3T=.FALSE.
952:          MLPB3NEWT=.FALSE.952:          MLPB3NEWT=.FALSE.
953:          MLPVB3T=.FALSE. 
954:          NOREGBIAS=.FALSE. 
955:          MLPNEWREG=.FALSE.953:          MLPNEWREG=.FALSE.
956:          MLPPROB=.FALSE.954:          MLPPROB=.FALSE.
957:          MLPDONE=.FALSE.955:          MLPDONE=.FALSE.
958:          MLPNORM=.FALSE.956:          MLPNORM=.FALSE.
959:          MLPLAMBDA=0.0D0957:          MLPLAMBDA=0.0D0
960:          MLPDATSTART=1958:          MLPDATSTART=1
961: !959: !
962: ! ML quadratic function960: ! ML quadratic function
963: !961: !
964:          MLQT=.FALSE.962:          MLQT=.FALSE.
4054:             IF (NMLP.NE.NATOMS) THEN4052:             IF (NMLP.NE.NATOMS) THEN
4055:                PRINT '(A,2I8)', 'keywords> ERROR *** NATOMS,NMLP=',NATOMS,NMLP4053:                PRINT '(A,2I8)', 'keywords> ERROR *** NATOMS,NMLP=',NATOMS,NMLP
4056:                STOP4054:                STOP
4057:             ENDIF4055:             ENDIF
4058:    4056:    
4059:             LUNIT=GETUNIT()4057:             LUNIT=GETUNIT()
4060:             OPEN(LUNIT,FILE='MLPdata',STATUS='OLD')4058:             OPEN(LUNIT,FILE='MLPdata',STATUS='OLD')
4061:             ALLOCATE(MLPDAT(MLPDATA,MLPIN),MLPOUTCOME(MLPDATA),MLPMEAN(MLPIN))4059:             ALLOCATE(MLPDAT(MLPDATA,MLPIN),MLPOUTCOME(MLPDATA),MLPMEAN(MLPIN))
4062:             MLPMEAN(1:MLPIN)=0.0D04060:             MLPMEAN(1:MLPIN)=0.0D0
4063:             DO J1=1,MLPDATA4061:             DO J1=1,MLPDATA
4064:                READ(LUNIT,*) MLPOUTCOME(J1),(DUMMY,J2=1,MLPSTART-1),MLPDAT(J1,1:MLPIN) 
4065:                MLPOUTCOME(J1)=MLPOUTCOME(J1)+1 ! to shift the range from 0 to from 1 
4066:                DO J2=1,MLPIN 
4067:                   MLPMEAN(J2)=MLPMEAN(J2)+ABS(MLPDAT(J1,J2)) 
4068:                ENDDO 
4069:             ENDDO 
4070:             CLOSE(LUNIT) 
4071:             IF (MLPNORM) THEN 
4072:                MLPMEAN(1:MLPIN)=MLPMEAN(1:MLPIN)/MLPDATA 
4073:                WRITE(*,'(A)') 'keyword> Rescaling inputs by mean absolute values:' 
4074:                WRITE(*,'(6G20.10)') MLPMEAN(1:MLPIN) 
4075:                DO J1=1,MLPIN 
4076:                   MLPDAT(1:MLPDATA,J1)=MLPDAT(1:MLPDATA,J1)/MLPMEAN(J1) 
4077:                ENDDO 
4078:             ENDIF 
4079:             DEALLOCATE(MLPMEAN) 
4080:             MLPDONE=.TRUE. 
4081:          ELSE IF (WORD.EQ.'NOREGBIAS') THEN 
4082:             NOREGBIAS=.TRUE. 
4083:          ELSE IF (WORD.EQ.'MLPVB3') THEN 
4084:             MLPVB3T=.TRUE. 
4085:             CALL READI(MLPIN)      ! number of inputs (data items after outcome) 
4086:             CALL READI(MLPSTART) ! starting position in data list, not counting outcome 
4087:             CALL READI(MLPHIDDEN) 
4088:             CALL READI(MLPOUT) 
4089:             CALL READI(MLPDATA) 
4090:             IF (NITEMS.GT.5) CALL READF(MLPLAMBDA) 
4091:             WRITE(*,'(A,5I8,G20.10)') ' keywords> MLP3 vector bias nodes and Nin, Ninstart, Nhidden, Nout, Ndata, lambda=', 
4092:      &                                MLPIN,MLPSTART,MLPHIDDEN,MLPOUT,MLPDATA,MLPLAMBDA 
4093:             NMLP=MLPHIDDEN*(MLPIN+MLPOUT)+MLPHIDDEN+MLPOUT 
4094:             IF (NMLP.NE.NATOMS) THEN 
4095:                PRINT '(A,2I8)', 'keywords> ERROR *** NATOMS,NMLP=',NATOMS,NMLP 
4096:                STOP 
4097:             ENDIF 
4098:     
4099:             LUNIT=GETUNIT() 
4100:             OPEN(LUNIT,FILE='MLPdata',STATUS='OLD') 
4101:             ALLOCATE(MLPDAT(MLPDATA,MLPIN),MLPOUTCOME(MLPDATA),MLPMEAN(MLPIN)) 
4102:             MLPMEAN(1:MLPIN)=0.0D0 
4103:             DO J1=1,MLPDATA 
4104:                READ(LUNIT,*) MLPOUTCOME(J1),(DUMMY,J2=1,MLPSTART-1),MLPDAT(J1,1:MLPIN)4062:                READ(LUNIT,*) MLPOUTCOME(J1),(DUMMY,J2=1,MLPSTART-1),MLPDAT(J1,1:MLPIN)
4105:                MLPOUTCOME(J1)=MLPOUTCOME(J1)+1 ! to shift the range from 0 to from 14063:                MLPOUTCOME(J1)=MLPOUTCOME(J1)+1 ! to shift the range from 0 to from 1
4106:                DO J2=1,MLPIN4064:                DO J2=1,MLPIN
4107:                   MLPMEAN(J2)=MLPMEAN(J2)+ABS(MLPDAT(J1,J2))4065:                   MLPMEAN(J2)=MLPMEAN(J2)+ABS(MLPDAT(J1,J2))
4108:                ENDDO4066:                ENDDO
4109:             ENDDO4067:             ENDDO
4110:             CLOSE(LUNIT)4068:             CLOSE(LUNIT)
4111:             IF (MLPNORM) THEN4069:             IF (MLPNORM) THEN
4112:                MLPMEAN(1:MLPIN)=MLPMEAN(1:MLPIN)/MLPDATA4070:                MLPMEAN(1:MLPIN)=MLPMEAN(1:MLPIN)/MLPDATA
4113:                WRITE(*,'(A)') 'keyword> Rescaling inputs by mean absolute values:'4071:                WRITE(*,'(A)') 'keyword> Rescaling inputs by mean absolute values:'


r31704/ljadd.f 2017-01-21 10:35:12.326250000 +0000 r31703/ljadd.f 2017-01-21 10:35:13.930250000 +0000
231:      2                 R8(N,N), G(N,N), XG(N,N),231:      2                 R8(N,N), G(N,N), XG(N,N),
232:      3                 R14(N,N), F(N,N), DUMMY, DUMMYX, DUMMYY, DUMMYZ, DIST, XMUL2232:      3                 R14(N,N), F(N,N), DUMMY, DUMMYX, DUMMYY, DUMMYZ, DIST, XMUL2
233: C233: C
234: C  Store distance matrices.234: C  Store distance matrices.
235: C235: C
236: !     WRITE(MYUNIT,'(A)') 'coords in LJADD2:'236: !     WRITE(MYUNIT,'(A)') 'coords in LJADD2:'
237: !     WRITE(MYUNIT,'(3G20.10)') X(1:3*N)237: !     WRITE(MYUNIT,'(3G20.10)') X(1:3*N)
238:       ENERGY=0.0D0238:       ENERGY=0.0D0
239:       IF (GTEST.AND.(.NOT.STEST)) THEN239:       IF (GTEST.AND.(.NOT.STEST)) THEN
240:          DO J1=1,N240:          DO J1=1,N
241:             MJ1=MOD(J1-1,NADDTARGET)+1241:             MJ1=MOD(J1,NADDTARGET)+1
242:             J3=3*J1242:             J3=3*J1
243:             XG(J1,J1)=0.0D0243:             XG(J1,J1)=0.0D0
244:             DO J2=J1+1,N244:             DO J2=J1+1,N
245:                MJ2=MOD(J2-1,NADDTARGET)+1245:                MJ2=MOD(J2,NADDTARGET)+1
246:                J4=3*J2246:                J4=3*J2
247:                DIST=(X(J3-2)-X(J4-2))**2+(X(J3-1)-X(J4-1))**2+(X(J3)-X(J4))**2247:                DIST=(X(J3-2)-X(J4-2))**2+(X(J3-1)-X(J4-1))**2+(X(J3)-X(J4))**2
248:                DIST=1.0D0/DIST248:                DIST=1.0D0/DIST
249:                R6=DIST**3249:                R6=DIST**3
250:                IF (MJ1.EQ.MJ2) THEN ! repulsive part only250:                IF (MJ1.EQ.MJ2) THEN ! repulsive part only
251:                   DUMMY=LJADDEPS(MJ2,MJ1)*R6*R6251:                   DUMMY=LJADDEPS(MJ2,MJ1)*R6*R6
252:                ELSE252:                ELSE
253:                   DUMMY=LJADDEPS(MJ2,MJ1)*R6*(R6-1.0D0)253:                   DUMMY=LJADDEPS(MJ2,MJ1)*R6*(R6-1.0D0)
254:                ENDIF254:                ENDIF
255:                ENERGY=ENERGY+DUMMY255:                ENERGY=ENERGY+DUMMY
258:                IF (MJ1.EQ.MJ2) THEN ! repulsive part only258:                IF (MJ1.EQ.MJ2) THEN ! repulsive part only
259:                   XG(J2,J1)=-LJADDEPS(MJ2,MJ1)*24.0D0*2.0D0*R6*DIST259:                   XG(J2,J1)=-LJADDEPS(MJ2,MJ1)*24.0D0*2.0D0*R6*DIST
260:                ELSE260:                ELSE
261:                   XG(J2,J1)=-LJADDEPS(MJ2,MJ1)*24.0D0*(2.0D0*R6-1.0D0)*DIST261:                   XG(J2,J1)=-LJADDEPS(MJ2,MJ1)*24.0D0*(2.0D0*R6-1.0D0)*DIST
262:                ENDIF262:                ENDIF
263:                XG(J1,J2)=XG(J2,J1)263:                XG(J1,J2)=XG(J2,J1)
264:             ENDDO264:             ENDDO
265:          ENDDO265:          ENDDO
266:       ELSEIF (GTEST) THEN266:       ELSEIF (GTEST) THEN
267:          DO J1=1,N267:          DO J1=1,N
268:             MJ1=MOD(J1-1,NADDTARGET)+1268:             MJ1=MOD(J1,NADDTARGET)+1
269:             XG(J1,J1)=0.0D0269:             XG(J1,J1)=0.0D0
270:             R2(J1,J1)=0.0D0270:             R2(J1,J1)=0.0D0
271:             R8(J1,J1)=0.0D0271:             R8(J1,J1)=0.0D0
272:             R14(J1,J1)=0.0D0272:             R14(J1,J1)=0.0D0
273:             DO J2=J1+1,N273:             DO J2=J1+1,N
274:                MJ2=MOD(J2-1,NADDTARGET)+1274:                MJ2=MOD(J2,NADDTARGET)+1
275:                R2(J2,J1)=(X(3*(J1-1)+1)-X(3*(J2-1)+1))**2275:                R2(J2,J1)=(X(3*(J1-1)+1)-X(3*(J2-1)+1))**2
276:      1                  +(X(3*(J1-1)+2)-X(3*(J2-1)+2))**2276:      1                  +(X(3*(J1-1)+2)-X(3*(J2-1)+2))**2
277:      2                  +(X(3*(J1-1)+3)-X(3*(J2-1)+3))**2277:      2                  +(X(3*(J1-1)+3)-X(3*(J2-1)+3))**2
278:                R2(J2,J1)=1.0D0/R2(J2,J1)278:                R2(J2,J1)=1.0D0/R2(J2,J1)
279:                R6=R2(J2,J1)**3279:                R6=R2(J2,J1)**3
280:                IF (MJ1.EQ.MJ2) THEN ! repulsive part only280:                IF (MJ1.EQ.MJ2) THEN ! repulsive part only
281:                   DUMMY=LJADDEPS(MJ2,MJ1)*R6*R6281:                   DUMMY=LJADDEPS(MJ2,MJ1)*R6*R6
282:                ELSE282:                ELSE
283:                   DUMMY=LJADDEPS(MJ2,MJ1)*R6*(R6-1.0D0)283:                   DUMMY=LJADDEPS(MJ2,MJ1)*R6*(R6-1.0D0)
284:                ENDIF284:                ENDIF
289:                IF (MJ1.EQ.MJ2) THEN ! repulsive part only289:                IF (MJ1.EQ.MJ2) THEN ! repulsive part only
290:                   XG(J2,J1)=-LJADDEPS(MJ2,MJ1)*24.0D0*2.0D0*R6*R2(J1,J2)*R6290:                   XG(J2,J1)=-LJADDEPS(MJ2,MJ1)*24.0D0*2.0D0*R6*R2(J1,J2)*R6
291:                ELSE291:                ELSE
292:                   XG(J2,J1)=-LJADDEPS(MJ2,MJ1)*24.0D0*(2.0D0*R6-1.0D0)*R2(J1,J2)*R6292:                   XG(J2,J1)=-LJADDEPS(MJ2,MJ1)*24.0D0*(2.0D0*R6-1.0D0)*R2(J1,J2)*R6
293:                ENDIF293:                ENDIF
294:                XG(J1,J2)=XG(J2,J1)294:                XG(J1,J2)=XG(J2,J1)
295:             ENDDO295:             ENDDO
296:          ENDDO296:          ENDDO
297:       ELSE297:       ELSE
298:          DO J1=1,N298:          DO J1=1,N
299:             MJ1=MOD(J1-1,NADDTARGET)+1299:             MJ1=MOD(J1,NADDTARGET)+1
300:             J3=3*(J1-1)300:             J3=3*(J1-1)
301:             DO J2=J1+1,N301:             DO J2=J1+1,N
302:                MJ2=MOD(J2-1,NADDTARGET)+1302:                MJ2=MOD(J2,NADDTARGET)+1
303:                J4=3*(J2-1)303:                J4=3*(J2-1)
304:                R2T=(X(J3+1)-X(J4+1))**2+(X(J3+2)-X(J4+2))**2+(X(J3+3)-X(J4+3))**2304:                R2T=(X(J3+1)-X(J4+1))**2+(X(J3+2)-X(J4+2))**2+(X(J3+3)-X(J4+3))**2
305:                R2T=1.0D0/R2T305:                R2T=1.0D0/R2T
306:                R6=R2T**3306:                R6=R2T**3
307:                IF (MJ1.EQ.MJ2) THEN ! repulsive part only307:                IF (MJ1.EQ.MJ2) THEN ! repulsive part only
308:                   ENERGY=ENERGY+LJADDEPS(MJ2,MJ1)*R6*R6308:                   ENERGY=ENERGY+LJADDEPS(MJ2,MJ1)*R6*R6
309:                ELSE309:                ELSE
310:                   ENERGY=ENERGY+LJADDEPS(MJ2,MJ1)*R6*(R6-1.0D0)310:                   ENERGY=ENERGY+LJADDEPS(MJ2,MJ1)*R6*(R6-1.0D0)
311:                ENDIF311:                ENDIF
312:             ENDDO312:             ENDDO
347:       IMPLICIT NONE347:       IMPLICIT NONE
348:       INTEGER N, J1, J2, J3, J4, J5, J6, MJ1, MJ2348:       INTEGER N, J1, J2, J3, J4, J5, J6, MJ1, MJ2
349:       DOUBLE PRECISION G(N,N), R14(N,N), R8(N,N),349:       DOUBLE PRECISION G(N,N), R14(N,N), R8(N,N),
350:      1                 R2(N,N), F(N,N),350:      1                 R2(N,N), F(N,N),
351:      2                 X(3*N),DUMMY351:      2                 X(3*N),DUMMY
352: 352: 
353: C353: C
354: C  Calculate the g tensor.354: C  Calculate the g tensor.
355: C355: C
356:       DO J1=1,N356:       DO J1=1,N
357:          MJ1=MOD(J1-1,NADDTARGET)+1357:          MJ1=MOD(J1,NADDTARGET)+1
358:          G(J1,J1)=0.0D0358:          G(J1,J1)=0.0D0
359:          DO J2=J1+1,N359:          DO J2=J1+1,N
360:             MJ2=MOD(J2-1,NADDTARGET)+1360:             MJ2=MOD(J2,NADDTARGET)+1
361:             IF (MJ1.EQ.MJ2) THEN ! repulsive part only361:             IF (MJ1.EQ.MJ2) THEN ! repulsive part only
362:                G(J2,J1)=-LJADDEPS(MJ2,MJ1)*24.0D0*2.0D0*R14(J2,J1)362:                G(J2,J1)=-LJADDEPS(MJ2,MJ1)*24.0D0*2.0D0*R14(J2,J1)
363:             ELSE363:             ELSE
364:                G(J2,J1)=-LJADDEPS(MJ2,MJ1)*24.0D0*(2.0D0*R14(J2,J1)-R8(J2,J1))364:                G(J2,J1)=-LJADDEPS(MJ2,MJ1)*24.0D0*(2.0D0*R14(J2,J1)-R8(J2,J1))
365:             ENDIF365:             ENDIF
366:             G(J1,J2)=G(J2,J1)366:             G(J1,J2)=G(J2,J1)
367:          ENDDO367:          ENDDO
368:       ENDDO368:       ENDDO
369: 369: 
370:       DO J1=1,N370:       DO J1=1,N
371:          MJ1=MOD(J1-1,NADDTARGET)+1371:          MJ1=MOD(J1,NADDTARGET)+1
372:          F(J1,J1)=0.0D0372:          F(J1,J1)=0.0D0
373:          DO J2=J1+1,N373:          DO J2=J1+1,N
374:             MJ2=MOD(J2-1,NADDTARGET)+1374:             MJ2=MOD(J2,NADDTARGET)+1
375:             IF (MJ1.EQ.MJ2) THEN ! repulsive part only375:             IF (MJ1.EQ.MJ2) THEN ! repulsive part only
376:                F(J2,J1)=LJADDEPS(MJ2,MJ1)*672.0D0*R14(J2,J1)376:                F(J2,J1)=LJADDEPS(MJ2,MJ1)*672.0D0*R14(J2,J1)
377:             ELSE377:             ELSE
378:                F(J2,J1)=LJADDEPS(MJ2,MJ1)*672.0D0*R14(J2,J1)-LJADDEPS(MJ2,MJ1)*192.0D0*R8(J2,J1)378:                F(J2,J1)=LJADDEPS(MJ2,MJ1)*672.0D0*R14(J2,J1)-LJADDEPS(MJ2,MJ1)*192.0D0*R8(J2,J1)
379:             ENDIF379:             ENDIF
380:             F(J1,J2)=F(J2,J1)380:             F(J1,J2)=F(J2,J1)
381:          ENDDO381:          ENDDO
382:       ENDDO382:       ENDDO
383: C383: C
384: C  Now do the hessian. First are the entirely diagonal terms.384: C  Now do the hessian. First are the entirely diagonal terms.


r31704/MLPVB3.f90 2017-01-21 10:35:10.942250000 +0000 r31703/MLPVB3.f90 2017-01-21 10:35:12.786250000 +0000
  1: SUBROUTINE MLPVB3(X,V,ENERGY,GTEST,SECT)  1: svn: E195012: Unable to find repository location for 'svn+ssh://svn.ch.private.cam.ac.uk/groups/wales/trunk/OPTIM/source/MLPVB3.f90' in revision 31703
  2: USE KEY 
  3: USE MODHESS 
  4: USE COMMONS, ONLY : DEBUG 
  5: IMPLICIT NONE 
  6: LOGICAL GTEST,SECT 
  7: DOUBLE PRECISION X(NMLP), V(NMLP), ENERGY, DUMMY1, DUMMY2, DUMMY3, DUMMY4, SUMINPUTS, AREA 
  8: DOUBLE PRECISION Y(MLPOUT), PROB(MLPOUT), PMLPOUTJ1, DMAX 
  9: DOUBLE PRECISION DYW1G(MLPHIDDEN), DPCW1BG(MLPOUT,MLPHIDDEN) 
 10: DOUBLE PRECISION DYW2G(MLPOUT,MLPHIDDEN,MLPIN), DPCW2BG(MLPHIDDEN,MLPIN), TANHSUM(MLPHIDDEN), SECH2(MLPHIDDEN) 
 11: DOUBLE PRECISION DYCDBHID(MLPOUT), DPCDBHID(MLPOUT), DPCDBHIDOUTJ1, D2PCDBHID2, D2YCDBHID2(MLPOUT) 
 12: DOUBLE PRECISION D2YCDBHIDDW2BG, D2PCDBHIDDW2BG 
 13: DOUBLE PRECISION DPCDBHIDDW1BG(MLPOUT,MLPHIDDEN) 
 14: DOUBLE PRECISION DPEDW2BG(MLPOUT,MLPHIDDEN,MLPIN), DPEDW1BG(MLPOUT,MLPOUT,MLPHIDDEN), DPEDWBHJ(MLPOUT,MLPHIDDEN) 
 15: DOUBLE PRECISION, PARAMETER :: BOUT=0.0D0 
 16: DOUBLE PRECISION, ALLOCATABLE :: PSAVE(:,:) 
 17: INTEGER MLPOUTJ1, MLPOFFSET, BOFFSET1, BOFFSET2, NREG 
 18: INTEGER GETUNIT, LUNIT, J1, J2, J3, J4, K4, K2, K3, J5, J6 
 19: CHARACTER(LEN=132) FNAME 
 20:  
 21: ! 
 22: ! Variables are ordered  
 23: ! w^2_{jk} at (j-1)*MLPIN+k 
 24: !   up to MLPHIDDEN*MLPIN, then 
 25: ! w^1_{ij} at MLPHIDDEN*MLPIN + (i-1)*MLPHIDDEN+j 
 26: !   up to MLPHIDDEN*MLPIN + MLPOUT*MLPHIDDEN 
 27: ! w^bh_j at MLPHIDDEN*(MLPIN+MLPOUT)+1 to MLPHIDDEN*(MLPIN+MLPOUT)+MLPHIDDEN 
 28: ! w^bo_i at MLPHIDDEN*(MLPIN+MLPOUT)+MLPHIDDEN+1 to MLPHIDDEN*(MLPIN+MLPOUT)+MLPHIDDEN+MLPOUT 
 29: ! 
 30:  
 31: IF (MLPPROB) THEN 
 32:    ALLOCATE(PSAVE(MLPDATA,MLPOUT)) 
 33: ENDIF 
 34: MLPOFFSET=MLPHIDDEN*MLPIN 
 35: BOFFSET1=MLPHIDDEN*(MLPIN+MLPOUT) 
 36: BOFFSET2=MLPHIDDEN*(MLPIN+MLPOUT+1) 
 37: ENERGY=0.0D0 
 38: V(1:NMLP)=0.0D0 
 39: IF (SECT) HESS(1:NMLP,1:NMLP)=0.0D0 
 40: DO J1=1,MLPDATA 
 41:    MLPOUTJ1=MLPOUTCOME(J1) 
 42:    DO J2=1,MLPHIDDEN 
 43:       DUMMY1=0.0D0 
 44:       DO J3=1,MLPIN 
 45:          DUMMY1=DUMMY1+X((J2-1)*MLPIN+J3)*MLPDAT(J1,J3) 
 46:       ENDDO 
 47:       DUMMY1=DUMMY1+X(BOFFSET1+J2) 
 48:       TANHSUM(J2)=TANH(DUMMY1)  
 49:       DYW1G(J2)=TANHSUM(J2) 
 50: !     PRINT *,'DUMMY1=',DUMMY1 
 51: !     IF (ABS(DUMMY1).GT.20.0D0) THEN 
 52: !        SECH2(J2)=0.0D0 
 53: !     ELSE 
 54:          SECH2(J2)=1.0D0/COSH(DUMMY1)**2  
 55: !     ENDIF 
 56:    ENDDO 
 57: !  DUMMY3=0.0D0 
 58:    DMAX=-1.0D100 
 59:    DO J4=1,MLPOUT 
 60:       DUMMY2=0.0D0 
 61:       DO J2=1,MLPHIDDEN 
 62:          DO J3=1,MLPIN 
 63:             DYW2G(J4,J2,J3)=(X( MLPOFFSET + (J4-1)*MLPHIDDEN + J2 )) * MLPDAT(J1,J3)*SECH2(J2) 
 64:          ENDDO 
 65:          DUMMY2=DUMMY2+X(MLPOFFSET+(J4-1)*MLPHIDDEN+J2)*TANHSUM(J2) 
 66:       ENDDO 
 67:       DUMMY2=DUMMY2+X(BOFFSET2+J4) 
 68: !     IF (DUMMY2.GT.50.0D0) DUMMY2=50.0D0 ! to prevent FPE 
 69:       IF (DUMMY2.GT.DMAX) DMAX=DUMMY2 
 70:       Y(J4)=DUMMY2 
 71: !     DUMMY3=DUMMY3+EXP(DUMMY2) 
 72:    ENDDO   
 73: ! 
 74: ! Factor DMAX out of the exponentials to prevent overflow. 
 75: ! 
 76:    DUMMY3=0.0D0 
 77:    DO J4=1,MLPOUT 
 78:       Y(J4)=Y(J4)-DMAX 
 79:       DUMMY3=DUMMY3+EXP(Y(J4)) 
 80:    ENDDO 
 81:    DO J4=1,MLPOUT 
 82:       PROB(J4)=EXP(Y(J4))/DUMMY3 
 83:    ENDDO 
 84:    PMLPOUTJ1=PROB(MLPOUTJ1) 
 85:    IF (MLPPROB) THEN 
 86: !     WRITE(*,'(A,I8,A)') 'MLP3> data point ',J1,' outputs and probabilities:' 
 87: !     WRITE(LUNIT,'(8G15.5)') Y(1:MLPOUT),PROB(1:MLPOUT), MLPOUTJ1 
 88: !     WRITE(*,'(100G20.10)') PROB(1:MLPOUT), MLPOUTJ1*1.0D0 
 89:       PSAVE(J1,1:MLPOUT)=PROB(1:MLPOUT) 
 90:    ENDIF 
 91:    ENERGY=ENERGY-LOG(PMLPOUTJ1) 
 92:    IF (GTEST) THEN 
 93: ! 
 94: ! We only need the probability derivative for the probability corresponding to the correct outcome for this data point 
 95: ! 
 96:       DPCW1BG(1:MLPOUT,1:MLPHIDDEN)=0.0D0 
 97:       DO J2=1,MLPHIDDEN 
 98:          DO J4=1,MLPOUT 
 99:             DPCW1BG(J4,J2)=DPCW1BG(J4,J2)-PMLPOUTJ1*PROB(J4)*DYW1G(J2)  
100:          ENDDO 
101:          DPCW1BG(MLPOUTJ1,J2)=DPCW1BG(MLPOUTJ1,J2)+PMLPOUTJ1*DYW1G(J2) ! d p_c / dw^1_beta gamma 
102:       ENDDO 
103:  
104:       DO J3=1,MLPIN 
105:          DO J2=1,MLPHIDDEN 
106:             DUMMY3=0.0D0 
107:             DO J4=1,MLPOUT 
108:                DUMMY3=DUMMY3+PROB(J4)*DYW2G(J4,J2,J3) 
109:             ENDDO 
110:             DPCW2BG(J2,J3)=PMLPOUTJ1*(DYW2G(MLPOUTJ1,J2,J3)-DUMMY3) 
111:          ENDDO 
112:       ENDDO 
113:       DO J4=1,MLPOUT 
114:          V(BOFFSET2+J4)=V(BOFFSET2+J4)+PROB(J4) 
115:       ENDDO 
116:       V(BOFFSET2+MLPOUTJ1)=V(BOFFSET2+MLPOUTJ1)-1.0D0 
117:  
118:       SUMINPUTS=0.0D0 
119:       DO J3=1,MLPIN 
120:          SUMINPUTS=SUMINPUTS+MLPDAT(J1,J3) 
121:       ENDDO 
122:  
123:       DO J4=1,MLPOUT 
124:          DUMMY2=0.0D0 
125:          DO J2=1,MLPHIDDEN 
126:             DUMMY2=DUMMY2+X(MLPOFFSET+(J4-1)*MLPHIDDEN+J2)*SECH2(J2) 
127:          ENDDO 
128:          DYCDBHID(J4)=DUMMY2*SUMINPUTS 
129:       ENDDO 
130:  
131:       DUMMY2=0.0D0 
132:       DO J4=1,MLPOUT 
133:          DUMMY2=DUMMY2+PROB(J4)*DYCDBHID(J4) 
134:       ENDDO 
135:  
136:       DPCDBHIDOUTJ1=PMLPOUTJ1*(DYCDBHID(MLPOUTJ1)-DUMMY2) 
137:           
138:       DO J4=1,MLPOUT 
139:          DPCDBHID(J4)=PROB(J4)*(DYCDBHID(J4)-DUMMY2) 
140:       ENDDO 
141:  
142:       DO J4=1,MLPOUT 
143:          DO J2=1,MLPHIDDEN 
144:             V(MLPOFFSET+(J4-1)*MLPHIDDEN+J2)=V(MLPOFFSET+(J4-1)*MLPHIDDEN+J2)-DPCW1BG(J4,J2)/PMLPOUTJ1 
145:          ENDDO 
146:       ENDDO 
147:  
148:       DO J2=1,MLPHIDDEN 
149:          DUMMY2=0.0D0 
150:          DO J4=1,MLPOUT 
151:             DUMMY2=DUMMY2-X(MLPOFFSET+(J4-1)*MLPHIDDEN+J2)*PROB(J4) 
152:          ENDDO 
153:          V(BOFFSET1+J2)=V(BOFFSET1+J2)-SECH2(J2)*(X(MLPOFFSET+(MLPOUTJ1-1)*MLPHIDDEN+J2)+DUMMY2) 
154:       ENDDO 
155:  
156:       DO J3=1,MLPIN 
157:          DO J2=1,MLPHIDDEN 
158:             V((J2-1)*MLPIN+J3)=V((J2-1)*MLPIN+J3)-DPCW2BG(J2,J3)/PMLPOUTJ1 
159:          ENDDO 
160:       ENDDO 
161:    ENDIF 
162:  
163:    IF (SECT) THEN 
164: ! 
165: ! This block w^1 with w^1 is locally symmetric 
166: ! 
167:       DO J4=1,MLPOUT ! J4 is beta  
168:          DO J2=1,MLPHIDDEN ! J2 is gamma 
169:             DO K4=1,J4 ! K4 is alpha 
170:                DO K2=1,MLPHIDDEN ! K2 is epsilon 
171:                   DUMMY1=0.0D0 
172:                   IF ((J4.EQ.MLPOUTJ1).AND.(K4.EQ.MLPOUTJ1)) DUMMY1=1.0D0 
173:                   IF (J4.EQ.MLPOUTJ1) DUMMY1=DUMMY1-PROB(K4) 
174:                   IF (K4.EQ.MLPOUTJ1) DUMMY1=DUMMY1-PROB(J4) 
175:                   IF (K4.EQ.J4) DUMMY1=DUMMY1-PROB(J4) 
176:                   DUMMY1=DUMMY1+2.0D0*PROB(J4)*PROB(K4) 
177:                   HESS(MLPOFFSET+(J4-1)*MLPHIDDEN+J2,MLPOFFSET+(K4-1)*MLPHIDDEN+K2)= & 
178:   &               HESS(MLPOFFSET+(J4-1)*MLPHIDDEN+J2,MLPOFFSET+(K4-1)*MLPHIDDEN+K2) & ! sum over data points 
179:   &               +DPCW1BG(J4,J2)*DPCW1BG(K4,K2)/PMLPOUTJ1**2 & 
180:   &               -DUMMY1*TANHSUM(J2)*TANHSUM(K2) 
181:                ENDDO 
182:             ENDDO 
183:          ENDDO 
184:       ENDDO 
185: ! 
186: ! Off-diagonal w^1 with w^2 blocks 
187: ! 
188:       DO J3=1,MLPOUT ! J3 is beta for w^1 outputs 
189:          DO J2=1,MLPHIDDEN ! J2 is gamma for w^1 hidden 
190:             DO K4=1,MLPHIDDEN ! K4 is alpha for w^2 hidden 
191:                DUMMY3=0.0D0 
192:                DO J5=1,MLPOUT 
193:                   DUMMY3=DUMMY3+PROB(J5)*X(MLPOFFSET + (J5-1)*MLPHIDDEN + K4)  
194:                ENDDO 
195:                DO K2=1,MLPIN ! K2 is epsilon for w^2 inputs 
196:                   DUMMY1=0.0D0 
197:                   IF (K4.EQ.J2) DUMMY1=DUMMY1-PMLPOUTJ1*PROB(J3)*MLPDAT(J1,K2)*SECH2(J2) 
198:                   IF ((K4.EQ.J2).AND.(J3.EQ.MLPOUTJ1)) DUMMY1=DUMMY1+PMLPOUTJ1*MLPDAT(J1,K2)*SECH2(J2) 
199:  
200:                   DUMMY2=TANHSUM(J2)*PMLPOUTJ1*MLPDAT(J1,K2)*SECH2(K4) & 
201:   &                      *(X(MLPOFFSET+(MLPOUTJ1-1)*MLPHIDDEN+K4)-DUMMY3) 
202:                   DUMMY1=DUMMY1-PROB(J3)*DUMMY2 
203:                   IF (MLPOUTJ1.EQ.J3) DUMMY1=DUMMY1+DUMMY2 
204:                   DUMMY1=DUMMY1-PMLPOUTJ1*PROB(J3)*MLPDAT(J1,K2)*SECH2(K4)*TANHSUM(J2) & 
205:   &                             *(X(MLPOFFSET + (J3-1)*MLPHIDDEN + K4)-DUMMY3) 
206:                   HESS(MLPOFFSET+(J3-1)*MLPHIDDEN+J2,(K4-1)*MLPIN+K2)= & 
207:   &               HESS(MLPOFFSET+(J3-1)*MLPHIDDEN+J2,(K4-1)*MLPIN+K2) & ! sum over data points 
208:   &               +DPCW1BG(J3,J2)*DPCW2BG(K4,K2)/PMLPOUTJ1**2 & 
209:   &               -DUMMY1/PMLPOUTJ1 
210:                ENDDO 
211:             ENDDO 
212:          ENDDO 
213:       ENDDO 
214: ! 
215: ! diagonal w^2 with w^2  
216: ! 
217:       DO J3=1,MLPIN ! J3 is gamma for w^2 inputs 
218:          DO J2=1,MLPHIDDEN ! J2 is beta for w^2 hidden 
219:             DUMMY2=0.0D0 
220:             DO J5=1,MLPOUT 
221:                DUMMY2=DUMMY2+PROB(J5)*X(MLPOFFSET + (J5-1)*MLPHIDDEN + J2)  
222:             ENDDO 
223:             DO K4=1,MLPIN ! K4 is epsilon for w^2 inputs 
224:                DO K2=1,J2 ! K2 is alpha for w^2 hidden 
225:                   DUMMY3=0.0D0 
226:                   DO J5=1,MLPOUT 
227:                      DUMMY3=DUMMY3+PROB(J5)*X(MLPOFFSET + (J5-1)*MLPHIDDEN + K2)  
228:                   ENDDO 
229:                   DUMMY4=0.0D0 
230:                   DO J5=1,MLPOUT 
231:                      DUMMY4=DUMMY4+PROB(J5)*X(MLPOFFSET + (J5-1)*MLPHIDDEN + K2)  & ! take out of loops 
232:   &                                        *X(MLPOFFSET + (J5-1)*MLPHIDDEN + J2) 
233:                   ENDDO 
234:                   DUMMY1=DPCW2BG(K2,K4)*MLPDAT(J1,J3)*SECH2(J2) & 
235:   &                      *(X(MLPOFFSET+(MLPOUTJ1-1)*MLPHIDDEN+J2)-DUMMY2) & 
236:   &                      -PMLPOUTJ1*MLPDAT(J1,J3)*SECH2(J2)*MLPDAT(J1,K4)*SECH2(K2) & 
237:   &                      *(DUMMY4-DUMMY2*DUMMY3) 
238:                   IF (K2.EQ.J2) DUMMY1=DUMMY1-2.0D0*PMLPOUTJ1*MLPDAT(J1,K4)*MLPDAT(J1,J3) & 
239:   &                                     *SECH2(J2)*TANHSUM(J2)*(X(MLPOFFSET + (MLPOUTJ1-1)*MLPHIDDEN + J2)-DUMMY2)  
240:  
241:                   HESS((J2-1)*MLPIN+J3,(K2-1)*MLPIN+K4)= & 
242:   &               HESS((J2-1)*MLPIN+J3,(K2-1)*MLPIN+K4) & ! sum over data points 
243:   &               +DPCW2BG(J2,J3)*DPCW2BG(K2,K4)/PMLPOUTJ1**2 - DUMMY1/PMLPOUTJ1 
244:                ENDDO 
245:             ENDDO 
246:          ENDDO 
247:       ENDDO 
248: ! 
249: ! w^2_{jk} at (j-1)*MLPIN+k 
250: !   up to MLPHIDDEN*MLPIN, then 
251: ! w^1_{ij} at MLPHIDDEN*MLPIN + (i-1)*MLPHIDDEN+j 
252: !   up to MLPHIDDEN*MLPIN + MLPOUT*MLPHIDDEN 
253: ! w^bh_j at MLPHIDDEN*(MLPIN+MLPOUT)+1 to MLPHIDDEN*(MLPIN+MLPOUT)+MLPHIDDEN 
254: ! w^bo_i at MLPHIDDEN*(MLPIN+MLPOUT)+MLPHIDDEN+1 to MLPHIDDEN*(MLPIN+MLPOUT)+MLPHIDDEN+MLPOUT 
255: ! 
256:  
257:       DO J5=1,MLPOUT ! eps 
258:          DO J3=1,MLPHIDDEN ! gamma   
259:             DO J2=1,MLPOUT ! beta 
260:                DPEDW1BG(J5,J2,J3)=-PROB(J5)*TANHSUM(J3)*PROB(J2) 
261:                IF (J2.EQ.J5) DPEDW1BG(J5,J2,J3)=DPEDW1BG(J5,J2,J3)+PROB(J5)*TANHSUM(J3) 
262:             ENDDO 
263:          ENDDO 
264:       ENDDO 
265:  
266:       DO J5=1,MLPOUT ! eps 
267:          DO J3=1,MLPIN ! gamma  
268:             DO J2=1,MLPHIDDEN ! beta 
269:                DUMMY3=0.0D0 
270:                DO J4=1,MLPOUT 
271:                   DUMMY3=DUMMY3+PROB(J4)*X(MLPOFFSET+(J4-1)*MLPHIDDEN+J2) 
272:                ENDDO 
273:                DPEDW2BG(J5,J2,J3)=PROB(J5)*SECH2(J2)*MLPDAT(J1,J3)*(X(MLPOFFSET+(J5-1)*MLPHIDDEN+J2)-DUMMY3) 
274:             ENDDO 
275:          ENDDO 
276:       ENDDO 
277:  
278:       DO J5=1,MLPOUT ! eps 
279:          DO J2=1,MLPHIDDEN ! j 
280:             DUMMY3=0.0D0 
281:             DO J4=1,MLPOUT 
282:                DUMMY3=DUMMY3+PROB(J4)*X(MLPOFFSET+(J4-1)*MLPHIDDEN+J2) 
283:             ENDDO 
284:             DPEDWBHJ(J5,J2)=PROB(J5)*SECH2(J2)*(X(MLPOFFSET+(J5-1)*MLPHIDDEN+J2)-DUMMY3) 
285:          ENDDO 
286:       ENDDO 
287: ! 
288: ! w^bo_i with w^1_bg (output hidden)  
289: ! 
290:       DO J5=1,MLPOUT 
291:          DO J3=1,MLPHIDDEN ! gamma 
292:             DO J2=1,MLPOUT ! beta 
293:                HESS(BOFFSET2+J5,MLPOFFSET+(J2-1)*MLPHIDDEN+J3)=HESS(BOFFSET2+J5,MLPOFFSET+(J2-1)*MLPHIDDEN+J3)+DPEDW1BG(J5,J2,J3) 
294:             ENDDO 
295:          ENDDO 
296:       ENDDO 
297: ! 
298: ! w^bo_i with w^2_bg (hidden input)  
299: ! 
300:       DO J5=1,MLPOUT ! i 
301:          DO J2=1,MLPHIDDEN ! beta 
302:             DO J3=1,MLPIN  ! gamma 
303:                HESS(BOFFSET2+J5,(J2-1)*MLPIN+J3)=HESS(BOFFSET2+J5,(J2-1)*MLPIN+J3)+DPEDW2BG(J5,J2,J3) 
304:             ENDDO 
305:          ENDDO 
306:       ENDDO 
307: ! 
308: ! w^bo_i with w^bo_j  
309: ! 
310:       DO J5=1,MLPOUT 
311:          DO J2=1,J5 
312:             HESS(BOFFSET2+J5,BOFFSET2+J2)=HESS(BOFFSET2+J5,BOFFSET2+J2)-PROB(J5)*PROB(J2) 
313:          ENDDO 
314:          HESS(BOFFSET2+J5,BOFFSET2+J5)=HESS(BOFFSET2+J5,BOFFSET2+J5)+PROB(J5) 
315:       ENDDO 
316: ! 
317: ! w^bo_i with w^bh_j  
318: ! 
319:       DO J2=1,MLPHIDDEN ! j 
320:          DUMMY1=0.0D0 
321:          DO J5=1,MLPOUT 
322:             DUMMY1=DUMMY1+PROB(J5)*X(MLPOFFSET+(J5-1)*MLPHIDDEN+J2) 
323:          ENDDO 
324:          DO J5=1,MLPOUT ! i 
325:             HESS(BOFFSET2+J5,BOFFSET1+J2)=HESS(BOFFSET2+J5,BOFFSET1+J2)+PROB(J5)*SECH2(J2)*(X(MLPOFFSET+(J5-1)*MLPHIDDEN+J2)-DUMMY1) 
326:          ENDDO 
327:       ENDDO 
328: ! 
329: ! w^bh_j with w^1_bg 
330: ! 
331:       DO J6=1,MLPOUT ! beta 
332:          DO J2=1,MLPHIDDEN ! j 
333:             DO J3=1,MLPHIDDEN ! gamma 
334:                DUMMY1=0.0D0 
335:                DO J5=1,MLPOUT 
336:                   DUMMY1=DUMMY1+X(MLPOFFSET+(J5-1)*MLPHIDDEN+J2)*DPEDW1BG(J5,J6,J3) 
337:                ENDDO 
338:                HESS(BOFFSET1+J2,MLPOFFSET+(J6-1)*MLPHIDDEN+J3)=HESS(BOFFSET1+J2,MLPOFFSET+(J6-1)*MLPHIDDEN+J3)+SECH2(J2)*DUMMY1 
339:                IF (J2.EQ.J3) THEN 
340:                   HESS(BOFFSET1+J2,MLPOFFSET+(J6-1)*MLPHIDDEN+J3)=HESS(BOFFSET1+J2,MLPOFFSET+(J6-1)*MLPHIDDEN+J3)+SECH2(J2)*PROB(J6) 
341:                   IF (J6.EQ.MLPOUTJ1) THEN 
342:                      HESS(BOFFSET1+J2,MLPOFFSET+(J6-1)*MLPHIDDEN+J3)=HESS(BOFFSET1+J2,MLPOFFSET+(J6-1)*MLPHIDDEN+J3)-SECH2(J2) 
343:                   ENDIF 
344:                ENDIF 
345:  
346: !              HESS(BOFFSET1+J2,MLPOFFSET+(J6-1)*MLPHIDDEN+J3)=HESS(BOFFSET1+J2,MLPOFFSET+(J6-1)*MLPHIDDEN+J3)+TANHSUM(J3)*DPEDWBHJ(J6,J2) 
347: !              IF (J6.EQ.MLPOUTJ1) THEN 
348: !                 IF (J3.EQ.J2) THEN 
349: !                    HESS(BOFFSET1+J2,MLPOFFSET+(J6-1)*MLPHIDDEN+J3)=HESS(BOFFSET1+J2,MLPOFFSET+(J6-1)*MLPHIDDEN+J3)-SECH2(J3)  
350: !                 ENDIF 
351: !              ENDIF 
352:             ENDDO 
353:          ENDDO 
354:       ENDDO 
355: ! 
356: ! w^bh_j hidden with w^bh k  
357: ! 
358:       DO J6=1,MLPHIDDEN ! j 
359:          DUMMY1=0.0D0 
360:          DO J4=1,MLPOUT 
361:             DUMMY1=DUMMY1+PROB(J4)*X(MLPOFFSET+(J4-1)*MLPHIDDEN+J6) 
362:          ENDDO 
363:          DO J2=1,MLPHIDDEN ! k 
364:             DUMMY2=0.0D0 
365:             DO J4=1,MLPOUT 
366:                DUMMY2=DUMMY2+X(MLPOFFSET+(J4-1)*MLPHIDDEN+J6)*DPEDWBHJ(J4,J2) 
367:             ENDDO 
368:             HESS(BOFFSET1+J6,BOFFSET1+J2)=HESS(BOFFSET1+J6,BOFFSET1+J2)+SECH2(J6)*DUMMY2 
369:             IF (J2.EQ.J6) THEN 
370:                HESS(BOFFSET1+J6,BOFFSET1+J2)=HESS(BOFFSET1+J6,BOFFSET1+J2)+2.0D0*SECH2(J6)*TANHSUM(J6)*(X(MLPOFFSET+(MLPOUTJ1-1)*MLPHIDDEN+J6)-DUMMY1) 
371:             ENDIF 
372:          ENDDO 
373:       ENDDO 
374: ! 
375: ! w^bh_j hidden with w^2 beta gamma (hidden in) 
376: ! 
377:       DO J6=1,MLPHIDDEN ! j 
378:          DUMMY1=0.0D0 
379:          DO J4=1,MLPOUT 
380:             DUMMY1=DUMMY1+PROB(J4)*X(MLPOFFSET+(J4-1)*MLPHIDDEN+J6) 
381:          ENDDO 
382:  
383:          DO J2=1,MLPHIDDEN ! J2 is beta 
384:             DO J3=1,MLPIN  ! J3 is gamma 
385:  
386:                DUMMY4=0.0D0 
387:                DO J4=1,MLPOUT 
388:                   DUMMY4=DUMMY4+X(MLPOFFSET+(J4-1)*MLPHIDDEN+J6)*DPEDW2BG(J4,J2,J3) 
389:                ENDDO 
390:  
391:                HESS(BOFFSET1+J6,(J2-1)*MLPIN+J3)=HESS(BOFFSET1+J6,(J2-1)*MLPIN+J3) + SECH2(J6)*DUMMY4 
392:                IF (J2.EQ.J6) THEN 
393:                   HESS(BOFFSET1+J6,(J2-1)*MLPIN+J3)=HESS(BOFFSET1+J6,(J2-1)*MLPIN+J3)+2.0D0*SECH2(J6)*TANHSUM(J6)*MLPDAT(J1,J3) & 
394:   &                  * (X(MLPOFFSET+(MLPOUTJ1-1)*MLPHIDDEN+J6)-DUMMY1) 
395:                ENDIF 
396:             ENDDO 
397:          ENDDO 
398:       ENDDO 
399:  
400:    ENDIF ! end of sect block 
401:  
402: ENDDO ! end of sum over data items counter J1 
403:  
404: NREG=NMLP 
405: ! 
406: ! If NOREGBIAS is true then exclude regularisation over bias weights, which start at 
407: ! variable BOFFSET1 + 1 
408: ! 
409: IF (NOREGBIAS) NREG=BOFFSET1 
410:  
411: DUMMY1=0.0D0 
412: DO J1=1,NREG 
413:    DUMMY1=DUMMY1+X(J1)**2 
414: ENDDO 
415:  
416: ENERGY=ENERGY/MLPDATA + MLPLAMBDA*DUMMY1 
417: ! IF (DEBUG) WRITE(*,'(A,G20.10)') 'MLP3> objective function=',ENERGY 
418:  
419: IF (MLPPROB) THEN 
420: !    CLOSE(LUNIT) 
421: !    WRITE(*,'(A)') ' MLPB3> predicted probabilities written to file probabilities' 
422: !    PRINT '(A)',' MLPB3> finished writing probabilities' 
423:    CALL ROC(PSAVE,AREA) 
424:    PRINT '(A,2G20.10)','energy, AUC=',ENERGY,AREA 
425: ENDIF 
426: IF (MLPPROB) DEALLOCATE(PSAVE) 
427:  
428: IF (GTEST) V(1:NMLP)=V(1:NMLP)/MLPDATA  
429: IF (GTEST) V(1:NREG)=V(1:NREG) + 2.0D0*MLPLAMBDA*X(1:NREG) 
430: ! 
431: ! Symmetrise Hessian here 
432: ! 
433: IF (SECT) HESS(1:NMLP,1:NMLP)=HESS(1:NMLP,1:NMLP)/MLPDATA 
434: IF (SECT) THEN 
435:    DO J1=1,NREG 
436:       HESS(J1,J1)=HESS(J1,J1)+2*MLPLAMBDA 
437: !     IF (MLPNEWREG) HESS(J1,J1)=HESS(J1,J1)+6.0D0*MLPLAMBDA/MAX(1.0D-100,X(J1)**4) 
438:    ENDDO 
439:    DO J1=1,NMLP 
440:       DO J2=1,J1-1 
441:          HESS(J2,J1)=HESS(J1,J2) 
442:       ENDDO 
443:    ENDDO 
444: ENDIF 
445:  
446: END SUBROUTINE MLPVB3 


r31704/MLQ.f90 2017-01-21 10:35:11.318250000 +0000 r31703/MLQ.f90 2017-01-21 10:35:13.006250000 +0000
 73:       Y(J4)=Y(J4)-DMAX 73:       Y(J4)=Y(J4)-DMAX
 74:       DUMMY3=DUMMY3+EXP(Y(J4)) 74:       DUMMY3=DUMMY3+EXP(Y(J4))
 75:    ENDDO 75:    ENDDO
 76:    DO J4=1,MLQOUT 76:    DO J4=1,MLQOUT
 77:       PROB(J4)=EXP(Y(J4))/DUMMY3 77:       PROB(J4)=EXP(Y(J4))/DUMMY3
 78:    ENDDO 78:    ENDDO
 79:    PMLQOUTJ1=PROB(MLQOUTJ1) ! this is p_c(d), the probability predicted for the known outcome for data point d 79:    PMLQOUTJ1=PROB(MLQOUTJ1) ! this is p_c(d), the probability predicted for the known outcome for data point d
 80:    IF (MLQPROB) THEN 80:    IF (MLQPROB) THEN
 81:       PSAVE(J1,1:MLQOUT)=PROB(1:MLQOUT) 81:       PSAVE(J1,1:MLQOUT)=PROB(1:MLQOUT)
 82:    ENDIF 82:    ENDIF
 83:    IF (PMLQOUTJ1.LE.0.0D0) THEN 
 84:        IF (DEBUG) WRITE(MYUNIT,'(A,G20.10)') 'WARNING PROB=',PROB(MLQOUTJ1) 
 85:        PMLQOUTJ1=1.0D-100 
 86: !      WRITE(MYUNIT,'(I10,G20.10)') (J4,PROB(J4),J4=1,MLQOUT) 
 87: !      STOP 
 88:    ENDIF 
 89:    ENERGY=ENERGY-LOG(PMLQOUTJ1) ! summing over contributions for data points, indexed by J1 83:    ENERGY=ENERGY-LOG(PMLQOUTJ1) ! summing over contributions for data points, indexed by J1
 90:  84: 
 91:    IF (GTEST) THEN 85:    IF (GTEST) THEN
 92: ! 86: !
 93: ! We only need the probability derivatives for the probability corresponding to the correct outcome for this data point 87: ! We only need the probability derivatives for the probability corresponding to the correct outcome for this data point
 94: ! 88: !
 95:       DO J3=1,NMLQ ! for each variable w 89:       DO J3=1,NMLQ ! for each variable w
 96:          DUMMY3=0.0D0 90:          DUMMY3=0.0D0
 97:          DO J4=1,MLQOUT ! the only non-zero terms in this sum will correspond to variables associated with the given output 91:          DO J4=1,MLQOUT ! the only non-zero terms in this sum will correspond to variables associated with the given output
 98:             DUMMY3=DUMMY3+PROB(J4)*DYDW(J3,J4) 92:             DUMMY3=DUMMY3+PROB(J4)*DYDW(J3,J4)


r31704/OPTIM.F 2017-01-21 10:35:11.586250000 +0000 r31703/OPTIM.F 2017-01-21 10:35:13.234250000 +0000
452:             WRITE(*,'(A,G20.10)') ' OPTIM> Using translational/rotational ev shift=',SHIFTV452:             WRITE(*,'(A,G20.10)') ' OPTIM> Using translational/rotational ev shift=',SHIFTV
453:          ELSE453:          ELSE
454:             WRITE(*,'(A)') ' OPTIM> z rotational ev shifting'454:             WRITE(*,'(A)') ' OPTIM> z rotational ev shifting'
455:          ENDIF455:          ENDIF
456:          SHIFTL(1)=SHIFTV456:          SHIFTL(1)=SHIFTV
457:          SHIFTL(2)=SHIFTV457:          SHIFTL(2)=SHIFTV
458:          SHIFTL(3)=SHIFTV458:          SHIFTL(3)=SHIFTV
459:          SHIFTL(4)=SHIFTV459:          SHIFTL(4)=SHIFTV
460:          SHIFTL(5)=SHIFTV460:          SHIFTL(5)=SHIFTV
461:          SHIFTL(6)=SHIFTV461:          SHIFTL(6)=SHIFTV
462:       ELSE IF (MIEFT.OR.MLP3T.OR.MLPB3T.OR.MLPVB3T.OR.NOTRANSROTT) THEN462:       ELSE IF (MIEFT.OR.MLP3T.OR.MLPB3T.OR.NOTRANSROTT) THEN
463:          WRITE(*,'(A,G20.10)') ' OPTIM> No eigenvalues shifted'463:          WRITE(*,'(A,G20.10)') ' OPTIM> No eigenvalues shifted'
464:       ELSE464:       ELSE
465:          WRITE(*,'(A,G20.10)') ' OPTIM> Using translational/rotational ev shift=',SHIFTV465:          WRITE(*,'(A,G20.10)') ' OPTIM> Using translational/rotational ev shift=',SHIFTV
466:          DO J1=1,6466:          DO J1=1,6
467:             SHIFTL(J1)=SHIFTV467:             SHIFTL(J1)=SHIFTV
468:          ENDDO468:          ENDDO
469:       ENDIF469:       ENDIF
470: 470: 
471:       !js850> 471:       !js850> 
472:       IF (ZSYM(NATOMS).EQ.'SS') THEN472:       IF (ZSYM(NATOMS).EQ.'SS') THEN


r31704/potential.f 2017-01-21 10:35:12.554250000 +0000 r31703/potential.f 2017-01-21 10:35:14.162250000 +0000
488:             CALL example1Dpotential(COORDS,VNEW,ENERGY,GTEST,STEST)488:             CALL example1Dpotential(COORDS,VNEW,ENERGY,GTEST,STEST)
489:             489:             
490:             IF (PTEST) THEN490:             IF (PTEST) THEN
491:                WRITE(*,10) ' potential> Energy for last cycle=',ENERGY,' '491:                WRITE(*,10) ' potential> Energy for last cycle=',ENERGY,' '
492:                WRITE(ESTRING,10) 'Energy for last cycle=',ENERGY,' '492:                WRITE(ESTRING,10) 'Energy for last cycle=',ENERGY,' '
493:             ENDIF493:             ENDIF
494: 494: 
495: 495: 
496:          ELSE IF (VARIABLES) THEN496:          ELSE IF (VARIABLES) THEN
497: 497: 
498:             IF (MLPVB3T) THEN498:             IF (MLPB3T) THEN
499:                CALL MLPVB3(COORDS, VNEW, ENERGY, GTEST, STEST)499:                CALL MLPB3(COORDS, VNEW, ENERGY, GTEST, STEST)
500: 500:             ELSEIF (MLP3T) THEN
 501:                CALL MLP3(COORDS, VNEW, ENERGY, GTEST, STEST)
 502:             ELSEIF (MLQT) THEN
 503:                CALL MLQ(COORDS, VNEW, ENERGY, GTEST, STEST)
501: !               DIFF=1.0D-4504: !               DIFF=1.0D-4
502: !               PRINT*,'analytic and numerical gradients:'505: !               PRINT*,'analytic and numerical gradients:'
503: !               IF (.NOT.(ALLOCATED(HESS))) ALLOCATE(HESS(NMLP,NMLP))506: !               IF (.NOT.(ALLOCATED(HESS))) ALLOCATE(HESS(NMLQ,NMLQ))
504: !               CALL MLPVB3(COORDS, VNEW, ENERGY, .TRUE., .TRUE.)507: !               CALL MLQ(COORDS, VNEW, ENERGY, .TRUE., .TRUE.)
505: !               PRINT '(A,I8)','SIZE(HESS)=',SIZE(HESS)508: !               PRINT '(A,I8)','SIZE(HESS)=',SIZE(HESS)
506: !               HESSDUM(1:NMLP,1:NMLP)=HESS(1:NMLP,1:NMLP)509: !               HESSDUM(1:NMLQ,1:NMLQ)=HESS(1:NMLQ,1:NMLQ)
507: !               DO J1=1,NATOMS510: !               DO J1=1,NATOMS
508: !                  COORDS(J1)=COORDS(J1)+DIFF511: !                  COORDS(J1)=COORDS(J1)+DIFF
509: !                  CALL MLPVB3(COORDS,VPLUS,EPLUS,.FALSE.,.FALSE.)512: !                  CALL MLQ(COORDS,VPLUS,EPLUS,.FALSE.,.FALSE.)
510: !                  COORDS(J1)=COORDS(J1)-2.0D0*DIFF513: !                  COORDS(J1)=COORDS(J1)-2.0D0*DIFF
511: !                  CALL MLPVB3(COORDS,VMINUS,EMINUS,.FALSE.,.FALSE.)514: !                  CALL MLQ(COORDS,VMINUS,EMINUS,.FALSE.,.FALSE.)
512: !                  COORDS(J1)=COORDS(J1)+DIFF515: !                  COORDS(J1)=COORDS(J1)+DIFF
513: !                  IF ((ABS(VNEW(J1)).NE.0.0D0).AND. 516: !                  IF ((ABS(VNEW(J1)).NE.0.0D0).AND. 
514: !     &               (ABS(100.0D0*(VNEW(J1)-(EPLUS-EMINUS)/(2.0D0*DIFF))/VNEW(J1)).GT.1.0D0)) THEN517: !    &                (ABS(100.0D0*(VNEW(J1)-(EPLUS-EMINUS)/(2.0D0*DIFF))/VNEW(J1)).GT.1.0D0)) THEN
515: !                     WRITE(*,'(A,I5,2F20.10,A)') 'anal and num ',J1,VNEW(J1),(EPLUS-EMINUS)/(2.0D0*DIFF),'   X'518: !                     WRITE(*,'(A,I5,2F20.10,A)') 'anal and num ',J1,VNEW(J1),(EPLUS-EMINUS)/(2.0D0*DIFF),'   X'
516: !                  ELSE519: !                  ELSE
517: !                     WRITE(*,'(A,I5,2F20.10)') 'anal and num ',J1,VNEW(J1),(EPLUS-EMINUS)/(2.0D0*DIFF)520: !                     WRITE(*,'(A,I5,2F20.10)') 'anal and num ',J1,VNEW(J1),(EPLUS-EMINUS)/(2.0D0*DIFF)
518: !                  ENDIF521: !                  ENDIF
519: !               ENDDO522: !               ENDDO
520: !               PRINT*,'analytic and numerical second derivatives:'523: !               PRINT*,'analytic and numerical second derivatives:'
521: !               DO J1=1,NATOMS524: !               DO J1=1,NATOMS
522: !                  COORDS(J1)=COORDS(J1)+DIFF525: !                  COORDS(J1)=COORDS(J1)+DIFF
523: !                  CALL MLPVB3(COORDS,VPLUS,EPLUS,.TRUE.,.FALSE.)526: !                  CALL MLQ(COORDS,VPLUS,EPLUS,.TRUE.,.FALSE.)
524: !                  COORDS(J1)=COORDS(J1)-2.0D0*DIFF527: !                  COORDS(J1)=COORDS(J1)-2.0D0*DIFF
525: !                  CALL MLPVB3(COORDS,VMINUS,EMINUS,.TRUE.,.FALSE.)528: !                  CALL MLQ(COORDS,VMINUS,EMINUS,.TRUE.,.FALSE.)
526: !                  COORDS(J1)=COORDS(J1)+DIFF529: !                  COORDS(J1)=COORDS(J1)+DIFF
527: !                  DO J2=1,NATOMS530: !                  DO J2=1,NATOMS
528: !                     IF ((ABS(HESS(J1,J2)).NE.0.0D0).AND. 531: !                     IF ((ABS(HESS(J1,J2)).NE.0.0D0).AND. 
529: !     &                  (ABS(100.0D0*(HESS(J1,J2)-(VPLUS(J2)-VMINUS(J2))/(2.0D0*DIFF))/HESS(J1,J2)).GT.1.0D0)) THEN532: !    &                   (ABS(100.0D0*(HESS(J1,J2)-(VPLUS(J2)-VMINUS(J2))/(2.0D0*DIFF))/HESS(J1,J2)).GT.1.0D0)) THEN
530: !                     WRITE(*,'(2I5,2G20.10,A)') J1,J2,HESS(J1,J2),(VPLUS(J2)-VMINUS(J2))/(2.0D0*DIFF),'   X'533: !                     WRITE(*,'(2I5,2G20.10,A)') J1,J2,HESS(J1,J2),(VPLUS(J2)-VMINUS(J2))/(2.0D0*DIFF),'   X'
531: !                     ELSE534: !                     ELSE
532: !                        WRITE(*,'(2I5,2G20.10,A)') J1,J2,HESS(J1,J2),(VPLUS(J2)-VMINUS(J2))/(2.0D0*DIFF)535: !                        WRITE(*,'(2I5,2G20.10,A)') J1,J2,HESS(J1,J2),(VPLUS(J2)-VMINUS(J2))/(2.0D0*DIFF)
533: !                     ENDIF536: !                     ENDIF
534: !                  ENDDO537: !                  ENDDO
535: !                ENDDO538: !               ENDDO
536: !                STOP539: !               STOP
537:             ELSEIF (MLPB3T) THEN 
538:                CALL MLPB3(COORDS, VNEW, ENERGY, GTEST, STEST) 
539:             ELSEIF (MLP3T) THEN 
540:                CALL MLP3(COORDS, VNEW, ENERGY, GTEST, STEST) 
541:             ELSEIF (MLQT) THEN 
542:                CALL MLQ(COORDS, VNEW, ENERGY, GTEST, STEST) 
543:             ELSE540:             ELSE
544:                CALL FUNCTIONAL( COORDS, VNEW, ENERGY, GTEST, STEST)541:                CALL FUNCTIONAL( COORDS, VNEW, ENERGY, GTEST, STEST)
545:             ENDIF542:             ENDIF
546:             IF (PTEST.OR.MLPPROB.OR.MLQPROB) THEN543:             IF (PTEST.OR.MLPPROB.OR.MLQPROB) THEN
547:                 WRITE(*,10) ' potential> Energy for last cycle=',ENERGY,'         '544:                 WRITE(*,10) ' potential> Energy for last cycle=',ENERGY,'         '
548:                 WRITE(ESTRING,10) 'Energy for last cycle=',ENERGY,'         '545:                 WRITE(ESTRING,10) 'Energy for last cycle=',ENERGY,'         '
549:             ENDIF546:             ENDIF
550: !           IF (MLPPROB) STOP547: !           IF (MLPPROB) STOP
551:             ! CALL CTEST(NATOMS, COORDS, VNEW, ENERGY, GTEST, STEST)548:             ! CALL CTEST(NATOMS, COORDS, VNEW, ENERGY, GTEST, STEST)
552:             ! CALL TWODFUNC(COORDS,VNEW,ENERGY,GTEST,STEST)549:             ! CALL TWODFUNC(COORDS,VNEW,ENERGY,GTEST,STEST)


legend
Lines Added 
Lines changed
 Lines Removed

hdiff - version: 2.1.0