LAMMPS WWW Site - LAMMPS Documentation - LAMMPS Mailing List Archives
Re: [lammps-users] TIP4P-2005 water model
[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: [lammps-users] TIP4P-2005 water model

From: Axel Kohlmeyer <akohlmey@...24...>
Date: Tue, 15 Aug 2017 08:43:21 -0400

On Fri, Aug 11, 2017 at 3:49 AM, 284237308@...1204... <284237308@...1204...> wrote:
Hello guys,
    I'm trying my water model with TIP4P-2005. Here is my script.
    However, when I start, it just stays running while nothing happens on the screen no matter how long I wait.
    So, is there anything wrong about my script?

​you have a force constant of 0 for both bond and angle interaction. that would be fine for as long as you *also* use fix shake, but that is not the case during minimization (because fix shake only works for MD). so there is nothing keeping your hydrogens from getting extremely close to the oxygens and thus leading to the explosion of potential energy and ultimately the divergence you see. GI-GO.



#Pure water

units           real
dimension       3
boundary        p p p 
atom_style      full


#Force field parameters
#H=1, O=2 

set             type 1 charge 0.5564
set             type 2 charge -1.1128              

#TIP4P-2005 water
pair_style      lj/cut/tip4p/long 2 1 1 1 0.1546 13.0 8.5
pair_modify     tail yes mix arithmetic
bond_style      harmonic
angle_style     harmonic
kspace_style    pppm/tip4p 1.0e-5

pair_coeff      1 1 0.0 0.0         #H-H
pair_coeff      2 2 0.1852 3.1589   #O-O
pair_coeff      1 2 0.0 0.0         #O-H
bond_coeff      1 0.0 0.9572        #H-O
angle_coeff     1 0.0 104.52        #H-O-H

neighbor        3.0 bin
neigh_modify    every 1 delay 1 check yes

group           water type 1 2
group           oxgen type 2


compute         tp water temp 
compute         peratom water stress/atom NULL
compute         pr water reduce sum c_peratom[1] c_peratom[2] c_peratom[3]
compute         pr3 water reduce sum c_peratom[4] c_peratom[5] c_peratom[6]

variable        T equal 300.0
variable        P equal 1.0
variable        dt equal 0.5
variable        ci equal 5                              #correlation interval
variable        si equal 2                              #sample interval
variable        ti equal ${ci}*${si}                    #total interval
variable        th equal 1000                           #thermo interval
variable        di equal 1000                           #dump interval
variable        Tdamp equal ${dt}*100
variable        Pdamp equal ${dt}*1000
variable        den equal density


timestep        ${dt}
thermo          ${th}
thermo_style    custom step c_tp ke pe etotal vol press v_den 

minimize        1.0e-4 1.0e-6 100 1000

fix             shakeh2o water shake 0.0001 200 0 b 1 a 1

fix             l all langevin ${T} ${T} ${Tdamp} 19930830
fix             e all nve 

dump            1 all atom ${di} langevin-nve.lammpstrj

run             10000

unfix           l
unfix           e 
undump          1

fix             1 all npt temp ${T} ${T} ${Tdamp} iso ${P} ${P} ${Pdamp}

dump            1 all atom ${di} npt.lammpstrj

run             10000

undump          1

Check out the vibrant tech community on one of the world's most
engaging tech sites,!
lammps-users mailing list

Dr. Axel Kohlmeyer  akohlmey@...24...
College of Science & Technology, Temple University, Philadelphia PA, USA
International Centre for Theoretical Physics, Trieste. Italy.