Dear,
I was trying to use the special bonds command in LAMMPS to neglect some of the intramolecular coulombic and LJ interactions in a cluster of molecules. I defined all possible lj/cut/coul/long pairs in a .nonbonding file and I scripted a .data file containing all information about bond, angle and dihedral parameters from the molecular FF. At the beginning, I thought it could be useful to introduce in the .data file the “Special Bond Counts” and “Special Bonds” sections to define which intramolecular interactions were to be neglected together with:
special_bonds lj 0.5 1 0 coul 0.833 1 0
in the input file to scale target 1-4 interactions and consider and exclude some of them. However, it seems that this kind of formatting if only valid for molecule templates. I guess I can generate such type of file and then provide lammps the structure of my molecular cluster in a valid format but I was wondering if there is any other way of achieving what I am trying because I also have a slab of a metallic oxide that I would like to introduce in my simulations.
Thank you!
Since you are using long-range Coulomb, you need more than a setting to exclude pair-wise interactions, you also need to remove the contribution from the kspace style.
The molecule file settings are only meant to be use in cases where the system does not yet allow to rebuild the special bond information. So those are not “permanent”.
In general, you need a very, very good justification to do what you are mentioning, since in most cases, you are likely to introduce some inconsistency and then whatever you may gain from the modifications, it may not be useful, since you may lose even more through the inconsistencies.
My main concern is that I am limited due to the QM-derivation procedure itself. I would have three types of pair interactions between 13 types of atoms:
buck/coul/long for atom types 1 and 2 (each atom is read a single molecule)
For the remaining 11 types, I would have some intramolecular lj/cut and lj/cut /coul/long for the intermolecular interactions between them but also with atom types 1 and 2.
I thought of hybrid/molecular setting up lj/cut and then lj/cut/coul/long but this cannot be implemented taking also into account the buckingham potential between types 1 and 2… am I right?
I scripted all pairs defining every single atom as a single atom type but then I am limited on the amount of memory in my machines. Can it be done by overriding some pair interactions with several “run 0” commands or creating intermediate .data files?
@manuperezesc I don’t know enough about your research and the details of the model (and I am not interested in spending the time needed learning it), but there are some general guidelines that I can offer from what I learned from reading and responding to questions about LAMMPS for 15+ years.
- I am always wary of elaborate hybrid potential setups. They are almost always an inferior solution to a fully consistent treatment. Have you heard of the saying the the cure may be worse than the disease? I suspect your worries about assigning different special_bonds settings are misguided since the correction it may offer is rather less than the error introduced by using a hybrid potential setup. If you have two different phases, it is often ok to use a different model for each phase, but that leaves the challenge of finding a suitable way to describe the interaction between those two models. Since your model includes Coulomb interactions this can be particularly difficult if you take into account my second point.
- You mention a metal oxide. Those are particularly difficult to model and particularly, you have to spend great attention to the surface(s). There is often significant reconstruction; there is a high risk of creating unphysical polar surfaces, and the charge distribution is usually different from the bulk and thus some form of polarization handling like charge equilibration is usually needed. This is often incompatible with how (fixed mean field) charges are determined for conventional molecular force fields. So the challenge becomes the same as it is for, say QM/MM.
- I am not a great fan of ReaxFF since its parameterizations are usually inferior to those for non-reactive models and I often see people use ReaxFF not to model reactions, but because they are too lazy (or too inexperienced?) to do the atom typing procedure required by conventional classical force fields). But in this case, it offers the unique opportunity to generate a single parameter set that can handle both, oxides (including charge equilibration) and molecular systems. Only, the generation of the parameter set tends to be a bit of a black art and will be rather specific to the model at hand.