Kinematics¶

Overview¶

The kinematics catalog supplies the small, composable primitives that encode non-mechanical contributions to a material’s deformation — the strains or volumetric changes that arise from temperature excursions, phase transformations, density changes during sintering, or freezing/melting of a contained fluid. Pair them with an elastic or inelastic constitutive model, which subtracts the non-mechanical contribution from the total kinematic measure to recover the mechanical part.

The catalog carves the space along two axes:

Output measure. Eigenstrain primitives return a symmetric rank-2 tensor eigenstrain : SR2 suitable for the small-strain additive split \(\boldsymbol{\varepsilon} = \boldsymbol{\varepsilon}^{\text{mech}} + \boldsymbol{\varepsilon}^{*}\). Deformation Jacobian primitives return a jacobian : Scalar \(J^{*}\) suitable for the finite-strain multiplicative split \(\boldsymbol{F} = \boldsymbol{F}^{\text{mech}} \boldsymbol{F}^{*}\) (combined with VolumeAdjustDeformationGradient to lift the scalar Jacobian into a deformation gradient correction).
Driving quantity. Temperature, volume change, phase fraction, or a coupled swelling/phase-change state.

For the small-strain workflow these eigenstrains are typically wired into an SR2LinearCombination that subtracts every non-mechanical strain from the total strain before the elastic constitutive call. For the finite-strain workflow the deformation Jacobian primitives are multiplied together and then fed through VolumeAdjustDeformationGradient to produce the volume-corrected deformation gradient.

Math¶

Eigenstrain (small-strain additive split)¶

For an isotropic eigenstrain driven by a scalar field \(q\) (temperature, volume, phase fraction), NEML2 adopts the standard cumulative form

\[ \boldsymbol{\varepsilon}^{*}(q) = \varepsilon^{*}_v(q)\,\boldsymbol{I}, \]

where \(\varepsilon^{*}_v\) is the (linear) volumetric eigenstrain and \(\boldsymbol{I}\) is the rank-2 identity. The three concrete leaves in the catalog differ only in how \(\varepsilon^{*}_v\) depends on its driver:

(21)¶\[\begin{align} \text{Thermal:}\quad &\varepsilon^{*}_v = \alpha\,(T - T_0), \\ \text{Volume change:}\quad &\varepsilon^{*}_v = \left(\dfrac{V}{V_0}\right)^{1/3} - 1, \\ \text{Phase transformation:}\quad &\varepsilon^{*}_v = \Delta V \, f. \end{align}\]

Here \(\alpha\) is the coefficient of thermal expansion, \(T_0\) the stress-free reference temperature, \(V_0\) the reference volume, \(\Delta V\) the volume fraction change between phases A and B, and \(f \in [0,1]\) the current phase fraction. Each variant exposes the same eigenstrain output name, which lets a downstream model treat them interchangeably.

The mechanical (elastic) strain is then recovered by the additive split

\[ \boldsymbol{\varepsilon}^{\text{el}} = \boldsymbol{\varepsilon} - \boldsymbol{\varepsilon}^{p} - \boldsymbol{\varepsilon}^{*}, \]

with any plastic strain \(\boldsymbol{\varepsilon}^{p}\) also subtracted in inelastic flows.

Deformation Jacobian (finite-strain multiplicative split)¶

In the finite-strain setting NEML2 represents the same physics through a scalar deformation Jacobian \(J^{*}\) — the determinant of the volumetric correction to be removed from the total deformation gradient:

(22)¶\[\begin{align} \text{Thermal:}\quad &J^{*} = 1 + \alpha\,(T - T_0), \\ \text{Swelling + phase change:}\quad &J^{*} = 1 + \alpha\,c\,\phi^{f} + (1-c)\,\phi^{f}\,\Delta\Omega. \end{align}\]

For the coupled swelling/phase-change case, \(\phi^{f}\) is the fluid volume fraction that participates in the swelling, \(c \in [0,1]\) is the phase fraction (0 for fully solid, 1 for fully liquid), \(\alpha\) is the swelling coefficient, and \(\Delta\Omega\) is the relative difference in reference volume between the two phases.

Once the per-mechanism Jacobians are multiplied together into a single \(J\), VolumeAdjustDeformationGradient lifts the scalar back to a rank-2 correction on the deformation gradient,

\[ \boldsymbol{F}^{\text{mech}} = J^{-1/3}\,\boldsymbol{F}, \]

so that the mechanical part carries only the deviatoric (shape-changing) piece while the volumetric piece \(J\) is removed before the elastic call.

Example: thermal eigenstrain in a free-sintering model¶

The free_sintering regression scenario composes ThermalEigenstrain into a full GTN poroplastic free-sintering model. The eigenstrain enters through the elastic-strain combination, so the elastic constitutive call sees only the mechanical strain after both plastic and thermal contributions are removed:

Listing 25 tests/regression/solid_mechanics/viscoplasticity/free_sintering/model.i¶

# neml2
# Native port of tests/regression/solid_mechanics/viscoplasticity/free_sintering/model.i.
# Free sintering at zero applied strain with thermal eigenstrain + Olevsky-Skorohod
# sintering stress driving GTN poroplasticity. Batch (100, 10): 100 time steps,
# 10 surface-tension values (gamma in [0, 150]).
[Tensors]
  # end_time = LogspaceScalar(3, 3, 10) -> 10 copies of 10^3 = 1000.
  [end_time]
    type = Python
    expr = 'Scalar(torch.logspace(3.0, 3.0, 10, dtype=torch.float64))'
  []
  # times = LinspaceScalar(0, end_time, 100) -> shape (100, 10).
  [times]
    type = Python
    expr = 'Scalar(end_time.data.unsqueeze(0) * torch.linspace(0.0, 1.0, 100, dtype=torch.float64).unsqueeze(-1))'
  []
  # start_temperature = LinspaceScalar(300, 300, 10) -> 10 copies of 300.
  [start_temperature]
    type = Python
    expr = 'Scalar.linspace(300.0, 300.0, 10)'
  []
  # end_temperature = LinspaceScalar(1800, 1800, 10) -> 10 copies of 1800.
  [end_temperature]
    type = Python
    expr = 'Scalar.linspace(1800.0, 1800.0, 10)'
  []
  # temperatures = LinspaceScalar(start_temperature, end_temperature, 100) -> (100, 10).
  [temperatures]
    type = Python
    expr = 'Scalar(start_temperature.data.unsqueeze(0) + (end_temperature.data - start_temperature.data).unsqueeze(0) * torch.linspace(0.0, 1.0, 100, dtype=torch.float64).unsqueeze(-1))'
  []
  # max_strain = FillSR2(exx=0, eyy=0, ezz=0) batched (10,) -> all zeros.
  # Diagonal 3-arg fill places exx/eyy/ezz on slots 0/1/2 with no Mandel scaling.
  [max_strain]
    type = Python
    expr = 'SR2(torch.zeros((10, 6), dtype=torch.float64))'
  []
  # strains = LinspaceSR2(0, max_strain, 100) -> (100, 10, 6), all zeros.
  [strains]
    type = Python
    expr = 'SR2(torch.zeros((100, 10, 6), dtype=torch.float64))'
  []
  # f0 = FullScalar(0.36, (10,))
  [f0]
    type = Python
    expr = 'Scalar.full(10, fill_value=0.36)'
  []
  # gamma = LinspaceScalar(0, 150, 10) -> shape (10,)
  [gamma]
    type = Python
    expr = 'Scalar.linspace(0.0, 150.0, 10)'
  []
[]

[Drivers]
  [driver]
    type = TransientDriver
    model = 'model'
    prescribed_time = 'times'
    force_SR2_names = 'E'
    force_SR2_values = 'strains'
    force_Scalar_names = 'temperature'
    force_Scalar_values = 'temperatures'
    ic_Scalar_names = 'void_fraction'
    ic_Scalar_values = 'f0'
    save_as = 'result.pt'
  []
  [regression]
    type = TransientRegression
    driver = 'driver'
    reference = 'gold/result.pt'
  []
[]

[Models]
  [isoharden]
    type = VoceIsotropicHardening
    saturated_hardening = 5
    saturation_rate = 1.2
  []
  [sintering_stress]
    type = OlevskySinteringStress
    surface_tension = 'gamma'
    particle_radius = 3e-4
  []
  [eigenstrain]
    type = ThermalEigenstrain
    reference_temperature = 300
    CTE = 1e-6
  []
  [elastic_strain]
    type = SR2LinearCombination
    from = 'E plastic_strain eigenstrain'
    to = 'elastic_strain'
    weights = '1 -1 -1'
  []
  [elasticity]
    type = LinearIsotropicElasticity
    coefficients = '3e4 0.3'
    coefficient_types = 'YOUNGS_MODULUS POISSONS_RATIO'
    strain = 'elastic_strain'
  []
  [mandel_stress]
    type = IsotropicMandelStress
    cauchy_stress = 'stress'
  []
  [j2]
    type = SR2Invariant
    invariant_type = 'VONMISES'
    tensor = 'mandel_stress'
    invariant = 'flow_invariant'
  []
  [sh]
    type = SR2Invariant
    invariant_type = 'I1'
    tensor = 'mandel_stress'
    invariant = 'hydrostatic_stress'
  []
  [sp]
    type = ScalarLinearCombination
    from = 'hydrostatic_stress sintering_stress'
    to = 'poro_invariant'
    weights = '1 -1'
  []
  [q1]
    type = ArrheniusParameter
    temperature = 'temperature'
    reference_value = 8000
    activation_energy = 5e4
    ideal_gas_constant = 8.314
  []
  [yield_surface]
    type = GTNYieldFunction
    yield_stress = 60.0
    q1 = 'q1'
    q2 = 0.01
    q3 = 1.57
    isotropic_hardening = 'isotropic_hardening'
  []
  [flow]
    type = ComposedModel
    models = 'j2 sh sp yield_surface'
  []
  [flow_rate]
    type = PerzynaPlasticFlowRate
    reference_stress = 500
    exponent = 2
  []
  [normality]
    type = Normality
    model = 'flow'
    function = 'yield_function'
    from = 'mandel_stress isotropic_hardening'
    to = 'flow_direction isotropic_hardening_direction'
  []
  [Eprate]
    type = AssociativePlasticFlow
  []
  [eprate]
    type = AssociativeIsotropicPlasticHardening
  []
  [voidrate]
    type = GursonCavitation
  []
  [integrate_Ep]
    type = SR2BackwardEulerTimeIntegration
    variable = 'plastic_strain'
  []
  [integrate_ep]
    type = ScalarBackwardEulerTimeIntegration
    variable = 'equivalent_plastic_strain'
  []
  [integrate_void]
    type = ScalarBackwardEulerTimeIntegration
    variable = 'void_fraction'
  []
  [surface]
    type = ComposedModel
    models = 'isoharden sintering_stress elastic_strain elasticity mandel_stress flow flow_rate normality Eprate eprate voidrate integrate_Ep integrate_ep integrate_void'
  []
[]

[EquationSystems]
  [eq_sys]
    type = NonlinearSystem
    model = 'surface'
    unknowns = 'plastic_strain equivalent_plastic_strain void_fraction'
    residuals = 'plastic_strain_residual equivalent_plastic_strain_residual void_fraction_residual'
  []
[]

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

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

Explanation¶

The kinematic piece of the model is concentrated in two [Models] blocks.

[eigenstrain] declares a ThermalEigenstrain with reference_temperature = 300 and CTE = 1e-6. At every batch entry its temperature input is sourced from the temperature force prescribed by the driver, and it produces eigenstrain : SR2 — the cumulative thermal strain \(\alpha(T - T_0)\,\boldsymbol{I}\) relative to the stress-free 300 K state.

[elastic_strain] is an SR2LinearCombination that implements the additive split

\[ \boldsymbol{\varepsilon}^{\text{el}} = 1\cdot\boldsymbol{\varepsilon} - 1\cdot\boldsymbol{\varepsilon}^{p} - 1\cdot\boldsymbol{\varepsilon}^{*}_T, \]

with from = 'E plastic_strain eigenstrain' and weights = '1 -1 -1'. The output elastic_strain is then the strain that the LinearIsotropicElasticity block consumes — so by the time stress is evaluated the thermal contribution has already been peeled off. This is the typical wiring pattern for eigenstrains in NEML2: declare the eigenstrain model, then subtract its eigenstrain output from the total strain inside an SR2LinearCombination whose result is consumed by elasticity.

To extend this composition to additional eigenstrain sources — say adding a phase-transformation contribution via PhaseTransformationEigenstrain — declare the second eigenstrain block, give its output a distinct variable name, and add it to the from / weights lists of the same SR2LinearCombination. No other part of the model needs to change because every eigenstrain leaf in the catalog publishes the same eigenstrain : SR2 surface.

For finite-strain workflows the analogous pattern uses ThermalDeformationJacobian or SwellingAndPhaseChangeDeformationJacobian, multiplied together via a ScalarMultiplication block, then handed to VolumeAdjustDeformationGradient to produce the mechanical deformation gradient that the elastic model consumes. On the small-strain side, VolumeChangeEigenstrain plays the corresponding role when a volume-change eigenstrain is needed.

Note

All eigenstrain leaves are cumulative — they return the total non-mechanical strain relative to the reference state, not an increment. This is why subtraction at the elastic-strain step yields the correct mechanical strain regardless of how the temperature or phase history evolves through time.