Line data Source code
1 : //* This file is part of the MOOSE framework
2 : //* https://mooseframework.inl.gov
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 : #include "TwoParameterPlasticityStressUpdate.h"
11 :
12 : #include "Conversion.h" // for stringify
13 : #include "libmesh/utility.h" // for Utility::pow
14 :
15 : InputParameters
16 1780 : TwoParameterPlasticityStressUpdate::validParams()
17 : {
18 1780 : InputParameters params = MultiParameterPlasticityStressUpdate::validParams();
19 1780 : params.addClassDescription("Return-map and Jacobian algorithms for (P, Q) plastic models");
20 1780 : return params;
21 0 : }
22 :
23 1338 : TwoParameterPlasticityStressUpdate::TwoParameterPlasticityStressUpdate(
24 1338 : const InputParameters & parameters, unsigned num_yf, unsigned num_intnl)
25 : : MultiParameterPlasticityStressUpdate(parameters, _num_pq, num_yf, num_intnl),
26 1338 : _p_trial(0.0),
27 1338 : _q_trial(0.0),
28 1338 : _Epp(0.0),
29 1338 : _Eqq(0.0),
30 1338 : _dgaE_dpt(0.0),
31 1338 : _dgaE_dqt(0.0),
32 1338 : _dp_dpt(0.0),
33 1338 : _dq_dpt(0.0),
34 1338 : _dp_dqt(0.0),
35 1338 : _dq_dqt(0.0)
36 : {
37 1338 : }
38 :
39 : void
40 1804132 : TwoParameterPlasticityStressUpdate::yieldFunctionValuesV(const std::vector<Real> & stress_params,
41 : const std::vector<Real> & intnl,
42 : std::vector<Real> & yf) const
43 : {
44 1804132 : const Real p = stress_params[0];
45 1804132 : const Real q = stress_params[1];
46 1804132 : yieldFunctionValues(p, q, intnl, yf);
47 1804132 : }
48 :
49 : void
50 507396 : TwoParameterPlasticityStressUpdate::setEffectiveElasticity(const RankFourTensor & Eijkl)
51 : {
52 507396 : setEppEqq(Eijkl, _Epp, _Eqq);
53 507396 : _Eij[0][0] = _Epp;
54 507396 : _Eij[1][0] = _Eij[0][1] = 0.0;
55 507396 : _Eij[1][1] = _Eqq;
56 507396 : _En = _Epp;
57 507396 : _Cij[0][0] = 1.0 / _Epp;
58 507396 : _Cij[1][0] = _Cij[0][1] = 0.0;
59 507396 : _Cij[1][1] = 1.0 / _Eqq;
60 507396 : }
61 :
62 : void
63 507396 : TwoParameterPlasticityStressUpdate::preReturnMapV(const std::vector<Real> & trial_stress_params,
64 : const RankTwoTensor & stress_trial,
65 : const std::vector<Real> & intnl_old,
66 : const std::vector<Real> & yf,
67 : const RankFourTensor & Eijkl)
68 : {
69 507396 : const Real p_trial = trial_stress_params[0];
70 507396 : const Real q_trial = trial_stress_params[1];
71 507396 : preReturnMap(p_trial, q_trial, stress_trial, intnl_old, yf, Eijkl);
72 507396 : }
73 :
74 : void
75 0 : TwoParameterPlasticityStressUpdate::preReturnMap(Real /*p_trial*/,
76 : Real /*q_trial*/,
77 : const RankTwoTensor & /*stress_trial*/,
78 : const std::vector<Real> & /*intnl_old*/,
79 : const std::vector<Real> & /*yf*/,
80 : const RankFourTensor & /*Eijkl*/)
81 : {
82 0 : return;
83 : }
84 :
85 : void
86 0 : TwoParameterPlasticityStressUpdate::initializeVars(Real p_trial,
87 : Real q_trial,
88 : const std::vector<Real> & intnl_old,
89 : Real & p,
90 : Real & q,
91 : Real & gaE,
92 : std::vector<Real> & intnl) const
93 : {
94 0 : p = p_trial;
95 0 : q = q_trial;
96 0 : gaE = 0.0;
97 0 : std::copy(intnl_old.begin(), intnl_old.end(), intnl.begin());
98 0 : }
99 :
100 : void
101 1636862 : TwoParameterPlasticityStressUpdate::computeAllQV(const std::vector<Real> & stress_params,
102 : const std::vector<Real> & intnl,
103 : std::vector<yieldAndFlow> & all_q) const
104 : {
105 1636862 : const Real p = stress_params[0];
106 1636862 : const Real q = stress_params[1];
107 1636862 : computeAllQ(p, q, intnl, all_q);
108 1636862 : }
109 :
110 : void
111 507396 : TwoParameterPlasticityStressUpdate::initializeVarsV(const std::vector<Real> & trial_stress_params,
112 : const std::vector<Real> & intnl_old,
113 : std::vector<Real> & stress_params,
114 : Real & gaE,
115 : std::vector<Real> & intnl) const
116 : {
117 507396 : const Real p_trial = trial_stress_params[0];
118 507396 : const Real q_trial = trial_stress_params[1];
119 : Real p;
120 : Real q;
121 507396 : initializeVars(p_trial, q_trial, intnl_old, p, q, gaE, intnl);
122 507396 : stress_params[0] = p;
123 507396 : stress_params[1] = q;
124 507396 : }
125 :
126 : void
127 1636778 : TwoParameterPlasticityStressUpdate::setIntnlValuesV(const std::vector<Real> & trial_stress_params,
128 : const std::vector<Real> & current_stress_params,
129 : const std::vector<Real> & intnl_old,
130 : std::vector<Real> & intnl) const
131 : {
132 1636778 : const Real p_trial = trial_stress_params[0];
133 1636778 : const Real q_trial = trial_stress_params[1];
134 1636778 : const Real p = current_stress_params[0];
135 1636778 : const Real q = current_stress_params[1];
136 1636778 : setIntnlValues(p_trial, q_trial, p, q, intnl_old, intnl);
137 1636778 : }
138 :
139 : void
140 1204678 : TwoParameterPlasticityStressUpdate::setIntnlDerivativesV(
141 : const std::vector<Real> & trial_stress_params,
142 : const std::vector<Real> & current_stress_params,
143 : const std::vector<Real> & intnl_old,
144 : std::vector<std::vector<Real>> & dintnl) const
145 : {
146 1204678 : const Real p_trial = trial_stress_params[0];
147 1204678 : const Real q_trial = trial_stress_params[1];
148 1204678 : const Real p = current_stress_params[0];
149 1204678 : const Real q = current_stress_params[1];
150 1204678 : setIntnlDerivatives(p_trial, q_trial, p, q, intnl_old, dintnl);
151 1204678 : }
152 :
153 : void
154 1630916 : TwoParameterPlasticityStressUpdate::computeStressParams(const RankTwoTensor & stress,
155 : std::vector<Real> & stress_params) const
156 : {
157 : Real p;
158 : Real q;
159 1630916 : computePQ(stress, p, q);
160 1630916 : stress_params[0] = p;
161 1630916 : stress_params[1] = q;
162 1630916 : }
163 :
164 : void
165 96544 : TwoParameterPlasticityStressUpdate::consistentTangentOperatorV(
166 : const RankTwoTensor & stress_trial,
167 : const std::vector<Real> & trial_stress_params,
168 : const RankTwoTensor & stress,
169 : const std::vector<Real> & stress_params,
170 : Real gaE,
171 : const yieldAndFlow & smoothed_q,
172 : const RankFourTensor & Eijkl,
173 : bool compute_full_tangent_operator,
174 : const std::vector<std::vector<Real>> & dvar_dtrial,
175 : RankFourTensor & cto)
176 : {
177 96544 : const Real p_trial = trial_stress_params[0];
178 96544 : const Real q_trial = trial_stress_params[1];
179 96544 : const Real p = stress_params[0];
180 96544 : const Real q = stress_params[1];
181 96544 : _dp_dpt = dvar_dtrial[0][0];
182 96544 : _dp_dqt = dvar_dtrial[0][1];
183 96544 : _dq_dpt = dvar_dtrial[1][0];
184 96544 : _dq_dqt = dvar_dtrial[1][1];
185 96544 : _dgaE_dpt = dvar_dtrial[2][0];
186 96544 : _dgaE_dqt = dvar_dtrial[2][1];
187 96544 : consistentTangentOperator(stress_trial,
188 : p_trial,
189 : q_trial,
190 : stress,
191 : p,
192 : q,
193 : gaE,
194 : smoothed_q,
195 : Eijkl,
196 : compute_full_tangent_operator,
197 : cto);
198 96544 : }
199 :
200 : void
201 0 : TwoParameterPlasticityStressUpdate::consistentTangentOperator(
202 : const RankTwoTensor & stress_trial,
203 : Real /*p_trial*/,
204 : Real /*q_trial*/,
205 : const RankTwoTensor & stress,
206 : Real /*p*/,
207 : Real /*q*/,
208 : Real gaE,
209 : const yieldAndFlow & smoothed_q,
210 : const RankFourTensor & elasticity_tensor,
211 : bool compute_full_tangent_operator,
212 : RankFourTensor & cto) const
213 : {
214 0 : cto = elasticity_tensor;
215 0 : if (!compute_full_tangent_operator)
216 0 : return;
217 :
218 0 : const RankTwoTensor dpdsig = dpdstress(stress);
219 0 : const RankTwoTensor dpdsig_trial = dpdstress(stress_trial);
220 0 : const RankTwoTensor dqdsig = dqdstress(stress);
221 0 : const RankTwoTensor dqdsig_trial = dqdstress(stress_trial);
222 :
223 0 : const RankTwoTensor s1 = elasticity_tensor * ((1.0 / _Epp) * (1.0 - _dp_dpt) * dpdsig +
224 0 : (1.0 / _Eqq) * (-_dq_dpt) * dqdsig);
225 0 : const RankTwoTensor s2 = elasticity_tensor * ((1.0 / _Epp) * (-_dp_dqt) * dpdsig +
226 0 : (1.0 / _Eqq) * (1.0 - _dq_dqt) * dqdsig);
227 0 : const RankTwoTensor t1 = elasticity_tensor * dpdsig_trial;
228 0 : const RankTwoTensor t2 = elasticity_tensor * dqdsig_trial;
229 :
230 0 : for (unsigned i = 0; i < _tensor_dimensionality; ++i)
231 0 : for (unsigned j = 0; j < _tensor_dimensionality; ++j)
232 0 : for (unsigned k = 0; k < _tensor_dimensionality; ++k)
233 0 : for (unsigned l = 0; l < _tensor_dimensionality; ++l)
234 0 : cto(i, j, k, l) -= s1(i, j) * t1(k, l) + s2(i, j) * t2(k, l);
235 :
236 0 : const RankFourTensor d2pdsig2 = d2pdstress2(stress);
237 0 : const RankFourTensor d2qdsig2 = d2qdstress2(stress);
238 :
239 0 : const RankFourTensor Tijab = elasticity_tensor * (gaE / _Epp) *
240 0 : (smoothed_q.dg[0] * d2pdsig2 + smoothed_q.dg[1] * d2qdsig2);
241 :
242 0 : RankFourTensor inv = RankFourTensor(RankFourTensor::initIdentityFour) + Tijab;
243 : try
244 : {
245 0 : inv = inv.transposeMajor().invSymm();
246 : }
247 0 : catch (const MooseException & e)
248 : {
249 : // Cannot form the inverse, so probably at some degenerate place in stress space.
250 : // Just return with the "best estimate" of the cto.
251 0 : mooseWarning("TwoParameterPlasticityStressUpdate: Cannot invert 1+T in consistent tangent "
252 : "operator computation at quadpoint ",
253 0 : _qp,
254 : " of element ",
255 0 : _current_elem->id());
256 : return;
257 0 : }
258 :
259 0 : cto = (cto.transposeMajor() * inv).transposeMajor();
260 : }
261 :
262 : void
263 507312 : TwoParameterPlasticityStressUpdate::setStressAfterReturnV(const RankTwoTensor & stress_trial,
264 : const std::vector<Real> & stress_params,
265 : Real gaE,
266 : const std::vector<Real> & intnl,
267 : const yieldAndFlow & smoothed_q,
268 : const RankFourTensor & Eijkl,
269 : RankTwoTensor & stress) const
270 : {
271 507312 : const Real p_ok = stress_params[0];
272 507312 : const Real q_ok = stress_params[1];
273 507312 : setStressAfterReturn(stress_trial, p_ok, q_ok, gaE, intnl, smoothed_q, Eijkl, stress);
274 507312 : }
275 :
276 : void
277 507312 : TwoParameterPlasticityStressUpdate::setInelasticStrainIncrementAfterReturn(
278 : const RankTwoTensor & /*stress_trial*/,
279 : Real gaE,
280 : const yieldAndFlow & smoothed_q,
281 : const RankFourTensor & /*elasticity_tensor*/,
282 : const RankTwoTensor & returned_stress,
283 : RankTwoTensor & inelastic_strain_increment) const
284 : {
285 507312 : inelastic_strain_increment = (gaE / _Epp) * (smoothed_q.dg[0] * dpdstress(returned_stress) +
286 507312 : smoothed_q.dg[1] * dqdstress(returned_stress));
287 507312 : }
288 :
289 : std::vector<RankTwoTensor>
290 0 : TwoParameterPlasticityStressUpdate::dstress_param_dstress(const RankTwoTensor & stress) const
291 : {
292 0 : std::vector<RankTwoTensor> dsp(_num_pq, RankTwoTensor());
293 0 : dsp[0] = dpdstress(stress);
294 0 : dsp[1] = dqdstress(stress);
295 0 : return dsp;
296 : }
297 :
298 : std::vector<RankFourTensor>
299 0 : TwoParameterPlasticityStressUpdate::d2stress_param_dstress(const RankTwoTensor & stress) const
300 : {
301 0 : std::vector<RankFourTensor> d2(_num_pq, RankFourTensor());
302 0 : d2[0] = d2pdstress2(stress);
303 0 : d2[1] = d2qdstress2(stress);
304 0 : return d2;
305 : }
|
