Dear Lammps users
Im using reaxFF to model the interaction of a protein with a nanotube. and im hoping to run a md simulation at 300 K
this is a part of my input script
velocity all create 300 87287 rot yes dist gaussian
min_style sd
minimize 0 1e-8 2000 10000
fix 1 all nvt 300.0 300.0 100.0
compute 1 all temp
dump 2 all atom 1000 dumpfile.lammpsstrj.*
thermo 100
thermo_style step epair temp
run 1000
undump 2
but the temperature of the system drops gradually as I keep on running the simulation. it starts from 300 K. what might be the reason for this.
best
Milinda Samaraweera
University of Connecticut
Department of Chemistry
55 N Eagleville road
unit 3060
Storrs CT
USA
Dear Lammps users
Im using reaxFF to model the interaction of a protein with a nanotube. and
hmmm... that is not the kind of system where
i would reax expect to be very applicable.
im hoping to run a md simulation at 300 K
this is a part of my input script
how many more times do i and others have to say
that a fragment of an input is *useless* in order to
debug a complex problem? if it is not easily possible
to reproduce a run and experiment with parameters
and settings, nobody will be able to provide accurate
help. all that is left is guessing.
velocity all create 300 87287 rot yes dist gaussian
min_style sd
minimize 0 1e-8 2000 10000
fix 1 all nvt 300.0 300.0 100.0
compute 1 all temp
dump 2 all atom 1000 dumpfile.lammpsstrj.*
thermo 100
thermo_style step epair temp
run 1000
undump 2
but the temperature of the system drops gradually as I keep on running the
simulation. it starts from 300 K. what might be the reason for this.
what? after 1000 steps?
why should it not drop?
you had run a minimization before.
1000 MD steps is *nothingness*.
axel.
Not much to add to Axel's reply, except that you have to think about
equipartitioning when using the velocity command! A velocity command
with an argument of 300K just adds kinetic energy. This will move into
potential energy and lower the temperature (to about half of the set
value (if you are close to the harmonic approximation)) even without a
thermostat after a handful of simulation steps.
Hi Dr. Axel
my full input is
#output stuff
echo both
newton on
#Initialization
dimension 3
processors 0 0 1
boundary s s s
units real
atom_style full
log output.txt
read_data 1PristineCNT-AI-REBO
#potential information
pair_style reax
OPLS considers 1-4 interactions with 50%.
#Force field parameters
pair_coeff * * ffield.reax 1 2 4 1 1 2 2 4 3 1 4 2 1
#Bond coefficients
#Angle Coefficients
#Dihedral Coefficients
velocity all create 300 87287 rot yes dist gaussian
min_style sd
minimize 0 1e-8 2000 10000
fix 1 all nve
compute 1 all temp
dump 2 all atom 1000 dumpfile.lammpsstrj.*
thermo 100
thermo_style step epair temp
run 10000
undump 2
the temperature keeps on decreasing from 300 K, do you think the input script has a bug?
Milinda Samaraweera
University of Connecticut
Department of Chemistry
55 N Eagleville road
unit 3060
Storrs CT
USA
Hi Dr. Axel
my full input please disregard the earlier e-mail
#output stuff
echo both
newton on
#Initialization
dimension 3
processors 0 0 1
boundary s s s
units real
atom_style full
log output.txt
read_data 1PristineCNT-AI-REBO
#potential information
pair_style reax
OPLS considers 1-4 interactions with 50%.
#Force field parameters
pair_coeff * * ffield.reax 1 2 4 1 1 2 2 4 3 1 4 2 1
#Bond coefficients
#Angle Coefficients
#Dihedral Coefficients
velocity all create 300 87287 rot yes dist gaussian
min_style sd
minimize 0 1e-8 2000 10000
fix 1 all nvt 300.0 300.0 100.0
compute 1 all temp
dump 2 all atom 1000 dumpfile.lammpsstrj.*
thermo 100
thermo_style step epair temp
run 10000
undump 2
the temperature keeps on decreasing from 300 K, do you think the input script has a bug?
Milinda Samaraweera
University of Connecticut
Department of Chemistry
55 N Eagleville road
unit 3060
Storrs CT
USA
More to consider:
The default timestep in lammps is 1.0fs for real units but when using
reaxff (depending on the system) you might need to go down to ~ 0.1 fs
in order to properly capture reactions.
See the reaxff papers. Your pair_coeff section seem bogus. Why so many
entries if you only have 4 atom types?
Carlos
Just to follow up on Carlos' comment. For a system of CHNOS I have not been able to use a timestep above 0.1 fs, the temperature is not stable.
The pair coefficient line must match the order of Masses in your data file. If the order of atom type masses is - C, H, N, O then the pair_coeff line would read * * field.reax 1 2 4 3
If the order is CHON then the pair_coeff would end with 1 2 3 4.
This is an easy mistake (I have done it many times) and it can be hard to detect until your simulation has run for some time - basically you are assigning the wrong atom parameters. In your case I think the you would have assigned the atom types as C, H, N and C (1 2 4 1) with the rest of the numbers on the line being ignored (on the presumption that there were only for elements in the Masses table).
Kevin