www.mooseframework.org
ComputeSmearedCrackingStress.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 /****************************************************************/
9 #include "StressUpdateBase.h"
10 
11 template <>
12 InputParameters
14 {
15  InputParameters params = validParams<ComputeMultipleInelasticStress>();
16  params.addClassDescription("Compute stress using a fixed smeared cracking model");
17  params.addParam<std::string>(
18  "cracking_release",
19  "abrupt",
20  "The cracking release type. Choices are abrupt (default) and exponential.");
21  params.addParam<Real>(
22  "cracking_residual_stress",
23  0.0,
24  "The fraction of the cracking stress allowed to be maintained following a crack.");
25  params.addRequiredCoupledVar(
26  "cracking_stress",
27  "The stress threshold beyond which cracking occurs. Negative values prevent cracking.");
28  params.addParam<std::vector<unsigned int>>(
29  "active_crack_planes", "Planes on which cracks are allowed (0,1,2 -> x,z,theta in RZ)");
30  params.addParam<unsigned int>(
31  "max_cracks", 3, "The maximum number of cracks allowed at a material point.");
32  params.addParam<Real>("cracking_neg_fraction",
33  0,
34  "The fraction of the cracking strain at which "
35  "a transitition begins during decreasing "
36  "strain to the original stiffness.");
37  params.addParam<Real>("cracking_beta",
38  1.0,
39  "Coefficient used to control the softening in the exponential model. "
40  "When set to 1, the initial softening slope is equal to the negative "
41  "of the Young's modulus. Smaller numbers scale down that slope.");
42  params.addParam<Real>(
43  "max_stress_correction",
44  1.0,
45  "Maximum permitted correction to the predicted stress as a ratio of the "
46  "stress change to the predicted stress from the previous step's damage level. "
47  "Values less than 1 will improve robustness, but not be as accurate.");
48 
49  params.addRangeCheckedParam<Real>(
50  "shear_retention_factor",
51  0.0,
52  "shear_retention_factor>=0 & shear_retention_factor<=1.0",
53  "Fraction of original shear stiffness to be retained after cracking");
54  params.set<std::vector<MaterialName>>("inelastic_models") = {};
55 
56  return params;
57 }
58 
59 namespace
60 {
62 getCrackingModel(const std::string & name)
63 {
64  std::string n(name);
65  std::transform(n.begin(), n.end(), n.begin(), ::tolower);
67  if (n == "abrupt")
69  else if (n == "exponential")
71  else if (n == "power")
74  mooseError("Unknown cracking model");
75  return cm;
76 }
77 }
78 
80  : ComputeMultipleInelasticStress(parameters),
81  _cracking_release(getCrackingModel(getParam<std::string>("cracking_release"))),
82  _cracking_residual_stress(getParam<Real>("cracking_residual_stress")),
83  _cracking_stress(coupledValue("cracking_stress")),
84  _active_crack_planes(3, 1),
85  _max_cracks(getParam<unsigned int>("max_cracks")),
86  _cracking_neg_fraction(getParam<Real>("cracking_neg_fraction")),
87  _cracking_beta(getParam<Real>("cracking_beta")),
88  _shear_retention_factor(getParam<Real>("shear_retention_factor")),
89  _max_stress_correction(getParam<Real>("max_stress_correction")),
90  _crack_flags(declareProperty<RealVectorValue>(_base_name + "crack_flags")),
91  _crack_flags_old(getMaterialPropertyOld<RealVectorValue>(_base_name + "crack_flags")),
92  _crack_count(NULL),
93  _crack_count_old(NULL),
94  _crack_rotation(declareProperty<RankTwoTensor>(_base_name + "crack_rotation")),
95  _crack_rotation_old(getMaterialPropertyOld<RankTwoTensor>(_base_name + "crack_rotation")),
96  _crack_strain(declareProperty<RealVectorValue>(_base_name + "crack_strain")),
97  _crack_strain_old(getMaterialPropertyOld<RealVectorValue>(_base_name + "crack_strain")),
98  _crack_max_strain(declareProperty<RealVectorValue>(_base_name + "crack_max_strain")),
99  _crack_max_strain_old(getMaterialPropertyOld<RealVectorValue>(_base_name + "crack_max_strain"))
100 {
102  {
103  _crack_count = &declareProperty<RealVectorValue>(_base_name + "crack_count");
104  _crack_count_old = &getMaterialPropertyOld<RealVectorValue>(_base_name + "crack_count");
105  }
106 
107  if (parameters.isParamValid("active_crack_planes"))
108  {
109  const std::vector<unsigned int> & planes =
110  getParam<std::vector<unsigned>>("active_crack_planes");
111  for (unsigned int i = 0; i < 3; ++i)
112  _active_crack_planes[i] = 0;
113 
114  for (unsigned int i = 0; i < planes.size(); ++i)
115  {
116  if (planes[i] > 2)
117  mooseError("Active planes must be 0, 1, or 2");
118  _active_crack_planes[planes[i]] = 1;
119  }
120  }
121  if (_cracking_residual_stress < 0 || _cracking_residual_stress > 1)
122  mooseError("cracking_residual_stress must be between 0 and 1");
123  if (parameters.isParamSetByUser("cracking_neg_fraction") &&
124  (_cracking_neg_fraction <= 0 || _cracking_neg_fraction > 1))
125  mooseError("cracking_neg_fraction must be > zero and <= 1");
126 }
127 
128 void
130 {
132 
133  _crack_flags[_qp](0) = _crack_flags[_qp](1) = _crack_flags[_qp](2) = 1;
134 
135  if (_crack_count)
136  (*_crack_count)[_qp] = 0;
137 
138  _crack_strain[_qp] = 0;
139  _crack_max_strain[_qp](0) = 0;
140  _crack_rotation[_qp].Identity();
141 }
142 
143 void
145 {
147  _is_finite_strain = hasBlockMaterialProperty<RankTwoTensor>(_base_name + "strain_increment");
148 
150  mooseError("ComputeSmearedCrackingStress requires that the elasticity tensor be "
151  "guaranteed isotropic");
152 }
153 
154 void
156 {
157  bool force_elasticity_rotation = false;
158 
159  if (!previouslyCracked())
161  else
162  {
164 
165  // Propagate behavior from the (now inactive) inelastic models
167  for (auto model : _models)
168  {
169  model->setQp(_qp);
170  model->propagateQpStatefulProperties();
171  }
172 
173  // Since the other inelastic models are inactive, they will not limit the time step
174  _matl_timestep_limit[_qp] = std::numeric_limits<Real>::max();
175 
176  // update elasticity tensor with old cracking status: crack_flags_old and crack_orientation_old
178 
179  // Calculate stress in intermediate configuration
181  // InitialStress Deprecation: remove these lines
185 
187  force_elasticity_rotation = true;
188  }
189 
190  // compute crack status and adjust stress
192 
194  {
195  finiteStrainRotation(force_elasticity_rotation);
196  _crack_rotation[_qp] =
197  _rotation_increment[_qp] * _crack_rotation[_qp] * _rotation_increment[_qp].transpose();
198  }
199 }
200 
201 void
203 {
204  const Real youngs_modulus =
206 
207  bool cracking_locally_active = false;
208 
209  const Real cracking_stress = _cracking_stress[_qp];
210 
211  if (cracking_stress > 0)
212  {
213  RealVectorValue crack_flags_local(1.0, 1.0, 1.0);
214  const RankTwoTensor & R = _crack_rotation_old[_qp];
215  RankTwoTensor ePrime(_elastic_strain_old[_qp]);
216  ePrime.rotate(R.transpose());
217 
218  for (unsigned int i = 0; i < 3; ++i)
219  {
220  // Update elasticity tensor based on crack status of the end of last time step
221  if (_crack_flags_old[_qp](i) < 1.0)
222  {
223  if (_cracking_neg_fraction == 0 && MooseUtils::absoluteFuzzyLessThan(ePrime(i, i), 0.0))
224  crack_flags_local(i) = 1.0;
225  else if (_cracking_neg_fraction > 0 &&
226  ePrime(i, i) < _crack_strain_old[_qp](i) * _cracking_neg_fraction &&
227  ePrime(i, i) > -_crack_strain_old[_qp](i) * _cracking_neg_fraction)
228  {
229  const Real etr = _cracking_neg_fraction * _crack_strain_old[_qp](i);
230  const Real Eo = cracking_stress / _crack_strain_old[_qp](i);
231  const Real Ec = Eo * _crack_flags_old[_qp](i);
232  const Real a = (Ec - Eo) / (4 * etr);
233  const Real b = (Ec + Eo) / 2;
234  // Compute the ratio of the current transition stiffness to the original stiffness
235  crack_flags_local(i) = (2.0 * a * etr + b) / Eo;
236  cracking_locally_active = true;
237  }
238  else
239  {
240  crack_flags_local(i) = _crack_flags_old[_qp](i);
241  cracking_locally_active = true;
242  }
243  }
244  }
245 
246  if (cracking_locally_active)
247  {
248  // Update the elasticity tensor in the crack coordinate system
249  const RealVectorValue & c = crack_flags_local;
250 
251  const bool c0_coupled = MooseUtils::absoluteFuzzyEqual(c(0), 1.0);
252  const bool c1_coupled = MooseUtils::absoluteFuzzyEqual(c(1), 1.0);
253  const bool c2_coupled = MooseUtils::absoluteFuzzyEqual(c(2), 1.0);
254 
255  const Real c01 = (c0_coupled && c1_coupled ? 1.0 : 0.0);
256  const Real c02 = (c0_coupled && c2_coupled ? 1.0 : 0.0);
257  const Real c12 = (c1_coupled && c2_coupled ? 1.0 : 0.0);
258 
259  const Real c01_shear_retention = (c0_coupled && c1_coupled ? 1.0 : _shear_retention_factor);
260  const Real c02_shear_retention = (c0_coupled && c2_coupled ? 1.0 : _shear_retention_factor);
261  const Real c12_shear_retention = (c1_coupled && c2_coupled ? 1.0 : _shear_retention_factor);
262 
263  // Filling with 9 components is sufficient because these are the only nonzero entries
264  // for isotropic or orthotropic materials.
265  std::vector<Real> local_elastic(9);
266 
267  local_elastic[0] = (c0_coupled ? _elasticity_tensor[_qp](0, 0, 0, 0) : c(0) * youngs_modulus);
268  local_elastic[1] = _elasticity_tensor[_qp](0, 0, 1, 1) * c01;
269  local_elastic[2] = _elasticity_tensor[_qp](0, 0, 2, 2) * c02;
270  local_elastic[3] = (c1_coupled ? _elasticity_tensor[_qp](1, 1, 1, 1) : c(1) * youngs_modulus);
271  local_elastic[4] = _elasticity_tensor[_qp](1, 1, 2, 2) * c12;
272  local_elastic[5] = (c2_coupled ? _elasticity_tensor[_qp](2, 2, 2, 2) : c(2) * youngs_modulus);
273  local_elastic[6] = _elasticity_tensor[_qp](1, 2, 1, 2) * c12_shear_retention;
274  local_elastic[7] = _elasticity_tensor[_qp](0, 2, 0, 2) * c02_shear_retention;
275  local_elastic[8] = _elasticity_tensor[_qp](0, 1, 0, 1) * c01_shear_retention;
276 
277  _local_elasticity_tensor.fillFromInputVector(local_elastic, RankFourTensor::symmetric9);
278 
279  // Rotate the modified elasticity tensor back into global coordinates
280  _local_elasticity_tensor.rotate(R);
281  }
282  }
283  if (!cracking_locally_active)
285 }
286 
287 void
289 {
290  const Real youngs_modulus =
292  const Real cracking_alpha = -youngs_modulus;
293 
294  Real cracking_stress = _cracking_stress[_qp];
295 
296  if (cracking_stress > 0)
297  {
298  // Initializing crack states
300 
301  for (unsigned i = 0; i < 3; ++i)
302  {
303  _crack_max_strain[_qp](i) = _crack_max_strain_old[_qp](i);
304  _crack_strain[_qp](i) = _crack_strain_old[_qp](i);
305  _crack_flags[_qp](i) = _crack_flags_old[_qp](i);
306  }
307 
308  // Compute crack orientations: updated _crack_rotation[_qp] based on current strain
309  RealVectorValue principal_strain;
310  computeCrackStrainAndOrientation(principal_strain);
311 
312  for (unsigned i = 0; i < 3; ++i)
313  {
314  if (principal_strain(i) > _crack_max_strain[_qp](i))
315  _crack_max_strain[_qp](i) = principal_strain(i);
316  }
317 
318  // Check for new cracks.
319  // Rotate stress to cracked orientation.
320  const RankTwoTensor & R = _crack_rotation[_qp];
321  RankTwoTensor sigmaPrime(_stress[_qp]);
322  sigmaPrime.rotate(R.transpose()); // stress in crack coordinates
323 
324  unsigned int num_cracks = 0;
325  for (unsigned int i = 0; i < 3; ++i)
326  {
327  if (_crack_flags_old[_qp](i) < 1)
328  ++num_cracks;
329  }
330 
331  bool cracked(false);
332  RealVectorValue sigma;
333  for (unsigned int i = 0; i < 3; ++i)
334  {
335  sigma(i) = sigmaPrime(i, i);
336 
337  Real crackFactor(1);
338 
340  (*_crack_count)[_qp](i) = (*_crack_count_old)[_qp](i);
341 
342  if ((_cracking_release == CR_POWER && sigma(i) > cracking_stress &&
343  _active_crack_planes[i] == 1))
344  {
345  cracked = true;
346  ++((*_crack_count)[_qp](i));
347 
348  // Assume Poisson's ratio drops to zero for this direction. Stiffness is then Young's
349  // modulus.
350  if ((*_crack_count_old)[_qp](i) == 0)
351  {
352  ++num_cracks;
353  _crack_strain[_qp](i) = cracking_stress / youngs_modulus;
354  }
355  // Compute stress, factor....
356  _crack_flags[_qp](i) *= 1. / 3.;
357 
358  if (_crack_max_strain[_qp](i) < _crack_strain[_qp](i))
359  _crack_max_strain[_qp](i) = _crack_strain[_qp](i);
360 
361  sigma(i) = _crack_flags[_qp](i) * youngs_modulus * principal_strain(i);
362  }
363  else if (_cracking_release != CR_POWER && _crack_flags_old[_qp](i) == 1 &&
364  sigma(i) > cracking_stress && num_cracks < _max_cracks &&
365  _active_crack_planes[i] == 1)
366  {
367  // A new crack
368  cracked = true;
369  ++num_cracks;
370 
371  // Assume Poisson's ratio drops to zero for this direction. Stiffness is then Young's
372  // modulus.
373  _crack_strain[_qp](i) = cracking_stress / youngs_modulus;
374 
375  if (_crack_max_strain[_qp](i) < _crack_strain[_qp](i))
376  _crack_max_strain[_qp](i) = _crack_strain[_qp](i);
377 
378  crackFactor = computeCrackFactor(
379  i, sigma(i), _crack_flags[_qp](i), cracking_stress, cracking_alpha, youngs_modulus);
380 
381  _crack_flags[_qp](i) = crackFactor;
382  }
383 
384  else if (_cracking_release != CR_POWER && _crack_flags_old[_qp](i) < 1 &&
385  std::abs(principal_strain(i) - _crack_max_strain[_qp](i)) < 1e-10)
386  {
387  // Previously cracked,
388  // Crack opening
389  cracked = true;
390  crackFactor = computeCrackFactor(
391  i, sigma(i), _crack_flags[_qp](i), cracking_stress, cracking_alpha, youngs_modulus);
392  _crack_flags[_qp](i) = crackFactor;
393  }
394 
395  else if (_cracking_neg_fraction > 0 &&
396  _crack_strain[_qp](i) * _cracking_neg_fraction > principal_strain(i) &&
397  -_crack_strain[_qp](i) * _cracking_neg_fraction < principal_strain(i))
398  {
399  cracked = true;
400  const Real etr = _cracking_neg_fraction * _crack_strain[_qp](i);
401  const Real Eo = cracking_stress / _crack_strain[_qp](i);
402  const Real Ec = Eo * _crack_flags_old[_qp](i);
403  const Real a = (Ec - Eo) / (4 * etr);
404  const Real b = (Ec + Eo) / 2;
405  const Real c = (Ec - Eo) * etr / 4;
406  sigma(i) = (a * principal_strain(i) + b) * principal_strain(i) + c;
407  }
408  }
409 
410  if (cracked)
411  applyCracksToTensor(_stress[_qp], sigma);
412  }
413 }
414 
415 void
417 {
418  // The rotation tensor is ordered such that known dirs appear last in the list of
419  // columns. So, if one dir is known, it corresponds with the last column in the
420  // rotation tensor.
421  //
422  // This convention is based on the eigen routine returning eigen values in
423  // ascending order.
424  const unsigned int numKnownDirs = getNumKnownCrackDirs();
425 
426  if (numKnownDirs == 0)
427  {
428  std::vector<Real> eigval(3, 0.0);
429  RankTwoTensor eigvec;
430 
431  _elastic_strain[_qp].symmetricEigenvaluesEigenvectors(eigval, eigvec);
432 
433  // If the elastic strain is beyond the cracking strain, save the eigen vectors as
434  // the rotation tensor.
435  _crack_rotation[_qp] = eigvec;
436 
437  principal_strain(0) = eigval[0];
438  principal_strain(1) = eigval[1];
439  principal_strain(2) = eigval[2];
440  }
441  else if (numKnownDirs == 1)
442  {
443  // This is easily the most complicated case.
444  // 1. Rotate the elastic strain to the orientation associated with the known
445  // crack.
446  // 2. Extract the upper 2x2 diagonal block into a separate tensor.
447  // 3. Run the eigen solver on the result.
448  // 4. Update the rotation tensor to reflect the effect of the 2 eigenvectors.
449 
450  // 1.
451  const RankTwoTensor & R = _crack_rotation[_qp];
452  RankTwoTensor ePrime(_elastic_strain[_qp]);
453  ePrime.rotate(R.transpose()); // elastic strain in principal coordinates
454 
455  // 2.
456  ColumnMajorMatrix e2x2(2, 2);
457  e2x2(0, 0) = ePrime(0, 0);
458  e2x2(1, 0) = ePrime(1, 0);
459  e2x2(0, 1) = ePrime(0, 1);
460  e2x2(1, 1) = ePrime(1, 1);
461 
462  // 3.
463  ColumnMajorMatrix e_val2x1(2, 1);
464  ColumnMajorMatrix e_vec2x2(2, 2);
465  e2x2.eigen(e_val2x1, e_vec2x2);
466 
467  // 4.
468  ColumnMajorMatrix e_vec(3, 3);
469  e_vec(0, 0) = e_vec2x2(0, 0);
470  e_vec(1, 0) = e_vec2x2(1, 0);
471  e_vec(2, 0) = 0;
472  e_vec(0, 1) = e_vec2x2(0, 1);
473  e_vec(1, 1) = e_vec2x2(1, 1);
474  e_vec(2, 1) = 0;
475  e_vec(2, 0) = 0;
476  e_vec(2, 1) = 0;
477  e_vec(2, 2) = 1;
478 
479  RankTwoTensor eigvec(e_vec(0, 0),
480  e_vec(1, 0),
481  e_vec(2, 0),
482  e_vec(0, 1),
483  e_vec(1, 1),
484  e_vec(2, 1),
485  e_vec(0, 2),
486  e_vec(1, 2),
487  e_vec(2, 2));
488 
489  _crack_rotation[_qp] = _crack_rotation_old[_qp] * eigvec; // Roe implementation
490 
491  principal_strain(0) = e_val2x1(0, 0);
492  principal_strain(1) = e_val2x1(1, 0);
493  principal_strain(2) = ePrime(2, 2);
494  }
495  else if (numKnownDirs == 2 || numKnownDirs == 3)
496  {
497  // Rotate to cracked orientation and pick off the strains in the rotated
498  // coordinate directions.
499  const RankTwoTensor & R = _crack_rotation[_qp];
500  RankTwoTensor ePrime(_elastic_strain[_qp]);
501  ePrime.rotate(R.transpose()); // elastic strain in principal coordinates
502 
503  principal_strain(0) = ePrime(0, 0);
504  principal_strain(1) = ePrime(1, 1);
505  principal_strain(2) = ePrime(2, 2);
506  }
507  else
508  mooseError("Invalid number of known crack directions");
509 }
510 
511 unsigned int
513 {
514  const unsigned fromElement = 0;
515  unsigned int retVal = 0;
516  for (unsigned int i = 0; i < 3 - fromElement; ++i)
517  retVal += (_crack_flags_old[_qp](i) < 1);
518  return retVal + fromElement;
519 }
520 
521 Real
523  Real & sigma,
524  Real & flag_value,
525  const Real cracking_stress,
526  const Real cracking_alpha,
527  const Real youngs_modulus)
528 {
530  {
531  if (_crack_max_strain[_qp](i) < _crack_strain[_qp](i))
532  {
533  mooseError("Max strain less than crack strain: ",
534  i,
535  " ",
536  sigma,
537  ", ",
538  _crack_max_strain[_qp](i),
539  ", ",
540  _crack_strain[_qp](i),
541  ", ",
542  _elastic_strain[_qp](0, 0),
543  ", ",
544  _elastic_strain[_qp](1, 1),
545  ", ",
546  _elastic_strain[_qp](2, 2),
547  ", ",
548  _elastic_strain[_qp](0, 1),
549  ", ",
550  _elastic_strain[_qp](0, 2),
551  ", ",
552  _elastic_strain[_qp](1, 2));
553  }
554  const Real crackMaxStrain = _crack_max_strain[_qp](i);
555  // Compute stress that follows exponental curve
556  sigma = cracking_stress * (_cracking_residual_stress +
557  (1.0 - _cracking_residual_stress) *
558  std::exp(cracking_alpha * _cracking_beta / cracking_stress *
559  (crackMaxStrain - _crack_strain[_qp](i))));
560  // Compute ratio of current stiffness to original stiffness
561  flag_value = sigma * _crack_strain[_qp](i) / (crackMaxStrain * cracking_stress);
562  }
563  else
564  {
565  if (_cracking_residual_stress == 0)
566  {
567  const Real tiny = 1e-16;
568  flag_value = tiny;
569  sigma = tiny * _crack_strain[_qp](i) * youngs_modulus;
570  }
571  else
572  {
573  sigma = _cracking_residual_stress * cracking_stress;
574  flag_value = sigma / (_crack_max_strain[_qp](i) * youngs_modulus);
575  }
576  }
577  if (flag_value < 0)
578  {
579  std::stringstream err;
580  err << "Negative crack flag found: " << i << " " << flag_value << ", "
581  << _crack_max_strain[_qp](i) << ", " << _crack_strain[_qp](i) << ", " << std::endl;
582  mooseError(err.str());
583  }
584  return flag_value;
585 }
586 
587 void
589  const RealVectorValue & sigma)
590 {
591  // Get transformation matrix
592  const RankTwoTensor & R = _crack_rotation[_qp];
593  // Rotate to crack frame
594  tensor.rotate(R.transpose());
595 
596  // Reset stress if cracked
597  for (unsigned int i = 0; i < 3; ++i)
598  if (_crack_flags[_qp](i) < 1)
599  {
600  const Real stress_correction_ratio = (tensor(i, i) - sigma(i)) / tensor(i, i);
601  if (stress_correction_ratio > _max_stress_correction)
602  tensor(i, i) *= (1.0 - _max_stress_correction);
603  else if (stress_correction_ratio < -_max_stress_correction)
604  tensor(i, i) *= (1.0 + _max_stress_correction);
605  else
606  tensor(i, i) = sigma(i);
607  }
608 
609  // Rotate back to global frame
610  tensor.rotate(R);
611 }
612 
613 bool
615 {
616  for (unsigned int i = 0; i < 3; ++i)
617  if (_crack_flags_old[_qp](i) < 1.0)
618  return true;
619  return false;
620 }
MaterialProperty< RankFourTensor > & _Jacobian_mult
derivative of stress w.r.t. strain (_dstress_dstrain)
const std::string _elasticity_tensor_name
ComputeMultipleInelasticStress computes the stress, the consistent tangent operator (or an approximat...
const MaterialProperty< RankTwoTensor > & _rotation_increment
MaterialProperty< RealVectorValue > & _crack_max_strain
const MaterialProperty< RankTwoTensor > & _elastic_strain_old
std::vector< unsigned int > _active_crack_planes
const MaterialProperty< RankTwoTensor > & _crack_rotation_old
const MaterialProperty< RealVectorValue > & _crack_flags_old
MaterialProperty< Real > & _matl_timestep_limit
const MaterialProperty< RankTwoTensor > & _strain_increment
MaterialProperty< RealVectorValue > & _crack_strain
const MaterialProperty< RealVectorValue > & _crack_max_strain_old
MaterialProperty< RankTwoTensor > & _crack_rotation
virtual Real computeCrackFactor(int i, Real &sigma, Real &flag_value, const Real cracking_stress, const Real cracking_alpha, const Real youngs_modulus)
const MaterialProperty< RealVectorValue > & _crack_strain_old
MaterialProperty< RealVectorValue > * _crack_count
void applyCracksToTensor(RankTwoTensor &tensor, const RealVectorValue &sigma)
MaterialProperty< RankTwoTensor > & _stress
ComputeSmearedCrackingStress(const InputParameters &parameters)
const MaterialProperty< RankTwoTensor > & _inelastic_strain_old
old value of inelastic strain
virtual void computeQpStressIntermediateConfiguration()
Compute the stress for the current QP, but do not rotate tensors from the intermediate configuration ...
InputParameters validParams< ComputeSmearedCrackingStress >()
virtual unsigned int getNumKnownCrackDirs() const
const MaterialProperty< RealVectorValue > * _crack_count_old
MaterialProperty< RealVectorValue > & _crack_flags
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Rank-4 and Rank-2 elasticity and elastic strain tensors.
void computeCrackStrainAndOrientation(RealVectorValue &principal_strain)
void addQpInitialStress()
InitialStress Deprecation: remove this method.
Real getIsotropicYoungsModulus(const RankFourTensor &elasticity_tensor)
Get the Young&#39;s modulus for an isotropic elasticity tensor param elasticity_tensor the tensor (must b...
const std::string _base_name
const CRACKING_RELEASE _cracking_release
Input parameters for smeared crack models.
virtual void finiteStrainRotation(const bool force_elasticity_rotation=false)
Rotate _elastic_strain, _stress, _inelastic_strain, and _Jacobian_mult to the new configuration...
std::vector< StressUpdateBase * > _models
The user supplied list of inelastic models to use in the simulation.
const bool _perform_finite_strain_rotations
after updateQpState, rotate the stress, elastic_strain, inelastic_strain and Jacobian_mult using _rot...
MaterialProperty< RankTwoTensor > & _elastic_strain
bool hasGuaranteedMaterialProperty(const MaterialPropertyName &prop, Guarantee guarantee)
InputParameters validParams< ComputeMultipleInelasticStress >()
virtual void rotateQpInitialStress()
InitialStress Deprecation: remove this method.
MaterialProperty< RankTwoTensor > & _inelastic_strain
The sum of the inelastic strains that come from the plastic models.