hdiff output

r33036/commons.f90 2017-07-19 13:30:12.269781491 +0100 r33035/commons.f90 2017-07-19 13:30:13.473797563 +0100
109:      &        LJSITECOORDST, VGW, ACKLANDT, G46, DF1T, PULLT, LOCALSAMPLET, CSMT, A9INTET, INTERESTORE, COLDFUSION, &109:      &        LJSITECOORDST, VGW, ACKLANDT, G46, DF1T, PULLT, LOCALSAMPLET, CSMT, A9INTET, INTERESTORE, COLDFUSION, &
110:      &        CSMGUIDET, MULTISITEPYT, CHAPERONINT, AVOIDRESEEDT, OHCELLT, UNFREEZEFINALQ, PERCOLATET, PERCT, PERCACCEPTED, PERCCOMPMARKOV, PERCGROUPT, &110:      &        CSMGUIDET, MULTISITEPYT, CHAPERONINT, AVOIDRESEEDT, OHCELLT, UNFREEZEFINALQ, PERCOLATET, PERCT, PERCACCEPTED, PERCCOMPMARKOV, PERCGROUPT, &
111:      &        GENALT, MINDENSITYT, RESTRICTREGION, RESTRICTREGIONTEST, RESTRICTCYL, ACK1, ACK2, HARMONICF, PERCGROUPRESEEDT, &111:      &        GENALT, MINDENSITYT, RESTRICTREGION, RESTRICTREGIONTEST, RESTRICTCYL, ACK1, ACK2, HARMONICF, PERCGROUPRESEEDT, &
112:      &        HARMONICDONTMOVE, DUMPUNIQUE, FREEZESAVE, TBP, RBSYMT, PTMCDUMPSTRUCT, PTMCDUMPENERT, PYCOLDFUSION, MONITORT,&112:      &        HARMONICDONTMOVE, DUMPUNIQUE, FREEZESAVE, TBP, RBSYMT, PTMCDUMPSTRUCT, PTMCDUMPENERT, PYCOLDFUSION, MONITORT,&
113:      &        CHARMMDFTBT, PERMINVOPT, BLOCKMOVET, MAXERISE_SET, PYT, BINARY_EXAB, CHIROT, POLYT, SANDBOXT, &113:      &        CHARMMDFTBT, PERMINVOPT, BLOCKMOVET, MAXERISE_SET, PYT, BINARY_EXAB, CHIROT, POLYT, SANDBOXT, &
114:      &        RESERVOIRT, DISTOPT, ONEDAPBCT, ONEDPBCT, TWODAPBCT, TWODPBCT, THREEDAPBCT, THREEDPBCT, RATIOT, &114:      &        RESERVOIRT, DISTOPT, ONEDAPBCT, ONEDPBCT, TWODAPBCT, TWODPBCT, THREEDAPBCT, THREEDPBCT, RATIOT, &
115:      &        PTRANDOM, PTINTERVAL, PTSINGLE, PTSETS, CHEMSHIFT, CHEMSHIFT2, CSH, DEBUGss2029, UNIFORMMOVE, RANSEEDT, &115:      &        PTRANDOM, PTINTERVAL, PTSINGLE, PTSETS, CHEMSHIFT, CHEMSHIFT2, CSH, DEBUGss2029, UNIFORMMOVE, RANSEEDT, &
116:      &        TTM3T, NOINVERSION, RIGIDCONTOURT, UPDATERIGIDREFT, HYBRIDMINT, COMPRESSRIGIDT, MWFILMT, &116:      &        TTM3T, NOINVERSION, RIGIDCONTOURT, UPDATERIGIDREFT, HYBRIDMINT, COMPRESSRIGIDT, MWFILMT, &
117:      &        SUPPRESST, MFETT, POLIRT, QUIPT, SWPOTT, MWPOTT, REPMATCHT, GLJT, MLJT, LJ_VAR_RADT, READMASST, SPECMASST, NEWTSALLIST, &117:      &        SUPPRESST, MFETT, POLIRT, QUIPT, SWPOTT, MWPOTT, REPMATCHT, GLJT, MLJT, LJ_VAR_RADT, READMASST, SPECMASST, NEWTSALLIST, &
118:      &        PHI4MODELT, CUDAT, CUDATIMET, AMBER12T, ENERGY_DECOMPT, NEWMOVEST, DUMPMINT, MBPOLT, MOLECULART, GCBHT, SEMIGRAND_MUT, USEROT, &118:      &        PHI4MODELT, CUDAT, CUDATIMET, AMBER12T, ENERGY_DECOMPT, NEWMOVEST, DUMPMINT, MBPOLT, MOLECULART, GCBHT, SEMIGRAND_MUT, USEROT, &
119:      &        SAVEMULTIMINONLY, GRADPROBLEMT, INTLJT, CONDATT, QCIPERMCHECK, RIGIDMBPOLT, &119:      &        SAVEMULTIMINONLY, GRADPROBLEMT, INTLJT, CONDATT, QCIPERMCHECK, &
120:      &        INTCONSTRAINTT, INTFREEZET, CHECKCONINT, CONCUTABST, CONCUTFRACT, INTERPCOSTFUNCTION, &120:      &        INTCONSTRAINTT, INTFREEZET, CHECKCONINT, CONCUTABST, CONCUTFRACT, INTERPCOSTFUNCTION, &
121:      &        RBAAT, FREEZENODEST, DUMPINTEOS, DUMPINTXYZ, QCIPOTT, QCIPOT2T, INTSPRINGACTIVET, LPERMDIST, LOCALPERMDIST, QCIRADSHIFTT, &121:      &        RBAAT, FREEZENODEST, DUMPINTEOS, DUMPINTXYZ, QCIPOTT, QCIPOT2T, INTSPRINGACTIVET, LPERMDIST, LOCALPERMDIST, QCIRADSHIFTT, &
122:      &        MLP3T, MKTRAPT, MLPB3T, MLPB3NEWT, MULTIPOTT, QCIAMBERT, MLPNEWREG, DJWRBT, STEALTHYT, LJADDT, QCINOREPINT, RIGIDMDT, &122:      &        MLP3T, MKTRAPT, MLPB3T, MLPB3NEWT, MULTIPOTT, QCIAMBERT, MLPNEWREG, DJWRBT, STEALTHYT, LJADDT, QCINOREPINT, RIGIDMDT, &
123:      &        DUMPMQT, MLQT, MLQPROB, LJADD2T, MLPVB3T, NOREGBIAS, PYADDT, PYADD2T, LJADD3T, REORDERADDT,  LJADD4T, &123:      &        DUMPMQT, MLQT, MLQPROB, LJADD2T, MLPVB3T, NOREGBIAS, PYADDT, PYADD2T, LJADD3T, REORDERADDT,  LJADD4T, &
124:      &        SQNMT, SQNM_DEBUGT, SQNM_BIOT, BENZRIGIDEWALDT, ORTHO, EWALDT, WATERMETHANET, CLATHRATET124:      &        SQNMT, SQNM_DEBUGT, SQNM_BIOT, BENZRIGIDEWALDT, ORTHO, EWALDT, WATERMETHANET
125: !125: !
126:       DOUBLE PRECISION, ALLOCATABLE :: SEMIGRAND_MU(:)126:       DOUBLE PRECISION, ALLOCATABLE :: SEMIGRAND_MU(:)
127:       DOUBLE PRECISION, ALLOCATABLE :: ATMASS(:)127:       DOUBLE PRECISION, ALLOCATABLE :: ATMASS(:)
128:       DOUBLE PRECISION, ALLOCATABLE :: SPECMASS(:)128:       DOUBLE PRECISION, ALLOCATABLE :: SPECMASS(:)
129: 129: 
130: ! dj337: Ewald summation variables130: ! dj337: Ewald summation variables
131:       DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:,:) :: RERHOARRAY, IMRHOARRAY131:       DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:,:) :: RERHOARRAY, IMRHOARRAY
132: 132: 
133: ! csw34> FREEZEGROUP variables133: ! csw34> FREEZEGROUP variables
134: !134: !
135:       INTEGER :: GROUPCENTRE135:       INTEGER :: GROUPCENTRE
136:       DOUBLE PRECISION :: GROUPRADIUS136:       DOUBLE PRECISION :: GROUPRADIUS
137:       CHARACTER (LEN=2) :: FREEZEGROUPTYPE137:       CHARACTER (LEN=2) :: FREEZEGROUPTYPE
138:       LOGICAL :: FREEZEGROUPT138:       LOGICAL :: FREEZEGROUPT
139: ! END139: ! END
140: 140: 
141: !141: !
142: ! csw34> DONTMOVE variables142: ! csw34> DONTMOVE variables
143: 143: !
144:       INTEGER :: NDONTMOVE, DONTMOVECENTRE144:       INTEGER :: NDONTMOVE, DONTMOVECENTRE
145:       DOUBLE PRECISION :: DONTMOVERADIUS145:       DOUBLE PRECISION :: DONTMOVERADIUS
146:       CHARACTER (LEN=2) :: DONTMOVEGROUPTYPE146:       CHARACTER (LEN=2) :: DONTMOVEGROUPTYPE
147:       LOGICAL :: DONTMOVET, DONTMOVEGROUPT, DONTMOVEREST, DONTMOVEALL, DOMOVEREST147:       LOGICAL :: DONTMOVET, DONTMOVEGROUPT, DONTMOVEREST, DONTMOVEALL, DOMOVEREST
148:       LOGICAL, ALLOCATABLE :: DONTMOVE(:),DONTMOVERES(:)148:       LOGICAL, ALLOCATABLE :: DONTMOVE(:),DONTMOVERES(:)
149:       INTEGER, ALLOCATABLE :: DUMPXYZUNIT(:), DUMPVUNIT(:)149:       INTEGER, ALLOCATABLE :: DUMPXYZUNIT(:), DUMPVUNIT(:)
150: !150: !
151: ! csw34> PAIRDIST variables151: ! csw34> PAIRDIST variables
152: !152: !
153:       INTEGER :: NPAIRS153:       INTEGER :: NPAIRS
357: !ds656> Stress tensor357: !ds656> Stress tensor
358:       LOGICAL :: STRESST358:       LOGICAL :: STRESST
359:       INTEGER :: STRESS_MODE359:       INTEGER :: STRESS_MODE
360:       DOUBLE PRECISION, ALLOCATABLE :: STRESS(:,:,:)360:       DOUBLE PRECISION, ALLOCATABLE :: STRESS(:,:,:)
361: 361: 
362: !ds656> A saw-tooth temperature protocol362: !ds656> A saw-tooth temperature protocol
363:       LOGICAL :: SAWTOOTH363:       LOGICAL :: SAWTOOTH
364:       INTEGER :: SAWTOOTH_NREJMAX364:       INTEGER :: SAWTOOTH_NREJMAX
365:       DOUBLE PRECISION :: SAWTOOTH_TMAX, SAWTOOTH_TFAC, &365:       DOUBLE PRECISION :: SAWTOOTH_TMAX, SAWTOOTH_TFAC, &
366:            SAWTOOTH_SFAC, SAWTOOTH_SFAC2366:            SAWTOOTH_SFAC, SAWTOOTH_SFAC2
367: !cv320> Variable for clathrates367: 
368:       INTEGER :: NWATER 
369: !ds656> Dump current Markov state at regular intervals368: !ds656> Dump current Markov state at regular intervals
370:       LOGICAL :: DUMP_MARKOV369:       LOGICAL :: DUMP_MARKOV
371:       INTEGER :: DUMP_MARKOV_NWAIT, DUMP_MARKOV_NFREQ370:       INTEGER :: DUMP_MARKOV_NWAIT, DUMP_MARKOV_NFREQ
372: 371: 
373: !   allocatables372: !   allocatables
374: 373: 
375:       INTEGER,ALLOCATABLE,DIMENSION(:) :: MOVABLEATOMLIST         ! a list containing the movable atom indices374:       INTEGER,ALLOCATABLE,DIMENSION(:) :: MOVABLEATOMLIST         ! a list containing the movable atom indices
376:       LOGICAL,ALLOCATABLE,DIMENSION(:) :: MOVABLEATOMLISTLOGICAL  ! is atom i movable?375:       LOGICAL,ALLOCATABLE,DIMENSION(:) :: MOVABLEATOMLISTLOGICAL  ! is atom i movable?
377:       INTEGER,ALLOCATABLE,DIMENSION(:) :: ATOMSINBLOCK            ! for BLOCKMOVE, to split movableatoms into separate blocks376:       INTEGER,ALLOCATABLE,DIMENSION(:) :: ATOMSINBLOCK            ! for BLOCKMOVE, to split movableatoms into separate blocks
378:       INTEGER,ALLOCATABLE,DIMENSION(:) :: NSPECIES(:), NSPECIES_INI(:)             ! for multicomponent systems377:       INTEGER,ALLOCATABLE,DIMENSION(:) :: NSPECIES(:), NSPECIES_INI(:)             ! for multicomponent systems


r33036/genrigid.f90 2017-07-19 13:30:12.489784426 +0100 r33035/genrigid.f90 2017-07-19 13:30:13.693800503 +0100
833:       NLATTICECOORDS=6833:       NLATTICECOORDS=6
834:   ENDIF834:   ENDIF
835: 835: 
836:   GTEST = .TRUE.836:   GTEST = .TRUE.
837:   GR(:) = 0.0D0837:   GR(:) = 0.0D0
838:   838:   
839:   DO J1 = 1, NRIGIDBODY839:   DO J1 = 1, NRIGIDBODY
840:      840:      
841:      PI = XR(3*NRIGIDBODY+3*J1-2 : 3*NRIGIDBODY+3*J1)841:      PI = XR(3*NRIGIDBODY+3*J1-2 : 3*NRIGIDBODY+3*J1)
842:      CALL RMDRVT(PI, RMI, DRMI1, DRMI2, DRMI3, GTEST)842:      CALL RMDRVT(PI, RMI, DRMI1, DRMI2, DRMI3, GTEST)
843: 843:      
844:      DO J2 = 1, NSITEPERBODY(J1)844:      DO J2 = 1, NSITEPERBODY(J1)
845:         J9 = RIGIDGROUPS(J2, J1)845:         J9 = RIGIDGROUPS(J2, J1)
846: 846: 
847: ! hk286 > translation847: ! hk286 > translation
848:         GR(3*J1-2:3*J1) = GR(3*J1-2:3*J1) + G(3*J9-2:3*J9)848:         GR(3*J1-2:3*J1) = GR(3*J1-2:3*J1) + G(3*J9-2:3*J9)
849:         849:         
850: ! hk286 > rotation850: ! hk286 > rotation
851:         DR1(:) = MATMUL(DRMI1,SITESRIGIDBODY(J2,:,J1))851:         DR1(:) = MATMUL(DRMI1,SITESRIGIDBODY(J2,:,J1))
852:         DR2(:) = MATMUL(DRMI2,SITESRIGIDBODY(J2,:,J1))852:         DR2(:) = MATMUL(DRMI2,SITESRIGIDBODY(J2,:,J1))
853:         DR3(:) = MATMUL(DRMI3,SITESRIGIDBODY(J2,:,J1))853:         DR3(:) = MATMUL(DRMI3,SITESRIGIDBODY(J2,:,J1))


r33036/keywords.f 2017-07-19 13:30:12.813788767 +0100 r33035/keywords.f 2017-07-19 13:30:13.921803543 +0100
1110:       MWPOTT=.FALSE.1110:       MWPOTT=.FALSE.
1111:       SWPOTT=.FALSE.1111:       SWPOTT=.FALSE.
1112: ! jdf43> SUPPRESS1112: ! jdf43> SUPPRESS
1113:       SUPPRESST=.FALSE.1113:       SUPPRESST=.FALSE.
1114: ! jdf43> MFET1114: ! jdf43> MFET
1115:       MFETT=.FALSE.1115:       MFETT=.FALSE.
1116: ! jdf43> POLIR1116: ! jdf43> POLIR
1117:       POLIRT=.FALSE.1117:       POLIRT=.FALSE.
1118: ! jdf43> MBPOL1118: ! jdf43> MBPOL
1119:       MBPOLT=.FALSE.1119:       MBPOLT=.FALSE.
1120:       RIGIDMBPOLT=.FALSE. 
1121:       MOLECULART=.FALSE.1120:       MOLECULART=.FALSE.
1122: ! jdf43>1121: ! jdf43>
1123:       REPMATCHT=.FALSE.1122:       REPMATCHT=.FALSE.
1124: 1123: 
1125:       UNIFORMMOVE=.FALSE.1124:       UNIFORMMOVE=.FALSE.
1126:       ORBITTOL=1.0D-31125:       ORBITTOL=1.0D-3
1127:       NOINVERSION=.FALSE.1126:       NOINVERSION=.FALSE.
1128: !cv320> rigid body watermethane1127: !cv320> rigid body watermethane
1129:       WATERMETHANET=.FALSE.1128:       WATERMETHANET=.FALSE.
1130:       CLATHRATET= .FALSE. 
1131: ! ds656> Parallelised generalised basin-hopping1129: ! ds656> Parallelised generalised basin-hopping
1132:       GBHT=.FALSE.1130:       GBHT=.FALSE.
1133: 1131: 
1134: ! General mixed LJ systems1132: ! General mixed LJ systems
1135:       GLJT=.FALSE.1133:       GLJT=.FALSE.
1136: ! ds656> Multicomponent LJ system (different implementation to GLJ!)1134: ! ds656> Multicomponent LJ system (different implementation to GLJ!)
1137:       MLJT=.FALSE.1135:       MLJT=.FALSE.
1138: ! jwrm2> LJ with each atom having a different radius. Works with RIGIDINIT.1136: ! jwrm2> LJ with each atom having a different radius. Works with RIGIDINIT.
1139:       LJ_VAR_RADT = .FALSE.1137:       LJ_VAR_RADT = .FALSE.
1140: 1138: 
3676:          CALL DEFPAHARIGID()3674:          CALL DEFPAHARIGID()
3677:          NCARBON = 63675:          NCARBON = 6
3678:          CALL DEFBENZENERIGID()3676:          CALL DEFBENZENERIGID()
3679: 3677: 
3680: ! ----------------------------------------------------------------------------------------3678: ! ----------------------------------------------------------------------------------------
3681: 3679: 
3682:       ELSE IF (WORD.EQ.'FAL') THEN3680:       ELSE IF (WORD.EQ.'FAL') THEN
3683:          FAL=.TRUE.3681:          FAL=.TRUE.
3684:       ELSE IF (WORD.EQ.'WATERMETHANE') THEN3682:       ELSE IF (WORD.EQ.'WATERMETHANE') THEN
3685:          WATERMETHANET=.TRUE.3683:          WATERMETHANET=.TRUE.
3686:       ELSE IF (WORD.EQ.'CLATHRATE') THEN 
3687:          CLATHRATET=.TRUE. 
3688:          CALL READI(NWATER) 
3689:          CALL init_ttm3f(NWATER) 
3690:       ELSE IF (WORD == 'FEBH') THEN3684:       ELSE IF (WORD == 'FEBH') THEN
3691:          CALL READF(FETEMP)3685:          CALL READF(FETEMP)
3692:          FEBHT = .TRUE.3686:          FEBHT = .TRUE.
3693:          FE_FILE_UNIT = GETUNIT()3687:          FE_FILE_UNIT = GETUNIT()
3694:          OPEN(UNIT = FE_FILE_UNIT, FILE = 'free_energy', STATUS = 'REPLACE')3688:          OPEN(UNIT = FE_FILE_UNIT, FILE = 'free_energy', STATUS = 'REPLACE')
3695:          WRITE(FE_FILE_UNIT, '(6A20)') '       Quench       ', '  Potential energy  ',3689:          WRITE(FE_FILE_UNIT, '(6A20)') '       Quench       ', '  Potential energy  ',
3696:      &   '   Harmonic term    ', '    Free energy     ', '   Markov energy    ', '        Time        '3690:      &   '   Harmonic term    ', '    Free energy     ', '   Markov energy    ', '        Time        '
3697:          WRITE(FE_FILE_UNIT, '(6A20)') ' ------------------ ', ' ------------------ ',3691:          WRITE(FE_FILE_UNIT, '(6A20)') ' ------------------ ', ' ------------------ ',
3698:      &   ' ------------------ ', ' ------------------ ', ' ------------------ ', ' ------------------ '3692:      &   ' ------------------ ', ' ------------------ ', ' ------------------ ', ' ------------------ '
3699:          IF (NITEMS .GT. 2) THEN3693:          IF (NITEMS .GT. 2) THEN
7658:       ELSE IF (WORD.EQ.'NEWTSALLIS') THEN7652:       ELSE IF (WORD.EQ.'NEWTSALLIS') THEN
7659:          NEWTSALLIST=.TRUE.7653:          NEWTSALLIST=.TRUE.
7660:          IF (NITEMS.GT.1) THEN7654:          IF (NITEMS.GT.1) THEN
7661:             CALL READF(QTSALLIS)7655:             CALL READF(QTSALLIS)
7662:          ENDIF7656:          ENDIF
7663: 7657: 
7664: !7658: !
7665: !  Xantheas' TTM3-F water potential7659: !  Xantheas' TTM3-F water potential
7666: !7660: !
7667:       ELSE IF (WORD.EQ.'TTM3') THEN7661:       ELSE IF (WORD.EQ.'TTM3') THEN
7668:          TTM3T=.TRUE.7662:         TTM3T=.TRUE.
7669:          CALL READI(NWATER)7663: 
7670:          CALL init_ttm3f(NWATER) 
7671: !! ns644> Adding dihedral angle potential via TWIST keyword:7664: !! ns644> Adding dihedral angle potential via TWIST keyword:
7672:       ELSE IF (WORD.EQ.'TWIST') THEN7665:       ELSE IF (WORD.EQ.'TWIST') THEN
7673:          TWISTT=.TRUE.7666:          TWISTT=.TRUE.
7674:          INQUIRE(FILE='twistgroups',EXIST=YESNO)7667:          INQUIRE(FILE='twistgroups',EXIST=YESNO)
7675:          IF (NITEMS.EQ.7) THEN7668:          IF (NITEMS.EQ.7) THEN
7676:             NTWISTGROUPS = 17669:             NTWISTGROUPS = 1
7677:             ALLOCATE(TWIST_K(1))7670:             ALLOCATE(TWIST_K(1))
7678:             ALLOCATE(TWIST_THETA0(1))7671:             ALLOCATE(TWIST_THETA0(1))
7679:             ALLOCATE(TWIST_ATOMS(1,1))7672:             ALLOCATE(TWIST_ATOMS(1,1))
7680:             CALL READI(TWIST_ATOMS(1,1))7673:             CALL READI(TWIST_ATOMS(1,1))
7917:          POLIRT=.TRUE.7910:          POLIRT=.TRUE.
7918: 7911: 
7919:       ELSE IF (WORD.EQ.'MBPOL') THEN7912:       ELSE IF (WORD.EQ.'MBPOL') THEN
7920:          CALL MBPOLINIT7913:          CALL MBPOLINIT
7921:          MBPOLT=.TRUE.7914:          MBPOLT=.TRUE.
7922:          MOLECULART=.TRUE.7915:          MOLECULART=.TRUE.
7923:          ALLOCATE(TYPECH(NATOMSALLOC))7916:          ALLOCATE(TYPECH(NATOMSALLOC))
7924:          TYPECH(1::3)(1:1)="O"7917:          TYPECH(1::3)(1:1)="O"
7925:          TYPECH(2::3)(1:1)="H"7918:          TYPECH(2::3)(1:1)="H"
7926:          TYPECH(3::3)(1:1)="H"7919:          TYPECH(3::3)(1:1)="H"
7927:       ELSE IF (WORD.EQ.'RIGIDMBPOL') THEN 
7928:          CALL MBPOLINIT 
7929:          RIGIDMBPOLT=.TRUE. 
7930:          RIGID= .TRUE. 
7931:          MOLECULART=.TRUE. 
7932:          ALLOCATE(TYPECH(NATOMSALLOC)) 
7933:          TYPECH(1::3)(1:1)="O" 
7934:          TYPECH(2::3)(1:1)="H" 
7935:          TYPECH(3::3)(1:1)="H" 
7936: ! cs675>7920: ! cs675>
7937:       ELSE IF (WORD.EQ.'ORBITALS') THEN7921:       ELSE IF (WORD.EQ.'ORBITALS') THEN
7938:          ORBITALS=.TRUE.7922:          ORBITALS=.TRUE.
7939:          CALL READI(ORBVAREXPONENT)7923:          CALL READI(ORBVAREXPONENT)
7940:          CALL ORBITALS_INIT(NORBS, NROTS)7924:          CALL ORBITALS_INIT(NORBS, NROTS)
7941:          IF (NROTS.NE.NATOMS) THEN7925:          IF (NROTS.NE.NATOMS) THEN
7942:             WRITE(MYUNIT,'(A,2I8)') 'keywords> ERROR *** NATOMS,NROTS=',NATOMS,NROTS7926:             WRITE(MYUNIT,'(A,2I8)') 'keywords> ERROR *** NATOMS,NROTS=',NATOMS,NROTS
7943:             STOP7927:             STOP
7944:          ENDIF7928:          ENDIF
7945:       ELSE7929:       ELSE


r33036/potential.f90 2017-07-19 13:30:13.029791635 +0100 r33035/potential.f90 2017-07-19 13:30:14.145806533 +0100
1157:    ELSE IF (NCAPT) THEN1157:    ELSE IF (NCAPT) THEN
1158:       CALL NEWCAPSID (X, GRAD, EREAL, GRADT)1158:       CALL NEWCAPSID (X, GRAD, EREAL, GRADT)
1159: 1159: 
1160:    ELSE IF (NPAHT) THEN1160:    ELSE IF (NPAHT) THEN
1161:       CALL NEWPAH (X, GRAD, EREAL, GRADT)1161:       CALL NEWPAH (X, GRAD, EREAL, GRADT)
1162: 1162: 
1163:    ELSE IF (NTIPT) THEN1163:    ELSE IF (NTIPT) THEN
1164:       CALL NEWTIP (X, GRAD, EREAL, GRADT)1164:       CALL NEWTIP (X, GRAD, EREAL, GRADT)
1165: 1165: 
1166:    ELSE IF (TTM3T) THEN ! Xantheas' TTM3-F water potential1166:    ELSE IF (TTM3T) THEN ! Xantheas' TTM3-F water potential
1167:       ! EREAL=0.0D01167:       EREAL=0.0D0
1168:       ! GRAD(1:3*NATOMS)=0.0D01168:       GRAD(1:3*NATOMS)=0.0D0
1169:      CALL TTM3FCALL(NWATER, X, EREAL, GRAD)1169: !      CALL TTM3FCALL(NATOMS/3, X, EREAL, GRAD)
1170:       IF (RIGIDINIT .AND. .NOT. ATOMRIGIDCOORDT) THEN 
1171:          XRIGIDCOORDS(1:DEGFREEDOMS)=X(1:DEGFREEDOMS) 
1172:          CALL TRANSFORMRIGIDTOC(1, NRIGIDBODY, X, XRIGIDCOORDS) 
1173:       END IF 
1174:      CALL TTM3FCALL(NATOMS/3, X, EREAL, GRAD) 
1175:       IF (RIGIDINIT .AND. .NOT. ATOMRIGIDCOORDT) THEN 
1176:          X(DEGFREEDOMS+1:3*NATOMS)=0.0D0 
1177:          X(1:DEGFREEDOMS)=XRIGIDCOORDS(1:DEGFREEDOMS) 
1178:          IF (GRADT) THEN 
1179:             GRADATOMS(1:3*NATOMS) = GRAD(1:3*NATOMS) ! Needed for AACONVERGENCE 
1180:             CALL TRANSFORMGRAD(GRADATOMS, XRIGIDCOORDS, XRIGIDGRAD) 
1181:             GRAD(DEGFREEDOMS+1:3*NATOMS)=0.0D0 
1182:             GRAD(1:DEGFREEDOMS)=XRIGIDGRAD(1:DEGFREEDOMS) 
1183:          END IF 
1184:       END IF 
1185: 1170: 
1186:    ELSE IF (PATCHY) THEN1171:    ELSE IF (PATCHY) THEN
1187:       CALL PATCHYPOT (X, GRAD, EREAL, GRADT)1172:       CALL PATCHYPOT (X, GRAD, EREAL, GRADT)
1188: 1173: 
1189:    ELSE IF (ASAOOS) THEN ! this is not anisotropic1174:    ELSE IF (ASAOOS) THEN ! this is not anisotropic
1190:       CALL ASAOOSPOT (X, GRAD, EREAL, GRADT)1175:       CALL ASAOOSPOT (X, GRAD, EREAL, GRADT)
1191: 1176: 
1192:    ELSE IF (GBT) THEN1177:    ELSE IF (GBT) THEN
1193:       CALL GBCALAMITIC (X, GRAD, EREAL, GRADT)1178:       CALL GBCALAMITIC (X, GRAD, EREAL, GRADT)
1194: 1179: 
1291: 1276: 
1292:    ELSE IF (MWFILMT) THEN1277:    ELSE IF (MWFILMT) THEN
1293:       CALL MWFILM (X, NATOMS, GRAD, EREAL, GRADT, BOXLX, BOXLY, BOXLZ, LAT)1278:       CALL MWFILM (X, NATOMS, GRAD, EREAL, GRADT, BOXLX, BOXLY, BOXLZ, LAT)
1294: 1279: 
1295:    ELSE IF (POLIRT) THEN1280:    ELSE IF (POLIRT) THEN
1296:       CALL POLIR (X, GRAD, EREAL, GRADT)1281:       CALL POLIR (X, GRAD, EREAL, GRADT)
1297: 1282: 
1298:    ELSE IF (MBPOLT) THEN1283:    ELSE IF (MBPOLT) THEN
1299:       CALL MBPOL (X, EREAL, GRAD, GRADT)1284:       CALL MBPOL (X, EREAL, GRAD, GRADT)
1300: 1285: 
1301:    ELSE IF (RIGIDMBPOLT) THEN 
1302:       IF (RIGIDINIT .AND. .NOT. ATOMRIGIDCOORDT) THEN 
1303:          XRIGIDCOORDS(1:DEGFREEDOMS)=X(1:DEGFREEDOMS) 
1304:          CALL TRANSFORMRIGIDTOC(1, NRIGIDBODY, X, XRIGIDCOORDS) 
1305:       END IF 
1306:       CALL MBPOL (X, EREAL, GRAD, .TRUE.) 
1307:       IF (RIGIDINIT .AND. .NOT. ATOMRIGIDCOORDT) THEN 
1308:          X(DEGFREEDOMS+1:3*NATOMS)=0.0D0 
1309:          X(1:DEGFREEDOMS)=XRIGIDCOORDS(1:DEGFREEDOMS) 
1310:          ! IF (GRADT) THEN 
1311:             GRADATOMS(1:3*NATOMS) = GRAD(1:3*NATOMS) ! Needed for AACONVERGENCE 
1312:             CALL TRANSFORMGRAD(GRADATOMS, XRIGIDCOORDS, XRIGIDGRAD) 
1313:             GRAD(DEGFREEDOMS+1:3*NATOMS)=0.0D0 
1314:             GRAD(1:DEGFREEDOMS)=XRIGIDGRAD(1:DEGFREEDOMS) 
1315:          ! END IF 
1316:       do i=1, natoms 
1317:          write(*,*) grad(3*i+1: 3*i+3) 
1318:       end do 
1319:       ! stop 
1320:       END IF 
1321:  
1322:    ELSE IF (WATERDCT) THEN1286:    ELSE IF (WATERDCT) THEN
1323:       CALL WATERPDC (X, GRAD, EREAL, GRADT)1287:       CALL WATERPDC (X, GRAD, EREAL, GRADT)
1324: 1288: 
1325:    ELSE IF (WATERKZT) THEN1289:    ELSE IF (WATERKZT) THEN
1326:       CALL WATERPKZ (X, GRAD, EREAL, GRADT)1290:       CALL WATERPKZ (X, GRAD, EREAL, GRADT)
1327: 1291: 
1328:    ELSE IF (WATERMETHANET) THEN1292:    ELSE IF (WATERMETHANET) THEN
1329:       IF (.NOT. RIGIDINIT) THEN1293:       IF (.NOT. RIGIDINIT) THEN
1330:          WRITE(*,*) "WATER-METHANE POTENTIAL ONLY IMPLEMENTED FOR GENRIGID FRAMEWORK"1294:          WRITE(*,*) "WATER-METHANE POTENTIAL ONLY IMPLEMENTED FOR GENRIGID FRAMEWORK"
1331:          STOP1295:          STOP
1339:          X(DEGFREEDOMS+1:3*NATOMS)=0.0D01303:          X(DEGFREEDOMS+1:3*NATOMS)=0.0D0
1340:          X(1:DEGFREEDOMS)=XRIGIDCOORDS(1:DEGFREEDOMS)1304:          X(1:DEGFREEDOMS)=XRIGIDCOORDS(1:DEGFREEDOMS)
1341:          IF (GRADT) THEN1305:          IF (GRADT) THEN
1342:             GRADATOMS(1:3*NATOMS) = GRAD(1:3*NATOMS) ! Needed for AACONVERGENCE1306:             GRADATOMS(1:3*NATOMS) = GRAD(1:3*NATOMS) ! Needed for AACONVERGENCE
1343:             CALL TRANSFORMGRAD(GRADATOMS, XRIGIDCOORDS, XRIGIDGRAD)1307:             CALL TRANSFORMGRAD(GRADATOMS, XRIGIDCOORDS, XRIGIDGRAD)
1344:             GRAD(DEGFREEDOMS+1:3*NATOMS)=0.0D01308:             GRAD(DEGFREEDOMS+1:3*NATOMS)=0.0D0
1345:             GRAD(1:DEGFREEDOMS)=XRIGIDGRAD(1:DEGFREEDOMS)1309:             GRAD(1:DEGFREEDOMS)=XRIGIDGRAD(1:DEGFREEDOMS)
1346:          END IF1310:          END IF
1347:      END IF1311:      END IF
1348: 1312: 
1349:    ELSE IF (CLATHRATET) THEN 
1350:       IF (RIGIDINIT .AND. .NOT. ATOMRIGIDCOORDT) THEN 
1351:          XRIGIDCOORDS(1:DEGFREEDOMS)=X(1:DEGFREEDOMS) 
1352:          CALL TRANSFORMRIGIDTOC(1, NRIGIDBODY, X, XRIGIDCOORDS) 
1353:       END IF 
1354:       CALL clath_pot(X, GRAD, EREAL, GRADT, NWATER) 
1355:       IF (RIGIDINIT .AND. .NOT. ATOMRIGIDCOORDT) THEN 
1356:          X(DEGFREEDOMS+1:3*NATOMS)=0.0D0 
1357:          X(1:DEGFREEDOMS)=XRIGIDCOORDS(1:DEGFREEDOMS) 
1358:          IF (GRADT) THEN 
1359:             GRADATOMS(1:3*NATOMS) = GRAD(1:3*NATOMS) ! Needed for AACONVERGENCE 
1360:             CALL TRANSFORMGRAD(GRADATOMS, XRIGIDCOORDS, XRIGIDGRAD) 
1361:             GRAD(DEGFREEDOMS+1:3*NATOMS)=0.0D0 
1362:             GRAD(1:DEGFREEDOMS)=XRIGIDGRAD(1:DEGFREEDOMS) 
1363:          END IF 
1364:       END IF 
1365: 1313: 
1366:    ELSE IF (GAYBERNET) THEN1314:    ELSE IF (GAYBERNET) THEN
1367:       CALL GAYBERNE(X, GRAD, EREAL, GRADT,.FALSE.)1315:       CALL GAYBERNE(X, GRAD, EREAL, GRADT,.FALSE.)
1368: 1316: 
1369:    ELSE IF (PYGPERIODICT .OR. PYBINARYT) THEN1317:    ELSE IF (PYGPERIODICT .OR. PYBINARYT) THEN
1370:       IF (PYADDT) THEN1318:       IF (PYADDT) THEN
1371:          CALL PYGPERIODICADD(X, GRAD, EREAL, GRADT)1319:          CALL PYGPERIODICADD(X, GRAD, EREAL, GRADT)
1372:       ELSEIF (PYADD2T) THEN1320:       ELSEIF (PYADD2T) THEN
1373:          CALL PYGPERIODICADD2(X, GRAD, EREAL, GRADT)1321:          CALL PYGPERIODICADD2(X, GRAD, EREAL, GRADT)
1374:       ELSE1322:       ELSE


r33036/watermethane.f90 2017-07-19 13:30:13.249794573 +0100 r33035/watermethane.f90 2017-07-19 13:30:14.365809472 +0100
 35: double precision, parameter::    methanecharge(methanedof)= (/0.279901d0,0.279901d0,0.279901d0,0.279901d0, 3.590472d0, & 35: double precision, parameter::    methanecharge(methanedof)= (/0.279901d0,0.279901d0,0.279901d0,0.279901d0, 3.590472d0, &
 36:      -1.177519d0,-1.177519d0,-1.177519d0,-1.177519d0/) 36:      -1.177519d0,-1.177519d0,-1.177519d0,-1.177519d0/)
 37: !Parameters in Angstroms and kcal/mol 37: !Parameters in Angstroms and kcal/mol
 38: double precision::               betaang(waterdof, methanedof), Aang0(waterdof, methanedof), Aang1(waterdof, methanedof), AangM(waterdof, methanedof) 38: double precision::               betaang(waterdof, methanedof), Aang0(waterdof, methanedof), Aang1(waterdof, methanedof), AangM(waterdof, methanedof)
 39: double precision::               Cang6(waterdof, methanedof), Cang8(waterdof, methanedof), Cang10(waterdof, methanedof), deltaang6(waterdof, methanedof) 39: double precision::               Cang6(waterdof, methanedof), Cang8(waterdof, methanedof), Cang10(waterdof, methanedof), deltaang6(waterdof, methanedof)
 40: double precision::               deltaang8(waterdof, methanedof) 40: double precision::               deltaang8(waterdof, methanedof)
 41: !Parameters in atomic units 41: !Parameters in atomic units
 42: double precision::               beta(waterdof, methanedof), A0(waterdof, methanedof), A1(waterdof, methanedof), AM(waterdof, methanedof) 42: double precision::               beta(waterdof, methanedof), A0(waterdof, methanedof), A1(waterdof, methanedof), AM(waterdof, methanedof)
 43: double precision::               C6(waterdof, methanedof), C8(waterdof, methanedof), C10(waterdof, methanedof), delta6(waterdof, methanedof) 43: double precision::               C6(waterdof, methanedof), C8(waterdof, methanedof), C10(waterdof, methanedof), delta6(waterdof, methanedof)
 44: double precision::               delta8(waterdof, methanedof) 44: double precision::               delta8(waterdof, methanedof)
 45: !Unit conversions 
 46: double precision, parameter::    autoang= 0.529177d0, autokcal= 627.503d0 
 47: !-------------------------------------------------------------------------------------------------------------------- 45: !--------------------------------------------------------------------------------------------------------------------
 48: !Input of data 46: !Input of data
 49: !H1 and H2: 47: !H1 and H2:
 50: data betaang(1,:)/2.84808454d0,2.84808454d0,2.84808454d0,2.84808454d0,2.7971225d0,2.75581866d0,2.75581866d0,2.75581866d0,2.75581866d0/ 48: data betaang(1,:)/2.84808454d0,2.84808454d0,2.84808454d0,2.84808454d0,2.7971225d0,2.75581866d0,2.75581866d0,2.75581866d0,2.75581866d0/
 51: data betaang(2,:)/2.84808454d0,2.84808454d0,2.84808454d0,2.84808454d0,2.7971225d0,2.75581866d0,2.75581866d0,2.75581866d0,2.75581866d0/ 49: data betaang(2,:)/2.84808454d0,2.84808454d0,2.84808454d0,2.84808454d0,2.7971225d0,2.75581866d0,2.75581866d0,2.75581866d0,2.75581866d0/
 52: !Q: 50: !Q:
 53: data betaang(3,:)/2.86928398d0,2.86928398d0,2.86928398d0,2.86928398d0,2.3463075d0,2.31474866d0,2.31474866d0,2.31474866d0,2.31474866d0/ 51: data betaang(3,:)/2.86928398d0,2.86928398d0,2.86928398d0,2.86928398d0,2.3463075d0,2.31474866d0,2.31474866d0,2.31474866d0,2.31474866d0/
 54: !D1 and D2: 52: !D1 and D2:
 55: data betaang(4,:)/5.71995231d0,5.71995231d0,5.71995231d0,5.71995231d0,3.16999754d0,2.35594058d0,2.35594058d0,2.35594058d0,2.35594058d0/ 53: data betaang(4,:)/5.71995231d0,5.71995231d0,5.71995231d0,5.71995231d0,3.16999754d0,2.35594058d0,2.35594058d0,2.35594058d0,2.35594058d0/
 56: data betaang(5,:)/5.71995231d0,5.71995231d0,5.71995231d0,5.71995231d0,3.16999754d0,2.35594058d0,2.35594058d0,2.35594058d0,2.35594058d0/ 54: data betaang(5,:)/5.71995231d0,5.71995231d0,5.71995231d0,5.71995231d0,3.16999754d0,2.35594058d0,2.35594058d0,2.35594058d0,2.35594058d0/
149: !Q:147: !Q:
150: data deltaang8(3,:)/3.643903d0,3.643903d0,3.643903d0,3.643903d0,4.138380d0,4.368576d0,4.368576d0,4.368576d0,4.368576d0/148: data deltaang8(3,:)/3.643903d0,3.643903d0,3.643903d0,3.643903d0,4.138380d0,4.368576d0,4.368576d0,4.368576d0,4.368576d0/
151: !D1 and D2:149: !D1 and D2:
152: data deltaang8(4,:)/4.415080d0,4.415080d0,4.415080d0,4.415080d0,3.962102d0,3.741763d0,3.741763d0,3.741763d0,3.741763d0/150: data deltaang8(4,:)/4.415080d0,4.415080d0,4.415080d0,4.415080d0,3.962102d0,3.741763d0,3.741763d0,3.741763d0,3.741763d0/
153: data deltaang8(5,:)/4.415080d0,4.415080d0,4.415080d0,4.415080d0,3.962102d0,3.741763d0,3.741763d0,3.741763d0,3.741763d0/151: data deltaang8(5,:)/4.415080d0,4.415080d0,4.415080d0,4.415080d0,3.962102d0,3.741763d0,3.741763d0,3.741763d0,3.741763d0/
154: !T1 and T2:152: !T1 and T2:
155: data deltaang8(6,:)/methanedof*0.0d0/153: data deltaang8(6,:)/methanedof*0.0d0/
156: data deltaang8(7,:)/methanedof*0.0d0/154: data deltaang8(7,:)/methanedof*0.0d0/
157: 155: 
158: contains156: contains
159: !------------------------------------------------------------------------------------------- 
160: subroutine clath_pot(x, grad, ereal, gradt, nwater) 
161: implicit none 
162:  
163: double precision::     x(nwater*24+3*methanedof), grad(nwater*24+3*methanedof), ereal, ewater, ewatermeth 
164: double precision::     xwatermeth(51), tempgrad(51), xwater(nwater*9), gradwater(nwater*9) 
165: integer::              i,j,k, nwater 
166: logical::             gradt 
167: 157: 
168: ereal=0.0d0 
169: grad(:)=0.0d0 
170: do i=1, nwater 
171:    xwatermeth(:)=0.0d0 
172:    do j=1, waterdof+1 
173:       do k=1, 3 
174:          xwatermeth(3*(j-1) + k)= x(3*((waterdof+1)*(i-1) + j -1) +k)/autoang 
175:       end do 
176:    end do 
177:    xwatermeth(25:51)= x(nwater*(waterdof+1)*3 + 1: nwater*(waterdof+1)*3 + methanedof*3)/autoang 
178:    call wmrb(xwatermeth, tempgrad, ewatermeth, gradt) 
179:    ereal= ereal+ ewatermeth*autokcal 
180:    do j=1, waterdof+1 
181:       do k=1, 3 
182:          grad(3*((waterdof+1)*(i-1) + j -1) +k)= grad(3*((waterdof+1)*(i-1) + j -1) +k)& 
183:               +tempgrad(3*(j-1) + k)*autokcal/autoang 
184:       end do 
185:    end do 
186:    grad(nwater*(waterdof+1)*3 + 1: nwater*(waterdof+1)*3 + methanedof*3)= & 
187:         grad(nwater*(waterdof+1)*3 + 1: nwater*(waterdof+1)*3 + methanedof*3)+& 
188:         tempgrad(25:51)*autokcal/autoang 
189:    do k=1,3 
190:       xwater(3*(3*(i-1)+1 -1) + k)= x(3*((waterdof+1)*(i-1) + 8 -1) +k) 
191:       xwater(3*(3*(i-1)+2 -1) + k)= x(3*((waterdof+1)*(i-1) + 1 -1) +k) 
192:       xwater(3*(3*(i-1)+3 -1) + k)= x(3*((waterdof+1)*(i-1) + 2 -1) +k) 
193:    end do 
194: end do 
195:  
196: if (gradt) then 
197:   call mbpolenergygradient(nwater,ewater,xwater,gradwater) 
198: else 
199:   call mbpolenergy(nwater,ewater,xwater) 
200: endif 
201: do i=1,nwater 
202:    do k=1,3 
203:       grad(3*((waterdof+1)*(i-1) + 8 -1) +k)= grad(3*((waterdof+1)*(i-1) + 8 -1) +k) + & 
204:            gradwater(3*(3*(i-1)+1 -1) + k) 
205:       grad(3*((waterdof+1)*(i-1) + 1 -1) +k)= grad(3*((waterdof+1)*(i-1) + 1 -1) +k) + & 
206:            gradwater(3*(3*(i-1)+2 -1) + k) 
207:       grad(3*((waterdof+1)*(i-1) + 2 -1) +k)= grad(3*((waterdof+1)*(i-1) + 2 -1) +k) + & 
208:            gradwater(3*(3*(i-1)+3 -1) + k) 
209:    end do 
210: end do 
211:  
212:  
213: ereal= ereal+ewater 
214: return 
215: end subroutine clath_pot 
216: !-------------------------------------------------------------------------------------------158: !-------------------------------------------------------------------------------------------
217: subroutine wmrb(x, grad, ereal, gradt)159: subroutine wmrb(x, grad, ereal, gradt)
218: !function that calculates the potential from the water-methane interaction160: !function that calculates the potential from the water-methane interaction
219: !for rigid body161: !for rigid body
220: !atoms/sites need to be: H H Q D D T T O/ H H H H C M M M M162: !atoms/sites need to be: H H Q D D T T O/ H H H H C M M M M
221: implicit none163: implicit none
222: double precision::    x(51), grad(51), ereal164: double precision::    x(51), grad(51), ereal
223: double precision::    r12, rk165: double precision::    r12, rk
224: logical::             gradt166: logical::             gradt
225: integer::             i,j,k, totdof167: integer::             i,j,k, totdof


legend
Lines Added 
Lines changed
 Lines Removed

hdiff - version: 2.1.0