The example input scripts included in the LAMMPS distribution and highlighted in Section_example also show how to @@ -2415,14 +2417,91 @@ to specify this command explicitly.
In polarizable force fields the charge distributions in molecules and +materials respond to their electrostatic environements. Polarizable +systems can be simulated in LAMMPS using three methods: +
+The fluctuating charge method calculates instantaneous charges on +interacting atoms based on the electronegativity equalization +principle. It is implemented in the fix qeq which is +available in several variants. It is a relatively efficient technique +since no additional particles are introduced. This method allows for +charge transfer between molecules or atom groups. However, because the +charges are located at the interaction sites, off-plane components of +polarization cannot be represented in planar molecules or atom groups. +
+The two other methods share the same basic idea: polarizable atoms are +split into one core atom and one satellite particle (called shell or +Drude particle) attached to it by a harmonic spring. Both atoms bear +a charge and they represent collectively an induced electric dipole. +These techniques are computationally more expensive than the QEq +method because of additional particles and bonds. These two +charge-on-spring methods differ in certain features, with the +core-shell model being normally used for ionic/crystalline materials, +whereas the so-called Drude model is normally used for molecular +systems and fluid states. +
+The core-shell model is applicable to crystalline materials where the +high symmetry around each site leads to stable trajectories of the +core-shell pairs. However, bonded atoms in molecules can be so close +that a core would interact too strongly or even capture the Drude +particle of a neighbor. The Drude dipole model is relatively more +complex in order to remediate this and other issues. Specifically, the +Drude model includes specific thermostating of the core-Drude pairs +and short-range damping of the induced dipoles. +
+The three polarization methods can be implemented through a +self-consistent calculation of charges or induced dipoles at each +timestep. In the fluctuating charge scheme this is done by the matrix +inversion method in fix qeq/point, but for core-shell +or Drude-dipoles the relaxed-dipoles technique would require an slow +iterative procedure. These self-consistent solutions yield accurate +trajectories since the additional degrees of freedom representing +polarization are massless. An alternative is to attribute a mass to +the additional degrees of freedom and perform time integration using +an extended Lagrangian technique. For the fluctuating charge scheme +this is done by fix qeq/dynamic, and for the +charge-on-spring models by the methods outlined in the next two +sections. The assignment of masses to the additional degrees of +freedom can lead to unphysical trajectories if care is not exerted in +choosing the parameters of the poarizable models and the simulation +conditions. +
+In the core-shell model the vibration of the shells is kept faster +than the ionic vibrations to mimic the fast response of the +polarizable electrons. But in molecular systems thermalizing the +core-Drude pairs at temperatures comparable to the rest of the +simulation leads to several problems (kinetic energy transfer, too +short a timestep, etc.) In order to avoid these problems the relative +motion of the Drude particles with respect to their cores is kept +"cold" so the vibration of the core-Drude pairs is very slow, +approaching the self-consistent regime. In both models the +temperature is regulated using the velocities of the center of mass of +core+shell (or Drude) pairs, but in the Drude model the actual +relative core-Drude particle motion is thermostated separately as +well. +
+The adiabatic core-shell model by Mitchell and Finchham is a simple method for adding polarizability to a system. In order to mimic the electron shell of -an ion, a ghost atom is attached to it. This way the ions are split -into a core and a shell where the latter is meant to react to the -electrostatic environment inducing polarizability. +an ion, a satellite particle is attached to it. This way the ions are +split into a core and a shell where the latter is meant to react to +the electrostatic environment inducing polarizability.
Technically, shells are attached to the cores by a spring force f = k*r where k is a parametrized spring constant and r is the distance @@ -2600,6 +2679,70 @@ fix ave_chunk all ave/time 10 1 10 c_cstherm file chunk.dump mode vector
The thermalized Drude model, similarly to the core-shell +model, representes induced dipoles by a pair of charges (the core atom +and the Drude particle) connected by a harmonic spring. The Drude +model has a number of features aimed at its use in molecular systems +(Lamoureux and Roux): +
+A detailed tutorial covering the usage of Drude induced dipoles in +LAMMPS is available here. +
+As with the core-shell model, the cores and Drude particles should +appear in the data file as standard atoms. The same holds for the +springs between them, which are described by standard harmonic bonds. +The nature of the atoms (core, Drude particle or non-polarizable) is +specified via the fix drude command. The special +list of neighbors is automatically refactored to account for the +equivalence of core and Drude particles as regards special 1-2 to 1-4 +screening. It may be necessary to use the extra keyword of the +special_bonds command. If using fix +shake, make sure no Drude particle is in this fix +group. +
+There are two ways to thermostat the Drude particles at a low +temperature: use either fix langevin/drude +for a Langevin thermostat, or fix +drude/transform for a Nose-Hoover +thermostat. The former requires use of the command comm_modify vel +yes. The latter requires two separate integration +fixes like nvt or npt. The correct temperatures of the reduced +degrees of freedom can be calculated using the compute +temp/drude. This requires also to use the +command comm_modify vel yes. +
+Short-range damping of the induced dipole interactions can be achieved +using Thole functions through the the pair style +thole in pair_style hybrid/overlay +with a Coulomb pair style. It may be useful to use coul/long/cs or +similar from the CORESHELL package if the core and Drude particle come +too close, which can cause numerical issues. +
+(Mitchell and Finchham) Mitchell, Finchham, J Phys Condensed Matter, 5, 1031-1038 (1993).
+ + +(Lamoureux and Roux) G. Lamoureux, B. Roux, J. Chem. Phys 119, 3025 (2003) +