The radial return mapping method, introduced by Simo and Taylor (1985), uses a von Mises yield surface to determine the increment of plastic strain necessary to return the stress state to the yield surface after a trial stress increment takes the computed stress state across the yield surface. Because the von Mises yield surface in the deviatoric stress space has the shape of a circle, the plastic correction stress is always directed towards the center of the yield surface circle.

In addition to the Sumo and Hughes (1998) textbook, Introduction to Computational Plasticity by Dunnu and Petrinic (2004) (ISBN: 0-19-856826-6) is an excellent reference for users working with the Recompute Radial Return materials; several of the isotropic plasticity and creep effective plastic strain increment algorithms are taken from Dunne and Petrinic.

The Recompute Radial Return Approach in a Nutshell

The recompute materials are not called by MOOSE directly but instead only by other materials using the computeProperties method. For the RecomputeRadialReturn materials, this calling material is ComputeReturnMappingStress(Doxygen). Separating the call to the recompute materials from MOOSE allows us to iteratively call the recompute materials as is required to achieve convergence. (For more details on the recompute materials within MOOSE, see this page on Recompute Materials).

Why use the Radial Recompute Materials? For isotropic materials the radial return approach offers distinct advantages:

  • Faster simulation run times: The isotropic material iteration algorithm uses single variable Reals to compute and converge the inelastic strain instead of inverting the full Rank-4 elasticity tensor required in more complicated anisotropic algorithms.
  • Easy to understand: The return mapping algorithm implemented in RecomputeRadialReturn (Doxygen) is the classical radial return method based on the von Mises yield criterion. To write a new radial recompute material, simply inherit from RecomputeRadialReturn and overwrite the virtual methods to calculate the effective plastic strain increment.
  • Applicable to a variety of models: The radial return method provides the flexibility to include creep, plasticity, and damage within a single simulation. The ComputeReturnMappingStress class calls each individual creep and plasticity model to iterate separately over the inelastic strain increment before checking for the convergence of the combined total radial return stress increment required to return the stress state to the yield surface.

Radial Return Algorithm Overview

A trial stress is shown outside of the devatoric yield surface and the radial return stress which is normal to the yield surface.

In this numerical approach, a trial stress is calculated at the start of each simulation time increment; the trial stress calculation assumed all of the new strain increment is elastic strain: $$$$\sigma_{trial} = C_{ijkl} \left( \Delta \epsilon_{assumed-elastic} + \epsilon_{elastic}^{old} \right)$$$$ The algorithms checks to see if the trial stress state is outside of the yield surface, as shown in the figure to the right. If the stress state is outside of the yield surface, the algorithm recomputes the scalar effective inelastic strain required to return the stress state to the yield surface. This approach is given the name Radial Return because the yield surface used is the von Mises yield surface: in the devitoric stress space , this yield surface has the shape of a circle, and the scalar inelastic strain is assumed to always be directed at the circle center.

'Recompute' Iterations on the Effective Plastic Strain Increment

The recompute radial return materials each individually calculate, using the Newton Method, the amount of effective inelastic strain required to return the stress state to the yield surface. $$$$\Delta p^{(t+1)} = \Delta p^t + d \Delta p$$$$ where the change in the iterative effective inelastic strain is defined as the yield surface over the derivative of the yield surface with respect to the inelastic strain increment. In the case of isotropic linear hardening plasticity, with the hardening function $$$r = hp$$$, the effective plastic strain increment has the form: $$$$d \Delta p = \frac{\sigma^{trial}_{effective} - 3 G \Delta p - r - \sigma_{yield}}{3G + h}$$$$ where G is the isotropic shear modulus, and $$$\sigma^{trial}_{effective}$$$ is the scalar von Mises trial stress.

Once convergence has been reached on the scalar inelastic strain increment, the full inelastic strain tensor is calculated. $$$$\Delta \epsilon^{inelastic}_{ij} = \frac{3}{2} \Delta p^{(t+1)} \frac{dev(\sigma^{trial}_{ij})}{\sigma^{trial}_{effective}}$$$$

The elastic strain is calculated by subtracting the return mapping inelastic strain increment tensor from the mechanical strain tensor. Mechanical strain is considered as the sum of the elastic and inelastic (plastic, creep, ect) strains. $$$$\epsilon_{total} = \epsilon_{mechanical} + \epsilon_{eigenstrain} = \left( \epsilon_{elastic} + \epsilon_{inelastic} \right) + \epsilon_{eigenstrain} = \epsilon_{elastic} + \left( \epsilon_{plastic} + \epsilon_{creep} + \epsilon_{damage} \right) + \epsilon_{eigenstrain}$$$$

The final inelastic strain is returned from the recompute radial material, and ComputeReturnMappingStress computes the stress, with a return mapping stress increment following elasticity theory for finite strains. The final stress is calculated from the elastic strain increment. $$$$\sigma^{new}_{ij} = C_{ijkl} \left( \Delta \epsilon^{elastic}_{kl} + \epsilon^{old-elastic}_{kl} \right)$$$$

When more than one radial recompute material is included in the simulation, as in Combined Power Law Creep and Linear Strain Hardening, ComputeReturnMappingStress will iterate over the change in the calculated stress until the return stress has reached a stable value.

Users can print out any of these strains and stresses using the RankTwoAux as described in Visualizing Tensors.

Recompute Material Methods

New radial return models must inherit from RecomputeRadialReturn(Doxygen) and must overwrite the six virtual methods.

  • initQpStatefulProperties: Set the initial values for all new material properties that are not initialized by an input parameter; generally the material properties initialized in this method are all set to zero.
  • computeStressInitialize: Calculate the initial trial stress state, the yield surface value, and any hardening or softening parameters at the start of the simulation time increment.
  • computeResidual: In each iteration over the inelastic strain increment, calculate the value of the effective scalar trial stress subtracted by the yield surface function.
  • computeDerivative: In each iteration over the inelastic strain increment, calculate the derivative of the yield surface function with respect to the inelastic strain increment.

  • iterationFinalize: Store the value of the inelastic strain increment at the end of each iteration.

  • computeStressFinalize: Update the stress after convergence on the inelastic strain increment has been reached.

Additionally, new radial return methods must also overwrite a single method from the MOOSE Material class.

  • resetQpProperties: Set the material property used in the iteration, usually $$$\Delta p$$$, to zero at the start the iteration. This method is necessary to avoid incorrect material property values.

More details on how to write the equivalent yield surface equation for a creep model are given in Dunne and Petrinic.