hdiff output

r31687/keywords.f 2017-01-21 10:35:54.334250000 +0000 r31686/keywords.f 2017-01-21 10:35:54.818250000 +0000
537:       GLJY = .FALSE.537:       GLJY = .FALSE.
538:       GLJ_EXP = 6538:       GLJ_EXP = 6
539:       YUK_A = 0.0D0539:       YUK_A = 0.0D0
540:       YUK_XI = 1.0D0540:       YUK_XI = 1.0D0
541:       541:       
542: !     ds656> Molecular Dynamics542: !     ds656> Molecular Dynamics
543:       MDT = .FALSE.543:       MDT = .FALSE.
544:       MD_NSTEPS = 0544:       MD_NSTEPS = 0
545:       MD_NWAIT = 1545:       MD_NWAIT = 1
546:       MD_NFREQ = 100546:       MD_NFREQ = 100
547:       MD_GAMMA = 0.0547:       MD_GAMMA = 0.1
548:       MD_TSTEP = 0.002548:       MD_TSTEP = 0.01
549: 549: 
550:       CHRMMT=.FALSE.550:       CHRMMT=.FALSE.
551:       CHARMMTYPE=1551:       CHARMMTYPE=1
552:       CHARMMDFTBT=.FALSE.552:       CHARMMDFTBT=.FALSE.
553:       ACESOLV=.FALSE.553:       ACESOLV=.FALSE.
554:       ACEUPSTEP=50554:       ACEUPSTEP=50
555:       CHRIGIDTRANST=.FALSE.555:       CHRIGIDTRANST=.FALSE.
556:       CHRIGIDROTT=.FALSE.556:       CHRIGIDROTT=.FALSE.
557:       CHNMIN=0.D0557:       CHNMIN=0.D0
558:       CHNMAX=HUGE(1.0D0)558:       CHNMAX=HUGE(1.0D0)


r31687/md.f90 2017-01-21 10:35:54.574250000 +0000 r31686/md.f90 2017-01-21 10:35:55.042250000 +0000
 11: !  but WITHOUT ANY WARRANTY; without even the implied warranty of 11: !  but WITHOUT ANY WARRANTY; without even the implied warranty of
 12: !  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the 12: !  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 13: !  GNU General Public License for more details. 13: !  GNU General Public License for more details.
 14: ! 14: !
 15: !  You should have received a copy of the GNU General Public License 15: !  You should have received a copy of the GNU General Public License
 16: !  along with this program; if not, write to the Free Software 16: !  along with this program; if not, write to the Free Software
 17: !  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA 17: !  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
 18: ! 18: !
 19: MODULE MOLECULAR_DYNAMICS 19: MODULE MOLECULAR_DYNAMICS
 20:   ! 20:   !
 21:   USE COMMONS, ONLY : NATOMS, MYUNIT, ATMASS, PRTFRQ, DUMPXYZUNIT, & 21:   USE COMMONS, ONLY : NATOMS, MYUNIT, ATMASS, PRTFRQ, DUMPXYZUNIT
 22:        TSTART 
 23:   !  22:   ! 
 24:   IMPLICIT NONE 23:   IMPLICIT NONE
 25:   ! 24:   !
 26:   ! Parameters set in the 'keywords' routine are to be saved! 25:   ! Parameters set in the 'keywords' routine are to be saved!
 27:   LOGICAL, SAVE :: MDT 26:   LOGICAL, SAVE :: MDT
 28:   INTEGER, SAVE :: MD_NWAIT, MD_NFREQ, MD_UNIT, MD_NSTEPS 27:   INTEGER, SAVE :: MD_NWAIT, MD_NFREQ, MD_UNIT, MD_NSTEPS
 29:   DOUBLE PRECISION, SAVE :: MD_GAMMA, MD_TSTEP 28:   DOUBLE PRECISION, SAVE :: MD_GAMMA, MD_TSTEP
 30:   ! 29:   !
 31:   ! MD-specific variables 30:   ! MD-specific variables
 32:   DOUBLE PRECISION :: EPOT, EKIN 31:   DOUBLE PRECISION :: EPOT, EKIN
 33:   DOUBLE PRECISION, ALLOCATABLE :: FORCE(:), VEL(:)  32:   DOUBLE PRECISION, ALLOCATABLE :: FORCE(:), VEL(:) 
 34:   ! 33:   !
 35:   DOUBLE PRECISION, PARAMETER :: PI=3.141592652D0 
 36:   ! 
 37: CONTAINS 34: CONTAINS
 38:   ! 35:   !
 39:   !============================================================ 36:   !============================================================
 40:   ! 37:   !
 41:   SUBROUTINE MDRUN(POS,NSTEPS,TEMP,DUMP) 38:   SUBROUTINE MDRUN(POS,NSTEPS,TEMP,DUMP)
 42:     ! 39:     !
 43:     IMPLICIT NONE 40:     IMPLICIT NONE
 44:     ! 41:     !
 45:     LOGICAL, INTENT(IN) :: DUMP 42:     LOGICAL, INTENT(IN) :: DUMP
 46:     INTEGER, INTENT(IN) :: NSTEPS 43:     INTEGER, INTENT(IN) :: NSTEPS
 47:     DOUBLE PRECISION, INTENT(IN) :: TEMP ! k_B*T in units of energy 44:     DOUBLE PRECISION, INTENT(IN) :: TEMP ! k_B*T in units of energy
 48:     DOUBLE PRECISION, INTENT(INOUT) :: POS(3*NATOMS) 45:     DOUBLE PRECISION, INTENT(INOUT) :: POS(3*NATOMS)
 49:     ! 46:     !
 50:     INTEGER :: I, J, K 47:     INTEGER :: I, J, K
 51:     DOUBLE PRECISION :: TIME 
 52:     ! 48:     !
 53:     IF (DUMP) WRITE(MYUNIT, '(A)') "mdrun> Start of MD run." 49:     IF (DUMP) WRITE(MYUNIT, '(A)') "mdrun> Start of MD run."
 54:         50:        
 55:     ! 51:     !
 56:     ALLOCATE(FORCE(3*NATOMS), VEL(3*NATOMS)) 52:     ALLOCATE(FORCE(3*NATOMS), VEL(3*NATOMS))
 57:     ! 53:     !
 58:     CALL INIVEL(POS, TEMP)  54:     CALL INIVEL(POS, TEMP) 
 59:     CALL POTENTIAL(POS,FORCE,EPOT,.TRUE.,.FALSE.) ! get gradient 55:     CALL POTENTIAL(POS,FORCE,EPOT,.TRUE.,.FALSE.) ! get gradient
 60:     FORCE(:) = -FORCE(:) ! get force     56:     FORCE(:) = -FORCE(:) ! get force    
 61:     ! 57:     !
 67:        ! 63:        !
 68:        CALL MDSTEP(POS, TEMP) 64:        CALL MDSTEP(POS, TEMP)
 69:        ! 65:        !
 70:        IF(DUMP) THEN ! Dump data to files 66:        IF(DUMP) THEN ! Dump data to files
 71:           ! 67:           !
 72:           IF(MOD(I,PRTFRQ).EQ.0) THEN ! Basic data 68:           IF(MOD(I,PRTFRQ).EQ.0) THEN ! Basic data
 73:              WRITE(MYUNIT, '(A,I8,A,F18.8,A,F14.8)') & 69:              WRITE(MYUNIT, '(A,I8,A,F18.8,A,F14.8)') &
 74:                   'mdrun> Step ', I,' EPOT= ',EPOT,' EKIN= ',EKIN 70:                   'mdrun> Step ', I,' EPOT= ',EPOT,' EKIN= ',EKIN
 75:           ENDIF 71:           ENDIF
 76:           ! 72:           !
 77:           IF(I == MD_NWAIT) THEN 73:           IF(I > MD_NWAIT .AND. &
 78:              CALL MYCPU_TIME(TIME) 
 79:              WRITE(MYUNIT, '(A,F11.1)') & 
 80:                   "mdrun> Finished MD equilibration t= ",TIME-TSTART 
 81:           ELSEIF(I > MD_NWAIT .AND. & 
 82:                MOD(I-MD_NWAIT,MD_NFREQ)==0) THEN ! XYZ dump 74:                MOD(I-MD_NWAIT,MD_NFREQ)==0) THEN ! XYZ dump
 83:              WRITE(DUMPXYZUNIT(1),'(I4)') NATOMS 75:              WRITE(DUMPXYZUNIT(1),'(I4)') NATOMS
 84:              WRITE(DUMPXYZUNIT(1),'(A,I9,A,F20.10)')  & 76:              WRITE(DUMPXYZUNIT(1),'(A,I9,A,F20.10)')  &
 85:                   'MD step ',I,'  EPOT= ', EPOT 77:                   'MD step ',I,'  EPOT= ', EPOT
 86:              DO J=1,NATOMS 78:              DO J=1,NATOMS
 87:                 WRITE(DUMPXYZUNIT(1),'(A,3(1X,F10.5))') & 79:                 WRITE(DUMPXYZUNIT(1),'(A,3(1X,F10.5))') &
 88:                      'X ', (POS(3*(J-1)+K), K=1,3) 80:                      'X ', (POS(3*(J-1)+K), K=1,3)
 89:              ENDDO              81:              ENDDO             
 90:           ENDIF 82:           ENDIF
 91:           ! 83:           !
 92:        ENDIF 84:        ENDIF
 93:        ! 85:        !
 94:     ENDDO 86:     ENDDO
 95:     ! 87:     !
 96:     DEALLOCATE(FORCE, VEL) 88:     DEALLOCATE(FORCE, VEL)
 97:     ! 89:     !
 98:     CALL MYCPU_TIME(TIME) 90:     IF (DUMP) WRITE(MYUNIT, '(A)') "mdrun> End of MD run."
 99:     IF (DUMP) WRITE(MYUNIT, '(A,F11.1)') & 
100:          "mdrun> Finished MD run t= ", TIME-TSTART 
101:     ! 91:     !
102:     RETURN 92:     RETURN
103:     ! 93:     !
104:   END SUBROUTINE MDRUN 94:   END SUBROUTINE MDRUN
105:   ! 95:   !
106:   !============================================================ 96:   !============================================================
107:   ! 97:   !
108:   SUBROUTINE MDSTEP(POS, TEMP) 98:   SUBROUTINE MDSTEP(POS, TEMP)
109:     ! 99:     !
110:     IMPLICIT NONE100:     IMPLICIT NONE
111:     !101:     !
112:     DOUBLE PRECISION, INTENT(IN) :: TEMP102:     DOUBLE PRECISION, INTENT(IN) :: TEMP
113:     DOUBLE PRECISION, INTENT(INOUT) :: POS(3*NATOMS)103:     DOUBLE PRECISION, INTENT(INOUT) :: POS(3*NATOMS)
114:     104:     
115:     INTEGER :: I,J,K105:     INTEGER :: I,J,K
116:     DOUBLE PRECISION :: HDT, GFRIC, NOISE(NATOMS), &106:     DOUBLE PRECISION :: HDT, HDT2, GFRIC, NOISE(NATOMS), &
117:          DPRAND,R1,R2,NR,NR2(3*NATOMS)107:          DPRAND
118:     !108:     !
119:     HDT=0.5D0*MD_TSTEP109:     HDT=0.5D0*MD_TSTEP
 110:     HDT2=HDT*MD_TSTEP
120:     GFRIC=1.0D0-MD_GAMMA*HDT111:     GFRIC=1.0D0-MD_GAMMA*HDT
121:     NOISE(:) = DSQRT(TEMP*MD_GAMMA*MD_TSTEP/ATMASS(:))112:     NOISE(:) = DSQRT(6.0D0*TEMP*ATMASS(:)*MD_GAMMA/MD_TSTEP)
122:     !113:     !
123:     K=0114:     K=0
124:     DO I=1,NATOMS115:     DO I=1,NATOMS
125:        DO J=1,3116:        DO J=1,3
126:           K=K+1117:           K=K+1
127:           ! Generate norman random numbers118:           POS(K) = POS(K) + VEL(K)*MD_TSTEP + FORCE(K)*HDT2/ATMASS(I)
128:           R1=DSQRT(-2.0D0*LOG(DPRAND()))119:           VEL(K) = GFRIC*VEL(K) + FORCE(K)*HDT/ATMASS(I)
129:           R2=2.0*PI*DPRAND() 
130:           NR=R1*COS(R2) 
131:           NR2(K)=R1*SIN(R2)           
132:           ! Update mid-step velocity 
133:           VEL(K) = GFRIC*VEL(K) + FORCE(K)*HDT/ATMASS(I) & 
134:                + NR*NOISE(I) 
135:           ! Update full-step position 
136:           POS(K) = POS(K) + VEL(K)*MD_TSTEP 
137:        ENDDO120:        ENDDO
138:     ENDDO121:     ENDDO
139:     !122:     !
140:     CALL POTENTIAL(POS,FORCE,EPOT,.TRUE.,.FALSE.) ! GET GRADIENT123:     CALL POTENTIAL(POS,FORCE,EPOT,.TRUE.,.FALSE.) ! GET GRADIENT
141:     FORCE(:)= -FORCE(:) 
142:     !124:     !
143:     K=0125:     K=0
144:     EKIN=0.0D0126:     EKIN=0.0D0
145:     DO I=1,NATOMS127:     DO I=1,NATOMS
146:        DO J=1,3128:        DO J=1,3
147:           K=K+1129:           K=K+1
148:           ! Update full-step velocity130:           FORCE(K)= -FORCE(K) + (2.0D0*DPRAND()-1.0D0)*NOISE(I)
149:           VEL(K) = GFRIC*VEL(K) + FORCE(K)*HDT/ATMASS(I) &131:           VEL(K) = GFRIC*VEL(K) + FORCE(K)*HDT/ATMASS(I)
150:                + NR2(K)*NOISE(I) 
151:           EKIN = EKIN + ATMASS(I)*VEL(K)*VEL(K)132:           EKIN = EKIN + ATMASS(I)*VEL(K)*VEL(K)
152:        ENDDO133:        ENDDO
153:     ENDDO134:     ENDDO
154:     EKIN=0.5D0*EKIN135:     EKIN=0.5D0*EKIN
155:     !136:     !
156:     RETURN137:     RETURN
157:     !138:     !
158:   END SUBROUTINE MDSTEP139:   END SUBROUTINE MDSTEP
159:   !140:   !
160:   !============================================================141:   !============================================================
162:   SUBROUTINE INIVEL(POS, TEMP)143:   SUBROUTINE INIVEL(POS, TEMP)
163:     !144:     !
164:     IMPLICIT NONE145:     IMPLICIT NONE
165:     !146:     !
166:     DOUBLE PRECISION, INTENT (IN) :: TEMP, POS(3*NATOMS)147:     DOUBLE PRECISION, INTENT (IN) :: TEMP, POS(3*NATOMS)
167:     !148:     !
168:     INTEGER :: I, J, K, L, IPIVOT(3)149:     INTEGER :: I, J, K, L, IPIVOT(3)
169:     DOUBLE PRECISION :: COM(3), COM_MOM(3), ANG_MOM(3), &150:     DOUBLE PRECISION :: COM(3), COM_MOM(3), ANG_MOM(3), &
170:          MOI(3,3), MOI_INV(3,3), DPRAND, SCALE, TMASS, &151:          MOI(3,3), MOI_INV(3,3), DPRAND, SCALE, TMASS, &
171:          P(3), R(3), W(3), R2, X(NATOMS,3)152:          P(3), R(3), W(3), R2, X(NATOMS,3)
 153:     DOUBLE PRECISION, PARAMETER :: PI=3.141592652D0
172:     !154:     !
173:     COM(:) = 0.0D0155:     COM(:) = 0.0D0
174:     COM_MOM(:) = 0.0D0156:     COM_MOM(:) = 0.0D0
175:     TMASS = 0.0D0157:     TMASS = 0.0D0
176:     !158:     !
177:     ! On first pass randomly generate velocities159:     ! On first pass randomly generate velocities
178:     K=0    160:     K=0    
179:     DO I=1,NATOMS161:     DO I=1,NATOMS
180:        TMASS = TMASS + ATMASS(I) 162:        TMASS = TMASS + ATMASS(I) 
181:        DO J=1,3163:        DO J=1,3


legend
Lines Added 
Lines changed
 Lines Removed

hdiff - version: 2.1.0