Plasticity¶

Overview¶

Plasticity models describe the (typically irreversible and dissipative) history-dependent deformation of solids. NEML2’s macroscale plasticity catalog under neml2/models/solid_mechanics/plasticity/ provides small composable ingredients that you wire together with the Syntax catalog primitives in common/ to assemble a return-map step:

Stress measures — IsotropicMandelStress, MandelStress, and SR2Invariant (from common/) project the Cauchy/Mandel stress onto the effective stress used by the yield surface.
Yield surfaces — YieldFunction for classical \(J_2\)-style envelopes, GTNYieldFunction and GursonCavitation for porous metal plasticity.
Hardening models — isotropic (LinearIsotropicHardening, VoceIsotropicHardening, SlopeSaturationVoceIsotropicHardening, PowerLawIsotropicHardeningStaticRecovery), kinematic (LinearKinematicHardening, ChabochePlasticHardening, FredrickArmstrongPlasticHardening, PowerLawKinematicHardeningStaticRecovery), and rate-temperature couplings built from the Kocks-Mecking family (KocksMeckingActivationEnergy, KocksMeckingFlowSwitch, KocksMeckingFlowViscosity, KocksMeckingIntercept, KocksMeckingRateSensitivity, KocksMeckingYieldStress).
Flow directions — Normality takes the gradient of the yield function with respect to its variational arguments, or AssociativeJ2FlowDirection provides a closed-form \(J_2\) direction.
Flow rules and rates — AssociativePlasticFlow and the hardening-variable rates AssociativeIsotropicPlasticHardening / AssociativeKinematicPlasticHardening multiply the consistency parameter by the corresponding associative direction. For rate-dependent flow, PerzynaPlasticFlowRate gives the consistency parameter itself as an explicit function of the yield function.

These pieces are typically glued together by ComposedModel and closed by either FBComplementarity (consistent plasticity) or a Perzyna-style rate equation (viscoplasticity), then wrapped in ImplicitUpdate for the implicit return map. A closed-form radial-return path is also available for the J2-elastic-perfectly-plastic case via LinearIsotropicElasticJ2TrialStressUpdate; see the rate_independent_plasticity/radial_return/ regression scenario.

Math¶

Plastic flow is governed by the Karush-Kuhn-Tucker (KKT) loading/unloading conditions

\[ f^p \leq 0, \quad \dot{\gamma} \geq 0, \quad \dot{\gamma}\, f^p = 0, \]

where \(f^p\) is the yield function and \(\gamma\) is the consistency parameter.

Consistent (rate-independent) plasticity¶

NEML2 enforces the KKT conditions by recasting them as the smooth Fischer-Burmeister complementarity residual, which the nonlinear solver drives to its convergence tolerance:

\[ r = \dot{\gamma} - f^p - \sqrt{\dot{\gamma}^{\,2} + {f^p}^{\,2}}, \]

implemented by FBComplementarity with a_inequality = 'LE'.

Note

“Consistent” plasticity is often called rate-independent. Rate sensitivity can still enter through the yield function or hardening laws expressed in terms of internal-variable rates, so the “consistent” label is the more precise one.

Viscoplasticity¶

The Perzyna regularization replaces the complementarity with an explicit power-law overstress relation,

\[ \dot{\gamma} = \left( \frac{\langle f^p \rangle}{\eta} \right)^{n}, \]

where \(\eta\) is the reference stress, \(n\) is the rate exponent, and \(\langle \cdot \rangle\) denotes the Macaulay bracket. This is the role of PerzynaPlasticFlowRate.

Effective stress and yield function¶

For classical \(J_2\) plasticity the effective stress is the von Mises norm of the deviatoric (over-)stress,

(16)¶\[\begin{align} \bar{\sigma} &= \sqrt{3 J_2}, \\ J_2 &= \tfrac{1}{2}\,\mathrm{dev}\,\boldsymbol{\Xi}\,:\, \mathrm{dev}\,\boldsymbol{\Xi}, \\ \boldsymbol{\Xi} &= \boldsymbol{\sigma} - \sum_i \boldsymbol{X}_i. \end{align}\]

The yield function combines the effective stress with the yield strength \(\sigma_y\), isotropic hardening \(k(\bar{\varepsilon}^p)\), and (implicitly, through \(\boldsymbol{\Xi}\)) the kinematic back stresses \(\boldsymbol{X}_i\):

\[ f^p = \bar{\sigma} - \sigma_y - k(\bar{\varepsilon}^p). \]

Associative flow¶

Associative flow rules derive the rates of the plastic strain and internal variables from the maximum-dissipation principle,

(17)¶\[\begin{align} \dot{\boldsymbol{\varepsilon}}^p &= \dot{\gamma}\, \frac{\partial f^p}{\partial \boldsymbol{\sigma}}, \\ \dot{\bar{\varepsilon}}^p &= -\dot{\gamma}\, \frac{\partial f^p}{\partial k}, \\ \dot{\boldsymbol{K}}^p &= \dot{\gamma}\, \frac{\partial f^p}{\partial \boldsymbol{X}}. \end{align}\]

Normality computes these gradients by symbolically differentiating a sub-model (the yield function) with respect to a list of input variables; AssociativePlasticFlow, AssociativeIsotropicPlasticHardening, and AssociativeKinematicPlasticHardening then scale the directions by \(\dot{\gamma}\).

Example: rate-independent \(J_2\) plasticity with mixed hardening¶

The input file below assembles a fully implicit return-map for elastic-plastic behaviour with both isotropic and kinematic hardening: linear-elastic stress-strain, VoceIsotropicHardening isotropic hardening, LinearKinematicHardening back stress, associative \(J_2\) flow, and the Fischer-Burmeister consistency condition.

Listing 23 tests/regression/solid_mechanics/rate_independent_plasticity/isokinharden/model.i¶

# neml2
[Tensors]
  [end_time]
    type = Python
    expr = 'Scalar(torch.logspace(-1.0, 5.0, 20, dtype=torch.float64))'
  []
  [times]
    type = Python
    expr = 'Scalar(end_time.data.unsqueeze(0) * torch.linspace(0.0, 1.0, 100, dtype=torch.float64).unsqueeze(-1))'
  []
  [max_strain]
    type = Python
    expr = 'SR2.fill(0.1, -0.05, -0.05, 0.0, 0.0, 0.0).dynamic_batch.expand(20)'
  []
  [strains]
    type = Python
    expr = 'SR2(max_strain.data.unsqueeze(0) * torch.linspace(0.0, 1.0, 100, dtype=torch.float64).reshape(100, 1, 1))'
  []
[]

[Drivers]
  [driver]
    type = TransientDriver
    model = 'model'
    prescribed_time = 'times'
    force_SR2_names = 'E'
    force_SR2_values = 'strains'
    save_as = 'result.pt'
  []
  [regression]
    type = TransientRegression
    driver = 'driver'
    reference = 'gold/result.pt'
  []
[]

[Models]
  [isoharden]
    type = VoceIsotropicHardening
    saturated_hardening = 100
    saturation_rate = 1.2
  []
  [kinharden]
    type = LinearKinematicHardening
    hardening_modulus = 1000
    back_stress = 'X'
  []
  [elastic_strain]
    type = SR2LinearCombination
    from = 'E plastic_strain'
    to = 'elastic_strain'
    weights = '1 -1'
  []
  [elasticity]
    type = LinearIsotropicElasticity
    coefficients = '1e5 0.3'
    coefficient_types = 'YOUNGS_MODULUS POISSONS_RATIO'
    strain = 'elastic_strain'
  []
  [mandel_stress]
    type = IsotropicMandelStress
    cauchy_stress = 'stress'
  []
  [overstress]
    type = SR2LinearCombination
    to = 'O'
    from = 'mandel_stress X'
    weights = '1 -1'
  []
  [vonmises]
    type = SR2Invariant
    invariant_type = 'VONMISES'
    tensor = 'O'
    invariant = 'effective_stress'
  []
  [yield_surface]
    type = YieldFunction
    yield_stress = 1000
    isotropic_hardening = 'isotropic_hardening'
  []
  [flow]
    type = ComposedModel
    models = 'overstress vonmises yield_surface'
  []
  [normality]
    type = Normality
    model = 'flow'
    function = 'yield_function'
    from = 'mandel_stress X isotropic_hardening'
    to = 'flow_direction kinematic_hardening_direction isotropic_hardening_direction'
  []
  [eprate]
    type = AssociativeIsotropicPlasticHardening
  []
  [Kprate]
    type = AssociativeKinematicPlasticHardening
  []
  [Eprate]
    type = AssociativePlasticFlow
  []
  [integrate_ep]
    type = ScalarBackwardEulerTimeIntegration
    variable = 'equivalent_plastic_strain'
  []
  [integrate_Kp]
    type = SR2BackwardEulerTimeIntegration
    variable = 'kinematic_plastic_strain'
  []
  [integrate_Ep]
    type = SR2BackwardEulerTimeIntegration
    variable = 'plastic_strain'
  []
  [consistency]
    type = FBComplementarity
    a = 'yield_function'
    a_inequality = 'LE'
    b = 'flow_rate'
  []
  [surface]
    type = ComposedModel
    models = 'isoharden kinharden elastic_strain elasticity
              mandel_stress overstress vonmises
              yield_surface normality eprate Kprate Eprate
              consistency integrate_ep integrate_Kp integrate_Ep'
  []
[]

[EquationSystems]
  [eq_sys]
    type = NonlinearSystem
    model = 'surface'
    unknowns = 'plastic_strain kinematic_plastic_strain equivalent_plastic_strain flow_rate'
    residuals = 'plastic_strain_residual kinematic_plastic_strain_residual equivalent_plastic_strain_residual complementarity'
  []
[]

[Solvers]
  [newton]
    type = Newton
    linear_solver = 'lu'
  []
  [lu]
    type = DenseLU
  []
[]

[Models]
  [predictor]
    type = ConstantExtrapolationPredictor
    unknowns_SR2 = 'plastic_strain kinematic_plastic_strain'
    unknowns_Scalar = 'equivalent_plastic_strain flow_rate'
  []
  [return_map]
    type = ImplicitUpdate
    equation_system = 'eq_sys'
    solver = 'newton'
    predictor = 'predictor'
  []
  [model]
    type = ComposedModel
    models = 'return_map elastic_strain elasticity'
    additional_outputs = 'plastic_strain kinematic_plastic_strain equivalent_plastic_strain'
  []
[]

How the pieces wire together¶

Hardening laws — [isoharden] (VoceIsotropicHardening) maps equivalent_plastic_strain to isotropic_hardening; [kinharden] (LinearKinematicHardening) maps kinematic_plastic_strain to the back-stress contribution X.
Stress side — [elastic_strain] (SR2LinearCombination) subtracts the plastic strain from the total strain, [elasticity] (LinearIsotropicElasticity) gives the Cauchy stress, and [mandel_stress] (IsotropicMandelStress) lifts it to the work-conjugate Mandel stress consumed by the yield surface.
Yield surface — [overstress] builds \(\boldsymbol{\Xi} = \boldsymbol{\sigma} - \boldsymbol{X}\) as O, [vonmises] (SR2Invariant) produces the effective stress, and [yield_surface] (YieldFunction) assembles \(f^p\) from the effective stress, yield strength, and isotropic hardening.
Flow — [flow] is a ComposedModel glueing the overstress, invariant, and yield steps so that Normality can differentiate the resulting yield_function with respect to mandel_stress, X, and isotropic_hardening, producing the three associative directions. [Eprate], [Kprate], [eprate] scale each direction by the consistency parameter flow_rate.
Time integration — three time-integration blocks (SR2BackwardEulerTimeIntegration and ScalarBackwardEulerTimeIntegration) turn the rates into the residuals plastic_strain_residual, kinematic_plastic_strain_residual, equivalent_plastic_strain_residual.
Consistency — [consistency] (FBComplementarity) closes the system with the Fischer-Burmeister residual between yield_function and flow_rate.
Return map — [surface] composes every residual-producing piece into one model. NonlinearSystem solves for plastic_strain, kinematic_plastic_strain, equivalent_plastic_strain, flow_rate using Newton + DenseLU. [return_map] (ImplicitUpdate) wraps the equation system, warm-started by [predictor] (ConstantExtrapolationPredictor). The top-level [model] re-evaluates elastic_strain and elasticity with the converged plastic strain to expose the final stress.

Switching to viscoplasticity¶

Swapping the FB closure for a direct Perzyna flow-rate evaluation is the basic move from rate-independent to viscoplastic; depending on the variation you may also integrate stress in rate form rather than plastic strain. The reference recipe lives at tests/regression/solid_mechanics/viscoplasticity/perfect/model.i.

Variations shipped in the test tree¶

The tests/regression/solid_mechanics/ tree carries ready-to-read compositions for the most common dialects:

rate_independent_plasticity/perfect/ — perfectly plastic, no hardening.
rate_independent_plasticity/isoharden/ — isotropic hardening only.
rate_independent_plasticity/kinharden/ — kinematic hardening only.
rate_independent_plasticity/isokinharden/ — combined (used above).
rate_independent_plasticity/radial_return/ — closed-form \(J_2\) radial-return via LinearIsotropicElasticJ2TrialStressUpdate.
rate_independent_plasticity/gurson/ — porous metal plasticity using GTNYieldFunction and GursonCavitation.
viscoplasticity/perfect/ — Perzyna viscoplastic counterpart.
recovery/ — adds isotropic/kinematic static recovery using PowerLawIsotropicHardeningStaticRecovery and PowerLawKinematicHardeningStaticRecovery.
km_flow/ — Kocks-Mecking temperature- and rate-dependent flow.

Plasticity¶

Overview¶

Math¶

Consistent (rate-independent) plasticity¶

Viscoplasticity¶

Effective stress and yield function¶

Associative flow¶

Example: rate-independent \(J_2\) plasticity with mixed hardening¶

How the pieces wire together¶

Switching to viscoplasticity¶

Variations shipped in the test tree¶

See also¶