calculate strain command

Purpose

Calculates the strain of a collection of particles

Synatx

calculate strain keyword/value pairs

The following keyword/value pairs are allowed:

Keywords	Description
id	user-assigned name for the calculate strain command
region	region in which the strain tensor is calculated
particle_group	particle group for which the strain tensor is calculated
kernel_radius	radius of the smoothing kernel default: 4*maximum radius; range: (0,∞); units: [length]
kernel_type	type of the smoothing kernel, must be one of: `Gaussian`, `Wendland` default: Gaussian;

Examples

calculate strain id cs
calculate strain id cs kernel_radius 1e-2 kernel_type Wendland region outlet
# calculate average xx strain:
calculate average id cs_xx quantity id_strainTensor_cs.xx

Description

The strain tensor is calculated for each particle according to Zhang et al. The formula is given by

$\epsilon_{i,ab} = \frac{1}{2 rho_i} \sum_{j,k} m_j m_k \phi(r_{ij}) dt (v_{jk,a} \nabla\phi(r_{ik},b) + v_{jk,b} \nabla\phi(r_{ik},a))$

where most of the variables are given as above and additionally

 is the a-th component of the velocity difference between i and j
 is the a-th component of the gradient of phi with respect to 
dt is the time-step size

The following three kernel_type values are implemented at the moment:

Gaussian - Gaussian kernel

$w(r) = a_g exp(-q^2) \quad (q = 3 r / kernel\_radius)$

Wendland - Quintic radial polynomial

$w(r) = a_w (1-q/2)^4 (1+2q) \quad (q = 2 r / kernel\_radius)$

Note that all kernels are equal to zero if $r$ is greater than kernel_radius (this implies a cut-off for the Gaussian). The constants $a$ (different for each kernel) are chosen such that the integral of $w$ over the ball of radius kernel_radius is equal to one. In case of the top hat kernel $a_t$ is equal to the volume of this sphere.

The command can be restricted to work on a specific group of particles or a specific region by using the particle_group or region keywords, respectively.

The strain tensor is written to a per-particle property which can be accessed using id_strainTensor_ID, where ID is the id of the calculate command. It has nine components which can be accessed using .xx, .xy, .xz, .yx, .yy, .yz, .zx, .zy and .zz for the respective components.

References

(Zhang) Zhang, Behringer, Goldhirsch; Coarse-Graining of a Physical Granular System, Progress of Theoretical Physics Supplement (2010)