www.mooseframework.org
TensorMechanicsPlasticJ2.C
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 /****************************************************************/
8 
9 template <>
10 InputParameters
12 {
13  InputParameters params = validParams<TensorMechanicsPlasticModel>();
14  params.addRequiredParam<UserObjectName>(
15  "yield_strength",
16  "A TensorMechanicsHardening UserObject that defines hardening of the yield strength");
17  params.addRangeCheckedParam<unsigned>(
18  "max_iterations", 10, "max_iterations>0", "Maximum iterations for custom J2 return map");
19  params.addParam<bool>("use_custom_returnMap",
20  true,
21  "Whether to use the custom returnMap "
22  "algorithm. Set to true if you are using "
23  "isotropic elasticity.");
24  params.addParam<bool>("use_custom_cto",
25  true,
26  "Whether to use the custom consistent tangent "
27  "operator computations. Set to true if you are "
28  "using isotropic elasticity.");
29  params.addClassDescription("J2 plasticity, associative, with hardening");
30 
31  return params;
32 }
33 
34 TensorMechanicsPlasticJ2::TensorMechanicsPlasticJ2(const InputParameters & parameters)
35  : TensorMechanicsPlasticModel(parameters),
36  _strength(getUserObject<TensorMechanicsHardeningModel>("yield_strength")),
37  _max_iters(getParam<unsigned>("max_iterations")),
38  _use_custom_returnMap(getParam<bool>("use_custom_returnMap")),
39  _use_custom_cto(getParam<bool>("use_custom_cto"))
40 {
41 }
42 
43 Real
44 TensorMechanicsPlasticJ2::yieldFunction(const RankTwoTensor & stress, Real intnl) const
45 {
46  return std::sqrt(3.0 * stress.secondInvariant()) - yieldStrength(intnl);
47 }
48 
49 RankTwoTensor
50 TensorMechanicsPlasticJ2::dyieldFunction_dstress(const RankTwoTensor & stress, Real /*intnl*/) const
51 {
52  Real sII = stress.secondInvariant();
53  if (sII == 0.0)
54  return RankTwoTensor();
55  else
56  return 0.5 * std::sqrt(3.0 / sII) * stress.dsecondInvariant();
57 }
58 
59 Real
60 TensorMechanicsPlasticJ2::dyieldFunction_dintnl(const RankTwoTensor & /*stress*/, Real intnl) const
61 {
62  return -dyieldStrength(intnl);
63 }
64 
65 RankTwoTensor
66 TensorMechanicsPlasticJ2::flowPotential(const RankTwoTensor & stress, Real intnl) const
67 {
68  return dyieldFunction_dstress(stress, intnl);
69 }
70 
71 RankFourTensor
72 TensorMechanicsPlasticJ2::dflowPotential_dstress(const RankTwoTensor & stress, Real /*intnl*/) const
73 {
74  Real sII = stress.secondInvariant();
75  if (sII == 0)
76  return RankFourTensor();
77 
78  RankFourTensor dfp = 0.5 * std::sqrt(3.0 / sII) * stress.d2secondInvariant();
79  Real pre = -0.25 * std::sqrt(3.0) * std::pow(sII, -1.5);
80  RankTwoTensor dII = stress.dsecondInvariant();
81  for (unsigned i = 0; i < 3; ++i)
82  for (unsigned j = 0; j < 3; ++j)
83  for (unsigned k = 0; k < 3; ++k)
84  for (unsigned l = 0; l < 3; ++l)
85  dfp(i, j, k, l) += pre * dII(i, j) * dII(k, l);
86  return dfp;
87 }
88 
89 RankTwoTensor
90 TensorMechanicsPlasticJ2::dflowPotential_dintnl(const RankTwoTensor & /*stress*/,
91  Real /*intnl*/) const
92 {
93  return RankTwoTensor();
94 }
95 
96 Real
98 {
99  return _strength.value(intnl);
100 }
101 
102 Real
104 {
105  return _strength.derivative(intnl);
106 }
107 
108 std::string
110 {
111  return "J2";
112 }
113 
114 bool
115 TensorMechanicsPlasticJ2::returnMap(const RankTwoTensor & trial_stress,
116  Real intnl_old,
117  const RankFourTensor & E_ijkl,
118  Real ep_plastic_tolerance,
119  RankTwoTensor & returned_stress,
120  Real & returned_intnl,
121  std::vector<Real> & dpm,
122  RankTwoTensor & delta_dp,
123  std::vector<Real> & yf,
124  bool & trial_stress_inadmissible) const
125 {
126  if (!(_use_custom_returnMap))
127  return TensorMechanicsPlasticModel::returnMap(trial_stress,
128  intnl_old,
129  E_ijkl,
130  ep_plastic_tolerance,
131  returned_stress,
132  returned_intnl,
133  dpm,
134  delta_dp,
135  yf,
136  trial_stress_inadmissible);
137 
138  yf.resize(1);
139 
140  Real yf_orig = yieldFunction(trial_stress, intnl_old);
141 
142  yf[0] = yf_orig;
143 
144  if (yf_orig < _f_tol)
145  {
146  // the trial_stress is admissible
147  trial_stress_inadmissible = false;
148  return true;
149  }
150 
151  trial_stress_inadmissible = true;
152  Real mu = E_ijkl(0, 1, 0, 1);
153 
154  // Perform a Newton-Raphson to find dpm when
155  // residual = 3*mu*dpm - trial_equivalent_stress + yieldStrength(intnl_old + dpm) = 0
156  Real trial_equivalent_stress = yf_orig + yieldStrength(intnl_old);
157  Real residual;
158  Real jac;
159  dpm[0] = 0;
160  unsigned int iter = 0;
161  do
162  {
163  residual = 3.0 * mu * dpm[0] - trial_equivalent_stress + yieldStrength(intnl_old + dpm[0]);
164  jac = 3.0 * mu + dyieldStrength(intnl_old + dpm[0]);
165  dpm[0] += -residual / jac;
166  if (iter > _max_iters) // not converging
167  return false;
168  iter++;
169  } while (residual * residual > _f_tol * _f_tol);
170 
171  // set the returned values
172  yf[0] = 0;
173  returned_intnl = intnl_old + dpm[0];
174  RankTwoTensor nn = 1.5 * trial_stress.deviatoric() /
175  trial_equivalent_stress; // = dyieldFunction_dstress(trial_stress, intnl_old) =
176  // the normal to the yield surface, at the trial
177  // stress
178  returned_stress = 2.0 / 3.0 * nn * yieldStrength(returned_intnl);
179  returned_stress.addIa(1.0 / 3.0 * trial_stress.trace());
180  delta_dp = nn * dpm[0];
181 
182  return true;
183 }
184 
185 RankFourTensor
186 TensorMechanicsPlasticJ2::consistentTangentOperator(const RankTwoTensor & trial_stress,
187  Real intnl_old,
188  const RankTwoTensor & stress,
189  Real intnl,
190  const RankFourTensor & E_ijkl,
191  const std::vector<Real> & cumulative_pm) const
192 {
193  if (!_use_custom_cto)
195  trial_stress, intnl_old, stress, intnl, E_ijkl, cumulative_pm);
196 
197  Real mu = E_ijkl(0, 1, 0, 1);
198 
199  Real h = 3 * mu + dyieldStrength(intnl);
200  RankTwoTensor sij = stress.deviatoric();
201  Real sII = stress.secondInvariant();
202  Real equivalent_stress = std::sqrt(3.0 * sII);
203  Real zeta = cumulative_pm[0] / (1.0 + 3.0 * mu * cumulative_pm[0] / equivalent_stress);
204 
205  return E_ijkl - 3.0 * mu * mu / sII / h * sij.outerProduct(sij) -
206  4.0 * mu * mu * zeta * dflowPotential_dstress(stress, intnl);
207 }
208 
209 bool
211 {
212  return _use_custom_returnMap;
213 }
214 
215 bool
217 {
218  return _use_custom_cto;
219 }
virtual RankFourTensor consistentTangentOperator(const RankTwoTensor &trial_stress, Real intnl_old, const RankTwoTensor &stress, Real intnl, const RankFourTensor &E_ijkl, const std::vector< Real > &cumulative_pm) const
Calculates a custom consistent tangent operator.
virtual bool returnMap(const RankTwoTensor &trial_stress, Real intnl_old, const RankFourTensor &E_ijkl, Real ep_plastic_tolerance, RankTwoTensor &returned_stress, Real &returned_intnl, std::vector< Real > &dpm, RankTwoTensor &delta_dp, std::vector< Real > &yf, bool &trial_stress_inadmissible) const override
Performs a custom return-map.
virtual RankFourTensor consistentTangentOperator(const RankTwoTensor &trial_stress, Real intnl_old, const RankTwoTensor &stress, Real intnl, const RankFourTensor &E_ijkl, const std::vector< Real > &cumulative_pm) const override
Calculates a custom consistent tangent operator.
virtual RankTwoTensor flowPotential(const RankTwoTensor &stress, Real intnl) const override
The flow potential.
InputParameters validParams< TensorMechanicsPlasticJ2 >()
virtual Real derivative(Real intnl) const
virtual bool returnMap(const RankTwoTensor &trial_stress, Real intnl_old, const RankFourTensor &E_ijkl, Real ep_plastic_tolerance, RankTwoTensor &returned_stress, Real &returned_intnl, std::vector< Real > &dpm, RankTwoTensor &delta_dp, std::vector< Real > &yf, bool &trial_stress_inadmissible) const
Performs a custom return-map.
Real dyieldFunction_dintnl(const RankTwoTensor &stress, Real intnl) const override
The derivative of yield function with respect to the internal parameter.
virtual bool useCustomReturnMap() const override
Returns false. You will want to override this in your derived class if you write a custom returnMap f...
TensorMechanicsPlasticJ2(const InputParameters &parameters)
const unsigned _max_iters
max iters for custom return map loop
virtual Real dyieldStrength(Real intnl) const
d(yieldStrength)/d(intnl)
virtual Real yieldStrength(Real intnl) const
YieldStrength.
virtual RankTwoTensor dyieldFunction_dstress(const RankTwoTensor &stress, Real intnl) const override
The derivative of yield function with respect to stress.
RankTwoTensor dflowPotential_dintnl(const RankTwoTensor &stress, Real intnl) const override
The derivative of the flow potential with respect to the internal parameter.
const bool _use_custom_cto
Whether to use the custom consistent tangent operator calculation.
virtual std::string modelName() const override
InputParameters validParams< TensorMechanicsPlasticModel >()
virtual Real value(Real intnl) const
ExpressionBuilder::EBTerm pow(const ExpressionBuilder::EBTerm &left, T exponent)
const bool _use_custom_returnMap
Whether to use the custom return-map algorithm.
const Real _f_tol
Tolerance on yield function.
virtual Real yieldFunction(const RankTwoTensor &stress, Real intnl) const override
The following functions are what you should override when building single-plasticity models...
const TensorMechanicsHardeningModel & _strength
yield strength, from user input
Plastic Model base class The virtual functions written below must be over-ridden in derived classes t...
virtual RankFourTensor dflowPotential_dstress(const RankTwoTensor &stress, Real intnl) const override
The derivative of the flow potential with respect to stress.
virtual bool useCustomCTO() const override
Returns false. You will want to override this in your derived class if you write a custom consistent ...