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)
 
virtual void initQpStatefulProperties () override
 
void subtractEigenstrainIncrementFromStrain (RankTwoTensor &strain)
 

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
 
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 15 of file ComputeFiniteStrain.h.

Member Enumeration Documentation

Enumerator
TaylorExpansion 
EigenSolution 

Definition at line 31 of file ComputeFiniteStrain.h.

32  {
33  TaylorExpansion,
34  EigenSolution
35  };

Constructor & Destructor Documentation

ComputeFiniteStrain::ComputeFiniteStrain ( const InputParameters &  parameters)

Definition at line 33 of file ComputeFiniteStrain.C.

34  : ComputeIncrementalStrainBase(parameters),
35  _Fhat(_fe_problem.getMaxQps()),
36  _decomposition_method(getParam<MooseEnum>("decomposition_method").getEnum<DecompMethod>())
37 {
38 }
const DecompMethod _decomposition_method
ComputeIncrementalStrainBase(const InputParameters &parameters)
std::vector< RankTwoTensor > _Fhat

Member Function Documentation

void ComputeFiniteStrain::computeProperties ( )
virtual

Reimplemented in ComputeRSphericalFiniteStrain, Compute2DFiniteStrain, and Compute1DFiniteStrain.

Definition at line 41 of file ComputeFiniteStrain.C.

42 {
43  RankTwoTensor ave_Fhat;
44  Real ave_dfgrd_det = 0.0;
45  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
46  {
47  // Deformation gradient
48  RankTwoTensor A((*_grad_disp[0])[_qp],
49  (*_grad_disp[1])[_qp],
50  (*_grad_disp[2])[_qp]); // Deformation gradient
51  RankTwoTensor Fbar((*_grad_disp_old[0])[_qp],
52  (*_grad_disp_old[1])[_qp],
53  (*_grad_disp_old[2])[_qp]); // Old Deformation gradient
54 
55  _deformation_gradient[_qp] = A;
56  _deformation_gradient[_qp].addIa(1.0); // Gauss point deformation gradient
57 
58  // A = gradU - gradUold
59  A -= Fbar;
60 
61  // Fbar = ( I + gradUold)
62  Fbar.addIa(1.0);
63 
64  // Incremental deformation gradient _Fhat = I + A Fbar^-1
65  _Fhat[_qp] = A * Fbar.inverse();
66  _Fhat[_qp].addIa(1.0);
67 
69  {
70  // Calculate average _Fhat for volumetric locking correction
71  ave_Fhat += _Fhat[_qp] * _JxW[_qp] * _coord[_qp];
72 
73  // Average deformation gradient
74  ave_dfgrd_det += _deformation_gradient[_qp].det() * _JxW[_qp] * _coord[_qp];
75  }
76  }
77 
79  {
80  // needed for volumetric locking correction
81  ave_Fhat /= _current_elem_volume;
82  // average deformation gradient
83  ave_dfgrd_det /= _current_elem_volume;
84  }
85 
86  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
87  {
89  {
90  // Finalize volumetric locking correction
91  _Fhat[_qp] *= std::cbrt(ave_Fhat.det() / _Fhat[_qp].det());
92  _deformation_gradient[_qp] *= std::cbrt(ave_dfgrd_det / _deformation_gradient[_qp].det());
93  }
94 
96  }
97 }
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
void ComputeFiniteStrain::computeQpIncrements ( RankTwoTensor &  e,
RankTwoTensor &  r 
)
protectedvirtual

Definition at line 129 of file ComputeFiniteStrain.C.

Referenced by computeQpStrain().

131 {
132  switch (_decomposition_method)
133  {
135  {
136  // inverse of _Fhat
137  RankTwoTensor invFhat(_Fhat[_qp].inverse());
138 
139  // A = I - _Fhat^-1
140  RankTwoTensor A(RankTwoTensor::initIdentity);
141  A -= invFhat;
142 
143  // Cinv - I = A A^T - A - A^T;
144  RankTwoTensor Cinv_I = A * A.transpose() - A - A.transpose();
145 
146  // strain rate D from Taylor expansion, Chat = (-1/2(Chat^-1 - I) + 1/4*(Chat^-1 - I)^2 + ...
147  total_strain_increment = -Cinv_I * 0.5 + Cinv_I * Cinv_I * 0.25;
148 
149  const Real a[3] = {invFhat(1, 2) - invFhat(2, 1),
150  invFhat(2, 0) - invFhat(0, 2),
151  invFhat(0, 1) - invFhat(1, 0)};
152 
153  Real q = (a[0] * a[0] + a[1] * a[1] + a[2] * a[2]) / 4.0;
154  Real trFhatinv_1 = invFhat.trace() - 1.0;
155  const Real p = trFhatinv_1 * trFhatinv_1 / 4.0;
156 
157  // cos theta_a
158  const Real C1 =
159  std::sqrt(p + 3.0 * Utility::pow<2>(p) * (1.0 - (p + q)) / Utility::pow<2>(p + q) -
160  2.0 * Utility::pow<3>(p) * (1.0 - (p + q)) / Utility::pow<3>(p + q));
161 
162  Real C2;
163  if (q > 0.01)
164  // (1-cos theta_a)/4q
165  C2 = (1.0 - C1) / (4.0 * q);
166  else
167  // alternate form for small q
168  C2 = 0.125 + q * 0.03125 * (Utility::pow<2>(p) - 12.0 * (p - 1.0)) / Utility::pow<2>(p) +
169  Utility::pow<2>(q) * (p - 2.0) * (Utility::pow<2>(p) - 10.0 * p + 32.0) /
170  Utility::pow<3>(p) +
171  Utility::pow<3>(q) * (1104.0 - 992.0 * p + 376.0 * Utility::pow<2>(p) -
172  72.0 * Utility::pow<3>(p) + 5.0 * Utility::pow<4>(p)) /
173  (512.0 * Utility::pow<4>(p));
174  const Real 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  RankTwoTensor 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<Real> e_value(3);
200  RankTwoTensor e_vector, N1, N2, N3;
201 
202  RankTwoTensor Chat = _Fhat[_qp].transpose() * _Fhat[_qp];
203  Chat.symmetricEigenvaluesEigenvectors(e_value, e_vector);
204 
205  const Real lambda1 = std::sqrt(e_value[0]);
206  const Real lambda2 = std::sqrt(e_value[1]);
207  const Real 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  RankTwoTensor Uhat = N1 * lambda1 + N2 * lambda2 + N3 * lambda3;
214  RankTwoTensor 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("ComputeFiniteStrain Error: Pass valid decomposition type: TaylorExpansion or "
225  "EigenSolution.");
226  }
227 }
const DecompMethod _decomposition_method
std::vector< RankTwoTensor > _Fhat
void ComputeFiniteStrain::computeQpStrain ( )
protectedvirtual

Definition at line 100 of file ComputeFiniteStrain.C.

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

101 {
102  RankTwoTensor total_strain_increment;
103 
104  // two ways to calculate these increments: TaylorExpansion(default) or EigenSolution
105  computeQpIncrements(total_strain_increment, _rotation_increment[_qp]);
106 
107  _strain_increment[_qp] = total_strain_increment;
108 
109  // Remove the eigenstrain increment
111 
112  if (_dt > 0)
113  _strain_rate[_qp] = _strain_increment[_qp] / _dt;
114  else
115  _strain_rate[_qp].zero();
116 
117  // Update strain in intermediate configuration
119  _total_strain[_qp] = _total_strain_old[_qp] + total_strain_increment;
120 
121  // Rotate strain to current configuration
122  _mechanical_strain[_qp] =
123  _rotation_increment[_qp] * _mechanical_strain[_qp] * _rotation_increment[_qp].transpose();
124  _total_strain[_qp] =
125  _rotation_increment[_qp] * _total_strain[_qp] * _rotation_increment[_qp].transpose();
126 }
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
MaterialProperty< RankTwoTensor > & _rotation_increment
MaterialProperty< RankTwoTensor > & _total_strain
MooseEnum ComputeFiniteStrain::decompositionType ( )
static

Definition at line 15 of file ComputeFiniteStrain.C.

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

16 {
17  return MooseEnum("TaylorExpansion EigenSolution", "TaylorExpansion");
18 }
void ComputeIncrementalStrainBase::initQpStatefulProperties ( )
overrideprotectedvirtualinherited

Reimplemented from ComputeStrainBase.

Reimplemented in ComputeCosseratIncrementalSmallStrain.

Definition at line 44 of file ComputeIncrementalStrainBase.C.

Referenced by ComputeCosseratIncrementalSmallStrain::initQpStatefulProperties(), and ComputeIncrementalStrainBase::~ComputeIncrementalStrainBase().

45 {
46  _mechanical_strain[_qp].zero();
47  _total_strain[_qp].zero();
48  _deformation_gradient[_qp].zero();
49  _deformation_gradient[_qp].addIa(1.0);
50 
51  // Note that for some models (small strain), the rotation increment is
52  // never updated. Because we always have stateful properties, this method
53  // always gets called, so we can rely on this getting set here without
54  // setting it again when properties get computed.
55  _rotation_increment[_qp].zero();
56  _rotation_increment[_qp].addIa(1.0);
57 }
MaterialProperty< RankTwoTensor > & _deformation_gradient
MaterialProperty< RankTwoTensor > & _mechanical_strain
MaterialProperty< RankTwoTensor > & _rotation_increment
MaterialProperty< RankTwoTensor > & _total_strain
void ComputeIncrementalStrainBase::subtractEigenstrainIncrementFromStrain ( RankTwoTensor &  strain)
protectedinherited

Definition at line 60 of file ComputeIncrementalStrainBase.C.

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

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

Member Data Documentation

std::string ComputeStrainBase::_base_name
protectedinherited

Definition at line 33 of file ComputeStrainBase.h.

Referenced by ComputeStrainBase::ComputeStrainBase().

const Real& ComputeStrainBase::_current_elem_volume
protectedinherited
const DecompMethod ComputeFiniteStrain::_decomposition_method
private

Definition at line 37 of file ComputeFiniteStrain.h.

Referenced by computeQpIncrements().

MaterialProperty<RankTwoTensor>& ComputeIncrementalStrainBase::_deformation_gradient
protectedinherited
std::vector<const VariableValue *> ComputeStrainBase::_disp
protectedinherited
std::vector<MaterialPropertyName> ComputeStrainBase::_eigenstrain_names
protectedinherited
std::vector<const MaterialProperty<RankTwoTensor> *> ComputeStrainBase::_eigenstrains
protectedinherited
std::vector<const MaterialProperty<RankTwoTensor> *> ComputeIncrementalStrainBase::_eigenstrains_old
protectedinherited
std::vector<RankTwoTensor> ComputeFiniteStrain::_Fhat
protected
std::vector<const VariableGradient *> ComputeStrainBase::_grad_disp
protectedinherited
std::vector<const VariableGradient *> ComputeIncrementalStrainBase::_grad_disp_old
protectedinherited
MaterialProperty<RankTwoTensor>& ComputeStrainBase::_mechanical_strain
protectedinherited
const MaterialProperty<RankTwoTensor>& ComputeIncrementalStrainBase::_mechanical_strain_old
protectedinherited
unsigned int ComputeStrainBase::_ndisp
protectedinherited
MaterialProperty<RankTwoTensor>& ComputeIncrementalStrainBase::_rotation_increment
protectedinherited
MaterialProperty<RankTwoTensor>& ComputeIncrementalStrainBase::_strain_increment
protectedinherited
MaterialProperty<RankTwoTensor>& ComputeIncrementalStrainBase::_strain_rate
protectedinherited
MaterialProperty<RankTwoTensor>& ComputeStrainBase::_total_strain
protectedinherited
const MaterialProperty<RankTwoTensor>& ComputeIncrementalStrainBase::_total_strain_old
protectedinherited
bool ComputeStrainBase::_volumetric_locking_correction
protectedinherited

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