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)
 

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 initialSetup () override
 
virtual void computeQpStress () override
 
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, _comsistent_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 rotateQpInitialStress ()
 InitialStress Deprecation: remove this method. More...
 
virtual void computeQpProperties () override
 
void addQpInitialStress ()
 InitialStress Deprecation: remove this method. More...
 
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...
 
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
 
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 MaterialProperty< RankTwoTensor > & _rotation_increment
 
const MaterialProperty< RankTwoTensor > & _stress_old
 
const std::string _base_name
 
const std::string _elasticity_tensor_name
 
const MaterialProperty< RankTwoTensor > & _mechanical_strain
 
MaterialProperty< RankTwoTensor > & _stress
 
MaterialProperty< RankTwoTensor > & _elastic_strain
 
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 bool _store_stress_old
 Parameter which decides whether to store old stress. This is required for HHT time integration and Rayleigh damping. More...
 
const bool _initial_stress_provided
 Whether initial stress was provided. InitialStress Deprecation: remove this. More...
 
MaterialProperty< RankTwoTensor > * _initial_stress
 Initial stress, if provided. InitialStress Deprecation: remove this. More...
 
const MaterialProperty< RankTwoTensor > * _initial_stress_old
 Old value of initial stress, which is needed to correctly implement finite-strain rotations. InitialStress Deprecation: remove this. 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 _output_iteration_info
 
const MaterialProperty< RankFourTensor > & _elasticity_tensor
 Rank-4 and Rank-2 elasticity and elastic strain tensors. More...
 
const MaterialProperty< RankTwoTensor > & _elastic_strain_old
 
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 30 of file ComputeMultipleInelasticStress.h.

Member Enumeration Documentation

what sort of Tangent operator to calculate

Enumerator
elastic 
nonlinear 

Definition at line 134 of file ComputeMultipleInelasticStress.h.

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

Constructor & Destructor Documentation

ComputeMultipleInelasticStress::ComputeMultipleInelasticStress ( const InputParameters &  parameters)

Definition at line 68 of file ComputeMultipleInelasticStress.C.

70  _max_iterations(parameters.get<unsigned int>("max_iterations")),
71  _relative_tolerance(parameters.get<Real>("relative_tolerance")),
72  _absolute_tolerance(parameters.get<Real>("absolute_tolerance")),
73  _output_iteration_info(getParam<bool>("output_iteration_info")),
74  _perform_finite_strain_rotations(getParam<bool>("perform_finite_strain_rotations")),
75  _elasticity_tensor(getMaterialPropertyByName<RankFourTensor>(_base_name + "elasticity_tensor")),
76  _elastic_strain_old(getMaterialPropertyOld<RankTwoTensor>(_base_name + "elastic_strain")),
77  _strain_increment(getMaterialProperty<RankTwoTensor>(_base_name + "strain_increment")),
78  _inelastic_strain(declareProperty<RankTwoTensor>(_base_name + "combined_inelastic_strain")),
80  getMaterialPropertyOld<RankTwoTensor>(_base_name + "combined_inelastic_strain")),
81  _tangent_operator_type(getParam<MooseEnum>("tangent_operator").getEnum<TangentOperatorEnum>()),
82  _num_models(getParam<std::vector<MaterialName>>("inelastic_models").size()),
83  _inelastic_weights(isParamValid("combined_inelastic_strain_weights")
84  ? getParam<std::vector<Real>>("combined_inelastic_strain_weights")
85  : std::vector<Real>(_num_models, true)),
87  _cycle_models(getParam<bool>("cycle_models")),
88  _matl_timestep_limit(declareProperty<Real>("matl_timestep_limit"))
89 {
90  if (_inelastic_weights.size() != _num_models)
91  mooseError(
92  "ComputeMultipleInelasticStress: combined_inelastic_strain_weights must contain the same "
93  "number of entries as inelastic_models ",
94  _inelastic_weights.size(),
95  " ",
96  _num_models);
97 }
std::vector< RankFourTensor > _consistent_tangent_operator
the consistent tangent operators computed by each plastic model
const MaterialProperty< RankTwoTensor > & _elastic_strain_old
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
const bool _cycle_models
whether to cycle through the models, using only one model per timestep
const unsigned _num_models
number of plastic models
ComputeFiniteStrainElasticStress(const InputParameters &parameters)
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Rank-4 and Rank-2 elasticity and elastic strain tensors.
const std::string _base_name
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

void ComputeStressBase::addQpInitialStress ( )
protectedinherited

InitialStress Deprecation: remove this method.

Adds initial stress, if it is provided, to _stress[_qp]. This function should NOT be used if you calculate stress using

stress = stress_old + elasticity * strain_increment

because stress_old will already include initial stress. However this function SHOULD be used if your code uses

stress = elasticity * (elastic_strain_old + strain_increment) or stress = elasticity * elastic_strain

since in these cases the elastic_strain and elastic_strain_old will not include any contribution from initial stress.

Definition at line 112 of file ComputeStressBase.C.

Referenced by ComputeLinearElasticStress::computeQpStress(), ComputeCosseratLinearElasticStress::computeQpStress(), ComputeFiniteStrainElasticStress::computeQpStress(), ComputeSmearedCrackingStress::computeQpStress(), computeQpStressIntermediateConfiguration(), ComputeStressBase::initQpStatefulProperties(), updateQpState(), and updateQpStateSingleModel().

113 {
115  _stress[_qp] += (*_initial_stress)[_qp];
116 }
MaterialProperty< RankTwoTensor > & _stress
const bool _initial_stress_provided
Whether initial stress was provided. InitialStress Deprecation: remove this.
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 364 of file ComputeMultipleInelasticStress.C.

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

368 {
369  _models[model_number]->updateState(elastic_strain_increment,
370  inelastic_strain_increment,
371  _rotation_increment[_qp],
372  _stress[_qp],
373  _stress_old[_qp],
374  _elasticity_tensor[_qp],
375  _elastic_strain_old[_qp],
377  consistent_tangent_operator);
378 }
const MaterialProperty< RankTwoTensor > & _rotation_increment
const MaterialProperty< RankTwoTensor > & _elastic_strain_old
MaterialProperty< RankTwoTensor > & _stress
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Rank-4 and Rank-2 elasticity and elastic strain tensors.
enum ComputeMultipleInelasticStress::TangentOperatorEnum _tangent_operator_type
const MaterialProperty< RankTwoTensor > & _stress_old
std::vector< StressUpdateBase * > _models
The user supplied list of inelastic models to use in the simulation.
void ComputeMultipleInelasticStress::computeQpJacobianMult ( )
protectedvirtual

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

Reimplemented in ComputeMultipleInelasticCosseratStress.

Definition at line 312 of file ComputeMultipleInelasticStress.C.

Referenced by updateQpState().

313 {
316  else
317  {
318  const RankFourTensor E_inv = _elasticity_tensor[_qp].invSymm();
320  for (unsigned i_rmm = 1; i_rmm < _num_models; ++i_rmm)
321  _Jacobian_mult[_qp] = _consistent_tangent_operator[i_rmm] * E_inv * _Jacobian_mult[_qp];
322  }
323 }
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 unsigned _num_models
number of plastic models
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Rank-4 and Rank-2 elasticity and elastic strain tensors.
enum ComputeMultipleInelasticStress::TangentOperatorEnum _tangent_operator_type
void ComputeStressBase::computeQpProperties ( )
overrideprotectedvirtualinherited

Reimplemented in ComputeBirchMurnaghanEquationOfStress, and ComputeStressEosBase.

Definition at line 99 of file ComputeStressBase.C.

100 {
101  // InitialStress Deprecation: remove the following 2 lines
103  (*_initial_stress)[_qp] = (*_initial_stress_old)[_qp];
104 
105  computeQpStress();
106 
107  // Add in extra stress
108  _stress[_qp] += _extra_stress[_qp];
109 }
virtual void computeQpStress()=0
MaterialProperty< RankTwoTensor > & _stress
const bool _initial_stress_provided
Whether initial stress was provided. InitialStress Deprecation: remove this.
const MaterialProperty< RankTwoTensor > & _extra_stress
Extra stress tensor.
void ComputeMultipleInelasticStress::computeQpStress ( )
overrideprotectedvirtual

Reimplemented from ComputeFiniteStrainElasticStress.

Reimplemented in ComputeMultipleInelasticCosseratStress, and ComputeSmearedCrackingStress.

Definition at line 130 of file ComputeMultipleInelasticStress.C.

Referenced by ComputeMultipleInelasticCosseratStress::computeQpStress().

131 {
135 }
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...
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 138 of file ComputeMultipleInelasticStress.C.

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

139 {
140  RankTwoTensor elastic_strain_increment;
141  RankTwoTensor combined_inelastic_strain_increment;
142 
143  if (_num_models == 0)
144  {
146 
147  // If the elasticity tensor values have changed and the tensor is isotropic,
148  // use the old strain to calculate the old stress
150  {
151  _stress[_qp] = _elasticity_tensor[_qp] * (_elastic_strain_old[_qp] + _strain_increment[_qp]);
152  // InitialStress Deprecation: remove these lines
156  }
157  else
158  _stress[_qp] = _stress_old[_qp] + _elasticity_tensor[_qp] * _strain_increment[_qp];
159 
160  if (_fe_problem.currentlyComputingJacobian())
162  }
163  else
164  {
165  if (_num_models == 1 || _cycle_models)
166  updateQpStateSingleModel((_t_step - 1) % _num_models,
167  elastic_strain_increment,
168  combined_inelastic_strain_increment);
169  else
170  updateQpState(elastic_strain_increment, combined_inelastic_strain_increment);
171 
172  _elastic_strain[_qp] = _elastic_strain_old[_qp] + elastic_strain_increment;
173  _inelastic_strain[_qp] = _inelastic_strain_old[_qp] + combined_inelastic_strain_increment;
174  }
175 }
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
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...
MaterialProperty< RankTwoTensor > & _stress
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< RankFourTensor > & _elasticity_tensor
Rank-4 and Rank-2 elasticity and elastic strain tensors.
void addQpInitialStress()
InitialStress Deprecation: remove this method.
const MaterialProperty< RankTwoTensor > & _stress_old
const bool _perform_finite_strain_rotations
after updateQpState, rotate the stress, elastic_strain, inelastic_strain and Jacobian_mult using _rot...
MaterialProperty< RankTwoTensor > & _elastic_strain
virtual void rotateQpInitialStress()
InitialStress Deprecation: remove this method.
MaterialProperty< RankTwoTensor > & _inelastic_strain
The sum of the inelastic strains that come from the plastic models.
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 178 of file ComputeMultipleInelasticStress.C.

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

179 {
180  _elastic_strain[_qp] =
181  _rotation_increment[_qp] * _elastic_strain[_qp] * _rotation_increment[_qp].transpose();
182  _stress[_qp] = _rotation_increment[_qp] * _stress[_qp] * _rotation_increment[_qp].transpose();
183  _inelastic_strain[_qp] =
184  _rotation_increment[_qp] * _inelastic_strain[_qp] * _rotation_increment[_qp].transpose();
185  if (force_elasticity_rotation ||
188  _Jacobian_mult[_qp].rotate(_rotation_increment[_qp]);
189 }
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
MaterialProperty< RankTwoTensor > & _stress
const unsigned _num_models
number of plastic models
enum ComputeMultipleInelasticStress::TangentOperatorEnum _tangent_operator_type
MaterialProperty< RankTwoTensor > & _elastic_strain
MaterialProperty< RankTwoTensor > & _inelastic_strain
The sum of the inelastic strains that come from the plastic models.
bool GuaranteeConsumer::hasGuaranteedMaterialProperty ( const MaterialPropertyName &  prop,
Guarantee  guarantee 
)
protectedinherited

Definition at line 26 of file GuaranteeConsumer.C.

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

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

Reimplemented in ComputeSmearedCrackingStress.

Definition at line 107 of file ComputeMultipleInelasticStress.C.

Referenced by ComputeSmearedCrackingStress::initialSetup().

108 {
111 
112  std::vector<MaterialName> models = getParam<std::vector<MaterialName>>("inelastic_models");
113  for (unsigned int i = 0; i < _num_models; ++i)
114  {
115  StressUpdateBase * rrr = dynamic_cast<StressUpdateBase *>(&getMaterialByName(models[i]));
116  if (rrr)
117  {
118  _models.push_back(rrr);
120  mooseError("Model " + models[i] +
121  " requires an isotropic elasticity tensor, but the one supplied is not "
122  "guaranteed isotropic");
123  }
124  else
125  mooseError("Model " + models[i] + " is not compatible with ComputeMultipleInelasticStress");
126  }
127 }
bool _is_elasticity_tensor_guaranteed_isotropic
is the elasticity tensor guaranteed to be isotropic?
const std::string _elasticity_tensor_name
StressUpdateBase is a material that is not called by MOOSE because of the compute=false flag set in t...
virtual bool requiresIsotropicTensor()=0
Does the model require the elasticity tensor to be isotropic?
const unsigned _num_models
number of plastic models
std::vector< StressUpdateBase * > _models
The user supplied list of inelastic models to use in the simulation.
bool hasGuaranteedMaterialProperty(const MaterialPropertyName &prop, Guarantee guarantee)
void ComputeMultipleInelasticStress::initQpStatefulProperties ( )
overrideprotectedvirtual

Reimplemented from ComputeFiniteStrainElasticStress.

Reimplemented in ComputeMultipleInelasticCosseratStress, and ComputeSmearedCrackingStress.

Definition at line 100 of file ComputeMultipleInelasticStress.C.

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

101 {
103  _inelastic_strain[_qp].zero();
104 }
virtual void initQpStatefulProperties() override
MaterialProperty< RankTwoTensor > & _inelastic_strain
The sum of the inelastic strains that come from the plastic models.
void ComputeFiniteStrainElasticStress::rotateQpInitialStress ( )
protectedvirtualinherited

InitialStress Deprecation: remove this method.

Rotates initial_stress via rotation_increment. In large-strain scenarios this must be used before addQpInitialStress

Definition at line 69 of file ComputeFiniteStrainElasticStress.C.

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

70 {
72  (*_initial_stress)[_qp] = _rotation_increment[_qp] * (*_initial_stress_old)[_qp] *
73  _rotation_increment[_qp].transpose();
74 }
const MaterialProperty< RankTwoTensor > & _rotation_increment
const bool _initial_stress_provided
Whether initial stress was provided. InitialStress Deprecation: remove this.
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 192 of file ComputeMultipleInelasticStress.C.

Referenced by computeQpStressIntermediateConfiguration().

194 {
195  if (_output_iteration_info == true)
196  {
197  _console << std::endl
198  << "iteration output for ComputeMultipleInelasticStress solve:"
199  << " time=" << _t << " int_pt=" << _qp << std::endl;
200  }
201  Real l2norm_delta_stress, first_l2norm_delta_stress;
202  unsigned int counter = 0;
203 
204  std::vector<RankTwoTensor> inelastic_strain_increment;
205  inelastic_strain_increment.resize(_num_models);
206 
207  for (unsigned i_rmm = 0; i_rmm < _models.size(); ++i_rmm)
208  inelastic_strain_increment[i_rmm].zero();
209 
210  RankTwoTensor stress_max, stress_min;
211 
212  do
213  {
214  for (unsigned i_rmm = 0; i_rmm < _num_models; ++i_rmm)
215  {
216  _models[i_rmm]->setQp(_qp);
217 
218  // initially assume the strain is completely elastic
219  elastic_strain_increment = _strain_increment[_qp];
220  // and subtract off all inelastic strain increments calculated so far
221  // except the one that we're about to calculate
222  for (unsigned j_rmm = 0; j_rmm < _num_models; ++j_rmm)
223  if (i_rmm != j_rmm)
224  elastic_strain_increment -= inelastic_strain_increment[j_rmm];
225 
226  // form the trial stress, with the check for changed elasticity constants
228  {
229  _stress[_qp] =
230  _elasticity_tensor[_qp] * (_elastic_strain_old[_qp] + elastic_strain_increment);
231  // InitialStress Deprecation: remove these lines
235  }
236  else
237  _stress[_qp] = _stress_old[_qp] + _elasticity_tensor[_qp] * elastic_strain_increment;
238 
239  // given a trial stress (_stress[_qp]) and a strain increment (elastic_strain_increment)
240  // let the i^th model produce an admissible stress (as _stress[_qp]), and decompose
241  // the strain increment into an elastic part (elastic_strain_increment) and an
242  // inelastic part (inelastic_strain_increment[i_rmm])
244  elastic_strain_increment,
245  inelastic_strain_increment[i_rmm],
247 
248  if (i_rmm == 0)
249  {
250  stress_max = _stress[_qp];
251  stress_min = _stress[_qp];
252  }
253  else
254  {
255  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
256  {
257  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
258  {
259  if (_stress[_qp](i, j) > stress_max(i, j))
260  stress_max(i, j) = _stress[_qp](i, j);
261  else if (stress_min(i, j) > _stress[_qp](i, j))
262  stress_min(i, j) = _stress[_qp](i, j);
263  }
264  }
265  }
266  }
267 
268  // now check convergence in the stress:
269  // once the change in stress is within tolerance after each recompute material
270  // consider the stress to be converged
271  RankTwoTensor delta_stress(stress_max - stress_min);
272  l2norm_delta_stress = delta_stress.L2norm();
273  if (counter == 0)
274  first_l2norm_delta_stress = l2norm_delta_stress;
275 
276  if (_output_iteration_info == true)
277  {
278  _console << "stress iteration number = " << counter << "\n"
279  << " relative l2 norm delta stress = "
280  << (0 == first_l2norm_delta_stress ? 0
281  : l2norm_delta_stress / first_l2norm_delta_stress)
282  << "\n"
283  << " stress convergence relative tolerance = " << _relative_tolerance << "\n"
284  << " absolute l2 norm delta stress = " << l2norm_delta_stress << "\n"
285  << " stress convergence absolute tolerance = " << _absolute_tolerance << std::endl;
286  }
287  ++counter;
288  } while (counter < _max_iterations && l2norm_delta_stress > _absolute_tolerance &&
289  (l2norm_delta_stress / first_l2norm_delta_stress) > _relative_tolerance &&
290  _num_models != 1);
291 
292  if (counter == _max_iterations && l2norm_delta_stress > _absolute_tolerance &&
293  (l2norm_delta_stress / first_l2norm_delta_stress) > _relative_tolerance)
294  throw MooseException("Max stress iteration hit during ComputeMultipleInelasticStress solve!");
295 
296  combined_inelastic_strain_increment.zero();
297  for (unsigned i_rmm = 0; i_rmm < _num_models; ++i_rmm)
298  combined_inelastic_strain_increment +=
299  _inelastic_weights[i_rmm] * inelastic_strain_increment[i_rmm];
300 
301  if (_fe_problem.currentlyComputingJacobian())
303 
304  _matl_timestep_limit[_qp] = 0.0;
305  for (unsigned i_rmm = 0; i_rmm < _num_models; ++i_rmm)
306  _matl_timestep_limit[_qp] += 1.0 / _models[i_rmm]->computeTimeStepLimit();
307 
308  _matl_timestep_limit[_qp] = 1.0 / _matl_timestep_limit[_qp];
309 }
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
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)
MaterialProperty< RankTwoTensor > & _stress
const unsigned _num_models
number of plastic models
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Rank-4 and Rank-2 elasticity and elastic strain tensors.
void addQpInitialStress()
InitialStress Deprecation: remove this method.
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
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...
virtual void rotateQpInitialStress()
InitialStress Deprecation: remove this method.
virtual void computeQpJacobianMult()
Using _elasticity_tensor[_qp] and the consistent tangent operators, _comsistent_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...
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 326 of file ComputeMultipleInelasticStress.C.

Referenced by computeQpStressIntermediateConfiguration().

330 {
331  for (auto model : _models)
332  model->setQp(_qp);
333 
334  elastic_strain_increment = _strain_increment[_qp];
335 
336  // If the elasticity tensor values have changed and the tensor is isotropic,
337  // use the old strain to calculate the old stress
339  {
340  _stress[_qp] = _elasticity_tensor[_qp] * (_elastic_strain_old[_qp] + elastic_strain_increment);
341  // InitialStress Deprecation: remove these lines
345  }
346  else
347  _stress[_qp] = _stress_old[_qp] + _elasticity_tensor[_qp] * elastic_strain_increment;
348 
349  computeAdmissibleState(model_number,
350  elastic_strain_increment,
351  combined_inelastic_strain_increment,
352  _Jacobian_mult[_qp]);
353 
354  _matl_timestep_limit[_qp] = _models[0]->computeTimeStepLimit();
355 
356  /* propagate internal variables, etc, to this timestep for those inelastic models where
357  * "updateState" is not called */
358  for (unsigned i_rmm = 0; i_rmm < _num_models; ++i_rmm)
359  if (i_rmm != model_number)
360  _models[i_rmm]->propagateQpStatefulProperties();
361 }
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
MaterialProperty< Real > & _matl_timestep_limit
const MaterialProperty< RankTwoTensor > & _strain_increment
MaterialProperty< RankTwoTensor > & _stress
const unsigned _num_models
number of plastic models
const MaterialProperty< RankFourTensor > & _elasticity_tensor
Rank-4 and Rank-2 elasticity and elastic strain tensors.
void addQpInitialStress()
InitialStress Deprecation: remove this method.
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
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...
virtual void rotateQpInitialStress()
InitialStress Deprecation: remove this method.

Member Data Documentation

const Real ComputeMultipleInelasticStress::_absolute_tolerance
protected

Definition at line 114 of file ComputeMultipleInelasticStress.h.

Referenced by updateQpState().

const std::string ComputeStressBase::_base_name
protectedinherited
std::vector<RankFourTensor> ComputeMultipleInelasticStress::_consistent_tangent_operator
protected

the consistent tangent operators computed by each plastic model

Definition at line 143 of file ComputeMultipleInelasticStress.h.

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

const bool ComputeMultipleInelasticStress::_cycle_models
protected

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

Definition at line 146 of file ComputeMultipleInelasticStress.h.

Referenced by computeQpStressIntermediateConfiguration().

MaterialProperty<RankTwoTensor>& ComputeStressBase::_elastic_strain
protectedinherited
const MaterialProperty<RankTwoTensor>& ComputeMultipleInelasticStress::_elastic_strain_old
protected
const MaterialProperty<RankFourTensor>& ComputeMultipleInelasticStress::_elasticity_tensor
protected
const std::string ComputeStressBase::_elasticity_tensor_name
protectedinherited
const MaterialProperty<RankTwoTensor>& ComputeStressBase::_extra_stress
protectedinherited
MaterialProperty<RankTwoTensor>& ComputeMultipleInelasticStress::_inelastic_strain
protected

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

Definition at line 128 of file ComputeMultipleInelasticStress.h.

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

const MaterialProperty<RankTwoTensor>& ComputeMultipleInelasticStress::_inelastic_strain_old
protected

old value of inelastic strain

Definition at line 131 of file ComputeMultipleInelasticStress.h.

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

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

_inelastic_strain = sum_i (_inelastic_weights_i * inelastic_strain_from_model_i)

Definition at line 140 of file ComputeMultipleInelasticStress.h.

Referenced by ComputeMultipleInelasticStress(), and updateQpState().

MaterialProperty<RankTwoTensor>* ComputeStressBase::_initial_stress
protectedinherited

Initial stress, if provided. InitialStress Deprecation: remove this.

Definition at line 74 of file ComputeStressBase.h.

Referenced by ComputeStressBase::initQpStatefulProperties().

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

initial stress components

Definition at line 62 of file ComputeStressBase.h.

Referenced by ComputeStressBase::ComputeStressBase(), and ComputeStressBase::initQpStatefulProperties().

const MaterialProperty<RankTwoTensor>* ComputeStressBase::_initial_stress_old
protectedinherited

Old value of initial stress, which is needed to correctly implement finite-strain rotations. InitialStress Deprecation: remove this.

Definition at line 77 of file ComputeStressBase.h.

const bool ComputeStressBase::_initial_stress_provided
protectedinherited
bool ComputeMultipleInelasticStress::_is_elasticity_tensor_guaranteed_isotropic
protected

is the elasticity tensor guaranteed to be isotropic?

Definition at line 160 of file ComputeMultipleInelasticStress.h.

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

MaterialProperty<RankFourTensor>& ComputeStressBase::_Jacobian_mult
protectedinherited
MaterialProperty<Real>& ComputeMultipleInelasticStress::_matl_timestep_limit
protected
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 112 of file ComputeMultipleInelasticStress.h.

Referenced by updateQpState().

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_mechanical_strain
protectedinherited
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 157 of file ComputeMultipleInelasticStress.h.

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

const unsigned ComputeMultipleInelasticStress::_num_models
protected
const bool ComputeMultipleInelasticStress::_output_iteration_info
protected

Definition at line 115 of file ComputeMultipleInelasticStress.h.

Referenced by updateQpState().

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 119 of file ComputeMultipleInelasticStress.h.

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

const Real ComputeMultipleInelasticStress::_relative_tolerance
protected

Definition at line 113 of file ComputeMultipleInelasticStress.h.

Referenced by updateQpState().

const MaterialProperty<RankTwoTensor>& ComputeFiniteStrainElasticStress::_rotation_increment
protectedinherited
const bool ComputeStressBase::_store_stress_old
protectedinherited

Parameter which decides whether to store old stress. This is required for HHT time integration and Rayleigh damping.

Definition at line 68 of file ComputeStressBase.h.

const MaterialProperty<RankTwoTensor>& ComputeMultipleInelasticStress::_strain_increment
protected
MaterialProperty<RankTwoTensor>& ComputeStressBase::_stress
protectedinherited

Definition at line 53 of file ComputeStressBase.h.

Referenced by ComputeStressBase::addQpInitialStress(), ComputeMultipleInelasticCosseratStress::computeAdmissibleState(), computeAdmissibleState(), ComputeStressBase::computeQpProperties(), ComputeStressEosBase::computeQpProperties(), ComputeBirchMurnaghanEquationOfStress::computeQpProperties(), ComputeLinearElasticStress::computeQpStress(), ComputeStrainIncrementBasedStress::computeQpStress(), ComputeLinearElasticPFFractureStress::computeQpStress(), ComputeCosseratLinearElasticStress::computeQpStress(), ComputeFiniteStrainElasticStress::computeQpStress(), ComputeIsotropicLinearElasticPFFractureStress::computeQpStress(), ComputeFiniteStrainElasticStressBirchMurnaghan::computeQpStress(), FiniteStrainPlasticMaterial::computeQpStress(), ComputeMultiPlasticityStress::computeQpStress(), ComputeSmearedCrackingStress::computeQpStress(), ComputeLinearViscoelasticStress::computeQpStress(), computeQpStressIntermediateConfiguration(), ComputeSmearedCrackingStress::crackingStressRotation(), finiteStrainRotation(), ComputeStressBase::initQpStatefulProperties(), FiniteStrainCrystalPlasticity::initQpStatefulProperties(), FiniteStrainUObasedCP::initQpStatefulProperties(), FiniteStrainHyperElasticViscoPlastic::initQpStatefulProperties(), ComputeMultiPlasticityStress::postReturnMap(), FiniteStrainUObasedCP::postSolveQp(), FiniteStrainHyperElasticViscoPlastic::postSolveQp(), FiniteStrainCrystalPlasticity::postSolveQp(), updateQpState(), updateQpStateSingleModel(), and LinearIsoElasticPFDamage::updateVar().

const MaterialProperty<RankTwoTensor>& ComputeFiniteStrainElasticStress::_stress_old
protectedinherited
enum ComputeMultipleInelasticStress::TangentOperatorEnum ComputeMultipleInelasticStress::_tangent_operator_type
protected

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