- Generated by
1.10.0
|
NEML2 1.4.0
|
The solid mechanics physics module is a collection objects serving as building blocks for composing material models for solids. Each category of the material model is explained below, with both the mathematical formulations and example input files.
Elasticity models describe the relationship between stress \( \boldsymbol{\sigma} \) and strain \( \boldsymbol{\varepsilon} \) without any history-dependent (internal) state variables. In general, the stress-strain relation can be written as
\[ \boldsymbol{\sigma} = \mathbb{C} : \boldsymbol{\varepsilon} \]
where \( \mathbb{C} \) is the fourth-order elasticity tensor. For linear isotropic elasticity, this relation can be simplified as
\[ \boldsymbol{\sigma} = 3 K \operatorname{vol} \boldsymbol{\varepsilon} + 2 G \operatorname{dev} \boldsymbol{\varepsilon} \]
where \( K \) is the bulk modulus, and \( G \) is the shear modulus.
Below is an example input file that defines a linear elasticity model.
Generally speaking, plasticity models describe (oftentimes irreversible and dissipative) history-dependent deformation of solid materials. The plastic deformation is governed by the plastic loading/unloading conditions (or more generally the Karush-Kuhn-Tucker conditions):
\begin{align*} f^p \leq 0, \quad \dot{\gamma} \geq 0, \quad \dot{\gamma}f^p = 0, \\ \end{align*}
where \( f^p \) is the yield function, and \( \gamma \) is the consistency parameter.
Consistent plasticity refers to the family of (macroscale) plasticity models that solve the plastic loading/unloading conditions (or the KKT conditions) exactly (up to machine precision).
Consistent plasticity models are sometimes considered rate-independent. But that is a misnomer as rate sensitivity can be baked into the yield function definition in terms of the rates of the internal variables.
Residual associated with the KKT conditions can be written as the Fischer-Burmeister complementarity condition
\begin{align*} r = \dot{\gamma} - f^p - \sqrt{{\dot{\gamma}}^2 + {f^p}^2}. \end{align*}
This complementarity condition is implemented by the object RateIndependentPlasticFlowConstraint. A complete example input file for consistent plasticity is shown below, and the composition and possible modifications are explained in the following subsections.
Viscoplasticity models regularize the KKT conditions by introducing approximations to the constraints. A widely adopted approximation is the Perzyna model where rate sensitivity is baked into the approximation following a power-law relation:
\begin{align*} \dot{\gamma} = \left( \dfrac{\left< f^p \right>}{\eta} \right)^n, \end{align*}
where \( \eta \) is the reference stress and \( n \) is the power-law exponent.
The Perzyna model is implemented by the object PerzynaPlasticFlowRate. A complete example input file for viscoplasticity is shown below, and the composition and possible modifications are explained in the following subsections.
The effective stress is a measure of stress describing how the plastic deformation "flows". For example, the widely-used \( J_2 \) plasticity uses the von Mises stress as the stress measure, i.e.,
\begin{align*} \bar{\sigma} &= \sqrt{3 J_2}, \\ J_2 &= \frac{1}{2} \operatorname{dev} \boldsymbol{\sigma} : \operatorname{dev} \boldsymbol{\sigma}. \end{align*}
Commonly used stress measures are defined using SR2Invariant.
For perfectly plastic materials, the yield function only depends on the effective stress and a constant yield stress, i.e., the envelope does not shrink or expand depending on the loading history.
\begin{align*} f^p &= \bar{\sigma} - \sigma_y. \end{align*}
Below is an example input file defining a perfectly plastic yield function with \( J_2 \) flow.
The equivalent plastic strain \( \bar{\varepsilon}^p \) is a scalar-valued internal variable that can be introduced to control the shape of the yield function. The isotropic strain hardening \( k \) is controlled by the accumulated equivalent plastic strain, and enters the yield function as
\begin{align*} f^p &= \bar{\sigma} - \sigma_y - k(\bar{\varepsilon}^p). \end{align*}
Below is an example input file defining a yield function with \( J_2 \) flow and linear isotropic hardening.
The kinematic plastic strain \( \boldsymbol{K}^p \) is a tensor-valued internal variable that can be introduced to control the shape of the yield function. The kinematic hardening \( X \) is controlled by the accumulated kinematic plastic strain, and the effective stress is defined in terms of the over stress. In case of \( J_2 \), the effective stress can be rewritten as
\begin{align*} \bar{\sigma} &= \sqrt{3 J_2}, \\ J_2 &= \frac{1}{2} \operatorname{dev} \boldsymbol{\Xi} : \operatorname{dev} \boldsymbol{\Xi}, \\ \boldsymbol{\Xi} &= \boldsymbol{\sigma} - \boldsymbol{X}. \end{align*}
Below is an example input file defining a yield function with \( J_2 \) flow and linear kinematic hardening.
An alternative way of introducing hardening is through back stresses. Instead of modeling the accumulation of kinematic plastic strain, back stress models directly describe the evolution of back stress. An example input file defining a yield function with \( J_2 \) flow and two back stresses is shown below.
Isotropic hardening, kinematic hardening, and back stresses are all optional and can be "mixed" in the definition of a yield function. The example input file below shows a yield function with \( J_2 \) flow, isotropic hardening, kinematic hardening, and two back stresses.
Flow rules are required to map from the consistency parameter \( \gamma \) to various internal variables describing the state of hardening, such as the equivalent plastic strain \( \bar{\varepsilon}^p \), the kinematic plastic strain \( \boldsymbol{K}^p \), and the back stress \( \boldsymbol{X} \).
Associative flow rules define flow directions variationally according to the principle of maximum dissipation, i.e.,
\begin{align*} \dot{\boldsymbol{\varepsilon}}^p &= \dot{\gamma} \dfrac{\partial f^p}{\partial \boldsymbol{\sigma}}, \\ \dot{\bar{\varepsilon}}^p &= -\dot{\gamma} \dfrac{\partial f^p}{\partial k}, \\ \dot{\boldsymbol{K}}^p &= \dot{\gamma} \dfrac{\partial f^p}{\partial \boldsymbol{X}}. \end{align*}
The example input file below defines associative \( J_2 \) flow rules
In the above example, a model named "normality" is used to compute the associative flow directions, and the rates of the internal variables are mapped using the rate of the consistency parameter and each of the associative flow direction. The cross-referenced model named "flow" (omitted in the example) is the composition of models defining the yield function \( f^p \) in terms of the variational arguments \( \boldsymbol{\sigma} \), \( k \), and \( \boldsymbol{X} \).
NEML2 models can be setup to take as input strain and provide stress or take as input stress and provide strain. With some additional work, a model that "natively" works in either strain or stress control can be modified to work under mixed stress/strain control.
Mixed control here means that some components of the strain tensor and some components of the stress tensor are provided as input and the model must solve for the conjugate, missing components of stress or strain.
Mathematically, at each time step consider a tensor of applied, mixed strain or stress conditions \( f_{ij} \) and a control signal \( c_{ij} \). When \( c_{ij} < h \) for some threshold \( h \) we consider the corresponding component of the input \( f_{ij} \) to be a strain value and the model must solve for the corresponding value of stress. If \( c_{ij} \ge h \) we consider the corresponding component of the input \( f_{ij} \) to be a stress value and the model must solve for the corresponding value of strain.
Modifying a model for mixed control only requires a few additional objects. The first maps from the mixed input and a conjugate mixed state vector into the actual model stress and strain input axes:
The second modification renames the resulting residual to map to the correct name for the mixed state:
In this example the "base" model is setup for strain control, so that the residual is formed on stress. For "base" stress control the only change would be the name of the from parameter in this object.
These two modifications will allow the model to be run in mixed control. One additional modification clarifies the output by mapping the mixed input and state back to stress and strain tensors
NEML2 adopts an incremental rate-form view of crystal plasticity. The fundemental kinematics derive from the rate expansion of the elastic-plastic multiplactive split:
\begin{align*} F = F^e F^p \end{align*}
where the spatial velocity gradient is then
\begin{align*} l = \dot{F} F^{-1} = \dot{F}^e F^{e-1} + F^e \dot{F}^{p} {F}^{p-1} F^{e-1} \end{align*}
The plastic deformation \( \bar{l}^p = \dot{F}^{p} {F}^{p-1} \) defines the crystal plasticity kinemtics and NEML2 assumes that the elastic stretch is small ( \( F^e = \left(I + \varepsilon \right) R^e \)) so that spatial velocity gradient becomes
\begin{align*} l = \dot{\varepsilon} + \Omega^e - \Omega^e \varepsilon + \varepsilon \Omega^e + l^p + \varepsilon l^p - l^p \varepsilon \end{align*}
defining \( l^p = R^e \bar{l}^p R^{eT} \) as the constitutive plastic velocity gradient rotated into the current configuration and \( \Omega^e = \dot{R}^e R^{eT} \) as the elastic spin and assuming that
The first assumption is accurate for metal plasticity, the second assumption is more questionable if the material deforms at a fast strain rate.
Define the current orientation of a crystal as the composition of its initial rotation from the crystal system to the lab frame and the elastic rotation, i.e.
\begin{align*} Q = R^e Q_0 \end{align*}
and note with this defintion we can rewrite the spin equation:
\begin{align*} \Omega^e = \dot{R}^e R^{eT} = \dot{Q} Q_0^T Q_0 Q^T = \dot{Q} Q^T \end{align*}
With this definition and the choice of kinematics above we can derive evolution equations for our fundemental constitutive quantities, the elastic stretch \( \varepsilon \) and the orientation $Q$ by splitting the spatial velocity gradient into symmetric and skew parts and rearranging the resulting equations:
\begin{align*} \dot{\varepsilon}= d -d^p-\varepsilon w + w \varepsilon \\ \dot{Q} = \left(w - w^p - \varepsilon d^p + d^p \varepsilon\right) Q . \end{align*}
where \(l = d + w\) and \(l^p = d^p + w^p\)
For most (or all) choices of crystal plasticity constitutive models we also need to define the Cauchy stress as:
\begin{align*} \sigma = C : \varepsilon \end{align*}
with \( C \) a (generally) anisotropic crystal elasticity tensor rotated into the current configuration.
The crystal plasticity examples in NEML2 integrate the elastic strain (and the constitutive internal variables) using a backward Euler integration rule and integrate the crystal orientation using either an implicit or explicit exponential rule. These integrate rules can either be coupled or decoupled, i.e. integrated together in a fully implicit manner or first integrate the strain and internal variables and then sequentially integrating the rotations.
A full constitutive model must then define the plastic deformation $l^p$ and whatever internal variables are used in this definition. A wide variety of choices are possible, but the examples use the basic assumption of Asaro:
\begin{align*} l^p = \sum_{i=1}^{n_{slip}} \dot{\gamma}_i Q \left(d_i \otimes n_i \right) Q^T \end{align*}
where now \( \dot{\gamma}_i \), the slip rate on each system, is the constitutive chioce. NEML provides a variety of options for defining these slip rates in terms of internal hardening variables and the results shear stress
\begin{align*} \tau_i = \sigma : Q \operatorname{sym}\left(d_i \otimes n_i \right) Q^T \end{align*}
Ancillary classes automatically generate lists of slip and twin systems from the crystal sytem, so the user does not need to manually provide these themselves.
NEML2 uses modified Rodrigues parameters to define orientations internally. These can be converted to Euler angles, quaternions, etc. for output.
1.10.0