hdiff output

r29875/geopt.f 2016-01-31 09:30:10.586285941 +0000 r29874/geopt.f 2016-01-31 09:30:11.374296353 +0000
836:                   ELSE836:                   ELSE
837:                      CALL DSYEV('V','U',NOPT,HESS,NOPT,EVALUES,TEMPA,9*NATOMS,INFO)837:                      CALL DSYEV('V','U',NOPT,HESS,NOPT,EVALUES,TEMPA,9*NATOMS,INFO)
838:                   ENDIF838:                   ENDIF
839:                ELSE 839:                ELSE 
840: ! hk286 - optional number of ENDHESS not yet implemented840: ! hk286 - optional number of ENDHESS not yet implemented
841:                   IF (RIGIDINIT) THEN841:                   IF (RIGIDINIT) THEN
842:                      RBAANORMALMODET = .TRUE.842:                      RBAANORMALMODET = .TRUE.
843:                      CALL GENRIGID_EIGENVALUES(Q, ATMASS, EVALUES, INFO)843:                      CALL GENRIGID_EIGENVALUES(Q, ATMASS, EVALUES, INFO)
844:                      RBAANORMALMODET = .FALSE.844:                      RBAANORMALMODET = .FALSE.
845:                   ELSE845:                   ELSE
846:                      IF ((NENDHESS.GE.NOPT).OR.(VARIABLES.AND.(NENDHESS.GE.NATOMS))) THEN846:                      IF ((NENDHESS.GE.3*NATOMS).OR.(VARIABLES.AND.(NENDHESS.GE.NATOMS))) THEN
847: C                       CALL DSYEV('N','U',NOPT,HESS,NOPT,EVALUES,TEMPA,9*NATOMS,INFO)847: C                       CALL DSYEV('N','U',NOPT,HESS,NOPT,EVALUES,TEMPA,9*NATOMS,INFO)
848: C  csw34> changed the call used here to be the same as if NENDHESS < NOPT below848: C  csw34> changed the call used here to be the same as if NENDHESS < NOPT below
849: C         as this leads to a 20% speed increase!849: C         as this leads to a 20% speed increase!
850:                         ABSTOL=DLAMCH('Safe  minimum')850:                         ABSTOL=DLAMCH('Safe  minimum')
851:                         ALLOCATE(ZWK(1,1)) ! not referenced for job type 'N'851:                         ALLOCATE(ZWK(1,1)) ! not referenced for job type 'N'
852:                         IF (RBAAT) THEN ! hk286852:                         IF (RBAAT) THEN ! hk286
853:                            RBAANORMALMODET = .TRUE.853:                            RBAANORMALMODET = .TRUE.
854:                            CALL NRMLMD (Q, EVALUES, .TRUE.)854:                            CALL NRMLMD (Q, EVALUES, .TRUE.)
855:                            RBAANORMALMODET = .FALSE.855:                            RBAANORMALMODET = .FALSE.
856:                         ELSE856:                         ELSE
875:                ENDIF875:                ENDIF
876:             ENDIF876:             ENDIF
877: C877: C
878: C  MASSWT2 and MASSWTRP do not mass weight Q and VNEW, but MASSWT does. Need to undo this878: C  MASSWT2 and MASSWTRP do not mass weight Q and VNEW, but MASSWT does. Need to undo this
879: C  if DUMPDATAT is .TRUE. for comparison with pathsample (which uses unit masses).879: C  if DUMPDATAT is .TRUE. for comparison with pathsample (which uses unit masses).
880: C  Probably best to always undo it, in case we need non-mass-weighted Q somewhere880: C  Probably best to always undo it, in case we need non-mass-weighted Q somewhere
881: C  further down.881: C  further down.
882: C882: C
883: C           IF (DUMPDATAT.AND.(.NOT.(CHRMMT.OR.RINGPOLYMERT))) THEN883: C           IF (DUMPDATAT.AND.(.NOT.(CHRMMT.OR.RINGPOLYMERT))) THEN
884:             IF (.NOT. (CHRMMT .OR. RBAAT)) THEN ! hk286884:             IF (.NOT. (CHRMMT .OR. RBAAT)) THEN ! hk286
885:                IF (.NOT.(RIGIDINIT.OR.VARIABLES)) THEN ! jmc49 bugfix, because the mass weighting is not done if rigidinit 885:                IF (.NOT.RIGIDINIT) THEN ! jmc49 bugfix, because the mass weighting is not done if rigidinit 
886:                   DO J1=1,NATOMS886:                   DO J1=1,NATOMS
887:                      AMASS=1/SQRT(ATMASS(J1))887:                      AMASS=1/SQRT(ATMASS(J1))
888:                      J3=3*J1888:                      J3=3*J1
889:                      Q(J3-2)=AMASS*Q(J3-2)889:                      Q(J3-2)=AMASS*Q(J3-2)
890:                      Q(J3-1)=AMASS*Q(J3-1)890:                      Q(J3-1)=AMASS*Q(J3-1)
891:                      Q(J3)=AMASS*Q(J3)891:                      Q(J3)=AMASS*Q(J3)
892:                   ENDDO892:                   ENDDO
893:                ENDIF893:                ENDIF
894:             ENDIF894:             ENDIF
895: 895: 
1061:                ENDIF 1061:                ENDIF 
1062:             ENDDO1062:             ENDDO
1063:          ELSE1063:          ELSE
1064:             IF (NENDHESS-NEXMODES.GT.0) THEN1064:             IF (NENDHESS-NEXMODES.GT.0) THEN
1065:                MINFRQ2=LOG(EVALUES(NENDHESS-NEXMODES))1065:                MINFRQ2=LOG(EVALUES(NENDHESS-NEXMODES))
1066:             ELSE1066:             ELSE
1067:                MINFRQ2=1.0D01067:                MINFRQ2=1.0D0
1068:             ENDIF1068:             ENDIF
1069:             EWARN=.FALSE.1069:             EWARN=.FALSE.
1070:             DO I1=1,NENDHESS - NEXMODES1070:             DO I1=1,NENDHESS - NEXMODES
1071:                IF (I1.GT.1) THEN1071:                IF (I1.GT.J1) THEN
1072:                   IF (EVALUES(I1-1).NE.0.0D0) THEN1072:                   IF (EVALUES(I1-1).NE.0.0D0) THEN
1073:                      IF (ABS(EVALUES(I1)/EVALUES(I1-1)).LT.1.0D-2) THEN1073:                      IF (ABS(EVALUES(I1)/EVALUES(I1-1)).LT.1.0D-2) THEN
1074:                         PRINT '(A,G20.10,A,G20.10)',' geopt> WARNING - decrease in magnitude of eigenvalues from ',EVALUES(I1-1),1074:                         PRINT '(A,G20.10,A,G20.10)',' geopt> WARNING - decrease in magnitude of eigenvalues from ',EVALUES(I1-1),
1075:      &                       ' to ',EVALUES(I1)1075:      &                       ' to ',EVALUES(I1)
1076:                         PRINT '(A)',' geopt> WARNING - this could indicate a stationary point of the wrong index'1076:                         PRINT '(A)',' geopt> WARNING - this could indicate a stationary point of the wrong index'
1077:                         EWARN=.TRUE.1077:                         EWARN=.TRUE.
1078:                      ENDIF1078:                      ENDIF
1079:                   ENDIF1079:                   ENDIF
1080:                ENDIF1080:                ENDIF
1081:                IF (EVALUES(I1).GT.0.0D0) THEN1081:                IF (EVALUES(I1).GT.0.0D0) THEN
1350:                 IF ((AMBERT.OR.NABT).AND.MWVECTORS) THEN1350:                 IF ((AMBERT.OR.NABT).AND.MWVECTORS) THEN
1351:                         CALL A9DUMPMODES(DIAG,ZT,NOPT,3*NATOMS)1351:                         CALL A9DUMPMODES(DIAG,ZT,NOPT,3*NATOMS)
1352:                 ENDIF1352:                 ENDIF
1353:               IF(NORMALMODET) CALL VISNORMODES(Q)1353:               IF(NORMALMODET) CALL VISNORMODES(Q)
1354:          ENDIF1354:          ENDIF
1355: 1355: 
1356:          PRINT*1356:          PRINT*
1357:          WRITE(*,'(A)') ' geopt>                          **** CONVERGED ****'1357:          WRITE(*,'(A)') ' geopt>                          **** CONVERGED ****'
1358:          PRINT*1358:          PRINT*
1359:          CALL FLUSH(6)1359:          CALL FLUSH(6)
1360:  
1361: !1360: !
1362: ! Print TIP frequencies in cm-1 for CoM/Euler angle representation for debugging angle-axis.1361: ! Print TIP frequencies in cm-1 for CoM/Euler angle representation for debugging angle-axis.
1363: !1362: !
1364:          IF (DEBUG.AND.(ZSYM(1)(1:1).EQ.'W')) THEN1363:          IF (DEBUG.AND.(ZSYM(1)(1:1).EQ.'W')) THEN
1365:             IF (ZSYM(1)(1:2).EQ.'W4') IPOT=41364:             IF (ZSYM(1)(1:2).EQ.'W4') IPOT=4
1366:             IF (ZSYM(1)(1:2).EQ.'W3') IPOT=31365:             IF (ZSYM(1)(1:2).EQ.'W3') IPOT=3
1367:             IF (ZSYM(1)(1:2).EQ.'W2') IPOT=21366:             IF (ZSYM(1)(1:2).EQ.'W2') IPOT=2
1368:             IF (ZSYM(1)(1:2).EQ.'W1') IPOT=11367:             IF (ZSYM(1)(1:2).EQ.'W1') IPOT=1
1369:             CALL H2OMODES(NATOMS/2,IPOT,Q,DIAG)1368:             CALL H2OMODES(NATOMS/2,IPOT,Q,DIAG)
1370:             PRINT '(A,I6,A)',' geopt> TIP normal mode frequencies in wavenumbers'1369:             PRINT '(A,I6,A)',' geopt> TIP normal mode frequencies in wavenumbers'


r29875/keywords.f 2016-01-31 09:30:10.786288583 +0000 r29874/keywords.f 2016-01-31 09:30:11.574298996 +0000
3689: !3689: !
3690:       ELSE IF ((WORD.EQ.'MLP3').OR.(WORD.EQ.'MLPB3')) THEN3690:       ELSE IF ((WORD.EQ.'MLP3').OR.(WORD.EQ.'MLPB3')) THEN
3691:          MLP3T=.TRUE.3691:          MLP3T=.TRUE.
3692:          IF (WORD.EQ.'MLPB3') MLPB3T=.TRUE.3692:          IF (WORD.EQ.'MLPB3') MLPB3T=.TRUE.
3693:          CALL READI(MLPIN)3693:          CALL READI(MLPIN)
3694:          CALL READI(MLPHIDDEN)3694:          CALL READI(MLPHIDDEN)
3695:          CALL READI(MLPOUT)3695:          CALL READI(MLPOUT)
3696:          CALL READI(MLPDATA)3696:          CALL READI(MLPDATA)
3697:          IF (NITEMS.GT.5) CALL READF(MLPLAMBDA)3697:          IF (NITEMS.GT.5) CALL READF(MLPLAMBDA)
3698:          IF (WORD.EQ.'MLPB3') THEN3698:          IF (WORD.EQ.'MLPB3') THEN
3699:             WRITE(*,'(A,4I8,G20.10)') 'MLPB3 potential with bias node and Nin, Nhidden, Nout, Ndata, lambda=',3699:             WRITE(*,'(A,4I8,G20.10)') 'MLPB3 potential with two bias nodes and Nin, Nhidden, Nout, Ndata, lambda=',
3700:      &                                   MLPIN,MLPHIDDEN,MLPOUT,MLPDATA,MLPLAMBDA3700:      &                                   MLPIN,MLPHIDDEN,MLPOUT,MLPDATA,MLPLAMBDA
3701:             NMLP=MLPHIDDEN*(MLPIN+MLPOUT)+13701:             NMLP=MLPHIDDEN*(MLPIN+MLPOUT)+2
3702:          ELSE3702:          ELSE
3703:             WRITE(*,'(A,4I8,G20.10)') 'MLP3 potential with Nin, Nhidden, Nout, Ndata, lambda=',3703:             WRITE(*,'(A,4I8,G20.10)') 'MLP3 potential with Nin, Nhidden, Nout, Ndata, lambda=',
3704:      &                                   MLPIN,MLPHIDDEN,MLPOUT,MLPDATA,MLPLAMBDA3704:      &                                   MLPIN,MLPHIDDEN,MLPOUT,MLPDATA,MLPLAMBDA
3705:             NMLP=MLPHIDDEN*(MLPIN+MLPOUT)3705:             NMLP=MLPHIDDEN*(MLPIN+MLPOUT)
3706:          ENDIF3706:          ENDIF
3707:          IF (NMLP.NE.NATOMS) THEN3707:          IF (NMLP.NE.NATOMS) THEN
3708:             PRINT '(A,2I8)', 'keywords> ERROR *** NATOMS,NMLP=',NATOMS,NMLP3708:             PRINT '(A,2I8)', 'keywords> ERROR *** NATOMS,NMLP=',NATOMS,NMLP
3709:             STOP3709:             STOP
3710:          ENDIF3710:          ENDIF
3711:          LUNIT=GETUNIT()3711:          LUNIT=GETUNIT()


r29875/MLPB3.f90 2016-01-31 09:30:10.390283350 +0000 r29874/MLPB3.f90 2016-01-31 09:30:11.182293816 +0000
  4: USE COMMONS, ONLY : DEBUG  4: USE COMMONS, ONLY : DEBUG
  5: IMPLICIT NONE  5: IMPLICIT NONE
  6: LOGICAL GTEST,SECT  6: LOGICAL GTEST,SECT
  7: DOUBLE PRECISION X(NMLP), V(NMLP), ENERGY, DUMMY1, DUMMY2, DUMMY3, DUMMY4, SUMINPUTS  7: DOUBLE PRECISION X(NMLP), V(NMLP), ENERGY, DUMMY1, DUMMY2, DUMMY3, DUMMY4, SUMINPUTS
  8: DOUBLE PRECISION Y(MLPOUT), PROB(MLPOUT), PMLPOUTJ1, DMAX  8: DOUBLE PRECISION Y(MLPOUT), PROB(MLPOUT), PMLPOUTJ1, DMAX
  9: DOUBLE PRECISION DYW1G(MLPHIDDEN), DPCW1BG(MLPOUT,MLPHIDDEN)  9: DOUBLE PRECISION DYW1G(MLPHIDDEN), DPCW1BG(MLPOUT,MLPHIDDEN)
 10: DOUBLE PRECISION DYW2G(MLPOUT,MLPHIDDEN,MLPIN), DPCW2BG(MLPHIDDEN,MLPIN), TANHSUM(MLPHIDDEN), SECH2(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) 11: DOUBLE PRECISION DYCDBHID(MLPOUT), DPCDBHID(MLPOUT), DPCDBHIDOUTJ1, D2PCDBHID2, D2YCDBHID2(MLPOUT)
 12: DOUBLE PRECISION DPCDW2BG(MLPOUT,MLPHIDDEN,MLPIN), D2YCDBHIDDW2BG, D2PCDBHIDDW2BG 12: DOUBLE PRECISION DPCDW2BG(MLPOUT,MLPHIDDEN,MLPIN), D2YCDBHIDDW2BG, D2PCDBHIDDW2BG
 13: DOUBLE PRECISION DPCDBHIDDW1BG(MLPOUT,MLPHIDDEN) 13: DOUBLE PRECISION DPCDBHIDDW1BG(MLPOUT,MLPHIDDEN)
 14: DOUBLE PRECISION, PARAMETER :: BOUT=0.0D0 
 15: INTEGER MLPOUTJ1, MLPOFFSET, BOFFSET 14: INTEGER MLPOUTJ1, MLPOFFSET, BOFFSET
 16: INTEGER GETUNIT, LUNIT, J1, J2, J3, J4, K4, K2, K3, J5 15: INTEGER GETUNIT, LUNIT, J1, J2, J3, J4, K4, K2, K3, J5
 17:  16: 
 18: ! 17: !
 19: ! Variables are ordered  18: ! Variables are ordered 
 20: ! w^2_{jk} at (j-1)*MLPIN+k 19: ! w^2_{jk} at (j-1)*MLPIN+k
 21: !   up to MLPHIDDEN*MLPIN, then 20: !   up to MLPHIDDEN*MLPIN, then
 22: ! w^1_{ij} at MLPHIDDEN*MLPIN + (i-1)*MLPHIDDEN+j 21: ! w^1_{ij} at MLPHIDDEN*MLPIN + (i-1)*MLPHIDDEN+j
 23: !   up to MLPHIDDEN*MLPIN + MLPOUT*MLPHIDDEN 22: !   up to MLPHIDDEN*MLPIN + MLPOUT*MLPHIDDEN
 24: ! bhidden at MLPHIDDEN*(MLPIN+MLPOUT)+1 23: ! bhidden at MLPHIDDEN*(MLPIN+MLPOUT)+1
 25: ! bout    at MLPHIDDEN*(MLPIN+MLPOUT)+2 24: ! bout    at MLPHIDDEN*(MLPIN+MLPOUT)+2
 26: ! 25: !
 27: ! bout literally does nothing - all first and second derivatives are zero. 
 28: !      may as well remove it. 
 29: ! 
 30:  26: 
 31: MLPOFFSET=MLPHIDDEN*MLPIN 27: MLPOFFSET=MLPHIDDEN*MLPIN
 32: BOFFSET=MLPHIDDEN*(MLPIN+MLPOUT) 28: BOFFSET=MLPHIDDEN*(MLPIN+MLPOUT)
 33: ENERGY=0.0D0 29: ENERGY=0.0D0
 34: V(1:NMLP)=0.0D0 30: V(1:NMLP)=0.0D0
 35: IF (SECT) HESS(1:NMLP,1:NMLP)=0.0D0 31: IF (SECT) HESS(1:NMLP,1:NMLP)=0.0D0
 36: DO J1=1,MLPDATA 32: DO J1=1,MLPDATA
 37:    MLPOUTJ1=MLPOUTCOME(J1) 33:    MLPOUTJ1=MLPOUTCOME(J1)
 38:    DO J2=1,MLPHIDDEN 34:    DO J2=1,MLPHIDDEN
 39:       DUMMY1=0.0D0 35:       DUMMY1=0.0D0
 48: !     ELSE 44: !     ELSE
 49:          SECH2(J2)=1.0D0/COSH(DUMMY1)**2  45:          SECH2(J2)=1.0D0/COSH(DUMMY1)**2 
 50: !     ENDIF 46: !     ENDIF
 51:    ENDDO 47:    ENDDO
 52: !  DUMMY3=0.0D0 48: !  DUMMY3=0.0D0
 53:    DMAX=-1.0D100 49:    DMAX=-1.0D100
 54:    DO J4=1,MLPOUT 50:    DO J4=1,MLPOUT
 55:       DUMMY2=0.0D0 51:       DUMMY2=0.0D0
 56:       DO J2=1,MLPHIDDEN 52:       DO J2=1,MLPHIDDEN
 57:          DO J3=1,MLPIN 53:          DO J3=1,MLPIN
 58: !           DYW2G(J4,J2,J3)=(X( MLPOFFSET + (J4-1)*MLPHIDDEN + J2 )+X(BOFFSET+2)) * MLPDAT(J1,J3)*SECH2(J2) 54:             DYW2G(J4,J2,J3)=(X( MLPOFFSET + (J4-1)*MLPHIDDEN + J2 )+X(BOFFSET+2)) * MLPDAT(J1,J3)*SECH2(J2)
 59:             DYW2G(J4,J2,J3)=(X( MLPOFFSET + (J4-1)*MLPHIDDEN + J2 )) * MLPDAT(J1,J3)*SECH2(J2) 
 60:          ENDDO 55:          ENDDO
 61: !        DUMMY2=DUMMY2+(X(MLPOFFSET+(J4-1)*MLPHIDDEN+J2)+X(BOFFSET+2))*TANHSUM(J2) 56:          DUMMY2=DUMMY2+(X(MLPOFFSET+(J4-1)*MLPHIDDEN+J2)+X(BOFFSET+2))*TANHSUM(J2)
 62:          DUMMY2=DUMMY2+X(MLPOFFSET+(J4-1)*MLPHIDDEN+J2)*TANHSUM(J2) 
 63:       ENDDO 57:       ENDDO
 64: !     IF (DUMMY2.GT.50.0D0) DUMMY2=50.0D0 ! to prevent FPE 58: !     IF (DUMMY2.GT.50.0D0) DUMMY2=50.0D0 ! to prevent FPE
 65:       IF (DUMMY2.GT.DMAX) DMAX=DUMMY2 59:       IF (DUMMY2.GT.DMAX) DMAX=DUMMY2
 66:       Y(J4)=DUMMY2 60:       Y(J4)=DUMMY2
 67: !     DUMMY3=DUMMY3+EXP(DUMMY2) 61: !     DUMMY3=DUMMY3+EXP(DUMMY2)
 68:    ENDDO   62:    ENDDO  
 69: ! 63: !
 70: ! Factor DMAX out of the exponentials to prevent overflow. 64: ! Factor DMAX out of the exponentials to prevent overflow.
 71: ! 65: !
 72:    DUMMY3=0.0D0 66:    DUMMY3=0.0D0
 98:       DO J3=1,MLPIN 92:       DO J3=1,MLPIN
 99:          DO J2=1,MLPHIDDEN 93:          DO J2=1,MLPHIDDEN
100:             DUMMY3=0.0D0 94:             DUMMY3=0.0D0
101:             DO J4=1,MLPOUT 95:             DO J4=1,MLPOUT
102:                DUMMY3=DUMMY3+PROB(J4)*DYW2G(J4,J2,J3) 96:                DUMMY3=DUMMY3+PROB(J4)*DYW2G(J4,J2,J3)
103:             ENDDO 97:             ENDDO
104:             DPCW2BG(J2,J3)=PMLPOUTJ1*(DYW2G(MLPOUTJ1,J2,J3)-DUMMY3) 98:             DPCW2BG(J2,J3)=PMLPOUTJ1*(DYW2G(MLPOUTJ1,J2,J3)-DUMMY3)
105:          ENDDO 99:          ENDDO
106:       ENDDO100:       ENDDO
107: 101: 
108: !     V(BOFFSET+2)=0.0D0102:       V(BOFFSET+2)=0.0D0
109: 103: 
110:       SUMINPUTS=0.0D0104:       SUMINPUTS=0.0D0
111:       DO J3=1,MLPIN105:       DO J3=1,MLPIN
112:          SUMINPUTS=SUMINPUTS+MLPDAT(J1,J3)106:          SUMINPUTS=SUMINPUTS+MLPDAT(J1,J3)
113:       ENDDO107:       ENDDO
114: 108: 
115:       DO J4=1,MLPOUT109:       DO J4=1,MLPOUT
116:          DUMMY2=0.0D0110:          DUMMY2=0.0D0
117:          DO J2=1,MLPHIDDEN111:          DO J2=1,MLPHIDDEN
118: !           DUMMY2=DUMMY2+(X(MLPOFFSET+(J4-1)*MLPHIDDEN+J2)+X(BOFFSET+2))*SECH2(J2)112:             DUMMY2=DUMMY2+(X(MLPOFFSET+(J4-1)*MLPHIDDEN+J2)+X(BOFFSET+2))*SECH2(J2)
119:             DUMMY2=DUMMY2+X(MLPOFFSET+(J4-1)*MLPHIDDEN+J2)*SECH2(J2) 
120:          ENDDO113:          ENDDO
121:          DYCDBHID(J4)=DUMMY2*SUMINPUTS114:          DYCDBHID(J4)=DUMMY2*SUMINPUTS
122:       ENDDO115:       ENDDO
123: 116: 
124:       DUMMY2=0.0D0117:       DUMMY2=0.0D0
125:       DO J4=1,MLPOUT118:       DO J4=1,MLPOUT
126:          DUMMY2=DUMMY2+PROB(J4)*DYCDBHID(J4)119:          DUMMY2=DUMMY2+PROB(J4)*DYCDBHID(J4)
127:       ENDDO120:       ENDDO
128:       V(BOFFSET+1)=V(BOFFSET+1)-(DYCDBHID(MLPOUTJ1)-DUMMY2)121:       V(BOFFSET+1)=V(BOFFSET+1)-(DYCDBHID(MLPOUTJ1)-DUMMY2)
129:       DPCDBHIDOUTJ1=PMLPOUTJ1*(DYCDBHID(MLPOUTJ1)-DUMMY2)122:       DPCDBHIDOUTJ1=PMLPOUTJ1*(DYCDBHID(MLPOUTJ1)-DUMMY2)
169:   &         HESS(BOFFSET+1,MLPOFFSET+(J4-1)*MLPHIDDEN+J2) + DPCDBHIDOUTJ1*DPCW1BG(J4,J2)/PMLPOUTJ1**2  &162:   &         HESS(BOFFSET+1,MLPOFFSET+(J4-1)*MLPHIDDEN+J2) + DPCDBHIDOUTJ1*DPCW1BG(J4,J2)/PMLPOUTJ1**2  &
170:   &                                -DPCDBHIDDW1BG(J4,J2)/PMLPOUTJ1163:   &                                -DPCDBHIDDW1BG(J4,J2)/PMLPOUTJ1
171:          ENDDO164:          ENDDO
172:       ENDDO165:       ENDDO
173: 166: 
174: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!167: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
175: 168: 
176:       DO J4=1,MLPOUT ! J4 is input epsilon169:       DO J4=1,MLPOUT ! J4 is input epsilon
177:          DUMMY1=0.0D0170:          DUMMY1=0.0D0
178:          DO J2=1,MLPHIDDEN ! J2 is gamma171:          DO J2=1,MLPHIDDEN ! J2 is gamma
179: !           DUMMY1=DUMMY1+SECH2(J2)*TANHSUM(J2)*(X(MLPHIDDEN*MLPIN + (J4-1)*MLPHIDDEN+J2)+X(BOFFSET+2))172:             DUMMY1=DUMMY1+SECH2(J2)*TANHSUM(J2)*(X(MLPHIDDEN*MLPIN + (J4-1)*MLPHIDDEN+J2)+X(BOFFSET+2))
180:             DUMMY1=DUMMY1+SECH2(J2)*TANHSUM(J2)*X(MLPHIDDEN*MLPIN + (J4-1)*MLPHIDDEN+J2) 
181:          ENDDO173:          ENDDO
182:          D2YCDBHID2(J4)=-2.0D0*SUMINPUTS**2*DUMMY1174:          D2YCDBHID2(J4)=-2.0D0*SUMINPUTS**2*DUMMY1
183:       ENDDO175:       ENDDO
184: 176: 
185:       DUMMY1=0.0D0177:       DUMMY1=0.0D0
186:       DUMMY3=0.0D0178:       DUMMY3=0.0D0
187:       DO J4=1,MLPOUT179:       DO J4=1,MLPOUT
188:          DUMMY1=DUMMY1+PROB(J4)*D2YCDBHID2(J4)+DPCDBHID(J4)*DYCDBHID(J4)180:          DUMMY1=DUMMY1+PROB(J4)*D2YCDBHID2(J4)+DPCDBHID(J4)*DYCDBHID(J4)
189:          DUMMY3=DUMMY3+PROB(J4)*DYCDBHID(J4)181:          DUMMY3=DUMMY3+PROB(J4)*DYCDBHID(J4)
190:       ENDDO182:       ENDDO
199: !191: !
200:       DUMMY1=0.0D0192:       DUMMY1=0.0D0
201:       DO J4=1,MLPOUT193:       DO J4=1,MLPOUT
202:          DUMMY1=DUMMY1+PROB(J4)*DYCDBHID(J4)194:          DUMMY1=DUMMY1+PROB(J4)*DYCDBHID(J4)
203:       ENDDO195:       ENDDO
204: 196: 
205:       DO J2=1,MLPHIDDEN ! J2 is beta197:       DO J2=1,MLPHIDDEN ! J2 is beta
206:          DO J3=1,MLPIN  ! J3 is gamma198:          DO J3=1,MLPIN  ! J3 is gamma
207:             DUMMY4=0.0D0199:             DUMMY4=0.0D0
208:             DO J4=1,MLPOUT200:             DO J4=1,MLPOUT
209: !              DUMMY4=DUMMY4+PROB(J4)*(X(MLPOFFSET+(J4-1)*MLPHIDDEN+J2)+X(BOFFSET+2))201:                DUMMY4=DUMMY4+PROB(J4)*(X(MLPOFFSET+(J4-1)*MLPHIDDEN+J2)+X(BOFFSET+2))
210:                DUMMY4=DUMMY4+PROB(J4)*X(MLPOFFSET+(J4-1)*MLPHIDDEN+J2) 
211:             ENDDO202:             ENDDO
212:             DO J4=1,MLPOUT203:             DO J4=1,MLPOUT
213: !              DPCDW2BG(J4,J2,J3)=PROB(J4)*MLPDAT(J1,J3)*SECH2(J2)*(X(MLPOFFSET+(J4-1)*MLPHIDDEN+J2)+X(BOFFSET+2)-DUMMY4)204:                DPCDW2BG(J4,J2,J3)=PROB(J4)*MLPDAT(J1,J3)*SECH2(J2)*(X(MLPOFFSET+(J4-1)*MLPHIDDEN+J2)+X(BOFFSET+2)-DUMMY4)
214:                DPCDW2BG(J4,J2,J3)=PROB(J4)*MLPDAT(J1,J3)*SECH2(J2)*(X(MLPOFFSET+(J4-1)*MLPHIDDEN+J2)-DUMMY4) 
215:             ENDDO205:             ENDDO
216:             DUMMY2=0.0D0206:             DUMMY2=0.0D0
217:             DUMMY3=0.0D0207:             DUMMY3=0.0D0
218:             DO J4=1,MLPOUT208:             DO J4=1,MLPOUT
219:                DUMMY2=DUMMY2+DYCDBHID(J4)*DPCDW2BG(J4,J2,J3)209:                DUMMY2=DUMMY2+DYCDBHID(J4)*DPCDW2BG(J4,J2,J3)
220: !              DUMMY3=DUMMY3-2.0D0*PROB(J4)*MLPDAT(J1,J3)*SUMINPUTS*(X(MLPOFFSET+(J4-1)*MLPHIDDEN+J2)+X(BOFFSET+2))*SECH2(J2)*TANHSUM(J2)210:                DUMMY3=DUMMY3-2.0D0*PROB(J4)*MLPDAT(J1,J3)*SUMINPUTS*(X(MLPOFFSET+(J4-1)*MLPHIDDEN+J2)+X(BOFFSET+2))*SECH2(J2)*TANHSUM(J2)
221:                DUMMY3=DUMMY3-2.0D0*PROB(J4)*MLPDAT(J1,J3)*SUMINPUTS*X(MLPOFFSET+(J4-1)*MLPHIDDEN+J2)*SECH2(J2)*TANHSUM(J2) 
222:             ENDDO211:             ENDDO
223: 212: 
224: !           D2YCDBHIDDW2BG=-2.0D0*MLPDAT(J1,J3)*SUMINPUTS*(X(MLPOFFSET+(MLPOUTJ1-1)*MLPHIDDEN+J2)+X(BOFFSET+2))*SECH2(J2)*TANHSUM(J2)213:             D2YCDBHIDDW2BG=-2.0D0*MLPDAT(J1,J3)*SUMINPUTS*(X(MLPOFFSET+(MLPOUTJ1-1)*MLPHIDDEN+J2)+X(BOFFSET+2))*SECH2(J2)*TANHSUM(J2)
225:             D2YCDBHIDDW2BG=-2.0D0*MLPDAT(J1,J3)*SUMINPUTS*X(MLPOFFSET+(MLPOUTJ1-1)*MLPHIDDEN+J2)*SECH2(J2)*TANHSUM(J2) 
226: 214: 
227:             D2PCDBHIDDW2BG=DPCDW2BG(MLPOUTJ1,J2,J3)*(DYCDBHID(MLPOUTJ1)-DUMMY1)+PMLPOUTJ1*(D2YCDBHIDDW2BG-DUMMY2-DUMMY3)215:             D2PCDBHIDDW2BG=DPCDW2BG(MLPOUTJ1,J2,J3)*(DYCDBHID(MLPOUTJ1)-DUMMY1)+PMLPOUTJ1*(D2YCDBHIDDW2BG-DUMMY2-DUMMY3)
228: 216: 
229:             HESS(BOFFSET+1,(J2-1)*MLPIN+J3)= &217:             HESS(BOFFSET+1,(J2-1)*MLPIN+J3)= &
230:   &         HESS(BOFFSET+1,(J2-1)*MLPIN+J3) + DPCDBHID(MLPOUTJ1)*DPCDW2BG(MLPOUTJ1,J2,J3)/PMLPOUTJ1**2 - D2PCDBHIDDW2BG/PMLPOUTJ1218:   &         HESS(BOFFSET+1,(J2-1)*MLPIN+J3) + DPCDBHID(MLPOUTJ1)*DPCDW2BG(MLPOUTJ1,J2,J3)/PMLPOUTJ1**2 - D2PCDBHIDDW2BG/PMLPOUTJ1
231:          ENDDO219:          ENDDO
232:       ENDDO220:       ENDDO
233: !221: !
234: ! This block w^1 with w^1 is locally symmetric222: ! This block w^1 with w^1 is locally symmetric
235: !223: !
311:   &               HESS((J2-1)*MLPIN+J3,(K2-1)*MLPIN+K4) & ! sum over data points299:   &               HESS((J2-1)*MLPIN+J3,(K2-1)*MLPIN+K4) & ! sum over data points
312:   &               +DPCW2BG(J2,J3)*DPCW2BG(K2,K4)/PMLPOUTJ1**2 - DUMMY1/PMLPOUTJ1300:   &               +DPCW2BG(J2,J3)*DPCW2BG(K2,K4)/PMLPOUTJ1**2 - DUMMY1/PMLPOUTJ1
313:                ENDDO301:                ENDDO
314:             ENDDO302:             ENDDO
315:          ENDDO303:          ENDDO
316:       ENDDO304:       ENDDO
317:    ENDIF305:    ENDIF
318: ENDDO306: ENDDO
319: 307: 
320: DUMMY1=0.0D0308: DUMMY1=0.0D0
321: DO J1=1,NMLP-1309: DO J1=1,NMLP-2
322:    DUMMY1=DUMMY1+X(J1)**2310:    DUMMY1=DUMMY1+X(J1)**2
323: ENDDO311: ENDDO
324: 312: 
325: ENERGY=ENERGY/MLPDATA + MLPLAMBDA*DUMMY1313: ENERGY=ENERGY/MLPDATA + MLPLAMBDA*DUMMY1
326: ! IF (DEBUG) WRITE(*,'(A,G20.10)') 'MLP3> objective function=',ENERGY314: ! IF (DEBUG) WRITE(*,'(A,G20.10)') 'MLP3> objective function=',ENERGY
327: 315: 
328: IF (GTEST) V(1:NMLP)=V(1:NMLP)/MLPDATA 316: IF (GTEST) V(1:NMLP)=V(1:NMLP)/MLPDATA 
329: IF (GTEST) V(1:NMLP-1)=V(1:NMLP-1) + 2.0D0*MLPLAMBDA*X(1:NMLP-1)317: IF (GTEST) V(1:NMLP-2)=V(1:NMLP-2) + 2.0D0*MLPLAMBDA*X(1:NMLP-2)
330: !318: !
331: ! Symmetrise Hessian here319: ! Symmetrise Hessian here
332: !320: !
333: IF (SECT) HESS(1:NMLP,1:NMLP)=HESS(1:NMLP,1:NMLP)/MLPDATA321: IF (SECT) HESS(1:NMLP,1:NMLP)=HESS(1:NMLP,1:NMLP)/MLPDATA
334: IF (SECT) THEN322: IF (SECT) THEN
335:    DO J1=1,NMLP-1323:    DO J1=1,NMLP-2
336:       HESS(J1,J1)=HESS(J1,J1)+2*MLPLAMBDA324:       HESS(J1,J1)=HESS(J1,J1)+2*MLPLAMBDA
337:    ENDDO325:    ENDDO
338:    DO J1=1,NMLP326:    DO J1=1,NMLP
339:       DO J2=1,J1-1327:       DO J2=1,J1-1
340:          HESS(J2,J1)=HESS(J1,J2)328:          HESS(J2,J1)=HESS(J1,J2)
341:       ENDDO329:       ENDDO
342:    ENDDO330:    ENDDO
343: ENDIF331: ENDIF
344: 332: 
345: END SUBROUTINE MLPB3333: END SUBROUTINE MLPB3


r29875/potential.f 2016-01-31 09:30:10.982291173 +0000 r29874/potential.f 2016-01-31 09:30:11.770301587 +0000
478:             478:             
479:             IF (PTEST) THEN479:             IF (PTEST) THEN
480:                WRITE(*,10) ' potential> Energy for last cycle=',ENERGY,' '480:                WRITE(*,10) ' potential> Energy for last cycle=',ENERGY,' '
481:                WRITE(ESTRING,10) 'Energy for last cycle=',ENERGY,' '481:                WRITE(ESTRING,10) 'Energy for last cycle=',ENERGY,' '
482:             ENDIF482:             ENDIF
483: 483: 
484: 484: 
485:          ELSE IF (VARIABLES) THEN485:          ELSE IF (VARIABLES) THEN
486: 486: 
487:             IF (MLPB3T) THEN487:             IF (MLPB3T) THEN
488:                CALL MLPB3(COORDS, VNEW, ENERGY, GTEST, STEST)488: !              CALL MLPB3(COORDS, VNEW, ENERGY, GTEST, STEST)
489: !               DIFF=1.0D-4489:                DIFF=1.0D-4
490: !               PRINT*,'analytic and numerical gradients:'490:                PRINT*,'analytic and numerical gradients:'
491: !               IF (.NOT.(ALLOCATED(HESS))) ALLOCATE(HESS(NMLP,NMLP))491:                IF (.NOT.(ALLOCATED(HESS))) ALLOCATE(HESS(NMLP,NMLP))
492: !               CALL MLPB3(COORDS, VNEW, ENERGY, .TRUE., .TRUE.)492:                CALL MLPB3(COORDS, VNEW, ENERGY, .TRUE., .TRUE.)
493: !               PRINT '(A,I8)','SIZE(HESS)=',SIZE(HESS)493:                PRINT '(A,I8)','SIZE(HESS)=',SIZE(HESS)
494: !               HESSDUM(1:NMLP,1:NMLP)=HESS(1:NMLP,1:NMLP)494:                HESSDUM(1:NMLP,1:NMLP)=HESS(1:NMLP,1:NMLP)
495: !               DO J1=1,NATOMS495:                DO J1=1,NATOMS
496: !                  COORDS(J1)=COORDS(J1)+DIFF496:                   COORDS(J1)=COORDS(J1)+DIFF
497: !                  CALL MLPB3(COORDS,VPLUS,EPLUS,.FALSE.,.FALSE.)497:                   CALL MLPB3(COORDS,VPLUS,EPLUS,.FALSE.,.FALSE.)
498: !                  COORDS(J1)=COORDS(J1)-2.0D0*DIFF498:                   COORDS(J1)=COORDS(J1)-2.0D0*DIFF
499: !                  CALL MLPB3(COORDS,VMINUS,EMINUS,.FALSE.,.FALSE.)499:                   CALL MLPB3(COORDS,VMINUS,EMINUS,.FALSE.,.FALSE.)
500: !                  COORDS(J1)=COORDS(J1)+DIFF500:                   COORDS(J1)=COORDS(J1)+DIFF
501: !                  WRITE(*,'(I5,2F20.10)') J1,VNEW(J1),(EPLUS-EMINUS)/(2.0D0*DIFF)501:                   WRITE(*,'(I5,2F20.10)') J1,VNEW(J1),(EPLUS-EMINUS)/(2.0D0*DIFF)
502: !               ENDDO502:                ENDDO
503: !               PRINT*,'analytic and numerical second derivatives:'503:                PRINT*,'analytic and numerical second derivatives:'
504: !               DO J1=1,NATOMS504:                DO J1=1,NATOMS
505: !                  COORDS(J1)=COORDS(J1)+DIFF505:                   COORDS(J1)=COORDS(J1)+DIFF
506: !                  CALL MLPB3(COORDS,VPLUS,EPLUS,.TRUE.,.FALSE.)506:                   CALL MLPB3(COORDS,VPLUS,EPLUS,.TRUE.,.FALSE.)
507: !                  COORDS(J1)=COORDS(J1)-2.0D0*DIFF507:                   COORDS(J1)=COORDS(J1)-2.0D0*DIFF
508: !                  CALL MLPB3(COORDS,VMINUS,EMINUS,.TRUE.,.FALSE.)508:                   CALL MLPB3(COORDS,VMINUS,EMINUS,.TRUE.,.FALSE.)
509: !                  COORDS(J1)=COORDS(J1)+DIFF509:                   COORDS(J1)=COORDS(J1)+DIFF
510: !                  DO J2=1,NATOMS510:                   DO J2=1,NATOMS
511: !                     IF ((ABS(HESS(J1,J2)).NE.0.0D0).AND. 511:                      IF ((ABS(HESS(J1,J2)).NE.0.0D0).AND. 
512: !     &                   (ABS(100.0D0*(HESS(J1,J2)-(VPLUS(J2)-VMINUS(J2))/(2.0D0*DIFF))/HESS(J1,J2)).GT.1.0D-3)) THEN512:      &                   (ABS(100.0D0*(HESS(J1,J2)-(VPLUS(J2)-VMINUS(J2))/(2.0D0*DIFF))/HESS(J1,J2)).GT.1.0D-3)) THEN
513: !                     WRITE(*,'(2I5,2F20.10,A)') J1,J2,HESS(J1,J2),(VPLUS(J2)-VMINUS(J2))/(2.0D0*DIFF),'   X'513:                      WRITE(*,'(2I5,2F20.10,A)') J1,J2,HESS(J1,J2),(VPLUS(J2)-VMINUS(J2))/(2.0D0*DIFF),'   X'
514: !                     ELSE514:                      ELSE
515: !                        WRITE(*,'(2I5,2F20.10,A)') J1,J2,HESS(J1,J2),(VPLUS(J2)-VMINUS(J2))/(2.0D0*DIFF)515:                         WRITE(*,'(2I5,2F20.10,A)') J1,J2,HESS(J1,J2),(VPLUS(J2)-VMINUS(J2))/(2.0D0*DIFF)
516: !                     ENDIF516:                      ENDIF
517: !                  ENDDO517:                   ENDDO
518: !               ENDDO518:                ENDDO
519: !               STOP519:                STOP
520:             ELSEIF (MLP3T) THEN520:             ELSEIF (MLP3T) THEN
521:                CALL MLP3(COORDS, VNEW, ENERGY, GTEST, STEST)521:                CALL MLP3(COORDS, VNEW, ENERGY, GTEST, STEST)
522:             ELSE522:             ELSE
523:                CALL FUNCTIONAL( COORDS, VNEW, ENERGY, GTEST, STEST)523:                CALL FUNCTIONAL( COORDS, VNEW, ENERGY, GTEST, STEST)
524:             ENDIF524:             ENDIF
525:             IF (PTEST) THEN525:             IF (PTEST) THEN
526:                 WRITE(*,10) ' potential> Energy for last cycle=',ENERGY,'         '526:                 WRITE(*,10) ' potential> Energy for last cycle=',ENERGY,'         '
527:                 WRITE(ESTRING,10) 'Energy for last cycle=',ENERGY,'         '527:                 WRITE(ESTRING,10) 'Energy for last cycle=',ENERGY,'         '
528:             ENDIF528:             ENDIF
529:             ! CALL CTEST(NATOMS, COORDS, VNEW, ENERGY, GTEST, STEST)529:             ! CALL CTEST(NATOMS, COORDS, VNEW, ENERGY, GTEST, STEST)


legend
Lines Added 
Lines changed
 Lines Removed

hdiff - version: 2.1.0