Hello LAMMPS users,
I am using LAMMPS to do MD simulations of C70 fullerene, with the aim currently being to see thermal vibrations in the molecule. I will later compress the molecule to see how it behaves.
The config file I am using is:
This type of bad dynamics resulted from bad structure, potential, or setup. Please search the mailing list archive on “bad dynamics” and you would find plenty of discussions. In addition, mixing Tersoff with a bond style would not a good model.
Ray
Hello LAMMPS users,
[...]
# ---------- Define Interatomic Potential ---------------------
#bond_style harmonic
#bond_coeff 1 86.84 1.4
pair_style tersoff
pair_coeff * * SiC.tersoff C
neighbor 2.0 bin
neigh_modify delay 10 check yes
[...]
Where the data file is:
#####################################
70 atoms
105 bonds
1 atom types
1 bond types
[...]
#################################
This simulations causes the atoms to collapse and is not giving any useful
results for thermal equilibration of the fullerene.I tried to add another bonded harmonic interaction but it did not improve
the situation much.Any suggestions to where I am going wrong will be highly appreciated.
when doing simulations with manybody potentials like tersoff, you must
not have any bonds defined in the data file, as bonded interactions
will be computed implicitly through the pair style, and the bonds will
result in particularly those pairs being excluded from the neighbor
lists. as a temporary workaround, you can try adding the following
command:
special_bonds lj/coul 1.0 1.0 1.0
axel.
Hello,
Thanks for the quick replies.
The problem, as pointed out by Dr. Kohlmeyer, was with defining bonds with along with a Tersoff Potential. The thermal equilibration as well as compressive stress analysis are giving feasible results now.
I had another question, while compressing the fullerene, how do I analyse if any bond breaks? On visualising the trajectory file, none of the bond breaks even if I compress the fullerene completely.
Regards
P. Mathur
Hello LAMMPS users,
[…]
---------- Define Interatomic Potential ---------------------
#bond_style harmonic
#bond_coeff 1 86.84 1.4
pair_style tersoff
pair_coeff * * SiC.tersoff C
neighbor 2.0 bin
neigh_modify delay 10 check yes
[…]
Where the data file is:
#####################################
70 atoms
105 bonds
1 atom types
1 bond types
[…]
#################################
This simulations causes the atoms to collapse and is not giving any useful
results for thermal equilibration of the fullerene.I tried to add another bonded harmonic interaction but it did not improve
the situation much.Any suggestions to where I am going wrong will be highly appreciated.
when doing simulations with manybody potentials like tersoff, you must
not have any bonds defined in the data file, as bonded interactions
will be computed implicitly through the pair style, and the bonds will
result in particularly those pairs being excluded from the neighbor
lists. as a temporary workaround, you can try adding the following
command:
special_bonds lj/coul 1.0 1.0 1.0
axel.
Hello,
Thanks for the quick replies.The problem, as pointed out by Dr. Kohlmeyer, was with defining bonds with
along with a Tersoff Potential. The thermal equilibration as well as
compressive stress analysis are giving feasible results now.I had another question, while compressing the fullerene, how do I analyse if
any bond breaks? On visualising the trajectory file, none of the bond breaks
even if I compress the fullerene completely.
this has little to do with the simulation, but is a technical and a
conceptual matter.
the technical aspect: visualization programs for molecules usually
have no concept of chemistry and generally assume that atoms that are
close to each other are bound and atoms that are significantly apart
are not bound. thus they usually use a simple distance cutoff
criterion with a large margin. programs that are primarily used for
visualizing classical MD trajectories and large system, especially of
biological relevance (e.g. VMD) also assume by default that bonds
don't change (they cannot in simple class 1 force fields) and thus
maintain the bonding pattern from the first frame. in short, you
cannot depend on the visualization program to tell you reliably where
a bond is and when exactly it breaks.
the conceptual aspect: it is actually not so simple to define
precisely when two atoms are bound with a chemical bond and when not.
it is definitely not a on/off relation, but a gradual process and it
is also not described by only a distance criterion. from inside the
simulation software, one can define a "bond order" that quantifies the
"amount of bonding", but few visualization tools can be interfaced to
utilize these properties.
so you have to figure out first, what suitable parameters in your
specific case are, that indicate a bond being (fully/partially?)
formed or (fully/partially?) broken, and *only then* you can start
thinking about how to implement and visualize that with your
visualization tool of your choice.
axel.