www.mooseframework.org
FiniteStrainPlasticMaterial.h
Go to the documentation of this file.
1 /****************************************************************/
2 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
3 /* */
4 /* All contents are licensed under LGPL V2.1 */
5 /* See LICENSE for full restrictions */
6 /****************************************************************/
7 // Original class author: A.M. Jokisaari, O. Heinonen
8 
9 #ifndef FINITESTRAINPLASTICMATERIAL_H
10 #define FINITESTRAINPLASTICMATERIAL_H
11 
12 #include "ComputeStressBase.h"
13 
15 
16 template <>
18 
28 {
29 public:
30  FiniteStrainPlasticMaterial(const InputParameters & parameters);
31 
32 protected:
33  virtual void computeQpStress();
34  virtual void initQpStatefulProperties();
35 
36  std::vector<Real> _yield_stress_vector;
37  MaterialProperty<RankTwoTensor> & _plastic_strain;
38  const MaterialProperty<RankTwoTensor> & _plastic_strain_old;
39  MaterialProperty<Real> & _eqv_plastic_strain;
40  const MaterialProperty<Real> & _eqv_plastic_strain_old;
41  const MaterialProperty<RankTwoTensor> & _stress_old;
42  const MaterialProperty<RankTwoTensor> & _strain_increment;
43  const MaterialProperty<RankTwoTensor> & _rotation_increment;
44  const MaterialProperty<RankFourTensor> & _elasticity_tensor;
45  Real _rtol;
46  Real _ftol;
47  Real _eptol;
48 
49  // outer and mixed product of the delta function tensors
50  RankFourTensor _deltaOuter, _deltaMixed;
51 
69  virtual void returnMap(const RankTwoTensor & sig_old,
70  const Real eqvpstrain_old,
71  const RankTwoTensor & plastic_strain_old,
72  const RankTwoTensor & delta_d,
73  const RankFourTensor & E_ijkl,
74  RankTwoTensor & sig,
75  Real & eqvpstrain,
76  RankTwoTensor & plastic_strain);
77 
84  virtual Real yieldFunction(const RankTwoTensor & stress, const Real yield_stress);
85 
89  virtual RankTwoTensor dyieldFunction_dstress(const RankTwoTensor & stress);
90 
94  virtual Real dyieldFunction_dinternal(const Real equivalent_plastic_strain);
95 
101  virtual RankTwoTensor flowPotential(const RankTwoTensor & stress);
102 
106  virtual Real internalPotential();
107 
111  Real getSigEqv(const RankTwoTensor & stress);
112 
121  virtual void getJac(const RankTwoTensor & sig,
122  const RankFourTensor & E_ijkl,
123  Real flow_incr,
124  RankFourTensor & dresid_dsig);
125 
130  Real getYieldStress(const Real equivalent_plastic_strain);
131 
135  Real getdYieldStressdPlasticStrain(const Real equivalent_plastic_strain);
136 };
137 
138 #endif // FINITESTRAINPLASTICMATERIAL_H
FiniteStrainPlasticMaterial implements rate-independent associative J2 plasticity with isotropic hard...
ComputeStressBase is the base class for stress tensors.
const MaterialProperty< RankTwoTensor > & _strain_increment
MaterialProperty< Real > & _eqv_plastic_strain
virtual Real dyieldFunction_dinternal(const Real equivalent_plastic_strain)
Derivative of yieldFunction with respect to the equivalent 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.
virtual RankTwoTensor flowPotential(const RankTwoTensor &stress)
Flow potential, which in this case is just dyieldFunction_dstress because we are doing associative fl...
const MaterialProperty< RankFourTensor > & _elasticity_tensor
const MaterialProperty< Real > & _eqv_plastic_strain_old
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.
FiniteStrainPlasticMaterial(const InputParameters &parameters)
InputParameters validParams< FiniteStrainPlasticMaterial >()
const MaterialProperty< RankTwoTensor > & _plastic_strain_old
Real getSigEqv(const RankTwoTensor &stress)
Equivalent stress.
virtual RankTwoTensor dyieldFunction_dstress(const RankTwoTensor &stress)
Derivative of yieldFunction with respect to the stress.
virtual Real internalPotential()
The internal potential.
const MaterialProperty< RankTwoTensor > & _rotation_increment
Real getdYieldStressdPlasticStrain(const Real equivalent_plastic_strain)
d(yieldstress)/d(equivalent plastic strain)
MaterialProperty< RankTwoTensor > & _plastic_strain
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}...
const MaterialProperty< RankTwoTensor > & _stress_old