www.mooseframework.org
ComputeMultipleInelasticStress.h
Go to the documentation of this file.
1 /****************************************************************/
2 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
3 /* */
4 /* All contents are licensed under LGPL V2.1 */
5 /* See LICENSE for full restrictions */
6 /****************************************************************/
7 #ifndef COMPUTEMULTIPLEINELASTICSTRESS_H
8 #define COMPUTEMULTIPLEINELASTICSTRESS_H
9 
11 
12 class StressUpdateBase;
13 
31 {
32 public:
33  ComputeMultipleInelasticStress(const InputParameters & parameters);
34 
35 protected:
36  virtual void initQpStatefulProperties() override;
37  virtual void initialSetup() override;
38 
39  virtual void computeQpStress() override;
40 
46 
53  virtual void finiteStrainRotation(const bool force_elasticity_rotation = false);
54 
67  virtual void updateQpState(RankTwoTensor & elastic_strain_increment,
68  RankTwoTensor & combined_inelastic_strain_increment);
69 
80  virtual void updateQpStateSingleModel(unsigned model_number,
81  RankTwoTensor & elastic_strain_increment,
82  RankTwoTensor & combined_inelastic_strain_increment);
83 
89  virtual void computeQpJacobianMult();
90 
106  virtual void computeAdmissibleState(unsigned model_number,
107  RankTwoTensor & elastic_strain_increment,
108  RankTwoTensor & inelastic_strain_increment,
109  RankFourTensor & consistent_tangent_operator);
110 
112  const unsigned int _max_iterations;
117 
120 
122  const MaterialProperty<RankFourTensor> & _elasticity_tensor;
123  const MaterialProperty<RankTwoTensor> & _elastic_strain_old;
124  const MaterialProperty<RankTwoTensor> & _strain_increment;
126 
128  MaterialProperty<RankTwoTensor> & _inelastic_strain;
129 
131  const MaterialProperty<RankTwoTensor> & _inelastic_strain_old;
132 
135 
137  const unsigned _num_models;
138 
140  const std::vector<Real> _inelastic_weights;
141 
143  std::vector<RankFourTensor> _consistent_tangent_operator;
144 
146  const bool _cycle_models;
147 
148  MaterialProperty<Real> & _matl_timestep_limit;
149 
157  std::vector<StressUpdateBase *> _models;
158 
161 };
162 
163 #endif // COMPUTEMULTIPLEINELASTICSTRESS_H
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
ComputeMultipleInelasticStress computes the stress, the consistent tangent operator (or an approximat...
const MaterialProperty< RankTwoTensor > & _elastic_strain_old
MaterialProperty< Real > & _matl_timestep_limit
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...
StressUpdateBase is a material that is not called by MOOSE because of the compute=false flag set in t...
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
virtual void computeQpStressIntermediateConfiguration()
Compute the stress for the current QP, but do not rotate tensors from the intermediate configuration ...
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...
TangentOperatorEnum
what sort of Tangent operator to calculate
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.
ComputeFiniteStrainElasticStress computes the stress following elasticity theory for finite strains...
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...
enum ComputeMultipleInelasticStress::TangentOperatorEnum _tangent_operator_type
virtual void finiteStrainRotation(const bool force_elasticity_rotation=false)
Rotate _elastic_strain, _stress, _inelastic_strain, and _Jacobian_mult to the new configuration...
std::vector< StressUpdateBase * > _models
The user supplied list of inelastic models to use in the simulation.
const bool _perform_finite_strain_rotations
after updateQpState, rotate the stress, elastic_strain, inelastic_strain and Jacobian_mult using _rot...
virtual void computeQpJacobianMult()
Using _elasticity_tensor[_qp] and the consistent tangent operators, _comsistent_tangent_operator[...] computed by the inelastic models, compute _Jacobian_mult[_qp].
ComputeMultipleInelasticStress(const InputParameters &parameters)
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...