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

ADComputeRSphericalFiniteStrain defines a strain increment and a rotation increment for finite strains in 1D spherical symmetry geometries. More...

#include <ADComputeRSphericalFiniteStrain.h>

Inheritance diagram for ADComputeRSphericalFiniteStrain< compute_stage >:
[legend]

Public Member Functions

 ADComputeRSphericalFiniteStrain (const InputParameters &parameters)
 
virtual void initialSetup ()
 
virtual void computeProperties ()
 Computes the current and old deformation gradients with the assumptions for 1D spherical symmetry geometries: \( \epsilon_{\theta} = \epsilon_{\phi} = \frac{u_r}{r} \). More...
 

Static Public Member Functions

static MooseEnum decompositionType ()
 

Protected Member Functions

virtual void computeQpStrain ()
 
virtual void computeQpIncrements (ADRankTwoTensor &e, ADRankTwoTensor &r)
 
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

const VariableValue & _disp_old_0
 the old value of the first component of the displacements vector More...
 
 usingComputeFiniteStrainMembers
 
std::vector< ADRankTwoTensor > _Fhat
 
 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 ADComputeRSphericalFiniteStrain< compute_stage >

ADComputeRSphericalFiniteStrain defines a strain increment and a rotation increment for finite strains in 1D spherical symmetry geometries.

The strains in the polar and azimuthal directions are functions of the radial displacement.

Definition at line 15 of file ADComputeRSphericalFiniteStrain.h.

Constructor & Destructor Documentation

◆ ADComputeRSphericalFiniteStrain()

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

Definition at line 26 of file ADComputeRSphericalFiniteStrain.C.

29  _disp_old_0(coupledValueOld("displacements", 0))
30 {
31 }
const VariableValue & _disp_old_0
the old value of the first component of the displacements vector
ADComputeFiniteStrain defines a strain increment and rotation increment, for finite strains...

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 ADComputeRSphericalFiniteStrain< compute_stage >::computeProperties ( )
virtual

Computes the current and old deformation gradients with the assumptions for 1D spherical symmetry geometries: \( \epsilon_{\theta} = \epsilon_{\phi} = \frac{u_r}{r} \).

Reimplemented from ADComputeFiniteStrain< compute_stage >.

Definition at line 47 of file ADComputeRSphericalFiniteStrain.C.

48 {
49  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
50  {
51  // Deformation gradient calculation in cylindrical coordinates
52  ADRankTwoTensor A; // Deformation gradient
53  RankTwoTensor Fbar; // Old Deformation gradient
54 
55  // Step through calculating the current and old deformation gradients
56  // Only diagonal components are nonzero because this is a 1D material
57  // Note: x_disp is the radial displacement
58  A(0, 0) = (*_grad_disp[0])[_qp](0);
59  Fbar(0, 0) = (*_grad_disp_old[0])[_qp](0);
60 
61  // The polar and azimuthal strains are functions of radial displacement
62  if (!MooseUtils::relativeFuzzyEqual(_q_point[_qp](0), 0.0))
63  {
64  A(1, 1) = (*_disp[0])[_qp] / _q_point[_qp](0);
65  Fbar(1, 1) = _disp_old_0[_qp] / _q_point[_qp](0);
66  }
67 
68  // The polar and azimuthal strains are equivalent in this 1D problem
69  A(2, 2) = A(1, 1);
70  Fbar(2, 2) = Fbar(1, 1);
71 
72  // very nearly A = gradU - gradUold, adapted to cylindrical coords
73  A -= Fbar;
74 
75  // Fbar = ( I + gradUold)
76  Fbar.addIa(1.0);
77 
78  // Incremental deformation gradient _Fhat = I + A Fbar^-1
79  _Fhat[_qp] = A * Fbar.inverse();
80  _Fhat[_qp].addIa(1.0);
81 
83  }
84 }
std::vector< const ADVariableValue * > _disp
std::vector< const ADVariableGradient * > _grad_disp
const VariableValue & _disp_old_0
the old value of the first component of the displacements vector
std::vector< const VariableGradient * > _grad_disp_old
std::vector< ADRankTwoTensor > _Fhat

◆ computeQpIncrements()

template<ComputeStage compute_stage>
void ADComputeFiniteStrain< compute_stage >::computeQpIncrements ( ADRankTwoTensor &  e,
ADRankTwoTensor &  r 
)
protectedvirtualinherited

Definition at line 122 of file ADComputeFiniteStrain.C.

124 {
125  switch (_decomposition_method)
126  {
128  {
129  // inverse of _Fhat
130  ADRankTwoTensor invFhat;
131  static const RankTwoTensor zero;
132  if (_Fhat[_qp] == zero)
133  invFhat.zero();
134  else
135  invFhat = _Fhat[_qp].inverse();
136 
137  // A = I - _Fhat^-1
138  ADRankTwoTensor A(RankTwoTensorType<compute_stage>::type::initIdentity);
139  A -= invFhat;
140 
141  // Cinv - I = A A^T - A - A^T;
142  ADRankTwoTensor Cinv_I = A * A.transpose() - A - A.transpose();
143 
144  // strain rate D from Taylor expansion, Chat = (-1/2(Chat^-1 - I) + 1/4*(Chat^-1 - I)^2 + ...
145  total_strain_increment = -Cinv_I * 0.5 + Cinv_I * Cinv_I * 0.25;
146 
147  const ADReal a[3] = {invFhat(1, 2) - invFhat(2, 1),
148  invFhat(2, 0) - invFhat(0, 2),
149  invFhat(0, 1) - invFhat(1, 0)};
150 
151  const auto q = (a[0] * a[0] + a[1] * a[1] + a[2] * a[2]) / 4.0;
152  const auto trFhatinv_1 = invFhat.trace() - 1.0;
153  const auto p = trFhatinv_1 * trFhatinv_1 / 4.0;
154 
155  // cos theta_a
156  const auto C1 =
157  std::sqrt(p + 3.0 * Utility::pow<2>(p) * (1.0 - (p + q)) / Utility::pow<2>(p + q) -
158  2.0 * Utility::pow<3>(p) * (1.0 - (p + q)) / Utility::pow<3>(p + q));
159 
160  ADReal C2;
161  if (q > 0.01)
162  // (1-cos theta_a)/4q
163  C2 = (1.0 - C1) / (4.0 * q);
164  else
165  // alternate form for small q
166  C2 = 0.125 + q * 0.03125 * (Utility::pow<2>(p) - 12.0 * (p - 1.0)) / Utility::pow<2>(p) +
167  Utility::pow<2>(q) * (p - 2.0) * (Utility::pow<2>(p) - 10.0 * p + 32.0) /
168  Utility::pow<3>(p) +
169  Utility::pow<3>(q) *
170  (1104.0 - 992.0 * p + 376.0 * Utility::pow<2>(p) - 72.0 * Utility::pow<3>(p) +
171  5.0 * Utility::pow<4>(p)) /
172  (512.0 * Utility::pow<4>(p));
173 
174  const auto C3 =
175  0.5 * std::sqrt((p * q * (3.0 - q) + Utility::pow<3>(p) + Utility::pow<2>(q)) /
176  Utility::pow<3>(p + q)); // sin theta_a/(2 sqrt(q))
177 
178  // Calculate incremental rotation. Note that this value is the transpose of that from Rashid,
179  // 93, so we transpose it before storing
180  ADRankTwoTensor R_incr;
181  R_incr.addIa(C1);
182  for (unsigned int i = 0; i < 3; ++i)
183  for (unsigned int j = 0; j < 3; ++j)
184  R_incr(i, j) += C2 * a[i] * a[j];
185 
186  R_incr(0, 1) += C3 * a[2];
187  R_incr(0, 2) -= C3 * a[1];
188  R_incr(1, 0) -= C3 * a[2];
189  R_incr(1, 2) += C3 * a[0];
190  R_incr(2, 0) += C3 * a[1];
191  R_incr(2, 1) -= C3 * a[0];
192 
193  rotation_increment = R_incr.transpose();
194  break;
195  }
196 
198  {
199  std::vector<ADReal> e_value(3);
200  ADRankTwoTensor e_vector, N1, N2, N3;
201 
202  const auto Chat = _Fhat[_qp].transpose() * _Fhat[_qp];
203  Chat.symmetricEigenvaluesEigenvectors(e_value, e_vector);
204 
205  const auto lambda1 = std::sqrt(e_value[0]);
206  const auto lambda2 = std::sqrt(e_value[1]);
207  const auto lambda3 = std::sqrt(e_value[2]);
208 
209  N1.vectorOuterProduct(e_vector.column(0), e_vector.column(0));
210  N2.vectorOuterProduct(e_vector.column(1), e_vector.column(1));
211  N3.vectorOuterProduct(e_vector.column(2), e_vector.column(2));
212 
213  const auto Uhat = N1 * lambda1 + N2 * lambda2 + N3 * lambda3;
214  ADRankTwoTensor invUhat(Uhat.inverse());
215 
216  rotation_increment = _Fhat[_qp] * invUhat;
217 
218  total_strain_increment =
219  N1 * std::log(lambda1) + N2 * std::log(lambda2) + N3 * std::log(lambda3);
220  break;
221  }
222 
223  default:
224  mooseError("ADComputeFiniteStrain Error: Pass valid decomposition type: TaylorExpansion or "
225  "EigenSolution.");
226  }
227 }
const DecompMethod _decomposition_method
std::vector< ADRankTwoTensor > _Fhat

◆ computeQpStrain()

template<ComputeStage compute_stage>
void ADComputeFiniteStrain< compute_stage >::computeQpStrain ( )
protectedvirtualinherited

Definition at line 89 of file ADComputeFiniteStrain.C.

90 {
91  ADRankTwoTensor total_strain_increment;
92 
93  // two ways to calculate these increments: TaylorExpansion(default) or EigenSolution
94  computeQpIncrements(total_strain_increment, _rotation_increment[_qp]);
95 
96  _strain_increment[_qp] = total_strain_increment;
97 
98  // Remove the eigenstrain increment
99  subtractEigenstrainIncrementFromStrain(_strain_increment[_qp]);
100 
101  if (_dt > 0)
102  _strain_rate[_qp] = _strain_increment[_qp] / _dt;
103  else
104  _strain_rate[_qp].zero();
105 
106  // Update strain in intermediate configuration
107  _mechanical_strain[_qp] = _mechanical_strain_old[_qp] + _strain_increment[_qp];
108  _total_strain[_qp] = _total_strain_old[_qp] + total_strain_increment;
109 
110  // Rotate strain to current configuration
111  _mechanical_strain[_qp] =
112  _rotation_increment[_qp] * _mechanical_strain[_qp] * _rotation_increment[_qp].transpose();
113  _total_strain[_qp] =
114  _rotation_increment[_qp] * _total_strain[_qp] * _rotation_increment[_qp].transpose();
115 
116  if (_global_strain)
117  _total_strain[_qp] += (*_global_strain)[_qp];
118 }
void subtractEigenstrainIncrementFromStrain(ADRankTwoTensor &strain)
const MaterialProperty< RankTwoTensor > & _mechanical_strain_old
const MaterialProperty< RankTwoTensor > & _total_strain_old
virtual void computeQpIncrements(ADRankTwoTensor &e, ADRankTwoTensor &r)

◆ decompositionType()

template<ComputeStage compute_stage>
MooseEnum ADComputeFiniteStrain< compute_stage >::decompositionType ( )
staticinherited

Definition at line 17 of file ADComputeFiniteStrain.C.

18 {
19  return MooseEnum("TaylorExpansion EigenSolution", "TaylorExpansion");
20 }

◆ 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 ADComputeRSphericalFiniteStrain< compute_stage >::initialSetup ( )
virtual

Definition at line 35 of file ADComputeRSphericalFiniteStrain.C.

36 {
38 
39  const auto & subdomainIDs = _mesh.meshSubdomains();
40  for (auto subdomainID : subdomainIDs)
41  if (_fe_problem.getCoordSystem(subdomainID) != Moose::COORD_RSPHERICAL)
42  mooseError("The coordinate system must be set to RSPHERICAL for 1D R spherical simulations.");
43 }

◆ 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 69 of file ADComputeStrainBase.h.

◆ _disp

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

Definition at line 55 of file ADComputeStrainBase.h.

◆ _disp_old_0

template<ComputeStage compute_stage>
const VariableValue& ADComputeRSphericalFiniteStrain< compute_stage >::_disp_old_0
protected

the old value of the first component of the displacements vector

Definition at line 39 of file ADComputeRSphericalFiniteStrain.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

◆ _Fhat

template<ComputeStage compute_stage>
std::vector<ADRankTwoTensor> ADComputeFiniteStrain< compute_stage >::_Fhat
protectedinherited

Definition at line 42 of file ADComputeFiniteStrain.h.

◆ _grad_disp

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

Definition at line 56 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 47 of file ADComputeIncrementalStrainBase.h.

◆ _mechanical_strain_old

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

Definition at line 53 of file ADComputeIncrementalStrainBase.h.

◆ _ndisp

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

Coupled displacement variables.

Definition at line 54 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 54 of file ADComputeIncrementalStrainBase.h.

◆ _volumetric_locking_correction

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

◆ usingComputeFiniteStrainMembers

template<ComputeStage compute_stage>
ADComputeRSphericalFiniteStrain< compute_stage >::usingComputeFiniteStrainMembers
protected

Definition at line 41 of file ADComputeRSphericalFiniteStrain.h.

◆ usingComputeIncrementalStrainBaseMembers

template<ComputeStage compute_stage>
ADComputeFiniteStrain< compute_stage >::usingComputeIncrementalStrainBaseMembers
protectedinherited

Definition at line 54 of file ADComputeFiniteStrain.h.

◆ usingComputeStrainBaseMembers

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

Definition at line 58 of file ADComputeIncrementalStrainBase.h.

◆ usingMaterialMembers

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

Definition at line 71 of file ADComputeStrainBase.h.


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