ICODE-ADC/lammps
4
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Proofing RHEO package

Browse Source
This commit is contained in:
jtclemm
2024-06-28 17:07:23 -06:00
parent 3c0eaf6870
commit 1326592f23
40 changed files with 3797 additions and 355 deletions
Download Patch File Download Diff File Expand all files Collapse all files

1
doc/src/Commands_fix.rst
View File

@ -28,6 +28,7 @@ OPT.
* :doc:`adapt <fix_adapt>`
* :doc:`adapt/fep <fix_adapt_fep>`
* :doc:`addforce <fix_addforce>`
* :doc:`add/heat <fix_add_heat>`
* :doc:`addtorque <fix_addtorque>`
* :doc:`alchemy <fix_alchemy>`
* :doc:`amoeba/bitorsion <fix_amoeba_bitorsion>`

38
doc/src/Howto_rheo.rst
View File

@ -2,8 +2,9 @@ Reproducing hydrodynamics and elastic objects (RHEO)
====================================================
The RHEO package is a hybrid implementation of smoothed particle
hydrodynamics (SPH) for fluid flow, coupled to the :doc:`BPM package <Howto_bpm>` to model
solid elements. RHEO combines these methods to enable mesh-free modeling of multiphase material systems. The SPH solver supports many advanced options
hydrodynamics (SPH) for fluid flow, coupled to the :doc:`BPM package <Howto_bpm>`
to model solid elements. RHEO combines these methods to enable mesh-free modeling
of multiphase material systems. The SPH solver supports many advanced options
including reproducing kernels, particle shifting, free surface identification,
and solid surface reconstruction. To model fluid-solid systems, the status of
particles can dynamically change between a fluid and solid state, e.g. during
@ -29,9 +30,10 @@ of reproducing kernels). In conjunction to fix rheo, one must specify an
instance of :doc:`fix rheo/pressure <fix_rheo_pressure>` and
:doc:`fix rheo/viscosity <fix_rheo_viscosity>` to define a pressure equation
of state and viscosity model, respectively. Optionally, one can model
a heat equation with :doc:`fix rheo/thermal`, which also allows the user
to specify equations for a particle's thermal conductivity, specific heat,
latent heat, and melting temperature. The ordering of these fixes in an an input script matters. Fix rheo must be defined prior to all
a heat equation with :doc:`fix rheo/thermal <fix_rhe0_thermal>`, which also
allows the user to specify equations for a particle's thermal conductivity,
specific heat, latent heat, and melting temperature. The ordering of these
fixes in an an input script matters. Fix rheo must be defined prior to all
other RHEO fixes.
Typically, RHEO requires atom style rheo. In addition to typical atom
@ -45,9 +47,9 @@ affect particles. Instead, one should use the *sph/e* attribute.
The status variable uses bitmasking to track various properties of a particle
such as its current state of matter (fluid or solid) and its location relative
to a surface. Some of these properties (and others) can be accessed using
:doc:`compute rheo/property/atom <fix_rheo_property_atom>`. The *status* attribute
in :doc:`the set command <set>` only allows control over the first bit which sets
the state of matter, 0 is fluid and 1 is solid.
:doc:`compute rheo/property/atom <compute_rheo_property_atom>`. The *status*
attribute in :doc:`the set command <set>` only allows control over the first bit
which sets the state of matter, 0 is fluid and 1 is solid.
Fluid interactions, including pressure forces, viscous forces, and heat exchange,
are calculated using :doc:`pair rheo <pair_rheo>`. Unlike typical pair styles,
@ -59,18 +61,18 @@ hydrodynamic forces are only calculated if a fluid particle is involved.
To model elastic objects, there are currently two mechanisms in RHEO, one designed
for bulk solid bodies and the other for thin shells. Both mechanisms rely on
introducing bonded forces between particles and therefore require a hybrid of atom style bond and rheo
(or rheo/thermal).
introducing bonded forces between particles and therefore require a hybrid of atom
style bond and rheo (or rheo/thermal).
To create an elastic solid body, one has to (a) change the status of constituent
particles to solid (e.g. with the :doc:`set <set>` command), (b) create bpm
bonds between the particles (see the :doc:`bpm howto <Howto_bpm>` page for
more details), and (c) use :doc:`pair rheo/solid <pair_rheo_solid>` to
apply repulsive contact forces between distinct solid bodies. Akin to pair rheo,
looks at a particles fluid/solid status to determine whether to apply forces.
However, unlike pair rheo, pair rheo/solid does obey special bond settings such
that contact forces do not have to be calculated between two bonded solid particles
in the same elastic body.
pair rheo/solid considers a particles fluid/solid phase to determine whether to
apply forces. However, unlike pair rheo, pair rheo/solid does obey special bond
settings such that contact forces do not have to be calculated between two bonded
solid particles in the same elastic body.
In systems with thermal evolution, fix rheo/thermal can optionally set a
melting/solidification temperature allowing particles to dynamically swap their
@ -79,9 +81,9 @@ critical temperature, respectively. Using the *react* option, one can specify a
bond length and a bond type. Then, when solidifying, particles will search their
local neighbors and automatically create bonds with any neighboring solid particles
in range. For BPM bond styles, bonds will then use the immediate position of the two
particles to calculate a reference state. When melting, particles will then delete
any bonds of the specified type when reverting to a fluid state. Special bonds are
updated as bonds are created/broken.
particles to calculate a reference state. When melting, particles will delete any
bonds of the specified type when reverting to a fluid state. Special bonds are updated
as bonds are created/broken.
The other option for elastic objects is an elastic shell that is nominally much
thinner than a particle diameter, e.g. a oxide skin which gradually forms over time
@ -106,7 +108,7 @@ criteria for creating/deleting a bond or altering force calculations).
.. _howto_rheo_palermo:
**(Palermo)** Palermo, Clemmer, Wolf, O'Connor, in preparation.
**(Palermo)** Palermo, Wolf, Clemmer, O'Connor, in preparation.
.. _howto_rheo_clemmer:

2
doc/src/atom_style.rst
View File

@ -194,7 +194,7 @@ the Additional Information section below.
- :ref:`RHEO <PKG-RHEO>`
- solid and fluid RHEO particles
* - *rheo/thermal*
- *atomic* + rho, status, temperature
- *atomic* + rho, status, energy, temperature
- :ref:`RHEO <PKG-RHEO>`
- RHEO particles with temperature
* - *sph*

47
doc/src/bond_rheo_shell.rst
View File

@ -44,12 +44,11 @@ The *rheo/shell* bond style is designed to work with
:doc:`fix rheo/oxidation <fix_rheo_oxidation>` which creates candidate
bonds between eligible surface or near-surface particles. When a bond
is first created, it computes no forces and starts a timer. Forces are
not computed until the timer reaches the specified bond formation time
and the bond is fully enabled. If the two particles move outside of the
maximum bond distance or move into the bulk before the timer reaches
the formation time, the bond automatically deletes itself. Not that
this deletion does not generate any broken bond data saved to a
*store/local* fix.
not computed until the timer reaches the specified bond formation time,
*t/form*, and the bond is enabled and applies forces. If the two particles
move outside of the maximum bond distance or move into the bulk before
the timer reaches *t/form*, the bond automatically deletes itself. This
deletion is not recorded as a broken bond in the optional *store/local* fix.
Before bonds are enabled, they are still treated as regular bonds by
all other parts of LAMMPS. This means they are written to data files
@ -60,17 +59,17 @@ To only count enabled bonds, use the *nbond/shell* attribute in
When enabled, the bond then computes forces based on deviations from
the initial reference state of the two atoms much like a BPM style
bond (as further discussed in the :doc:`BPM howto page <Howto_bpm>`).
The reference state is stored by each bond when it is first enabled
the setup of a run. Data is then preserved across run commands and is
written to :doc:`binary restart files <restart>` such that restarting
the system will not reset the reference state of a bond or the timer.
The reference state is stored by each bond when it is first enabled.
Data is then preserved across run commands and is written to
:doc:`binary restart files <restart>` such that restarting the system
will not reset the reference state of a bond or the timer.
This bond style is based on a model described in
:ref:`(Clemmer) <howto_rheo_clemmer>`. The force has a magnitude of
:ref:`(Clemmer) <rheo_clemmer>`. The force has a magnitude of
.. math::
F = 2 k (r - r_0) + \frac{2 * k}{r_0^2 \epsilon_c^2} (r - r_0)^3
F = 2 k (r - r_0) + \frac{2 k}{r_0^2 \epsilon_c^2} (r - r_0)^3
where :math:`k` is a stiffness, :math:`r` is the current distance
and :math:`r_0` is the initial distance between the two particles, and
@ -78,17 +77,15 @@ and :math:`r_0` is the initial distance between the two particles, and
is done by setting the bond type to 0 such that forces are no longer
computed.
An additional damping force is applied to the bonded
particles. This forces is proportional to the difference in the
normal velocity of particles using a similar construction as
dissipative particle dynamics :ref:`(Groot) <Groot4>`:
A damping force proportional to the difference in the normal velocity
of particles is also applied to bonded particles:
.. math::
F_D = - \gamma w (\hat{r} \bullet \vec{v})
where :math:`\gamma` is the damping strength, :math:`\hat{r}` is the
radial normal vector, and :math:`\vec{v}` is the velocity difference
displacement normal vector, and :math:`\vec{v}` is the velocity difference
between the two particles.
The following coefficients must be defined for each bond type via the
@ -103,7 +100,7 @@ the data file or restart files read by the :doc:`read_data
Unlike other BPM-style bonds, this bond style does not update special
bond settings when bonds are created or deleted. This bond style also
does not enforce specific :doc:`special_bonds <special_bonds>` settings.
This behavior is purposeful such :doc:`RHEO pair forces <pair_rheo>`
This behavior is purposeful such :doc:`RHEO pair <pair_rheo>` forces
and heat flows are still calculated.
If the *store/local* keyword is used, an internal fix will track bonds that
@ -162,10 +159,10 @@ the specified attribute.
The single() function of this bond style returns 0.0 for the energy
of a bonded interaction, since energy is not conserved in these
dissipative potentials. The single() function also calculates an
extra bond quantity, the initial distance :math:`r_0`. This
extra quantity can be accessed by the
:doc:`compute bond/local <compute_bond_local>` command as *b1*\ .
dissipative potentials. The single() function also calculates two
extra bond quantities, the initial distance :math:`r_0` and a time.
These extra quantities can be accessed by the
:doc:`compute bond/local <compute_bond_local>` command as *b1* and *b2*\ .
Restrictions
""""""""""""
@ -186,10 +183,6 @@ NA
----------
.. _howto_rheo_clemmer:
.. _rheo_clemmer:
**(Clemmer)** Clemmer, Pierce, O'Connor, Nevins, Jones, Lechman, Tencer, Appl. Math. Model., 130, 310-326 (2024).
.. _Groot4:
**(Groot)** Groot and Warren, J Chem Phys, 107, 4423-35 (1997).

1
doc/src/bond_style.rst
View File

@ -105,6 +105,7 @@ accelerated styles exist.
* :doc:`oxdna2/fene <bond_oxdna>` - same as oxdna but used with different pair styles
* :doc:`oxrna2/fene <bond_oxdna>` - modified FENE bond suitable for RNA modeling
* :doc:`quartic <bond_quartic>` - breakable quartic bond
* :doc:`rheo/shell <bond_rheo_shell>` - shell bond for oxidation modeling in RHEO
* :doc:`special <bond_special>` - enable special bond exclusions for 1-5 pairs and beyond
* :doc:`table <bond_table>` - tabulated by bond length

43
doc/src/compute_rheo_property_atom.rst
View File

@ -1,7 +1,7 @@
.. index:: compute rheo/property/atom
compute rheo/property/atom command
=============================
==================================
Syntax
""""""
@ -17,11 +17,10 @@ Syntax
.. parsed-literal::
possible attributes = phase, surface, surface/r,
surface/divr, surface/n/:math:`\alpha`, coordination,
shift/v/:math:`\alpha`, energy, temperature, heatflow,
surface/divr, surface/n/a, coordination,
shift/v/a, energy, temperature, heatflow,
conductivity, cv, viscosity, pressure, rho,
grad/v/:math:`\alpha \beta`, stress/v/:math:`\alpha \beta`,
stress/t/:math:`\alpha \beta`, nbond/shell
grad/v/ab, stress/v/ab, stress/t/ab, nbond/shell
.. parsed-literal::
@ -29,9 +28,9 @@ Syntax
*surface* = atom surface status
*surface/r* = atom distance from the surface
*surface/divr* = divergence of position at atom position
*surface/n/:math:`\alpha`* = surface normal vector
*surface/n/a* = a-component of surface normal vector
*coordination* = coordination number
*shift/v/:math:`\alpha`* = atom shifting velocity
*shift/v/a* = a-component of atom shifting velocity
*energy* = atom energy
*temperature* = atom temperature
*heatflow* = atom heat flow
@ -40,9 +39,9 @@ Syntax
*viscosity* = atom viscosity
*pressure* = atom pressure
*rho* = atom density
*grad/v/:math:`\alpha \beta`* = atom velocity gradient
*stress/v/:math:`\alpha \beta`* = atom viscous stress
*stress/t/:math:`\alpha \beta`* = atom total stress (pressure and viscous)
*grad/v/ab* = ab-component of atom velocity gradient tensor
*stress/v/ab* = ab-component of atom viscous stress tensor
*stress/t/ab* = ab-component of atom total stress tensor (pressure and viscous)
*nbond/shell* = number of oxide bonds
Examples
@ -68,12 +67,12 @@ computes as inputs. See for example, the
:doc:`fix ave/chunk <fix_ave_chunk>`, and
:doc:`atom-style variable <variable>` commands.
For vector attributes, e.g. *shift/v/:math:`\alpha`*, one must specify
For vector attributes, e.g. *shift/v/*:math:`\alpha`, one must specify
:math:`\alpha` as the *x*, *y*, or *z* component, e.g. *shift/v/x*.
Alternatively, a wild card \* will include all components, *x* and *y* in
2D or *x*, *y*, and *z* in 3D.
For tensor attributes, e.g. *grad/v/:math:`\alpha \beta`*, one must specify
For tensor attributes, e.g. *grad/v/*:math:`\alpha \beta`, one must specify
both :math:`\alpha` and :math:`\beta` as *x*, *y*, or *z*, e.g. *grad/v/xy*.
Alternatively, a wild card \* will include all components. In 2D, this
includes *xx*, *xy*, *yx*, and *yy*. In 3D, this includes *xx*, *xy*, *xz*,
@ -93,11 +92,11 @@ the *interface/reconstruct* option of :doc:`fix rheo <fix_rheo>`. Bulk
particles have a value of 0, surface particles have a value of 1, and
splash particles have a value of 2. The *surface/r* property is the
distance from the surface, up to the kernel cutoff length. Surface particles
have a value of 0. The *surface/n* properties are the components of the
surface normal vector.
have a value of 0. The *surface/n/*:math:`\alpha` properties are the
components of the surface normal vector.
The *shift/v/:math:`\alpha`* properties are the components of the shifting velocity
produced by the *shift* option of :doc:`fix rheo <fix_rheo>`.
The *shift/v/*:math:`\alpha` properties are the components of the shifting
velocity produced by the *shift* option of :doc:`fix rheo <fix_rheo>`.
The *nbond/shell* property is the number of shell bonds that have been
activated from :doc:`bond style rheo/shell <bond_rheo_shell>`.
@ -110,13 +109,13 @@ Output info
"""""""""""
This compute calculates a per-atom vector or per-atom array depending
on the number of input values. If a single input is specified, a
per-atom vector is produced. If two or more inputs are specified, a
on the number of input values. Generally, if a single input is specified,
a per-atom vector is produced. If two or more inputs are specified, a
per-atom array is produced where the number of columns = the number of
inputs. If a wild card \* is used for a vector or tensor, then the number
of inputs is considered to be incremented by the dimensiod or dimension
squared, respectively. The vector or array can be accessed by any command
that uses per-atom values from a compute as input. See the
inputs. However, if a wild card \* is used for a vector or tensor, then
the number of inputs is considered to be incremented by the dimension or
the dimension squared, respectively. The vector or array can be accessed
by any command that uses per-atom values from a compute as input. See the
:doc:`Howto output <Howto_output>` page for an overview of LAMMPS output
options.

6
doc/src/fix.rst
View File

@ -193,6 +193,7 @@ accelerated styles exist.
* :doc:`adapt <fix_adapt>` - change a simulation parameter over time
* :doc:`adapt/fep <fix_adapt_fep>` - enhanced version of fix adapt
* :doc:`addforce <fix_addforce>` - add a force to each atom
* :doc:`add/heat <fix_add_heat>` - add a heat flux to each atom
* :doc:`addtorque <fix_addtorque>` - add a torque to a group of atoms
* :doc:`alchemy <fix_alchemy>` - perform an "alchemical transformation" between two partitions
* :doc:`amoeba/bitorsion <fix_amoeba_bitorsion>` - torsion/torsion terms in AMOEBA force field
@ -371,8 +372,9 @@ accelerated styles exist.
* :doc:`restrain <fix_restrain>` - constrain a bond, angle, dihedral
* :doc:`rheo <fix_rheo>` - integrator for the RHEO package
* :doc:`rheo/thermal <fix_rheo_thermal>` - thermal integrator for the RHEO package
* :doc:`rheo/pressure <fix_rheo_pressure>` - pressure derivation for the RHEO package
* :doc:`rheo/viscosity <fix_rheo_pressure>` - viscosity derivation for the RHEO package
* :doc:`rheo/oxidation <fix_rheo_oxidation>` - create oxidation bonds for the RHEO package
* :doc:`rheo/pressure <fix_rheo_pressure>` - pressure calculation for the RHEO package
* :doc:`rheo/viscosity <fix_rheo_pressure>` - viscosity calculation for the RHEO package
* :doc:`rhok <fix_rhok>` - add bias potential for long-range ordered systems
* :doc:`rigid <fix_rigid>` - constrain one or more clusters of atoms to move as a rigid body with NVE integration
* :doc:`rigid/meso <fix_rigid_meso>` - constrain clusters of mesoscopic SPH/SDPD particles to move as a rigid body

36
doc/src/fix_add_heat.rst
View File

@ -16,14 +16,14 @@ Syntax
.. parsed-literal::
*constant* args = rate
rate = rate of heat flow (energy/time units)
*linear* args = t_target k
t_target = target temperature (temperature units)
k = prefactor (energy/(time*temperature) units)
*quartic* args = t_target k
t_target = target temperature (temperature units)
k = prefactor (energy/(time*temperature^4) units)
*constant* args = *rate*
*rate* = rate of heat flow (energy/time units)
*linear* args = :math:`T_{target}` *k*
:math:`T_{target}` = target temperature (temperature units)
*k* = prefactor (energy/(time*temperature) units)
*quartic* args = :math:`T_{target}` *k*
:math:`T_{target}` = target temperature (temperature units)
*k* = prefactor (energy/(time*temperature^4) units)
* zero or more keyword/value pairs may be appended to args
* keyword = *overwrite*
@ -45,7 +45,9 @@ Examples
Description
"""""""""""
This fix adds heat to particles every timestep.
This fix adds heat to particles with the temperature attribute every timestep.
Note that this is an internal temperature of a particle intended for use with
non-atomistic models like the discrete element method.
For the *constant* style, heat is added at the specified rate. For the *linear* style,
heat is added at a rate of :math:`k (T_{target} - T)` where :math:`k` is the
@ -62,22 +64,22 @@ determine the rate of heat added.
Equal-style variables can specify formulas with various mathematical
functions and include :doc:`thermo_style <thermo_style>` command
keywords for the simulation box parameters, time step, and elapsed time.
Thus, it is easy to specify time-dependent heating.
keywords for the simulation box parameters, time step, and elapsed time
to specify time-dependent heating.
Atom-style variables can specify the same formulas as equal-style
variables but can also include per-atom values, such as atom
coordinates. Thus, it is easy to specify a spatially-dependent heating
field with optional time-dependence as well.
coordinates to specify spatially-dependent heating.
If the *overwrite* keyword is set to *yes*, this fix will effectively set
the total heat flow on a particle, overwriting contributions from pair
styles or other fixes.
If the *overwrite* keyword is set to *yes*, this fix will set the total
heat flow on a particle every timestep, overwriting contributions from pair
styles or other fixes. If *overwrite* is *no*, this fix will add heat on
top of other contributions.
----------
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.

74
doc/src/fix_rheo.rst
View File

@ -1,7 +1,7 @@
.. index:: fix rheo
fix rheo command
===============
================
Syntax
""""""
@ -38,8 +38,8 @@ Examples
.. code-block:: LAMMPS
fix 1 all rheo 1.0 quintic thermal density 0.1 speed/sound 10.0
fix 1 all rheo 1.0 RK1 shift surface/detection coordination 40
fix 1 all rheo 3.0 quintic 0 thermal density 0.1 0.1 speed/sound 10.0 1.0
fix 1 all rheo 3.0 RK1 10 shift surface/detection coordination 40
Description
"""""""""""
@ -48,14 +48,15 @@ Description
Perform time integration for RHEO particles, updating positions, velocities,
and densities. For a detailed breakdown of the integration timestep and
numerical details, see :ref:`(Palermo) <howto_rheo_palermo>`. For an
numerical details, see :ref:`(Palermo) <rheo_palermo>`. For an
overview of other features available in the RHEO package, see
:doc:`the RHEO howto <Howto_rheo>`.
The type of kernel is specified using *kstyle* and the cutoff is *cut*. Four
kernels are currently available. The *quintic* kernel is a standard quintic
spline function commonly used in SPH. The other options, *RK0*, *RK1*, and
*RK2*, are zeroth, first, and second order reproducing. To generate a reproducing kernel, a particle must have sufficient neighbors inside the
*RK2*, are zeroth, first, and second order reproducing. To generate a
reproducing kernel, a particle must have sufficient neighbors inside the
kernel cutoff distance (a coordination number) to accurately calculate
moments. This threshold is set by *zmin*. If reproducing kernels are
requested but a particle has fewer neighbors, then it will revert to a
@ -81,28 +82,27 @@ surfaces.
A modified form of Fickian particle shifting can be enabled with the
*shift* keyword. This effectively shifts particle positions to generate a
more uniform spatial distribution. Shifting currently does consider the
more uniform spatial distribution. Shifting currently does not consider the
type of a particle and therefore may be inappropriate in systems consisting
of multiple materials.
of multiple fluid phases.
In systems with free surfaces, the *surface/detection* keyword can be used
to classify the location of particles as being within the bulk fluid, on a
free surface, or isolated from other particles in a splash or droplet.
Shifting is then disabled in the direction away from the free surface to
prevent it from diffusing particles away from the bulk fluid. Surface
detection can also be used to control surface-nucleated effects like
oxidation when used in combination with
:doc:`fix rheo/oxidation <fix_rheo_oxidation>`. Surface detection is not
Shifting is then disabled in the normal direction away from the free surface
to prevent particles from difusing away. Surface detection can also be used
to control surface-nucleated effects like oxidation when used in combination
with :doc:`fix rheo/oxidation <fix_rheo_oxidation>`. Surface detection is not
performed on solid bodies.
The *surface/detection* keyword takes three arguments: *sdstyle*, *limit*,
and *limi/splash*. The first, *sdstyle*, specifies whether surface particles
and *limit/splash*. The first, *sdstyle*, specifies whether surface particles
are identified using a coordination number (*coordination*) or the divergence
of the local particle positions (*divergence*). The threshold value for a
surface particle for either of these criteria is set by the numerical value
of *limit*. Additionally, if a particle's coordination number is too low,
i.e. if it has separated off from the bulk in a droplet, it is not possible
to define surfaces and a particle is classified as a splash. The coordination
to define surfaces and the particle is classified as a splash. The coordination
threshold for this classification is set by the numerical value of
*limit/splash*.
@ -112,23 +112,23 @@ a kernel summation of the masses of neighboring particles by specifying the *rho
keyword.
The *self/mass* keyword modifies the behavior of the density summation in *rho/sum*.
Typically, the density :math:`\rho` of a particle is calculated as the sum
Typically, the density :math:`\rho` of a particle is calculated as the sum over neighbors
.. math::
\rho_i = \Sum_{j} W_{ij} M_j
\rho_i = \sum_{j} W_{ij} M_j
where the summation is over neighbors, :math:`W_{ij}` is the kernel, and :math:`M_j`
is the mass of particle :math:`j`. The *self/mass* keyword augments this expression
by replacing :math:`M_j` with :math:`M_i`. This may be useful in simulations of
multiple fluid phases with large differences in density, :ref:`(Hu) <fix_rheo_hu>`.
where :math:`W_{ij}` is the kernel, and :math:`M_j` is the mass of particle :math:`j`.
The *self/mass* keyword augments this expression by replacing :math:`M_j` with
:math:`M_i`. This may be useful in simulations of multiple fluid phases with large
differences in density, :ref:`(Hu) <fix_rheo_hu>`.
The *density* keyword is used to specify the equilbrium density of each of the N
particle types. It must be followed by N numerical values specifying each
type's equilibrium density *rho0*.
The *density* keyword is used to specify the equilibrium density of each of the N
particle types. It must be followed by N numerical values specifying each type's
equilibrium density *rho0*.
The *speed/sound* keyword is used to specify the speed of sound of each of the
N particle types. It must be followed by N numerical values specifying each
type's speed of sound *cs*.
N particle types. It must be followed by N numerical values specifying each type's
speed of sound *cs*.
Restart, fix_modify, output, run start/stop, minimize info
"""""""""""""""""""""""""""""""""""""""""""""""""""""""""""
@ -144,18 +144,16 @@ the :doc:`run <run>` command. This fix is not invoked during
Restrictions
""""""""""""
This fix must be used with atom style rheo or rheo/thermal.
This fix must be used in conjuction with
:doc:`fix rheo/pressure <fix_rheo_pressure>`. and
:doc:`fix rheo/viscosity <fix_rheo_viscosity>`, If the *thermal*
setting is used, there must also be an instance of
:doc:`fix rheo/thermal <fix_rheo_thermal>`. The fix group must be
set to all. Only one instance of fix rheo may be defined and it
must be defined prior to all other RHEO fixes.
This fix must be used with atom style rheo or rheo/thermal. This fix must
be used in conjuction with :doc:`fix rheo/pressure <fix_rheo_pressure>`.
and :doc:`fix rheo/viscosity <fix_rheo_viscosity>`. If the *thermal* setting
is used, there must also be an instance of
:doc:`fix rheo/thermal <fix_rheo_thermal>`. The fix group must be set to all.
Only one instance of fix rheo may be defined and it must be defined prior
to all other RHEO fixes in the input script.
This fix is part of the RHEO package. It is only enabled if
LAMMPS was built with that package. See the
:doc:`Build package <Build_package>` page for more info.
This fix is part of the RHEO package. It is only enabled if LAMMPS was built
with that package. See the :doc:`Build package <Build_package>` page for more info.
Related commands
""""""""""""""""
@ -173,9 +171,9 @@ Default
----------
.. _howto_rheo_palermo:
.. _rheo_palermo:
**(Palermo)** Palermo, Clemmer, Wolf, O'Connor, in preparation.
**(Palermo)** Palermo, Wolf, Clemmer, O'Connor, in preparation.
.. _fix_rheo_hu:

21
doc/src/fix_rheo_oxidation.rst
View File

@ -36,21 +36,21 @@ for use with bond style :doc:`bond rheo/shell <bond_rheo_shell>`.
Every timestep, particles check neighbors within a distance of *cut*.
This distance must be smaller than the kernel length defined in
:doc:`fix rheo <fix_rheo>`. Bonds of type *btype* are created between
pairs of particles that satisfy one of two conditions. First, if both
particles are on the fluid surface, or within a distance of *rsurf*
from the surface. Secondly, if one particle is on the fluid surface
and the other bond is solid. This process is further described in
:ref:`(Clemmer) <howto_rheo_clemmer>`.
a fluid particle and either a fluid or solid neighbor. The fluid particles
must also be on the fluid surface, or within a distance of *rsurf* from
the surface. This process is further described in
:ref:`(Clemmer) <howto_rheo_clemmer2>`.
If used in conjunction with solid bodies, such as those generated
by the *react* option of :doc:`fix rheo/thermal <fix_rheo_thermal>`,
it is recommended that one uses a :doc:`hybrid bond style <bond_hybrid>`
it is recommended to use a :doc:`hybrid bond style <bond_hybrid>`
with different bond types for solid and oxide bonds.
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
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
@ -59,11 +59,12 @@ the :doc:`run <run>` command. This fix is not invoked during :doc:`energy minim
Restrictions
""""""""""""
This fix must be used with an bond style :doc:`rheo/shell <bond_rheo_shell>`
This fix must be used with the bond style :doc:`rheo/shell <bond_rheo_shell>`
and :doc:`fix rheo <fix_rheo>` with surface detection enabled.
This fix is part of the RHEO package. It is only enabled if
LAMMPS was built with that package. See the :doc:`Build package <Build_package>` page for more info.
LAMMPS was built with that package. See the :doc:`Build package <Build_package>`
page for more info.
Related commands
""""""""""""""""
@ -79,6 +80,6 @@ none
----------
.. _howto_rheo_clemmer:
.. _howto_rheo_clemmer2:
**(Clemmer)** Clemmer, Pierce, O'Connor, Nevins, Jones, Lechman, Tencer, Appl. Math. Model., 130, 310-326 (2024).

14
doc/src/fix_rheo_pressure.rst
View File

@ -1,7 +1,7 @@
.. index:: fix rheo/pressure
fix rheo/pressure command
===============
=========================
Syntax
""""""
@ -20,7 +20,7 @@ Syntax
*linear* args = none
*taitwater* args = none
*cubic* args = cubic term prefactor :math:`A_3` (pressure/density\^2)
*cubic* args = cubic prefactor :math:`A_3` (pressure/density\^2)
Examples
""""""""
@ -36,8 +36,8 @@ Description
.. versionadded:: TBD
This fix defines a pressure equation of state for RHEO particles. One can
define different equations of state for different atom types, but an
equation must be specified for every atom type.
define different equations of state for different atom types. An equation
must be specified for every atom type.
One first defines the atom *types*. A wild-card asterisk can be used in place
of or in conjunction with the *types* argument to set the coefficients for
@ -74,7 +74,8 @@ Style *taitwater* is Tait's equation of state:
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
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
@ -89,7 +90,8 @@ conjuction with :doc:`fix rheo <fix_rheo>`. The fix group must be
set to all. Only one instance of fix rheo/pressure can be defined.
This fix is part of the RHEO package. It is only enabled if
LAMMPS was built with that package. See the :doc:`Build package <Build_package>` page for more info.
LAMMPS was built with that package. See the :doc:`Build package <Build_package>`
page for more info.
Related commands
""""""""""""""""

20
doc/src/fix_rheo_thermal.rst
View File

@ -1,7 +1,7 @@
.. index:: fix rheo/thermal
fix rheo/thermal command
===============
========================
Syntax
""""""
@ -50,16 +50,16 @@ Description
.. versionadded:: TBD
This fix performs time integration of temperature evolution for atom style
rheo/thermal. In addition, it defines multiple thermal properties of
particles and handles melting/solidification, if applicable. For more details
on phase transitions in RHEO, see :doc:`the RHEO howto <Howto_rheo>`.
This fix performs time integration of temperature for atom style rheo/thermal.
In addition, it defines multiple thermal properties of particles and handles
melting/solidification, if applicable. For more details on phase transitions
in RHEO, see :doc:`the RHEO howto <Howto_rheo>`.
Note that the temperature of a particle is always derived from the energy.
This implies the *temperature* attribute of :doc:`the set command <set>` does
not affect particles. Instead, one should use the *sph/e* attribute.
For each atom type, one can define attributes for the *conductivity*,
For each atom type, one can define expressions for the *conductivity*,
*specific/heat*, *latent/heat*, and critical temperature (*Tfreeze*).
The conductivity and specific heat must be defined for all atom types.
The latent heat and critical temperature are optional. However, a
@ -88,13 +88,14 @@ types that have a defined value of *Tfreeze*. When a fluid particle's
temperature drops below *Tfreeze*, bonds of type *btype* are created between
nearby solid particles within a distance of *cut*. The particle's status also
swaps to a solid state. When a solid particle's temperature rises above
*Tfreeze*, all bonds of type *btype* are broken and the particle's tatus swaps
*Tfreeze*, all bonds of type *btype* are broken and the particle's status swaps
to a fluid state.
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
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
@ -110,7 +111,8 @@ must be used in conjuction with :doc:`fix rheo <fix_rheo>` with the
instance of fix rheo/pressure can be defined.
This fix is part of the RHEO package. It is only enabled if
LAMMPS was built with that package. See the :doc:`Build package <Build_package>` page for more info.
LAMMPS was built with that package. See the :doc:`Build package <Build_package>`
page for more info.
Related commands
""""""""""""""""

12
doc/src/fix_rheo_viscosity.rst
View File

@ -1,7 +1,7 @@
.. index:: fix rheo/viscosity
fix rheo/viscosity command
===============
==========================
Syntax
""""""
@ -63,16 +63,16 @@ for the stress :math:`\tau`
\tau = \left(\frac{\tau_0}{\dot{\gamma}} + K \dot{\gamma}^{n - 1}\right) \dot{\gamma}, \tau \ge \tau_0
where :math:`\dot{\gamma}` is the strain rate and :math:`tau_0` is the critical
where :math:`\dot{\gamma}` is the strain rate and :math:`\tau_0` is the critical
yield stress, below which :math:`\dot{\gamma} = 0.0`. To avoid divergences, this
expression is regularized by defining a critical strain rate *gd0*. If the local
strain rate on a particle falls below this limit, a constant viscosity of *eta*
is assigned. This implies a value of
.. math::
\tau_0 = \eta * \dot{\gamma}_0 - K \dot{\gamma}_0^N
\tau_0 = \eta \dot{\gamma}_0 - K \dot{\gamma}_0^N
as further discussed in :ref:`(Palermo) <howto_rheo_palermo>`.
as further discussed in :ref:`(Palermo) <rheo_palermo2>`.
Restart, fix_modify, output, run start/stop, minimize info
@ -112,6 +112,6 @@ none
----------
.. _howto_rheo_palermo:
.. _rheo_palermo2:
**(Palermo)** Palermo, Clemmer, Wolf, O'Connor, in preparation.
**(Palermo)** Palermo, Wolf, Clemmer, O'Connor, in preparation.

36
doc/src/pair_rheo.rst
View File

@ -1,7 +1,7 @@
.. index:: pair_style rheo
pair_style rheo command
=========================
=======================
Syntax
""""""
@ -25,7 +25,7 @@ Examples
.. code-block:: LAMMPS
pair_style rheo 1.0 quintic rho/damp 1.0 artificial/visc 2.0
pair_style rheo 3.0 rho/damp 1.0 artificial/visc 2.0
pair_coeff * *
Description
@ -41,22 +41,23 @@ heat exchanged between particles.
The *artificial/viscosity* keyword is used to specify the magnitude
:math:`\zeta` of an optional artificial viscosity contribution to forces.
This factor can help stabilize simulations by smoothing out small length
scale variations in velocity fields. Artificial viscous forces are only
exchanged by fluid particles unless interfaces are not reconstructed in
fix rheo, in which fluid particles will also exchange artificial viscous
forces with solid particles to improve stability.
scale variations in velocity fields. Artificial viscous forces typically
are only exchanged by fluid particles. However, if interfaces are not
reconstructed in fix rheo, fluid particles will also exchange artificial
viscous forces with solid particles to improve stability.
The *rho/damp* keyword is used to specify the magnitude :math:`\xi` of
an optional pairwise damping term between the density of particles. This
factor can help stabilize simulations by smoothing out small length
scale variations in density fields.
scale variations in density fields. However, in systems that develop
a density gradient in equilibrium (e.g. in a hydrostatic column underlying
gravity), this option may be inappropriate.
If particles have different viscosities or conductivities, the
*harmonic/means* keyword changes how they are averaged before calculating
pairwise forces or heat exchanges. By default, an arithmetic averaged is
used, however, a harmonic mean may improve stability in multiphase systems
with large disparities in viscosities. This keyword has no effect on
results if viscosities and conductivities are constant.
used, however, a harmonic mean may improve stability in systems with multiple
fluid phases with large disparities in viscosities.
No coefficients are defined for each pair of atoms types via the
:doc:`pair_coeff <pair_coeff>` command as in the examples
@ -70,16 +71,20 @@ Mixing, shift, table, tail correction, restart, rRESPA info
This style does not support the :doc:`pair_modify <pair_modify>`
shift, table, and tail options.
This style does not write information to :doc:`binary restart files <restart>`. Thus, you need to re-specify the pair_style and
pair_coeff commands in an input script that reads a restart file.
This style does not write information to :doc:`binary restart files <restart>`.
Thus, you need to re-specify the pair_style and pair_coeff commands in an input
script that reads a restart file.
This style can only be used via the *pair* keyword of the :doc:`run_style respa <run_style>` command. It does not support the *inner*, *middle*, *outer* keywords.
This style can only be used via the *pair* keyword of the
:doc:`run_style respa <run_style>` command. It does not support the *inner*,
*middle*, *outer* keywords.
Restrictions
""""""""""""
This fix is part of the RHEO package. It is only enabled if
LAMMPS was built with that package. See the :doc:`Build package <Build_package>` page for more info.
LAMMPS was built with that package. See the
:doc:`Build package <Build_package>` page for more info.
Related commands
""""""""""""""""
@ -93,4 +98,5 @@ Related commands
Default
"""""""
Density damping and artificial viscous forces are not calculated. Arithmetic means are used for mixing particle properties.
Density damping and artificial viscous forces are not calculated.
Arithmetic means are used for mixing particle properties.

6
doc/src/pair_rheo_solid.rst
View File

@ -1,7 +1,7 @@
.. index:: pair_style rheo/solid
pair_style rheo/solid command
=========================
=============================
Syntax
""""""
@ -44,7 +44,7 @@ normal velocity of particles
F_D = - \gamma w (\hat{r} \bullet \vec{v})
where :math:`\gamma` is the damping strength, :math:`\hat{r}` is the
radial normal vector, :math:`\vec{v}` is the velocity difference
displacement normal vector, :math:`\vec{v}` is the velocity difference
between the two particles, and :math:`w` is a smoothing factor.
This smoothing factor is constructed such that damping forces go to zero
as particles come out of contact to avoid discontinuities. It is
@ -95,7 +95,7 @@ This pair style can only be used via the *pair* keyword of the
Restrictions
""""""""""""
This pair style is part of the BPM package. It is only enabled if
This pair style is part of the RHEO package. It is only enabled if
LAMMPS was built with that package. See the :doc:`Build package
<Build_package>` page for more info.

2
doc/src/pair_style.rst
View File

@ -337,6 +337,8 @@ accelerated styles exist.
* :doc:`reaxff <pair_reaxff>` - ReaxFF potential
* :doc:`rebo <pair_airebo>` - Second generation REBO potential of Brenner
* :doc:`rebomos <pair_rebomos>` - REBOMoS potential for MoS2
* :doc:`rheo <pair_rheo>` - fluid interactions in RHEO package
* :doc:`rheo/solid <pair_rheo_solid>` - solid interactions in RHEO package
* :doc:`resquared <pair_resquared>` - Everaers RE-Squared ellipsoidal potential
* :doc:`saip/metal <pair_saip_metal>` - Interlayer potential for hetero-junctions formed with hexagonal 2D materials and metal surfaces
* :doc:`sdpd/taitwater/isothermal <pair_sdpd_taitwater_isothermal>` - Smoothed dissipative particle dynamics for water at isothermal conditions

2
examples/rheo/balloon/in.rheo.balloon
View File

@ -71,5 +71,5 @@ compute nbond all nbond/atom
thermo 200
thermo_style custom step time ke press atoms
dump 1 all custom 200 atomDump id type x y vx vy fx fy c_phase c_nbond c_rho
#dump 1 all custom 200 atomDump id type x y vx vy fx fy c_phase c_nbond c_rho
run 30000

382
examples/rheo/balloon/log.17Apr2024.balloon.g++.4 Normal file
View File

@ -0,0 +1,382 @@
LAMMPS (17 Apr 2024 - Development - patch_5May2020-18508-g3c0eaf6870-modified)
# ------ 2D water balloon ------ #
dimension 2
units lj
atom_style hybrid rheo bond
boundary m m p
comm_modify vel yes
newton off
region box block -40 40 0 80 -0.01 0.01 units box
create_box 1 box bond/types 1 extra/bond/per/atom 15 extra/special/per/atom 50
Created orthogonal box = (-40 0 -0.01) to (40 80 0.01)
2 by 2 by 1 MPI processor grid
region fluid sphere -10 40 0 30 units box side in
lattice hex 1.0
Lattice spacing in x,y,z = 1.0745699 1.8612097 1.0745699
create_atoms 1 region fluid
Created 2830 atoms
using lattice units in orthogonal box = (-40 0 -0.01) to (40 80 0.01)
create_atoms CPU = 0.001 seconds
region shell sphere -10 40 0 27 units box side out
group shell region shell
544 atoms in group shell
set group shell rheo/status 1
Setting atom values ...
544 settings made for rheo/status
set group all vx 0.005 vy -0.04
Setting atom values ...
2830 settings made for vx
2830 settings made for vy
# ------ Model parameters ------#
variable cut equal 3.0
variable n equal 1.0
variable rho0 equal 1.0
variable cs equal 1.0
variable mp equal ${rho0}/${n}
variable mp equal 1/${n}
variable mp equal 1/1
variable zeta equal 0.05
variable kappa equal 0.01*${rho0}/${mp}
variable kappa equal 0.01*1/${mp}
variable kappa equal 0.01*1/1
variable dt_max equal 0.1*${cut}/${cs}/3
variable dt_max equal 0.1*3/${cs}/3
variable dt_max equal 0.1*3/1/3
variable eta equal 0.05
variable Cv equal 1.0
variable L equal 1.0
variable Tf equal 1.0
mass * ${mp}
mass * 1
timestep 0.1
pair_style hybrid/overlay rheo ${cut} artificial/visc ${zeta} rheo/solid
pair_style hybrid/overlay rheo 3 artificial/visc ${zeta} rheo/solid
pair_style hybrid/overlay rheo 3 artificial/visc 0.05 rheo/solid
pair_coeff * * rheo
pair_coeff * * rheo/solid 1.0 1.0 1.0
special_bonds lj 0.0 1.0 1.0 coul 0.0 1.0 1.0
Finding 1-2 1-3 1-4 neighbors ...
special bond factors lj: 0 1 1
special bond factors coul: 0 1 1
0 = max # of 1-2 neighbors
101 = max # of special neighbors
special bonds CPU = 0.000 seconds
create_bonds many shell shell 1 0 1.5
Generated 0 of 0 mixed pair_coeff terms from geometric mixing rule
Neighbor list info ...
update: every = 1 steps, delay = 0 steps, check = yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 3.3
ghost atom cutoff = 3.3
binsize = 1.65, bins = 49 49 1
3 neighbor lists, perpetual/occasional/extra = 2 1 0
(1) command create_bonds, occasional
attributes: full, newton off
pair build: full/bin
stencil: full/bin/2d
bin: standard
(2) pair rheo, perpetual
attributes: half, newton off
pair build: half/bin/newtoff
stencil: full/bin/2d
bin: standard
(3) pair rheo/solid, perpetual, trim from (2)
attributes: half, newton off, cut 1.3
pair build: trim
stencil: none
bin: none
Added 1263 bonds, new total = 1263
Finding 1-2 1-3 1-4 neighbors ...
special bond factors lj: 0 1 1
special bond factors coul: 0 1 1
6 = max # of 1-2 neighbors
101 = max # of special neighbors
special bonds CPU = 0.000 seconds
special_bonds lj 0.0 1.0 1.0 coul 1.0 1.0 1.0
bond_style bpm/spring
bond_coeff 1 1.0 1.0 1.0
# A lower critical strain allows the balloon to pop
#bond_coeff 1 1.0 0.05 1.0
# ------ Drop balloon ------#
fix 1 all rheo ${cut} quintic 0 shift surface/detection coordination 22 8
fix 1 all rheo 3 quintic 0 shift surface/detection coordination 22 8
fix 2 all rheo/viscosity * constant ${eta}
fix 2 all rheo/viscosity * constant 0.05
fix 3 all rheo/pressure * linear
fix 4 all wall/harmonic ylo EDGE 2.0 1.0 1.0
fix 5 all enforce2d
compute rho all rheo/property/atom rho
compute phase all rheo/property/atom phase
compute nbond all nbond/atom
# ------ Output & Run ------ #
thermo 200
thermo_style custom step time ke press atoms
dump 1 all custom 200 atomDump id type x y vx vy fx fy c_phase c_nbond c_rho
run 30000
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Your simulation uses code contributions which should be cited:
- BPM bond style: doi:10.1039/D3SM01373A
@Article{Clemmer2024,
author = {Clemmer, Joel T. and Monti, Joseph M. and Lechman, Jeremy B.},
title = {A soft departure from jamming: the compaction of deformable
granular matter under high pressures},
journal = {Soft Matter},
year = 2024,
volume = 20,
number = 8,
pages = {1702--1718}
}
- @article{PalermoInPrep,
journal = {in prep},
title = {RHEO: A Hybrid Mesh-Free Model Framework for Dynamic Multi-Phase Flows},
year = {2024},
author = {Eric T. Palermo and Ki T. Wolf and Joel T. Clemmer and Thomas C. O'Connor},
}
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Generated 0 of 0 mixed pair_coeff terms from geometric mixing rule
Neighbor list info ...
update: every = 1 steps, delay = 0 steps, check = yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 3.3
ghost atom cutoff = 3.3
binsize = 1.65, bins = 49 49 1
6 neighbor lists, perpetual/occasional/extra = 6 0 0
(1) pair rheo, perpetual, half/full from (3)
attributes: half, newton off
pair build: halffull/newtoff
stencil: none
bin: none
(2) pair rheo/solid, perpetual, trim from (1)
attributes: half, newton off, cut 1.3
pair build: trim
stencil: none
bin: none
(3) compute RHEO/KERNEL, perpetual
attributes: full, newton off
pair build: full/bin
stencil: full/bin/2d
bin: standard
(4) compute RHEO/GRAD, perpetual, copy from (1)
attributes: half, newton off
pair build: copy
stencil: none
bin: none
(5) compute RHEO/VSHIFT, perpetual, copy from (1)
attributes: half, newton off
pair build: copy
stencil: none
bin: none
(6) compute RHEO/SURFACE, perpetual, copy from (1)
attributes: half, newton off
pair build: copy
stencil: none
bin: none
Per MPI rank memory allocation (min/avg/max) = 17.63 | 17.64 | 17.65 Mbytes
Step Time KinEng Press Atoms
0 0 0.0008125 0.00035927734 2830
200 20 0.0008125 0.00035927734 2830
400 40 0.0008125 0.00035927734 2830
600 60 0.0008125 0.00035927734 2830
800 80 0.0008125 0.00035927734 2830
1000 100 0.0008125 0.00035927734 2830
1200 120 0.0008125 0.00035927734 2830
1400 140 0.0008125 0.00035927734 2830
1600 160 0.0008125 0.00035927734 2830
1800 180 0.0008125 0.00035927734 2830
2000 200 0.0008125 0.00035927734 2830
2200 220 0.0008125 0.00035927734 2830
2400 240 0.00079033569 0.00043037861 2830
2600 260 0.0007549229 0.00045188383 2830
2800 280 0.00072808836 0.00031695003 2830
3000 300 0.0007017958 1.6121754e-05 2830
3200 320 0.00067479047 -0.00015725514 2830
3400 340 0.00064762254 -0.00023361314 2830
3600 360 0.00061960255 -0.00033837679 2830
3800 380 0.0005857206 -0.00051770716 2830
4000 400 0.00055061733 -0.00070309251 2830
4200 420 0.00051884719 -0.0008247795 2830
4400 440 0.00049022236 -0.00099918413 2830
4600 460 0.00046060011 -0.0010923159 2830
4800 480 0.00042900173 -0.0011524571 2830
5000 500 0.00039751503 -0.0012586358 2830
5200 520 0.00036620054 -0.0013973543 2830
5400 540 0.00033130023 -0.0015185231 2830
5600 560 0.00030565892 -0.0016159836 2830
5800 580 0.00028209836 -0.0016925198 2830
6000 600 0.00024695044 -0.0017796892 2830
6200 620 0.00021190635 -0.0018706272 2830
6400 640 0.0001947093 -0.0019146643 2830
6600 660 0.00018903936 -0.0019146199 2830
6800 680 0.00017753371 -0.0019390155 2830
7000 700 0.00015170593 -0.0020247472 2830
7200 720 0.00011509692 -0.0021222209 2830
7400 740 7.9861785e-05 -0.0022033181 2830
7600 760 6.1350463e-05 -0.0022511971 2830
7800 780 6.5269523e-05 -0.0022222806 2830
8000 800 8.5709569e-05 -0.0021089664 2830
8200 820 0.00011746348 -0.0019351493 2830
8400 840 0.00015698134 -0.0017079928 2830
8600 860 0.00019758065 -0.0014618965 2830
8800 880 0.00023338199 -0.0012365832 2830
9000 900 0.00026282353 -0.0010348527 2830
9200 920 0.00028604776 -0.00085287884 2830
9400 940 0.00030388767 -0.000681122 2830
9600 960 0.000317589 -0.00052203521 2830
9800 980 0.00032716728 -0.00037501187 2830
10000 1000 0.00033270692 -0.00025576132 2830
10200 1020 0.00033485986 -0.00016554207 2830
10400 1040 0.00033476763 -9.8525417e-05 2830
10600 1060 0.00033351922 -5.1166347e-05 2830
10800 1080 0.00033161645 -2.0773965e-05 2830
11000 1100 0.00032913022 2.2384421e-07 2830
11200 1120 0.00032618376 1.2304773e-05 2830
11400 1140 0.00032310409 1.3725982e-05 2830
11600 1160 0.0003202128 9.0431945e-06 2830
11800 1180 0.00031760386 -5.3537879e-07 2830
12000 1200 0.00031518884 -1.331708e-05 2830
12200 1220 0.00031283958 -3.0838612e-05 2830
12400 1240 0.0003104901 -5.0038548e-05 2830
12600 1260 0.00030811597 -6.9699925e-05 2830
12800 1280 0.00030555782 -8.9972287e-05 2830
13000 1300 0.00030256671 -0.00011712941 2830
13200 1320 0.00029907961 -0.00015495826 2830
13400 1340 0.00029504656 -0.00020292633 2830
13600 1360 0.0002905184 -0.00024892421 2830
13800 1380 0.00028564542 -0.000295085 2830
14000 1400 0.00028073246 -0.00034571956 2830
14200 1420 0.00027611457 -0.00039341977 2830
14400 1440 0.00027217382 -0.0004281012 2830
14600 1460 0.00026919129 -0.00045342545 2830
14800 1480 0.00026727674 -0.00047323419 2830
15000 1500 0.0002663482 -0.00048423944 2830
15200 1520 0.00026616663 -0.0004816085 2830
15400 1540 0.00026634862 -0.00047573486 2830
15600 1560 0.0002664314 -0.00046803192 2830
15800 1580 0.00026603348 -0.00045753668 2830
16000 1600 0.00026511015 -0.00044676105 2830
16200 1620 0.00026373403 -0.00044075794 2830
16400 1640 0.00026217342 -0.00043684036 2830
16600 1660 0.0002607038 -0.00042774771 2830
16800 1680 0.00025951097 -0.00041603026 2830
17000 1700 0.00025869088 -0.00040302996 2830
17200 1720 0.00025825588 -0.00038415247 2830
17400 1740 0.00025818373 -0.00035742127 2830
17600 1760 0.00025843381 -0.00032854722 2830
17800 1780 0.00025897836 -0.00029821183 2830
18000 1800 0.00025981472 -0.00026108907 2830
18200 1820 0.00026095775 -0.00021731058 2830
18400 1840 0.00026239688 -0.00017030825 2830
18600 1860 0.00026404432 -0.00011868778 2830
18800 1880 0.00026574247 -5.9556286e-05 2830
19000 1900 0.00026729563 2.3014881e-06 2830
19200 1920 0.00026852418 6.2100169e-05 2830
19400 1940 0.00026929086 0.00012090325 2830
19600 1960 0.0002695407 0.00017904223 2830
19800 1980 0.00026929677 0.00023112254 2830
20000 2000 0.00026863577 0.0002756697 2830
20200 2020 0.00026765699 0.0003158399 2830
20400 2040 0.00026646841 0.00035200747 2830
20600 2060 0.00026516938 0.00038018442 2830
20800 2080 0.00026383495 0.00040179111 2830
21000 2100 0.00026252489 0.00042030921 2830
21200 2120 0.00026128616 0.00043466976 2830
21400 2140 0.00026014896 0.00044221445 2830
21600 2160 0.00025912325 0.00044531883 2830
21800 2180 0.00025821515 0.00044661709 2830
22000 2200 0.00025742576 0.00044409089 2830
22200 2220 0.00025674938 0.00043634999 2830
22400 2240 0.00025617111 0.00042630344 2830
22600 2260 0.0002556791 0.00041561603 2830
22800 2280 0.00025525963 0.00040166735 2830
23000 2300 0.00025489538 0.00038430419 2830
23200 2320 0.00025456861 0.0003669402 2830
23400 2340 0.00025426747 0.00034972373 2830
23600 2360 0.00025398353 0.0003302242 2830
23800 2380 0.00025370842 0.00030993088 2830
24000 2400 0.00025344084 0.00029143258 2830
24200 2420 0.00025318683 0.00027421708 2830
24400 2440 0.0002529591 0.00025603123 2830
24600 2460 0.0002527713 0.00023950245 2830
24800 2480 0.00025264228 0.00022644812 2830
25000 2500 0.00025259021 0.00021540748 2830
25200 2520 0.00025262892 0.00020544201 2830
25400 2540 0.00025276229 0.00019845807 2830
25600 2560 0.0002529876 0.00019449958 2830
25800 2580 0.00025329374 0.00019082606 2830
26000 2600 0.00025366066 0.00018700064 2830
26200 2620 0.00025406164 0.00018426061 2830
26400 2640 0.00025446737 0.00018098339 2830
26600 2660 0.00025484714 0.00017471869 2830
26800 2680 0.00025516604 0.00016565557 2830
27000 2700 0.00025538911 0.00015493626 2830
27200 2720 0.00025548177 0.00014075592 2830
27400 2740 0.00025541168 0.00012205573 2830
27600 2760 0.00025514889 0.00010039772 2830
27800 2780 0.00025467547 7.7069215e-05 2830
28000 2800 0.0002539915 5.1158042e-05 2830
28200 2820 0.00025312083 2.3468384e-05 2830
28400 2840 0.00025211323 -3.2184465e-06 2830
28600 2860 0.00025104366 -2.7726301e-05 2830
28800 2880 0.00025000263 -5.0202987e-05 2830
29000 2900 0.00024907814 -6.9244776e-05 2830
29200 2920 0.00024833815 -8.2874516e-05 2830
29400 2940 0.0002478155 -9.1854992e-05 2830
29600 2960 0.00024750313 -9.766055e-05 2830
29800 2980 0.00024735538 -9.9681291e-05 2830
30000 3000 0.00024730191 -9.818759e-05 2830
Loop time of 177.982 on 4 procs for 30000 steps with 2830 atoms
Performance: 1456330.235 tau/day, 168.557 timesteps/s, 477.016 katom-step/s
99.7% CPU use with 4 MPI tasks x no OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 22.913 | 27.061 | 34.594 | 87.2 | 15.20
Bond | 0.22386 | 0.26159 | 0.30792 | 6.0 | 0.15
Neigh | 0.84412 | 0.84509 | 0.8462 | 0.1 | 0.47
Comm | 0.50015 | 0.55579 | 0.60346 | 5.2 | 0.31
Output | 0.65854 | 0.69412 | 0.72473 | 2.8 | 0.39
Modify | 133.13 | 136 | 137.38 | 14.5 | 76.41
Other | | 12.57 | | | 7.06
Nlocal: 707.5 ave 1576 max 53 min
Histogram: 2 0 0 0 0 0 1 0 0 1
Nghost: 164.75 ave 239 max 94 min
Histogram: 1 0 1 0 0 0 0 1 0 1
Neighs: 12307.8 ave 27380 max 983 min
Histogram: 2 0 0 0 0 0 1 0 0 1
FullNghs: 23517 ave 53040 max 1502 min
Histogram: 2 0 0 0 0 0 1 0 0 1
Total # of neighbors = 94068
Ave neighs/atom = 33.239576
Ave special neighs/atom = 0.89257951
Neighbor list builds = 783
Dangerous builds = 0
Total wall time: 0:02:58

1694
examples/rheo/dam-break/log.17Apr2024.dam.g++.4 Normal file
View File

File diff suppressed because it is too large Load Diff

3
examples/rheo/ice-cubes/in.rheo.ice.cubes
View File

@ -77,5 +77,6 @@ compute nbond all nbond/atom
thermo 200
thermo_style custom step time ke press atoms
dump 1 all custom 200 atomDump id type x y vx vy fx fy c_phase c_temp c_eng c_nbond c_rho
#dump 1 all custom 200 atomDump id type x y vx vy fx fy c_phase c_temp c_eng c_nbond c_rho
run 30000

379
examples/rheo/ice-cubes/log.17Apr2024.ice.g++.4 Normal file
View File

@ -0,0 +1,379 @@
LAMMPS (17 Apr 2024 - Development - patch_5May2020-18508-g3c0eaf6870-modified)
# ------ 2D Ice Cube Pour ------ #
dimension 2
units lj
atom_style hybrid rheo/thermal bond
boundary m m p
comm_modify vel yes
newton off
special_bonds lj 0.0 1.0 1.0 coul 1.0 1.0 1.0
region box block -25 25 0 100 -0.01 0.01 units box
create_box 1 box bond/types 1 extra/bond/per/atom 15 extra/special/per/atom 50
Created orthogonal box = (-25 0 -0.01) to (25 100 0.01)
2 by 2 by 1 MPI processor grid
region fluid block $(xlo+1) $(xhi-1) $(ylo+1) $(ylo+30) EDGE EDGE units box
region fluid block -24 $(xhi-1) $(ylo+1) $(ylo+30) EDGE EDGE units box
region fluid block -24 24 $(ylo+1) $(ylo+30) EDGE EDGE units box
region fluid block -24 24 1 $(ylo+30) EDGE EDGE units box
region fluid block -24 24 1 30 EDGE EDGE units box
lattice sq 1.0
Lattice spacing in x,y,z = 1 1 1
create_atoms 1 region fluid
Created 1470 atoms
using lattice units in orthogonal box = (-25 0 -0.01) to (25 100 0.01)
create_atoms CPU = 0.001 seconds
set group all sph/e 8.0
Setting atom values ...
1470 settings made for sph/e
# ------ Model parameters ------#
variable cut equal 3.0
variable n equal 1.0
variable rho0 equal 1.0
variable cs equal 1.0
variable mp equal ${rho0}/${n}
variable mp equal 1/${n}
variable mp equal 1/1
variable zeta equal 0.05
variable kappa equal 0.01*${rho0}/${mp}
variable kappa equal 0.01*1/${mp}
variable kappa equal 0.01*1/1
variable dt_max equal 0.1*${cut}/${cs}/3
variable dt_max equal 0.1*3/${cs}/3
variable dt_max equal 0.1*3/1/3
variable eta equal 0.05
variable Cv equal 1.0
variable L equal 1.0
variable Tf equal 1.0
mass * ${mp}
mass * 1
timestep 0.1
pair_style hybrid/overlay rheo ${cut} artificial/visc ${zeta} rheo/solid
pair_style hybrid/overlay rheo 3 artificial/visc ${zeta} rheo/solid
pair_style hybrid/overlay rheo 3 artificial/visc 0.05 rheo/solid
pair_coeff * * rheo
pair_coeff * * rheo/solid 1.0 1.0 1.0
bond_style bpm/spring
bond_coeff 1 1.0 1.0 1.0
# ------ Pour particles ------#
molecule my_mol "square.mol"
Read molecule template my_mol:
#Made with create_mol.py
1 molecules
0 fragments
100 atoms with max type 1
342 bonds with max type 1
0 angles with max type 0
0 dihedrals with max type 0
0 impropers with max type 0
# Wall region extends far enough in z to avoid contact
region wall block EDGE EDGE EDGE EDGE -5 5 side in open 4 units box
region drop block -16 16 70 90 EDGE EDGE side in units box
fix 1 all rheo ${cut} quintic 0 thermal shift surface/detection coordination 22 8
fix 1 all rheo 3 quintic 0 thermal shift surface/detection coordination 22 8
fix 2 all rheo/viscosity * constant ${eta}
fix 2 all rheo/viscosity * constant 0.05
fix 3 all rheo/pressure * linear
fix 4 all rheo/thermal conductivity * constant ${kappa} specific/heat * constant ${Cv} Tfreeze * constant ${Tf} latent/heat * constant ${L} react 1.5 1
fix 4 all rheo/thermal conductivity * constant 0.01 specific/heat * constant ${Cv} Tfreeze * constant ${Tf} latent/heat * constant ${L} react 1.5 1
fix 4 all rheo/thermal conductivity * constant 0.01 specific/heat * constant 1 Tfreeze * constant ${Tf} latent/heat * constant ${L} react 1.5 1
fix 4 all rheo/thermal conductivity * constant 0.01 specific/heat * constant 1 Tfreeze * constant 1 latent/heat * constant ${L} react 1.5 1
fix 4 all rheo/thermal conductivity * constant 0.01 specific/heat * constant 1 Tfreeze * constant 1 latent/heat * constant 1 react 1.5 1
fix 5 all wall/region wall harmonic 1.0 1.0 1.0
fix 6 all gravity 5e-4 vector 0 -1 0
fix 7 all deposit 8 0 1000 37241459 mol my_mol region drop near 2.0 vy -0.02 -0.02
WARNING: Molecule attributes do not match system attributes (../molecule.cpp:1881)
fix 8 all enforce2d
compute rho all rheo/property/atom rho
compute phase all rheo/property/atom phase
compute temp all rheo/property/atom temperature
compute eng all rheo/property/atom energy
compute nbond all nbond/atom
# ------ Output & Run ------ #
thermo 200
thermo_style custom step time ke press atoms
dump 1 all custom 200 atomDump id type x y vx vy fx fy c_phase c_temp c_eng c_nbond c_rho
run 30000
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Your simulation uses code contributions which should be cited:
- BPM bond style: doi:10.1039/D3SM01373A
@Article{Clemmer2024,
author = {Clemmer, Joel T. and Monti, Joseph M. and Lechman, Jeremy B.},
title = {A soft departure from jamming: the compaction of deformable
granular matter under high pressures},
journal = {Soft Matter},
year = 2024,
volume = 20,
number = 8,
pages = {1702--1718}
}
- @article{PalermoInPrep,
journal = {in prep},
title = {RHEO: A Hybrid Mesh-Free Model Framework for Dynamic Multi-Phase Flows},
year = {2024},
author = {Eric T. Palermo and Ki T. Wolf and Joel T. Clemmer and Thomas C. O'Connor},
}
- @article{ApplMathModel.130.310,
title = {A hybrid smoothed-particle hydrodynamics model of oxide skins on molten aluminum},
journal = {Applied Mathematical Modelling},
volume = {130},
pages = {310-326},
year = {2024},
issn = {0307-904X},
doi = {https://doi.org/10.1016/j.apm.2024.02.027},
author = {Joel T. Clemmer and Flint Pierce and Thomas C. O'Connor and Thomas D. Nevins and Elizabeth M.C. Jones and Jeremy B. Lechman and John Tencer},
}
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Generated 0 of 0 mixed pair_coeff terms from geometric mixing rule
Neighbor list info ...
update: every = 1 steps, delay = 0 steps, check = yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 3.3
ghost atom cutoff = 3.3
binsize = 1.65, bins = 31 61 1
7 neighbor lists, perpetual/occasional/extra = 6 1 0
(1) pair rheo, perpetual, half/full from (3)
attributes: half, newton off
pair build: halffull/newtoff
stencil: none
bin: none
(2) pair rheo/solid, perpetual, trim from (4)
attributes: half, newton off, cut 1.3
pair build: trim
stencil: none
bin: none
(3) compute RHEO/KERNEL, perpetual
attributes: full, newton off
pair build: full/bin
stencil: full/bin/2d
bin: standard
(4) compute RHEO/GRAD, perpetual, copy from (1)
attributes: half, newton off
pair build: copy
stencil: none
bin: none
(5) compute RHEO/VSHIFT, perpetual, copy from (1)
attributes: half, newton off
pair build: copy
stencil: none
bin: none
(6) compute RHEO/SURFACE, perpetual, copy from (1)
attributes: half, newton off
pair build: copy
stencil: none
bin: none
(7) fix rheo/thermal, occasional, trim from (4)
attributes: half, newton off, cut 3
pair build: trim
stencil: none
bin: none
Per MPI rank memory allocation (min/avg/max) = 15.53 | 15.61 | 15.69 Mbytes
Step Time KinEng Press Atoms
0 0 0 0 1470
200 20 5.6002982e-05 3.4434234e-05 1570
400 40 8.2173099e-05 8.6171768e-05 1570
600 60 8.019018e-05 0.00010750355 1570
800 80 0.00013866953 0.00010265608 1570
1000 100 0.00018965028 8.1985605e-05 1570
1200 120 0.00022033242 7.4736443e-05 1670
1400 140 0.00030767062 0.00011264333 1670
1600 160 0.00040770127 0.00018779992 1670
1800 180 0.00047476332 0.00023153009 1670
2000 200 0.00059116774 0.00027200445 1670
2200 220 0.0007151733 0.0002919963 1770
2400 240 0.00083392135 0.00029757889 1770
2600 260 0.00099653466 0.00036547269 1770
2800 280 0.0011964069 0.00045983458 1770
3000 300 0.0013716953 0.00055013647 1770
3200 320 0.0015174096 0.00064203572 1870
3400 340 0.0016539743 0.00086671622 1870
3600 360 0.0015887858 0.00066353749 1870
3800 380 0.0016451684 0.00070551483 1870
4000 400 0.0017330971 0.00080722283 1870
4200 420 0.001812193 0.00073573903 1970
4400 440 0.001755871 0.0010621909 1970
4600 460 0.0016190772 0.00072913706 1970
4800 480 0.0015741931 0.00073524088 1970
5000 500 0.0016488815 0.00088684275 1970
5200 520 0.0017213288 0.00077042378 2070
5400 540 0.0018509598 0.0010219434 2070
5600 560 0.0020251064 0.00083182483 2070
5800 580 0.0022473255 0.00095076144 2070
6000 600 0.0024843519 0.0011247014 2070
6200 620 0.0022282321 0.0018105932 2170
6400 640 0.0020289063 0.0014158497 2170
6600 660 0.002145241 0.0011359383 2170
6800 680 0.0024313937 0.0016475504 2170
7000 700 0.0021000599 0.0020983745 2170
7200 720 0.0019137235 0.0010439152 2270
7400 740 0.0018801367 0.00095436448 2270
7600 760 0.0017979449 0.0011184039 2270
7800 780 0.0018005205 0.0009243205 2270
8000 800 0.0017827073 0.0013671228 2270
8200 820 0.0018387108 0.0015426012 2270
8400 840 0.0016000788 0.0016751514 2270
8600 860 0.0013954964 0.0016884335 2270
8800 880 0.0013283728 0.0012399398 2270
9000 900 0.001389385 0.0012968496 2270
9200 920 0.0012295438 0.0012995821 2270
9400 940 0.0010522655 0.00082245528 2270
9600 960 0.00097085496 0.00053833131 2270
9800 980 0.0009398987 0.00063467387 2270
10000 1000 0.00092710392 0.00059494446 2270
10200 1020 0.00095545471 0.00074560644 2270
10400 1040 0.0009645841 0.00085429807 2270
10600 1060 0.00064037148 0.0017222246 2270
10800 1080 0.00046790978 0.00088204234 2270
11000 1100 0.00030106229 0.00074950209 2270
11200 1120 0.00027746016 0.00052831745 2270
11400 1140 0.0002533348 0.0006272715 2270
11600 1160 0.00021825085 0.00029691552 2270
11800 1180 0.0001451308 0.00015037478 2270
12000 1200 0.0001314823 0.00017227174 2270
12200 1220 0.00013693632 0.00017791384 2270
12400 1240 0.00014987347 0.0002286677 2270
12600 1260 0.00015092598 0.0003698436 2270
12800 1280 0.0001291653 0.00047229532 2270
13000 1300 0.00011949988 0.00049560375 2270
13200 1320 0.00011694665 0.00057542084 2270
13400 1340 9.6164519e-05 0.00062714755 2270
13600 1360 8.4517591e-05 0.00044156913 2270
13800 1380 0.00019140516 0.0003264745 2270
14000 1400 0.00013868599 0.00037753497 2270
14200 1420 9.3701636e-05 0.00031517848 2270
14400 1440 6.7389077e-05 0.0002946861 2270
14600 1460 5.3640086e-05 0.00026650711 2270
14800 1480 4.2699992e-05 0.00023789279 2270
15000 1500 5.3012016e-05 0.00019933234 2270
15200 1520 5.8834197e-05 0.00022407007 2270
15400 1540 5.0899982e-05 0.00029695531 2270
15600 1560 3.0476742e-05 0.00039119066 2270
15800 1580 1.6633264e-05 0.00033770401 2270
16000 1600 1.098906e-05 0.00036684894 2270
16200 1620 1.464848e-05 0.00036449759 2270
16400 1640 1.9598429e-05 0.00021056689 2270
16600 1660 1.2644955e-05 0.00020781781 2270
16800 1680 8.8428553e-06 0.000165 2270
17000 1700 8.8971439e-06 0.00012266475 2270
17200 1720 1.7032781e-05 0.00019873443 2270
17400 1740 1.9448563e-05 0.00025661663 2270
17600 1760 1.3714713e-05 0.000324022 2270
17800 1780 9.1326468e-06 0.00031392513 2270
18000 1800 9.2464802e-06 0.00029729527 2270
18200 1820 1.5553042e-05 0.00027488475 2270
18400 1840 1.4132933e-05 0.00019565459 2270
18600 1860 9.4734832e-06 0.00016716988 2270
18800 1880 5.5115145e-06 0.00013728033 2270
19000 1900 8.268812e-06 0.00015119605 2270
19200 1920 1.2470136e-05 0.00020222131 2270
19400 1940 9.9387775e-06 0.00024503373 2270
19600 1960 5.4241999e-06 0.00026921858 2270
19800 1980 2.7987348e-06 0.00026201267 2270
20000 2000 6.272538e-06 0.00025626323 2270
20200 2020 8.0157781e-06 0.000220139 2270
20400 2040 6.1652093e-06 0.00017089058 2270
20600 2060 2.9967592e-06 0.00014582864 2270
20800 2080 3.016678e-06 0.000148629 2270
21000 2100 7.287645e-06 0.00016486102 2270
21200 2120 8.6905277e-06 0.00020276916 2270
21400 2140 6.8453018e-06 0.00023156153 2270
21600 2160 3.3853799e-06 0.0002432462 2270
21800 2180 4.1241209e-06 0.00022829024 2270
22000 2200 7.0802396e-06 0.00020784823 2270
22200 2220 7.3361691e-06 0.00018114134 2270
22400 2240 5.0764593e-06 0.00014351106 2270
22600 2260 2.7487537e-06 0.00012919872 2270
22800 2280 4.620167e-06 0.00013746218 2270
23000 2300 6.9819357e-06 0.00015985102 2270
23200 2320 6.8923916e-06 0.00018713045 2270
23400 2340 4.1795088e-06 0.00019846682 2270
23600 2360 2.2871028e-06 0.00021068421 2270
23800 2380 3.862046e-06 0.00019553306 2270
24000 2400 5.2448555e-06 0.00017398041 2270
24200 2420 4.7565441e-06 0.00015008142 2270
24400 2440 2.2952135e-06 0.00012747106 2270
24600 2460 2.1575617e-06 0.00012516996 2270
24800 2480 4.1777868e-06 0.0001331902 2270
25000 2500 5.5679133e-06 0.00015504562 2270
25200 2520 4.5758741e-06 0.00017146032 2270
25400 2540 2.3403277e-06 0.00017611666 2270
25600 2560 2.7029302e-06 0.00016850788 2270
25800 2580 4.3601102e-06 0.00015884642 2270
26000 2600 5.2244249e-06 0.00013793898 2270
26200 2620 3.4577672e-06 0.00012395875 2270
26400 2640 2.361577e-06 0.00011600057 2270
26600 2660 2.8515644e-06 0.00011277063 2270
26800 2680 4.0851213e-06 0.0001290832 2270
27000 2700 4.2579644e-06 0.0001476495 2270
27200 2720 2.6593858e-06 0.00015977745 2270
27400 2740 1.990115e-06 0.00015612787 2270
27600 2760 2.6756835e-06 0.00014913772 2270
27800 2780 3.9032806e-06 0.00014014763 2270
28000 2800 3.2729446e-06 0.00012216846 2270
28200 2820 1.9357278e-06 0.00011078621 2270
28400 2840 1.7094832e-06 0.00010910509 2270
28600 2860 2.8731406e-06 0.00011179644 2270
28800 2880 3.7062354e-06 0.00012254091 2270
29000 2900 2.7844262e-06 0.00013060331 2270
29200 2920 1.7680655e-06 0.00013797514 2270
29400 2940 1.706873e-06 0.0001350685 2270
29600 2960 2.8764562e-06 0.00012428508 2270
29800 2980 3.1502029e-06 0.00011456718 2270
30000 3000 2.1833409e-06 0.00010317469 2270
Loop time of 165.611 on 4 procs for 30000 steps with 2270 atoms
Performance: 1565111.240 tau/day, 181.147 timesteps/s, 411.204 katom-step/s
99.7% CPU use with 4 MPI tasks x no OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 0.63183 | 21.226 | 42.266 | 444.6 | 12.82
Bond | 0.095073 | 0.17799 | 0.27877 | 17.0 | 0.11
Neigh | 2.0745 | 2.0781 | 2.0822 | 0.2 | 1.25
Comm | 0.32024 | 0.38703 | 0.45564 | 8.1 | 0.23
Output | 0.60459 | 0.76798 | 0.93724 | 18.6 | 0.46
Modify | 119.85 | 140.76 | 161.36 | 172.2 | 85.00
Other | | 0.2124 | | | 0.13
Nlocal: 567.5 ave 1139 max 0 min
Histogram: 2 0 0 0 0 0 0 0 0 2
Nghost: 75.5 ave 152 max 0 min
Histogram: 2 0 0 0 0 0 0 0 0 2
Neighs: 9238.25 ave 18490 max 0 min
Histogram: 2 0 0 0 0 0 0 0 0 2
FullNghs: 17945 ave 35917 max 0 min
Histogram: 2 0 0 0 0 0 0 0 0 2
Total # of neighbors = 71780
Ave neighs/atom = 31.621145
Ave special neighs/atom = 0.22026432
Neighbor list builds = 2071
Dangerous builds = 0
Total wall time: 0:02:45

3
examples/rheo/oxidation/in.rheo.oxidation
View File

@ -98,5 +98,6 @@ compute nbond_solid all nbond/atom bond/type 1
thermo 200
thermo_style custom step time ke press atoms
dump 1 all custom 200 atomDump id type x y vx vy fx fy c_phase c_temp c_eng c_nbond_solid c_nbond_shell c_rho c_surf c_status
#dump 1 all custom 200 atomDump id type x y vx vy fx fy c_phase c_temp c_eng c_nbond_solid c_nbond_shell c_rho c_surf c_status
run 40000

488
examples/rheo/oxidation/log.17Apr2024.oxidation.g++.4 Normal file
View File

@ -0,0 +1,488 @@
LAMMPS (17 Apr 2024 - Development - patch_5May2020-18508-g3c0eaf6870-modified)
# ------ 2D oxidizing bar ------ #
dimension 2
units lj
atom_style hybrid rheo/thermal bond
boundary m m p
comm_modify vel yes
newton off
region box block -60 60 0 80 -0.01 0.01 units box
create_box 3 box bond/types 2 extra/bond/per/atom 20 extra/special/per/atom 50
Created orthogonal box = (-60 0 -0.01) to (60 80 0.01)
2 by 2 by 1 MPI processor grid
region lbar block -15 0 3 80 EDGE EDGE units box
region rbar block 0 15 3 80 EDGE EDGE units box
region bar union 2 lbar rbar
region floor block EDGE EDGE EDGE 3.0 EDGE EDGE units box
lattice hex 1.0
Lattice spacing in x,y,z = 1.0745699 1.8612097 1.0745699
create_atoms 1 region bar
Created 2255 atoms
using lattice units in orthogonal box = (-60 0 -0.01) to (60 80 0.01)
create_atoms CPU = 0.001 seconds
create_atoms 3 region floor
Created 446 atoms
using lattice units in orthogonal box = (-60 0 -0.01) to (60 80 0.01)
create_atoms CPU = 0.000 seconds
set region rbar type 2
Setting atom values ...
1148 settings made for type
group bar type 1 2
2255 atoms in group bar
group rbar type 2
1148 atoms in group rbar
group floor type 3
446 atoms in group floor
set group all sph/e 0.0
Setting atom values ...
2701 settings made for sph/e
set group all rheo/status 1
Setting atom values ...
2701 settings made for rheo/status
# ------ Model parameters ------#
variable cut equal 3.0
variable n equal 1.0
variable rho0 equal 1.0
variable cs equal 1.0
variable mp equal ${rho0}/${n}
variable mp equal 1/${n}
variable mp equal 1/1
variable zeta equal 0.05
variable kappa equal 0.1*${rho0}/${mp}
variable kappa equal 0.1*1/${mp}
variable kappa equal 0.1*1/1
variable dt_max equal 0.1*${cut}/${cs}/3
variable dt_max equal 0.1*3/${cs}/3
variable dt_max equal 0.1*3/1/3
variable eta equal 0.05
variable Cv equal 1.0
variable L equal 0.1
variable Tf equal 1.0
mass * ${mp}
mass * 1
timestep 0.1
pair_style hybrid/overlay rheo ${cut} artificial/visc ${zeta} rheo/solid
pair_style hybrid/overlay rheo 3 artificial/visc ${zeta} rheo/solid
pair_style hybrid/overlay rheo 3 artificial/visc 0.05 rheo/solid
pair_coeff * * rheo
pair_coeff * * rheo/solid 1.0 1.0 1.0
special_bonds lj 0.0 1.0 1.0 coul 0.0 1.0 1.0
Finding 1-2 1-3 1-4 neighbors ...
special bond factors lj: 0 1 1
special bond factors coul: 0 1 1
0 = max # of 1-2 neighbors
101 = max # of special neighbors
special bonds CPU = 0.000 seconds
create_bonds many bar bar 1 0 1.5
Generated 0 of 3 mixed pair_coeff terms from geometric mixing rule
Neighbor list info ...
update: every = 1 steps, delay = 0 steps, check = yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 3.3
ghost atom cutoff = 3.3
binsize = 1.65, bins = 73 49 1
3 neighbor lists, perpetual/occasional/extra = 2 1 0
(1) command create_bonds, occasional
attributes: full, newton off
pair build: full/bin
stencil: full/bin/2d
bin: standard
(2) pair rheo, perpetual
attributes: half, newton off
pair build: half/bin/newtoff
stencil: full/bin/2d
bin: standard
(3) pair rheo/solid, perpetual, trim from (2)
attributes: half, newton off, cut 1.3
pair build: trim
stencil: none
bin: none
Added 6547 bonds, new total = 6547
Finding 1-2 1-3 1-4 neighbors ...
special bond factors lj: 0 1 1
special bond factors coul: 0 1 1
6 = max # of 1-2 neighbors
101 = max # of special neighbors
special bonds CPU = 0.000 seconds
special_bonds lj 0.0 1.0 1.0 coul 1.0 1.0 1.0
bond_style hybrid bpm/spring rheo/shell t/form 100
bond_coeff 1 bpm/spring 1.0 1.0 1.0
bond_coeff 2 rheo/shell 0.2 0.2 0.1
# ------ Apply dynamics ------#
# Note: surface detection is not performed on solid bodies, so cannot use surface property
compute coord all rheo/property/atom coordination
variable surf atom c_coord<22
group surf dynamic all var surf every 10
dynamic group surf defined
fix 1 all rheo ${cut} quintic 0 thermal shift surface/detection coordination 22 8
fix 1 all rheo 3 quintic 0 thermal shift surface/detection coordination 22 8
fix 2 all rheo/viscosity * constant ${eta}
fix 2 all rheo/viscosity * constant 0.05
fix 3 all rheo/pressure * linear
fix 4 all rheo/thermal conductivity * constant ${kappa} specific/heat * constant ${Cv} Tfreeze * constant ${Tf} latent/heat * constant ${L} react 1.5 1
fix 4 all rheo/thermal conductivity * constant 0.1 specific/heat * constant ${Cv} Tfreeze * constant ${Tf} latent/heat * constant ${L} react 1.5 1
fix 4 all rheo/thermal conductivity * constant 0.1 specific/heat * constant 1 Tfreeze * constant ${Tf} latent/heat * constant ${L} react 1.5 1
fix 4 all rheo/thermal conductivity * constant 0.1 specific/heat * constant 1 Tfreeze * constant 1 latent/heat * constant ${L} react 1.5 1
fix 4 all rheo/thermal conductivity * constant 0.1 specific/heat * constant 1 Tfreeze * constant 1 latent/heat * constant 0.1 react 1.5 1
fix 5 rbar rheo/oxidation 1.5 2 1.0
fix 6 all wall/harmonic ylo EDGE 2.0 1.0 1.0
fix 7 all gravity 5e-5 vector 0 -1 0
fix 8 floor setforce 0.0 0.0 0.0
fix 9 surf add/heat linear 1.1 0.05
fix 10 floor add/heat constant 0 overwrite yes # fix the temperature of the floor
fix 11 all enforce2d
compute surf all rheo/property/atom surface
compute rho all rheo/property/atom rho
compute phase all rheo/property/atom phase
compute status all rheo/property/atom status
compute temp all rheo/property/atom temperature
compute eng all rheo/property/atom energy
compute nbond_shell all rheo/property/atom nbond/shell
compute nbond_solid all nbond/atom bond/type 1
# ------ Output & Run ------ #
thermo 200
thermo_style custom step time ke press atoms
dump 1 all custom 200 atomDump id type x y vx vy fx fy c_phase c_temp c_eng c_nbond_solid c_nbond_shell c_rho c_surf c_status
run 40000
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Your simulation uses code contributions which should be cited:
- BPM bond style: doi:10.1039/D3SM01373A
@Article{Clemmer2024,
author = {Clemmer, Joel T. and Monti, Joseph M. and Lechman, Jeremy B.},
title = {A soft departure from jamming: the compaction of deformable
granular matter under high pressures},
journal = {Soft Matter},
year = 2024,
volume = 20,
number = 8,
pages = {1702--1718}
}
- @article{PalermoInPrep,
journal = {in prep},
title = {RHEO: A Hybrid Mesh-Free Model Framework for Dynamic Multi-Phase Flows},
year = {2024},
author = {Eric T. Palermo and Ki T. Wolf and Joel T. Clemmer and Thomas C. O'Connor},
}
- @article{ApplMathModel.130.310,
title = {A hybrid smoothed-particle hydrodynamics model of oxide skins on molten aluminum},
journal = {Applied Mathematical Modelling},
volume = {130},
pages = {310-326},
year = {2024},
issn = {0307-904X},
doi = {https://doi.org/10.1016/j.apm.2024.02.027},
author = {Joel T. Clemmer and Flint Pierce and Thomas C. O'Connor and Thomas D. Nevins and Elizabeth M.C. Jones and Jeremy B. Lechman and John Tencer},
}
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Generated 0 of 3 mixed pair_coeff terms from geometric mixing rule
Neighbor list info ...
update: every = 1 steps, delay = 0 steps, check = yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 3.3
ghost atom cutoff = 3.3
binsize = 1.65, bins = 73 49 1
8 neighbor lists, perpetual/occasional/extra = 7 1 0
(1) pair rheo, perpetual, half/full from (3)
attributes: half, newton off
pair build: halffull/newtoff
stencil: none
bin: none
(2) pair rheo/solid, perpetual, trim from (4)
attributes: half, newton off, cut 1.3
pair build: trim
stencil: none
bin: none
(3) compute RHEO/KERNEL, perpetual
attributes: full, newton off
pair build: full/bin
stencil: full/bin/2d
bin: standard
(4) compute RHEO/GRAD, perpetual, copy from (1)
attributes: half, newton off
pair build: copy
stencil: none
bin: none
(5) compute RHEO/VSHIFT, perpetual, copy from (1)
attributes: half, newton off
pair build: copy
stencil: none
bin: none
(6) compute RHEO/SURFACE, perpetual, copy from (1)
attributes: half, newton off
pair build: copy
stencil: none
bin: none
(7) fix rheo/thermal, occasional, trim from (4)
attributes: half, newton off, cut 3
pair build: trim
stencil: none
bin: none
(8) fix rheo/oxidation, perpetual, trim from (3)
attributes: full, newton off, cut 1.8
pair build: trim
stencil: none
bin: none
Per MPI rank memory allocation (min/avg/max) = 25.96 | 25.96 | 25.96 Mbytes
Step Time KinEng Press Atoms
0 0 0 0 2701
200 20 4.1743799e-07 1.1743617e-07 2701
400 40 1.6697519e-06 4.6974469e-07 2701
600 60 3.7127333e-06 1.0646825e-05 2701
800 80 4.6683656e-06 0.00015182605 2701
1000 100 4.7368707e-06 0.00028128761 2701
1200 120 3.4384322e-06 0.00045913378 2701
1400 140 1.4119866e-06 0.00055627091 2701
1600 160 4.4114517e-07 0.00058247308 2701
1800 180 4.8289229e-07 0.0005510948 2701
2000 200 1.8494183e-06 0.00048386222 2701
2200 220 3.3319816e-06 0.00037903264 2701
2400 240 3.8128922e-06 0.00024115906 2701
2600 260 3.1943401e-06 9.727407e-05 2701
2800 280 1.6172816e-06 -2.632162e-05 2701
3000 300 3.6100709e-07 -8.5761867e-05 2701
3200 320 1.4745502e-07 -5.9204127e-05 2701
3400 340 8.3369782e-07 8.8312464e-07 2701
3600 360 2.0484052e-06 5.8521477e-05 2701
3800 380 3.1639387e-06 0.0001685663 2701
4000 400 3.1692907e-06 0.00026875988 2701
4200 420 2.391933e-06 0.00038621787 2701
4400 440 1.1964404e-06 0.00048901286 2701
4600 460 4.0508824e-07 0.00051863639 2701
4800 480 5.4908507e-07 0.00049263754 2701
5000 500 1.3139665e-06 0.00041984264 2701
5200 520 2.1939161e-06 0.00033095351 2701
5400 540 2.3687031e-06 0.00022422981 2701
5600 560 1.8280882e-06 0.00011544328 2701
5800 580 8.8610517e-07 2.9307791e-05 2701
6000 600 2.0989359e-07 -1.7340941e-05 2701
6200 620 2.8658301e-07 -8.1237835e-06 2701
6400 640 9.7636239e-07 4.3755922e-05 2701
6600 660 1.891303e-06 0.0001185719 2701
6800 680 2.4149904e-06 0.00020830273 2701
7000 700 2.3174953e-06 0.00030114767 2701
7200 720 1.7918612e-06 0.00037821537 2701
7400 740 1.2114987e-06 0.0004233475 2701
7600 760 9.9661553e-07 0.00042958263 2701
7800 780 1.1552559e-06 0.00039944618 2701
8000 800 1.5249138e-06 0.00034034478 2701
8200 820 1.7453861e-06 0.00026826463 2701
8400 840 1.6259021e-06 0.00019131768 2701
8600 860 1.2612805e-06 0.0001162957 2701
8800 880 8.6964518e-07 7.1771506e-05 2701
9000 900 7.6892472e-07 5.6170687e-05 2701
9200 920 1.0780045e-06 7.1925995e-05 2701
9400 940 1.6514902e-06 0.00011635293 2701
9600 960 2.1891377e-06 0.00017599885 2701
9800 980 2.4551701e-06 0.00024127934 2701
10000 1000 2.4277051e-06 0.00029918622 2701
10200 1020 2.2655987e-06 0.00034067996 2701
10400 1040 2.1767207e-06 0.00035598133 2701
10600 1060 2.2796719e-06 0.00034359076 2701
10800 1080 2.4884225e-06 0.00030749714 2701
11000 1100 2.6387215e-06 0.00025725198 2701
11200 1120 2.5968908e-06 0.00020170699 2701
11400 1140 2.4108931e-06 0.00015185858 2701
11600 1160 2.2375166e-06 0.00011800349 2701
11800 1180 2.2407196e-06 0.00010646971 2701
12000 1200 2.4845263e-06 0.00011817498 2701
12200 1220 2.8733204e-06 0.00015013186 2701
12400 1240 3.2437087e-06 0.00019211975 2701
12600 1260 3.4732728e-06 0.00023620276 2701
12800 1280 3.5836611e-06 0.00027352269 2701
13000 1300 3.6592211e-06 0.00029533734 2701
13200 1320 3.782506e-06 0.00030032559 2701
13400 1340 3.9807086e-06 0.00028395722 2701
13600 1360 4.2023176e-06 0.00025390325 2701
13800 1380 4.3559781e-06 0.00021794236 2701
14000 1400 4.4273371e-06 0.00018026034 2701
14200 1420 4.49867e-06 0.0001526569 2701
14400 1440 4.6591574e-06 0.00013707051 2701
14600 1460 4.9589583e-06 0.00013803875 2701
14800 1480 5.3859375e-06 0.00015455425 2701
15000 1500 5.8639557e-06 0.00017954785 2701
15200 1520 6.3075561e-06 0.0002084257 2701
15400 1540 6.7022179e-06 0.0002347669 2701
15600 1560 7.0789688e-06 0.00025020766 2701
15800 1580 7.4734777e-06 0.00025394845 2701
16000 1600 7.8884743e-06 0.00024571725 2701
16200 1620 8.3224059e-06 0.00022706648 2701
16400 1640 8.7337783e-06 0.00020320706 2701
16600 1660 9.1454649e-06 0.00017824346 2701
16800 1680 9.5948793e-06 0.00015961835 2701
17000 1700 1.0106407e-05 0.00015135471 2701
17200 1720 1.0707273e-05 0.00015166884 2701
17400 1740 1.1392597e-05 0.0001645916 2701
17600 1760 1.2118829e-05 0.00018119729 2701
17800 1780 1.2846056e-05 0.0002003616 2701
18000 1800 1.3555288e-05 0.00021585952 2701
18200 1820 1.4301024e-05 0.00022290158 2701
18400 1840 1.5089217e-05 0.00021970192 2701
18600 1860 1.5902351e-05 0.00020911128 2701
18800 1880 1.6753175e-05 0.00019278718 2701
19000 1900 1.7602996e-05 0.00017584076 2701
19200 1920 1.8479378e-05 0.00016206226 2701
19400 1940 1.9421603e-05 0.00015575677 2701
19600 1960 2.0477421e-05 0.00015687558 2701
19800 1980 2.1617288e-05 0.00016424998 2701
20000 2000 2.2814347e-05 0.00017466664 2701
20200 2020 2.4029097e-05 0.00018647149 2701
20400 2040 2.5255953e-05 0.00019516077 2701
20600 2060 2.649418e-05 0.00019906384 2701
20800 2080 2.7755897e-05 0.00019630586 2701
21000 2100 2.9067854e-05 0.00018674721 2701
21200 2120 3.0396477e-05 0.0001758048 2701
21400 2140 3.1759719e-05 0.00016782801 2701
21600 2160 3.3193597e-05 0.00016324138 2701
21800 2180 3.4729384e-05 0.00016124274 2701
22000 2200 3.6367594e-05 0.00016437457 2701
22200 2220 3.8095131e-05 0.00017015573 2701
22400 2240 3.9867003e-05 0.00017649465 2701
22600 2260 4.169511e-05 0.00018111374 2701
22800 2280 4.3566134e-05 0.00018104136 2701
23000 2300 4.5461538e-05 0.00017822707 2701
23200 2320 4.7377333e-05 0.00017285066 2701
23400 2340 4.9354403e-05 0.00016826524 2701
23600 2360 5.1399791e-05 0.00016517913 2701
23800 2380 5.3510931e-05 0.00016299649 2701
24000 2400 5.5681048e-05 0.00016256674 2701
24200 2420 5.7902429e-05 0.00016513449 2701
24400 2440 6.0216049e-05 0.00016895109 2701
24600 2460 6.270982e-05 0.00016946227 2701
24800 2480 6.5390117e-05 0.00016589426 2701
25000 2500 6.8121899e-05 0.00016241676 2701
25200 2520 7.0947331e-05 0.00015624292 2701
25400 2540 7.4304148e-05 0.0001449537 2701
25600 2560 7.7745077e-05 0.00013179658 2701
25800 2580 8.0739829e-05 0.00013098838 2701
26000 2600 8.3827874e-05 0.00014278841 2701
26200 2620 8.7060677e-05 0.00015381649 2701
26400 2640 9.0266508e-05 0.00016130999 2701
26600 2660 9.3339049e-05 0.00016908268 2701
26800 2680 9.6347013e-05 0.00016771087 2701
27000 2700 9.9294711e-05 0.00016577315 2701
27200 2720 0.00010230007 0.0001670893 2701
27400 2740 0.00010547172 0.00016569077 2701
27600 2760 0.00010872426 0.00016506303 2701
27800 2780 0.00011201844 0.00016482702 2701
28000 2800 0.00011532129 0.00016694886 2701
28200 2820 0.00011869854 0.00016163005 2701
28400 2840 0.00012209747 0.00015339281 2701
28600 2860 0.00012549322 0.00014765883 2701
28800 2880 0.00012898685 0.00014241765 2701
29000 2900 0.00013259039 0.00014215724 2701
29200 2920 0.00013628209 0.00014881155 2701
29400 2940 0.00014001213 0.00015671333 2701
29600 2960 0.00014379216 0.00016446215 2701
29800 2980 0.00014764687 0.0001639602 2701
30000 3000 0.00015142301 0.00015664816 2701
30200 3020 0.00015496407 0.00015545099 2701
30400 3040 0.00015797338 0.00015368625 2701
30600 3060 0.00016042141 0.00015679918 2701
30800 3080 0.00016244716 0.00016093678 2701
31000 3100 0.00016202247 0.00016066954 2701
31200 3120 0.0001613312 0.00015932059 2701
31400 3140 0.00016274961 0.00015988567 2701
31600 3160 0.00016541518 0.00015724809 2701
31800 3180 0.00016809362 0.00015498827 2701
32000 3200 0.00017067801 0.00014830489 2701
32200 3220 0.00017333906 0.00014371345 2701
32400 3240 0.0001759011 0.00014421259 2701
32600 3260 0.00017849952 0.00014228443 2701
32800 3280 0.00017801812 0.00014117391 2701
33000 3300 0.00017718857 0.00014644675 2701
33200 3320 0.00017833666 0.0001291286 2701
33400 3340 0.000178576 0.00014878558 2701
33600 3360 0.00017846711 0.00013905481 2701
33800 3380 0.00017822937 0.00015535996 2701
34000 3400 0.00017899663 0.00016094303 2701
34200 3420 0.00017924661 0.00015017553 2701
34400 3440 0.00018024855 0.00014723549 2701
34600 3460 0.00018143865 0.00013903131 2701
34800 3480 0.00018258173 0.00013722112 2701
35000 3500 0.00018404873 0.00014675949 2701
35200 3520 0.00018538521 0.00015108242 2701
35400 3540 0.00018669649 0.00014564852 2701
35600 3560 0.00018814608 0.00013762161 2701
35800 3580 0.00018967415 0.00014602307 2701
36000 3600 0.00019146735 0.000126909 2701
36200 3620 0.00019414036 0.00012384379 2701
36400 3640 0.00019613057 0.00011059573 2701
36600 3660 0.00019897104 0.00013621801 2701
36800 3680 0.00020169688 0.00013665462 2701
37000 3700 0.00020447655 0.00013929258 2701
37200 3720 0.00020711105 0.0001363895 2701
37400 3740 0.00021077854 0.00013610672 2701
37600 3760 0.00021303084 0.00015051235 2701
37800 3780 0.00021619561 0.00012664801 2701
38000 3800 0.0002194018 0.00012808247 2701
38200 3820 0.00022242646 0.0001360174 2701
38400 3840 0.00022531568 0.00013311221 2701
38600 3860 0.00022821731 0.00013523939 2701
38800 3880 0.000231228 0.00014090695 2701
39000 3900 0.00023404038 0.00013661835 2701
39200 3920 0.00023755044 0.00013659469 2701
39400 3940 0.00024009059 0.00012097907 2701
39600 3960 0.0002432098 9.7877876e-05 2701
39800 3980 0.00024475294 0.0001164688 2701
40000 4000 0.00024171274 0.00012432219 2701
Loop time of 192.659 on 4 procs for 40000 steps with 2701 atoms
Performance: 1793840.118 tau/day, 207.620 timesteps/s, 560.783 katom-step/s
99.6% CPU use with 4 MPI tasks x no OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 16.881 | 24.402 | 30.74 | 114.6 | 12.67
Bond | 1.1126 | 1.8917 | 2.6935 | 43.3 | 0.98
Neigh | 35.387 | 35.508 | 35.625 | 1.5 | 18.43
Comm | 1.5499 | 1.6694 | 1.8006 | 7.4 | 0.87
Output | 0.99755 | 1.0072 | 1.0165 | 0.8 | 0.52
Modify | 120.6 | 127.43 | 135.54 | 54.8 | 66.14
Other | | 0.7553 | | | 0.39
Nlocal: 675.25 ave 1373 max 7 min
Histogram: 2 0 0 0 0 0 0 0 0 2
Nghost: 103 ave 163 max 50 min
Histogram: 2 0 0 0 0 0 0 0 1 1
Neighs: 10509 ave 21592 max 126 min
Histogram: 2 0 0 0 0 0 0 0 0 2
FullNghs: 20367 ave 41981 max 141 min
Histogram: 2 0 0 0 0 0 0 0 0 2
Total # of neighbors = 81468
Ave neighs/atom = 30.162162
Ave special neighs/atom = 1.6593854
Neighbor list builds = 39932
Dangerous builds = 0
Total wall time: 0:03:12

2
examples/rheo/poiseuille/in.rheo.poiseuille
View File

@ -69,7 +69,7 @@ compute rho all rheo/property/atom rho
thermo 200
thermo_style custom step time ke press
dump 1 all custom 200 atomDump id type x y vx vy fx fy c_rho
#dump 1 all custom 200 atomDump id type x y vx vy fx fy c_rho
run 20000

288
examples/rheo/poiseuille/log.17Apr2024.poiseuille.g++.4 Normal file
View File

@ -0,0 +1,288 @@
LAMMPS (17 Apr 2024 - Development - patch_5May2020-18508-g3c0eaf6870-modified)
# ------ 2D Poiseuille flow ------ #
dimension 2
units lj
atom_style rheo
boundary p p p
comm_modify vel yes
# ------ Create simulation box ------ #
variable n equal 1.0
variable cut equal 3.0
region box block 0 20 -10 10 -0.01 0.01
create_box 2 box
Created orthogonal box = (0 -10 -0.01) to (20 10 0.01)
2 by 2 by 1 MPI processor grid
lattice sq ${n}
lattice sq 1
Lattice spacing in x,y,z = 1 1 1
region inner block INF INF -7.5 7.5 INF INF units box
region walls block INF INF -7.5 7.5 INF INF units box side out
create_atoms 2 region walls
Created 100 atoms
using lattice units in orthogonal box = (0 -10 -0.01) to (20 10 0.01)
create_atoms CPU = 0.000 seconds
create_atoms 1 region inner
Created 300 atoms
using lattice units in orthogonal box = (0 -10 -0.01) to (20 10 0.01)
create_atoms CPU = 0.000 seconds
group fluid type 1
300 atoms in group fluid
group rig type 2
100 atoms in group rig
displace_atoms fluid random 0.1 0.1 0 135414 units box
Displacing atoms ...
# ------ Model parameters ------ #
variable rho0 equal 1.0
variable cs equal 1.0
variable mp equal ${rho0}/${n}
variable mp equal 1/${n}
variable mp equal 1/1
variable zeta equal 1.0
variable kappa equal 1.0*${rho0}/${mp}
variable kappa equal 1.0*1/${mp}
variable kappa equal 1.0*1/1
variable fext equal 1e-4/${n}
variable fext equal 1e-4/1
variable dt_max equal 0.1*${cut}/${cs}/3
variable dt_max equal 0.1*3/${cs}/3
variable dt_max equal 0.1*3/1/3
variable Dr equal 0.05*${cut}*${cs}
variable Dr equal 0.05*3*${cs}
variable Dr equal 0.05*3*1
variable eta equal 0.1
variable gd0 equal 5e-4
variable npow equal 0.5
variable K equal 0.001
mass * ${mp}
mass * 1
set group all rheo/rho ${rho0}
set group all rheo/rho 1
Setting atom values ...
400 settings made for rheo/rho
set group all rheo/status 0
Setting atom values ...
400 settings made for rheo/status
set group rig rheo/status 1
Setting atom values ...
100 settings made for rheo/status
timestep ${dt_max}
timestep 0.1
pair_style rheo ${cut} artificial/visc ${zeta} rho/damp ${Dr}
pair_style rheo 3 artificial/visc ${zeta} rho/damp ${Dr}
pair_style rheo 3 artificial/visc 1 rho/damp ${Dr}
pair_style rheo 3 artificial/visc 1 rho/damp 0.15
pair_coeff * *
# ------ Fixes & computes ------ #
fix 1 all rheo ${cut} quintic 0 shift
fix 1 all rheo 3 quintic 0 shift
fix 2 all rheo/viscosity * constant ${eta}
fix 2 all rheo/viscosity * constant 0.1
#fix 2 all rheo/viscosity * power ${eta} ${gd0} ${K} ${npow}
fix 3 all rheo/pressure * linear
fix 4 rig setforce 0.0 0.0 0.0
fix 5 fluid addforce ${fext} 0.0 0.0
fix 5 fluid addforce 0.0001 0.0 0.0
fix 6 all enforce2d
compute rho all rheo/property/atom rho
# ------ Output & Run ------ #
thermo 200
thermo_style custom step time ke press
dump 1 all custom 200 atomDump id type x y vx vy fx fy c_rho
run 20000
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Your simulation uses code contributions which should be cited:
- @article{PalermoInPrep,
journal = {in prep},
title = {RHEO: A Hybrid Mesh-Free Model Framework for Dynamic Multi-Phase Flows},
year = {2024},
author = {Eric T. Palermo and Ki T. Wolf and Joel T. Clemmer and Thomas C. O'Connor},
}
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Generated 0 of 1 mixed pair_coeff terms from geometric mixing rule
Neighbor list info ...
update: every = 1 steps, delay = 0 steps, check = yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 3.3
ghost atom cutoff = 3.3
binsize = 1.65, bins = 13 13 1
4 neighbor lists, perpetual/occasional/extra = 4 0 0
(1) pair rheo, perpetual, half/full from (2)
attributes: half, newton on
pair build: halffull/newton
stencil: none
bin: none
(2) compute RHEO/KERNEL, perpetual
attributes: full, newton on
pair build: full/bin/atomonly
stencil: full/bin/2d
bin: standard
(3) compute RHEO/GRAD, perpetual, copy from (1)
attributes: half, newton on
pair build: copy
stencil: none
bin: none
(4) compute RHEO/VSHIFT, perpetual, copy from (1)
attributes: half, newton on
pair build: copy
stencil: none
bin: none
Per MPI rank memory allocation (min/avg/max) = 5.693 | 5.693 | 5.693 Mbytes
Step Time KinEng Press
0 0 0 0
200 20 1.2220462e-06 3.7383146e-05
400 40 4.345762e-06 7.5866885e-05
600 60 8.8559433e-06 0.00011353743
800 80 1.4370506e-05 0.00015135634
1000 100 2.0576198e-05 0.00018903722
1200 120 2.721926e-05 0.00022533997
1400 140 3.4099653e-05 0.00026016069
1600 160 4.1064175e-05 0.00029445207
1800 180 4.8001225e-05 0.00032893763
2000 200 5.4832849e-05 0.00036402396
2200 220 6.1508431e-05 0.00039945249
2400 240 6.8000141e-05 0.00043534411
2600 260 7.430136e-05 0.00046943441
2800 280 8.0415328e-05 0.00049807225
3000 300 8.6335032e-05 0.00051815375
3200 320 9.2021626e-05 0.00052618224
3400 340 9.7387936e-05 0.00051877918
3600 360 0.00010231526 0.00048650828
3800 380 0.00010676617 0.00044578079
4000 400 0.00011080098 0.00044777126
4200 420 0.00011448127 0.00047047629
4400 440 0.00011787852 0.00050280249
4600 460 0.00012106805 0.0005397213
4800 480 0.00012412056 0.00057885539
5000 500 0.0001271078 0.00061396896
5200 520 0.00013006637 0.00063981812
5400 540 0.00013295039 0.00065094073
5600 560 0.00013561487 0.00063918847
5800 580 0.00013791796 0.00059087656
6000 600 0.00013983422 0.00052171998
6200 620 0.00014144833 0.00050658002
6400 640 0.00014286538 0.0005248626
6600 660 0.00014417734 0.00055826606
6800 680 0.00014546931 0.00060063748
7000 700 0.00014682553 0.00064421411
7200 720 0.0001482833 0.00068252242
7400 740 0.00014977996 0.00070671308
7600 760 0.00015114829 0.00069774026
7800 780 0.0001522719 0.00064408311
8000 800 0.00015312897 0.00055977044
8200 820 0.00015375669 0.0005225573
8400 840 0.00015425683 0.00053833691
8600 860 0.00015471278 0.00057447427
8800 880 0.0001552059 0.00061980921
9000 900 0.00015581593 0.0006659836
9200 920 0.0001565564 0.00070813532
9400 940 0.00015733573 0.00073378551
9600 960 0.00015802107 0.00071560835
9800 980 0.00015855339 0.00065636189
10000 1000 0.00015890743 0.0005699855
10200 1020 0.00015908095 0.00053138971
10400 1040 0.00015915523 0.00054790708
10600 1060 0.00015921254 0.00058899454
10800 1080 0.00015934193 0.00063964906
11000 1100 0.00015959891 0.00069241358
11200 1120 0.0001599636 0.00073734651
11400 1140 0.00016036526 0.00074477329
11600 1160 0.00016075471 0.00071047555
11800 1180 0.00016109516 0.00064173183
12000 1200 0.00016131524 0.00055500553
12200 1220 0.00016136366 0.0005290215
12400 1240 0.0001613025 0.00055124296
12600 1260 0.00016123023 0.00059758627
12800 1280 0.00016123043 0.00065488735
13000 1300 0.00016132935 0.0007140876
13200 1320 0.00016152165 0.00074795629
13400 1340 0.00016180372 0.00074730778
13600 1360 0.00016216585 0.00071370995
13800 1380 0.0001625339 0.00065176323
14000 1400 0.00016274999 0.00057515371
14200 1420 0.00016271295 0.00055878258
14400 1440 0.00016249768 0.00058448193
14600 1460 0.00016223675 0.00063096229
14800 1480 0.00016201846 0.00068639548
15000 1500 0.00016190593 0.00072444357
15200 1520 0.00016194466 0.00073830636
15400 1540 0.00016216164 0.00072773256
15600 1560 0.00016253174 0.00069215481
15800 1580 0.00016290895 0.00063239408
16000 1600 0.00016306463 0.00057466273
16200 1620 0.00016292218 0.00057951567
16400 1640 0.00016261117 0.00061504156
16600 1660 0.00016225906 0.00066066637
16800 1680 0.00016197993 0.00069751908
17000 1700 0.0001618568 0.00072202303
17200 1720 0.00016194264 0.00073255034
17400 1740 0.00016225911 0.0007231031
17600 1760 0.00016270465 0.00068931224
17800 1780 0.00016304053 0.00062934836
18000 1800 0.00016302624 0.00058060272
18200 1820 0.00016274847 0.00058859513
18400 1840 0.00016236893 0.00061804803
18600 1860 0.00016202777 0.00065393237
18800 1880 0.0001618184 0.00068747094
19000 1900 0.0001618044 0.00071352541
19200 1920 0.00016204402 0.00072351769
19400 1940 0.00016249999 0.00071330322
19600 1960 0.00016297924 0.00067984167
19800 1980 0.00016317435 0.00061634142
20000 2000 0.00016301186 0.00057234115
Loop time of 15.6198 on 4 procs for 20000 steps with 400 atoms
Performance: 11062881.511 tau/day, 1280.426 timesteps/s, 512.170 katom-step/s
99.7% CPU use with 4 MPI tasks x no OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 2.1979 | 2.4473 | 2.6992 | 15.7 | 15.67
Neigh | 0.024709 | 0.027006 | 0.029223 | 1.3 | 0.17
Comm | 0.4657 | 0.71686 | 0.9662 | 29.0 | 4.59
Output | 0.033698 | 0.036781 | 0.039359 | 1.1 | 0.24
Modify | 12.306 | 12.313 | 12.319 | 0.2 | 78.83
Other | | 0.07916 | | | 0.51
Nlocal: 100 ave 107 max 93 min
Histogram: 1 0 0 1 0 0 1 0 0 1
Nghost: 185.5 ave 192 max 179 min
Histogram: 1 0 0 1 0 0 1 0 0 1
Neighs: 1712 ave 1848 max 1598 min
Histogram: 1 0 1 0 0 1 0 0 0 1
FullNghs: 3424 ave 3682 max 3174 min
Histogram: 1 0 1 0 0 0 1 0 0 1
Total # of neighbors = 13696
Ave neighs/atom = 34.24
Neighbor list builds = 331
Dangerous builds = 0
Total wall time: 0:00:15

224
examples/rheo/taylor-green/log.17Apr2024.taylor.green.g++.4 Normal file
View File

@ -0,0 +1,224 @@
LAMMPS (17 Apr 2024 - Development - patch_5May2020-18508-g3c0eaf6870-modified)
# ------ 2D Taylor Green vortex ------ #
dimension 2
units lj
atom_style rheo
boundary p p p
comm_modify vel yes
newton off
# ------ Create simulation box ------ #
variable n equal 1.0
variable cut equal 3.0
region box block 0 40 0 40 -0.01 0.01
create_box 1 box
Created orthogonal box = (0 0 -0.01) to (40 40 0.01)
2 by 2 by 1 MPI processor grid
lattice sq ${n}
lattice sq 1
Lattice spacing in x,y,z = 1 1 1
create_atoms 1 region box
Created 1600 atoms
using lattice units in orthogonal box = (0 0 -0.01) to (40 40 0.01)
create_atoms CPU = 0.001 seconds
displace_atoms all random 0.1 0.1 0 135414 units box
Displacing atoms ...
# ------ Model parameters ------ #
variable rho0 equal 1.0
variable mp equal ${rho0}/${n}
variable mp equal 1/${n}
variable mp equal 1/1
variable cs equal 1.0
variable eta equal 0.05
variable zeta equal 1
variable dt_max equal 0.1*${cut}/${cs}/3
variable dt_max equal 0.1*3/${cs}/3
variable dt_max equal 0.1*3/1/3
variable Dr equal 0.1*${cut}*${cs}
variable Dr equal 0.1*3*${cs}
variable Dr equal 0.1*3*1
mass * ${mp}
mass * 1
set group all rheo/rho ${rho0}
set group all rheo/rho 1
Setting atom values ...
1600 settings made for rheo/rho
set group all rheo/status 0
Setting atom values ...
1600 settings made for rheo/status
variable u0 equal 0.05
variable uy atom ${u0}*sin(2*PI*x/lx)*cos(2*PI*y/ly)
variable uy atom 0.05*sin(2*PI*x/lx)*cos(2*PI*y/ly)
variable ux atom -${u0}*sin(2*PI*y/ly)*cos(2*PI*x/ly)
variable ux atom -0.05*sin(2*PI*y/ly)*cos(2*PI*x/ly)
variable d0 atom ${rho0}-${u0}*${u0}*${rho0}*0.25*(cos(4*PI*x/lx)+cos(4*PI*y/ly))/${cs}/${cs}
variable d0 atom 1-${u0}*${u0}*${rho0}*0.25*(cos(4*PI*x/lx)+cos(4*PI*y/ly))/${cs}/${cs}
variable d0 atom 1-0.05*${u0}*${rho0}*0.25*(cos(4*PI*x/lx)+cos(4*PI*y/ly))/${cs}/${cs}
variable d0 atom 1-0.05*0.05*${rho0}*0.25*(cos(4*PI*x/lx)+cos(4*PI*y/ly))/${cs}/${cs}
variable d0 atom 1-0.05*0.05*1*0.25*(cos(4*PI*x/lx)+cos(4*PI*y/ly))/${cs}/${cs}
variable d0 atom 1-0.05*0.05*1*0.25*(cos(4*PI*x/lx)+cos(4*PI*y/ly))/1/${cs}
variable d0 atom 1-0.05*0.05*1*0.25*(cos(4*PI*x/lx)+cos(4*PI*y/ly))/1/1
velocity all set v_ux v_uy 0.0 units box
timestep ${dt_max}
timestep 0.1
pair_style rheo ${cut} artificial/visc ${zeta} rho/damp ${Dr}
pair_style rheo 3 artificial/visc ${zeta} rho/damp ${Dr}
pair_style rheo 3 artificial/visc 1 rho/damp ${Dr}
pair_style rheo 3 artificial/visc 1 rho/damp 0.3
pair_coeff * *
# ------ Fixes & computes ------ #
fix 1 all rheo ${cut} RK1 8 shift
fix 1 all rheo 3 RK1 8 shift
fix 2 all rheo/viscosity * constant ${eta}
fix 2 all rheo/viscosity * constant 0.05
fix 3 all rheo/pressure * linear
fix 4 all enforce2d
compute rho all rheo/property/atom rho
# ------ Output & Run ------ #
thermo 200
thermo_style custom step time ke press
dump 1 all custom 200 atomDump id type x y vx vy fx fy c_rho
run 10000
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Your simulation uses code contributions which should be cited:
- @article{PalermoInPrep,
journal = {in prep},
title = {RHEO: A Hybrid Mesh-Free Model Framework for Dynamic Multi-Phase Flows},
year = {2024},
author = {Eric T. Palermo and Ki T. Wolf and Joel T. Clemmer and Thomas C. O'Connor},
}
CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE-CITE
Generated 0 of 0 mixed pair_coeff terms from geometric mixing rule
Neighbor list info ...
update: every = 1 steps, delay = 0 steps, check = yes
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 3.3
ghost atom cutoff = 3.3
binsize = 1.65, bins = 25 25 1
4 neighbor lists, perpetual/occasional/extra = 4 0 0
(1) pair rheo, perpetual, half/full from (2)
attributes: half, newton off
pair build: halffull/newtoff
stencil: none
bin: none
(2) compute RHEO/KERNEL, perpetual
attributes: full, newton off
pair build: full/bin/atomonly
stencil: full/bin/2d
bin: standard
(3) compute RHEO/GRAD, perpetual, copy from (1)
attributes: half, newton off
pair build: copy
stencil: none
bin: none
(4) compute RHEO/VSHIFT, perpetual, copy from (1)
attributes: half, newton off
pair build: copy
stencil: none
bin: none
Per MPI rank memory allocation (min/avg/max) = 6.835 | 6.835 | 6.835 Mbytes
Step Time KinEng Press
0 0 0.00062497276 0.00062607301
200 20 0.00056200647 0.00056633785
400 40 0.00050570968 0.00051098771
600 60 0.00045586684 0.00046081672
800 80 0.00041124523 0.00041549607
1000 100 0.00037065341 0.00037412741
1200 120 0.00033391585 0.00033580899
1400 140 0.00030078316 0.00030057307
1600 160 0.00027093231 0.00026842603
1800 180 0.00024403239 0.00023839026
2000 200 0.0002197865 0.00021148941
2200 220 0.0001979269 0.00018659386
2400 240 0.00017822267 0.00016430442
2600 260 0.00016047141 0.00014408514
2800 280 0.00014448504 0.00012574125
3000 300 0.00013009159 0.00010869938
3200 320 0.00011713578 9.414951e-05
3400 340 0.00010547564 8.1900579e-05
3600 360 9.4982139e-05 7.1285649e-05
3800 380 8.5538983e-05 6.1571123e-05
4000 400 7.7040171e-05 5.3462572e-05
4200 420 6.9390317e-05 4.6338308e-05
4400 440 6.2503763e-05 3.9697323e-05
4600 460 5.6303766e-05 3.4234465e-05
4800 480 5.0721595e-05 3.0841338e-05
5000 500 4.5695301e-05 2.7788566e-05
5200 520 4.1169161e-05 2.5744409e-05
5400 540 3.7093059e-05 2.3912739e-05
5600 560 3.3421819e-05 2.2494185e-05
5800 580 3.0114735e-05 2.1594384e-05
6000 600 2.7135224e-05 2.1164421e-05
6200 620 2.4450446e-05 2.0979349e-05
6400 640 2.2030925e-05 2.0858567e-05
6600 660 1.9850196e-05 2.098115e-05
6800 680 1.7884553e-05 2.1134827e-05
7000 700 1.6112763e-05 2.1242242e-05
7200 720 1.4515783e-05 2.1312763e-05
7400 740 1.3076456e-05 2.1370947e-05
7600 760 1.1779327e-05 2.1332126e-05
7800 780 1.0610469e-05 2.1156562e-05
8000 800 9.5573298e-06 2.0898126e-05
8200 820 8.6085799e-06 2.0517958e-05
8400 840 7.7539888e-06 1.9841551e-05
8600 860 6.9843033e-06 1.9114769e-05
8800 880 6.2911575e-06 1.8362959e-05
9000 900 5.6669785e-06 1.7473404e-05
9200 920 5.1049208e-06 1.6452745e-05
9400 940 4.5987908e-06 1.5578629e-05
9600 960 4.1429972e-06 1.4427274e-05
9800 980 3.7324962e-06 1.3169649e-05
10000 1000 3.3627455e-06 1.1938723e-05
Loop time of 38.2006 on 4 procs for 10000 steps with 1600 atoms
Performance: 2261743.875 tau/day, 261.776 timesteps/s, 418.841 katom-step/s
99.7% CPU use with 4 MPI tasks x no OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 8.2958 | 8.7273 | 9.3582 | 15.2 | 22.85
Neigh | 0.034282 | 0.035689 | 0.037115 | 0.7 | 0.09
Comm | 0.16788 | 0.17018 | 0.17278 | 0.4 | 0.45
Output | 0.066977 | 0.06882 | 0.071704 | 0.7 | 0.18
Modify | 28.483 | 28.793 | 28.962 | 3.6 | 75.37
Other | | 0.4053 | | | 1.06
Nlocal: 400 ave 402 max 399 min
Histogram: 2 0 0 1 0 0 0 0 0 1
Nghost: 307.25 ave 308 max 305 min
Histogram: 1 0 0 0 0 0 0 0 0 3
Neighs: 7618.25 ave 7697 max 7564 min
Histogram: 1 0 1 1 0 0 0 0 0 1
FullNghs: 13343 ave 13497 max 13258 min
Histogram: 1 1 1 0 0 0 0 0 0 1
Total # of neighbors = 53372
Ave neighs/atom = 33.3575
Neighbor list builds = 123
Dangerous builds = 0
Total wall time: 0:00:38

2
src/BPM/bond_bpm_spring.cpp
View File

@ -23,8 +23,6 @@
#include "modify.h"
#include "neighbor.h"
#include "update.h"
#include <cmath>
#include <cstring>

2
src/BPM/compute_nbond_atom.cpp
View File

@ -30,7 +30,7 @@ ComputeNBondAtom::ComputeNBondAtom(LAMMPS *_lmp, int narg, char **arg) :
if (narg < 3) utils::missing_cmd_args(FLERR, "compute nbond/atom", error);
if (atom->avec->bonds_allow == 0)
error->all(FLERR,"Compute nbond/atom used when bonds are not allowed");
error->all(FLERR,"Compute nbond/atom used in system without bonds");
btype = -1;
int iarg = 3;

2
src/BPM/fix_update_special_bonds.cpp
View File

@ -169,7 +169,7 @@ void FixUpdateSpecialBonds::pre_force(int /*vflag*/)
tagint *tag = atom->tag;
// In theory could communicate a list of broken bonds to neighboring processors here
// to remove restriction that users use Newton bond off
// to remove restriction on Newton bond off
for (int ilist = 0; ilist < neighbor->nlist; ilist++) {
list = neighbor->lists[ilist];

5
src/RHEO/bond_rheo_shell.cpp
View File

@ -388,9 +388,8 @@ void BondRHEOShell::settings(int narg, char **arg)
for (std::size_t i = 0; i < leftover_iarg.size(); i++) {
iarg = leftover_iarg[i];
if (strcmp(arg[iarg], "t/form") == 0) {
if (iarg + 1 > narg) error->all(FLERR, "Illegal bond rheo/shell command, missing option for t/form");
if (iarg + 1 > narg) utils::missing_cmd_args(FLERR, "bond rheo/shell t/form", error);
tform = utils::numeric(FLERR, arg[iarg + 1], false, lmp);
if (tform < 0.0) error->all(FLERR, "Illegal bond rheo/shell value for t/form, {}", tform);
i += 1;
} else {
error->all(FLERR, "Illegal bond rheo/shell command, invalid argument {}", arg[iarg]);
@ -398,7 +397,7 @@ void BondRHEOShell::settings(int narg, char **arg)
}
if (tform < 0.0)
error->all(FLERR, "Illegal bond rheo/shell command, must specify formation time");
error->all(FLERR, "Illegal bond rheo/shell command, must specify positive formation time");
}

150
src/RHEO/compute_rheo_kernel.cpp
View File

@ -110,27 +110,27 @@ void ComputeRHEOKernel::init()
compute_interface = fix_rheo->compute_interface;
zmin = fix_rheo->zmin_kernel;
h = fix_rheo->h;
hsq = h * h;
hinv = 1.0 / h;
hsqinv = hinv * hinv;
cut = fix_rheo->cut;
cutsq = cut * cut;
cutinv = 1.0 / cut;
cutsqinv = cutinv * cutinv;
if (kernel_style != WENDLANDC4) {
if (dim == 3) {
pre_w = 1.0 / (120.0 * MY_PI) * 27.0 * hsqinv * hinv;
pre_wp = pre_w * 3.0 * hinv;
pre_w = 1.0 / (120.0 * MY_PI) * 27.0 * cutsqinv * cutinv;
pre_wp = pre_w * 3.0 * cutinv;
} else {
pre_w = 7.0 / (478.0 * MY_PI) * 9 * hsqinv;
pre_wp = pre_w * 3.0 * hinv;
pre_w = 7.0 / (478.0 * MY_PI) * 9 * cutsqinv;
pre_wp = pre_w * 3.0 * cutinv;
}
} else {
if (dim == 3) {
pre_w = 495.0 / (32.0 * MY_PI * hsq * h);
pre_wp = pre_w * hinv;
pre_w = 495.0 / (32.0 * MY_PI * cutsq * cut);
pre_wp = pre_w * cutinv;
} else {
pre_w = 9.0 / (MY_PI * hsq);
pre_wp = pre_w * hinv;
pre_w = 9.0 / (MY_PI * cutsq);
pre_wp = pre_w * cutinv;
}
}
@ -246,7 +246,7 @@ double ComputeRHEOKernel::calc_dw(int i, int j, double delx, double dely, double
double ComputeRHEOKernel::calc_w_quintic(int i, int j, double delx, double dely, double delz, double r)
{
double w, tmp1, tmp2, tmp3, tmp1sq, tmp2sq, tmp3sq, s;
s = r * 3.0 * hinv;
s = r * 3.0 * cutinv;
if (s > 3.0) {
w = 0.0;
@ -284,7 +284,7 @@ double ComputeRHEOKernel::calc_dw_quintic(int i, int j, double delx, double dely
double *mass = atom->mass;
int *type = atom->type;
s = r * 3.0 * hinv;
s = r * 3.0 * cutinv;
if (s > 3.0) {
wp = 0.0;
@ -323,7 +323,7 @@ double ComputeRHEOKernel::calc_dw_quintic(int i, int j, double delx, double dely
double ComputeRHEOKernel::calc_w_wendlandc4(int i, int j, double delx, double dely, double delz, double r)
{
double w, tmp6, s;
s = r * hinv;
s = r * cutinv;
if (s > 1.0) {
w = 0.0;
@ -349,7 +349,7 @@ double ComputeRHEOKernel::calc_dw_wendlandc4(int i, int j, double delx, double d
double *mass = atom->mass;
int *type = atom->type;
s = r * hinv;
s = r * cutinv;
if (s > 1.0) {
wp = 0.0;
@ -403,13 +403,13 @@ double ComputeRHEOKernel::calc_w_rk1(int i, int j, double delx, double dely, dou
if (dim == 2) {
H[0] = 1.0;
H[1] = dx[0] * hinv;
H[2] = dx[1] * hinv;
H[1] = dx[0] * cutinv;
H[2] = dx[1] * cutinv;
} else {
H[0] = 1.0;
H[1] = dx[0] * hinv;
H[2] = dx[1] * hinv;
H[3] = dx[2] * hinv;
H[1] = dx[0] * cutinv;
H[2] = dx[1] * cutinv;
H[3] = dx[2] * cutinv;
}
Wij = 0;
for (b = 0; b < Mdim; b++) {
@ -444,22 +444,22 @@ double ComputeRHEOKernel::calc_w_rk2(int i, int j, double delx, double dely, dou
if (dim == 2) {
H[0] = 1.0;
H[1] = dx[0] * hinv;
H[2] = dx[1] * hinv;
H[3] = 0.5 * dx[0] * dx[0] * hsqinv;
H[4] = 0.5 * dx[1] * dx[1] * hsqinv;
H[5] = dx[0] * dx[1] * hsqinv;
H[1] = dx[0] * cutinv;
H[2] = dx[1] * cutinv;
H[3] = 0.5 * dx[0] * dx[0] * cutsqinv;
H[4] = 0.5 * dx[1] * dx[1] * cutsqinv;
H[5] = dx[0] * dx[1] * cutsqinv;
} else {
H[0] = 1.0;
H[1] = dx[0] * hinv;
H[2] = dx[1] * hinv;
H[3] = dx[2] * hinv;
H[4] = 0.5 * dx[0] * dx[0] * hsqinv;
H[5] = 0.5 * dx[1] * dx[1] * hsqinv;
H[6] = 0.5 * dx[2] * dx[2] * hsqinv;
H[7] = dx[0] * dx[1] * hsqinv;
H[8] = dx[0] * dx[2] * hsqinv;
H[9] = dx[1] * dx[2] * hsqinv;
H[1] = dx[0] * cutinv;
H[2] = dx[1] * cutinv;
H[3] = dx[2] * cutinv;
H[4] = 0.5 * dx[0] * dx[0] * cutsqinv;
H[5] = 0.5 * dx[1] * dx[1] * cutsqinv;
H[6] = 0.5 * dx[2] * dx[2] * cutsqinv;
H[7] = dx[0] * dx[1] * cutsqinv;
H[8] = dx[0] * dx[2] * cutsqinv;
H[9] = dx[1] * dx[2] * cutsqinv;
}
Wij = 0;
for (b = 0; b < Mdim; b++) {
@ -496,13 +496,13 @@ void ComputeRHEOKernel::calc_dw_rk1(int i, int j, double delx, double dely, doub
//Populate correction basis
if (dim == 2) {
H[0] = 1.0;
H[1] = dx[0] * hinv;
H[2] = dx[1] * hinv;
H[1] = dx[0] * cutinv;
H[2] = dx[1] * cutinv;
} else {
H[0] = 1.0;
H[1] = dx[0] * hinv;
H[2] = dx[1] * hinv;
H[3] = dx[2] * hinv;
H[1] = dx[0] * cutinv;
H[2] = dx[1] * cutinv;
H[3] = dx[2] * cutinv;
}
// dWij[] = dWx dWy (dWz)
@ -533,22 +533,22 @@ void ComputeRHEOKernel::calc_dw_rk2(int i, int j, double delx, double dely, doub
//Populate correction basis
if (dim == 2) {
H[0] = 1.0;
H[1] = dx[0] * hinv;
H[2] = dx[1] * hinv;
H[3] = 0.5 * dx[0] * dx[0] * hsqinv;
H[4] = 0.5 * dx[1] * dx[1] * hsqinv;
H[5] = dx[0] * dx[1] * hsqinv;
H[1] = dx[0] * cutinv;
H[2] = dx[1] * cutinv;
H[3] = 0.5 * dx[0] * dx[0] * cutsqinv;
H[4] = 0.5 * dx[1] * dx[1] * cutsqinv;
H[5] = dx[0] * dx[1] * cutsqinv;
} else {
H[0] = 1.0;
H[1] = dx[0] * hinv;
H[2] = dx[1] * hinv;
H[3] = dx[2] * hinv;
H[4] = 0.5 * dx[0] * dx[0] * hsqinv;
H[5] = 0.5 * dx[1] * dx[1] * hsqinv;
H[6] = 0.5 * dx[2] * dx[2] * hsqinv;
H[7] = dx[0] * dx[1] * hsqinv;
H[8] = dx[0] * dx[2] * hsqinv;
H[9] = dx[1] * dx[2] * hsqinv;
H[1] = dx[0] * cutinv;
H[2] = dx[1] * cutinv;
H[3] = dx[2] * cutinv;
H[4] = 0.5 * dx[0] * dx[0] * cutsqinv;
H[5] = 0.5 * dx[1] * dx[1] * cutsqinv;
H[6] = 0.5 * dx[2] * dx[2] * cutsqinv;
H[7] = dx[0] * dx[1] * cutsqinv;
H[8] = dx[0] * dx[2] * cutsqinv;
H[9] = dx[1] * dx[2] * cutsqinv;
}
// dWij[] = dWx dWy (dWz)
@ -621,7 +621,7 @@ void ComputeRHEOKernel::compute_peratom()
dx[2] = ztmp - x[j][2];
rsq = lensq3(dx);
if (rsq < hsq) {
if (rsq < cutsq) {
r = sqrt(rsq);
w = calc_w_quintic(i, j, dx[0], dx[1], dx[2], r);
rhoj = rho[j];
@ -639,7 +639,7 @@ void ComputeRHEOKernel::compute_peratom()
}
} else if (correction_order > 0) {
// Moment matrix M and polynomial basis vector H (1d for gsl compatibility)
// Moment matrix M and polynomial basis vector cut (1d for gsl compatibility)
double H[Mdim], M[Mdim * Mdim];
for (ii = 0; ii < inum; ii++) {
@ -652,7 +652,7 @@ void ComputeRHEOKernel::compute_peratom()
jnum = numneigh[i];
itype = type[i];
// Zero upper-triangle M and H (will be symmetric):
// Zero upper-triangle M and cut (will be symmetric):
for (a = 0; a < Mdim; a++) {
for (b = a; b < Mdim; b++) {
M[a * Mdim + b] = 0;
@ -669,7 +669,7 @@ void ComputeRHEOKernel::compute_peratom()
rsq = lensq3(dx);
if (rsq < hsq) {
if (rsq < cutsq) {
r = sqrt(rsq);
w = calc_w_quintic(i, j, dx[0], dx[1], dx[2], r);
@ -683,25 +683,25 @@ void ComputeRHEOKernel::compute_peratom()
//Populate the H-vector of polynomials (2D)
if (dim == 2) {
H[0] = 1.0;
H[1] = dx[0] * hinv;
H[2] = dx[1] * hinv;
H[1] = dx[0] * cutinv;
H[2] = dx[1] * cutinv;
if (kernel_style == RK2) {
H[3] = 0.5 * dx[0] * dx[0] * hsqinv;
H[4] = 0.5 * dx[1] * dx[1] * hsqinv;
H[5] = dx[0] * dx[1] * hsqinv;
H[3] = 0.5 * dx[0] * dx[0] * cutsqinv;
H[4] = 0.5 * dx[1] * dx[1] * cutsqinv;
H[5] = dx[0] * dx[1] * cutsqinv;
}
} else {
H[0] = 1.0;
H[1] = dx[0] * hinv;
H[2] = dx[1] * hinv;
H[3] = dx[2] * hinv;
H[1] = dx[0] * cutinv;
H[2] = dx[1] * cutinv;
H[3] = dx[2] * cutinv;
if (kernel_style == RK2) {
H[4] = 0.5 * dx[0] * dx[0] * hsqinv;
H[5] = 0.5 * dx[1] * dx[1] * hsqinv;
H[6] = 0.5 * dx[2] * dx[2] * hsqinv;
H[7] = dx[0] * dx[1] * hsqinv;
H[8] = dx[0] * dx[2] * hsqinv;
H[9] = dx[1] * dx[2] * hsqinv;
H[4] = 0.5 * dx[0] * dx[0] * cutsqinv;
H[5] = 0.5 * dx[1] * dx[1] * cutsqinv;
H[6] = 0.5 * dx[2] * dx[2] * cutsqinv;
H[7] = dx[0] * dx[1] * cutsqinv;
H[8] = dx[0] * dx[2] * cutsqinv;
H[9] = dx[1] * dx[2] * cutsqinv;
}
}
@ -765,12 +765,12 @@ void ComputeRHEOKernel::compute_peratom()
C[i][0][a] = M[a * Mdim + 0]; // all rows of column 0
for (b = 0; b < dim; b++) {
//First derivatives
C[i][1 + b][a] = -M[a * Mdim + b + 1] * hinv;
C[i][1 + b][a] = -M[a * Mdim + b + 1] * cutinv;
// columns 1-2 (2D) or 1-3 (3D)
//Second derivatives
if (kernel_style == RK2)
C[i][1 + dim + b][a] = M[a * Mdim + b + 1 + dim] * hsqinv;
C[i][1 + dim + b][a] = M[a * Mdim + b + 1 + dim] * cutsqinv;
// columns 3-4 (2D) or 4-6 (3D)
}
}
@ -822,7 +822,7 @@ void ComputeRHEOKernel::compute_coordination()
dx[2] = ztmp - x[j][2];
rsq = lensq3(dx);
if (rsq < hsq)
if (rsq < cutsq)
coordination[i] += 1;
}
}

2
src/RHEO/compute_rheo_kernel.h
View File

@ -59,7 +59,7 @@ class ComputeRHEOKernel : public Compute {
int corrections_calculated;
int kernel_style, zmin, dim, Mdim, ncor;
int nmax_store;
double h, hsq, hinv, hsqinv, pre_w, pre_wp;
double cut, cutsq, cutinv, cutsqinv, pre_w, pre_wp;
double ***C;
double *C0;

4
src/RHEO/compute_rheo_surface.cpp
View File

@ -181,7 +181,7 @@ void ComputeRHEOSurface::compute_peratom()
gradC[i][a] += dWij[a] * Volj;
}
if (j < nlocal || newton) {
if ((j < nlocal) || newton) {
for (a = 0; a < dim; a++){
divr[j] += dWji[a] * dx[a] * Voli;
gradC[j][a] += dWji[a] * Voli;
@ -287,7 +287,7 @@ void ComputeRHEOSurface::compute_peratom()
}
}
// clear normal vectors for non surface particles
// clear normal vectors for non-surface particles
for (i = 0; i < nall; i++) {
if (mask[i] & groupbit) {

47
src/RHEO/fix_rheo.cpp
View File

@ -39,7 +39,12 @@ using namespace RHEO_NS;
using namespace FixConst;
static const char cite_rheo[] =
"TBD\n\n";
"@article{PalermoInPrep,\n"
" journal = {in prep},\n"
" title = {RHEO: A Hybrid Mesh-Free Model Framework for Dynamic Multi-Phase Flows},\n"
" year = {2024},\n"
" author = {Eric T. Palermo and Ki T. Wolf and Joel T. Clemmer and Thomas C. O'Connor},\n"
"}\n\n";
/* ---------------------------------------------------------------------- */
@ -84,10 +89,9 @@ FixRHEO::FixRHEO(LAMMPS *lmp, int narg, char **arg) :
error->all(FLERR, "fix rheo command requires atom_style with status");
if (narg < 5)
error->all(FLERR, "Insufficient arguments for fix rheo command");
utils::missing_cmd_args(FLERR, "fix rheo", error);
h = utils::numeric(FLERR, arg[3], false, lmp);
cut = h;
cut = utils::numeric(FLERR, arg[3], false, lmp);
if (strcmp(arg[4], "quintic") == 0) {
kernel_style = QUINTIC;
} else if (strcmp(arg[4], "wendland/c4") == 0) {
@ -109,7 +113,7 @@ FixRHEO::FixRHEO(LAMMPS *lmp, int narg, char **arg) :
thermal_flag = 1;
} else if (strcmp(arg[iarg], "surface/detection") == 0) {
surface_flag = 1;
if(iarg + 3 >= narg) error->all(FLERR, "Illegal surface/detection option in fix rheo");
if(iarg + 3 >= narg) utils::missing_cmd_args(FLERR, "fix rheo surface/detection", error);
if (strcmp(arg[iarg + 1], "coordination") == 0) {
surface_style = COORDINATION;
zmin_surface = utils::inumeric(FLERR, arg[iarg + 2], false, lmp);
@ -130,12 +134,12 @@ FixRHEO::FixRHEO(LAMMPS *lmp, int narg, char **arg) :
} else if (strcmp(arg[iarg], "self/mass") == 0) {
self_mass_flag = 1;
} else if (strcmp(arg[iarg], "density") == 0) {
if (iarg + n >= narg) error->all(FLERR, "Illegal rho0 option in fix rheo");
if (iarg + n >= narg) utils::missing_cmd_args(FLERR, "fix rheo density", error);
for (i = 1; i <= n; i++)
rho0[i] = utils::numeric(FLERR, arg[iarg + i], false, lmp);
iarg += n;
} else if (strcmp(arg[iarg], "speed/sound") == 0) {
if (iarg + n >= narg) error->all(FLERR, "Illegal csq option in fix rheo");
if (iarg + n >= narg) utils::missing_cmd_args(FLERR, "fix rheo speed/sound", error);
for (i = 1; i <= n; i++) {
csq[i] = utils::numeric(FLERR, arg[iarg + i], false, lmp);
csq[i] *= csq[i];
@ -281,7 +285,7 @@ void FixRHEO::setup(int /*vflag*/)
if (!pressure_fix_defined)
error->all(FLERR, "Missing fix rheo/pressure");
if((!thermal_fix_defined) && thermal_flag)
if(thermal_flag && !thermal_fix_defined)
error->all(FLERR, "Missing fix rheo/thermal");
// Reset to zero for future runs
@ -290,33 +294,6 @@ void FixRHEO::setup(int /*vflag*/)
pressure_fix_defined = 0;
oxidation_fix_defined = 0;
// Check fixes cover all atoms (may still fail if atoms are created)
// FixRHEOPressure currently requires group all
auto visc_fixes = modify->get_fix_by_style("rheo/viscosity");
auto therm_fixes = modify->get_fix_by_style("rheo/thermal");
int *mask = atom->mask;
int v_coverage_flag = 1;
int t_coverage_flag = 1;
int covered;
for (int i = 0; i < atom->nlocal; i++) {
covered = 0;
for (auto fix : visc_fixes)
if (mask[i] & fix->groupbit) covered = 1;
if (!covered) v_coverage_flag = 0;
if (thermal_flag) {
covered = 0;
for (auto fix : therm_fixes)
if (mask[i] & fix->groupbit) covered = 1;
if (!covered) t_coverage_flag = 0;
}
}
if (!v_coverage_flag)
error->one(FLERR, "Fix rheo/viscosity does not fully cover all atoms");
if (!t_coverage_flag)
error->one(FLERR, "Fix rheo/thermal does not fully cover all atoms");
if (rhosum_flag)
compute_rhosum->compute_peratom();
}

2
src/RHEO/fix_rheo.h
View File

@ -39,7 +39,7 @@ class FixRHEO : public Fix {
void reset_dt() override;
// Model parameters
double h, cut;
double cut;
double *rho0, *csq;
int self_mass_flag;
int zmin_kernel, zmin_surface, zmin_splash;

6
src/RHEO/fix_rheo_oxidation.cpp
View File

@ -99,10 +99,10 @@ void FixRHEOOxidation::init()
if (fixes.size() == 0) error->all(FLERR, "Need to define fix rheo to use fix rheo/oxidation");
fix_rheo = dynamic_cast<FixRHEO *>(fixes[0]);
if (cut > fix_rheo->h)
if (cut > fix_rheo->cut)
error->all(FLERR, "Bonding length exceeds kernel cutoff");
if (rsurf >= fix_rheo->h)
if (rsurf >= fix_rheo->cut)
error->all(FLERR, "Surface distance must be less than kernel cutoff");
if (!force->bond) error->all(FLERR, "Must define a bond style with fix rheo/oxidation");
@ -233,7 +233,7 @@ void FixRHEOOxidation::post_integrate()
// Add bonds to owned atoms
// If newton bond off, add to both, otherwise add to whichever has a smaller tag
if (!newton_bond || tagi < tagj) {
if (!newton_bond || (tagi < tagj)) {
if (num_bond[i] == atom->bond_per_atom)
error->one(FLERR,"New bond exceeded bonds per atom in fix rheo/oxidation for atom {}", tagi);
bond_type[i][num_bond[i]] = btype;

13
src/RHEO/fix_rheo_thermal.cpp
View File

@ -71,6 +71,7 @@ FixRHEOThermal::FixRHEOThermal(LAMMPS *lmp, int narg, char **arg) :
cut_bond = 0;
comm_forward = 0;
// Currently can only have one instance of fix rheo/thermal
if (igroup != 0)
error->all(FLERR,"fix rheo/thermal command requires group all");
@ -249,7 +250,7 @@ void FixRHEOThermal::init()
auto fixes = modify->get_fix_by_style("^rheo$");
if (fixes.size() == 0) error->all(FLERR, "Need to define fix rheo to use fix rheo/viscosity");
fix_rheo = dynamic_cast<FixRHEO *>(fixes[0]);
cut_kernel = fix_rheo->h;
cut_kernel = fix_rheo->cut;
if (cut_bond > cut_kernel)
error->all(FLERR, "Bonding length exceeds kernel cutoff");
@ -561,7 +562,7 @@ void FixRHEOThermal::break_bonds()
// Update special unless two owned atoms melt simultaneously then
// only update for atom with lower tag
if (fix_update_special_bonds) {
if (i < nlocal && j < nlocal && melti && meltj) {
if ((i < nlocal) && (j < nlocal) && melti && meltj) {
if (tag[i] < tag[j]) {
fix_update_special_bonds->add_broken_bond(i, j);
}
@ -590,7 +591,7 @@ void FixRHEOThermal::break_bonds()
// Delete bonds for non-melted local atoms (shifting)
if (i < nlocal && !melti) {
for (m = 0; m < num_bond[i]; m++) {
if (bond_atom[i][m] == tag[j] && bond_type[i][m] == btype) {
if ((bond_atom[i][m] == tag[j]) && (bond_type[i][m] == btype)) {
nmax = num_bond[i] - 1;
bond_type[i][m] = bond_type[i][nmax];
bond_atom[i][m] = bond_atom[i][nmax];
@ -609,7 +610,7 @@ void FixRHEOThermal::break_bonds()
if (j < nlocal && !meltj) {
for (m = 0; m < num_bond[j]; m++) {
if (bond_atom[j][m] == tag[i] && bond_type[j][m] == btype) {
if ((bond_atom[j][m] == tag[i]) && (bond_type[j][m] == btype)) {
nmax = num_bond[j] - 1;
bond_type[j][m] = bond_type[j][nmax];
bond_atom[j][m] = bond_atom[j][nmax];
@ -692,7 +693,7 @@ void FixRHEOThermal::create_bonds()
// Add bonds to owned atoms
// If newton bond off, add to both, otherwise add to whichever has a smaller tag
if (i < nlocal && (!newton_bond || tag[i] < tag[j])) {
if ((i < nlocal) && (!newton_bond || (tag[i] < tag[j]))) {
if (num_bond[i] == atom->bond_per_atom)
error->one(FLERR,"New bond exceeded bonds per atom in fix rheo/thermal");
bond_type[i][num_bond[i]] = btype;
@ -700,7 +701,7 @@ void FixRHEOThermal::create_bonds()
num_bond[i]++;
}
if (j < nlocal && (!newton_bond || tag[j] < tag[i])) {
if ((j < nlocal) && (!newton_bond || (tag[j] < tag[i]))) {
if (num_bond[j] == atom->bond_per_atom)
error->one(FLERR,"New bond exceeded bonds per atom in fix rheo/thermal");
bond_type[j][num_bond[j]] = btype;

1
src/RHEO/fix_rheo_viscosity.cpp
View File

@ -46,6 +46,7 @@ FixRHEOViscosity::FixRHEOViscosity(LAMMPS *lmp, int narg, char **arg) :
constant_flag = 0;
evolve_flag = 0;
// Currently can only have one instance of fix rheo/viscosity
if (igroup != 0)
error->all(FLERR,"fix rheo/viscosity command requires group all");

34
src/RHEO/pair_rheo.cpp
View File

@ -33,7 +33,6 @@
#include "modify.h"
#include "neighbor.h"
#include "neigh_list.h"
#include "update.h"
#include "utils.h"
#include <cmath>
@ -112,7 +111,6 @@ void PairRHEO::compute(int eflag, int vflag)
int *type = atom->type;
int *status = atom->status;
tagint *tag = atom->tag;
double fnorm, ftang[3];
double **fp_store, *chi;
if (compute_interface) {
@ -169,7 +167,7 @@ void PairRHEO::compute(int eflag, int vflag)
rsq = lensq3(dx);
jtype = type[j];
if (rsq < hsq) {
if (rsq < cutksq) {
r = sqrt(rsq);
rinv = 1 / r;
@ -219,16 +217,16 @@ void PairRHEO::compute(int eflag, int vflag)
rhoj = compute_interface->correct_rho(j, i);
Pj = fix_pressure->calc_pressure(rhoj, jtype);
if ((chi[j] > 0.9) && (r < (h * 0.5)))
fmag = (chi[j] - 0.9) * (h * 0.5 - r) * rho0j * csq_ave * h * rinv;
if ((chi[j] > 0.9) && (r < (cutk * 0.5)))
fmag = (chi[j] - 0.9) * (cutk * 0.5 - r) * rho0j * csq_ave * cutk * rinv;
} else if ((!fluidi) && fluidj) {
compute_interface->correct_v(vi, vj, i, j);
rhoi = compute_interface->correct_rho(i, j);
Pi = fix_pressure->calc_pressure(rhoi, itype);
if (chi[i] > 0.9 && r < (h * 0.5))
fmag = (chi[i] - 0.9) * (h * 0.5 - r) * rho0i * csq_ave * h * rinv;
if (chi[i] > 0.9 && r < (cutk * 0.5))
fmag = (chi[i] - 0.9) * (cutk * 0.5 - r) * rho0i * csq_ave * cutk * rinv;
} else if ((!fluidi) && (!fluidj)) {
rhoi = rho0i;
@ -281,8 +279,8 @@ void PairRHEO::compute(int eflag, int vflag)
for (b = 0; b < dim; b++)
du[a] -= 0.5 * (gradv[i][a * dim + b] + gradv[j][a * dim + b]) * dx[b];
mu = dot3(du, dx) * hinv3;
mu /= (rsq * hinv3 * hinv3 + EPSILON);
mu = dot3(du, dx) * cutkinv3;
mu /= (rsq * cutkinv3 * cutkinv3 + EPSILON);
mu = MIN(0.0, mu);
q = av * (-2.0 * cs_ave * mu + mu * mu);
fp_prefactor += voli * volj * q * (rhoj + rhoi);
@ -406,17 +404,17 @@ void PairRHEO::settings(int narg, char **arg)
{
if (narg < 1) error->all(FLERR, "Illegal pair_style command");
h = utils::numeric(FLERR, arg[0], false, lmp);
cutk = utils::numeric(FLERR, arg[0], false, lmp);
int iarg = 1;
while (iarg < narg) {
if (strcmp(arg[iarg], "rho/damp") == 0) {
if (iarg + 1 >= narg) error->all(FLERR, "Illegal pair_style command");
if (iarg + 1 >= narg) utils::missing_cmd_args(FLERR, "pair rheo rho/damp", error);
rho_damp_flag = 1;
rho_damp = utils::numeric(FLERR, arg[iarg + 1], false, lmp);
iarg++;
} else if (strcmp(arg[iarg], "artificial/visc") == 0) {
if (iarg + 1 >= narg) error->all(FLERR, "Illegal pair_style command");
if (iarg + 1 >= narg) utils::missing_cmd_args(FLERR, "pair rheo artificial/visc", error);
artificial_visc_flag = 1;
av = utils::numeric(FLERR, arg[iarg + 1], false, lmp);
iarg++;
@ -477,12 +475,12 @@ void PairRHEO::setup()
csq = fix_rheo->csq;
rho0 = fix_rheo->rho0;
if (h != fix_rheo->h)
error->all(FLERR, "Pair rheo cutoff {} does not agree with fix rheo cutoff {}", h, fix_rheo->h);
if (cutk != fix_rheo->cut)
error->all(FLERR, "Pair rheo cutoff {} does not agree with fix rheo cutoff {}", cutk, fix_rheo->cut);
hsq = h * h;
hinv = 1.0 / h;
hinv3 = hinv * 3.0;
cutksq = cutk * cutk;
cutkinv = 1.0 / cutk;
cutkinv3 = cutkinv * 3.0;
laplacian_order = -1;
int n = atom->ntypes;
@ -512,7 +510,7 @@ double PairRHEO::init_one(int i, int j)
if (setflag[i][j] == 0)
error->all(FLERR, "All pair rheo coeffs are not set");
return h;
return cutk;
}
/* ---------------------------------------------------------------------- */

4
src/RHEO/pair_rheo.h
View File

@ -37,8 +37,8 @@ class PairRHEO : public Pair {
void unpack_reverse_comm(int, int *, double *) override;
protected:
double h, *csq, *rho0; // From fix RHEO
double *cs, hsq, hinv, hinv3, av, rho_damp;
double cutk, *csq, *rho0; // From fix RHEO
double *cs, cutksq, cutkinv, cutkinv3, av, rho_damp;
int laplacian_order;
int artificial_visc_flag;