www.mooseframework.org
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
ADComputeIncrementalSmallStrain< compute_stage > Class Template Reference

ADComputeIncrementalSmallStrain defines a strain increment and rotation increment (=1), for small strains. More...

#include <ADComputeIncrementalSmallStrain.h>

Inheritance diagram for ADComputeIncrementalSmallStrain< compute_stage >:
[legend]

Public Member Functions

 ADComputeIncrementalSmallStrain (const InputParameters &parameters)
 
virtual void computeProperties () override
 
void initialSetup () override
 

Protected Member Functions

virtual void computeTotalStrainIncrement (ADRankTwoTensor &total_strain_increment)
 Computes the current and old deformation gradients and passes back the total strain increment tensor. More...
 
virtual void initQpStatefulProperties () override
 
void subtractEigenstrainIncrementFromStrain (ADRankTwoTensor &strain)
 
 ADMaterialProperty (RankTwoTensor) &_strain_rate
 
 ADMaterialProperty (RankTwoTensor) &_strain_increment
 
 ADMaterialProperty (RankTwoTensor) &_rotation_increment
 
const ADMaterialProperty (RankTwoTensor) *_global_strain
 
virtual void displacementIntegrityCheck ()
 

Protected Attributes

 usingComputeIncrementalStrainBaseMembers
 
std::vector< const VariableGradient * > _grad_disp_old
 
const MaterialProperty< RankTwoTensor > & _mechanical_strain_old
 
const MaterialProperty< RankTwoTensor > & _total_strain_old
 
std::vector< const MaterialProperty< RankTwoTensor > * > _eigenstrains_old
 
 usingComputeStrainBaseMembers
 
const unsigned int _ndisp
 Coupled displacement variables. More...
 
std::vector< const ADVariableValue * > _disp
 
std::vector< const ADVariableGradient * > _grad_disp
 
std::string _base_name
 
std::vector< MaterialPropertyName > _eigenstrain_names
 
std::vector< const ADMaterialProperty(RankTwoTensor) * > _eigenstrains
 
const bool _volumetric_locking_correction
 
const Real & _current_elem_volume
 
 usingMaterialMembers
 

Detailed Description

template<ComputeStage compute_stage>
class ADComputeIncrementalSmallStrain< compute_stage >

ADComputeIncrementalSmallStrain defines a strain increment and rotation increment (=1), for small strains.

Definition at line 20 of file ADComputeIncrementalSmallStrain.h.

Constructor & Destructor Documentation

◆ ADComputeIncrementalSmallStrain()

template<ComputeStage compute_stage>
ADComputeIncrementalSmallStrain< compute_stage >::ADComputeIncrementalSmallStrain ( const InputParameters &  parameters)

Definition at line 21 of file ADComputeIncrementalSmallStrain.C.

24 {
25 }
ADComputeIncrementalStrainBase is the base class for strain tensors using incremental formulations...

Member Function Documentation

◆ ADMaterialProperty() [1/4]

template<ComputeStage compute_stage>
ADComputeIncrementalStrainBase< compute_stage >::ADMaterialProperty ( RankTwoTensor  ) &
protectedinherited

◆ ADMaterialProperty() [2/4]

template<ComputeStage compute_stage>
ADComputeIncrementalStrainBase< compute_stage >::ADMaterialProperty ( RankTwoTensor  ) &
protectedinherited

◆ ADMaterialProperty() [3/4]

template<ComputeStage compute_stage>
ADComputeIncrementalStrainBase< compute_stage >::ADMaterialProperty ( RankTwoTensor  ) &
protectedinherited

◆ ADMaterialProperty() [4/4]

template<ComputeStage compute_stage>
const ADComputeStrainBase< compute_stage >::ADMaterialProperty ( RankTwoTensor  )
protectedinherited

◆ computeProperties()

template<ComputeStage compute_stage>
void ADComputeIncrementalSmallStrain< compute_stage >::computeProperties ( )
overridevirtual

Definition at line 29 of file ADComputeIncrementalSmallStrain.C.

30 {
31  ADReal volumetric_strain = 0.0;
32  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
33  {
34  ADRankTwoTensor total_strain_increment;
35  computeTotalStrainIncrement(total_strain_increment);
36 
37  _strain_increment[_qp] = total_strain_increment;
38 
40  volumetric_strain += total_strain_increment.trace() * _JxW[_qp] * _coord[_qp];
41  }
42 
44  volumetric_strain /= _current_elem_volume;
45 
46  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
47  {
49  {
50  const auto correction = (volumetric_strain - _strain_increment[_qp].trace()) / 3.0;
51  _strain_increment[_qp](0, 0) += correction;
52  _strain_increment[_qp](1, 1) += correction;
53  _strain_increment[_qp](2, 2) += correction;
54  }
55 
56  _total_strain[_qp] = _strain_increment[_qp] + _total_strain_old[_qp];
57 
58  // Remove the Eigen strain increment
59  subtractEigenstrainIncrementFromStrain(_strain_increment[_qp]);
60 
61  // strain rate
62  if (_dt > 0)
63  _strain_rate[_qp] = _strain_increment[_qp] / _dt;
64  else
65  _strain_rate[_qp].zero();
66 
67  // Update strain in intermediate configuration: rotations are not needed
68  _mechanical_strain[_qp] = _strain_increment[_qp] + _mechanical_strain_old[_qp];
69 
70  // incremental small strain does not include rotation
71  _rotation_increment[_qp].setToIdentity();
72  }
73 
74  copyDualNumbersToValues();
75 }
const bool _volumetric_locking_correction
virtual void computeTotalStrainIncrement(ADRankTwoTensor &total_strain_increment)
Computes the current and old deformation gradients and passes back the total strain increment tensor...
void subtractEigenstrainIncrementFromStrain(ADRankTwoTensor &strain)
const Real & _current_elem_volume
const MaterialProperty< RankTwoTensor > & _mechanical_strain_old
const MaterialProperty< RankTwoTensor > & _total_strain_old

◆ computeTotalStrainIncrement()

template<ComputeStage compute_stage>
void ADComputeIncrementalSmallStrain< compute_stage >::computeTotalStrainIncrement ( ADRankTwoTensor &  total_strain_increment)
protectedvirtual

Computes the current and old deformation gradients and passes back the total strain increment tensor.

Reimplemented in ADCompute2DIncrementalStrain< compute_stage >, ADCompute1DIncrementalStrain< compute_stage >, and ADComputeRSphericalIncrementalStrain< compute_stage >.

Definition at line 79 of file ADComputeIncrementalSmallStrain.C.

81 {
82  // Deformation gradient
83  ADRankTwoTensor A(
84  (*_grad_disp[0])[_qp], (*_grad_disp[1])[_qp], (*_grad_disp[2])[_qp]); // Deformation gradient
85  RankTwoTensor Fbar((*_grad_disp_old[0])[_qp],
86  (*_grad_disp_old[1])[_qp],
87  (*_grad_disp_old[2])[_qp]); // Old Deformation gradient
88 
89  A -= Fbar; // A = grad_disp - grad_disp_old
90 
91  total_strain_increment = 0.5 * (A + A.transpose());
92 }
std::vector< const ADVariableGradient * > _grad_disp
std::vector< const VariableGradient * > _grad_disp_old

◆ displacementIntegrityCheck()

template<ComputeStage compute_stage>
void ADComputeStrainBase< compute_stage >::displacementIntegrityCheck ( )
protectedvirtualinherited

Reimplemented in ADCompute2DFiniteStrain< compute_stage >, ADCompute2DIncrementalStrain< compute_stage >, and ADCompute2DSmallStrain< compute_stage >.

Definition at line 87 of file ADComputeStrainBase.C.

88 {
89  // Checking for consistency between mesh size and length of the provided displacements vector
90  if (_ndisp != _mesh.dimension())
91  paramError(
92  "displacements",
93  "The number of variables supplied in 'displacements' must match the mesh dimension.");
94 }
const unsigned int _ndisp
Coupled displacement variables.

◆ initialSetup()

template<ComputeStage compute_stage>
void ADComputeIncrementalStrainBase< compute_stage >::initialSetup ( )
overrideinherited

Definition at line 34 of file ADComputeIncrementalStrainBase.C.

Referenced by ADComputeAxisymmetricRZIncrementalStrain< compute_stage >::initialSetup(), ADComputeRSphericalFiniteStrain< compute_stage >::initialSetup(), and ADComputeRSphericalIncrementalStrain< compute_stage >::initialSetup().

35 {
37  for (unsigned int i = 0; i < 3; ++i)
38  {
39  if (_fe_problem.isTransient() && i < _ndisp)
40  _grad_disp_old[i] = &coupledGradientOld("displacements", i);
41  else
42  _grad_disp_old[i] = &_grad_zero;
43  }
44 }
const unsigned int _ndisp
Coupled displacement variables.
void initialSetup() override
std::vector< const VariableGradient * > _grad_disp_old

◆ initQpStatefulProperties()

template<ComputeStage compute_stage>
void ADComputeIncrementalStrainBase< compute_stage >::initQpStatefulProperties ( )
overrideprotectedvirtualinherited

Reimplemented from ADComputeStrainBase< compute_stage >.

Definition at line 48 of file ADComputeIncrementalStrainBase.C.

49 {
50  _mechanical_strain[_qp].zero();
51  _total_strain[_qp].zero();
52 }

◆ subtractEigenstrainIncrementFromStrain()

template<ComputeStage compute_stage>
void ADComputeIncrementalStrainBase< compute_stage >::subtractEigenstrainIncrementFromStrain ( ADRankTwoTensor &  strain)
protectedinherited

Definition at line 56 of file ADComputeIncrementalStrainBase.C.

58 {
59  for (unsigned int i = 0; i < _eigenstrains.size(); ++i)
60  {
61  strain -= (*_eigenstrains[i])[_qp];
62  strain += (*_eigenstrains_old[i])[_qp];
63  }
64 }
std::vector< const ADMaterialProperty(RankTwoTensor) * > _eigenstrains
std::vector< const MaterialProperty< RankTwoTensor > * > _eigenstrains_old

Member Data Documentation

◆ _base_name

template<ComputeStage compute_stage>
std::string ADComputeStrainBase< compute_stage >::_base_name
protectedinherited

◆ _current_elem_volume

template<ComputeStage compute_stage>
const Real& ADComputeStrainBase< compute_stage >::_current_elem_volume
protectedinherited

Definition at line 70 of file ADComputeStrainBase.h.

◆ _disp

template<ComputeStage compute_stage>
std::vector<const ADVariableValue *> ADComputeStrainBase< compute_stage >::_disp
protectedinherited

Definition at line 56 of file ADComputeStrainBase.h.

◆ _eigenstrain_names

template<ComputeStage compute_stage>
std::vector<MaterialPropertyName> ADComputeStrainBase< compute_stage >::_eigenstrain_names
protectedinherited

◆ _eigenstrains

template<ComputeStage compute_stage>
std::vector<const ADMaterialProperty(RankTwoTensor) *> ADComputeStrainBase< compute_stage >::_eigenstrains
protectedinherited

◆ _eigenstrains_old

template<ComputeStage compute_stage>
std::vector<const MaterialProperty<RankTwoTensor> *> ADComputeIncrementalStrainBase< compute_stage >::_eigenstrains_old
protectedinherited

◆ _grad_disp

template<ComputeStage compute_stage>
std::vector<const ADVariableGradient *> ADComputeStrainBase< compute_stage >::_grad_disp
protectedinherited

Definition at line 57 of file ADComputeStrainBase.h.

◆ _grad_disp_old

template<ComputeStage compute_stage>
std::vector<const VariableGradient *> ADComputeIncrementalStrainBase< compute_stage >::_grad_disp_old
protectedinherited

Definition at line 48 of file ADComputeIncrementalStrainBase.h.

◆ _mechanical_strain_old

template<ComputeStage compute_stage>
const MaterialProperty<RankTwoTensor>& ADComputeIncrementalStrainBase< compute_stage >::_mechanical_strain_old
protectedinherited

Definition at line 54 of file ADComputeIncrementalStrainBase.h.

◆ _ndisp

template<ComputeStage compute_stage>
const unsigned int ADComputeStrainBase< compute_stage >::_ndisp
protectedinherited

Coupled displacement variables.

Definition at line 55 of file ADComputeStrainBase.h.

Referenced by ADComputeStrainBase< compute_stage >::ADComputeStrainBase().

◆ _total_strain_old

template<ComputeStage compute_stage>
const MaterialProperty<RankTwoTensor>& ADComputeIncrementalStrainBase< compute_stage >::_total_strain_old
protectedinherited

Definition at line 55 of file ADComputeIncrementalStrainBase.h.

◆ _volumetric_locking_correction

template<ComputeStage compute_stage>
const bool ADComputeStrainBase< compute_stage >::_volumetric_locking_correction
protectedinherited

◆ usingComputeIncrementalStrainBaseMembers

template<ComputeStage compute_stage>
ADComputeIncrementalSmallStrain< compute_stage >::usingComputeIncrementalStrainBaseMembers
protected

Definition at line 43 of file ADComputeIncrementalSmallStrain.h.

◆ usingComputeStrainBaseMembers

template<ComputeStage compute_stage>
ADComputeIncrementalStrainBase< compute_stage >::usingComputeStrainBaseMembers
protectedinherited

Definition at line 59 of file ADComputeIncrementalStrainBase.h.

◆ usingMaterialMembers

template<ComputeStage compute_stage>
ADComputeStrainBase< compute_stage >::usingMaterialMembers
protectedinherited

Definition at line 72 of file ADComputeStrainBase.h.


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