Merge remote-tracking branch 'upstream/develop' into numdiff-stress

doc/src/Developer_org.rst

@@ -225,7 +225,7 @@ follows:
   commands in an input script.
 
 - The Force class computes various forces between atoms.  The Pair
-  parent class is for non-bonded or pair-wise forces, which in LAMMPS
+  parent class is for non-bonded or pairwise forces, which in LAMMPS
   also includes many-body forces such as the Tersoff 3-body potential if
   those are computed by walking pairwise neighbor lists.  The Bond,
   Angle, Dihedral, Improper parent classes are styles for bonded

doc/src/Errors_warnings.rst

@@ -416,7 +416,7 @@ This will most likely cause errors in kinetic fluctuations.
    not defined for the specified atom style.
 
 *Molecule has bond topology but no special bond settings*
-   This means the bonded atoms will not be excluded in pair-wise
+   This means the bonded atoms will not be excluded in pairwise
    interactions.
 
 *Molecule template for create_atoms has multiple molecules*

doc/src/Modify_style.rst

@@ -305,6 +305,40 @@ you are uncertain, please ask.
   FILE pointers and only be done on MPI rank 0.  Use the :cpp:func:`utils::logmesg`
   convenience function where possible.
 
+- Usage of C++11 `virtual`, `override`, `final` keywords: Please follow the
+  `C++ Core Guideline C.128 <https://isocpp.github.io/CppCoreGuidelines/CppCoreGuidelines#Rh-override>`_.
+  That means, you should only use `virtual` to declare a new virtual
+  function, `override` to indicate you are overriding an existing virtual
+  function, and `final` to prevent any further overriding.
+
+- Trivial destructors: Prefer not writing destructors when they are empty and `default`.
+
+  .. code-block:: c++
+
+     // don't write destructors for A or B like this
+     class A : protected Pointers {
+      public:
+        A();
+        ~A() override {}
+     };
+
+     class B : protected Pointers {
+      public:
+        B();
+        ~B() override = default;
+     };
+
+     // instead, let the compiler create the implicit default destructor by not writing it
+     class A : protected Pointers {
+      public:
+        A();
+     };
+
+     class B : protected Pointers {
+      public:
+        B();
+     };
+
 - Header files, especially those defining a "style", should only use
   the absolute minimum number of include files and **must not** contain
   any ``using`` statements. Typically that would be only the header for

doc/src/Run_output.rst

@@ -106,7 +106,7 @@ individual ranks. Here is an example output for this section:
 ----------
 
 The third section above lists the number of owned atoms (Nlocal),
-ghost atoms (Nghost), and pair-wise neighbors stored per processor.
+ghost atoms (Nghost), and pairwise neighbors stored per processor.
 The max and min values give the spread of these values across
 processors with a 10-bin histogram showing the distribution. The total
 number of histogram counts is equal to the number of processors.
@@ -114,7 +114,7 @@ number of histogram counts is equal to the number of processors.
 ----------
 
 The last section gives aggregate statistics (across all processors)
-for pair-wise neighbors and special neighbors that LAMMPS keeps track
+for pairwise neighbors and special neighbors that LAMMPS keeps track
 of (see the :doc:`special_bonds <special_bonds>` command).  The number
 of times neighbor lists were rebuilt is tallied, as is the number of
 potentially *dangerous* rebuilds.  If atom movement triggered neighbor

doc/src/Speed_kokkos.rst

@@ -214,7 +214,7 @@ threads/task as Nt. The product of these two values should be N, i.e.
    The default for the :doc:`package kokkos <package>` command when
    running on KNL is to use "half" neighbor lists and set the Newton flag
    to "on" for both pairwise and bonded interactions. This will typically
-   be best for many-body potentials. For simpler pair-wise potentials, it
+   be best for many-body potentials. For simpler pairwise potentials, it
    may be faster to use a "full" neighbor list with Newton flag to "off".
    Use the "-pk kokkos" :doc:`command-line switch <Run_options>` to change
    the default :doc:`package kokkos <package>` options. See its page for

doc/src/balance.rst

@@ -383,7 +383,7 @@ multiple groups, its weight is the product of the weight factors.
 
 This weight style is useful in combination with pair style
 :doc:`hybrid <pair_hybrid>`, e.g. when combining a more costly many-body
-potential with a fast pair-wise potential.  It is also useful when
+potential with a fast pairwise potential.  It is also useful when
 using :doc:`run_style respa <run_style>` where some portions of the
 system have many bonded interactions and others none.  It assumes that
 the computational cost for each group remains constant over time.

doc/src/compute_pressure_cylinder.rst

@@ -61,7 +61,7 @@ Restrictions
 This compute currently calculates the pressure tensor contributions
 for pair styles only (i.e. no bond, angle, dihedral, etc. contributions
 and in the presence of bonded interactions, the result will be incorrect
-due to exclusions for special bonds)  and requires pair-wise force
+due to exclusions for special bonds)  and requires pairwise force
 calculations not available for most many-body pair styles. K-space
 calculations are also excluded. Note that this pressure compute outputs
 the configurational terms only; the kinetic contribution is not included

doc/src/dump_modify.rst

@@ -17,13 +17,14 @@ Syntax
 * one or more keyword/value pairs may be appended
 
 * these keywords apply to various dump styles
-* keyword = *append* or *at* or *buffer* or *delay* or *element* or *every* or *every/time* or *fileper* or *first* or *flush* or *format* or *header* or *image* or *label* or *maxfiles* or *nfile* or *pad* or *pbc* or *precision* or *region* or *refresh* or *scale* or *sfactor* or *sort* or *tfactor* or *thermo* or *thresh* or *time* or *units* or *unwrap*
+* keyword = *append* or *at* or *balance* or *buffer* or *delay* or *element* or *every* or *every/time* or *fileper* or *first* or *flush* or *format* or *header* or *image* or *label* or *maxfiles* or *nfile* or *pad* or *pbc* or *precision* or *region* or *refresh* or *scale* or *sfactor* or *sort* or *tfactor* or *thermo* or *thresh* or *time* or *units* or *unwrap*
 
   .. parsed-literal::
 
        *append* arg = *yes* or *no*
        *at* arg = N
          N = index of frame written upon first dump
+       *balance* arg = *yes* or *no*
        *buffer* arg = *yes* or *no*
        *delay* arg = Dstep
          Dstep = delay output until this timestep
@@ -667,6 +668,14 @@ keywords are set to non-default values (i.e. the number of dump file
 pieces is not equal to the number of procs), then sorting cannot be
 performed.
 
+In a parallel run, the per-processor dump file pieces can have
+significant imbalance in number of lines of per-atom info. The *balance*
+keyword determines whether the number of lines in each processor
+snapshot are balanced to be nearly the same. A balance value of *no*
+means no balancing will be done, while *yes* means balancing will be
+performed. This balancing preserves dump sorting order. For a serial
+run, this option is ignored since the output is already balanced.
+
 .. note::
 
    Unless it is required by the dump style, sorting dump file
@@ -832,6 +841,7 @@ Default
 The option defaults are
 
 * append = no
+* balance = no
 * buffer = yes for dump styles *atom*, *custom*, *loca*, and *xyz*
 * element = "C" for every atom type
 * every = whatever it was set to via the :doc:`dump <dump>` command

doc/src/fix_charge_regulation.rst

@@ -20,13 +20,13 @@ Syntax
 
   .. parsed-literal::
 
-     keyword = *pH*, *pKa*, *pKb*, *pIp*, *pIm*, *pKs*, *acid_type*, *base_type*, *lunit_nm*, *temp*, *tempfixid*, *nevery*, *nmc*, *xrd*, *seed*, *tag*, *group*, *onlysalt*, *pmcmoves*
+     keyword = *pH*, *pKa*, *pKb*, *pIp*, *pIm*, *pKs*, *acid_type*, *base_type*, *lunit_nm*, *temp*, *tempfixid*, *nevery*, *nmc*, *rxd*, *seed*, *tag*, *group*, *onlysalt*, *pmcmoves*
      *pH* value = pH of the solution (can be specified as an equal-style variable)
-     *pKa* value = acid dissociation constant
-     *pKb* value = base dissociation constant
-     *pIp* value = chemical potential of free cations
-     *pIm* value = chemical potential of free anions
-     *pKs* value = solution self-dissociation constant
+     *pKa* value = acid dissociation constant (in the -log10 representation)
+     *pKb* value = base dissociation constant (in the -log10 representation)
+     *pIp* value = activity (effective concentration) of free cations (in the -log10 representation)
+     *pIm* value = activity (effective concentration) of free anions (in the -log10 representation)
+     *pKs* value = solvent self-dissociation constant (in the -log10 representation)
      *acid_type* = atom type of acid groups
      *base_type*  = atom type of base groups
      *lunit_nm* value = unit length used by LAMMPS (# in the units of nanometers)
@@ -34,7 +34,7 @@ Syntax
      *tempfixid* value = fix ID of temperature thermostat
      *nevery* value = invoke this fix every nevery steps
      *nmc* value = number of charge regulation MC moves to attempt every nevery steps
-     *xrd* value = cutoff distance for acid/base reaction
+     *rxd* value = cutoff distance for acid/base reaction
      *seed* value = random # seed (positive integer)
      *tag* value = yes or no (yes: The code assign unique tags to inserted ions; no: The tag of all inserted ions is "0")
      *group* value = group-ID, inserted ions are assigned to group group-ID. Can be used multiple times to assign inserted ions to multiple groups.
@@ -47,7 +47,7 @@ Examples
 """"""""
 .. code-block:: LAMMPS
 
-    fix chareg all charge/regulation 1 2 acid_type 3 base_type 4 pKa 5 pKb 7 lb 1.0 nevery 200 nexchange 200 seed 123 tempfixid fT
+    fix chareg all charge/regulation 1 2 acid_type 3 base_type 4 pKa 5.0 pKb 6.0 pH 7.0 pIp 3.0 pIm 3.0 nevery 200 nmc 200 seed 123 tempfixid fT
 
     fix chareg all charge/regulation 1 2 pIp 3 pIm 3 onlysalt yes 2 -1 seed 123 tag yes temp 1.0
 
@@ -92,7 +92,11 @@ where the fix attempts to charge :math:`\mathrm{A}` (discharge
 :math:`\mathrm{A}^-`) to :math:`\mathrm{A}^-` (:math:`\mathrm{A}`) and
 insert (delete) a proton (atom type 2). Besides, the fix implements
 self-ionization reaction of water :math:`\emptyset \rightleftharpoons
-\mathrm{H}^++\mathrm{OH}^-`.  However, this approach is highly
+\mathrm{H}^++\mathrm{OH}^-`.
+
+
+
+However, this approach is highly
 inefficient at :math:`\mathrm{pH} \approx 7` when the concentration of
 both protons and hydroxyl ions is low, resulting in a relatively low
 acceptance rate of MC moves.
@@ -102,24 +106,31 @@ reactions, which can be easily achieved via
 
 .. code-block:: LAMMPS
 
-    fix acid_reaction all charge/regulation 4 5 acid_type 1 pH 7.0 pKa 5.0 pIp 2.0 pIm 2.0
+    fix acid_reaction2 all charge/regulation 4 5 acid_type 1 pH 7.0 pKa 5.0 pIp 2.0 pIm 2.0
 
-where particles of atom type 4 and 5 are the salt cations and anions,
-both at chemical potential pI=2.0, see :ref:`(Curk1) <Curk1>` and
+where particles of atom type 4 and 5 are the salt cations and anions, both at activity (effective concentration) of :math:`10^{-2}` mol/l, see :ref:`(Curk1) <Curk1>` and
 :ref:`(Landsgesell) <Landsgesell>` for more details.
 
-
- Similarly, we could have simultaneously added a base ionization reaction
- (:math:`\mathrm{B} \rightleftharpoons \mathrm{B}^++\mathrm{OH}^-`)
+We could have simultaneously added a base ionization reaction (:math:`\mathrm{B} \rightleftharpoons \mathrm{B}^++\mathrm{OH}^-`)
 
 .. code-block:: LAMMPS
 
-    fix base_reaction all charge/regulation 2 3 base_type 6 pH 7.0 pKb 6.0 pIp 7.0 pIm 7.0
+    fix acid_base_reaction all charge/regulation 2 3 acid_type 1 base_type 6 pH 7.0 pKa 5.0 pKb 6.0 pIp 7.0 pIm 7.0
 
 where the fix will attempt to charge :math:`\mathrm{B}` (discharge
 :math:`\mathrm{B}^+`) to :math:`\mathrm{B}^+` (:math:`\mathrm{B}`) and
-insert (delete) a hydroxyl ion :math:`\mathrm{OH}^-` of atom type 3.  If
-neither the acid or the base type is specified, for example,
+insert (delete) a hydroxyl ion :math:`\mathrm{OH}^-` of atom type 3.
+
+
+Dissociated ions and salt ions can be combined into a single particle type, which reduces the number of necessary MC moves and increases sampling performance, see :ref:`(Curk1) <Curk1>`. The :math:`\mathrm{H}^+` and monovalent salt cation (:math:`\mathrm{S}^+`) are combined into a single particle type, :math:`\mathrm{X}^+ = \{\mathrm{H}^+, \mathrm{S}^+\}`. In this case "pIp" refers to the effective concentration of the combined cation type :math:`\mathrm{X}^+` and its value is determined by :math:`10^{-\mathrm{pIp}} = 10^{-\mathrm{pH}} + 10^{-\mathrm{pSp}}`, where :math:`10^{-\mathrm{pSp}}` is the effective concentration of salt cations. For example, at pH=7 and pSp=6 we would find pIp~5.958 and the command that performs reactions with combined ions could read,
+
+.. code-block:: LAMMPS
+
+    fix acid_reaction_combined all charge/regulation 2 3 acid_type 1 pH 7.0 pKa 5.0 pIp 5.958 pIm 5.958
+
+
+
+If neither the acid or the base type is specified, for example,
 
 .. code-block:: LAMMPS
 
@@ -127,11 +138,11 @@ neither the acid or the base type is specified, for example,
 
 the fix simply inserts or deletes an ion pair of a free cation (atom
 type 4) and a free anion (atom type 5) as done in a conventional
-grand-canonical MC simulation.
+grand-canonical MC simulation. Multivalent ions can be inserted (deleted) by using the *onlysalt* keyword.
 
 
 The fix is compatible with LAMMPS sub-packages such as *molecule* or
-*rigid*. That said, the acid and base particles can be part of larger
+*rigid*. The acid and base particles can be part of larger
 molecules or rigid bodies. Free ions that are inserted to or deleted
 from the system must be defined as single particles (no bonded
 interactions allowed) and cannot be part of larger molecules or rigid
@@ -153,14 +164,14 @@ Langevin thermostat:
     fix fT all langevin 1.0 1.0 1.0 123
     fix_modify fT temp dtemp
 
-The chemical potential units (e.g. pH) are in the standard log10
+The units of  pH, pKa, pKb, pIp, pIm are considered to be in the standard -log10
 representation assuming reference concentration :math:`\rho_0 =
-\mathrm{mol}/\mathrm{l}`.  Therefore, to perform the internal unit
-conversion, the length (in nanometers) of the LAMMPS unit length must be
-specified via *lunit_nm* (default is set to the Bjerrum length in water
-at room temperature *lunit_nm* = 0.71nm). For example, in the dilute
-ideal solution limit, the concentration of free ions will be
-:math:`c_\mathrm{I} = 10^{-\mathrm{pIp}}\mathrm{mol}/\mathrm{l}`.
+\mathrm{mol}/\mathrm{l}`.  For example, in the dilute
+ideal solution limit, the concentration of free cations will be
+:math:`c_\mathrm{I} = 10^{-\mathrm{pIp}}\mathrm{mol}/\mathrm{l}`. To perform the internal unit
+conversion, the the value of the LAMMPS unit length must be
+specified in nanometers via *lunit_nm*. The default value is set to the Bjerrum length in water
+at room temperature (0.71 nm), *lunit_nm* = 0.71.
 
 The temperature used in MC acceptance probability is set by *temp*. This
 temperature should be the same as the temperature set by the molecular
@@ -171,10 +182,10 @@ thermostat fix-ID is *fT*. The inserted particles attain a random
 velocity corresponding to the specified temperature. Using *tempfixid*
 overrides any fixed temperature set by *temp*.
 
-The *xrd* keyword can be used to restrict the inserted/deleted
+The *rxd* keyword can be used to restrict the inserted/deleted
 counterions to a specific radial distance from an acid or base particle
 that is currently participating in a reaction. This can be used to
-simulate more realist reaction dynamics. If *xrd* = 0 or *xrd* > *L* /
+simulate more realist reaction dynamics. If *rxd* = 0 or *rxd* > *L* /
 2, where *L* is the smallest box dimension, the radial restriction is
 automatically turned off and free ion can be inserted or deleted
 anywhere in the simulation box.
@@ -258,18 +269,18 @@ Default
 
 pH = 7.0; pKa = 100.0; pKb = 100.0; pIp = 5.0; pIm = 5.0; pKs = 14.0;
 acid_type = -1; base_type = -1; lunit_nm = 0.71; temp = 1.0; nevery =
-100; nmc = 100; xrd = 0; seed = 0; tag = no; onlysalt = no, pmcmoves =
+100; nmc = 100; rxd = 0; seed = 0; tag = no; onlysalt = no, pmcmoves =
 [1/3, 1/3, 1/3], group-ID = all
 
 ----------
 
 .. _Curk1:
 
-**(Curk1)** T. Curk, J. Yuan, and E. Luijten, "Coarse-grained simulation of charge regulation using LAMMPS", preprint (2021).
+**(Curk1)** T. Curk, J. Yuan, and E. Luijten, "Accelerated simulation method for charge regulation effects", JCP 156 (2022).
 
 .. _Curk2:
 
-**(Curk2)** T. Curk and E. Luijten, "Charge-regulation effects in nanoparticle self-assembly", PRL (2021)
+**(Curk2)** T. Curk and E. Luijten, "Charge-regulation effects in nanoparticle self-assembly", PRL 126 (2021)
 
 .. _Landsgesell:
 

doc/src/package.rst

@@ -460,7 +460,7 @@ using *neigh/thread* *on*, a full neighbor list must also be used. Using
 is turned on by default only when there are 16K atoms or less owned by
 an MPI rank and when using a full neighbor list. Not all KOKKOS-enabled
 potentials support this keyword yet, and only thread over atoms. Many
-simple pair-wise potentials such as Lennard-Jones do support threading
+simple pairwise potentials such as Lennard-Jones do support threading
 over both atoms and neighbors.
 
 The *newton* keyword sets the Newton flags for pairwise and bonded

doc/src/pair_charmm.rst

@@ -119,7 +119,7 @@ name are the older, original LAMMPS implementations.  They compute the
 LJ and Coulombic interactions with an energy switching function (esw,
 shown in the formula below as S(r)), which ramps the energy smoothly
 to zero between the inner and outer cutoff.  This can cause
-irregularities in pair-wise forces (due to the discontinuous second
+irregularities in pairwise forces (due to the discontinuous second
 derivative of energy at the boundaries of the switching region), which
 in some cases can result in detectable artifacts in an MD simulation.
 

doc/src/pair_dpd.rst

@@ -50,7 +50,7 @@ Style *dpd* computes a force field for dissipative particle dynamics
 
 Style *dpd/tstat* invokes a DPD thermostat on pairwise interactions,
 which is equivalent to the non-conservative portion of the DPD force
-field.  This pair-wise thermostat can be used in conjunction with any
+field.  This pairwise thermostat can be used in conjunction with any
 :doc:`pair style <pair_style>`, and in leiu of per-particle thermostats
 like :doc:`fix langevin <fix_langevin>` or ensemble thermostats like
 Nose Hoover as implemented by :doc:`fix nvt <fix_nh>`.  To use

doc/src/pair_python.rst

@@ -164,7 +164,7 @@ Following the *LJCutMelt* example, here are the two functions:
 .. note::
 
    The evaluation of scripted python code will slow down the
-   computation pair-wise interactions quite significantly. However, this
+   computation pairwise interactions quite significantly. However, this
    can be largely worked around through using the python pair style not
    for the actual simulation, but to generate tabulated potentials on the
    fly using the :doc:`pair_write <pair_write>` command. Please see below

doc/src/pair_style.rst

@@ -154,10 +154,10 @@ accelerated styles exist.
 * :doc:`coul/wolf/cs <pair_cs>` - Coulomb via Wolf potential with core/shell adjustments
 * :doc:`dpd <pair_dpd>` - dissipative particle dynamics (DPD)
 * :doc:`dpd/ext <pair_dpd_ext>` - generalized force field for DPD
-* :doc:`dpd/ext/tstat <pair_dpd_ext>` - pair-wise DPD thermostatting  with generalized force field
+* :doc:`dpd/ext/tstat <pair_dpd_ext>` - pairwise DPD thermostatting  with generalized force field
 * :doc:`dpd/fdt <pair_dpd_fdt>` - DPD for constant temperature and pressure
 * :doc:`dpd/fdt/energy <pair_dpd_fdt>` - DPD for constant energy and enthalpy
-* :doc:`dpd/tstat <pair_dpd>` - pair-wise DPD thermostatting
+* :doc:`dpd/tstat <pair_dpd>` - pairwise DPD thermostatting
 * :doc:`dsmc <pair_dsmc>` - Direct Simulation Monte Carlo (DSMC)
 * :doc:`e3b <pair_e3b>` - Explicit-three body (E3B) water model
 * :doc:`drip <pair_drip>` - Dihedral-angle-corrected registry-dependent interlayer potential (DRIP)

doc/src/pair_sw.rst

@@ -202,7 +202,7 @@ elements are the same.  Thus the two-body parameters for Si
 interacting with C, comes from the SiCC entry.  The three-body
 parameters can in principle be specific to the three elements of the
 configuration. In the literature, however, the three-body parameters
-are usually defined by simple formulas involving two sets of pair-wise
+are usually defined by simple formulas involving two sets of pairwise
 parameters, corresponding to the ij and ik pairs, where i is the
 center atom. The user must ensure that the correct combining rule is
 used to calculate the values of the three-body parameters for

doc/src/pair_write.rst

@@ -76,8 +76,10 @@ must be set before this command can be invoked.
 Due to how the pairwise force is computed, an inner value > 0.0 must
 be specified even if the potential has a finite value at r = 0.0.
 
-For EAM potentials, the pair_write command only tabulates the
-pairwise portion of the potential, not the embedding portion.
+The *pair_write* command can only be used for pairwise additive
+interactions for which a `Pair::single()` function can be and has
+been implemented.  This excludes for example manybody potentials
+or TIP4P coulomb styles.
 
 Related commands
 """"""""""""""""

doc/src/run_style.rst

@@ -89,7 +89,7 @@ in its 3d FFTs.  In this scenario, splitting your P total processors
 into 2 subsets of processors, P1 in the first partition and P2 in the
 second partition, can enable your simulation to run faster.  This is
 because the long-range forces in PPPM can be calculated at the same
-time as pair-wise and bonded forces are being calculated, and the FFTs
+time as pairwise and bonded forces are being calculated, and the FFTs
 can actually speed up when running on fewer processors.
 
 To use this style, you must define 2 partitions where P1 is a multiple