www.mooseframework.org
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
FiniteStrainPlasticMaterial Class Reference

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 Member Functions

 FiniteStrainPlasticMaterial (const InputParameters &parameters)
 

Protected Member Functions

virtual void computeQpStress ()
 Compute the stress and store it in the _stress material property for the current quadrature point. More...
 
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. More...
 
virtual Real yieldFunction (const RankTwoTensor &stress, const Real yield_stress)
 Calculates the yield function. More...
 
virtual RankTwoTensor dyieldFunction_dstress (const RankTwoTensor &stress)
 Derivative of yieldFunction with respect to the stress. More...
 
virtual Real dyieldFunction_dinternal (const Real equivalent_plastic_strain)
 Derivative of yieldFunction with respect to the equivalent plastic strain. More...
 
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. More...
 
virtual Real internalPotential ()
 The internal potential. More...
 
Real getSigEqv (const RankTwoTensor &stress)
 Equivalent stress. More...
 
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. More...
 
Real getYieldStress (const Real equivalent_plastic_strain)
 yield stress as a function of equivalent plastic strain. More...
 
Real getdYieldStressdPlasticStrain (const Real equivalent_plastic_strain)
 d(yieldstress)/d(equivalent plastic strain) More...
 
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
 
const MaterialProperty< RankTwoTensor > & _strain_increment
 
const MaterialProperty< RankTwoTensor > & _rotation_increment
 
const std::string _elasticity_tensor_name
 Name of the elasticity tensor material property. More...
 
const MaterialProperty< RankFourTensor > & _elasticity_tensor
 Elasticity tensor material property. More...
 
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. More...
 
const MaterialProperty< RankTwoTensor > & _mechanical_strain
 Mechanical strain material property. More...
 
MaterialProperty< RankTwoTensor > & _stress
 Stress material property. More...
 
MaterialProperty< RankTwoTensor > & _elastic_strain
 Elastic strain material property. More...
 
const MaterialProperty< RankTwoTensor > & _extra_stress
 Extra stress tensor. More...
 
std::vector< Function * > _initial_stress_fcn
 initial stress components More...
 
MaterialProperty< RankFourTensor > & _Jacobian_mult
 derivative of stress w.r.t. strain (_dstress_dstrain) More...
 

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 30 of file FiniteStrainPlasticMaterial.h.

Constructor & Destructor Documentation

◆ FiniteStrainPlasticMaterial()

FiniteStrainPlasticMaterial::FiniteStrainPlasticMaterial ( const InputParameters &  parameters)

Definition at line 33 of file FiniteStrainPlasticMaterial.C.

34  : ComputeStressBase(parameters),
35  _yield_stress_vector(getParam<std::vector<Real>>("yield_stress")), // Read from input file
36  _plastic_strain(declareProperty<RankTwoTensor>(_base_name + "plastic_strain")),
37  _plastic_strain_old(getMaterialPropertyOld<RankTwoTensor>(_base_name + "plastic_strain")),
38  _eqv_plastic_strain(declareProperty<Real>(_base_name + "eqv_plastic_strain")),
39  _eqv_plastic_strain_old(getMaterialPropertyOld<Real>(_base_name + "eqv_plastic_strain")),
40  _stress_old(getMaterialPropertyOld<RankTwoTensor>(_base_name + "stress")),
41  _strain_increment(getMaterialProperty<RankTwoTensor>(_base_name + "strain_increment")),
42  _rotation_increment(getMaterialProperty<RankTwoTensor>(_base_name + "rotation_increment")),
43  _elasticity_tensor_name(_base_name + "elasticity_tensor"),
44  _elasticity_tensor(getMaterialPropertyByName<RankFourTensor>(_elasticity_tensor_name)),
45  _rtol(getParam<Real>("rtol")),
46  _ftol(getParam<Real>("ftol")),
47  _eptol(getParam<Real>("eptol")),
48  _deltaOuter(RankTwoTensor::Identity().outerProduct(RankTwoTensor::Identity())),
49  _deltaMixed(RankTwoTensor::Identity().mixedProductIkJl(RankTwoTensor::Identity()))
50 {
51 }
const MaterialProperty< RankTwoTensor > & _strain_increment
MaterialProperty< Real > & _eqv_plastic_strain
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Elasticity tensor material property.
const MaterialProperty< Real > & _eqv_plastic_strain_old
const std::string _elasticity_tensor_name
Name of the elasticity tensor material property.
const std::string _base_name
Base name prepended to all material property names to allow for multi-material systems.
ComputeStressBase(const InputParameters &parameters)
const MaterialProperty< RankTwoTensor > & _plastic_strain_old
const MaterialProperty< RankTwoTensor > & _rotation_increment
MaterialProperty< RankTwoTensor > & _plastic_strain
const MaterialProperty< RankTwoTensor > & _stress_old

Member Function Documentation

◆ computeQpProperties()

void ComputeStressBase::computeQpProperties ( )
overrideprotectedvirtualinherited

Definition at line 49 of file ComputeStressBase.C.

50 {
52 
53  // Add in extra stress
54  _stress[_qp] += _extra_stress[_qp];
55 }
virtual void computeQpStress()=0
Compute the stress and store it in the _stress material property for the current quadrature point...
MaterialProperty< RankTwoTensor > & _stress
Stress material property.
const MaterialProperty< RankTwoTensor > & _extra_stress
Extra stress tensor.

◆ computeQpStress()

void FiniteStrainPlasticMaterial::computeQpStress ( )
protectedvirtual

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

Implements ComputeStressBase.

Definition at line 62 of file FiniteStrainPlasticMaterial.C.

63 {
64 
65  // perform the return-mapping algorithm
69  _strain_increment[_qp],
70  _elasticity_tensor[_qp],
71  _stress[_qp],
73  _plastic_strain[_qp]);
74 
75  // Rotate the stress tensor to the current configuration
76  _stress[_qp] = _rotation_increment[_qp] * _stress[_qp] * _rotation_increment[_qp].transpose();
77 
78  // Rotate plastic strain tensor to the current configuration
79  _plastic_strain[_qp] =
80  _rotation_increment[_qp] * _plastic_strain[_qp] * _rotation_increment[_qp].transpose();
81 
82  // Calculate the elastic strain_increment
84 
86 }
MaterialProperty< RankFourTensor > & _Jacobian_mult
derivative of stress w.r.t. strain (_dstress_dstrain)
const MaterialProperty< RankTwoTensor > & _strain_increment
MaterialProperty< Real > & _eqv_plastic_strain
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.
MaterialProperty< RankTwoTensor > & _stress
Stress material property.
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Elasticity tensor material property.
const MaterialProperty< Real > & _eqv_plastic_strain_old
const MaterialProperty< RankTwoTensor > & _mechanical_strain
Mechanical strain material property.
const MaterialProperty< RankTwoTensor > & _plastic_strain_old
MaterialProperty< RankTwoTensor > & _elastic_strain
Elastic strain material property.
const MaterialProperty< RankTwoTensor > & _rotation_increment
MaterialProperty< RankTwoTensor > & _plastic_strain
const MaterialProperty< RankTwoTensor > & _stress_old

◆ 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 325 of file FiniteStrainPlasticMaterial.C.

Referenced by returnMap().

326 {
327  return -getdYieldStressdPlasticStrain(equivalent_plastic_strain);
328 }
Real getdYieldStressdPlasticStrain(const Real equivalent_plastic_strain)
d(yieldstress)/d(equivalent plastic strain)

◆ dyieldFunction_dstress()

RankTwoTensor FiniteStrainPlasticMaterial::dyieldFunction_dstress ( const RankTwoTensor stress)
protectedvirtual

Derivative of yieldFunction with respect to the stress.

Definition at line 317 of file FiniteStrainPlasticMaterial.C.

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

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

◆ 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 331 of file FiniteStrainPlasticMaterial.C.

Referenced by getJac(), and returnMap().

332 {
333  return dyieldFunction_dstress(sig); // this plasticity model assumes associative flow
334 }
virtual RankTwoTensor dyieldFunction_dstress(const RankTwoTensor &stress)
Derivative of yieldFunction with respect to the stress.

◆ getdYieldStressdPlasticStrain()

Real FiniteStrainPlasticMaterial::getdYieldStressdPlasticStrain ( const Real  equivalent_plastic_strain)
protected

d(yieldstress)/d(equivalent plastic strain)

Definition at line 415 of file FiniteStrainPlasticMaterial.C.

Referenced by dyieldFunction_dinternal().

416 {
417  unsigned nsize;
418 
419  nsize = _yield_stress_vector.size();
420 
421  if (_yield_stress_vector[0] > 0.0 || nsize % 2 > 0) // Error check for input inconsitency
422  mooseError("Error in yield stress input: Should be a vector with eqv plastic strain and yield "
423  "stress pair values.\n");
424 
425  unsigned int ind = 0;
426 
427  while (ind < nsize)
428  {
429  if (ind + 2 < nsize)
430  {
431  if (eqpe >= _yield_stress_vector[ind] && eqpe < _yield_stress_vector[ind + 2])
432  return (_yield_stress_vector[ind + 3] - _yield_stress_vector[ind + 1]) /
433  (_yield_stress_vector[ind + 2] - _yield_stress_vector[ind]);
434  }
435  else
436  return 0.0;
437 
438  ind += 2;
439  }
440 
441  return 0.0;
442 }

◆ 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
sigstress
E_ijklelasticity tensor (sig = E*(strain - plastic_strain))
flow_incrconsistency parameter
dresid_dsigthe required derivative (this is an output variable)

Definition at line 350 of file FiniteStrainPlasticMaterial.C.

Referenced by returnMap().

354 {
355  RankTwoTensor sig_dev, df_dsig, flow_dirn;
356  RankTwoTensor dfi_dft, dfi_dsig;
357  RankFourTensor dft_dsig, dfd_dft, dfd_dsig;
358  Real sig_eqv;
359  Real f1, f2, f3;
360  RankFourTensor temp;
361 
362  sig_dev = sig.deviatoric();
363  sig_eqv = getSigEqv(sig);
364  df_dsig = dyieldFunction_dstress(sig);
365  flow_dirn = flowPotential(sig);
366 
367  f1 = 3.0 / (2.0 * sig_eqv);
368  f2 = f1 / 3.0;
369  f3 = 9.0 / (4.0 * Utility::pow<3>(sig_eqv));
370 
371  dft_dsig = f1 * _deltaMixed - f2 * _deltaOuter - f3 * sig_dev.outerProduct(sig_dev);
372 
373  dfd_dsig = dft_dsig;
374  dresid_dsig = E_ijkl.invSymm() + dfd_dsig * flow_incr;
375 }
virtual RankTwoTensor flowPotential(const RankTwoTensor &stress)
Flow potential, which in this case is just dyieldFunction_dstress because we are doing associative fl...
Real getSigEqv(const RankTwoTensor &stress)
Equivalent stress.
virtual RankTwoTensor dyieldFunction_dstress(const RankTwoTensor &stress)
Derivative of yieldFunction with respect to the stress.

◆ getSigEqv()

Real FiniteStrainPlasticMaterial::getSigEqv ( const RankTwoTensor stress)
protected

Equivalent stress.

Definition at line 343 of file FiniteStrainPlasticMaterial.C.

Referenced by getJac(), and yieldFunction().

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

◆ 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 379 of file FiniteStrainPlasticMaterial.C.

Referenced by returnMap().

380 {
381  unsigned nsize;
382 
383  nsize = _yield_stress_vector.size();
384 
385  if (_yield_stress_vector[0] > 0.0 || nsize % 2 > 0) // Error check for input inconsitency
386  mooseError("Error in yield stress input: Should be a vector with eqv plastic strain and yield "
387  "stress pair values.\n");
388 
389  unsigned int ind = 0;
390  Real tol = 1e-8;
391 
392  while (ind < nsize)
393  {
394  if (std::abs(eqpe - _yield_stress_vector[ind]) < tol)
395  return _yield_stress_vector[ind + 1];
396 
397  if (ind + 2 < nsize)
398  {
399  if (eqpe > _yield_stress_vector[ind] && eqpe < _yield_stress_vector[ind + 2])
400  return _yield_stress_vector[ind + 1] +
401  (eqpe - _yield_stress_vector[ind]) /
402  (_yield_stress_vector[ind + 2] - _yield_stress_vector[ind]) *
403  (_yield_stress_vector[ind + 3] - _yield_stress_vector[ind + 1]);
404  }
405  else
406  return _yield_stress_vector[nsize - 1];
407 
408  ind += 2;
409  }
410 
411  return 0.0;
412 }
const double tol
Definition: Setup.h:19

◆ initQpStatefulProperties()

void FiniteStrainPlasticMaterial::initQpStatefulProperties ( )
protectedvirtual

Reimplemented from ComputeStressBase.

Definition at line 54 of file FiniteStrainPlasticMaterial.C.

55 {
57  _plastic_strain[_qp].zero();
58  _eqv_plastic_strain[_qp] = 0.0;
59 }
MaterialProperty< Real > & _eqv_plastic_strain
virtual void initQpStatefulProperties() override
MaterialProperty< RankTwoTensor > & _plastic_strain

◆ internalPotential()

Real FiniteStrainPlasticMaterial::internalPotential ( )
protectedvirtual

The internal potential.

For associative J2 plasticity this is just -1

Definition at line 337 of file FiniteStrainPlasticMaterial.C.

Referenced by returnMap().

338 {
339  return -1;
340 }

◆ 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_oldThe stress at the previous "time" step
eqvpstrain_oldThe equivalent plastic strain at the previous "time" step
plastic_strain_oldThe value of plastic strain at the previous "time" step
delta_dThe total strain increment for this "time" step
E_ijklThe elasticity tensor. If no plasiticity then sig_new = sig_old + E_ijkl*delta_d
sigThe stress after returning to the yield surface (this is an output variable)
eqvpstrainThe equivalent plastic strain after returning to the yield surface (this is an output variable)
plastic_strainThe 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 150 of file FiniteStrainPlasticMaterial.C.

Referenced by computeQpStress().

158 {
159  // the yield function, must be non-positive
160  // Newton-Raphson sets this to zero if trial stress enters inadmissible region
161  Real f;
162 
163  // the consistency parameter, must be non-negative
164  // change in plastic strain in this timestep = flow_incr*flow_potential
165  Real flow_incr = 0.0;
166 
167  // direction of flow defined by the potential
168  RankTwoTensor flow_dirn;
169 
170  // Newton-Raphson sets this zero
171  // resid_ij = flow_incr*flow_dirn_ij - (plastic_strain - plastic_strain_old)
172  RankTwoTensor resid;
173 
174  // Newton-Raphson sets this zero
175  // rep = -flow_incr*internalPotential - (eqvpstrain - eqvpstrain_old)
176  Real rep;
177 
178  // change in the stress (sig) in a Newton-Raphson iteration
179  RankTwoTensor ddsig;
180 
181  // change in the consistency parameter in a Newton-Raphson iteration
182  Real dflow_incr = 0.0;
183 
184  // change in equivalent plastic strain in one Newton-Raphson iteration
185  Real deqvpstrain = 0.0;
186 
187  // convenience variable that holds the change in plastic strain incurred during the return
188  // delta_dp = plastic_strain - plastic_strain_old
189  // delta_dp = E^{-1}*(trial_stress - sig), where trial_stress = E*(strain - plastic_strain_old)
190  RankTwoTensor delta_dp;
191 
192  // d(yieldFunction)/d(stress)
193  RankTwoTensor df_dsig;
194 
195  // d(resid_ij)/d(sigma_kl)
196  RankFourTensor dr_dsig;
197 
198  // 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
199  // zero otherwise
200  RankFourTensor dr_dsig_inv;
201 
202  // d(yieldFunction)/d(eqvpstrain)
203  Real fq;
204 
205  // yield stress at the start of this "time step" (ie, evaluated with
206  // eqvpstrain_old). It is held fixed during the Newton-Raphson return,
207  // even if eqvpstrain != eqvpstrain_old.
208  Real yield_stress;
209 
210  // measures of whether the Newton-Raphson process has converged
211  Real err1, err2, err3;
212 
213  // number of Newton-Raphson iterations performed
214  unsigned int iter = 0;
215 
216  // maximum number of Newton-Raphson iterations allowed
217  unsigned int maxiter = 100;
218 
219  // plastic loading occurs if yieldFunction > toly
220  Real toly = 1.0e-8;
221 
222  // Assume this strain increment does not induce any plasticity
223  // This is the elastic-predictor
224  sig = sig_old + E_ijkl * delta_d; // the trial stress
225  eqvpstrain = eqvpstrain_old;
226  plastic_strain = plastic_strain_old;
227 
228  yield_stress = getYieldStress(eqvpstrain); // yield stress at this equivalent plastic strain
229  if (yieldFunction(sig, yield_stress) > toly)
230  {
231  // the sig just calculated is inadmissable. We must return to the yield surface.
232  // This is done iteratively, using a Newton-Raphson process.
233 
234  delta_dp.zero();
235 
236  sig = sig_old + E_ijkl * delta_d; // this is the elastic predictor
237 
238  flow_dirn = flowPotential(sig);
239 
240  resid = flow_dirn * flow_incr - delta_dp; // Residual 1 - refer Hughes Simo
241  f = yieldFunction(sig, yield_stress);
242  rep = -eqvpstrain + eqvpstrain_old - flow_incr * internalPotential(); // Residual 3 rep=0
243 
244  err1 = resid.L2norm();
245  err2 = std::abs(f);
246  err3 = std::abs(rep);
247 
248  while ((err1 > _rtol || err2 > _ftol || err3 > _eptol) &&
249  iter < maxiter) // Stress update iteration (hardness fixed)
250  {
251  iter++;
252 
253  df_dsig = dyieldFunction_dstress(sig);
254  getJac(sig, E_ijkl, flow_incr, dr_dsig); // gets dr_dsig = d(resid_ij)/d(sig_kl)
255  fq = dyieldFunction_dinternal(eqvpstrain); // d(f)/d(eqvpstrain)
256 
268  dr_dsig_inv = dr_dsig.invSymm();
269 
277  dflow_incr = (f - df_dsig.doubleContraction(dr_dsig_inv * resid) + fq * rep) /
278  (df_dsig.doubleContraction(dr_dsig_inv * flow_dirn) - fq);
279  ddsig =
280  dr_dsig_inv *
281  (-resid -
282  flow_dirn * dflow_incr); // from solving the top row of linear system, given dflow_incr
283  deqvpstrain =
284  rep + dflow_incr; // from solving the bottom row of linear system, given dflow_incr
285 
286  // update the variables
287  flow_incr += dflow_incr;
288  delta_dp -= E_ijkl.invSymm() * ddsig;
289  sig += ddsig;
290  eqvpstrain += deqvpstrain;
291 
292  // evaluate the RHS equations ready for next Newton-Raphson iteration
293  flow_dirn = flowPotential(sig);
294  resid = flow_dirn * flow_incr - delta_dp;
295  f = yieldFunction(sig, yield_stress);
296  rep = -eqvpstrain + eqvpstrain_old - flow_incr * internalPotential();
297 
298  err1 = resid.L2norm();
299  err2 = std::abs(f);
300  err3 = std::abs(rep);
301  }
302 
303  if (iter >= maxiter)
304  mooseError("Constitutive failure");
305 
306  plastic_strain += delta_dp;
307  }
308 }
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 fl...
virtual Real yieldFunction(const RankTwoTensor &stress, const Real yield_stress)
Calculates the yield function.
Real getYieldStress(const Real equivalent_plastic_strain)
yield stress as a function of equivalent plastic strain.
virtual RankTwoTensor dyieldFunction_dstress(const RankTwoTensor &stress)
Derivative of yieldFunction with respect to the stress.
virtual Real internalPotential()
The internal potential.
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}...

◆ yieldFunction()

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

Calculates the yield function.

Parameters
stressthe stress at which to calculate the yield function
yield_stressthe current value of the yield stress
Returns
equivalentstress - yield_stress

Definition at line 311 of file FiniteStrainPlasticMaterial.C.

Referenced by returnMap().

312 {
313  return getSigEqv(stress) - yield_stress;
314 }
Real getSigEqv(const RankTwoTensor &stress)
Equivalent stress.

Member Data Documentation

◆ _base_name

const std::string ComputeStressBase::_base_name
protectedinherited

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

Definition at line 44 of file ComputeStressBase.h.

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

◆ _deltaMixed

RankFourTensor FiniteStrainPlasticMaterial::_deltaMixed
protected

Definition at line 55 of file FiniteStrainPlasticMaterial.h.

Referenced by getJac().

◆ _deltaOuter

RankFourTensor FiniteStrainPlasticMaterial::_deltaOuter
protected

Definition at line 55 of file FiniteStrainPlasticMaterial.h.

Referenced by getJac().

◆ _elastic_strain

MaterialProperty<RankTwoTensor>& ComputeStressBase::_elastic_strain
protectedinherited

◆ _elasticity_tensor

const MaterialProperty<RankFourTensor>& FiniteStrainPlasticMaterial::_elasticity_tensor
protected

Elasticity tensor material property.

Definition at line 49 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 47 of file FiniteStrainPlasticMaterial.h.

◆ _eptol

Real FiniteStrainPlasticMaterial::_eptol
protected

Definition at line 52 of file FiniteStrainPlasticMaterial.h.

Referenced by returnMap().

◆ _eqv_plastic_strain

MaterialProperty<Real>& FiniteStrainPlasticMaterial::_eqv_plastic_strain
protected

Definition at line 41 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 42 of file FiniteStrainPlasticMaterial.h.

Referenced by computeQpStress().

◆ _extra_stress

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_extra_stress
protectedinherited

Extra stress tensor.

Definition at line 54 of file ComputeStressBase.h.

Referenced by ComputeStressBase::computeQpProperties().

◆ _ftol

Real FiniteStrainPlasticMaterial::_ftol
protected

Definition at line 51 of file FiniteStrainPlasticMaterial.h.

Referenced by returnMap().

◆ _initial_stress_fcn

std::vector<Function *> ComputeStressBase::_initial_stress_fcn
protectedinherited

initial stress components

Definition at line 57 of file ComputeStressBase.h.

◆ _Jacobian_mult

MaterialProperty<RankFourTensor>& ComputeStressBase::_Jacobian_mult
protectedinherited

◆ _mechanical_strain

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_mechanical_strain
protectedinherited

◆ _plastic_strain

MaterialProperty<RankTwoTensor>& FiniteStrainPlasticMaterial::_plastic_strain
protected

Definition at line 39 of file FiniteStrainPlasticMaterial.h.

Referenced by computeQpStress(), and initQpStatefulProperties().

◆ _plastic_strain_old

const MaterialProperty<RankTwoTensor>& FiniteStrainPlasticMaterial::_plastic_strain_old
protected

Definition at line 40 of file FiniteStrainPlasticMaterial.h.

Referenced by computeQpStress().

◆ _rotation_increment

const MaterialProperty<RankTwoTensor>& FiniteStrainPlasticMaterial::_rotation_increment
protected

Definition at line 45 of file FiniteStrainPlasticMaterial.h.

Referenced by computeQpStress().

◆ _rtol

Real FiniteStrainPlasticMaterial::_rtol
protected

Definition at line 50 of file FiniteStrainPlasticMaterial.h.

Referenced by returnMap().

◆ _strain_increment

const MaterialProperty<RankTwoTensor>& FiniteStrainPlasticMaterial::_strain_increment
protected

Definition at line 44 of file FiniteStrainPlasticMaterial.h.

Referenced by computeQpStress().

◆ _stress

MaterialProperty<RankTwoTensor>& ComputeStressBase::_stress
protectedinherited

Stress material property.

Definition at line 49 of file ComputeStressBase.h.

Referenced by ComputeMultipleInelasticCosseratStress::computeAdmissibleState(), ComputeMultipleInelasticStress::computeAdmissibleState(), ComputeStressBase::computeQpProperties(), ComputeStrainIncrementBasedStress::computeQpStress(), ComputeLinearElasticStress::computeQpStress(), ComputeDamageStress::computeQpStress(), ComputeFiniteStrainElasticStress::computeQpStress(), ComputeCosseratLinearElasticStress::computeQpStress(), ComputeSmearedCrackingStress::computeQpStress(), computeQpStress(), ComputeLinearElasticPFFractureStress::computeQpStress(), ComputeMultiPlasticityStress::computeQpStress(), ComputeLinearViscoelasticStress::computeQpStress(), ComputeMultipleInelasticStress::computeQpStressIntermediateConfiguration(), ComputeLinearElasticPFFractureStress::computeStrainSpectral(), ComputeLinearElasticPFFractureStress::computeStrainVolDev(), ComputeLinearElasticPFFractureStress::computeStressSpectral(), ComputeMultipleInelasticStress::finiteStrainRotation(), ComputeStressBase::initQpStatefulProperties(), FiniteStrainCrystalPlasticity::initQpStatefulProperties(), FiniteStrainUObasedCP::initQpStatefulProperties(), FiniteStrainHyperElasticViscoPlastic::initQpStatefulProperties(), ComputeMultiPlasticityStress::postReturnMap(), FiniteStrainUObasedCP::postSolveQp(), FiniteStrainHyperElasticViscoPlastic::postSolveQp(), FiniteStrainCrystalPlasticity::postSolveQp(), ComputeSmearedCrackingStress::updateCrackingStateAndStress(), ComputeMultipleInelasticStress::updateQpState(), and ComputeMultipleInelasticStress::updateQpStateSingleModel().

◆ _stress_old

const MaterialProperty<RankTwoTensor>& FiniteStrainPlasticMaterial::_stress_old
protected

Definition at line 43 of file FiniteStrainPlasticMaterial.h.

Referenced by computeQpStress().

◆ _yield_stress_vector

std::vector<Real> FiniteStrainPlasticMaterial::_yield_stress_vector
protected

Definition at line 38 of file FiniteStrainPlasticMaterial.h.

Referenced by getdYieldStressdPlasticStrain(), and getYieldStress().


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