Merge remote-tracking branch 'github/develop' into neigh-request-refactor

This commit is contained in:
Axel Kohlmeyer
2022-03-18 18:28:35 -04:00
177 changed files with 13666 additions and 199758 deletions

View File

@ -97,8 +97,6 @@ OPT.
* :doc:`latte <fix_latte>`
* :doc:`lb/fluid <fix_lb_fluid>`
* :doc:`lb/momentum <fix_lb_momentum>`
* :doc:`lb/pc <fix_lb_pc>`
* :doc:`lb/rigid/pc/sphere <fix_lb_rigid_pc_sphere>`
* :doc:`lb/viscous <fix_lb_viscous>`
* :doc:`lineforce <fix_lineforce>`
* :doc:`manifoldforce <fix_manifoldforce>`

View File

@ -8,11 +8,20 @@ without recompiling LAMMPS. The functionality for this and the
Plugins use the operating system's capability to load dynamic shared
object (DSO) files in a way similar shared libraries and then reference
specific functions in those DSOs. Any DSO file with plugins has to include
an initialization function with a specific name, "lammpsplugin_init", that
has to follow specific rules described below. When loading the DSO with
the "plugin" command, this function is looked up and called and will then
register the contained plugin(s) with LAMMPS.
specific functions in those DSOs. Any DSO file with plugins has to
include an initialization function with a specific name,
"lammpsplugin_init", that has to follow specific rules described below.
When loading the DSO with the "plugin" command, this function is looked
up and called and will then register the contained plugin(s) with
LAMMPS.
When the environment variable ``LAMMPS_PLUGIN_PATH`` is set, then LAMMPS
will search the directory (or directories) listed in this path for files
with names that end in ``plugin.so`` (e.g. ``helloplugin.so``) and will
try to load the contained plugins automatically at start-up. For
plugins that are loaded this way, the behavior of LAMMPS should be
identical to a binary where the corresponding code was compiled in
statically as a package.
From the programmer perspective this can work because of the object
oriented design of LAMMPS where all pair style commands are derived from
@ -65,19 +74,18 @@ Members of ``lammpsplugin_t``
* - handle
- Pointer to the open DSO file handle
Only one of the three alternate creator entries can be used at a time
and which of those is determined by the style of plugin. The
"creator.v1" element is for factory functions of supported styles
computing forces (i.e. command, pair, bond, angle, dihedral, or
improper styles) and the function takes as single argument the pointer
to the LAMMPS instance. The factory function is cast to the
``lammpsplugin_factory1`` type before assignment. The "creator.v2"
element is for factory functions creating an instance of a fix, compute,
or region style and takes three arguments: a pointer to the LAMMPS
instance, an integer with the length of the argument list and a ``char
**`` pointer to the list of arguments. The factory function pointer
needs to be cast to the ``lammpsplugin_factory2`` type before
assignment.
Only one of the two alternate creator entries can be used at a time and
which of those is determined by the style of plugin. The "creator.v1"
element is for factory functions of supported styles computing forces
(i.e. pair, bond, angle, dihedral, or improper styles) or command styles
and the function takes as single argument the pointer to the LAMMPS
instance. The factory function is cast to the ``lammpsplugin_factory1``
type before assignment. The "creator.v2" element is for factory
functions creating an instance of a fix, compute, or region style and
takes three arguments: a pointer to the LAMMPS instance, an integer with
the length of the argument list and a ``char **`` pointer to the list of
arguments. The factory function pointer needs to be cast to the
``lammpsplugin_factory2`` type before assignment.
Pair style example
^^^^^^^^^^^^^^^^^^
@ -249,3 +257,8 @@ by ``#ifdef PAIR_CLASS`` is not needed, since the mapping of the class
name to the style name is done by the plugin registration function with
the information from the ``lammpsplugin_t`` struct. It may be included
in case the new code is intended to be later included in LAMMPS directly.
A plugin may be registered under an existing style name. In that case
the plugin will override the existing code. This can be used to modify
the behavior of existing styles or to debug new versions of them without
having to recompile/reinstall all of LAMMPS.

View File

@ -2154,6 +2154,11 @@ A :doc:`plugin <plugin>` command that can load and unload several
kind of styles in LAMMPS from shared object files at runtime without
having to recompile and relink LAMMPS.
When the environment variable ``LAMMPS_PLUGIN_PATH`` is set, then LAMMPS
will search the directory (or directories) listed in this path for files
with names that end in ``plugin.so`` (e.g. ``helloplugin.so``) and will
try to load the contained plugins automatically at start-up.
**Authors:** Axel Kohlmeyer (Temple U)
**Supporting info:**

View File

@ -33,7 +33,7 @@ Syntax
* R_1, R_2,... = list of cutoff radii, one for each type (distance units)
* w_1, w_2,... = list of neighbor weights, one for each type
* zero or more keyword/value pairs may be appended
* keyword = *rmin0* or *switchflag* or *bzeroflag* or *quadraticflag* or *chem* or *bnormflag* or *wselfallflag* or *bikflag*
* keyword = *rmin0* or *switchflag* or *bzeroflag* or *quadraticflag* or *chem* or *bnormflag* or *wselfallflag* or *bikflag* or *switchinnerflag*
.. parsed-literal::
@ -59,6 +59,9 @@ Syntax
*bikflag* value = *0* or *1* (only implemented for compute snap)
*0* = per-atom bispectrum descriptors are summed over atoms
*1* = per-atom bispectrum descriptors are not summed over atoms
*switchinnerflag* values = *rinnerlist* *drinnerlist*
*rinnerlist* = *ntypes* values of rinner (distance units)
*drinnerlist* = *ntypes* values of drinner (distance units)
Examples
""""""""
@ -70,6 +73,7 @@ Examples
compute vb all sna/atom 1.4 0.95 6 2.0 1.0
compute snap all snap 1.4 0.95 6 2.0 1.0
compute snap all snap 1.0 0.99363 6 3.81 3.83 1.0 0.93 chem 2 0 1
compute snap all snap 1.0 0.99363 6 3.81 3.83 1.0 0.93 switchinnerflag 1.1 1.3 0.5 0.6
Description
"""""""""""
@ -308,6 +312,26 @@ the resulting bispectrum rows are :math:`B_{i,k}` instead of just
:math:`B_k`. In this case, the entries in the final column for these rows
are set to zero.
The keyword *switchinnerflag* activates an additional radial switching
function similar to :math:`f_c(r)` above, but acting to switch off
smoothly contributions from neighbor atoms at short separation distances.
This is useful when SNAP is used in combination with a simple
repulsive potential. The keyword is followed by the *ntypes*
values for :math:`r_{inner}` and the *ntypes*
values for :math:`\Delta r_{inner}`. For a neighbor atom at
distance :math:`r`, its contribution is scaled by a multiplicative
factor :math:`f_{inner}(r)` defined as follows:
.. math::
= & 0, r \leq r_{inner} \\
f_{inner}(r) = & \frac{1}{2}(1 - \cos(\pi \frac{r-r_{inner}}{\Delta r_{inner}})), r_{inner} < r \leq r_{inner} + \Delta r_{inner} \\
= & 1, r > r_{inner} + \Delta r_{inner}
The values of :math:`r_{inner}` and :math:`\Delta r_{inner}` are
the arithmetic means of the values for the central atom of type I
and the neighbor atom of type J.
.. note::
If you have a bonded system, then the settings of :doc:`special_bonds

View File

@ -240,8 +240,6 @@ accelerated styles exist.
* :doc:`latte <fix_latte>` - wrapper on LATTE density-functional tight-binding code
* :doc:`lb/fluid <fix_lb_fluid>` -
* :doc:`lb/momentum <fix_lb_momentum>` -
* :doc:`lb/pc <fix_lb_pc>` -
* :doc:`lb/rigid/pc/sphere <fix_lb_rigid_pc_sphere>` -
* :doc:`lb/viscous <fix_lb_viscous>` -
* :doc:`lineforce <fix_lineforce>` - constrain atoms to move in a line
* :doc:`manifoldforce <fix_manifoldforce>` - restrain atoms to a manifold during minimization

View File

@ -451,7 +451,7 @@ well in parallel. Currently, molecule translations and rotations
are not supported with more than one MPI process.
It is still possible to do parallel molecule exchange without
translation and rotation moves by setting MC moves to zero
and/or by using mcmoves keyword with *Pmoltrans*=*Pmolrotate*=0.
and/or by using the *mcmoves* keyword with *Pmoltrans* = *Pmolrotate* = 0 .
When using fix gcmc in combination with fix shake or fix rigid,

View File

@ -12,55 +12,63 @@ Syntax
* ID, group-ID are documented in :doc:`fix <fix>` command
* lb/fluid = style name of this fix command
* nevery = update the lattice-Boltzmann fluid every this many timesteps
* LBtype = 1 to use the standard finite difference LB integrator,
2 to use the LB integrator of :ref:`Ollila et al. <Ollila>`
* nevery = update the lattice-Boltzmann fluid every this many timesteps (should normally be 1)
* viscosity = the fluid viscosity (units of mass/(time\*length)).
* density = the fluid density.
* zero or more keyword/value pairs may be appended
* keyword = *setArea* or *setGamma* or *scaleGamma* or *dx* or *dm* or *a0* or *noise* or *calcforce* or *trilinear* or *D3Q19* or *read_restart* or *write_restart* or *zwall_velocity* or *bodyforce* or *printfluid*
* keyword = *dx* or *dm* or *noise* or *stencil* or *read_restart* or *write_restart* or *zwall_velocity* or *pressurebcx* or *bodyforce* or *D3Q19* or *dumpxdmf* or *dof* or *scaleGamma* or *a0* or *npits* or *wp* or *sw*
.. parsed-literal::
*setArea* values = type node_area
type = atom type (1-N)
node_area = portion of the surface area of the composite object associated with the particular atom type (used when the force coupling constant is set by default).
*setGamma* values = gamma
gamma = user set value for the force coupling constant.
*scaleGamma* values = type gammaFactor
type = atom type (1-N)
gammaFactor = factor to scale the *setGamma* gamma value by, for the specified atom type.
*dx* values = dx_LB = the lattice spacing.
*dm* values = dm_LB = the lattice-Boltzmann mass unit.
*a0* values = a_0_real = the square of the speed of sound in the fluid.
*noise* values = Temperature seed
Temperature = fluid temperature.
seed = random number generator seed (positive integer)
*calcforce* values = N forcegroup-ID
N = output the force and torque every N timesteps
forcegroup-ID = ID of the particle group to calculate the force and torque of
*trilinear* values = none (used to switch from the default Peskin interpolation stencil to the trilinear stencil).
*D3Q19* values = none (used to switch from the default D3Q15, 15 velocity lattice, to the D3Q19, 19 velocity lattice).
*stencil* values = 2 (trilinear stencil, the default), 3 (3-point immersed boundary stencil), or 4 (4-point Keys' interpolation stencil)
*read_restart* values = restart file = name of the restart file to use to restart a fluid run.
*write_restart* values = N = write a restart file every N MD timesteps.
*zwall_velocity* values = velocity_bottom velocity_top = velocities along the y-direction of the bottom and top walls (located at z=zmin and z=zmax).
*pressurebcx* values = pgradav = imposes a pressure jump at the (periodic) x-boundary of pgradav*Lx*1000.
*bodyforce* values = bodyforcex bodyforcey bodyforcez = the x,y and z components of a constant body force added to the fluid.
*printfluid* values = N = print the fluid density and velocity at each grid point every N timesteps.
*D3Q19* values = none (used to switch from the default D3Q15, 15 velocity lattice, to the D3Q19, 19 velocity lattice).
*dumpxdmf* values = N file timeI
N = output the force and torque every N timesteps
file = output file name
timeI = 1 (use simulation time to index xdmf file), 0 (use output frame number to index xdmf file)
*dof* values = dof = specify the number of degrees of freedom for temperature calculation
*scaleGamma* values = type gammaFactor
type = atom type (1-N)
gammaFactor = factor to scale the *setGamma* gamma value by, for the specified atom type.
*a0* values = a_0_real = the square of the speed of sound in the fluid.
*npits* values = npits h_p l_p l_pp l_e
npits = number of pit regions
h_p = z-height of pit regions (floor to bottom of slit)
l_p = x-length of pit regions
l_pp = x-length of slit regions between consecutive pits
l_e = x-length of slit regions at ends
*wp* values = w_p = y-width of slit regions (defaults to full width if not present or if sw active)
*sw* values = none (turns on y-sidewalls (in xz plane) if npits option active)
Examples
""""""""
.. code-block:: LAMMPS
fix 1 all lb/fluid 1 2 1.0 1.0 setGamma 13.0 dx 4.0 dm 10.0 calcforce sphere1
fix 1 all lb/fluid 1 1 1.0 0.0009982071 setArea 1 1.144592082 dx 2.0 dm 0.3 trilinear noise 300.0 8979873
fix 1 all lb/fluid 1 1.0 0.0009982071 dx 1.2 dm 0.001
fix 1 all lb/fluid 1 1.0 0.0009982071 dx 1.2 dm 0.001 noise 300.0 2761
fix 1 all lb/fluid 1 1.0 1.0 dx 4.0 dm 10.0 dumpxdmf 500 fflow 0 pressurebcx 0.01 npits 2 20 40 5 0 wp 30
Description
"""""""""""
Implement a lattice-Boltzmann fluid on a uniform mesh covering the LAMMPS
simulation domain. The MD particles described by *group-ID* apply a velocity
dependent force to the fluid.
Implement a lattice-Boltzmann fluid on a uniform mesh covering the
LAMMPS simulation domain. Note that this fix was updated in 2021 and is
not backward compatible with the previous version. If you need the
previous version, please download an older version of LAMMPS. The MD
particles described by *group-ID* apply a velocity dependent force to
the fluid.
The lattice-Boltzmann algorithm solves for the fluid motion governed by
the Navier Stokes equations,
@ -86,28 +94,23 @@ respectively. Here, we have implemented
\sigma_{\alpha \beta} = -P_{\alpha \beta} = -\rho a_0 \delta_{\alpha \beta}
with :math:`a_0` set to :math:`\frac{1}{3} \frac{dx}{dt}^2` by default.
You should not normally need to change this default.
The algorithm involves tracking the time evolution of a set of partial
distribution functions which evolve according to a velocity
discretized version of the Boltzmann equation,
distribution functions which evolve according to a velocity discretized
version of the Boltzmann equation,
.. math::
\left(\partial_t + e_{i\alpha}\partial_{\alpha}\right)f_i = -\frac{1}{\tau}\left(f_i - f_i^{eq}\right) + W_i
where the first term on the right hand side represents a single time
relaxation towards the equilibrium distribution function, and :math:`\tau` is a
parameter physically related to the viscosity. On a technical note,
we have implemented a 15 velocity model (D3Q15) as default; however,
the user can switch to a 19 velocity model (D3Q19) through the use of
the *D3Q19* keyword. This fix provides the user with the choice of
two algorithms to solve this equation, through the specification of
the keyword *LBtype*\ . If *LBtype* is set equal to 1, the standard
finite difference LB integrator is used. If *LBtype* is set equal to
2, the algorithm of :ref:`Ollila et al. <Ollila>` is used.
Physical variables are then defined in terms of moments of the distribution
functions,
relaxation towards the equilibrium distribution function, and
:math:`\tau` is a parameter physically related to the viscosity. On a
technical note, we have implemented a 15 velocity model (D3Q15) as
default; however, the user can switch to a 19 velocity model (D3Q19)
through the use of the *D3Q19* keyword. Physical variables are then
defined in terms of moments of the distribution functions,
.. math::
@ -115,7 +118,7 @@ functions,
\rho u_{\alpha} = & \displaystyle\sum\limits_{i} f_i e_{i\alpha}
Full details of the lattice-Boltzmann algorithm used can be found in
:ref:`Mackay et al. <fluid-Mackay>`.
:ref:`Denniston et al. <fluid-Denniston>`.
The fluid is coupled to the MD particles described by *group-ID* through
a velocity dependent force. The contribution to the fluid force on a
@ -127,92 +130,66 @@ calculated as:
{\bf F}_{j \alpha} = \gamma \left({\bf v}_n - {\bf u}_f \right) \zeta_{j\alpha}
where :math:`\mathbf{v}_n` is the velocity of the MD particle,
:math:`\mathbf{u}_f` is the fluid
velocity interpolated to the particle location, and :math:`\gamma` is the force
coupling constant. :math:`\zeta` is a weight assigned to the grid point,
obtained by distributing the particle to the nearest lattice sites.
For this, the user has the choice between a trilinear stencil, which
provides a support of 8 lattice sites, or the immersed boundary method
Peskin stencil, which provides a support of 64 lattice sites. While
the Peskin stencil is seen to provide more stable results, the
trilinear stencil may be better suited for simulation of objects close
to walls, due to its smaller support. Therefore, by default, the
Peskin stencil is used; however the user may switch to the trilinear
stencil by specifying the keyword, *trilinear*\ .
:math:`\mathbf{u}_f` is the fluid velocity interpolated to the particle
location, and :math:`\gamma` is the force coupling constant. This
force, as with most forces in LAMMPS, and hence the velocities, are
calculated at the half-time step. :math:`\zeta` is a weight assigned to
the grid point, obtained by distributing the particle to the nearest
lattice sites.
By default, the force coupling constant, :math:`\gamma`, is calculated
The force coupling constant, :math:`\gamma`, is calculated
according to
.. math::
\gamma = \frac{2m_um_v}{m_u+m_v}\left(\frac{1}{\Delta t_{collision}}\right)
\gamma = \frac{2m_um_v}{m_u+m_v}\left(\frac{1}{\Delta t}\right)
Here, :math:`m_v` is the mass of the MD particle, :math:`m_u` is a
representative fluid mass at the particle location, and :math:`\Delta
t_{collision}` is a collision time, chosen such that
:math:`\frac{\tau}{\Delta t_{collision}} = 1` (see :ref:`Mackay and
Denniston <Mackay2>` for full details). In order to calculate :math:`m_u`,
the fluid density is interpolated to the MD particle location, and
multiplied by a volume, node_area * :math:`dx_{LB}`, where node_area
represents the portion of the surface area of the composite object
associated with a given MD particle. By default, node_area is set
equal to :math:`dx_{LB}^2`; however specific values for given atom types
can be set using the *setArea* keyword.
The user also has the option of specifying their own value for the
force coupling constant, for all the MD particles associated with the
fix, through the use of the *setGamma* keyword. This may be useful
when modelling porous particles. See :ref:`Mackay et al. <fluid-Mackay>` for a
detailed description of the method by which the user can choose an
appropriate :math:`\gamma` value.
representative fluid mass at the particle location, and :math:`\Delta t`
is the time step. The fluid mass :math:`m_u` that the MD particle
interacts with is calculated internally. This coupling is chosen to
constrain the particle and associated fluid velocity to match at the end
of the time step. As with other constraints, such as :doc:`shake
<fix_shake>`, this constraint can remove degrees of freedom from the
simulation which are accounted for internally in the algorithm.
.. note::
while this fix applies the force of the particles on the fluid,
it does not apply the force of the fluid to the particles. When the
force coupling constant is set using the default method, there is only
one option to include this hydrodynamic force on the particles, and
that is through the use of the :doc:`lb/viscous <fix_lb_viscous>` fix.
This fix adds the hydrodynamic force to the total force acting on the
particles, after which any of the built-in LAMMPS integrators can be
used to integrate the particle motion. However, if the user specifies
their own value for the force coupling constant, as mentioned in
:ref:`Mackay et al. <fluid-Mackay>`, the built-in LAMMPS integrators may prove to
be unstable. Therefore, we have included our own integrators
:doc:`fix lb/rigid/pc/sphere <fix_lb_rigid_pc_sphere>`, and
:doc:`fix lb/pc <fix_lb_pc>`, to solve for the particle motion in these
cases. These integrators should not be used with the
:doc:`lb/viscous <fix_lb_viscous>` fix, as they add hydrodynamic forces
to the particles directly. In addition, they can not be used if the
force coupling constant has been set the default way.
.. note::
if the force coupling constant is set using the default method,
and the :doc:`lb/viscous <fix_lb_viscous>` fix is NOT used to add the
While this fix applies the force of the particles on the fluid, it
does not apply the force of the fluid to the particles. There is
only one option to include this hydrodynamic force on the particles,
and that is through the use of the :doc:`lb/viscous <fix_lb_viscous>`
fix. This fix adds the hydrodynamic force to the total force acting
on the particles, after which any of the built-in LAMMPS integrators
can be used to integrate the particle motion. If the
:doc:`lb/viscous <fix_lb_viscous>` fix is NOT used to add the
hydrodynamic force to the total force acting on the particles, this
physically corresponds to a situation in which an infinitely massive
particle is moving through the fluid (since collisions between the
particle and the fluid do not act to change the particle's velocity).
Therefore, the user should set the mass of the particle to be
significantly larger than the mass of the fluid at the particle
location, in order to approximate an infinitely massive particle (see
the dragforce test run for an example).
In this case, setting *scaleGamma* to -1 for the corresponding
particle type will explicitly take this limit (of infinite particle
mass) in computing the force coupling for the fluid force.
----------
Inside the fix, parameters are scaled by the lattice-Boltzmann
Physical parameters describing the fluid are specified through
*viscosity* and *density*. These parameters should all be given in
terms of the mass, distance, and time units chosen for the main LAMMPS
run, as they are scaled by the LB timestep, lattice spacing, and mass
unit, inside the fix.
The *dx* keyword allows the user to specify a value for the LB grid
spacing and the *dm* keyword allows the user to specify the LB mass
unit. Inside the fix, parameters are scaled by the lattice-Boltzmann
timestep, :math:`dt_{LB}`, grid spacing, :math:`dx_{LB}`, and mass unit,
:math:`dm_{LB}`. :math:`dt_{LB}` is set equal to
:math:`\mathrm{nevery}\cdot dt_{MD}`, where :math:`dt_{MD}` is the MD timestep.
By default,
:math:`dm_{LB}` is set equal to 1.0, and :math:`dx_{LB}` is chosen so that
:math:`\frac{\tau}{dt} = \frac{3\eta dt}{\rho dx^2}` is approximately equal to 1.
However, the user has the option of specifying their own values for
:math:`dm_{LB}`, and :math:`dx_{LB}`, by using
the optional keywords *dm*, and *dx* respectively.
:math:`\mathrm{nevery}\cdot dt_{MD}`, where :math:`dt_{MD}` is the MD
timestep. By default, :math:`dm_{LB}` is set equal to 1.0, and
:math:`dx_{LB}` is chosen so that :math:`\frac{\tau}{dt} = \frac{3\eta
dt}{\rho dx^2}` is approximately equal to 1.
.. note::
.. note::
Care must be taken when choosing both a value for :math:`dx_{LB}`,
and a simulation domain size. This fix uses the same subdivision of
@ -223,74 +200,27 @@ the optional keywords *dm*, and *dx* respectively.
with equal lengths in all dimensions, and the default value for
:math:`dx_{LB}` is used, this will automatically be satisfied.
Physical parameters describing the fluid are specified through
*viscosity*, *density*, and *a0*\ . If the force coupling constant is
set the default way, the surface area associated with the MD particles
is specified using the *setArea* keyword. If the user chooses to
specify a value for the force coupling constant, this is set using the
*setGamma* keyword. These parameters should all be given in terms of
the mass, distance, and time units chosen for the main LAMMPS run, as
they are scaled by the LB timestep, lattice spacing, and mass unit,
inside the fix.
----------
The *setArea* keyword allows the user to associate a surface area with
a given atom type. For example if a spherical composite object of
radius R is represented as a spherical shell of N evenly distributed
MD particles, all of the same type, the surface area per particle
associated with that atom type should be set equal to :math:`\frac{4\pi R^2}{N}`.
This keyword should only be used if the force coupling constant,
:math:`\gamma`, is set the default way.
The *setGamma* keyword allows the user to specify their own value for
the force coupling constant, :math:`\gamma`, instead of using the default
value.
The *scaleGamma* keyword should be used in conjunction with the
*setGamma* keyword, when the user wishes to specify different :math:`\gamma`
values for different atom types. This keyword allows the user to
scale the *setGamma* :math:`\gamma` value by a factor, gammaFactor,
for a given atom type.
The *dx* keyword allows the user to specify a value for the LB grid
spacing.
The *dm* keyword allows the user to specify the LB mass unit.
If the *a0* keyword is used, the value specified is used for the
square of the speed of sound in the fluid. If this keyword is not
present, the speed of sound squared is set equal to
:math:`\frac{1}{3}\left(\frac{dx_{LB}}{dt_{LB}}\right)^2`.
Setting :math:`a0 > (\frac{dx_{LB}}{dt_{LB}})^2` is not allowed,
as this may lead to instabilities.
If the *noise* keyword is used, followed by a positive temperature
value, and a positive integer random number seed, a thermal
lattice-Boltzmann algorithm is used. If *LBtype* is set equal to 1
(i.e. the standard LB integrator is chosen), the thermal LB algorithm
of :ref:`Adhikari et al. <Adhikari>` is used; however if *LBtype* is set
equal to 2 both the LB integrator, and thermal LB algorithm described
in :ref:`Ollila et al. <Ollila>` are used.
value, and a positive integer random number seed, the thermal LB algorithm
of :ref:`Adhikari et al. <Adhikari>` is used.
If the *calcforce* keyword is used, both the fluid force and torque
acting on the specified particle group are printed to the screen every
N timesteps.
If the keyword *stencil* is used, the value sets the number of
interpolation points used in each direction. For this, the user has the
choice between a trilinear stencil (*stencil* 2), which provides a
support of 8 lattice sites, or the 3-point immersed boundary method
stencil (*stencil* 3), which provides a support of 27 lattice sites, or
the 4-point Keys' interpolation stencil (stencil 4), which provides a
support of 64 lattice sites. The trilinear stencil is the default as it
is better suited for simulation of objects close to walls or other
objects, due to its smaller support. The 3-point stencil provides
smoother motion of the lattice and is suitable for particles not likely
to be to close to walls or other objects.
If the keyword *trilinear* is used, the trilinear stencil is used to
interpolate the particle nodes onto the fluid mesh. By default, the
immersed boundary method, Peskin stencil is used. Both of these
interpolation methods are described in :ref:`Mackay et al. <fluid-Mackay>`.
If the keyword *D3Q19* is used, the 19 velocity (D3Q19) lattice is
used by the lattice-Boltzmann algorithm. By default, the 15 velocity
(D3Q15) lattice is used.
If the keyword *write_restart* is used, followed by a positive
integer, N, a binary restart file is printed every N LB timesteps.
This restart file only contains information about the fluid.
Therefore, a LAMMPS restart file should also be written in order to
print out full details of the simulation.
If the keyword *write_restart* is used, followed by a positive integer,
N, a binary restart file is printed every N LB timesteps. This restart
file only contains information about the fluid. Therefore, a LAMMPS
restart file should also be written in order to print out full details
of the simulation.
.. note::
@ -308,19 +238,100 @@ conditions in the z-direction. If fixed boundary conditions are
present in the z-direction, and this keyword is not used, the walls
are assumed to be stationary.
If the *pressurebcx* keyword is used, a pressure jump (implemented by a
step jump in density) is imposed at the (periodic) x-boundary. The
value set specifies what would be the resulting equilibrium average
pressure gradient in the x-direction if the system had a constant
cross-section (i.e. resistance to flow). It is converted to a pressure
jump by multiplication by the system size in the x-direction. As this
value should normally be quite small, it is also assumed to be scaled
by 1000.
If the *bodyforce* keyword is used, a constant body force is added to
the fluid, defined by it's x, y and z components.
If the *printfluid* keyword is used, followed by a positive integer, N,
the fluid densities and velocities at each lattice site are printed to the
screen every N timesteps.
If the keyword *D3Q19* is used, the 19 velocity (D3Q19) lattice is
used by the lattice-Boltzmann algorithm. By default, the 15 velocity
(D3Q15) lattice is used.
If the *dumpxdmf* keyword is used, followed by a positive integer, N,
and a file name, the fluid densities and velocities at each lattice site
are output to an xdmf file every N timesteps. This is a binary file
format that can be read by visualization packages such as `Paraview
<https://www.paraview.org/>`_ . The xdmf file format contains a time
index for each frame dump and the value timeI = 1 uses simulation time
while 0 uses the output frame number to index xdmf file. The later can
be useful if the :doc:`dump vtk <dump_vtk>` command is used to output
the particle positions at the same timesteps and you want to visualize
both the fluid and particle data together in `Paraview
<https://www.paraview.org/>`_ .
The *scaleGamma* keyword allows the user to scale the :math:`\gamma`
value by a factor, gammaFactor, for a given atom type. Setting
*scaleGamma* to -1 for the corresponding particle type will explicitly
take the limit of infinite particle mass in computing the force coupling
for the fluid force (see note above).
If the *a0* keyword is used, the value specified is used for the square
of the speed of sound in the fluid. If this keyword is not present, the
speed of sound squared is set equal to
:math:`\frac{1}{3}\left(\frac{dx_{LB}}{dt_{LB}}\right)^2`. Setting
:math:`a0 > (\frac{dx_{LB}}{dt_{LB}})^2` is not allowed, as this may
lead to instabilities. As the speed of sound should usually be much
larger than any fluid velocity of interest, its value does not normally
have a significant impact on the results. As such, it is usually best
to use the default for this option.
The *npits* keyword (followed by integer arguments: npits, h_p, l_p,
l_pp, l_e) sets the fluid domain to the pits geometry. These arguments
should only be used if you actually want something more complex than a
rectangular/cubic geometry. The npits value sets the number of pits
regions (arranged along x). The remaining arguments are sizes measured
in multiples of dx_lb: h_p is the z-height of the pit regions, l_p is
the x-length of the pit regions, l_pp is the length of the region
between consecutive pits (referred to as a "slit" region), and l_e is
the x-length of the slit regions at each end of the channel. The pit
geometry must fill the system in the x-direction but can be longer, in
which case it is truncated (which enables asymmetric entrance/exit end
sections). The additional *wp* keyword allows the width (in
y-direction) of the pit to be specified (the default is full width) and
the *sw* keyword indicates that there should be sidewalls in the
y-direction (default is periodic in y-direction). These parameters are
illustrated below::
Sideview (in xz plane) of pit geometry:
______________________________________________________________________
slit slit slit ^
|
<---le---><---------lp-------><---lpp---><-------lp--------><---le---> hs = (Nbz-1) - hp
|
__________ __________ __________ v
| | | | ^ z
| | | | | |
| pit | | pit | hp +-x
| | | | |
|__________________| |__________________| v
Endview (in yz plane) of pit geometry (no sw so wp is active):
_____________________
^
|
hs
|
_____________________ v
| | ^
| | | z
|<---wp--->| hp |
| | | +-y
|__________| v
----------
For further details, as well as descriptions and results of several
test runs, see :ref:`Mackay et al. <fluid-Mackay>`. Please include a citation to
this paper if the lb_fluid fix is used in work contributing to
published research.
For further details, as well as descriptions and results of several test
runs, see :ref:`Denniston et al. <fluid-Denniston>`. Please include a
citation to this paper if the lb_fluid fix is used in work contributing
to published research.
----------
@ -333,68 +344,77 @@ binary restart files, if requested, independent of the main LAMMPS
is written to the main LAMMPS :doc:`binary restart files <restart>`.
None of the :doc:`fix_modify <fix_modify>` options are relevant to this
fix. No global or per-atom quantities are stored by this fix for
access by various :doc:`output commands <Howto_output>`. No parameter
of this fix can be used with the *start/stop* keywords of the
:doc:`run <run>` command. This fix is not invoked during :doc:`energy minimization <minimize>`.
fix.
The fix computes a global scalar which can be accessed by various
:doc:`output commands <Howto_output>`. The scalar is the current
temperature of the group of particles described by *group-ID* along with
the fluid constrained to move with them. The temperature is computed via
the kinetic energy of the group and fluid constrained to move with them
and the total number of degrees of freedom (calculated internally). If
the particles are not integrated independently (such as via :doc:`fix
NVE <fix_nve>`) but have additional constraints imposed on them (such as
via integration using :doc:`fix rigid <fix_rigid>`) the degrees of
freedom removed from these additional constraints will not be properly
accounted for. In this case, the user can specify the total degrees of
freedom independently using the *dof* keyword.
The fix also computes a global array of values which can be accessed by
various :doc:`output commands <Howto_output>`. There are 5 entries in
the array. The first entry is the temperature of the fluid, the second
entry is the total mass of the fluid plus particles, the third through
fifth entries give the x, y, and z total momentum of the fluid plus
particles.
No parameter of this fix can be used with the *start/stop* keywords of
the :doc:`run <run>` command. This fix is not invoked during
:doc:`energy minimization <minimize>`.
Restrictions
""""""""""""
This fix is part of the LATBOLTZ package. It is only enabled if LAMMPS
was built with that package. See the :doc:`Build package <Build_package>` page for more info.
was built with that package. See the :doc:`Build package
<Build_package>` page for more info.
This fix can only be used with an orthogonal simulation domain.
Walls have only been implemented in the z-direction. Therefore, the
boundary conditions, as specified via the main LAMMPS boundary command
must be periodic for x and y, and either fixed or periodic for z.
Shrink-wrapped boundary conditions are not permitted with this fix.
The boundary conditions for the fluid are specified independently to the
particles. However, these should normally be specified consistently via
the main LAMMPS :doc:`boundary <boundary>` command (p p p, p p f, and p
f f are the only consistent possibilities). Shrink-wrapped boundary
conditions are not permitted with this fix.
This fix must be used before any of :doc:`fix lb/viscous <fix_lb_viscous>`, :doc:`fix lb/momentum <fix_lb_momentum>`, :doc:`fix lb/rigid/pc/sphere <fix_lb_rigid_pc_sphere>`, and/ or :doc:`fix lb/pc <fix_lb_pc>` , as the fluid needs to be initialized before
any of these routines try to access its properties. In addition, in
order for the hydrodynamic forces to be added to the particles, this
fix must be used in conjunction with the
:doc:`lb/viscous <fix_lb_viscous>` fix if the force coupling constant is
set by default, or either the :doc:`lb/viscous <fix_lb_viscous>` fix or
one of the :doc:`lb/rigid/pc/sphere <fix_lb_rigid_pc_sphere>` or
:doc:`lb/pc <fix_lb_pc>` integrators, if the user chooses to specify
their own value for the force coupling constant.
This fix must be used before any of :doc:`fix lb/viscous
<fix_lb_viscous>` and :doc:`fix lb/momentum <fix_lb_momentum>` as the
fluid needs to be initialized before any of these routines try to access
its properties. In addition, in order for the hydrodynamic forces to be
added to the particles, this fix must be used in conjunction with the
:doc:`lb/viscous <fix_lb_viscous>` fix.
This fix needs to be used in conjunction with a standard LAMMPS
integrator such as :doc:`fix NVE <fix_nve>` or :doc:`fix rigid
<fix_rigid>`.
Related commands
""""""""""""""""
:doc:`fix lb/viscous <fix_lb_viscous>`, :doc:`fix lb/momentum <fix_lb_momentum>`, :doc:`fix lb/rigid/pc/sphere <fix_lb_rigid_pc_sphere>`, :doc:`fix lb/pc <fix_lb_pc>`
:doc:`fix lb/viscous <fix_lb_viscous>`, :doc:`fix lb/momentum <fix_lb_momentum>`
Default
"""""""
By default, the force coupling constant is set according to
.. math::
\gamma = \frac{2m_um_v}{m_u+m_v}\left(\frac{1}{\Delta t_{collision}}\right)
and an area of :math:`dx_{LB}^2` per node, used to calculate the fluid mass at
the particle node location, is assumed.
*dx* is chosen such that :math:`\frac{\tau}{dt_{LB}} =
\frac{3\eta dt_{LB}}{\rho dx_{LB}^2}` is approximately equal to 1.
*dx* is chosen such that :math:`\frac{\tau}{dt_{LB}} = \frac{3\eta dt_{LB}}{\rho dx_{LB}^2}` is approximately equal to 1.
*dm* is set equal to 1.0.
*a0* is set equal to :math:`\frac{1}{3}\left(\frac{dx_{LB}}{dt_{LB}}\right)^2`.
The Peskin stencil is used as the default interpolation method.
The trilinear stencil is used as the default interpolation method.
The D3Q15 lattice is used for the lattice-Boltzmann algorithm.
If walls are present, they are assumed to be stationary.
----------
.. _Ollila:
.. _fluid-Denniston:
**(Ollila et al.)** Ollila, S.T.T., Denniston, C., Karttunen, M., and Ala-Nissila, T., Fluctuating lattice-Boltzmann model for complex fluids, J. Chem. Phys. 134 (2011) 064902.
.. _fluid-Mackay:
**(Mackay et al.)** Mackay, F. E., Ollila, S.T.T., and Denniston, C., Hydrodynamic Forces Implemented into LAMMPS through a lattice-Boltzmann fluid, Computer Physics Communications 184 (2013) 2021-2031.
**(Denniston et al.)** Denniston, C., Afrasiabian, N., Cole-Andre, M.G., Mackay, F. E., Ollila, S.T.T., and Whitehead, T., LAMMPS lb/fluid fix version 2: Improved Hydrodynamic Forces Implemented into LAMMPS through a lattice-Boltzmann fluid, Computer Physics Communications 275 (2022) `108318 <https://doi.org/10.1016/j.cpc.2022.108318>`_ .
.. _Mackay2:
@ -403,3 +423,4 @@ If walls are present, they are assumed to be stationary.
.. _Adhikari:
**(Adhikari et al.)** Adhikari, R., Stratford, K., Cates, M. E., and Wagner, A. J., Fluctuating lattice Boltzmann, Europhys. Lett. 71 (2005) 473-479.

View File

@ -38,7 +38,7 @@ lattice-Boltzmann fluid is present.
Zero the total linear momentum of the system, including both the atoms
specified by group-ID and the lattice-Boltzmann fluid every nevery
timesteps. This is accomplished by adjusting the particle velocities
timesteps. If there are no atoms specified by group-ID only the fluid momentum is affected. This is accomplished by adjusting the particle velocities
and the fluid velocities at each lattice site.
.. note::

View File

@ -1,65 +0,0 @@
.. index:: fix lb/pc
fix lb/pc command
=================
Syntax
""""""
.. parsed-literal::
fix ID group-ID lb/pc
* ID, group-ID are documented in the :doc:`fix <fix>` command
* lb/pc = style name of this fix command
Examples
""""""""
.. code-block:: LAMMPS
fix 1 all lb/pc
Description
"""""""""""
Update the positions and velocities of the individual particles
described by *group-ID*, experiencing velocity-dependent hydrodynamic
forces, using the integration algorithm described in :ref:`Mackay et al. <Mackay1>`. This integration algorithm should only be used if a
user-specified value for the force-coupling constant used in :doc:`fix lb/fluid <fix_lb_fluid>` has been set; do not use this integration
algorithm if the force coupling constant has been set by default.
Restart, fix_modify, output, run start/stop, minimize info
"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
No information about this fix is written to :doc:`binary restart files <restart>`. None of the :doc:`fix_modify <fix_modify>` options
are relevant to this fix. No global or per-atom quantities are stored
by this fix for access by various :doc:`output commands <Howto_output>`.
No parameter of this fix can be used with the *start/stop* keywords of
the :doc:`run <run>` command. This fix is not invoked during :doc:`energy minimization <minimize>`.
Restrictions
""""""""""""
This fix is part of the LATBOLTZ package. It is only enabled if LAMMPS
was built with that package. See the :doc:`Build package <Build_package>` page for more info.
Can only be used if a lattice-Boltzmann fluid has been created via the
:doc:`fix lb/fluid <fix_lb_fluid>` command, and must come after this
command.
Related commands
""""""""""""""""
:doc:`fix lb/fluid <fix_lb_fluid>` :doc:`fix lb/rigid/pc/sphere <fix_lb_rigid_pc_sphere>`
Default
"""""""
none.
----------
.. _Mackay1:
**(Mackay et al.)** Mackay, F. E., Ollila, S.T.T., and Denniston, C., Hydrodynamic Forces Implemented into LAMMPS through a lattice-Boltzmann fluid, Computer Physics Communications 184 (2013) 2021-2031.

View File

@ -1,167 +0,0 @@
.. index:: fix lb/rigid/pc/sphere
fix lb/rigid/pc/sphere command
==============================
Syntax
""""""
.. parsed-literal::
fix ID group-ID lb/rigid/pc/sphere bodystyle args keyword values ...
* ID, group-ID are documented in :doc:`fix <fix>` command
* lb/rigid/pc/sphere = style name of this fix command
* bodystyle = *single* or *molecule* or *group*
.. parsed-literal::
*single* args = none
*molecule* args = none
*group* args = N groupID1 groupID2 ...
N = # of groups
* zero or more keyword/value pairs may be appended
* keyword = *force* or *torque* or *innerNodes*
.. parsed-literal::
*force* values = M xflag yflag zflag
M = which rigid body from 1-Nbody (see asterisk form below)
xflag,yflag,zflag = off/on if component of center-of-mass force is active
*torque* values = M xflag yflag zflag
M = which rigid body from 1-Nbody (see asterisk form below)
xflag,yflag,zflag = off/on if component of center-of-mass torque is active
*innerNodes* values = innergroup-ID
innergroup-ID = ID of the atom group which does not experience a hydrodynamic force from the lattice-Boltzmann fluid
Examples
""""""""
.. code-block:: LAMMPS
fix 1 spheres lb/rigid/pc/sphere
fix 1 all lb/rigid/pc/sphere force 1 0 0 innerNodes ForceAtoms
Description
"""""""""""
This fix is based on the :doc:`fix rigid <fix_rigid>` command, and was
created to be used in place of that fix, to integrate the equations of
motion of spherical rigid bodies when a lattice-Boltzmann fluid is
present with a user-specified value of the force-coupling constant.
The fix uses the integration algorithm described in :ref:`Mackay et
al. <Mackay>` to update the positions, velocities, and orientations of
a set of spherical rigid bodies experiencing velocity dependent
hydrodynamic forces. The spherical bodies are assumed to rotate as
solid, uniform density spheres, with moments of inertia calculated
using the combined sum of the masses of all the constituent particles
(which are assumed to be point particles).
----------
By default, all of the atoms that this fix acts on experience a
hydrodynamic force due to the presence of the lattice-Boltzmann fluid.
However, the *innerNodes* keyword allows the user to specify atoms
belonging to a rigid object which do not interact with the
lattice-Boltzmann fluid (i.e. these atoms do not feel a hydrodynamic
force from the lattice-Boltzmann fluid). This can be used to
distinguish between atoms on the surface of a non-porous object, and
those on the inside.
This feature can be used, for example, when implementing a hard sphere
interaction between two spherical objects. Instead of interactions
occurring between the particles on the surfaces of the two spheres, it
is desirable simply to place an atom at the center of each sphere,
which does not contribute to the hydrodynamic force, and have these
central atoms interact with one another.
----------
Apart from the features described above, this fix is very similar to
the rigid fix (although it includes fewer optional arguments, and
assumes the constituent atoms are point particles); see
:doc:`fix rigid <fix_rigid>` for a complete documentation.
Restart, fix_modify, output, run start/stop, minimize info
"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
No information about the *rigid* and *rigid/nve* fixes are written to
:doc:`binary restart files <restart>`.
The :doc:`fix_modify <fix_modify>` *virial* option is supported by
this fix to add the contribution due to the added forces on atoms to
both the global pressure and per-atom stress of the system via the
:doc:`compute pressure <compute_pressure>` and :doc:`compute
stress/atom <compute_stress_atom>` commands. The former can be
accessed by :doc:`thermodynamic output <thermo_style>`. The default
setting for this fix is :doc:`fix_modify virial yes <fix_modify>`.
Similar to the :doc:`fix rigid <fix_rigid>` command: The rigid fix
computes a global scalar which can be accessed by various :doc:`output
commands <Howto_output>`. The scalar value calculated by these fixes
is "intensive". The scalar is the current temperature of the
collection of rigid bodies. This is averaged over all rigid bodies
and their translational and rotational degrees of freedom. The
translational energy of a rigid body is 1/2 m v\^2, where m = total
mass of the body and v = the velocity of its center of mass. The
rotational energy of a rigid body is 1/2 I w\^2, where I = the moment
of inertia tensor of the body and w = its angular velocity. Degrees
of freedom constrained by the *force* and *torque* keywords are
removed from this calculation.
All of these fixes compute a global array of values which can be
accessed by various :doc:`output commands <Howto_output>`. The number
of rows in the array is equal to the number of rigid bodies. The
number of columns is 15. Thus for each rigid body, 15 values are
stored: the xyz coords of the center of mass (COM), the xyz components
of the COM velocity, the xyz components of the force acting on the
COM, the xyz components of the torque acting on the COM, and the xyz
image flags of the COM, which have the same meaning as image flags for
atom positions (see the "dump" command). The force and torque values
in the array are not affected by the *force* and *torque* keywords in
the fix rigid command; they reflect values before any changes are made
by those keywords.
The ordering of the rigid bodies (by row in the array) is as follows.
For the *single* keyword there is just one rigid body. For the
*molecule* keyword, the bodies are ordered by ascending molecule ID.
For the *group* keyword, the list of group IDs determines the ordering
of bodies.
The array values calculated by these fixes are "intensive", meaning
they are independent of the number of atoms in the simulation.
No parameter of these fixes can be used with the *start/stop* keywords
of the :doc:`run <run>` command. These fixes are not invoked during
:doc:`energy minimization <minimize>`.
Restrictions
""""""""""""
This fix is part of the LATBOLTZ package. It is only enabled if LAMMPS
was built with that package. See the :doc:`Build package
<Build_package>` page for more info.
Can only be used if a lattice-Boltzmann fluid has been created via the
:doc:`fix lb/fluid <fix_lb_fluid>` command, and must come after this
command. Should only be used if the force coupling constant used in
:doc:`fix lb/fluid <fix_lb_fluid>` has been set by the user; this
integration fix cannot be used if the force coupling constant is set
by default.
Related commands
""""""""""""""""
:doc:`fix lb/fluid <fix_lb_fluid>`, :doc:`fix lb/pc <fix_lb_pc>`
Default
"""""""
The defaults are force \* on on on, and torque \* on on on.
----------
.. _Mackay:
**(Mackay et al.)** Mackay, F. E., Ollila, S.T.T., and Denniston, C., Hydrodynamic Forces Implemented into LAMMPS through a lattice-Boltzmann fluid, Computer Physics Communications 184 (2013) 2021-2031.

View File

@ -25,27 +25,14 @@ Description
This fix is similar to the :doc:`fix viscous <fix_viscous>` command, and
is to be used in place of that command when a lattice-Boltzmann fluid
is present, and the user wishes to integrate the particle motion using
one of the built in LAMMPS integrators.
is present using the :doc:`fix lb/fluid <fix_lb_fluid>`. This should be used in conjunction with one of the built-in LAMMPS integrators, such as :doc:`fix NVE <fix_nve>` or :doc:`fix rigid <fix_rigid>`.
This fix adds a force, F = - Gamma\*(velocity-fluid_velocity), to each
atom, where Gamma is the force coupling constant described in the :doc:`fix lb/fluid <fix_lb_fluid>` command (which applies an equal and
opposite force to the fluid).
.. note::
This fix should only be used in conjunction with one of the
built in LAMMPS integrators; it should not be used with the :doc:`fix lb/pc <fix_lb_pc>` or :doc:`fix lb/rigid/pc/sphere <fix_lb_rigid_pc_sphere>` integrators, which
already include the hydrodynamic forces. These latter fixes should
only be used if the force coupling constant has been set by the user
(instead of using the default value); if the default force coupling
value is used, then this fix provides the only method for adding the
hydrodynamic forces to the particles.
This fix adds a viscous force to each atom to cause it move with the same velocity as the fluid (an equal and opposite force is applied to the fluid via :doc:`fix lb/fluid <fix_lb_fluid>`). When :doc:`fix lb/fluid <fix_lb_fluid>` is called with the noise option, the atoms will also experience random forces which will thermalize them to the same temperature as the fluid. In this way, the combination of this fix with :doc:`fix lb/fluid <fix_lb_fluid>` and a LAMMPS integrator like :doc:`fix NVE <fix_nve>` is analogous to :doc:`fix langevin <fix_langevin>` except here the fluid is explicit. The temperature of the particles can be monitored via the scalar output of :doc:`fix lb/fluid <fix_lb_fluid>`.
----------
For further details, as well as descriptions and results of several
test runs, see :ref:`Mackay et al. <Mackay3>`. Please include a citation to
For details of this fix, as well as descriptions and results of several
test runs, see :ref:`Denniston et al. <fluid-Denniston2>`. Please include a citation to
this paper if this fix is used in work contributing to published
research.
@ -78,14 +65,11 @@ Can only be used if a lattice-Boltzmann fluid has been created via the
:doc:`fix lb/fluid <fix_lb_fluid>` command, and must come after this
command.
This fix should not be used if either the :doc:`fix lb/pc <fix_lb_pc>`
or :doc:`fix lb/rigid/pc/sphere <fix_lb_rigid_pc_sphere>` integrator is
used.
Related commands
""""""""""""""""
:doc:`fix lb/fluid <fix_lb_fluid>`, :doc:`fix lb/pc <fix_lb_pc>`, :doc:`fix lb/rigid/pc/sphere <fix_lb_rigid_pc_sphere>`
:doc:`fix lb/fluid <fix_lb_fluid>`
Default
"""""""
@ -94,6 +78,6 @@ none
----------
.. _Mackay3:
.. _fluid-Denniston2:
**(Mackay et al.)** Mackay, F. E., Ollila, S.T.T., and Denniston, C., Hydrodynamic Forces Implemented into LAMMPS through a lattice-Boltzmann fluid, Computer Physics Communications 184 (2013) 2021-2031.
**(Denniston et al.)** Denniston, C., Afrasiabian, N., Cole-Andre, M.G., Mackay, F. E., Ollila, S.T.T., and Whitehead, T., LAMMPS lb/fluid fix version 2: Improved Hydrodynamic Forces Implemented into LAMMPS through a lattice-Boltzmann fluid, Computer Physics Communications 275 (2022) `108318 <https://doi.org/10.1016/j.cpc.2022.108318>`_ .

View File

@ -21,7 +21,7 @@ Syntax
* ID, group-ID are documented in :doc:`fix <fix>` command
* style = *polarize/bem/gmres* or *polarize/bem/icc* or *polarize/functional*
* Nevery = this fixed is invoked every this many timesteps
* tolerance = the tolerance for the iterative solver to stop
* tolerance = the relative tolerance for the iterative solver to stop
Examples
@ -45,10 +45,45 @@ Description
"""""""""""
These fixes compute induced charges at the interface between two
impermeable media with different dielectric constants.
impermeable media with different dielectric constants. The interfaces
need to be discretized into vertices, each representing a boundary element.
The vertices are treated as if they were regular atoms or particles.
:doc:`atom_style dielectric <atom_style>` should be used since it defines
the additional properties of each interface particle such as
interface normal vectors, element areas, and local dielectric mismatch.
These fixes also require the use of :doc:`pair_style <pair_style>` and
:doc:`kspace_style <kspace_style>` with the *dielectric* suffix.
At every time step, given a configuration of the physical charges in the system
(such as atoms and charged particles) these fixes compute and update
the charge of the interface particles. The interfaces are allowed to move
during the simulation with appropriate time integrators (for example,
with :doc:`fix_rigid <fix_rigid>`).
There are some example scripts for using this fix
with LAMMPS in the examples/PACKAGES/dielectric directory.
Consider an interface between two media: one with dielectric constant
of 78 (water), the other of 4 (silica). The interface is discretized
into 2000 boundary elements, each represented by an interface particle. Suppose that
each interface particle has a normal unit vector pointing from the silica medium to water.
The dielectric difference along the normal vector is then 78 - 4 = 74,
the mean dielectric value is (78 + 4) / 2 = 41. Each boundary element
also has its area and the local mean curvature (which is used by these fixes
for computing a correction term in the local electric field).
To model charged interfaces, the interface particle will have a non-zero charge value,
coming from its area and surface charge density.
For non-interface particles such as atoms and charged particles,
the interface normal vectors, element area, and dielectric mismatch are
irrelevant. Their local dielectric value is used to rescale their actual charge
when computing the Coulombic interactions. For instance, for a cation carrying
a charge of +2 (in charge unit) in an implicit solvent with dielectric constant of 40
would have actual charge of +2, and a local dielectric constant value of 40.
It is assumed that the particles cannot pass through the interface during the simulation
so that its local dielectric constant value does not change.
There are some example scripts for using these fixes
with LAMMPS in the ``examples/PACKAGES/dielectric directory``. The README file
therein contains specific details on the system setup. Note that the example data files
show the additional fields (columns) needed for :doc:`atom_style dielectric <atom_style>`
beyond the conventional fields *id*, *mol*, *type*, *q*, *x*, *y*, and *z*.
----------
@ -75,12 +110,41 @@ as described in :ref:`(Barros) <Barros>` to solve :math:`\sigma_b`.
Fix *polarize/bem/icc* employs the successive over-relaxation algorithm
as described in :ref:`(Tyagi) <Tyagi>` to solve :math:`\sigma_b`.
Fix *polarize/functional* ...
The iterative solvers would terminate either when the maximum relative change
in the induced charges in consecutive iterations is below the set tolerance,
or when the number of iterations reaches *iter_max* (see below).
Fix *polarize/functional* employs the energy functional variation approach
as described in :ref:`(Jadhao) <Jadhao>` to solve :math:`\sigma_b`.
More details on the implementation of these fixes and their recommended use
are described in :ref:`(NguyenTD) <NguyenTD>`.
Restart, fix_modify, output, run start/stop, minimize info
""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
...
No information about this fix is written to :doc:`binary restart files <restart>`.
The :doc:`fix_modify <fix_modify>` command provides certain options to
control the induced charge solver and the initial values of the interface elements:
.. parsed-literal::
*itr_max* arg
arg = maximum number of iterations for convergence
*dielectrics* ediff emean epsilon area charge
ediff = dielectric difference
emean = dielectric mean
epsilon = local dielectric value
aree = element area
charge = real interface charge
*polarize/bem/gmres* or *polarize/bem/icc* compute a global 2-element vector
which can be accessed by various :doc:`output commands <Howto_output>`.
The first element is the number of iterations when the solver terminates
(of which the upperbound is set by *iter_max*). The second element is the RMS error.
Restrictions
""""""""""""
@ -94,12 +158,15 @@ KSPACE package is installed. See the :doc:`Build package
Related commands
""""""""""""""""
:doc:`pair_coeff <pair_coeff>`, :doc:`fix polarize <fix_polarize>`, :doc:`read_data <read_data>`,
:doc:`pair_style lj/cut/coul/long/dielectric <pair_dielectric>`,
:doc:`kspace_style pppm/dielectric <kspace_style>`,
:doc:`compute efield/atom <compute_efield_atom>`
Default
"""""""
None.
*iter_max* = 20
----------

View File

@ -135,7 +135,8 @@ with #) anywhere. Each non-blank non-comment line must contain one
keyword/value pair. The required keywords are *rcutfac* and
*twojmax*\ . Optional keywords are *rfac0*, *rmin0*,
*switchflag*, *bzeroflag*, *quadraticflag*, *chemflag*,
*bnormflag*, *wselfallflag*, *chunksize*, and *parallelthresh*\ .
*bnormflag*, *wselfallflag*, *switchinnerflag*,
*rinner*, *drinner*, *chunksize*, and *parallelthresh*\ .
The default values for these keywords are
@ -147,6 +148,7 @@ The default values for these keywords are
* *chemflag* = 0
* *bnormflag* = 0
* *wselfallflag* = 0
* *switchinnerflag* = 0
* *chunksize* = 32768
* *parallelthresh* = 8192
@ -189,6 +191,16 @@ corresponding *K*-vector of linear coefficients for element
which must equal the number of unique elements appearing in the LAMMPS
pair_coeff command, to avoid ambiguity in the number of coefficients.
The keyword *switchinnerflag* activates an additional switching function
that smoothly turns off contributions to the SNAP potential from neighbor
atoms at short separations. If *switchinnerflag* is set to 1 then
the additional keywords *rinner* and *drinner* must also be provided.
Each of these is followed by *nelements* values, where *nelements*
is the number of unique elements appearing in appearing in the LAMMPS
pair_coeff command. The element order should correspond to the order
in which elements first appear in the pair_coeff command reading from
left to right.
The keywords *chunksize* and *parallelthresh* are only applicable when
using the pair style *snap* with the KOKKOS package on GPUs and are
ignored otherwise. The *chunksize* keyword controls the number of atoms

View File

@ -56,6 +56,15 @@ styles and names.
The *clear* command will unload all currently loaded plugins.
.. admonition:: Automatic loading of plugins
:class: note
When the environment variable ``LAMMPS_PLUGIN_PATH`` is set, then
LAMMPS will search the directory (or directories) listed in this path
for files with names that end in ``plugin.so``
(e.g. ``helloplugin.so``) and will try to load the contained plugins
automatically at start-up.
Restrictions
""""""""""""