FiniteStrainPlasticMaterial implements rate-independent associative J2 plasticity with isotropic hardening in the finite-strain framework. More...

#include <FiniteStrainPlasticMaterial.h>

Inheritance diagram for FiniteStrainPlasticMaterial:

[legend]

Public Types
typedef DerivativeMaterialPropertyNameInterface::SymbolName	SymbolName

Public Member Functions
	FiniteStrainPlasticMaterial (const InputParameters &parameters)

const GenericMaterialProperty< U, is_ad > &	getDefaultMaterialProperty (const std::string &name)

const GenericMaterialProperty< U, is_ad > &	getDefaultMaterialPropertyByName (const std::string &name)

void	validateDerivativeMaterialPropertyBase (const std::string &base)

const MaterialPropertyName	derivativePropertyName (const MaterialPropertyName &base, const std::vector< SymbolName > &c) const

const MaterialPropertyName	derivativePropertyNameFirst (const MaterialPropertyName &base, const SymbolName &c1) const

const MaterialPropertyName	derivativePropertyNameSecond (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2) const

const MaterialPropertyName	derivativePropertyNameThird (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2, const SymbolName &c3) const

GenericMaterialProperty< U, is_ad > &	declarePropertyDerivative (const std::string &base, const std::vector< VariableName > &c)

GenericMaterialProperty< U, is_ad > &	declarePropertyDerivative (const std::string &base, const std::vector< SymbolName > &c)

GenericMaterialProperty< U, is_ad > &	declarePropertyDerivative (const std::string &base, const SymbolName &c1, const SymbolName &c2="", const SymbolName &c3="")

GenericMaterialProperty< U, is_ad > &	declarePropertyDerivative (const std::string &base, const std::vector< VariableName > &c)

GenericMaterialProperty< U, is_ad > &	declarePropertyDerivative (const std::string &base, const std::vector< SymbolName > &c)

GenericMaterialProperty< U, is_ad > &	declarePropertyDerivative (const std::string &base, const SymbolName &c1, const SymbolName &c2="", const SymbolName &c3="")

const GenericMaterialProperty< U, is_ad > &	getMaterialPropertyDerivative (const std::string &base, const std::vector< VariableName > &c)

const GenericMaterialProperty< U, is_ad > &	getMaterialPropertyDerivative (const std::string &base, const std::vector< SymbolName > &c)

const GenericMaterialProperty< U, is_ad > &	getMaterialPropertyDerivative (const std::string &base, const SymbolName &c1, const SymbolName &c2="", const SymbolName &c3="")

const GenericMaterialProperty< U, is_ad > &	getMaterialPropertyDerivative (const std::string &base, const SymbolName &c1, unsigned int v2, unsigned int v3=libMesh::invalid_uint)

const GenericMaterialProperty< U, is_ad > &	getMaterialPropertyDerivative (const std::string &base, unsigned int v1, unsigned int v2=libMesh::invalid_uint, unsigned int v3=libMesh::invalid_uint)

const GenericMaterialProperty< U, is_ad > &	getMaterialPropertyDerivative (const std::string &base, const std::vector< VariableName > &c)

const GenericMaterialProperty< U, is_ad > &	getMaterialPropertyDerivative (const std::string &base, const std::vector< SymbolName > &c)

const GenericMaterialProperty< U, is_ad > &	getMaterialPropertyDerivative (const std::string &base, const SymbolName &c1, const SymbolName &c2="", const SymbolName &c3="")

const GenericMaterialProperty< U, is_ad > &	getMaterialPropertyDerivative (const std::string &base, const SymbolName &c1, unsigned int v2, unsigned int v3=libMesh::invalid_uint)

const GenericMaterialProperty< U, is_ad > &	getMaterialPropertyDerivative (const std::string &base, unsigned int v1, unsigned int v2=libMesh::invalid_uint, unsigned int v3=libMesh::invalid_uint)

const GenericMaterialProperty< U, is_ad > &	getMaterialPropertyDerivativeByName (const MaterialPropertyName &base, const std::vector< VariableName > &c)

const GenericMaterialProperty< U, is_ad > &	getMaterialPropertyDerivativeByName (const MaterialPropertyName &base, const std::vector< SymbolName > &c)

const GenericMaterialProperty< U, is_ad > &	getMaterialPropertyDerivativeByName (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2="", const SymbolName &c3="")

const GenericMaterialProperty< U, is_ad > &	getMaterialPropertyDerivativeByName (const MaterialPropertyName &base, const std::vector< VariableName > &c)

const GenericMaterialProperty< U, is_ad > &	getMaterialPropertyDerivativeByName (const MaterialPropertyName &base, const std::vector< SymbolName > &c)

const GenericMaterialProperty< U, is_ad > &	getMaterialPropertyDerivativeByName (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2="", const SymbolName &c3="")

void	validateCoupling (const MaterialPropertyName &base, const std::vector< VariableName > &c, bool validate_aux=true)

void	validateCoupling (const MaterialPropertyName &base, const VariableName &c1="", const VariableName &c2="", const VariableName &c3="")

void	validateCoupling (const MaterialPropertyName &base, const std::vector< VariableName > &c, bool validate_aux=true)

void	validateCoupling (const MaterialPropertyName &base, const VariableName &c1="", const VariableName &c2="", const VariableName &c3="")

void	validateNonlinearCoupling (const MaterialPropertyName &base, const VariableName &c1="", const VariableName &c2="", const VariableName &c3="")

void	validateNonlinearCoupling (const MaterialPropertyName &base, const VariableName &c1="", const VariableName &c2="", const VariableName &c3="")

const MaterialPropertyName	propertyName (const MaterialPropertyName &base, const std::vector< SymbolName > &c) const

const MaterialPropertyName	propertyName (const MaterialPropertyName &base, const std::vector< SymbolName > &c) const

const MaterialPropertyName	propertyNameFirst (const MaterialPropertyName &base, const SymbolName &c1) const

const MaterialPropertyName	propertyNameFirst (const MaterialPropertyName &base, const SymbolName &c1) const

const MaterialPropertyName	propertyNameSecond (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2) const

const MaterialPropertyName	propertyNameSecond (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2) const

const MaterialPropertyName	propertyNameThird (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2, const SymbolName &c3) const

const MaterialPropertyName	propertyNameThird (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2, const SymbolName &c3) const

Static Public Member Functions
static InputParameters	validParams ()

Protected Member Functions
virtual void	computeQpStress ()
	Compute the stress and store it in the _stress material property for the current quadrature point.

virtual void	initQpStatefulProperties ()

virtual void	returnMap (const RankTwoTensor &sig_old, const Real eqvpstrain_old, const RankTwoTensor &plastic_strain_old, const RankTwoTensor &delta_d, const RankFourTensor &E_ijkl, RankTwoTensor &sig, Real &eqvpstrain, RankTwoTensor &plastic_strain)
	Implements the return map.

virtual Real	yieldFunction (const RankTwoTensor &stress, const Real yield_stress)
	Calculates the yield function.

virtual RankTwoTensor	dyieldFunction_dstress (const RankTwoTensor &stress)
	Derivative of yieldFunction with respect to the stress.

virtual Real	dyieldFunction_dinternal (const Real equivalent_plastic_strain)
	Derivative of yieldFunction with respect to the equivalent plastic strain.

virtual RankTwoTensor	flowPotential (const RankTwoTensor &stress)
	Flow potential, which in this case is just dyieldFunction_dstress because we are doing associative flow, and hence does not depend on the internal hardening parameter equivalent_plastic_strain.

virtual Real	internalPotential ()
	The internal potential.

Real	getSigEqv (const RankTwoTensor &stress)
	Equivalent stress.

virtual void	getJac (const RankTwoTensor &sig, const RankFourTensor &E_ijkl, Real flow_incr, RankFourTensor &dresid_dsig)
	Evaluates the derivative d(resid_ij)/d(sig_kl), where resid_ij = flow_incr*flowPotential_ij - (E^{-1}(trial_stress - sig))_ij.

Real	getYieldStress (const Real equivalent_plastic_strain)
	yield stress as a function of equivalent plastic strain.

Real	getdYieldStressdPlasticStrain (const Real equivalent_plastic_strain)
	d(yieldstress)/d(equivalent plastic strain)

	usingTensorIndices (i_, j_, k_, l_)
	make i,j,k,l available as tensor indices for mixedProduct

virtual void	computeQpProperties () override

Protected Attributes
std::vector< Real >	_yield_stress_vector

MaterialProperty< RankTwoTensor > &	_plastic_strain

const MaterialProperty< RankTwoTensor > &	_plastic_strain_old

MaterialProperty< Real > &	_eqv_plastic_strain

const MaterialProperty< Real > &	_eqv_plastic_strain_old

const MaterialProperty< RankTwoTensor > &	_stress_old
	The old stress tensor.

const MaterialProperty< RankTwoTensor > &	_strain_increment

const MaterialProperty< RankTwoTensor > &	_rotation_increment

const std::string	_elasticity_tensor_name
	Name of the elasticity tensor material property.

const MaterialProperty< RankFourTensor > &	_elasticity_tensor
	Elasticity tensor material property.

Real	_rtol

Real	_ftol

Real	_eptol

RankFourTensor	_deltaOuter

RankFourTensor	_deltaMixed

const std::string	_base_name
	Base name prepended to all material property names to allow for multi-material systems.

const MaterialProperty< RankTwoTensor > &	_mechanical_strain
	Mechanical strain material property.

MaterialProperty< RankTwoTensor > &	_stress
	Stress material property.

MaterialProperty< RankTwoTensor > &	_elastic_strain
	Elastic strain material property.

const MaterialProperty< RankTwoTensor > &	_extra_stress
	Extra stress tensor.

std::vector< const Function * >	_initial_stress_fcn
	initial stress components

MaterialProperty< RankFourTensor > &	_Jacobian_mult
	derivative of stress w.r.t. strain (_dstress_dstrain)

Private Member Functions
bool	haveMaterialProperty (const std::string &prop_name)

std::vector< VariableName >	buildVariableVector (const VariableName &c1, const VariableName &c2, const VariableName &c3)

void	validateCouplingHelper (const MaterialPropertyName &base, const std::vector< VariableName > &c, const System &system, std::vector< VariableName > &missing)

bool	isNotObjectVariable (const VariableName &name)

Private Attributes
FEProblemBase &	_dmi_fe_problem

Detailed Description

FiniteStrainPlasticMaterial implements rate-independent associative J2 plasticity with isotropic hardening in the finite-strain framework.

Yield function = sqrt(3*s_ij*s_ij/2) - K(equivalent plastic strain) where s_ij = stress_ij - delta_ij*trace(stress)/3 is the deviatoric stress and K is the yield stress, specified as a piecewise-linear function by the user. Integration is performed in an incremental manner using Newton Raphson.

Definition at line 24 of file FiniteStrainPlasticMaterial.h.

Constructor & Destructor Documentation

◆ FiniteStrainPlasticMaterial()

FiniteStrainPlasticMaterial::FiniteStrainPlasticMaterial ( const InputParameters & parameters )

Definition at line 32 of file FiniteStrainPlasticMaterial.C.

  : ComputeStressBase(parameters),
    _yield_stress_vector(getParam<std::vector<Real>>("yield_stress")), // Read from input file
    _plastic_strain(declareProperty<RankTwoTensor>(_base_name + "plastic_strain")),
    _plastic_strain_old(getMaterialPropertyOld<RankTwoTensor>(_base_name + "plastic_strain")),
    _eqv_plastic_strain(declareProperty<Real>(_base_name + "eqv_plastic_strain")),
    _eqv_plastic_strain_old(getMaterialPropertyOld<Real>(_base_name + "eqv_plastic_strain")),
    _stress_old(getMaterialPropertyOld<RankTwoTensor>(_base_name + "stress")),
    _strain_increment(getMaterialProperty<RankTwoTensor>(_base_name + "strain_increment")),
    _rotation_increment(getMaterialProperty<RankTwoTensor>(_base_name + "rotation_increment")),
    _elasticity_tensor_name(_base_name + "elasticity_tensor"),
    _elasticity_tensor(getMaterialPropertyByName<RankFourTensor>(_elasticity_tensor_name)),
    _rtol(getParam<Real>("rtol")),
    _ftol(getParam<Real>("ftol")),
    _eptol(getParam<Real>("eptol")),
    _deltaOuter(RankTwoTensor::Identity().times<i_, j_, k_, l_>(RankTwoTensor::Identity())),
    _deltaMixed(RankTwoTensor::Identity().times<i_, k_, j_, l_>(RankTwoTensor::Identity()))
{
}

Member Function Documentation

◆ computeQpProperties()

void ComputeGeneralStressBase::computeQpProperties ( )

overrideprotectedvirtualinherited

Definition at line 44 of file ComputeGeneralStressBase.C.

{
  computeQpStress();
 
  // Add in extra stress
  _stress[_qp] += _extra_stress[_qp];
}

◆ computeQpStress()

void FiniteStrainPlasticMaterial::computeQpStress ( )

protectedvirtual

Compute the stress and store it in the _stress material property for the current quadrature point.

Implements ComputeGeneralStressBase.

Definition at line 61 of file FiniteStrainPlasticMaterial.C.

{
 
  // perform the return-mapping algorithm
  returnMap(_stress_old[_qp],
            _eqv_plastic_strain_old[_qp],
            _plastic_strain_old[_qp],
            _strain_increment[_qp],
            _elasticity_tensor[_qp],
            _stress[_qp],
            _eqv_plastic_strain[_qp],
            _plastic_strain[_qp]);
 
  // Rotate the stress tensor to the current configuration
  _stress[_qp] = _rotation_increment[_qp] * _stress[_qp] * _rotation_increment[_qp].transpose();
 
  // Rotate plastic strain tensor to the current configuration
  _plastic_strain[_qp] =
      _rotation_increment[_qp] * _plastic_strain[_qp] * _rotation_increment[_qp].transpose();
 
  // Calculate the elastic strain_increment
  _elastic_strain[_qp] = _mechanical_strain[_qp] - _plastic_strain[_qp];
 
  _Jacobian_mult[_qp] = _elasticity_tensor[_qp];
}

◆ dyieldFunction_dinternal()

Real FiniteStrainPlasticMaterial::dyieldFunction_dinternal ( const Real equivalent_plastic_strain )

protectedvirtual

Derivative of yieldFunction with respect to the equivalent plastic strain.

Definition at line 324 of file FiniteStrainPlasticMaterial.C.

{
  return -getdYieldStressdPlasticStrain(equivalent_plastic_strain);
}

Referenced by returnMap().

◆ dyieldFunction_dstress()

RankTwoTensor FiniteStrainPlasticMaterial::dyieldFunction_dstress ( const RankTwoTensor & stress )

protectedvirtual

Derivative of yieldFunction with respect to the stress.

Definition at line 316 of file FiniteStrainPlasticMaterial.C.

{
  RankTwoTensor deriv = sig.dsecondInvariant();
  deriv *= std::sqrt(3.0 / sig.secondInvariant()) / 2.0;
  return deriv;
}

Referenced by flowPotential(), getJac(), and returnMap().

◆ flowPotential()

RankTwoTensor FiniteStrainPlasticMaterial::flowPotential ( const RankTwoTensor & stress )

protectedvirtual

Flow potential, which in this case is just dyieldFunction_dstress because we are doing associative flow, and hence does not depend on the internal hardening parameter equivalent_plastic_strain.

Definition at line 330 of file FiniteStrainPlasticMaterial.C.

{
  return dyieldFunction_dstress(sig); // this plasticity model assumes associative flow
}

Referenced by getJac(), and returnMap().

◆ getdYieldStressdPlasticStrain()

Real FiniteStrainPlasticMaterial::getdYieldStressdPlasticStrain ( const Real equivalent_plastic_strain )

protected

d(yieldstress)/d(equivalent plastic strain)

Definition at line 414 of file FiniteStrainPlasticMaterial.C.

{
  unsigned nsize;
 
  nsize = _yield_stress_vector.size();
 
  if (_yield_stress_vector[0] > 0.0 || nsize % 2 > 0) // Error check for input inconsitency
    mooseError("Error in yield stress input: Should be a vector with eqv plastic strain and yield "
               "stress pair values.\n");
 
  unsigned int ind = 0;
 
  while (ind < nsize)
  {
    if (ind + 2 < nsize)
    {
      if (eqpe >= _yield_stress_vector[ind] && eqpe < _yield_stress_vector[ind + 2])
        return (_yield_stress_vector[ind + 3] - _yield_stress_vector[ind + 1]) /
               (_yield_stress_vector[ind + 2] - _yield_stress_vector[ind]);
    }
    else
      return 0.0;
 
    ind += 2;
  }
 
  return 0.0;
}

Referenced by dyieldFunction_dinternal().

◆ getJac()

void FiniteStrainPlasticMaterial::getJac	(	const RankTwoTensor &	sig,
		const RankFourTensor &	E_ijkl,
		Real	flow_incr,
		RankFourTensor &	dresid_dsig
	)

protectedvirtual

Evaluates the derivative d(resid_ij)/d(sig_kl), where resid_ij = flow_incr*flowPotential_ij - (E^{-1}(trial_stress - sig))_ij.

Parameters

sig	stress
E_ijkl	elasticity tensor (sig = E*(strain - plastic_strain))
flow_incr	consistency parameter
dresid_dsig	the required derivative (this is an output variable)

Definition at line 349 of file FiniteStrainPlasticMaterial.C.

{
  RankTwoTensor sig_dev, df_dsig, flow_dirn;
  RankTwoTensor dfi_dft, dfi_dsig;
  RankFourTensor dft_dsig, dfd_dft, dfd_dsig;
  Real sig_eqv;
  Real f1, f2, f3;
  RankFourTensor temp;
 
  sig_dev = sig.deviatoric();
  sig_eqv = getSigEqv(sig);
  df_dsig = dyieldFunction_dstress(sig);
  flow_dirn = flowPotential(sig);
 
  f1 = 3.0 / (2.0 * sig_eqv);
  f2 = f1 / 3.0;
  f3 = 9.0 / (4.0 * Utility::pow<3>(sig_eqv));
 
  dft_dsig = f1 * _deltaMixed - f2 * _deltaOuter - f3 * sig_dev.outerProduct(sig_dev);
 
  dfd_dsig = dft_dsig;
  dresid_dsig = E_ijkl.invSymm() + dfd_dsig * flow_incr;
}

Referenced by returnMap().

◆ getSigEqv()

Real FiniteStrainPlasticMaterial::getSigEqv ( const RankTwoTensor & stress )

protected

Equivalent stress.

Definition at line 342 of file FiniteStrainPlasticMaterial.C.

{
  return std::sqrt(3 * stress.secondInvariant());
}

Referenced by getJac(), and yieldFunction().

◆ getYieldStress()

Real FiniteStrainPlasticMaterial::getYieldStress ( const Real equivalent_plastic_strain )

protected

yield stress as a function of equivalent plastic strain.

This is a piecewise linear function entered by the user in the yield_stress vector

Definition at line 378 of file FiniteStrainPlasticMaterial.C.

{
  unsigned nsize;
 
  nsize = _yield_stress_vector.size();
 
  if (_yield_stress_vector[0] > 0.0 || nsize % 2 > 0) // Error check for input inconsitency
    mooseError("Error in yield stress input: Should be a vector with eqv plastic strain and yield "
               "stress pair values.\n");
 
  unsigned int ind = 0;
  Real tol = 1e-8;
 
  while (ind < nsize)
  {
    if (std::abs(eqpe - _yield_stress_vector[ind]) < tol)
      return _yield_stress_vector[ind + 1];
 
    if (ind + 2 < nsize)
    {
      if (eqpe > _yield_stress_vector[ind] && eqpe < _yield_stress_vector[ind + 2])
        return _yield_stress_vector[ind + 1] +
               (eqpe - _yield_stress_vector[ind]) /
                   (_yield_stress_vector[ind + 2] - _yield_stress_vector[ind]) *
                   (_yield_stress_vector[ind + 3] - _yield_stress_vector[ind + 1]);
    }
    else
      return _yield_stress_vector[nsize - 1];
 
    ind += 2;
  }
 
  return 0.0;
}

Referenced by returnMap().

◆ initQpStatefulProperties()

void FiniteStrainPlasticMaterial::initQpStatefulProperties ( )

protectedvirtual

Reimplemented from ComputeGeneralStressBase.

Definition at line 53 of file FiniteStrainPlasticMaterial.C.

{
  ComputeStressBase::initQpStatefulProperties();
  _plastic_strain[_qp].zero();
  _eqv_plastic_strain[_qp] = 0.0;
}

◆ internalPotential()

Real FiniteStrainPlasticMaterial::internalPotential ( )

protectedvirtual

The internal potential.

For associative J2 plasticity this is just -1

Definition at line 336 of file FiniteStrainPlasticMaterial.C.

{
  return -1;
}

Referenced by returnMap().

◆ returnMap()

void FiniteStrainPlasticMaterial::returnMap	(	const RankTwoTensor &	sig_old,
		const Real	eqvpstrain_old,
		const RankTwoTensor &	plastic_strain_old,
		const RankTwoTensor &	delta_d,
		const RankFourTensor &	E_ijkl,
		RankTwoTensor &	sig,
		Real &	eqvpstrain,
		RankTwoTensor &	plastic_strain
	)

protectedvirtual

Implements the return map.

Implements the return-map algorithm via a Newton-Raphson process.

Parameters

sig_old	The stress at the previous "time" step
eqvpstrain_old	The equivalent plastic strain at the previous "time" step
plastic_strain_old	The value of plastic strain at the previous "time" step
delta_d	The total strain increment for this "time" step
E_ijkl	The elasticity tensor. If no plasticity then sig_new = sig_old + E_ijkl*delta_d
sig	The stress after returning to the yield surface (this is an output variable)
eqvpstrain	The equivalent plastic strain after returning to the yield surface (this is an output variable)
plastic_strain	The value of plastic strain after returning to the yield surface (this is an output variable) Note that this algorithm doesn't do any rotations. In order to find the final stress and plastic_strain, sig and plastic_strain must be rotated using _rotation_increment.

This idea is fully explained in Simo and Hughes "Computational Inelasticity" Springer 1997, for instance, as well as many other books on plasticity. The basic idea is as follows. Given: sig_old - the stress at the start of the "time step" plastic_strain_old - the plastic strain at the start of the "time step" eqvpstrain_old - equivalent plastic strain at the start of the "time step" (In general we would be given some number of internal parameters that the yield function depends upon.) delta_d - the prescribed strain increment for this "time step" we want to determine the following parameters at the end of this "time step": sig - the stress plastic_strain - the plastic strain eqvpstrain - the equivalent plastic strain (again, in general, we would have an arbitrary number of internal parameters).

To determine these parameters, introduce the "yield function", f the "consistency parameter", flow_incr the "flow potential", flow_dirn_ij the "internal potential", internalPotential (in general there are as many internalPotential functions as there are internal parameters). All three of f, flow_dirn_ij, and internalPotential, are functions of sig and eqvpstrain. To find sig, plastic_strain and eqvpstrain, we need to solve the following resid_ij = 0 f = 0 rep = 0 This is done by using Newton-Raphson. There are 8 equations here: six from resid_ij=0 (more generally there are nine but in this case resid is symmetric); one from f=0; one from rep=0 (more generally, for N internal parameters there are N of these equations).

resid_ij = flow_incr*flow_dirn_ij - (plastic_strain - plastic_strain_old)_ij = flow_incr*flow_dirn_ij - (E^{-1}(trial_stress - sig))_ij Here trial_stress = E*(strain - plastic_strain_old) sig = E*(strain - plastic_strain) Note: flow_dirn_ij is evaluated at sig and eqvpstrain (not the old values).

f is the yield function, evaluated at sig and eqvpstrain

rep = -flow_incr*internalPotential - (eqvpstrain - eqvpstrain_old) Here internalPotential are evaluated at sig and eqvpstrain (not the old values).

The Newton-Raphson procedure has sig, flow_incr, and eqvpstrain as its variables. Therefore we need the derivatives of resid_ij, f, and rep with respect to these parameters

In this associative J2 with isotropic hardening, things are a little more specialised. (1) f = sqrt(3*s_ij*s_ij/2) - K(eqvpstrain) (this is called "isotropic hardening") (2) associativity means that flow_dirn_ij = df/d(sig_ij) = s_ij*sqrt(3/2/(s_ij*s_ij)), and this means that flow_dirn_ij*flow_dirn_ij = 3/2, so when resid_ij=0, we get (plastic_strain_dot)_ij*(plastic_strain_dot)_ij = (3/2)*flow_incr^2, where plastic_strain_dot = plastic_strain - plastic_strain_old (3) The definition of equivalent plastic strain is through eqvpstrain_dot = sqrt(2*plastic_strain_dot_ij*plastic_strain_dot_ij/3), so using (2), we obtain eqvpstrain_dot = flow_incr, and this yields internalPotential = -1 in the "rep" equation.

The linear system is ( dr_dsig flow_dirn 0 )( ddsig ) ( - resid ) ( df_dsig 0 fq )( dflow_incr ) = ( - f ) ( 0 1 -1 )( deqvpstrain ) ( - rep ) The zeroes are: d(resid_ij)/d(eqvpstrain) = flow_dirn*d(df/d(sig_ij))/d(eqvpstrain) = 0 and df/d(flow_dirn) = 0 (this is always true, even for general hardening and non-associative) and d(rep)/d(sig_ij) = -flow_incr*d(internalPotential)/d(sig_ij) = 0

Because of the zeroes and ones, the linear system is not impossible to solve by hand. NOTE: andy believes there was originally a sign-error in the next line. The next line is unchanged, however andy's definition of fq is negative of the original definition of fq. andy can't see any difference in any tests!

Definition at line 149 of file FiniteStrainPlasticMaterial.C.

{
  // the yield function, must be non-positive
  // Newton-Raphson sets this to zero if trial stress enters inadmissible region
  Real f;
 
  // the consistency parameter, must be non-negative
  // change in plastic strain in this timestep = flow_incr*flow_potential
  Real flow_incr = 0.0;
 
  // direction of flow defined by the potential
  RankTwoTensor flow_dirn;
 
  // Newton-Raphson sets this zero
  // resid_ij = flow_incr*flow_dirn_ij - (plastic_strain - plastic_strain_old)
  RankTwoTensor resid;
 
  // Newton-Raphson sets this zero
  // rep = -flow_incr*internalPotential - (eqvpstrain - eqvpstrain_old)
  Real rep;
 
  // change in the stress (sig) in a Newton-Raphson iteration
  RankTwoTensor ddsig;
 
  // change in the consistency parameter in a Newton-Raphson iteration
  Real dflow_incr = 0.0;
 
  // change in equivalent plastic strain in one Newton-Raphson iteration
  Real deqvpstrain = 0.0;
 
  // convenience variable that holds the change in plastic strain incurred during the return
  // delta_dp = plastic_strain - plastic_strain_old
  // delta_dp = E^{-1}*(trial_stress - sig), where trial_stress = E*(strain - plastic_strain_old)
  RankTwoTensor delta_dp;
 
  // d(yieldFunction)/d(stress)
  RankTwoTensor df_dsig;
 
  // d(resid_ij)/d(sigma_kl)
  RankFourTensor dr_dsig;
 
  // dr_dsig_inv_ijkl*dr_dsig_klmn = 0.5*(de_ij de_jn + de_ij + de_jm), where de_ij = 1 if i=j, but
  // zero otherwise
  RankFourTensor dr_dsig_inv;
 
  // d(yieldFunction)/d(eqvpstrain)
  Real fq;
 
  // yield stress at the start of this "time step" (ie, evaluated with
  // eqvpstrain_old).  It is held fixed during the Newton-Raphson return,
  // even if eqvpstrain != eqvpstrain_old.
  Real yield_stress;
 
  // measures of whether the Newton-Raphson process has converged
  Real err1, err2, err3;
 
  // number of Newton-Raphson iterations performed
  unsigned int iter = 0;
 
  // maximum number of Newton-Raphson iterations allowed
  unsigned int maxiter = 100;
 
  // plastic loading occurs if yieldFunction > toly
  Real toly = 1.0e-8;
 
  // Assume this strain increment does not induce any plasticity
  // This is the elastic-predictor
  sig = sig_old + E_ijkl * delta_d; // the trial stress
  eqvpstrain = eqvpstrain_old;
  plastic_strain = plastic_strain_old;
 
  yield_stress = getYieldStress(eqvpstrain); // yield stress at this equivalent plastic strain
  if (yieldFunction(sig, yield_stress) > toly)
  {
    // the sig just calculated is inadmissable.  We must return to the yield surface.
    // This is done iteratively, using a Newton-Raphson process.
 
    delta_dp.zero();
 
    sig = sig_old + E_ijkl * delta_d; // this is the elastic predictor
 
    flow_dirn = flowPotential(sig);
 
    resid = flow_dirn * flow_incr - delta_dp; // Residual 1 - refer Hughes Simo
    f = yieldFunction(sig, yield_stress);
    rep = -eqvpstrain + eqvpstrain_old - flow_incr * internalPotential(); // Residual 3 rep=0
 
    err1 = resid.L2norm();
    err2 = std::abs(f);
    err3 = std::abs(rep);
 
    while ((err1 > _rtol || err2 > _ftol || err3 > _eptol) &&
           iter < maxiter) // Stress update iteration (hardness fixed)
    {
      iter++;
 
      df_dsig = dyieldFunction_dstress(sig);
      getJac(sig, E_ijkl, flow_incr, dr_dsig);   // gets dr_dsig = d(resid_ij)/d(sig_kl)
      fq = dyieldFunction_dinternal(eqvpstrain); // d(f)/d(eqvpstrain)
 
      dr_dsig_inv = dr_dsig.invSymm();
 
      dflow_incr = (f - df_dsig.doubleContraction(dr_dsig_inv * resid) + fq * rep) /
                   (df_dsig.doubleContraction(dr_dsig_inv * flow_dirn) - fq);
      ddsig =
          dr_dsig_inv *
          (-resid -
           flow_dirn * dflow_incr); // from solving the top row of linear system, given dflow_incr
      deqvpstrain =
          rep + dflow_incr; // from solving the bottom row of linear system, given dflow_incr
 
      // update the variables
      flow_incr += dflow_incr;
      delta_dp -= E_ijkl.invSymm() * ddsig;
      sig += ddsig;
      eqvpstrain += deqvpstrain;
 
      // evaluate the RHS equations ready for next Newton-Raphson iteration
      flow_dirn = flowPotential(sig);
      resid = flow_dirn * flow_incr - delta_dp;
      f = yieldFunction(sig, yield_stress);
      rep = -eqvpstrain + eqvpstrain_old - flow_incr * internalPotential();
 
      err1 = resid.L2norm();
      err2 = std::abs(f);
      err3 = std::abs(rep);
    }
 
    if (iter >= maxiter)
      mooseError("Constitutive failure");
 
    plastic_strain += delta_dp;
  }
}

Referenced by computeQpStress().

◆ usingTensorIndices()

FiniteStrainPlasticMaterial::usingTensorIndices	(	i_	,
		j_	,
		k_	,
		l_
	)

protected

make i,j,k,l available as tensor indices for mixedProduct

◆ validParams()

InputParameters FiniteStrainPlasticMaterial::validParams ( )

static

Definition at line 16 of file FiniteStrainPlasticMaterial.C.

{
  InputParameters params = ComputeStressBase::validParams();
 
  params.addRequiredParam<std::vector<Real>>(
      "yield_stress",
      "Input data as pairs of equivalent plastic strain and yield stress: Should "
      "start with equivalent plastic strain 0");
  params.addParam<Real>("rtol", 1e-8, "Plastic strain NR tolerance");
  params.addParam<Real>("ftol", 1e-4, "Consistency condition NR tolerance");
  params.addParam<Real>("eptol", 1e-7, "Equivalent plastic strain NR tolerance");
  params.addClassDescription("Associative J2 plasticity with isotropic hardening.");
 
  return params;
}

◆ yieldFunction()

Real FiniteStrainPlasticMaterial::yieldFunction	(	const RankTwoTensor &	stress,
		const Real	yield_stress
	)

protectedvirtual

Calculates the yield function.

Parameters

stress	the stress at which to calculate the yield function
yield_stress	the current value of the yield stress

Returns: equivalentstress - yield_stress

Definition at line 310 of file FiniteStrainPlasticMaterial.C.

{
  return getSigEqv(stress) - yield_stress;
}

Referenced by returnMap().

Member Data Documentation

◆ _base_name

const std::string ComputeGeneralStressBase::_base_name

protectedinherited

Base name prepended to all material property names to allow for multi-material systems.

Definition at line 43 of file ComputeGeneralStressBase.h.

Referenced by ComputeCosseratLinearElasticStress::initialSetup(), and ComputeLinearElasticStress::initialSetup().

◆ _deltaMixed

RankFourTensor FiniteStrainPlasticMaterial::_deltaMixed

protected

Definition at line 52 of file FiniteStrainPlasticMaterial.h.

Referenced by getJac().

◆ _deltaOuter

RankFourTensor FiniteStrainPlasticMaterial::_deltaOuter

protected

Definition at line 52 of file FiniteStrainPlasticMaterial.h.

Referenced by getJac().

◆ _elastic_strain

MaterialProperty<RankTwoTensor>& ComputeGeneralStressBase::_elastic_strain

protectedinherited

Elastic strain material property.

Definition at line 50 of file ComputeGeneralStressBase.h.

Referenced by ComputeSmearedCrackingStress::computeCrackStrainAndOrientation(), ComputeMultiPlasticityStress::computeQpStress(), computeQpStress(), ComputeCosseratLinearElasticStress::computeQpStress(), ComputeFiniteStrainElasticStress::computeQpStress(), ComputeLinearElasticStress::computeQpStress(), ComputeLinearViscoelasticStress::computeQpStress(), ComputeSmearedCrackingStress::computeQpStress(), ComputeMultipleInelasticStressBase::computeQpStressIntermediateConfiguration(), ComputeMultipleInelasticStressBase::finiteStrainRotation(), and ComputeGeneralStressBase::initQpStatefulProperties().

◆ _elasticity_tensor

const MaterialProperty<RankFourTensor>& FiniteStrainPlasticMaterial::_elasticity_tensor

protected

Elasticity tensor material property.

Definition at line 46 of file FiniteStrainPlasticMaterial.h.

Referenced by computeQpStress().

◆ _elasticity_tensor_name

const std::string FiniteStrainPlasticMaterial::_elasticity_tensor_name

protected

Name of the elasticity tensor material property.

Definition at line 44 of file FiniteStrainPlasticMaterial.h.

◆ _eptol

Real FiniteStrainPlasticMaterial::_eptol

protected

Definition at line 49 of file FiniteStrainPlasticMaterial.h.

Referenced by returnMap().

◆ _eqv_plastic_strain

MaterialProperty<Real>& FiniteStrainPlasticMaterial::_eqv_plastic_strain

protected

Definition at line 37 of file FiniteStrainPlasticMaterial.h.

Referenced by computeQpStress(), and initQpStatefulProperties().

◆ _eqv_plastic_strain_old

const MaterialProperty<Real>& FiniteStrainPlasticMaterial::_eqv_plastic_strain_old

protected

Definition at line 38 of file FiniteStrainPlasticMaterial.h.

Referenced by computeQpStress().

◆ _extra_stress

const MaterialProperty<RankTwoTensor>& ComputeGeneralStressBase::_extra_stress

protectedinherited

Extra stress tensor.

Definition at line 53 of file ComputeGeneralStressBase.h.

Referenced by ComputeGeneralStressBase::computeQpProperties().

◆ _ftol

Real FiniteStrainPlasticMaterial::_ftol

protected

Definition at line 48 of file FiniteStrainPlasticMaterial.h.

Referenced by returnMap().

◆ _initial_stress_fcn

std::vector<const Function *> ComputeGeneralStressBase::_initial_stress_fcn

protectedinherited

initial stress components

Definition at line 56 of file ComputeGeneralStressBase.h.

◆ _Jacobian_mult

MaterialProperty<RankFourTensor>& ComputeGeneralStressBase::_Jacobian_mult

protectedinherited

derivative of stress w.r.t. strain (_dstress_dstrain)

Definition at line 59 of file ComputeGeneralStressBase.h.

◆ _mechanical_strain

const MaterialProperty<RankTwoTensor>& ComputeGeneralStressBase::_mechanical_strain

protectedinherited

Mechanical strain material property.

Definition at line 46 of file ComputeGeneralStressBase.h.

Referenced by ComputeStrainIncrementBasedStress::computeQpStress(), computeQpStress(), ComputeCosseratLinearElasticStress::computeQpStress(), ComputeFiniteStrainElasticStress::computeQpStress(), ComputeLinearElasticPFFractureStress::computeQpStress(), ComputeLinearElasticStress::computeQpStress(), ComputeLinearViscoelasticStress::computeQpStress(), ComputeLinearElasticPFFractureStress::computeStrainSpectral(), ComputeLinearElasticPFFractureStress::computeStrainVolDev(), and ComputeLinearElasticPFFractureStress::computeStressSpectral().

◆ _plastic_strain

MaterialProperty<RankTwoTensor>& FiniteStrainPlasticMaterial::_plastic_strain

protected

Definition at line 35 of file FiniteStrainPlasticMaterial.h.

Referenced by computeQpStress(), and initQpStatefulProperties().

◆ _plastic_strain_old

const MaterialProperty<RankTwoTensor>& FiniteStrainPlasticMaterial::_plastic_strain_old

protected

Definition at line 36 of file FiniteStrainPlasticMaterial.h.

Referenced by computeQpStress().

◆ _rotation_increment

const MaterialProperty<RankTwoTensor>& FiniteStrainPlasticMaterial::_rotation_increment

protected

Definition at line 42 of file FiniteStrainPlasticMaterial.h.

Referenced by computeQpStress().

◆ _rtol

Real FiniteStrainPlasticMaterial::_rtol

protected

Definition at line 47 of file FiniteStrainPlasticMaterial.h.

Referenced by returnMap().

◆ _strain_increment

const MaterialProperty<RankTwoTensor>& FiniteStrainPlasticMaterial::_strain_increment

protected

Definition at line 41 of file FiniteStrainPlasticMaterial.h.

Referenced by computeQpStress().

◆ _stress

MaterialProperty<RankTwoTensor>& ComputeGeneralStressBase::_stress

protectedinherited

Stress material property.

Definition at line 48 of file ComputeGeneralStressBase.h.

◆ _stress_old

const MaterialProperty<RankTwoTensor>& FiniteStrainPlasticMaterial::_stress_old

protected

The old stress tensor.

Definition at line 40 of file FiniteStrainPlasticMaterial.h.

Referenced by computeQpStress().

◆ _yield_stress_vector

std::vector<Real> FiniteStrainPlasticMaterial::_yield_stress_vector

protected

Definition at line 34 of file FiniteStrainPlasticMaterial.h.

Referenced by getdYieldStressdPlasticStrain(), and getYieldStress().

The documentation for this class was generated from the following files:

solid_mechanics/include/materials/FiniteStrainPlasticMaterial.h
solid_mechanics/src/materials/FiniteStrainPlasticMaterial.C

Public Types

Public Member Functions

Static Public Member Functions

Protected Member Functions

Protected Attributes

Private Member Functions

Private Attributes

Detailed Description

Constructor & Destructor Documentation

◆ FiniteStrainPlasticMaterial()

Member Function Documentation

◆ computeQpProperties()

◆ computeQpStress()

◆ dyieldFunction_dinternal()

◆ dyieldFunction_dstress()

◆ flowPotential()

◆ getdYieldStressdPlasticStrain()

◆ getJac()

◆ getSigEqv()

◆ getYieldStress()

◆ initQpStatefulProperties()

◆ internalPotential()

◆ returnMap()

◆ usingTensorIndices()

◆ validParams()

◆ yieldFunction()

Member Data Documentation

◆ _base_name

◆ _deltaMixed

◆ _deltaOuter

◆ _elastic_strain

◆ _elasticity_tensor

◆ _elasticity_tensor_name

◆ _eptol

◆ _eqv_plastic_strain

◆ _eqv_plastic_strain_old

◆ _extra_stress

◆ _ftol

◆ _initial_stress_fcn

◆ _Jacobian_mult

◆ _mechanical_strain

◆ _plastic_strain

◆ _plastic_strain_old

◆ _rotation_increment

◆ _rtol

◆ _strain_increment

◆ _stress

◆ _stress_old

◆ _yield_stress_vector