OPLS-aa + Tersoff for fullerene(C60)-derivative PCBM (updated)

[Jason]

Hi lammps users,

To clarify from previous messages, I realize that my issue is that I want to apply bonding terms on top of the Tersoff potential. This is not allowed, according to:

http://lammps.sandia.gov/doc/special_bonds.html

I realize that there is good reason to exclude many-body pair styles and bonding terms, as discussed in previous messages:
http://lammps.sandia.gov/threads/msg38173.html
http://lammps.sandia.gov/threads/msg38174.html

However, I believe this model (Tersoff + bonds) is physically reasonable, and necessary, for the stability of the PCBM molecule I am trying to study (also attached):

https://raw.github.com/jasonlarkin/disorder/master/pcbm/pcbm_pdb/vmd_output/pcbm_cheung.png

The reason to use Tersoff is so I do not have to produce the (large) list of bonds for the 5- and 6-member rings of the C60 (atoms labeled CF in pcbm_cheung.png).

The reason to use bonds on top of the Tersoff potential is to secure the satellite molecule to the C60, where the required bonds are between, for example, CF-CT-CA, CF-CT-CT, CF-CF-CT-CA, etc.

I have tested this potential using the lattice/molecular dynamics program GULP, and the system seems to be stable.

Is it possible to implement Tersoff+bonds in some other way with the current version of LAMMPS?

If not, is there an easy way to modify the source code for my own purpose, and which routine(s) do I need to modify?

Thanks, and sorry for any previous confusion,

Jason

[akohlmey]

Hi lammps users,

To clarify from previous messages, I realize that my issue is that I want to
apply bonding terms on top of the Tersoff potential. This is not allowed,
according to:

LAMMPS Molecular Dynamics Simulator

IMPORTANT NOTE: These weighting factors are NOT used by pair styles that
compute many-body interactions, since the "bonds" that result from such
interactions are not permanent, but are created and broken dynamically as
atom conformations change. Examples of pair styles in this category are EAM,
MEAM, Stillinger-Weber, Tersoff, COMB, AIREBO, and ReaxFF. In fact, it
generally makes no sense to define permanent bonds between atoms that
interact via these potentials, though such bonds may exist elsewhere in your
system, e.g. when using the pair_style hybrid command. Thus LAMMPS ignores
special_bonds settings when manybody potentials are calculated.

it *is* allowed. however, having fully(!) excluded pair interactions
will result in unphysical results
for the manybody potential.

however, if you have non-zero scaling factors, those pair interactions
*will* show up and they will work as before.

now the tricky thing is if you have both, a regular force field that
requires 1-2, 1-3, 1-4 interactions to be excluded and a case where
you want to overlay bonds with a manybody potential.

the easy way out would be to simply use a tiny scaling factor that is
not zero, say 1.0e-30. that will keep all pairs in the neighbor lists,
but still practically exclude them for the regular potentials
(remember, the scaling factor is ignored in manybody potentials).

the second option would be to add those bonds not via defining bonds,
but using pair style list in combination with pair style hybrid
overlay.
this kind of purpose (adding bond-like pairwise interactions to a
system) is what it was written for.
the input for it can be written with a little bit of VMD scripting for example.

axel.

[sjplimp]

If I understand your model (and I didn’t the first couple time),
I think this works straight-forwardly. Say the C60 is type 1 atoms,
and the molecule is type 2-N atoms.

Use pair hybrid with tersoff only defined for interactions
between type 1 atoms. You can have extra bonds
defined between type 1 and type M atoms (M != 1),
or between any of the type 2-N atoms. You can
also have another pair style defined between type 2-N
atoms (just as you would for any molecule). And
you don’t need a pair style defined between type 1 and M (M != 1)
so just define “none” for those interactions in pair hybrid.

Steve