[lammps-users] Special Bonds 1-3 Interactions

LAMMPS LAMMPS Mailing List Mirror
1

Hello,

I notice on the Docs page for the special_bonds page that: “For purposes of computing weighted pairwise interactions, 1-3 and 1-4 interactions are not defined from the list of angles or dihedrals used by the simulation. Rather, they are inferred topologically from the set of bonds defined when the simulation is defined from a data or restart file.”

However, one can use the dihedral keyword in the special_bonds command to effectively force the 1-4 weighting factor to only be applied if atoms 1 and 4 are explicitly listed as the first and last atoms in a dihedral defined in a data file. Is there any way to do something similar for the 1-3 weighting factor?

I have a situation where a thiol-terminated organic molecule is bonded to a metallic surface, and while the molecule itself should retain the usual 1-3 and 1-4 weighting factors, an atom that is bonded to a metallic atom should interact with all other metallic atoms via a non-bonded potential. Essentially I need those atoms that are not explicitly listed as the first and last atoms in a bond angle definition to interact via a LJ interaction.

Any tips/thoughts?

Thanks,

William French

2

I have a situation where a thiol-terminated organic molecule is bonded
to a metallic surface, and while the molecule itself should retain the
usual 1-3 and 1-4 weighting factors, an atom that is bonded to a
metallic atom should interact with all other metallic atoms via a
non-bonded potential. Essentially I need those atoms that are not
explicitly listed as the first and last atoms in a bond angle
definition to interact via a LJ interaction.

william,

whether this is a problem or not, strongly depends on the
model that you are using to describe the metal surface.
does it have explicit bonds? if not, which is what i would
expect, then i don't quite see the problem.

cheers,
    axel.

3

Axel,

I am modeling the metal (gold) surface with a many-body potential. However, I also have a situation where a sulfur atom may bond to two adjacent Au atoms at a “bridge” site on the metallic surface. In this case, I model the bonded interaction with a Morse potential between each S-Au pair, and then define a “bond” between the two Au atoms so as to constrain the distance between these two atoms. This make it a sort of three-body case, I suppose.

Unfortunately, this leads to the situation I am currently having trouble addressing – in certain instances I have an Au-S-Au-S-Au configuration, where each S is bonded to two gold atoms. “Bonds” are also applied between Au1 and Au3 and Au3 and Au5, so based on the topology of this configuration LAMMPS does not compute nonbonded interactions for Au1-S4 or Au5-S2, though I would like for them to be. I am also seeing now that these same two pairs could be considered not only 1-3 pairs but also 1-4 pairs from their bonding topology. Eg, Au1-Au3-S4 and Au1-S2-Au3-S4.

I suppose if I could constrain the Au-Au distance within a bridge site in some way other than having to define a bond, this would also solve my problem.

Will

4

will,

Axel,

I am modeling the metal (gold) surface with a many-body potential. However,
I also have a situation where a sulfur atom may bond to two adjacent Au
atoms at a "bridge" site on the metallic surface. In this case, I model the
bonded interaction with a Morse potential between each S-Au pair, and then
define a "bond" between the two Au atoms so as to constrain the distance
between these two atoms. This make it a sort of three-body case, I suppose.

yes, very much so. i have some personal experience with this
kind of system: a couple of years ago, i did a very difficult and
time consuming DFT CPMD calculation that was used in combination
with other calculations on smaller systems to help interpreting
experimental results. http://dx.doi.org/10.1126/science.1158532

some people on the project decided to develop a many-body
potential based on the results, but i have not yet heard that they
were successful. it is a very tricky system due to the sulfur atoms'
ability to pull gold atoms out of the flat surface and the dynamic
nature of this system.

Unfortunately, this leads to the situation I am currently having trouble
addressing -- in certain instances I have an Au-S-Au-S-Au configuration,
where each S is bonded to two gold atoms. "Bonds" are also applied between
Au1 and Au3 and Au3 and Au5, so based on the topology of this configuration
LAMMPS does not compute nonbonded interactions for Au1-S4 or Au5-S2, though
I would like for them to be. I am also seeing now that these same two pairs
could be considered not only 1-3 pairs but also 1-4 pairs from their bonding
topology. Eg, Au1-Au3-S4 and Au1-S2-Au3-S4.

I suppose if I could constrain the Au-Au distance within a bridge site in
some way other than having to define a bond, this would also solve my
problem.

hmmm... what you could do is to define a non-bonded potential
to be added to your system via pair_style hybrid/overlay to add
the constraints. based on the fact that, according to the quantum
chemical calculations, your model cannot be overly accurate to
boot, it should be acceptable. you can just have multiple bond
types representing gold, map them all to your gold many-body
potential and then map them selectively with an appropriately
short cutoff to the non-bonded add-on potential. if you want it
to be harmonic, you can just use a tabulated potential.

that may take some effort to set up and will not allow dynamical
changes, but should sidestep the exclusions issue completely.

whether the total model will be useful. is hard to predict.
after my experience with the DFT calculation (i had to run
on a full rack of an IBM blue gene/l for over a month solid)
it will be extremely challenging to get anything with statistical
relevance out of higher level methods.

cheers,
     axel.

5

Axel,

I’ve considered using some sort of nonbonded constraint with a short cutoff between Au bridge atoms, but the problem is that for the Au-S-Au motif, LAMMPS considers these two Au atoms to be a 1-3 pair, and thus by default will not compute any nonbonded interaction between the two atoms. I think it would be nice if there was a option in the special_bonds command in which the user could explicitly tell LAMMPS the atoms he/she wants 1-3 interactions to be computed for (much like the option that already exists for dihedrals). I don’t know how easy this would be to code up, but in terms of giving the user more control over his/her simulation, it would certainly seem to be useful. Especially in cases where the user is hoping to create three-body or many-body terms, in which case determining the conventional 1-3 or 1-4 pairs based on bonding topology is not necessarily the best approach.

Will

6

However, one can use the dihedral keyword in the special_bonds command to effectively force the 1-4 weighting factor to >only be applied if atoms 1 and 4 are explicitly listed as the first and last atoms in a dihedral defined in a data file. Is there >any way to do something similar for the 1-3 weighting factor?

This could be added. We added the "dihedral" option b/c there was a specific
force field that did it this way that we wanted to match - can't
recall which one it was.

Steve

7

Steve,

I see the dihedral_trim() function in special.cpp that’s responsible for deleting 1-4 pairs if they are not explicit in the data file. I suppose if I want to do the same for 1-3 pairs then I can more or less use this function as my template. Can you think of any details that would make an angle-trim() function any different to implement in comparison to the dihedral_trim() function?

Thank you,

Will

8

Will,

Sounds reasonable to me — I don’t know of major differences or details you should be aware of.

Paul

9

Steve,

I have added the functionality to the special.cpp and related files so that angle terms that are not explicitly listed by a user in his/her data file will not be considered 1-3 pairs by LAMMPS, much like is already available for the dihedral terms. To do this, I created a new function in the special.cpp source file called angle_trim(), which is analogous to the dihedral_trim() function that was already implemented. The new feature works as expected in all the examples I’ve tested.

I will be happy to send you the updated files if you’re interested in including this feature in LAMMPS.

Best Regards,

Will

10

Will,

Very good. Thank you for implementing this. Please do send it along to Steve and me.

Thanks,

Paul