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