<li><pclass="first">style = <em>none</em> or <em>ewald</em> or <em>ewald/disp</em> or <em>ewald/omp</em> or <em>pppm</em> or <em>pppm/cg</em> or <em>pppm/disp</em> or <em>pppm/tip4p</em> or <em>pppm/stagger</em> or <em>pppm/disp/tip4p</em> or <em>pppm/gpu</em> or <em>pppm/kk</em> or <em>pppm/omp</em> or <em>pppm/cg/omp</em> or <em>pppm/tip4p/omp</em> or <em>msm</em> or <em>msm/cg</em> or <em>msm/omp</em> or <em>msm/cg/omp</em></p>
<preclass="literal-block">
<em>none</em> value = none
<em>ewald</em> value = accuracy
accuracy = desired relative error in forces
<em>ewald/disp</em> value = accuracy
accuracy = desired relative error in forces
<em>ewald/omp</em> value = accuracy
accuracy = desired relative error in forces
<em>pppm</em> value = accuracy
accuracy = desired relative error in forces
<em>pppm/cg</em> value = accuracy (smallq)
accuracy = desired relative error in forces
smallq = cutoff for charges to be considered (optional) (charge units)
<em>pppm/disp</em> value = accuracy
accuracy = desired relative error in forces
<em>pppm/tip4p</em> value = accuracy
accuracy = desired relative error in forces
<em>pppm/disp/tip4p</em> value = accuracy
accuracy = desired relative error in forces
<em>pppm/gpu</em> value = accuracy
accuracy = desired relative error in forces
<em>pppm/kk</em> value = accuracy
accuracy = desired relative error in forces
<em>pppm/omp</em> value = accuracy
accuracy = desired relative error in forces
<em>pppm/cg/omp</em> value = accuracy
accuracy = desired relative error in forces
<em>pppm/tip4p/omp</em> value = accuracy
accuracy = desired relative error in forces
<em>pppm/stagger</em> value = accuracy
accuracy = desired relative error in forces
<em>msm</em> value = accuracy
accuracy = desired relative error in forces
<em>msm/cg</em> value = accuracy (smallq)
accuracy = desired relative error in forces
smallq = cutoff for charges to be considered (optional) (charge units)
<em>msm/omp</em> value = accuracy
accuracy = desired relative error in forces
<em>msm/cg/omp</em> value = accuracy (smallq)
accuracy = desired relative error in forces
smallq = cutoff for charges to be considered (optional) (charge units)
</pre>
</li>
</ul>
</div>
<divclass="section"id="examples">
<h2>Examples</h2>
<preclass="literal-block">
kspace_style pppm 1.0e-4
kspace_style pppm/cg 1.0e-5 1.0e-6
kspace style msm 1.0e-4
kspace_style none
</pre>
</div>
<divclass="section"id="description">
<h2>Description</h2>
<p>Define a long-range solver for LAMMPS to use each timestep to compute
long-range Coulombic interactions or long-range 1/r^6 interactions.
Most of the long-range solvers perform their computation in K-space,
hence the name of this command.</p>
<p>When such a solver is used in conjunction with an appropriate pair
style, the cutoff for Coulombic or 1/r^N interactions is effectively
infinite. If the Coulombic case, this means each charge in the system
interacts with charges in an infinite array of periodic images of the
simulation domain.</p>
<p>Note that using a long-range solver requires use of a matching <aclass="reference internal"href="pair_style.html"><spanclass="doc">pair style</span></a> to perform consistent short-range pairwise
calculations. This means that the name of the pair style contains a
matching keyword to the name of the KSpace style, as in this table:</p>
<tableborder="1"class="docutils">
<colgroup>
<colwidth="49%"/>
<colwidth="51%"/>
</colgroup>
<tbodyvalign="top">
<trclass="row-odd"><td>Pair style</td>
<td>KSpace style</td>
</tr>
<trclass="row-even"><td>coul/long</td>
<td>ewald or pppm</td>
</tr>
<trclass="row-odd"><td>coul/msm</td>
<td>msm</td>
</tr>
<trclass="row-even"><td>lj/long or buck/long</td>
<td>disp (for dispersion)</td>
</tr>
<trclass="row-odd"><td>tip4p/long</td>
<td>tip4p</td>
</tr>
</tbody>
</table>
<hrclass="docutils"/>
<p>The <em>ewald</em> style performs a standard Ewald summation as described in
any solid-state physics text.</p>
<p>The <em>ewald/disp</em> style adds a long-range dispersion sum option for
1/r^6 potentials and is useful for simulation of interfaces
<aclass="reference internal"href="pair_lj_long.html#veld"><spanclass="std std-ref">(Veld)</span></a>. It also performs standard Coulombic Ewald summations,
but in a more efficient manner than the <em>ewald</em> style. The 1/r^6
capability means that Lennard-Jones or Buckingham potentials can be
used without a cutoff, i.e. they become full long-range potentials.
The <em>ewald/disp</em> style can also be used with point-dipoles
<aclass="reference internal"href="pair_dipole.html#toukmaji"><spanclass="std std-ref">(Toukmaji)</span></a> and is currently the only kspace solver in
LAMMPS with this capability.</p>
<hrclass="docutils"/>
<p>The <em>pppm</em> style invokes a particle-particle particle-mesh solver
<aclass="reference internal"href="#hockney"><spanclass="std std-ref">(Hockney)</span></a> which maps atom charge to a 3d mesh, uses 3d FFTs
to solve Poisson’s equation on the mesh, then interpolates electric
fields on the mesh points back to the atoms. It is closely related to
the particle-mesh Ewald technique (PME) <aclass="reference internal"href="#darden"><spanclass="std std-ref">(Darden)</span></a> used in
AMBER and CHARMM. The cost of traditional Ewald summation scales as
N^(3/2) where N is the number of atoms in the system. The PPPM solver
scales as Nlog(N) due to the FFTs, so it is almost always a faster
<pclass="last">All of the PPPM styles can be used with single-precision FFTs by
using the compiler switch -DFFT_SINGLE for the FFT_INC setting in your
lo-level Makefile. This setting also changes some of the PPPM
operations (e.g. mapping charge to mesh and interpolating electric
fields to particles) to be performed in single precision. This option
can speed-up long-range calulations, particularly in parallel or on
GPUs. The use of the -DFFT_SINGLE flag is discussed in <aclass="reference internal"href="Section_start.html#start-2-4"><spanclass="std std-ref">this section</span></a> of the manual. MSM does not
currently support the -DFFT_SINGLE compiler switch.</p>
</div>
<hrclass="docutils"/>
<p>The <em>msm</em> style invokes a multi-level summation method MSM solver,
<aclass="reference internal"href="#hardy"><spanclass="std std-ref">(Hardy)</span></a> or <aclass="reference internal"href="#hardy2"><spanclass="std std-ref">(Hardy2)</span></a>, which maps atom charge to a 3d
mesh, and uses a multi-level hierarchy of coarser and coarser meshes
on which direct coulomb solves are done. This method does not use
FFTs and scales as N. It may therefore be faster than the other
K-space solvers for relatively large problems when running on large
core counts. MSM can also be used for non-periodic boundary conditions and
for mixed periodic and non-periodic boundaries.</p>
<p>MSM is most competitive versus Ewald and PPPM when only relatively
low accuracy forces, about 1e-4 relative error or less accurate,
are needed. Note that use of a larger coulomb cutoff (i.e. 15
angstroms instead of 10 angstroms) provides better MSM accuracy for
both the real space and grid computed forces.</p>
<p>Currently calculation of the full pressure tensor in MSM is expensive.
Using the <aclass="reference internal"href="kspace_modify.html"><spanclass="doc">kspace_modify</span></a><em>pressure/scalar yes</em>
command provides a less expensive way to compute the scalar pressure
(Pxx + Pyy + Pzz)/3.0. The scalar pressure can be used, for example,
to run an isotropic barostat. If the full pressure tensor is needed,
then calculating the pressure at every timestep or using a fixed
pressure simulation with MSM will cause the code to run slower.</p>
<hrclass="docutils"/>
<p>The specified <em>accuracy</em> determines the relative RMS error in per-atom
forces calculated by the long-range solver. It is set as a
dimensionless number, relative to the force that two unit point
charges (e.g. 2 monovalent ions) exert on each other at a distance of
1 Angstrom. This reference value was chosen as representative of the
magnitude of electrostatic forces in atomic systems. Thus an accuracy
value of 1.0e-4 means that the RMS error will be a factor of 10000
smaller than the reference force.</p>
<p>The accuracy setting is used in conjunction with the pairwise cutoff
to determine the number of K-space vectors for style <em>ewald</em> or the
grid size for style <em>pppm</em> or <em>msm</em>.</p>
<p>Note that style <em>pppm</em> only computes the grid size at the beginning of
a simulation, so if the length or triclinic tilt of the simulation
cell increases dramatically during the course of the simulation, the
accuracy of the simulation may degrade. Likewise, if the
<aclass="reference internal"href="kspace_modify.html"><spanclass="doc">kspace_modify slab</span></a> option is used with
shrink-wrap boundaries in the z-dimension, and the box size changes
dramatically in z. For example, for a triclinic system with all three
tilt factors set to the maximum limit, the PPPM grid should be
increased roughly by a factor of 1.5 in the y direction and 2.0 in the
z direction as compared to the same system using a cubic orthogonal
simulation cell. One way to ensure the accuracy requirement is being
met is to run a short simulation at the maximum expected tilt or
length, note the required grid size, and then use the
<aclass="reference internal"href="kspace_modify.html"><spanclass="doc">kspace_modify</span></a><em>mesh</em> command to manually set the
PPPM grid size to this value.</p>
<p>RMS force errors in real space for <em>ewald</em> and <em>pppm</em> are estimated
using equation 18 of <aclass="reference internal"href="#kolafa"><spanclass="std std-ref">(Kolafa)</span></a>, which is also referenced as
equation 9 of <aclass="reference internal"href="#petersen"><spanclass="std std-ref">(Petersen)</span></a>. RMS force errors in K-space for
<em>ewald</em> are estimated using equation 11 of <aclass="reference internal"href="#petersen"><spanclass="std std-ref">(Petersen)</span></a>,
which is similar to equation 32 of <aclass="reference internal"href="#kolafa"><spanclass="std std-ref">(Kolafa)</span></a>. RMS force
errors in K-space for <em>pppm</em> are estimated using equation 38 of
<aclass="reference internal"href="#deserno"><spanclass="std std-ref">(Deserno)</span></a>. RMS force errors for <em>msm</em> are estimated
using ideas from chapter 3 of <aclass="reference internal"href="#hardy"><spanclass="std std-ref">(Hardy)</span></a>, with equation 3.197
of particular note. When using <em>msm</em> with non-periodic boundary
conditions, it is expected that the error estimation will be too
pessimistic. RMS force errors for dipoles when using <em>ewald/disp</em>
are estimated using equations 33 and 46 of <aclass="reference internal"href="pair_polymorphic.html#wang"><spanclass="std std-ref">(Wang)</span></a>.</p>
<p>See the <aclass="reference internal"href="kspace_modify.html"><spanclass="doc">kspace_modify</span></a> command for additional
options of the K-space solvers that can be set, including a <em>force</em>
option for setting an absoulte RMS error in forces, as opposed to a
relative RMS error.</p>
<hrclass="docutils"/>
<p>Styles with a <em>gpu</em>, <em>intel</em>, <em>kk</em>, <em>omp</em>, or <em>opt</em> suffix are
functionally the same as the corresponding style without the suffix.
They have been optimized to run faster, depending on your available
hardware, as discussed in <aclass="reference internal"href="Section_accelerate.html"><spanclass="doc">Section 5</span></a>
of the manual. The accelerated styles take the same arguments and
should produce the same results, except for round-off and precision
issues.</p>
<p>More specifically, the <em>pppm/gpu</em> style performs charge assignment and
force interpolation calculations on the GPU. These processes are
performed either in single or double precision, depending on whether
the -DFFT_SINGLE setting was specified in your lo-level Makefile, as
discussed above. The FFTs themselves are still calculated on the CPU.
If <em>pppm/gpu</em> is used with a GPU-enabled pair style, part of the PPPM
calculation can be performed concurrently on the GPU while other
calculations for non-bonded and bonded force calculation are performed
on the CPU.</p>
<p>The <em>pppm/kk</em> style also performs charge assignment and force
interpolation calculations on the GPU while the FFTs themselves are
calculated on the CPU in non-threaded mode.</p>
<p>These accelerated styles are part of the GPU, USER-INTEL,
KOKKOS, USER-OMP, and OPT packages respectively. They are only
enabled if LAMMPS was built with those packages. See the <aclass="reference internal"href="Section_start.html#start-3"><spanclass="std std-ref">Making LAMMPS</span></a> section for more info.</p>
<p>See <aclass="reference internal"href="Section_accelerate.html"><spanclass="doc">Section 5</span></a> of the manual for
more instructions on how to use the accelerated styles effectively.</p>
</div>
<divclass="section"id="restrictions">
<h2>Restrictions</h2>
<p>Note that the long-range electrostatic solvers in LAMMPS assume conducting
metal (tinfoil) boundary conditions for both charge and dipole
interactions. Vacuum boundary conditions are not currently supported.</p>
<p>The <em>ewald/disp</em>, <em>ewald</em>, <em>pppm</em>, and <em>msm</em> styles support
non-orthogonal (triclinic symmetry) simulation boxes. However,
triclinic simulation cells may not yet be supported by suffix versions
of these styles.</p>
<p>All of the kspace styles are part of the KSPACE package. They are
only enabled if LAMMPS was built with that package. See the <aclass="reference internal"href="Section_start.html#start-3"><spanclass="std std-ref">Making LAMMPS</span></a> section for more info. Note that
the KSPACE package is installed by default.</p>
<p>For MSM, a simulation must be 3d and one can use any combination of
periodic, non-periodic, or shrink-wrapped boundaries (specified using
the <aclass="reference internal"href="boundary.html"><spanclass="doc">boundary</span></a> command).</p>
<p>For Ewald and PPPM, a simulation must be 3d and periodic in all
dimensions. The only exception is if the slab option is set with
<aclass="reference internal"href="kspace_modify.html"><spanclass="doc">kspace_modify</span></a>, in which case the xy dimensions
must be periodic and the z dimension must be non-periodic.</p>
Built with <ahref="http://sphinx-doc.org/">Sphinx</a> using a <ahref="https://github.com/snide/sphinx_rtd_theme">theme</a> provided by <ahref="https://readthedocs.org">Read the Docs</a>.