Go-like interactions in LAMMPS

Hello,

Does LAMMPS have the capability to easily simulate Go-like interactions? By Go-like, I mean:

  1. Non-bonded interactions with parameters (epsilon, sigma) that are unique to atom i and atom j (atom ID not type).
  2. Bonded interactions (bonds, angles, dihedrals) with equilibrium distances and angles that are unique to the atom IDs of the constitutive atoms.

Obviously, one might specify a unique type for every bond and provide the corresponding coefficients; however, that seems rather crude.

Thank you,

Dan

Hi, Dan:

There isn’t a good way to attach a specific set of sigma and epsilon values to individual atoms; most LAMMPS pair potentials are geared to using atom “types” to determine the coefficients.

So you’d need to assign each atom to its own type, tinker with the software to store the pair potential as part of the atom type, or find some way to “import” this information and make it part of the pair potential calculation (the way parameter information is read in for the EAM or ReaxFF pair styles). This would require a decent amount of coding, although it is possible.

—AEI

I suggest you write a pair_go.cpp pair style that
computes the interactions you want. You'd still
need to define the param for each atom, which you
could do either with a new atom_style (also not hard
to create). Or with a fix that stored the additional
values/atom. The fix could read the definitions out
of the data file. There is a hook in read_data to enable
that. The fix could also store the per-atom bond/angle
params, used by a bond_go.cpp style, etc.

Steve

If the number of atoms/beads in your polymer is small (n < 10^2), then
I think it would be almost as efficient, (and a lot less work) to
assign each atom in your polymer a different type, and generate a text
file containing a long list of n*(n+1)/2 "pair_coeff" commands (one
for every pair of different atoms/beads in your polymer). As usual,
you can implement purely repulsive interactions either by setting the
distance cutoff to the location of the minima in the Lennard-Jones
interaction, or by taking the limit of small-epsilon & large-sigma to
approximate 1/r^12 repulsion. (Or use the pairs style described
below.)

  The only scenario I can think of where you would need to implement
your own pair_style would be if you need Go polymers containing
thousands of atoms or beads. (For reference, I think most
protein-domains are on the order of 100 amino acids or less.) (The
case of multiple polymers is discussed below.) I could be wrong about
all this. Try it the easy way, using existing LAMMPS pair styles and
see if it is fast enough for your needs.

Cheers!
Andrew

  --- multiple polymers ---
If you have multiple polymers in your solution,and you want to modify,
weaken, or prevent attractive interactions between different polymers,
then you can use "pair_style lj/charmm/coul/charmm/inter" which is
available here:
http://moltemplate.org/lammps_code/index.html
and documented here:
http://moltemplate.org/lammps_code/pair_lj_charmm_coul_charmm_inter.html
(does not work with CUDA)
In this case, purely-repulsive interactions (1/r^12) can be
implemented in this pair_style using the "eskl" option:
and setting the "K" coefficient to 0. (With this pair_style, you
don't have the freedom to customize each cutoff radius at the moment.)
Cheers again.

dan,

for a local project here i just coded up a pair style that allows to define individual interactions between pairs of atoms identified by their atom id.

so far, it can handle 12-6 LJ, Morse and harmonic bond-like interactions and reads the information simply from a file.
the benefit is that it requires (much) less storage. it was conceived to add bond-like interactions via pair_style hybrid/overlay, but it may address your issue 1) as well.

if you think this would be of use to you, please let me know and i can provide the code and more info.

for point 2) i don’t see an alternative. you’d need to store that kind of information in a rather similar way in any case.

axel.