LIGGGHTS(R)-PUBLIC WWW Site - LIGGGHTS(R)-PUBLIC Documentation - LIGGGHTS(R)-PUBLIC Commands

fix continuum/weighted command

Syntax:

fix ID group-ID continuum/weighted keyword value

Examples:

fix 1 all continuum/weighted kernel_radius 0.01 compute stress
fix 1 all continuum/weighted kernel_radius 0.2 kernel_type Wendland compute stress_strain 

Description:

If the compute keyword is set to either stress or stress_strain this fix calculates the complete stress tensor at each particle according to Goldhirsch. The formula is given by:

  sigma_{i,ab} = -1/2 sum_{j,k} f_{jk,a} r_{jk,b} integral_0^1 phi(r_i - r_j + s r_{jk}) ds
               - sum_j m_j v'_{ij,a} v'_{ij,b} phi(r_i - r_j)
               + sigma_wall_{i,ab} 

where

  v'_{ij} = v_j - V_i
  V_i = p_i / rho_i
  p_i = sum_j m_j v_j phi(r_i - r_j)
  rho_i = sum_j m_j phi(r_i - r_j)
  phi(r) = H(radius - |r|) / Omega(radius)
  Omega(radius) = 4/3 pi radius^3 

and

  v is the velocity vector
  f_{ij} force acting from j onto i
  r_{ij} vector from center of j to center of i
  m_i mass of particle i
  H is the Heavyside function 

In case solid boundaries are present the last term is given according to Weinhart et al. by

    sigma_wall_{i,ab} = - sum_{j,k} f_{jk,a} a_{jk,b} integral_0^1 phi(r_i - r_j + s a_{jk}) ds 

where

    a_{jk} = r_j - c_{jk} 

and c_{jk} is the contact point of particle j with wall k and the sum runs over all particles j and walls k.

If the compute keyword is set to either strain or stress_strain this fix calculates the incremental strain tensor at each particle according to Zhang et al. The formula is given by

  epsilon_i,ab = 1/(2 rho_i) sum_j,k m_j m_k phi(r_{ij}) dt (v_{jk,a} grad_phi(r_{ik},b) + v_{jk,b} grad_phi(r_{ik},a)) 

where most of the variables are given as above and additionally

    v_{ij,a} is the a-th component of the velocity difference between i and j
    grad_phi(r_{ij},a) is the a-th component of the gradient of phi with respect to r_i
    dt is the time-step size 

The following three kernel types are implemented at the moment:

Note that all kernels are equal to zero if r > 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.


Restart, fix_modify, output, run start/stop, minimize info:

No information about this fix is written to binary restart files. None of the fix_modify options are relevant to this fix.

The values can be dumped by using the f_stressTensor_[i] and/or f_strainTensor_[i] (0 <= i <= 8) values in dump commands

No parameter of this fix can be used with the start/stop keywords of the run command. This fix is not invoked during energy minimization.

Restrictions:

Strain computation does not work with the default TOP_HAT kernel as its derivative is zero.

In order to ensure that all particles in the kernel radius are considered make use of the neigh_modify command. In particular the contact_distance_factor which should be set such that

  2*min(radius)*contact_distance_factor >= kernel_radius 

Related commands:

compute, compute stress/atom, fix ave/atom, fix ave/histo, fix ave/time, fix ave/spatial, fix ave/euler

Default: none

References:

(Goldhirsch) Goldhirsch; Stress, stress asymmetry and couple stress: from discrete particles to continuous fields, Granular Matter (2010)

(Weinhart) Weinhart, Thornton, Luding, Bokhove; From discrete particles to continuum fields near a boundary (2012)

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