www.mooseframework.org
TwoParameterPlasticityStressUpdate.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 #include "Conversion.h" // for stringify
10 #include "libmesh/utility.h" // for Utility::pow
11 
12 template <>
13 InputParameters
15 {
16  InputParameters params = validParams<MultiParameterPlasticityStressUpdate>();
17  params.addClassDescription("Return-map and Jacobian algorithms for (P, Q) plastic models");
18  return params;
19 }
20 
22  const InputParameters & parameters, unsigned num_yf, unsigned num_intnl)
23  : MultiParameterPlasticityStressUpdate(parameters, _num_pq, num_yf, num_intnl),
24  _p_trial(0.0),
25  _q_trial(0.0),
26  _Epp(0.0),
27  _Eqq(0.0),
28  _dgaE_dpt(0.0),
29  _dgaE_dqt(0.0),
30  _dp_dpt(0.0),
31  _dq_dpt(0.0),
32  _dp_dqt(0.0),
33  _dq_dqt(0.0)
34 {
35 }
36 
37 void
38 TwoParameterPlasticityStressUpdate::yieldFunctionValuesV(const std::vector<Real> & stress_params,
39  const std::vector<Real> & intnl,
40  std::vector<Real> & yf) const
41 {
42  const Real p = stress_params[0];
43  const Real q = stress_params[1];
44  yieldFunctionValues(p, q, intnl, yf);
45 }
46 
47 void
49 {
50  setEppEqq(Eijkl, _Epp, _Eqq);
51  _Eij[0][0] = _Epp;
52  _Eij[1][0] = _Eij[0][1] = 0.0;
53  _Eij[1][1] = _Eqq;
54  _En = _Epp;
55  _Cij[0][0] = 1.0 / _Epp;
56  _Cij[1][0] = _Cij[0][1] = 0.0;
57  _Cij[1][1] = 1.0 / _Eqq;
58 }
59 
60 void
61 TwoParameterPlasticityStressUpdate::preReturnMapV(const std::vector<Real> & trial_stress_params,
62  const RankTwoTensor & stress_trial,
63  const std::vector<Real> & intnl_old,
64  const std::vector<Real> & yf,
65  const RankFourTensor & Eijkl)
66 {
67  const Real p_trial = trial_stress_params[0];
68  const Real q_trial = trial_stress_params[1];
69  preReturnMap(p_trial, q_trial, stress_trial, intnl_old, yf, Eijkl);
70 }
71 
72 void
74  Real /*q_trial*/,
75  const RankTwoTensor & /*stress_trial*/,
76  const std::vector<Real> & /*intnl_old*/,
77  const std::vector<Real> & /*yf*/,
78  const RankFourTensor & /*Eijkl*/)
79 {
80  return;
81 }
82 
83 void
85  Real q_trial,
86  const std::vector<Real> & intnl_old,
87  Real & p,
88  Real & q,
89  Real & gaE,
90  std::vector<Real> & intnl) const
91 {
92  p = p_trial;
93  q = q_trial;
94  gaE = 0.0;
95  std::copy(intnl_old.begin(), intnl_old.end(), intnl.begin());
96 }
97 
98 void
99 TwoParameterPlasticityStressUpdate::computeAllQV(const std::vector<Real> & stress_params,
100  const std::vector<Real> & intnl,
101  std::vector<yieldAndFlow> & all_q) const
102 {
103  const Real p = stress_params[0];
104  const Real q = stress_params[1];
105  computeAllQ(p, q, intnl, all_q);
106 }
107 
108 void
109 TwoParameterPlasticityStressUpdate::initializeVarsV(const std::vector<Real> & trial_stress_params,
110  const std::vector<Real> & intnl_old,
111  std::vector<Real> & stress_params,
112  Real & gaE,
113  std::vector<Real> & intnl) const
114 {
115  const Real p_trial = trial_stress_params[0];
116  const Real q_trial = trial_stress_params[1];
117  Real p;
118  Real q;
119  initializeVars(p_trial, q_trial, intnl_old, p, q, gaE, intnl);
120  stress_params[0] = p;
121  stress_params[1] = q;
122 }
123 
124 void
125 TwoParameterPlasticityStressUpdate::setIntnlValuesV(const std::vector<Real> & trial_stress_params,
126  const std::vector<Real> & current_stress_params,
127  const std::vector<Real> & intnl_old,
128  std::vector<Real> & intnl) const
129 {
130  const Real p_trial = trial_stress_params[0];
131  const Real q_trial = trial_stress_params[1];
132  const Real p = current_stress_params[0];
133  const Real q = current_stress_params[1];
134  setIntnlValues(p_trial, q_trial, p, q, intnl_old, intnl);
135 }
136 
137 void
139  const std::vector<Real> & trial_stress_params,
140  const std::vector<Real> & current_stress_params,
141  const std::vector<Real> & intnl_old,
142  std::vector<std::vector<Real>> & dintnl) const
143 {
144  const Real p_trial = trial_stress_params[0];
145  const Real q_trial = trial_stress_params[1];
146  const Real p = current_stress_params[0];
147  const Real q = current_stress_params[1];
148  setIntnlDerivatives(p_trial, q_trial, p, q, intnl_old, dintnl);
149 }
150 
151 void
153  std::vector<Real> & stress_params) const
154 {
155  Real p;
156  Real q;
157  computePQ(stress, p, q);
158  stress_params[0] = p;
159  stress_params[1] = q;
160 }
161 
162 void
164  const RankTwoTensor & stress_trial,
165  const std::vector<Real> & trial_stress_params,
166  const RankTwoTensor & stress,
167  const std::vector<Real> & stress_params,
168  Real gaE,
169  const yieldAndFlow & smoothed_q,
170  const RankFourTensor & Eijkl,
171  bool compute_full_tangent_operator,
172  const std::vector<std::vector<Real>> & dvar_dtrial,
173  RankFourTensor & cto)
174 {
175  const Real p_trial = trial_stress_params[0];
176  const Real q_trial = trial_stress_params[1];
177  const Real p = stress_params[0];
178  const Real q = stress_params[1];
179  _dp_dpt = dvar_dtrial[0][0];
180  _dp_dqt = dvar_dtrial[0][1];
181  _dq_dpt = dvar_dtrial[1][0];
182  _dq_dqt = dvar_dtrial[1][1];
183  _dgaE_dpt = dvar_dtrial[2][0];
184  _dgaE_dqt = dvar_dtrial[2][1];
185  consistentTangentOperator(stress_trial,
186  p_trial,
187  q_trial,
188  stress,
189  p,
190  q,
191  gaE,
192  smoothed_q,
193  Eijkl,
194  compute_full_tangent_operator,
195  cto);
196 }
197 
198 void
200  const RankTwoTensor & stress_trial,
201  Real /*p_trial*/,
202  Real /*q_trial*/,
203  const RankTwoTensor & stress,
204  Real /*p*/,
205  Real /*q*/,
206  Real gaE,
207  const yieldAndFlow & smoothed_q,
208  const RankFourTensor & elasticity_tensor,
209  bool compute_full_tangent_operator,
210  RankFourTensor & cto) const
211 {
212  cto = elasticity_tensor;
213  if (!compute_full_tangent_operator)
214  return;
215 
216  const RankTwoTensor dpdsig = dpdstress(stress);
217  const RankTwoTensor dpdsig_trial = dpdstress(stress_trial);
218  const RankTwoTensor dqdsig = dqdstress(stress);
219  const RankTwoTensor dqdsig_trial = dqdstress(stress_trial);
220 
221  const RankTwoTensor s1 = elasticity_tensor * ((1.0 / _Epp) * (1.0 - _dp_dpt) * dpdsig +
222  (1.0 / _Eqq) * (-_dq_dpt) * dqdsig);
223  const RankTwoTensor s2 = elasticity_tensor * ((1.0 / _Epp) * (-_dp_dqt) * dpdsig +
224  (1.0 / _Eqq) * (1.0 - _dq_dqt) * dqdsig);
225  const RankTwoTensor t1 = elasticity_tensor * dpdsig_trial;
226  const RankTwoTensor t2 = elasticity_tensor * dqdsig_trial;
227 
228  for (unsigned i = 0; i < _tensor_dimensionality; ++i)
229  for (unsigned j = 0; j < _tensor_dimensionality; ++j)
230  for (unsigned k = 0; k < _tensor_dimensionality; ++k)
231  for (unsigned l = 0; l < _tensor_dimensionality; ++l)
232  cto(i, j, k, l) -= s1(i, j) * t1(k, l) + s2(i, j) * t2(k, l);
233 
234  const RankFourTensor d2pdsig2 = d2pdstress2(stress);
235  const RankFourTensor d2qdsig2 = d2qdstress2(stress);
236 
237  const RankFourTensor Tijab = elasticity_tensor * (gaE / _Epp) *
238  (smoothed_q.dg[0] * d2pdsig2 + smoothed_q.dg[1] * d2qdsig2);
239 
240  RankFourTensor inv = RankFourTensor(RankFourTensor::initIdentityFour) + Tijab;
241  try
242  {
243  inv = inv.transposeMajor().invSymm();
244  }
245  catch (const MooseException & e)
246  {
247  // Cannot form the inverse, so probably at some degenerate place in stress space.
248  // Just return with the "best estimate" of the cto.
249  mooseWarning("TwoParameterPlasticityStressUpdate: Cannot invert 1+T in consistent tangent "
250  "operator computation at quadpoint ",
251  _qp,
252  " of element ",
253  _current_elem->id());
254  return;
255  }
256 
257  cto = (cto.transposeMajor() * inv).transposeMajor();
258 }
259 
260 void
262  const std::vector<Real> & stress_params,
263  Real gaE,
264  const std::vector<Real> & intnl,
265  const yieldAndFlow & smoothed_q,
266  const RankFourTensor & Eijkl,
267  RankTwoTensor & stress) const
268 {
269  const Real p_ok = stress_params[0];
270  const Real q_ok = stress_params[1];
271  setStressAfterReturn(stress_trial, p_ok, q_ok, gaE, intnl, smoothed_q, Eijkl, stress);
272 }
273 
274 void
276  const RankTwoTensor & /*stress_trial*/,
277  Real gaE,
278  const yieldAndFlow & smoothed_q,
279  const RankFourTensor & /*elasticity_tensor*/,
280  const RankTwoTensor & returned_stress,
281  RankTwoTensor & inelastic_strain_increment) const
282 {
283  inelastic_strain_increment = (gaE / _Epp) * (smoothed_q.dg[0] * dpdstress(returned_stress) +
284  smoothed_q.dg[1] * dqdstress(returned_stress));
285 }
286 
287 std::vector<RankTwoTensor>
289 {
290  std::vector<RankTwoTensor> dsp(_num_pq, RankTwoTensor());
291  dsp[0] = dpdstress(stress);
292  dsp[1] = dqdstress(stress);
293  return dsp;
294 }
295 
296 std::vector<RankFourTensor>
298 {
299  std::vector<RankFourTensor> d2(_num_pq, RankFourTensor());
300  d2[0] = d2pdstress2(stress);
301  d2[1] = d2qdstress2(stress);
302  return d2;
303 }
virtual void setIntnlValues(Real p_trial, Real q_trial, Real p, Real q, const std::vector< Real > &intnl_old, std::vector< Real > &intnl) const =0
Sets the internal parameters based on the trial values of p and q, their current values, and the old values of the internal parameters.
static constexpr int _num_pq
Number of variables = 2 = (p, q)
InputParameters validParams< TwoParameterPlasticityStressUpdate >()
InputParameters validParams< MultiParameterPlasticityStressUpdate >()
virtual void consistentTangentOperator(const RankTwoTensor &stress_trial, Real p_trial, Real q_trial, const RankTwoTensor &stress, Real p, Real q, Real gaE, const yieldAndFlow &smoothed_q, const RankFourTensor &Eijkl, bool compute_full_tangent_operator, RankFourTensor &cto) const
Calculates the consistent tangent operator.
Real _dp_dpt
derivative of Variable with respect to trial variable (used in consistent-tangent-operator calculatio...
virtual void initializeVars(Real p_trial, Real q_trial, const std::vector< Real > &intnl_old, Real &p, Real &q, Real &gaE, std::vector< Real > &intnl) const
Sets (p, q, gaE, intnl) at "good guesses" of the solution to the Return-Map algorithm.
Real _dq_dpt
derivative of Variable with respect to trial variable (used in consistent-tangent-operator calculatio...
virtual void setEffectiveElasticity(const RankFourTensor &Eijkl) override
Sets _Eij and _En and _Cij.
Struct designed to hold info about a single yield function and its derivatives, as well as the flow d...
virtual void setIntnlDerivativesV(const std::vector< Real > &trial_stress_params, const std::vector< Real > &current_stress_params, const std::vector< Real > &intnl, std::vector< std::vector< Real >> &dintnl) const override
Sets the derivatives of internal parameters, based on the trial values of stress_params, their current values, and the current values of the internal parameters.
std::vector< std::vector< Real > > _Cij
_Cij[i, j] * _Eij[j, k] = 1 iff j == k
Real _dgaE_dqt
derivative of Variable with respect to trial variable (used in consistent-tangent-operator calculatio...
void setStressAfterReturnV(const RankTwoTensor &stress_trial, const std::vector< Real > &stress_params, Real gaE, const std::vector< Real > &intnl, const yieldAndFlow &smoothed_q, const RankFourTensor &Eijkl, RankTwoTensor &stress) const override
Sets stress from the admissible parameters.
virtual void setEppEqq(const RankFourTensor &Eijkl, Real &Epp, Real &Eqq) const =0
Set Epp and Eqq based on the elasticity tensor Derived classes must override this.
Real _dgaE_dpt
derivative of Variable with respect to trial variable (used in consistent-tangent-operator calculatio...
virtual void computePQ(const RankTwoTensor &stress, Real &p, Real &q) const =0
Computes p and q, given stress.
TwoParameterPlasticityStressUpdate(const InputParameters &parameters, unsigned num_yf, unsigned num_intnl)
std::vector< std::vector< Real > > _Eij
E[i, j] in the system of equations to be solved.
Real _dq_dqt
derivative of Variable with respect to trial variable (used in consistent-tangent-operator calculatio...
virtual std::vector< RankTwoTensor > dstress_param_dstress(const RankTwoTensor &stress) const override
d(stress_param[i])/d(stress) at given stress
virtual void computeAllQ(Real p, Real q, const std::vector< Real > &intnl, std::vector< yieldAndFlow > &all_q) const =0
Completely fills all_q with correct values.
virtual void setIntnlDerivatives(Real p_trial, Real q_trial, Real p, Real q, const std::vector< Real > &intnl, std::vector< std::vector< Real >> &dintnl) const =0
Sets the derivatives of internal parameters, based on the trial values of p and q, their current values, and the old values of the internal parameters.
virtual void setStressAfterReturn(const RankTwoTensor &stress_trial, Real p_ok, Real q_ok, Real gaE, const std::vector< Real > &intnl, const yieldAndFlow &smoothed_q, const RankFourTensor &Eijkl, RankTwoTensor &stress) const =0
Sets stress from the admissible parameters.
void consistentTangentOperatorV(const RankTwoTensor &stress_trial, const std::vector< Real > &trial_stress_params, const RankTwoTensor &stress, const std::vector< Real > &stress_params, Real gaE, const yieldAndFlow &smoothed_q, const RankFourTensor &Eijkl, bool compute_full_tangent_operator, const std::vector< std::vector< Real >> &dvar_dtrial, RankFourTensor &cto) override
Calculates the consistent tangent operator.
void setIntnlValuesV(const std::vector< Real > &trial_stress_params, const std::vector< Real > &current_stress_params, const std::vector< Real > &intnl_old, std::vector< Real > &intnl) const override
Sets the internal parameters based on the trial values of stress_params, their current values...
Real _Epp
elasticity tensor in p direction
void setInelasticStrainIncrementAfterReturn(const RankTwoTensor &stress_trial, Real gaE, const yieldAndFlow &smoothed_q, const RankFourTensor &elasticity_tensor, const RankTwoTensor &returned_stress, RankTwoTensor &inelastic_strain_increment) const override
Sets inelastic strain increment from the returned configuration This is called after the return-map p...
virtual void computeStressParams(const RankTwoTensor &stress, std::vector< Real > &stress_params) const override
Computes stress_params, given stress.
virtual void yieldFunctionValues(Real p, Real q, const std::vector< Real > &intnl, std::vector< Real > &yf) const =0
Computes the values of the yield functions, given p, q and intnl parameters.
virtual RankFourTensor d2qdstress2(const RankTwoTensor &stress) const =0
d2(q)/d(stress)/d(stress) Derived classes must override this
void initializeVarsV(const std::vector< Real > &trial_stress_params, const std::vector< Real > &intnl_old, std::vector< Real > &stress_params, Real &gaE, std::vector< Real > &intnl) const override
Sets (stress_params, intnl) at "good guesses" of the solution to the Return-Map algorithm.
virtual std::vector< RankFourTensor > d2stress_param_dstress(const RankTwoTensor &stress) const override
d2(stress_param[i])/d(stress)/d(stress) at given stress
Real _Eqq
elasticity tensor in q direction
MultiParameterPlasticityStressUpdate performs the return-map algorithm and associated stress updates ...
void preReturnMapV(const std::vector< Real > &trial_stress_params, const RankTwoTensor &stress_trial, const std::vector< Real > &intnl_old, const std::vector< Real > &yf, const RankFourTensor &Eijkl) override
Derived classes may employ this function to record stuff or do other computations prior to the return...
virtual void preReturnMap(Real p_trial, Real q_trial, const RankTwoTensor &stress_trial, const std::vector< Real > &intnl_old, const std::vector< Real > &yf, const RankFourTensor &Eijkl)
Derived classes may employ this function to record stuff or do other computations prior to the return...
void yieldFunctionValuesV(const std::vector< Real > &stress_params, const std::vector< Real > &intnl, std::vector< Real > &yf) const override
Computes the values of the yield functions, given stress_params and intnl parameters.
void computeAllQV(const std::vector< Real > &stress_params, const std::vector< Real > &intnl, std::vector< yieldAndFlow > &all_q) const override
Completely fills all_q with correct values.
virtual RankTwoTensor dpdstress(const RankTwoTensor &stress) const =0
d(p)/d(stress) Derived classes must override this
Real _dp_dqt
derivative of Variable with respect to trial variable (used in consistent-tangent-operator calculatio...
static constexpr unsigned _tensor_dimensionality
Internal dimensionality of tensors (currently this is 3 throughout tensor_mechanics) ...
virtual RankFourTensor d2pdstress2(const RankTwoTensor &stress) const =0
d2(p)/d(stress)/d(stress) Derived classes must override this
virtual RankTwoTensor dqdstress(const RankTwoTensor &stress) const =0
d(q)/d(stress) Derived classes must override this