fix powder/update command

Purpose

Powder/update can be used to simulate extreme plastic compression of particles (e.g. tabletting, uniaxial compression)

Warning

GPU support for this command has not been tested and may not work as expected.

Syntax

fix ID group powder/update keyword value

ID, group are documented in fix command
powder/update = style name of this fix command
zero or more keyword/value pairs may be appended
keyword = alpha_relaxation_time or use_hardening or pressure_avg_time or slow_compression
keyword/value = optional keyword-value pair

alpha_relaxation_time keyword = time
  time = a limiting value that ensures that  (the compression ratio) does not change too much over consecutive time steps
use_hardening keyword = 'yes', 'no', 'base', 'exponent', 'exponential'
  no = do not use a hardening law
  'yes', 'base', 'exponent', 'exponential' = specific hardening model is activated
pressure_avg_time keyword = time
  time = weight in the sliding average of pressure computation
slow_compression keyword = 'yes' or 'no'
  yes, no = determines whether particles will compress instantaneously
elastic_dilation value = 'yes' or 'no'
  yes, no = determines whether an elastic component in the compression is present
non_linear value = 'yes' or 'no'
  yes, no = determines whether hardening is linear in compression density

Examples

fix powder all powder/update
fix powder all powder/update alpha_relaxation_time 0.5 pressure_avg_time 0.01 use_hardening yes elastic_dilation yes non_linear yes

Associated material properties

Material properties

alphaI_min ( $\alpha_{min}$ ): minimum compression density [ $\cdot$ ]
compressibility ( $c$ ): measure on how easy it is to compress a particle [ $\cdot$ ]
dilatability ( $d$ ): coefficient for determining how much a particle dilates [ $\cdot$ ]
hardening ( $h$ ): coefficient for modifying the compressibility for compacted particles. Only required if use_hardening is not no [ $\cdot$ ]
compression_exponent ( $n$ ): exponent for the non-linear hardening model (only required if non_linear yes) [ $\cdot$ ]

Description

This fix implements a powder compression model. The current compression state of each particle is given by a scalar $\alpha$ . This scalar can be imagined as a relative particle density.

All particles are initialized with an $\alpha$ equal to $\alpha_{min}$ which corresponds to an uncompressed state. If a particle is compressed its $\alpha$ will reach up to a maximum of 1 which corresponds to the fully compressed state.

During a simulation a simplistic instantaneous particle pressure is computed based on the contact forces, i.e.

$p_{this_i} = \frac{1}{4 r_i \pi Y_i} \sum_j f_{ij} \cdot n_{ij}$

where f are the forces and n the normals between two particles and Y is the Youngs Modulus. The pressure is computed as sliding average where the current particle is weighted by $dt/\tau_{p_{avg}}$ , where $\tau_{p_{avg}}$ is the pressure_avg_time parameter, i.e.

$w = \frac{dt}{\tau_{p_{avg}}}\\ p_i = p_{this_i} w + p_i (1-w)$

If the previous maximum pressure ( $p_{max_i}$ ) exceeds the current particle pressure, which was initialized as zero, then the particle is compressed according to the compressibility scalar and the pressure difference, i.e.

$\alpha^n += \alpha^{n-1} + c (p_i - p_{max_i})\\ p_{max_i} = p_i$

In case the pressure is lower than zero (e.g. when using cohesion) the particle is dilated according to the pressure difference and the dilatability * compressibility scalar, i.e.

$\alpha^n = \alpha^{n-1} + d\,c\,p_i\\ p_{max_i} = p_{max_i} - p$

where the change in $\alpha$ is negative as $p_i$ is.

If the keyword elastic_dilation is set, the dilation will also occur if the pressure is lower than the previously experienced $p_{max_i}$ . After such an elastic dilation the particle can be compressed again. However, in this case it will be first compressed with a rate equal to dilatability * compressibility until the all-time maximum $p_{max_i}$ is attained again. Then, the normal compression rate of compressibility is assumed.

Note that in the case of dilation, $\alpha$ is actually not set directly. Instead a target value of $\alpha$ is defined. The real value of $\alpha$ is then allowed to slowly approach this target value. The maximum rate of change is governed by the alpha_relaxation_time value that is specified with the fix. Within one time-step $\alpha$ is allowed to change by dt/alpha_relaxation_time, thus this parameter specifies the time scale over which $\alpha$ relaxes to its uncompressed minimum value. In order to apply the same algorithm to the compression case, the slow_compression needs to be set. This option will also use the alpha_relaxation_time as time scale to allow particles to compress over several time steps.

In case use_hardening is activated, the compressibility parameter is varied according to $\alpha$ . There are three different models available. For use_hardening yes or use_hardening base the model is

$c \rightarrow c\left(1 - h \omega(\alpha)\right), \omega(\alpha) = \frac{\alpha - \alpha_{min}}{1 - \alpha_{min}}$

with a hardening set to 0.5 this would imply that the compressibility is only half as large if the particles are fully compressed. If additionally non_linear is activated a compression exponent $n$ alters the equation above to read

$c \rightarrow c \left( 1 - h \left(\omega(\alpha)\right)^n\right)$

This allows for a better control over the loading curve.

The second hardening model can be activated using use_hardening exponent, where the linear model is given by

$c \rightarrow c^{1 - h \omega(\alpha)}$

and the non-linear version is

$c \rightarrow c^{1 - \left(h \omega(\alpha)\right)^n}.$

The non-linear model in particular is not recommended as the hardening parameter is rather useless as it is usually chosen to be close to one.

The third and final model can be activated with use_hardening exponential, where the linear model reads

$c \rightarrow c \exp\left(- 1 - h \omega(\alpha)\right),$

and the non-linear version is

$c \rightarrow c \exp\left(- 1 - \left(h \omega(\alpha)\right)^n\right).$

Restart, fix_modify, output, run start/stop

No information about this fix is written to binary restart files.

fix_modify cannot be used on the parameter of this fix.

This command adds a per-particle pressure vector, which can be accessed via id_pressure. Its components are average (index 1), max_lower (index 2) and max_lower (index 3). These can be accessed via id_pressure.component_name or id_pressure[index], and will be output automatically as a particle property via the output_settings command.

Restrictions

None.

Default

alpha_relaxation_time = 0.01, pressure_avg_time = 0.01, use_hardening = ‘no’, slow_compression = ‘no’