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

ComputeFiniteStrain defines a strain increment and rotation increment, for finite strains. More...

#include <ComputeFiniteStrain.h>

Inheritance diagram for ComputeFiniteStrain:
[legend]

Public Member Functions

 ComputeFiniteStrain (const InputParameters &parameters)
 
virtual void computeProperties ()
 

Static Public Member Functions

static MooseEnum decompositionType ()
 

Protected Member Functions

virtual void computeQpStrain ()
 
virtual void computeQpIncrements (RankTwoTensor &e, RankTwoTensor &r)
 
void initialSetup () override
 
virtual void initQpStatefulProperties () override
 
void subtractEigenstrainIncrementFromStrain (RankTwoTensor &strain)
 
virtual void displacementIntegrityCheck ()
 

Protected Attributes

std::vector< RankTwoTensor > _Fhat
 
std::vector< const VariableGradient * > _grad_disp_old
 
MaterialProperty< RankTwoTensor > & _strain_rate
 
MaterialProperty< RankTwoTensor > & _strain_increment
 
MaterialProperty< RankTwoTensor > & _rotation_increment
 
MaterialProperty< RankTwoTensor > & _deformation_gradient
 
const MaterialProperty< RankTwoTensor > & _mechanical_strain_old
 
const MaterialProperty< RankTwoTensor > & _total_strain_old
 
std::vector< const MaterialProperty< RankTwoTensor > * > _eigenstrains_old
 
unsigned int _ndisp
 Coupled displacement variables. More...
 
std::vector< const VariableValue * > _disp
 
std::vector< const VariableGradient * > _grad_disp
 
std::string _base_name
 
MaterialProperty< RankTwoTensor > & _mechanical_strain
 
MaterialProperty< RankTwoTensor > & _total_strain
 
std::vector< MaterialPropertyName > _eigenstrain_names
 
std::vector< const MaterialProperty< RankTwoTensor > * > _eigenstrains
 
const MaterialProperty< RankTwoTensor > * _global_strain
 
bool _volumetric_locking_correction
 
const Real & _current_elem_volume
 

Private Types

enum  DecompMethod { DecompMethod::TaylorExpansion, DecompMethod::EigenSolution }
 

Private Attributes

const DecompMethod _decomposition_method
 

Detailed Description

ComputeFiniteStrain defines a strain increment and rotation increment, for finite strains.

Definition at line 23 of file ComputeFiniteStrain.h.

Member Enumeration Documentation

◆ DecompMethod

Enumerator
TaylorExpansion 
EigenSolution 

Definition at line 39 of file ComputeFiniteStrain.h.

40  {
41  TaylorExpansion,
42  EigenSolution
43  };

Constructor & Destructor Documentation

◆ ComputeFiniteStrain()

ComputeFiniteStrain::ComputeFiniteStrain ( const InputParameters &  parameters)

Definition at line 37 of file ComputeFiniteStrain.C.

38  : ComputeIncrementalStrainBase(parameters),
39  _Fhat(_fe_problem.getMaxQps()),
40  _decomposition_method(getParam<MooseEnum>("decomposition_method").getEnum<DecompMethod>())
41 {
42 }
const DecompMethod _decomposition_method
ComputeIncrementalStrainBase(const InputParameters &parameters)
std::vector< RankTwoTensor > _Fhat

Member Function Documentation

◆ computeProperties()

void ComputeFiniteStrain::computeProperties ( )
virtual

Reimplemented in ComputeRSphericalFiniteStrain, Compute2DFiniteStrain, and Compute1DFiniteStrain.

Definition at line 45 of file ComputeFiniteStrain.C.

46 {
47  RankTwoTensor ave_Fhat;
48  Real ave_dfgrd_det = 0.0;
49  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
50  {
51  // Deformation gradient
52  RankTwoTensor A((*_grad_disp[0])[_qp],
53  (*_grad_disp[1])[_qp],
54  (*_grad_disp[2])[_qp]); // Deformation gradient
55  RankTwoTensor Fbar((*_grad_disp_old[0])[_qp],
56  (*_grad_disp_old[1])[_qp],
57  (*_grad_disp_old[2])[_qp]); // Old Deformation gradient
58 
59  _deformation_gradient[_qp] = A;
60  _deformation_gradient[_qp].addIa(1.0); // Gauss point deformation gradient
61 
62  // A = gradU - gradUold
63  A -= Fbar;
64 
65  // Fbar = ( I + gradUold)
66  Fbar.addIa(1.0);
67 
68  // Incremental deformation gradient _Fhat = I + A Fbar^-1
69  _Fhat[_qp] = A * Fbar.inverse();
70  _Fhat[_qp].addIa(1.0);
71 
73  {
74  // Calculate average _Fhat for volumetric locking correction
75  ave_Fhat += _Fhat[_qp] * _JxW[_qp] * _coord[_qp];
76 
77  // Average deformation gradient
78  ave_dfgrd_det += _deformation_gradient[_qp].det() * _JxW[_qp] * _coord[_qp];
79  }
80  }
81 
83  {
84  // needed for volumetric locking correction
85  ave_Fhat /= _current_elem_volume;
86  // average deformation gradient
87  ave_dfgrd_det /= _current_elem_volume;
88  }
89 
90  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
91  {
93  {
94  // Finalize volumetric locking correction
95  _Fhat[_qp] *= std::cbrt(ave_Fhat.det() / _Fhat[_qp].det());
96  _deformation_gradient[_qp] *= std::cbrt(ave_dfgrd_det / _deformation_gradient[_qp].det());
97  }
98 
100  }
101 }
const Real & _current_elem_volume
MaterialProperty< RankTwoTensor > & _deformation_gradient
virtual void computeQpStrain()
std::vector< const VariableGradient * > _grad_disp_old
std::vector< RankTwoTensor > _Fhat
std::vector< const VariableGradient * > _grad_disp

◆ computeQpIncrements()

void ComputeFiniteStrain::computeQpIncrements ( RankTwoTensor &  e,
RankTwoTensor &  r 
)
protectedvirtual

Definition at line 136 of file ComputeFiniteStrain.C.

Referenced by computeQpStrain().

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

◆ computeQpStrain()

void ComputeFiniteStrain::computeQpStrain ( )
protectedvirtual

Definition at line 104 of file ComputeFiniteStrain.C.

Referenced by computeProperties(), Compute1DFiniteStrain::computeProperties(), Compute2DFiniteStrain::computeProperties(), and ComputeRSphericalFiniteStrain::computeProperties().

105 {
106  RankTwoTensor total_strain_increment;
107 
108  // two ways to calculate these increments: TaylorExpansion(default) or EigenSolution
109  computeQpIncrements(total_strain_increment, _rotation_increment[_qp]);
110 
111  _strain_increment[_qp] = total_strain_increment;
112 
113  // Remove the eigenstrain increment
115 
116  if (_dt > 0)
117  _strain_rate[_qp] = _strain_increment[_qp] / _dt;
118  else
119  _strain_rate[_qp].zero();
120 
121  // Update strain in intermediate configuration
123  _total_strain[_qp] = _total_strain_old[_qp] + total_strain_increment;
124 
125  // Rotate strain to current configuration
126  _mechanical_strain[_qp] =
127  _rotation_increment[_qp] * _mechanical_strain[_qp] * _rotation_increment[_qp].transpose();
128  _total_strain[_qp] =
129  _rotation_increment[_qp] * _total_strain[_qp] * _rotation_increment[_qp].transpose();
130 
131  if (_global_strain)
132  _total_strain[_qp] += (*_global_strain)[_qp];
133 }
const MaterialProperty< RankTwoTensor > & _total_strain_old
MaterialProperty< RankTwoTensor > & _strain_increment
MaterialProperty< RankTwoTensor > & _mechanical_strain
virtual void computeQpIncrements(RankTwoTensor &e, RankTwoTensor &r)
MaterialProperty< RankTwoTensor > & _strain_rate
void subtractEigenstrainIncrementFromStrain(RankTwoTensor &strain)
const MaterialProperty< RankTwoTensor > & _mechanical_strain_old
const MaterialProperty< RankTwoTensor > * _global_strain
MaterialProperty< RankTwoTensor > & _rotation_increment
MaterialProperty< RankTwoTensor > & _total_strain

◆ decompositionType()

MooseEnum ComputeFiniteStrain::decompositionType ( )
static

Definition at line 17 of file ComputeFiniteStrain.C.

Referenced by validParams< ComputeFiniteStrain >(), and validParams< TensorMechanicsActionBase >().

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

◆ displacementIntegrityCheck()

void ComputeStrainBase::displacementIntegrityCheck ( )
protectedvirtualinherited

Reimplemented in Compute2DFiniteStrain, Compute2DIncrementalStrain, and Compute2DSmallStrain.

Definition at line 86 of file ComputeStrainBase.C.

Referenced by ComputeStrainBase::initialSetup().

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

◆ initialSetup()

void ComputeIncrementalStrainBase::initialSetup ( )
overrideprotectedinherited

Definition at line 37 of file ComputeIncrementalStrainBase.C.

Referenced by ComputeAxisymmetric1DFiniteStrain::initialSetup(), ComputeAxisymmetricRZIncrementalStrain::initialSetup(), ComputeRSphericalIncrementalStrain::initialSetup(), ComputeRSphericalFiniteStrain::initialSetup(), and ComputeAxisymmetric1DIncrementalStrain::initialSetup().

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

◆ initQpStatefulProperties()

void ComputeIncrementalStrainBase::initQpStatefulProperties ( )
overrideprotectedvirtualinherited

Reimplemented from ComputeStrainBase.

Reimplemented in ComputeCosseratIncrementalSmallStrain.

Definition at line 50 of file ComputeIncrementalStrainBase.C.

Referenced by ComputeCosseratIncrementalSmallStrain::initQpStatefulProperties().

51 {
52  _mechanical_strain[_qp].zero();
53  _total_strain[_qp].zero();
54  _deformation_gradient[_qp].zero();
55  _deformation_gradient[_qp].addIa(1.0);
56 
57  // Note that for some models (small strain), the rotation increment is
58  // never updated. Because we always have stateful properties, this method
59  // always gets called, so we can rely on this getting set here without
60  // setting it again when properties get computed.
61  _rotation_increment[_qp].zero();
62  _rotation_increment[_qp].addIa(1.0);
63 }
MaterialProperty< RankTwoTensor > & _deformation_gradient
MaterialProperty< RankTwoTensor > & _mechanical_strain
MaterialProperty< RankTwoTensor > & _rotation_increment
MaterialProperty< RankTwoTensor > & _total_strain

◆ subtractEigenstrainIncrementFromStrain()

void ComputeIncrementalStrainBase::subtractEigenstrainIncrementFromStrain ( RankTwoTensor &  strain)
protectedinherited

Definition at line 66 of file ComputeIncrementalStrainBase.C.

Referenced by ComputeIncrementalSmallStrain::computeProperties(), ComputeCosseratIncrementalSmallStrain::computeQpProperties(), and computeQpStrain().

67 {
68  for (unsigned int i = 0; i < _eigenstrains.size(); ++i)
69  {
70  strain -= (*_eigenstrains[i])[_qp];
71  strain += (*_eigenstrains_old[i])[_qp];
72  }
73 }
std::vector< const MaterialProperty< RankTwoTensor > * > _eigenstrains_old
std::vector< const MaterialProperty< RankTwoTensor > * > _eigenstrains

Member Data Documentation

◆ _base_name

std::string ComputeStrainBase::_base_name
protectedinherited

Definition at line 43 of file ComputeStrainBase.h.

Referenced by ComputeStrainBase::ComputeStrainBase().

◆ _current_elem_volume

const Real& ComputeStrainBase::_current_elem_volume
protectedinherited

◆ _decomposition_method

const DecompMethod ComputeFiniteStrain::_decomposition_method
private

Definition at line 45 of file ComputeFiniteStrain.h.

Referenced by computeQpIncrements().

◆ _deformation_gradient

MaterialProperty<RankTwoTensor>& ComputeIncrementalStrainBase::_deformation_gradient
protectedinherited

◆ _disp

std::vector<const VariableValue *> ComputeStrainBase::_disp
protectedinherited

◆ _eigenstrain_names

std::vector<MaterialPropertyName> ComputeStrainBase::_eigenstrain_names
protectedinherited

◆ _eigenstrains

std::vector<const MaterialProperty<RankTwoTensor> *> ComputeStrainBase::_eigenstrains
protectedinherited

◆ _eigenstrains_old

std::vector<const MaterialProperty<RankTwoTensor> *> ComputeIncrementalStrainBase::_eigenstrains_old
protectedinherited

◆ _Fhat

std::vector<RankTwoTensor> ComputeFiniteStrain::_Fhat
protected

◆ _global_strain

const MaterialProperty<RankTwoTensor>* ComputeStrainBase::_global_strain
protectedinherited

Definition at line 52 of file ComputeStrainBase.h.

Referenced by ComputeSmallStrain::computeProperties(), and computeQpStrain().

◆ _grad_disp

std::vector<const VariableGradient *> ComputeStrainBase::_grad_disp
protectedinherited

◆ _grad_disp_old

std::vector<const VariableGradient *> ComputeIncrementalStrainBase::_grad_disp_old
protectedinherited

◆ _mechanical_strain

MaterialProperty<RankTwoTensor>& ComputeStrainBase::_mechanical_strain
protectedinherited

◆ _mechanical_strain_old

const MaterialProperty<RankTwoTensor>& ComputeIncrementalStrainBase::_mechanical_strain_old
protectedinherited

◆ _ndisp

unsigned int ComputeStrainBase::_ndisp
protectedinherited

◆ _rotation_increment

MaterialProperty<RankTwoTensor>& ComputeIncrementalStrainBase::_rotation_increment
protectedinherited

◆ _strain_increment

MaterialProperty<RankTwoTensor>& ComputeIncrementalStrainBase::_strain_increment
protectedinherited

◆ _strain_rate

MaterialProperty<RankTwoTensor>& ComputeIncrementalStrainBase::_strain_rate
protectedinherited

◆ _total_strain

MaterialProperty<RankTwoTensor>& ComputeStrainBase::_total_strain
protectedinherited

◆ _total_strain_old

const MaterialProperty<RankTwoTensor>& ComputeIncrementalStrainBase::_total_strain_old
protectedinherited

◆ _volumetric_locking_correction

bool ComputeStrainBase::_volumetric_locking_correction
protectedinherited

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