www.mooseframework.org
FiniteStrainPlasticMaterial.h
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #pragma once
11 // Original class author: A.M. Jokisaari, O. Heinonen
12 
13 
14 #include "ComputeStressBase.h"
15 
17 
18 template <>
20 
30 {
31 public:
32  FiniteStrainPlasticMaterial(const InputParameters & parameters);
33 
34 protected:
35  virtual void computeQpStress();
36  virtual void initQpStatefulProperties();
37  std::vector<Real> _yield_stress_vector;
38  MaterialProperty<RankTwoTensor> & _plastic_strain;
39  const MaterialProperty<RankTwoTensor> & _plastic_strain_old;
40  MaterialProperty<Real> & _eqv_plastic_strain;
41  const MaterialProperty<Real> & _eqv_plastic_strain_old;
42  const MaterialProperty<RankTwoTensor> & _stress_old;
43  const MaterialProperty<RankTwoTensor> & _strain_increment;
44  const MaterialProperty<RankTwoTensor> & _rotation_increment;
46  const std::string _elasticity_tensor_name;
48  const MaterialProperty<RankFourTensor> & _elasticity_tensor;
49  Real _rtol;
50  Real _ftol;
51  Real _eptol;
52 
53  // outer and mixed product of the delta function tensors
55 
73  virtual void returnMap(const RankTwoTensor & sig_old,
74  const Real eqvpstrain_old,
75  const RankTwoTensor & plastic_strain_old,
76  const RankTwoTensor & delta_d,
77  const RankFourTensor & E_ijkl,
78  RankTwoTensor & sig,
79  Real & eqvpstrain,
80  RankTwoTensor & plastic_strain);
81 
88  virtual Real yieldFunction(const RankTwoTensor & stress, const Real yield_stress);
89 
93  virtual RankTwoTensor dyieldFunction_dstress(const RankTwoTensor & stress);
94 
98  virtual Real dyieldFunction_dinternal(const Real equivalent_plastic_strain);
99 
105  virtual RankTwoTensor flowPotential(const RankTwoTensor & stress);
106 
110  virtual Real internalPotential();
111 
115  Real getSigEqv(const RankTwoTensor & stress);
116 
125  virtual void getJac(const RankTwoTensor & sig,
126  const RankFourTensor & E_ijkl,
127  Real flow_incr,
128  RankFourTensor & dresid_dsig);
129 
134  Real getYieldStress(const Real equivalent_plastic_strain);
135 
139  Real getdYieldStressdPlasticStrain(const Real equivalent_plastic_strain);
140 };
141 
virtual void computeQpStress()
Compute the stress and store it in the _stress material property for the current quadrature point...
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
Elasticity tensor material property.
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.
const std::string _elasticity_tensor_name
Name of the elasticity tensor material property.
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