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 "HillCreepStressUpdate.h"
11 : #include "ElasticityTensorTools.h"
12 : #include "Function.h"
13 :
14 : registerMooseObject("SolidMechanicsApp", ADHillCreepStressUpdate);
15 : registerMooseObject("SolidMechanicsApp", HillCreepStressUpdate);
16 :
17 : template <bool is_ad>
18 : InputParameters
19 304 : HillCreepStressUpdateTempl<is_ad>::validParams()
20 : {
21 304 : InputParameters params = AnisotropicReturnCreepStressUpdateBaseTempl<is_ad>::validParams();
22 304 : params.addClassDescription(
23 : "This class uses the stress update material in a generalized radial return anisotropic power "
24 : "law creep "
25 : "model. This class can be used in conjunction with other creep and plasticity materials for "
26 : "more complex simulations.");
27 :
28 : // Linear strain hardening parameters
29 608 : params.addCoupledVar("temperature", "Coupled temperature");
30 608 : params.addRequiredParam<Real>("coefficient", "Leading coefficient in power-law equation");
31 608 : params.addRequiredParam<Real>("n_exponent", "Exponent on effective stress in power-law equation");
32 608 : params.addParam<Real>("m_exponent", 0.0, "Exponent on time in power-law equation");
33 608 : params.addRequiredParam<Real>("activation_energy", "Activation energy");
34 608 : params.addParam<Real>("gas_constant", 8.3143, "Universal gas constant");
35 608 : params.addParam<Real>("start_time", 0.0, "Start time (if not zero)");
36 608 : params.addParam<bool>("anisotropic_elasticity", false, "Enable using anisotropic elasticity");
37 608 : params.addParam<FunctionName>(
38 : "creep_prefactor", "Optional function to use as a scalar prefactor on the creep strain.");
39 304 : return params;
40 0 : }
41 :
42 : template <bool is_ad>
43 228 : HillCreepStressUpdateTempl<is_ad>::HillCreepStressUpdateTempl(const InputParameters & parameters)
44 : : AnisotropicReturnCreepStressUpdateBaseTempl<is_ad>(parameters),
45 228 : _has_temp(this->isParamValid("temperature")),
46 228 : _temperature(this->_has_temp ? this->coupledValue("temperature") : this->_zero),
47 456 : _coefficient(this->template getParam<Real>("coefficient")),
48 456 : _n_exponent(this->template getParam<Real>("n_exponent")),
49 456 : _m_exponent(this->template getParam<Real>("m_exponent")),
50 456 : _activation_energy(this->template getParam<Real>("activation_energy")),
51 456 : _gas_constant(this->template getParam<Real>("gas_constant")),
52 456 : _start_time(this->template getParam<Real>("start_time")),
53 228 : _exponential(1.0),
54 228 : _exp_time(1.0),
55 456 : _hill_constants(this->template getMaterialPropertyByName<std::vector<Real>>(this->_base_name +
56 : "hill_constants")),
57 456 : _hill_tensor(this->_use_transformation
58 612 : ? &this->template getMaterialPropertyByName<DenseMatrix<Real>>(
59 : this->_base_name + "hill_tensor")
60 : : nullptr),
61 228 : _C(6, 6),
62 228 : _elasticity_tensor_name(this->_base_name + "elasticity_tensor"),
63 228 : _elasticity_tensor(
64 228 : this->template getGenericMaterialProperty<RankFourTensor, is_ad>(_elasticity_tensor_name)),
65 456 : _anisotropic_elasticity(this->template getParam<bool>("anisotropic_elasticity")),
66 228 : _prefactor_function(
67 684 : this->isParamValid("creep_prefactor") ? &this->getFunction("creep_prefactor") : nullptr)
68 : {
69 228 : if (_start_time < this->_app.getStartTime() && (std::trunc(_m_exponent) != _m_exponent))
70 0 : this->paramError("start_time",
71 : "Start time must be equal to or greater than the Executioner start_time if a "
72 : "non-integer m_exponent is used");
73 228 : }
74 :
75 : template <bool is_ad>
76 : void
77 1814752 : HillCreepStressUpdateTempl<is_ad>::computeStressInitialize(
78 : const GenericDenseVector<is_ad> & /*stress_dev*/,
79 : const GenericDenseVector<is_ad> & /*stress*/,
80 : const GenericRankFourTensor<is_ad> & elasticity_tensor)
81 : {
82 1814752 : if (_has_temp)
83 0 : _exponential = std::exp(-_activation_energy / (_gas_constant * _temperature[_qp]));
84 :
85 2574528 : _two_shear_modulus = 2.0 * ElasticityTensorTools::getIsotropicShearModulus(elasticity_tensor);
86 :
87 1814752 : _exp_time = std::pow(_t - _start_time, _m_exponent);
88 1814752 : }
89 :
90 : template <bool is_ad>
91 : GenericReal<is_ad>
92 1781696 : HillCreepStressUpdateTempl<is_ad>::initialGuess(const GenericDenseVector<is_ad> & /*stress_dev*/)
93 : {
94 1781696 : return 0.0;
95 : }
96 :
97 : template <bool is_ad>
98 : GenericReal<is_ad>
99 9130708 : HillCreepStressUpdateTempl<is_ad>::computeResidual(
100 : const GenericDenseVector<is_ad> & /*effective_trial_stress*/,
101 : const GenericDenseVector<is_ad> & stress_new,
102 : const GenericReal<is_ad> & delta_gamma)
103 : {
104 : GenericReal<is_ad> qsigma_changed;
105 9130708 : GenericRankTwoTensor<is_ad> stress_changed;
106 :
107 : // Get the qsigma_changed and stress_changed (not used)
108 : // delta_gamma (scalar plastic strain) will change the stress (and qsigma as well)
109 9130708 : computeQsigmaChanged(qsigma_changed, stress_new, delta_gamma, stress_changed);
110 :
111 7900584 : GenericReal<is_ad> creep_rate =
112 9130708 : _coefficient * std::pow(qsigma_changed, _n_exponent) * _exponential * _exp_time;
113 :
114 9130708 : if (_prefactor_function)
115 0 : creep_rate *= _prefactor_function->value(_t, _q_point[_qp]);
116 :
117 9130708 : this->_inelastic_strain_rate[_qp] = MetaPhysicL::raw_value(creep_rate);
118 : // Return iteration difference between creep strain and inelastic strain multiplier
119 13081000 : return creep_rate * _dt - delta_gamma;
120 : }
121 :
122 : template <bool is_ad>
123 : Real
124 9130708 : HillCreepStressUpdateTempl<is_ad>::computeReferenceResidual(
125 : const GenericDenseVector<is_ad> & /*effective_trial_stress*/,
126 : const GenericDenseVector<is_ad> & /*stress_new*/,
127 : const GenericReal<is_ad> & /*residual*/,
128 : const GenericReal<is_ad> & /*scalar_effective_inelastic_strain*/)
129 : {
130 9130708 : return 1.0;
131 : }
132 :
133 : template <bool is_ad>
134 : GenericReal<is_ad>
135 7349012 : HillCreepStressUpdateTempl<is_ad>::computeDerivative(
136 : const GenericDenseVector<is_ad> & /*effective_trial_stress*/,
137 : const GenericDenseVector<is_ad> & stress_new,
138 : const GenericReal<is_ad> & delta_gamma)
139 : {
140 : GenericReal<is_ad> qsigma_changed;
141 7349012 : GenericRankTwoTensor<is_ad> stress_changed;
142 :
143 : // Get the qsigma_changed and stress_changed
144 : // delta_gamma (scalar plastic strain) will change the stress (and qsigma as well)
145 7349012 : computeQsigmaChanged(qsigma_changed, stress_new, delta_gamma, stress_changed);
146 :
147 7349012 : ElasticityTensorTools::toMooseVoigtNotation<is_ad>(_C, _elasticity_tensor[_qp]);
148 : const unsigned int dimension = _C.n();
149 7349012 : GenericDenseMatrix<is_ad> d_stress_d_inelasticStrainIncrement(dimension, dimension);
150 :
151 51443084 : for (unsigned int index_i = 0; index_i < dimension; index_i++)
152 308658504 : for (unsigned int index_j = 0; index_j < dimension; index_j++)
153 : {
154 264564432 : if (index_j < 3)
155 132282216 : d_stress_d_inelasticStrainIncrement(index_i, index_j) =
156 132282216 : -1.0 * MetaPhysicL::raw_value(_C(index_i, index_j));
157 : else
158 132282216 : d_stress_d_inelasticStrainIncrement(index_i, index_j) =
159 132282216 : -2.0 * MetaPhysicL::raw_value(_C(index_i, index_j));
160 : }
161 :
162 : // Hill constants, some constraints apply
163 7349012 : const Real & F = _hill_constants[_qp][0];
164 : const Real & G = _hill_constants[_qp][1];
165 : const Real & H = _hill_constants[_qp][2];
166 : const Real & L = _hill_constants[_qp][3];
167 : const Real & M = _hill_constants[_qp][4];
168 : const Real & N = _hill_constants[_qp][5];
169 :
170 7349012 : GenericDenseVector<is_ad> d_qsigma_d_inelasticStrainIncrement(6);
171 51443084 : for (unsigned int index_k = 0; index_k < 6; index_k++)
172 : {
173 19278456 : d_qsigma_d_inelasticStrainIncrement(index_k) =
174 24815616 : F * (stress_changed(1, 1) - stress_changed(2, 2)) *
175 24815616 : (d_stress_d_inelasticStrainIncrement(1, index_k) -
176 24815616 : d_stress_d_inelasticStrainIncrement(2, index_k)) +
177 24815616 : G * (stress_changed(2, 2) - stress_changed(0, 0)) *
178 24815616 : (d_stress_d_inelasticStrainIncrement(2, index_k) -
179 24815616 : d_stress_d_inelasticStrainIncrement(0, index_k)) +
180 24815616 : H * (stress_changed(0, 0) - stress_changed(1, 1)) *
181 24815616 : (d_stress_d_inelasticStrainIncrement(0, index_k) -
182 24815616 : d_stress_d_inelasticStrainIncrement(1, index_k)) +
183 44094072 : 2.0 * L * stress_changed(1, 2) * d_stress_d_inelasticStrainIncrement(4, index_k) +
184 44094072 : 2.0 * M * stress_changed(2, 0) * d_stress_d_inelasticStrainIncrement(5, index_k) +
185 44094072 : 2.0 * N * stress_changed(0, 1) * d_stress_d_inelasticStrainIncrement(3, index_k);
186 44094072 : d_qsigma_d_inelasticStrainIncrement(index_k) /= qsigma_changed;
187 : }
188 :
189 7349012 : GenericDenseVector<is_ad> d_qsigma_d_sigma(6);
190 :
191 7349012 : d_qsigma_d_sigma(0) = (H * (stress_changed(0, 0) - stress_changed(1, 1)) -
192 4135936 : G * (stress_changed(2, 2) - stress_changed(0, 0))) /
193 : qsigma_changed;
194 7349012 : d_qsigma_d_sigma(1) = (F * (stress_changed(1, 1) - stress_changed(2, 2)) -
195 4135936 : H * (stress_changed(0, 0) - stress_changed(1, 1))) /
196 : qsigma_changed;
197 7349012 : d_qsigma_d_sigma(2) = (G * (stress_changed(2, 2) - stress_changed(0, 0)) -
198 4135936 : F * (stress_changed(1, 1) - stress_changed(2, 2))) /
199 : qsigma_changed;
200 10562088 : d_qsigma_d_sigma(3) = 2.0 * N * stress_changed(0, 1) / qsigma_changed;
201 10562088 : d_qsigma_d_sigma(4) = 2.0 * L * stress_changed(1, 2) / qsigma_changed;
202 10562088 : d_qsigma_d_sigma(5) = 2.0 * M * stress_changed(0, 2) / qsigma_changed;
203 :
204 3213076 : GenericReal<is_ad> d_qsigma_d_deltaGamma =
205 4135936 : d_qsigma_d_inelasticStrainIncrement.dot(d_qsigma_d_sigma);
206 :
207 7349012 : GenericReal<is_ad> creep_rate_derivative = _coefficient * d_qsigma_d_deltaGamma * _n_exponent *
208 10562088 : std::pow(qsigma_changed, _n_exponent - 1.0) *
209 7349012 : _exponential * _exp_time;
210 :
211 7349012 : if (_prefactor_function)
212 0 : creep_rate_derivative *= _prefactor_function->value(_t, _q_point[_qp]);
213 :
214 17911100 : return (creep_rate_derivative * _dt - 1.0);
215 7349012 : }
216 :
217 : template <bool is_ad>
218 : void
219 1781696 : HillCreepStressUpdateTempl<is_ad>::computeStrainFinalize(
220 : GenericRankTwoTensor<is_ad> & inelasticStrainIncrement,
221 : const GenericRankTwoTensor<is_ad> & stress,
222 : const GenericDenseVector<is_ad> & stress_dev,
223 : const GenericReal<is_ad> & delta_gamma)
224 : {
225 : GenericReal<is_ad> qsigma_square;
226 1781696 : if (!this->_use_transformation)
227 : {
228 : // Hill constants, some constraints apply
229 522240 : const Real & F = _hill_constants[_qp][0];
230 : const Real & G = _hill_constants[_qp][1];
231 : const Real & H = _hill_constants[_qp][2];
232 : const Real & L = _hill_constants[_qp][3];
233 : const Real & M = _hill_constants[_qp][4];
234 : const Real & N = _hill_constants[_qp][5];
235 :
236 : // Equivalent deviatoric stress function.
237 522240 : qsigma_square = F * (stress(1, 1) - stress(2, 2)) * (stress(1, 1) - stress(2, 2));
238 522240 : qsigma_square += G * (stress(2, 2) - stress(0, 0)) * (stress(2, 2) - stress(0, 0));
239 522240 : qsigma_square += H * (stress(0, 0) - stress(1, 1)) * (stress(0, 0) - stress(1, 1));
240 522240 : qsigma_square += 2 * L * stress(1, 2) * stress(1, 2);
241 522240 : qsigma_square += 2 * M * stress(0, 2) * stress(0, 2);
242 522240 : qsigma_square += 2 * N * stress(0, 1) * stress(0, 1);
243 : }
244 : else
245 : {
246 1259456 : GenericDenseVector<is_ad> Ms(6);
247 1259456 : (*_hill_tensor)[_qp].vector_mult(Ms, stress_dev);
248 1259456 : qsigma_square = Ms.dot(stress_dev);
249 : }
250 :
251 1781696 : if (qsigma_square == 0)
252 : {
253 0 : inelasticStrainIncrement.zero();
254 :
255 0 : AnisotropicReturnCreepStressUpdateBaseTempl<is_ad>::computeStrainFinalize(
256 : inelasticStrainIncrement, stress, stress_dev, delta_gamma);
257 :
258 0 : this->_effective_inelastic_strain[_qp] = this->_effective_inelastic_strain_old[_qp];
259 :
260 0 : return;
261 : }
262 :
263 : // Use Hill-type flow rule to compute the time step inelastic increment.
264 1781696 : GenericReal<is_ad> prefactor = delta_gamma / std::sqrt(qsigma_square);
265 :
266 1781696 : if (!this->_use_transformation)
267 : {
268 : // Hill constants, some constraints apply
269 522240 : const Real & F = _hill_constants[_qp][0];
270 : const Real & G = _hill_constants[_qp][1];
271 : const Real & H = _hill_constants[_qp][2];
272 : const Real & L = _hill_constants[_qp][3];
273 : const Real & M = _hill_constants[_qp][4];
274 : const Real & N = _hill_constants[_qp][5];
275 :
276 522240 : inelasticStrainIncrement(0, 0) =
277 522240 : prefactor * (H * (stress(0, 0) - stress(1, 1)) - G * (stress(2, 2) - stress(0, 0)));
278 522240 : inelasticStrainIncrement(1, 1) =
279 522240 : prefactor * (F * (stress(1, 1) - stress(2, 2)) - H * (stress(0, 0) - stress(1, 1)));
280 522240 : inelasticStrainIncrement(2, 2) =
281 522240 : prefactor * (G * (stress(2, 2) - stress(0, 0)) - F * (stress(1, 1) - stress(2, 2)));
282 :
283 522240 : inelasticStrainIncrement(0, 1) = inelasticStrainIncrement(1, 0) =
284 522240 : prefactor * 2.0 * N * stress(0, 1);
285 522240 : inelasticStrainIncrement(0, 2) = inelasticStrainIncrement(2, 0) =
286 522240 : prefactor * 2.0 * M * stress(0, 2);
287 522240 : inelasticStrainIncrement(1, 2) = inelasticStrainIncrement(2, 1) =
288 522240 : prefactor * 2.0 * L * stress(1, 2);
289 : }
290 : else
291 : {
292 1259456 : GenericDenseVector<is_ad> inelastic_strain_increment(6);
293 1259456 : GenericDenseVector<is_ad> Ms(6);
294 1259456 : (*_hill_tensor)[_qp].vector_mult(Ms, stress_dev);
295 :
296 8816192 : for (unsigned int i = 0; i < 6; i++)
297 7556736 : inelastic_strain_increment(i) = Ms(i) * prefactor;
298 :
299 1259456 : inelasticStrainIncrement(0, 0) = inelastic_strain_increment(0);
300 1259456 : inelasticStrainIncrement(1, 1) = inelastic_strain_increment(1);
301 1259456 : inelasticStrainIncrement(2, 2) = inelastic_strain_increment(2);
302 :
303 1259456 : inelasticStrainIncrement(0, 1) = inelasticStrainIncrement(1, 0) = inelastic_strain_increment(3);
304 1259456 : inelasticStrainIncrement(1, 2) = inelasticStrainIncrement(2, 1) = inelastic_strain_increment(4);
305 1259456 : inelasticStrainIncrement(0, 2) = inelasticStrainIncrement(2, 0) = inelastic_strain_increment(5);
306 : }
307 :
308 1781696 : AnisotropicReturnCreepStressUpdateBaseTempl<is_ad>::computeStrainFinalize(
309 : inelasticStrainIncrement, stress, stress_dev, delta_gamma);
310 :
311 2518912 : this->_effective_inelastic_strain[_qp] = this->_effective_inelastic_strain_old[_qp] + delta_gamma;
312 : }
313 :
314 : template <bool is_ad>
315 : void
316 1814752 : HillCreepStressUpdateTempl<is_ad>::computeStressFinalize(
317 : const GenericRankTwoTensor<is_ad> & creepStrainIncrement,
318 : const GenericReal<is_ad> & /*delta_gamma*/,
319 : GenericRankTwoTensor<is_ad> & stress_new,
320 : const GenericDenseVector<is_ad> & /*stress_dev*/,
321 : const GenericRankTwoTensor<is_ad> & stress_old,
322 : const GenericRankFourTensor<is_ad> & elasticity_tensor)
323 : {
324 : // Class only valid for isotropic elasticity (for now)
325 1814752 : stress_new -= elasticity_tensor * creepStrainIncrement;
326 :
327 : // Compute the maximum time step allowed due to creep strain numerical integration error
328 2574528 : Real stress_dif = MetaPhysicL::raw_value(stress_new - stress_old).L2norm();
329 :
330 : // Get a representative value of the elasticity tensor
331 1814752 : Real elasticity_value =
332 : 1.0 / 3.0 *
333 1814752 : MetaPhysicL::raw_value((elasticity_tensor(0, 0, 0, 0) + elasticity_tensor(1, 1, 1, 1) +
334 : elasticity_tensor(2, 2, 2, 2)));
335 :
336 1814752 : if (std::abs(stress_dif) > TOLERANCE * TOLERANCE)
337 1562400 : this->_max_integration_error_time_step =
338 1562400 : _dt / (stress_dif / elasticity_value / this->_max_integration_error);
339 : else
340 252352 : this->_max_integration_error_time_step = std::numeric_limits<Real>::max();
341 1814752 : }
342 :
343 : template <bool is_ad>
344 : void
345 16479720 : HillCreepStressUpdateTempl<is_ad>::computeQsigmaChanged(
346 : GenericReal<is_ad> & qsigma_changed,
347 : const GenericDenseVector<is_ad> & stress_new,
348 : const GenericReal<is_ad> & delta_gamma,
349 : GenericRankTwoTensor<is_ad> & stress_changed)
350 : {
351 : GenericReal<is_ad> qsigma_square;
352 :
353 : // Hill constants, some constraints apply
354 16479720 : const Real & F = _hill_constants[_qp][0];
355 : const Real & G = _hill_constants[_qp][1];
356 : const Real & H = _hill_constants[_qp][2];
357 : const Real & L = _hill_constants[_qp][3];
358 : const Real & M = _hill_constants[_qp][4];
359 : const Real & N = _hill_constants[_qp][5];
360 :
361 16479720 : if (!this->_use_transformation)
362 : {
363 4658176 : qsigma_square = F * (stress_new(1) - stress_new(2)) * (stress_new(1) - stress_new(2));
364 4658176 : qsigma_square += G * (stress_new(2) - stress_new(0)) * (stress_new(2) - stress_new(0));
365 4658176 : qsigma_square += H * (stress_new(0) - stress_new(1)) * (stress_new(0) - stress_new(1));
366 4658176 : qsigma_square += 2 * L * stress_new(4) * stress_new(4);
367 4658176 : qsigma_square += 2 * M * stress_new(5) * stress_new(5);
368 4658176 : qsigma_square += 2 * N * stress_new(3) * stress_new(3);
369 : }
370 : else
371 : {
372 11821544 : GenericDenseVector<is_ad> Ms(6);
373 11821544 : (*_hill_tensor)[_qp].vector_mult(Ms, stress_new);
374 11821544 : qsigma_square = Ms.dot(stress_new);
375 : }
376 :
377 16479720 : GenericReal<is_ad> qsigma = std::sqrt(qsigma_square);
378 :
379 16479720 : if (!_anisotropic_elasticity)
380 23576208 : qsigma_changed = qsigma - 1.5 * _two_shear_modulus * delta_gamma;
381 : else
382 : {
383 33440 : GenericRankTwoTensor<is_ad> stress;
384 33440 : stress(0, 0) = stress_new(0);
385 33440 : stress(1, 1) = stress_new(1);
386 33440 : stress(2, 2) = stress_new(2);
387 33440 : stress(0, 1) = stress(1, 0) = stress_new(3);
388 33440 : stress(1, 2) = stress(2, 1) = stress_new(4);
389 33440 : stress(0, 2) = stress(2, 0) = stress_new(5);
390 :
391 33440 : GenericDenseVector<is_ad> b(6);
392 33440 : b(0) = H * (stress(0, 0) - stress(1, 1)) - G * (stress(2, 2) - stress(0, 0));
393 33440 : b(1) = F * (stress(1, 1) - stress(2, 2)) - H * (stress(0, 0) - stress(1, 1));
394 33440 : b(2) = G * (stress(2, 2) - stress(0, 0)) - F * (stress(1, 1) - stress(2, 2));
395 66880 : b(3) = 2.0 * N * stress(0, 1);
396 66880 : b(4) = 2.0 * L * stress(1, 2);
397 66880 : b(5) = 2.0 * M * stress(0, 2);
398 :
399 33440 : GenericRankTwoTensor<is_ad> inelasticStrainIncrement;
400 33440 : inelasticStrainIncrement(0, 0) = delta_gamma * b(0) / qsigma;
401 33440 : inelasticStrainIncrement(1, 1) = delta_gamma * b(1) / qsigma;
402 33440 : inelasticStrainIncrement(2, 2) = delta_gamma * b(2) / qsigma;
403 33440 : inelasticStrainIncrement(0, 1) = inelasticStrainIncrement(1, 0) = delta_gamma * b(3) / qsigma;
404 33440 : inelasticStrainIncrement(1, 2) = inelasticStrainIncrement(2, 1) = delta_gamma * b(4) / qsigma;
405 33440 : inelasticStrainIncrement(0, 2) = inelasticStrainIncrement(2, 0) = delta_gamma * b(5) / qsigma;
406 :
407 66880 : stress_changed = stress - _elasticity_tensor[_qp] * inelasticStrainIncrement;
408 :
409 : GenericReal<is_ad> qsigma_square_changed;
410 33440 : qsigma_square_changed = F * (stress_changed(1, 1) - stress_changed(2, 2)) *
411 : (stress_changed(1, 1) - stress_changed(2, 2));
412 33440 : qsigma_square_changed += G * (stress_changed(2, 2) - stress_changed(0, 0)) *
413 : (stress_changed(2, 2) - stress_changed(0, 0));
414 33440 : qsigma_square_changed += H * (stress_changed(0, 0) - stress_changed(1, 1)) *
415 : (stress_changed(0, 0) - stress_changed(1, 1));
416 66880 : qsigma_square_changed += 2 * L * stress_changed(1, 2) * stress_changed(1, 2);
417 66880 : qsigma_square_changed += 2 * M * stress_changed(0, 2) * stress_changed(0, 2);
418 66880 : qsigma_square_changed += 2 * N * stress_changed(0, 1) * stress_changed(0, 1);
419 33440 : qsigma_changed = std::sqrt(qsigma_square_changed);
420 : }
421 16479720 : }
422 :
423 : template class HillCreepStressUpdateTempl<false>;
424 : template class HillCreepStressUpdateTempl<true>;
|
