Modelling Thermalized Drude ions

Hello

So I want to simulate a system composed of water and Na+ and Cl- ions using drude oscillators. I follow the method described here: 8.5.7. Tutorial for Thermalized Drude oscillators in LAMMPS — LAMMPS documentation, using the SWM4-NDP water model. This paper describes how to model polarizable ions: https://doi.org/10.1021/ct900576a. The paper describes the “polarization catastrophe” which I think I have encountered regarding the Cl- ions. I guess the problem is that when the drude and core particle have high negative and positive atomic charges (which is the case for Cl-, Br-, I-), then the drude particle can get far away from the core particle and the simulation becomes unstable. As explained in the paper: “The simplest treatment is to introduce an additional anharmonic restoring force to prevent excessively large excursions of the Drude particle away from the atom”. For highly polarizable anions such as Br- and I- (and Cl- I think), they implement an anharmonic restoring force active only beyond a certain stretching distance:

U_hyp = K_hyp · (ΔR - ΔR_cut)^4 , where ΔR_cut is the cutoff distance and ΔR is the distance.

K_hyp is chosen to be a very high value (40 000 kcal/mol/Å^4) so that this potential effectively reduces the distance between the drude and core particle, also reducing the induced dipole and stabilizes the system. This has been implemented in CHARMM apparently.
So the question is: is this implemented in LAMMPS as well, (I haven’t found anything about it though), or is there an alternative way of doing this?

All usable features in LAMMPS are documented. This is a requirement for contributing to the distribution. So if it is not mentioned in the manual it is not available. You would need to modify the source to implement it.

From my memory, the Drude pairs in LAMMPS can be connected with any (suitable) bond_style available in LAMMPS. So you can look through that list and see if any suits (not bond_style quartic, funny enough) and use those to link the Drude pairs. Maybe bond_style table will suit your use.

Bond style table works as a charm! Didn’t know about that one. (It is not so easy to know all the features of LAMMPS, thankfully there is always help to get in this forum :)) Thanks a lot!

1 Like