www.mooseframework.org
Public Member Functions | Protected Types | Protected Member Functions | Protected Attributes | List of all members
ComputeMultipleInelasticStress Class Reference

ComputeMultipleInelasticStress computes the stress, the consistent tangent operator (or an approximation to it), and a decomposition of the strain into elastic and inelastic parts. More...

#include <ComputeMultipleInelasticStress.h>

Inheritance diagram for ComputeMultipleInelasticStress:
[legend]

Public Member Functions

 ComputeMultipleInelasticStress (const InputParameters &parameters)
 
virtual void initialSetup () override
 

Protected Types

enum  TangentOperatorEnum { TangentOperatorEnum::elastic, TangentOperatorEnum::nonlinear }
 what sort of Tangent operator to calculate More...
 

Protected Member Functions

virtual void initQpStatefulProperties () override
 
virtual void computeQpStress () override
 Compute the stress and store it in the _stress material property for the current quadrature point. More...
 
virtual void computeQpStressIntermediateConfiguration ()
 Compute the stress for the current QP, but do not rotate tensors from the intermediate configuration to the new configuration. More...
 
virtual void finiteStrainRotation (const bool force_elasticity_rotation=false)
 Rotate _elastic_strain, _stress, _inelastic_strain, and _Jacobian_mult to the new configuration. More...
 
virtual void updateQpState (RankTwoTensor &elastic_strain_increment, RankTwoTensor &combined_inelastic_strain_increment)
 Given the _strain_increment[_qp], iterate over all of the user-specified recompute materials in order to find an admissible stress (which is placed into _stress[_qp]) and set of inelastic strains, as well as the tangent operator (which is placed into _Jacobian_mult[_qp]). More...
 
virtual void updateQpStateSingleModel (unsigned model_number, RankTwoTensor &elastic_strain_increment, RankTwoTensor &combined_inelastic_strain_increment)
 An optimised version of updateQpState that gets used when the number of plastic models is unity, or when we're cycling through models Given the _strain_increment[_qp], find an admissible stress (which is put into _stress[_qp]) and inelastic strain, as well as the tangent operator (which is placed into _Jacobian_mult[_qp]) More...
 
virtual void computeQpJacobianMult ()
 Using _elasticity_tensor[_qp] and the consistent tangent operators, _consistent_tangent_operator[...] computed by the inelastic models, compute _Jacobian_mult[_qp]. More...
 
virtual void computeAdmissibleState (unsigned model_number, RankTwoTensor &elastic_strain_increment, RankTwoTensor &inelastic_strain_increment, RankFourTensor &consistent_tangent_operator)
 Given a trial stress (_stress[_qp]) and a strain increment (elastic_strain_increment) let the model_number model produce an admissible stress (gets placed back in _stress[_qp]), and decompose the strain increment into an elastic part (gets placed back into elastic_strain_increment) and an inelastic part (inelastic_strain_increment), as well as computing the consistent_tangent_operator. More...
 
virtual void computeQpProperties () override
 
bool hasGuaranteedMaterialProperty (const MaterialPropertyName &prop, Guarantee guarantee)
 

Protected Attributes

const bool _perform_finite_strain_rotations
 after updateQpState, rotate the stress, elastic_strain, inelastic_strain and Jacobian_mult using _rotation_increment More...
 
MaterialProperty< RankTwoTensor > & _inelastic_strain
 The sum of the inelastic strains that come from the plastic models. More...
 
const MaterialProperty< RankTwoTensor > & _inelastic_strain_old
 old value of inelastic strain More...
 
enum ComputeMultipleInelasticStress::TangentOperatorEnum _tangent_operator_type
 
const unsigned _num_models
 number of plastic models More...
 
std::vector< bool > _tangent_computation_flag
 Flags to compute tangent during updateState call. More...
 
TangentCalculationMethod _tangent_calculation_method
 Calculation method for the tangent modulus. More...
 
const std::vector< Real > _inelastic_weights
 _inelastic_strain = sum_i (_inelastic_weights_i * inelastic_strain_from_model_i) More...
 
std::vector< RankFourTensor_consistent_tangent_operator
 the consistent tangent operators computed by each plastic model More...
 
const bool _cycle_models
 whether to cycle through the models, using only one model per timestep More...
 
MaterialProperty< Real > & _matl_timestep_limit
 
const RankFourTensor _identity_symmetric_four
 Rank four symmetric identity tensor. More...
 
std::vector< StressUpdateBase * > _models
 The user supplied list of inelastic models to use in the simulation. More...
 
bool _is_elasticity_tensor_guaranteed_isotropic
 is the elasticity tensor guaranteed to be isotropic? More...
 
const std::string _elasticity_tensor_name
 Name of the elasticity tensor material property. More...
 
const MaterialProperty< RankFourTensor > & _elasticity_tensor
 Elasticity tensor material property. More...
 
const MaterialProperty< RankTwoTensor > & _rotation_increment
 Rotation increment material property. More...
 
const MaterialProperty< RankTwoTensor > & _stress_old
 Old state of the stress tensor material property. More...
 
const std::string _base_name
 Base name prepended to all material property names to allow for multi-material systems. More...
 
const MaterialProperty< RankTwoTensor > & _mechanical_strain
 Mechanical strain material property. More...
 
MaterialProperty< RankTwoTensor > & _stress
 Stress material property. More...
 
MaterialProperty< RankTwoTensor > & _elastic_strain
 Elastic strain material property. More...
 
const MaterialProperty< RankTwoTensor > & _extra_stress
 Extra stress tensor. More...
 
std::vector< Function * > _initial_stress_fcn
 initial stress components More...
 
MaterialProperty< RankFourTensor > & _Jacobian_mult
 derivative of stress w.r.t. strain (_dstress_dstrain) More...
 
const unsigned int _max_iterations
 Input parameters associated with the recompute iteration to return the stress state to the yield surface. More...
 
const Real _relative_tolerance
 
const Real _absolute_tolerance
 
const bool _internal_solve_full_iteration_history
 
const MaterialProperty< RankTwoTensor > & _elastic_strain_old
 Strain tensors. More...
 
const MaterialProperty< RankTwoTensor > & _strain_increment
 

Detailed Description

ComputeMultipleInelasticStress computes the stress, the consistent tangent operator (or an approximation to it), and a decomposition of the strain into elastic and inelastic parts.

By default finite strains are assumed.

The elastic strain is calculated by subtracting the computed inelastic strain increment tensor from the mechanical strain tensor. Mechanical strain is considered as the sum of the elastic and inelastic (plastic, creep, ect) strains.

This material is used to call the recompute iterative materials of a number of specified inelastic models that inherit from StressUpdateBase. It iterates over the specified inelastic models until the change in stress is within a user-specified tolerance, in order to produce the stress, the consistent tangent operator and the elastic and inelastic strains for the time increment.

Definition at line 38 of file ComputeMultipleInelasticStress.h.

Member Enumeration Documentation

◆ TangentOperatorEnum

what sort of Tangent operator to calculate

Enumerator
elastic 
nonlinear 

Definition at line 141 of file ComputeMultipleInelasticStress.h.

141 { elastic, nonlinear } _tangent_operator_type;
enum ComputeMultipleInelasticStress::TangentOperatorEnum _tangent_operator_type

Constructor & Destructor Documentation

◆ ComputeMultipleInelasticStress()

ComputeMultipleInelasticStress::ComputeMultipleInelasticStress ( const InputParameters &  parameters)

Definition at line 73 of file ComputeMultipleInelasticStress.C.

75  _max_iterations(parameters.get<unsigned int>("max_iterations")),
76  _relative_tolerance(parameters.get<Real>("relative_tolerance")),
77  _absolute_tolerance(parameters.get<Real>("absolute_tolerance")),
78  _internal_solve_full_iteration_history(getParam<bool>("internal_solve_full_iteration_history")),
79  _perform_finite_strain_rotations(getParam<bool>("perform_finite_strain_rotations")),
80  _elastic_strain_old(getMaterialPropertyOld<RankTwoTensor>(_base_name + "elastic_strain")),
81  _strain_increment(getMaterialProperty<RankTwoTensor>(_base_name + "strain_increment")),
82  _inelastic_strain(declareProperty<RankTwoTensor>(_base_name + "combined_inelastic_strain")),
84  getMaterialPropertyOld<RankTwoTensor>(_base_name + "combined_inelastic_strain")),
85  _tangent_operator_type(getParam<MooseEnum>("tangent_operator").getEnum<TangentOperatorEnum>()),
86  _num_models(getParam<std::vector<MaterialName>>("inelastic_models").size()),
89  _inelastic_weights(isParamValid("combined_inelastic_strain_weights")
90  ? getParam<std::vector<Real>>("combined_inelastic_strain_weights")
91  : std::vector<Real>(_num_models, true)),
93  _cycle_models(getParam<bool>("cycle_models")),
94  _matl_timestep_limit(declareProperty<Real>("matl_timestep_limit")),
95  _identity_symmetric_four(RankFourTensor::initIdentitySymmetricFour)
96 {
97  if (_inelastic_weights.size() != _num_models)
98  mooseError(
99  "ComputeMultipleInelasticStress: combined_inelastic_strain_weights must contain the same "
100  "number of entries as inelastic_models ",
101  _inelastic_weights.size(),
102  " ",
103  _num_models);
104 }
std::vector< RankFourTensor > _consistent_tangent_operator
the consistent tangent operators computed by each plastic model
const MaterialProperty< RankTwoTensor > & _elastic_strain_old
Strain tensors.
MaterialProperty< Real > & _matl_timestep_limit
const MaterialProperty< RankTwoTensor > & _strain_increment
const std::vector< Real > _inelastic_weights
_inelastic_strain = sum_i (_inelastic_weights_i * inelastic_strain_from_model_i)
const MaterialProperty< RankTwoTensor > & _inelastic_strain_old
old value of inelastic strain
std::vector< bool > _tangent_computation_flag
Flags to compute tangent during updateState call.
const RankFourTensor _identity_symmetric_four
Rank four symmetric identity tensor.
const bool _cycle_models
whether to cycle through the models, using only one model per timestep
TangentCalculationMethod _tangent_calculation_method
Calculation method for the tangent modulus.
const unsigned _num_models
number of plastic models
ComputeFiniteStrainElasticStress(const InputParameters &parameters)
const std::string _base_name
Base name prepended to all material property names to allow for multi-material systems.
enum ComputeMultipleInelasticStress::TangentOperatorEnum _tangent_operator_type
const bool _perform_finite_strain_rotations
after updateQpState, rotate the stress, elastic_strain, inelastic_strain and Jacobian_mult using _rot...
MaterialProperty< RankTwoTensor > & _inelastic_strain
The sum of the inelastic strains that come from the plastic models.
const unsigned int _max_iterations
Input parameters associated with the recompute iteration to return the stress state to the yield surf...

Member Function Documentation

◆ computeAdmissibleState()

void ComputeMultipleInelasticStress::computeAdmissibleState ( unsigned  model_number,
RankTwoTensor elastic_strain_increment,
RankTwoTensor inelastic_strain_increment,
RankFourTensor consistent_tangent_operator 
)
protectedvirtual

Given a trial stress (_stress[_qp]) and a strain increment (elastic_strain_increment) let the model_number model produce an admissible stress (gets placed back in _stress[_qp]), and decompose the strain increment into an elastic part (gets placed back into elastic_strain_increment) and an inelastic part (inelastic_strain_increment), as well as computing the consistent_tangent_operator.

Parameters
model_numberThe inelastic model to use
elastic_strain_incrementUpon input, this is the strain increment. Upon output, it is the elastic part of the strain increment
inelastic_strain_incrementThe inelastic strain increment corresponding to the supplied strain increment
consistent_tangent_operatorThe consistent tangent operator

Reimplemented in ComputeMultipleInelasticCosseratStress.

Definition at line 418 of file ComputeMultipleInelasticStress.C.

Referenced by ComputeMultipleInelasticCosseratStress::computeAdmissibleState(), updateQpState(), and updateQpStateSingleModel().

422 {
423  const bool jac = _fe_problem.currentlyComputingJacobian();
424  _models[model_number]->updateState(elastic_strain_increment,
425  inelastic_strain_increment,
426  _rotation_increment[_qp],
427  _stress[_qp],
428  _stress_old[_qp],
429  _elasticity_tensor[_qp],
430  _elastic_strain_old[_qp],
431  (jac && _tangent_computation_flag[model_number]),
432  consistent_tangent_operator);
433 
434  if (jac && !_tangent_computation_flag[model_number])
435  {
437  consistent_tangent_operator.zero();
438  else
439  consistent_tangent_operator = _elasticity_tensor[_qp];
440  }
441 }
const MaterialProperty< RankTwoTensor > & _rotation_increment
Rotation increment material property.
const MaterialProperty< RankTwoTensor > & _elastic_strain_old
Strain tensors.
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Elasticity tensor material property.
MaterialProperty< RankTwoTensor > & _stress
Stress material property.
std::vector< bool > _tangent_computation_flag
Flags to compute tangent during updateState call.
TangentCalculationMethod _tangent_calculation_method
Calculation method for the tangent modulus.
const MaterialProperty< RankTwoTensor > & _stress_old
Old state of the stress tensor material property.
std::vector< StressUpdateBase * > _models
The user supplied list of inelastic models to use in the simulation.

◆ computeQpJacobianMult()

void ComputeMultipleInelasticStress::computeQpJacobianMult ( )
protectedvirtual

Using _elasticity_tensor[_qp] and the consistent tangent operators, _consistent_tangent_operator[...] computed by the inelastic models, compute _Jacobian_mult[_qp].

Reimplemented in ComputeMultipleInelasticCosseratStress.

Definition at line 350 of file ComputeMultipleInelasticStress.C.

Referenced by updateQpState().

351 {
355  {
357  for (unsigned i_rmm = 0; i_rmm < _num_models; ++i_rmm)
358  A += _consistent_tangent_operator[i_rmm];
359  mooseAssert(A.isSymmetric(), "Tangent operator isn't symmetric");
360  _Jacobian_mult[_qp] = A.invSymm() * _elasticity_tensor[_qp];
361  }
362  else
363  {
364  const RankFourTensor E_inv = _elasticity_tensor[_qp].invSymm();
366  for (unsigned i_rmm = 1; i_rmm < _num_models; ++i_rmm)
367  _Jacobian_mult[_qp] = _consistent_tangent_operator[i_rmm] * E_inv * _Jacobian_mult[_qp];
368  }
369 }
MaterialProperty< RankFourTensor > & _Jacobian_mult
derivative of stress w.r.t. strain (_dstress_dstrain)
std::vector< RankFourTensor > _consistent_tangent_operator
the consistent tangent operators computed by each plastic model
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Elasticity tensor material property.
const RankFourTensor _identity_symmetric_four
Rank four symmetric identity tensor.
TangentCalculationMethod _tangent_calculation_method
Calculation method for the tangent modulus.
const unsigned _num_models
number of plastic models

◆ computeQpProperties()

void ComputeStressBase::computeQpProperties ( )
overrideprotectedvirtualinherited

Definition at line 49 of file ComputeStressBase.C.

50 {
52 
53  // Add in extra stress
54  _stress[_qp] += _extra_stress[_qp];
55 }
virtual void computeQpStress()=0
Compute the stress and store it in the _stress material property for the current quadrature point...
MaterialProperty< RankTwoTensor > & _stress
Stress material property.
const MaterialProperty< RankTwoTensor > & _extra_stress
Extra stress tensor.

◆ computeQpStress()

void ComputeMultipleInelasticStress::computeQpStress ( )
overrideprotectedvirtual

Compute the stress and store it in the _stress material property for the current quadrature point.

Reimplemented from ComputeFiniteStrainElasticStress.

Reimplemented in ComputeMultipleInelasticCosseratStress, and ComputeSmearedCrackingStress.

Definition at line 173 of file ComputeMultipleInelasticStress.C.

Referenced by ComputeMultipleInelasticCosseratStress::computeQpStress().

174 {
178 }
virtual void computeQpStressIntermediateConfiguration()
Compute the stress for the current QP, but do not rotate tensors from the intermediate configuration ...
virtual void finiteStrainRotation(const bool force_elasticity_rotation=false)
Rotate _elastic_strain, _stress, _inelastic_strain, and _Jacobian_mult to the new configuration...
const bool _perform_finite_strain_rotations
after updateQpState, rotate the stress, elastic_strain, inelastic_strain and Jacobian_mult using _rot...

◆ computeQpStressIntermediateConfiguration()

void ComputeMultipleInelasticStress::computeQpStressIntermediateConfiguration ( )
protectedvirtual

Compute the stress for the current QP, but do not rotate tensors from the intermediate configuration to the new configuration.

Definition at line 181 of file ComputeMultipleInelasticStress.C.

Referenced by ComputeSmearedCrackingStress::computeQpStress(), and computeQpStress().

182 {
183  RankTwoTensor elastic_strain_increment;
184  RankTwoTensor combined_inelastic_strain_increment;
185 
186  if (_num_models == 0)
187  {
189 
190  // If the elasticity tensor values have changed and the tensor is isotropic,
191  // use the old strain to calculate the old stress
194  else
195  _stress[_qp] = _stress_old[_qp] + _elasticity_tensor[_qp] * _strain_increment[_qp];
196 
197  if (_fe_problem.currentlyComputingJacobian())
199  }
200  else
201  {
202  if (_num_models == 1 || _cycle_models)
203  updateQpStateSingleModel((_t_step - 1) % _num_models,
204  elastic_strain_increment,
205  combined_inelastic_strain_increment);
206  else
207  updateQpState(elastic_strain_increment, combined_inelastic_strain_increment);
208 
209  _elastic_strain[_qp] = _elastic_strain_old[_qp] + elastic_strain_increment;
210  _inelastic_strain[_qp] = _inelastic_strain_old[_qp] + combined_inelastic_strain_increment;
211  }
212 }
bool _is_elasticity_tensor_guaranteed_isotropic
is the elasticity tensor guaranteed to be isotropic?
MaterialProperty< RankFourTensor > & _Jacobian_mult
derivative of stress w.r.t. strain (_dstress_dstrain)
const MaterialProperty< RankTwoTensor > & _elastic_strain_old
Strain tensors.
const MaterialProperty< RankTwoTensor > & _strain_increment
virtual void updateQpState(RankTwoTensor &elastic_strain_increment, RankTwoTensor &combined_inelastic_strain_increment)
Given the _strain_increment[_qp], iterate over all of the user-specified recompute materials in order...
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Elasticity tensor material property.
MaterialProperty< RankTwoTensor > & _stress
Stress material property.
const MaterialProperty< RankTwoTensor > & _inelastic_strain_old
old value of inelastic strain
virtual void updateQpStateSingleModel(unsigned model_number, RankTwoTensor &elastic_strain_increment, RankTwoTensor &combined_inelastic_strain_increment)
An optimised version of updateQpState that gets used when the number of plastic models is unity...
const bool _cycle_models
whether to cycle through the models, using only one model per timestep
const unsigned _num_models
number of plastic models
const MaterialProperty< RankTwoTensor > & _stress_old
Old state of the stress tensor material property.
const bool _perform_finite_strain_rotations
after updateQpState, rotate the stress, elastic_strain, inelastic_strain and Jacobian_mult using _rot...
MaterialProperty< RankTwoTensor > & _elastic_strain
Elastic strain material property.
MaterialProperty< RankTwoTensor > & _inelastic_strain
The sum of the inelastic strains that come from the plastic models.

◆ finiteStrainRotation()

void ComputeMultipleInelasticStress::finiteStrainRotation ( const bool  force_elasticity_rotation = false)
protectedvirtual

Rotate _elastic_strain, _stress, _inelastic_strain, and _Jacobian_mult to the new configuration.

Parameters
force_elasticity_rotationForce the elasticity tensor to be rotated, even if it is not deemed necessary.

Definition at line 215 of file ComputeMultipleInelasticStress.C.

Referenced by ComputeSmearedCrackingStress::computeQpStress(), and computeQpStress().

216 {
217  _elastic_strain[_qp] =
218  _rotation_increment[_qp] * _elastic_strain[_qp] * _rotation_increment[_qp].transpose();
219  _stress[_qp] = _rotation_increment[_qp] * _stress[_qp] * _rotation_increment[_qp].transpose();
220  _inelastic_strain[_qp] =
221  _rotation_increment[_qp] * _inelastic_strain[_qp] * _rotation_increment[_qp].transpose();
222  if (force_elasticity_rotation ||
225  _Jacobian_mult[_qp].rotate(_rotation_increment[_qp]);
226 }
bool _is_elasticity_tensor_guaranteed_isotropic
is the elasticity tensor guaranteed to be isotropic?
MaterialProperty< RankFourTensor > & _Jacobian_mult
derivative of stress w.r.t. strain (_dstress_dstrain)
const MaterialProperty< RankTwoTensor > & _rotation_increment
Rotation increment material property.
MaterialProperty< RankTwoTensor > & _stress
Stress material property.
TangentCalculationMethod _tangent_calculation_method
Calculation method for the tangent modulus.
const unsigned _num_models
number of plastic models
MaterialProperty< RankTwoTensor > & _elastic_strain
Elastic strain material property.
MaterialProperty< RankTwoTensor > & _inelastic_strain
The sum of the inelastic strains that come from the plastic models.

◆ hasGuaranteedMaterialProperty()

bool GuaranteeConsumer::hasGuaranteedMaterialProperty ( const MaterialPropertyName &  prop,
Guarantee  guarantee 
)
protectedinherited

Definition at line 28 of file GuaranteeConsumer.C.

Referenced by ComputeFiniteStrainElasticStress::initialSetup(), ComputeSmearedCrackingStress::initialSetup(), ComputeLinearElasticPFFractureStress::initialSetup(), and initialSetup().

30 {
31  if (!_gc_feproblem->startedInitialSetup())
32  mooseError("hasGuaranteedMaterialProperty() needs to be called in initialSetup()");
33 
34  // Reference to MaterialWarehouse for testing and retrieving block ids
35  const auto & warehouse = _gc_feproblem->getMaterialWarehouse();
36 
37  // Complete set of ids that this object is active
38  const auto & ids = (_gc_block_restrict && _gc_block_restrict->blockRestricted())
39  ? _gc_block_restrict->blockIDs()
40  : _gc_feproblem->mesh().meshSubdomains();
41 
42  // Loop over each id for this object
43  for (const auto & id : ids)
44  {
45  // If block materials exist, look if any issue the required guarantee
46  if (warehouse.hasActiveBlockObjects(id))
47  {
48  const std::vector<std::shared_ptr<Material>> & mats = warehouse.getActiveBlockObjects(id);
49  for (const auto & mat : mats)
50  {
51  const auto & mat_props = mat->getSuppliedItems();
52  if (mat_props.count(prop_name))
53  {
54  auto guarantee_mat = dynamic_cast<GuaranteeProvider *>(mat.get());
55  if (guarantee_mat && !guarantee_mat->hasGuarantee(prop_name, guarantee))
56  {
57  // we found at least one material on the set of block we operate on
58  // that does _not_ provide the requested guarantee
59  return false;
60  }
61  }
62  }
63  }
64  }
65 
66  return true;
67 }
Add-on class that provides the functionality to issue guarantees for declared material properties...
BlockRestrictable *const _gc_block_restrict
Access block restrictions of the object with this interface.
FEProblemBase *const _gc_feproblem
Reference to the FEProblemBase class.

◆ initialSetup()

void ComputeMultipleInelasticStress::initialSetup ( )
overridevirtual

Reimplemented in ComputeSmearedCrackingStress.

Definition at line 114 of file ComputeMultipleInelasticStress.C.

Referenced by ComputeSmearedCrackingStress::initialSetup().

115 {
118 
119  std::vector<MaterialName> models = getParam<std::vector<MaterialName>>("inelastic_models");
120 
121  for (unsigned int i = 0; i < _num_models; ++i)
122  {
123  StressUpdateBase * rrr = dynamic_cast<StressUpdateBase *>(&getMaterialByName(models[i]));
124 
125  if (rrr)
126  {
127  _models.push_back(rrr);
129  mooseError("Model " + models[i] +
130  " requires an isotropic elasticity tensor, but the one supplied is not "
131  "guaranteed isotropic");
132  }
133  else
134  mooseError("Model " + models[i] + " is not compatible with ComputeMultipleInelasticStress");
135  }
136 
137  // Check if tangent calculation methods are consistent. If all models have
138  // TangentOperatorEnum::ELASTIC or tangent_operator is set by the user as elasic, then the tangent
139  // is never calculated: J_tot = C. If PARTIAL and NONE models are present, utilize PARTIAL
140  // formulation: J_tot = (I + J_1 + ... J_N)^-1 C. If FULL and NONE models are present, utilize
141  // FULL formulation: J_tot = J_1 * C^-1 * J_2 * C^-1 * ... J_N * C. If PARTIAL and FULL models are
142  // present, error out.
143 
145  {
146  bool full_present = false;
147  bool partial_present = false;
148  for (unsigned int i = 0; i < _num_models; ++i)
149  {
150  if (_models[i]->getTangentCalculationMethod() == TangentCalculationMethod::FULL)
151  {
152  full_present = true;
153  _tangent_computation_flag[i] = true;
155  }
156  else if (_models[i]->getTangentCalculationMethod() == TangentCalculationMethod::PARTIAL)
157  {
158  partial_present = true;
159  _tangent_computation_flag[i] = true;
161  }
162  }
163  if (full_present && partial_present)
164  mooseError("In ",
165  _name,
166  ": Models that calculate the full tangent operator and the partial tangent "
167  "operator are being combined. Either set tangent_operator to elastic, implement "
168  "the corrent tangent formulations, or use different models.");
169  }
170 }
bool _is_elasticity_tensor_guaranteed_isotropic
is the elasticity tensor guaranteed to be isotropic?
StressUpdateBase is a material that is not called by MOOSE because of the compute=false flag set in t...
const std::string _elasticity_tensor_name
Name of the elasticity tensor material property.
std::vector< bool > _tangent_computation_flag
Flags to compute tangent during updateState call.
virtual bool requiresIsotropicTensor()=0
Does the model require the elasticity tensor to be isotropic?
TangentCalculationMethod _tangent_calculation_method
Calculation method for the tangent modulus.
const unsigned _num_models
number of plastic models
enum ComputeMultipleInelasticStress::TangentOperatorEnum _tangent_operator_type
std::vector< StressUpdateBase * > _models
The user supplied list of inelastic models to use in the simulation.
bool hasGuaranteedMaterialProperty(const MaterialPropertyName &prop, Guarantee guarantee)

◆ initQpStatefulProperties()

void ComputeMultipleInelasticStress::initQpStatefulProperties ( )
overrideprotectedvirtual

Reimplemented from ComputeStressBase.

Reimplemented in ComputeMultipleInelasticCosseratStress, and ComputeSmearedCrackingStress.

Definition at line 107 of file ComputeMultipleInelasticStress.C.

Referenced by ComputeSmearedCrackingStress::initQpStatefulProperties(), and ComputeMultipleInelasticCosseratStress::initQpStatefulProperties().

108 {
110  _inelastic_strain[_qp].zero();
111 }
virtual void initQpStatefulProperties() override
MaterialProperty< RankTwoTensor > & _inelastic_strain
The sum of the inelastic strains that come from the plastic models.

◆ updateQpState()

void ComputeMultipleInelasticStress::updateQpState ( RankTwoTensor elastic_strain_increment,
RankTwoTensor combined_inelastic_strain_increment 
)
protectedvirtual

Given the _strain_increment[_qp], iterate over all of the user-specified recompute materials in order to find an admissible stress (which is placed into _stress[_qp]) and set of inelastic strains, as well as the tangent operator (which is placed into _Jacobian_mult[_qp]).

Parameters
elastic_strain_incrementThe elastic part of _strain_increment[_qp] after the iterative process has converged
combined_inelastic_strain_incrementThe inelastic part of _strain_increment[_qp] after the iterative process has converged. This is a weighted sum of all the inelastic strains computed by all the plastic models during the Picard iterative scheme. The weights are dictated by the user using _inelastic_weights

Definition at line 229 of file ComputeMultipleInelasticStress.C.

Referenced by computeQpStressIntermediateConfiguration().

231 {
233  {
234  _console << std::endl
235  << "iteration output for ComputeMultipleInelasticStress solve:"
236  << " time=" << _t << " int_pt=" << _qp << std::endl;
237  }
238  Real l2norm_delta_stress;
239  Real first_l2norm_delta_stress = 1.0;
240  unsigned int counter = 0;
241 
242  std::vector<RankTwoTensor> inelastic_strain_increment;
243  inelastic_strain_increment.resize(_num_models);
244 
245  for (unsigned i_rmm = 0; i_rmm < _models.size(); ++i_rmm)
246  inelastic_strain_increment[i_rmm].zero();
247 
248  RankTwoTensor stress_max, stress_min;
249 
250  do
251  {
252  for (unsigned i_rmm = 0; i_rmm < _num_models; ++i_rmm)
253  {
254  _models[i_rmm]->setQp(_qp);
255 
256  // initially assume the strain is completely elastic
257  elastic_strain_increment = _strain_increment[_qp];
258  // and subtract off all inelastic strain increments calculated so far
259  // except the one that we're about to calculate
260  for (unsigned j_rmm = 0; j_rmm < _num_models; ++j_rmm)
261  if (i_rmm != j_rmm)
262  elastic_strain_increment -= inelastic_strain_increment[j_rmm];
263 
264  // form the trial stress, with the check for changed elasticity constants
266  _stress[_qp] =
267  _elasticity_tensor[_qp] * (_elastic_strain_old[_qp] + elastic_strain_increment);
268  else
269  _stress[_qp] = _stress_old[_qp] + _elasticity_tensor[_qp] * elastic_strain_increment;
270 
271  // given a trial stress (_stress[_qp]) and a strain increment (elastic_strain_increment)
272  // let the i^th model produce an admissible stress (as _stress[_qp]), and decompose
273  // the strain increment into an elastic part (elastic_strain_increment) and an
274  // inelastic part (inelastic_strain_increment[i_rmm])
276  elastic_strain_increment,
277  inelastic_strain_increment[i_rmm],
279 
280  if (i_rmm == 0)
281  {
282  stress_max = _stress[_qp];
283  stress_min = _stress[_qp];
284  }
285  else
286  {
287  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
288  {
289  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
290  {
291  if (_stress[_qp](i, j) > stress_max(i, j))
292  stress_max(i, j) = _stress[_qp](i, j);
293  else if (stress_min(i, j) > _stress[_qp](i, j))
294  stress_min(i, j) = _stress[_qp](i, j);
295  }
296  }
297  }
298  }
299 
300  // now check convergence in the stress:
301  // once the change in stress is within tolerance after each recompute material
302  // consider the stress to be converged
303  l2norm_delta_stress = (stress_max - stress_min).L2norm();
304  if (counter == 0 && l2norm_delta_stress > 0.0)
305  first_l2norm_delta_stress = l2norm_delta_stress;
306 
308  {
309  _console << "stress iteration number = " << counter << "\n"
310  << " relative l2 norm delta stress = "
311  << (0 == first_l2norm_delta_stress ? 0
312  : l2norm_delta_stress / first_l2norm_delta_stress)
313  << "\n"
314  << " stress convergence relative tolerance = " << _relative_tolerance << "\n"
315  << " absolute l2 norm delta stress = " << l2norm_delta_stress << "\n"
316  << " stress convergence absolute tolerance = " << _absolute_tolerance << std::endl;
317  }
318  ++counter;
319  } while (counter < _max_iterations && l2norm_delta_stress > _absolute_tolerance &&
320  (l2norm_delta_stress / first_l2norm_delta_stress) > _relative_tolerance &&
321  _num_models != 1);
322 
323  if (counter == _max_iterations && l2norm_delta_stress > _absolute_tolerance &&
324  (l2norm_delta_stress / first_l2norm_delta_stress) > _relative_tolerance)
325  throw MooseException("Max stress iteration hit during ComputeMultipleInelasticStress solve!");
326 
327  combined_inelastic_strain_increment.zero();
328  for (unsigned i_rmm = 0; i_rmm < _num_models; ++i_rmm)
329  combined_inelastic_strain_increment +=
330  _inelastic_weights[i_rmm] * inelastic_strain_increment[i_rmm];
331 
332  if (_fe_problem.currentlyComputingJacobian())
334 
335  _matl_timestep_limit[_qp] = 0.0;
336  for (unsigned i_rmm = 0; i_rmm < _num_models; ++i_rmm)
337  _matl_timestep_limit[_qp] += 1.0 / _models[i_rmm]->computeTimeStepLimit();
338 
339  if (MooseUtils::absoluteFuzzyEqual(_matl_timestep_limit[_qp], 0.0))
340  {
341  _matl_timestep_limit[_qp] = std::numeric_limits<Real>::max();
342  }
343  else
344  {
345  _matl_timestep_limit[_qp] = 1.0 / _matl_timestep_limit[_qp];
346  }
347 }
bool _is_elasticity_tensor_guaranteed_isotropic
is the elasticity tensor guaranteed to be isotropic?
std::vector< RankFourTensor > _consistent_tangent_operator
the consistent tangent operators computed by each plastic model
const MaterialProperty< RankTwoTensor > & _elastic_strain_old
Strain tensors.
MaterialProperty< Real > & _matl_timestep_limit
const MaterialProperty< RankTwoTensor > & _strain_increment
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Elasticity tensor material property.
const std::vector< Real > _inelastic_weights
_inelastic_strain = sum_i (_inelastic_weights_i * inelastic_strain_from_model_i)
MaterialProperty< RankTwoTensor > & _stress
Stress material property.
const unsigned _num_models
number of plastic models
virtual void computeAdmissibleState(unsigned model_number, RankTwoTensor &elastic_strain_increment, RankTwoTensor &inelastic_strain_increment, RankFourTensor &consistent_tangent_operator)
Given a trial stress (_stress[_qp]) and a strain increment (elastic_strain_increment) let the model_n...
const MaterialProperty< RankTwoTensor > & _stress_old
Old state of the stress tensor material property.
std::vector< StressUpdateBase * > _models
The user supplied list of inelastic models to use in the simulation.
Real L2norm(const RankTwoTensor &r2tensor)
const bool _perform_finite_strain_rotations
after updateQpState, rotate the stress, elastic_strain, inelastic_strain and Jacobian_mult using _rot...
virtual void computeQpJacobianMult()
Using _elasticity_tensor[_qp] and the consistent tangent operators, _consistent_tangent_operator[...] computed by the inelastic models, compute _Jacobian_mult[_qp].
static unsigned int counter
const unsigned int _max_iterations
Input parameters associated with the recompute iteration to return the stress state to the yield surf...

◆ updateQpStateSingleModel()

void ComputeMultipleInelasticStress::updateQpStateSingleModel ( unsigned  model_number,
RankTwoTensor elastic_strain_increment,
RankTwoTensor combined_inelastic_strain_increment 
)
protectedvirtual

An optimised version of updateQpState that gets used when the number of plastic models is unity, or when we're cycling through models Given the _strain_increment[_qp], find an admissible stress (which is put into _stress[_qp]) and inelastic strain, as well as the tangent operator (which is placed into _Jacobian_mult[_qp])

Parameters
model_numberUse this model number
elastic_strain_incrementThe elastic part of _strain_increment[_qp]
combined_inelastic_strain_incrementThe inelastic part of _strain_increment[_qp]

Definition at line 372 of file ComputeMultipleInelasticStress.C.

Referenced by computeQpStressIntermediateConfiguration().

376 {
377  for (auto model : _models)
378  model->setQp(_qp);
379 
380  elastic_strain_increment = _strain_increment[_qp];
381 
382  // If the elasticity tensor values have changed and the tensor is isotropic,
383  // use the old strain to calculate the old stress
385  _stress[_qp] = _elasticity_tensor[_qp] * (_elastic_strain_old[_qp] + elastic_strain_increment);
386  else
387  _stress[_qp] = _stress_old[_qp] + _elasticity_tensor[_qp] * elastic_strain_increment;
388 
389  computeAdmissibleState(model_number,
390  elastic_strain_increment,
391  combined_inelastic_strain_increment,
393 
394  if (_fe_problem.currentlyComputingJacobian())
395  {
399  {
401  mooseAssert(A.isSymmetric(), "Tangent operator isn't symmetric");
402  _Jacobian_mult[_qp] = A.invSymm() * _elasticity_tensor[_qp];
403  }
404  else
406  }
407 
408  _matl_timestep_limit[_qp] = _models[0]->computeTimeStepLimit();
409 
410  /* propagate internal variables, etc, to this timestep for those inelastic models where
411  * "updateState" is not called */
412  for (unsigned i_rmm = 0; i_rmm < _num_models; ++i_rmm)
413  if (i_rmm != model_number)
414  _models[i_rmm]->propagateQpStatefulProperties();
415 }
bool _is_elasticity_tensor_guaranteed_isotropic
is the elasticity tensor guaranteed to be isotropic?
MaterialProperty< RankFourTensor > & _Jacobian_mult
derivative of stress w.r.t. strain (_dstress_dstrain)
std::vector< RankFourTensor > _consistent_tangent_operator
the consistent tangent operators computed by each plastic model
const MaterialProperty< RankTwoTensor > & _elastic_strain_old
Strain tensors.
MaterialProperty< Real > & _matl_timestep_limit
const MaterialProperty< RankTwoTensor > & _strain_increment
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Elasticity tensor material property.
MaterialProperty< RankTwoTensor > & _stress
Stress material property.
const RankFourTensor _identity_symmetric_four
Rank four symmetric identity tensor.
TangentCalculationMethod _tangent_calculation_method
Calculation method for the tangent modulus.
const unsigned _num_models
number of plastic models
virtual void computeAdmissibleState(unsigned model_number, RankTwoTensor &elastic_strain_increment, RankTwoTensor &inelastic_strain_increment, RankFourTensor &consistent_tangent_operator)
Given a trial stress (_stress[_qp]) and a strain increment (elastic_strain_increment) let the model_n...
const MaterialProperty< RankTwoTensor > & _stress_old
Old state of the stress tensor material property.
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...

Member Data Documentation

◆ _absolute_tolerance

const Real ComputeMultipleInelasticStress::_absolute_tolerance
protected

Definition at line 122 of file ComputeMultipleInelasticStress.h.

Referenced by updateQpState().

◆ _base_name

const std::string ComputeStressBase::_base_name
protectedinherited

Base name prepended to all material property names to allow for multi-material systems.

Definition at line 44 of file ComputeStressBase.h.

Referenced by ComputeLinearElasticStress::initialSetup(), and ComputeCosseratLinearElasticStress::initialSetup().

◆ _consistent_tangent_operator

std::vector<RankFourTensor> ComputeMultipleInelasticStress::_consistent_tangent_operator
protected

the consistent tangent operators computed by each plastic model

Definition at line 156 of file ComputeMultipleInelasticStress.h.

Referenced by ComputeMultipleInelasticCosseratStress::computeQpJacobianMult(), computeQpJacobianMult(), updateQpState(), and updateQpStateSingleModel().

◆ _cycle_models

const bool ComputeMultipleInelasticStress::_cycle_models
protected

whether to cycle through the models, using only one model per timestep

Definition at line 159 of file ComputeMultipleInelasticStress.h.

Referenced by computeQpStressIntermediateConfiguration().

◆ _elastic_strain

MaterialProperty<RankTwoTensor>& ComputeStressBase::_elastic_strain
protectedinherited

◆ _elastic_strain_old

const MaterialProperty<RankTwoTensor>& ComputeMultipleInelasticStress::_elastic_strain_old
protected

◆ _elasticity_tensor

const MaterialProperty<RankFourTensor>& ComputeFiniteStrainElasticStress::_elasticity_tensor
protectedinherited

◆ _elasticity_tensor_name

const std::string ComputeFiniteStrainElasticStress::_elasticity_tensor_name
protectedinherited

Name of the elasticity tensor material property.

Definition at line 36 of file ComputeFiniteStrainElasticStress.h.

Referenced by ComputeFiniteStrainElasticStress::initialSetup(), ComputeSmearedCrackingStress::initialSetup(), and initialSetup().

◆ _extra_stress

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_extra_stress
protectedinherited

Extra stress tensor.

Definition at line 54 of file ComputeStressBase.h.

Referenced by ComputeStressBase::computeQpProperties().

◆ _identity_symmetric_four

const RankFourTensor ComputeMultipleInelasticStress::_identity_symmetric_four
protected

Rank four symmetric identity tensor.

Definition at line 166 of file ComputeMultipleInelasticStress.h.

Referenced by computeQpJacobianMult(), and updateQpStateSingleModel().

◆ _inelastic_strain

MaterialProperty<RankTwoTensor>& ComputeMultipleInelasticStress::_inelastic_strain
protected

The sum of the inelastic strains that come from the plastic models.

Definition at line 135 of file ComputeMultipleInelasticStress.h.

Referenced by ComputeSmearedCrackingStress::computeQpStress(), computeQpStressIntermediateConfiguration(), finiteStrainRotation(), and initQpStatefulProperties().

◆ _inelastic_strain_old

const MaterialProperty<RankTwoTensor>& ComputeMultipleInelasticStress::_inelastic_strain_old
protected

old value of inelastic strain

Definition at line 138 of file ComputeMultipleInelasticStress.h.

Referenced by ComputeSmearedCrackingStress::computeQpStress(), and computeQpStressIntermediateConfiguration().

◆ _inelastic_weights

const std::vector<Real> ComputeMultipleInelasticStress::_inelastic_weights
protected

_inelastic_strain = sum_i (_inelastic_weights_i * inelastic_strain_from_model_i)

Definition at line 153 of file ComputeMultipleInelasticStress.h.

Referenced by ComputeMultipleInelasticStress(), and updateQpState().

◆ _initial_stress_fcn

std::vector<Function *> ComputeStressBase::_initial_stress_fcn
protectedinherited

initial stress components

Definition at line 57 of file ComputeStressBase.h.

◆ _internal_solve_full_iteration_history

const bool ComputeMultipleInelasticStress::_internal_solve_full_iteration_history
protected

Definition at line 123 of file ComputeMultipleInelasticStress.h.

Referenced by updateQpState().

◆ _is_elasticity_tensor_guaranteed_isotropic

bool ComputeMultipleInelasticStress::_is_elasticity_tensor_guaranteed_isotropic
protected

is the elasticity tensor guaranteed to be isotropic?

Definition at line 178 of file ComputeMultipleInelasticStress.h.

Referenced by computeQpStressIntermediateConfiguration(), finiteStrainRotation(), initialSetup(), updateQpState(), and updateQpStateSingleModel().

◆ _Jacobian_mult

MaterialProperty<RankFourTensor>& ComputeStressBase::_Jacobian_mult
protectedinherited

◆ _matl_timestep_limit

MaterialProperty<Real>& ComputeMultipleInelasticStress::_matl_timestep_limit
protected

◆ _max_iterations

const unsigned int ComputeMultipleInelasticStress::_max_iterations
protected

Input parameters associated with the recompute iteration to return the stress state to the yield surface.

Definition at line 120 of file ComputeMultipleInelasticStress.h.

Referenced by updateQpState().

◆ _mechanical_strain

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_mechanical_strain
protectedinherited

◆ _models

std::vector<StressUpdateBase *> ComputeMultipleInelasticStress::_models
protected

The user supplied list of inelastic models to use in the simulation.

Users should take care to list creep models first and plasticity models last to allow for the case when a creep model relaxes the stress state inside of the yield surface in an iteration.

Definition at line 175 of file ComputeMultipleInelasticStress.h.

Referenced by computeAdmissibleState(), ComputeSmearedCrackingStress::computeQpStress(), initialSetup(), updateQpState(), and updateQpStateSingleModel().

◆ _num_models

const unsigned ComputeMultipleInelasticStress::_num_models
protected

◆ _perform_finite_strain_rotations

const bool ComputeMultipleInelasticStress::_perform_finite_strain_rotations
protected

after updateQpState, rotate the stress, elastic_strain, inelastic_strain and Jacobian_mult using _rotation_increment

Definition at line 127 of file ComputeMultipleInelasticStress.h.

Referenced by ComputeSmearedCrackingStress::computeQpStress(), ComputeMultipleInelasticCosseratStress::computeQpStress(), computeQpStress(), computeQpStressIntermediateConfiguration(), updateQpState(), and updateQpStateSingleModel().

◆ _relative_tolerance

const Real ComputeMultipleInelasticStress::_relative_tolerance
protected

Definition at line 121 of file ComputeMultipleInelasticStress.h.

Referenced by updateQpState().

◆ _rotation_increment

const MaterialProperty<RankTwoTensor>& ComputeFiniteStrainElasticStress::_rotation_increment
protectedinherited

◆ _strain_increment

const MaterialProperty<RankTwoTensor>& ComputeMultipleInelasticStress::_strain_increment
protected

◆ _stress

MaterialProperty<RankTwoTensor>& ComputeStressBase::_stress
protectedinherited

Stress material property.

Definition at line 49 of file ComputeStressBase.h.

Referenced by ComputeMultipleInelasticCosseratStress::computeAdmissibleState(), computeAdmissibleState(), ComputeStressBase::computeQpProperties(), ComputeStrainIncrementBasedStress::computeQpStress(), ComputeLinearElasticStress::computeQpStress(), ComputeCosseratLinearElasticStress::computeQpStress(), ComputeDamageStress::computeQpStress(), ComputeFiniteStrainElasticStress::computeQpStress(), ComputeSmearedCrackingStress::computeQpStress(), FiniteStrainPlasticMaterial::computeQpStress(), ComputeLinearElasticPFFractureStress::computeQpStress(), ComputeMultiPlasticityStress::computeQpStress(), ComputeLinearViscoelasticStress::computeQpStress(), computeQpStressIntermediateConfiguration(), ComputeLinearElasticPFFractureStress::computeStrainSpectral(), ComputeLinearElasticPFFractureStress::computeStrainVolDev(), ComputeLinearElasticPFFractureStress::computeStressSpectral(), finiteStrainRotation(), ComputeStressBase::initQpStatefulProperties(), FiniteStrainCrystalPlasticity::initQpStatefulProperties(), FiniteStrainUObasedCP::initQpStatefulProperties(), FiniteStrainHyperElasticViscoPlastic::initQpStatefulProperties(), ComputeMultiPlasticityStress::postReturnMap(), FiniteStrainUObasedCP::postSolveQp(), FiniteStrainHyperElasticViscoPlastic::postSolveQp(), FiniteStrainCrystalPlasticity::postSolveQp(), ComputeSmearedCrackingStress::updateCrackingStateAndStress(), updateQpState(), and updateQpStateSingleModel().

◆ _stress_old

const MaterialProperty<RankTwoTensor>& ComputeFiniteStrainElasticStress::_stress_old
protectedinherited

Old state of the stress tensor material property.

Definition at line 44 of file ComputeFiniteStrainElasticStress.h.

Referenced by computeAdmissibleState(), computeQpStressIntermediateConfiguration(), updateQpState(), and updateQpStateSingleModel().

◆ _tangent_calculation_method

TangentCalculationMethod ComputeMultipleInelasticStress::_tangent_calculation_method
protected

Calculation method for the tangent modulus.

Definition at line 150 of file ComputeMultipleInelasticStress.h.

Referenced by computeAdmissibleState(), computeQpJacobianMult(), finiteStrainRotation(), initialSetup(), and updateQpStateSingleModel().

◆ _tangent_computation_flag

std::vector<bool> ComputeMultipleInelasticStress::_tangent_computation_flag
protected

Flags to compute tangent during updateState call.

Definition at line 147 of file ComputeMultipleInelasticStress.h.

Referenced by computeAdmissibleState(), and initialSetup().

◆ _tangent_operator_type

enum ComputeMultipleInelasticStress::TangentOperatorEnum ComputeMultipleInelasticStress::_tangent_operator_type
protected

The documentation for this class was generated from the following files: