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

Material class for bond-associated correspondence material model for finite strain. More...

#include <ComputeFiniteStrainNOSPD.h>

Inheritance diagram for ComputeFiniteStrainNOSPD:
[legend]

Public Member Functions

 ComputeFiniteStrainNOSPD (const InputParameters &parameters)
 
virtual void initQpStatefulProperties () override
 

Static Public Member Functions

static MooseEnum decompositionType ()
 

Protected Member Functions

virtual void computeQpStrain () override
 Function to compute strain tensors. More...
 
virtual void computeQpFhat ()
 
void computeQpStrainRotationIncrements (RankTwoTensor &e, RankTwoTensor &r)
 Function to compute strain and rotational increments. More...
 
void subtractEigenstrainIncrementFromStrain (RankTwoTensor &strain)
 Function to compute the mechanical strain tensor by subtracting thermal strain from the total strain. More...
 
virtual void computeProperties () override
 
virtual void computeBondStretch () override
 
virtual void computeQpDeformationGradient ()
 Function to compute bond-associated deformation gradient. More...
 

Protected Attributes

std::vector< RankTwoTensor_Fhat
 'Incremental' deformation gradient More...
 
const bool _plane_strain
 Plane strain problem or not, this is only used for mechanical stretch calculation. More...
 
MaterialProperty< RankTwoTensor > & _strain_rate
 Material properties to store. More...
 
MaterialProperty< RankTwoTensor > & _strain_increment
 
MaterialProperty< RankTwoTensor > & _rotation_increment
 
const MaterialProperty< RankTwoTensor > & _deformation_gradient_old
 Material properties to fetch. More...
 
const MaterialProperty< RankTwoTensor > & _mechanical_strain_old
 
const MaterialProperty< RankTwoTensor > & _total_strain_old
 
std::vector< const MaterialProperty< RankTwoTensor > * > _eigenstrains_old
 
const MaterialProperty< RankFourTensor > & _Cijkl
 Material properties to fetch. More...
 
std::vector< MaterialPropertyName > _eigenstrain_names
 
std::vector< const MaterialProperty< RankTwoTensor > * > _eigenstrains
 
MaterialProperty< RankTwoTensor > & _shape2
 Material properties to store. More...
 
MaterialProperty< RankTwoTensor > & _deformation_gradient
 
MaterialProperty< RankTwoTensor > & _ddgraddu
 
MaterialProperty< RankTwoTensor > & _ddgraddv
 
MaterialProperty< RankTwoTensor > & _ddgraddw
 
MaterialProperty< RankTwoTensor > & _total_strain
 
MaterialProperty< RankTwoTensor > & _mechanical_strain
 
MaterialProperty< Real > & _multi
 
enum  DecompMethod { DecompMethod::TaylorExpansion, DecompMethod::EigenSolution }
 Method to decompose into rotation increment and strain increment. More...
 
const DecompMethod _decomposition_method
 

Detailed Description

Material class for bond-associated correspondence material model for finite strain.

Definition at line 23 of file ComputeFiniteStrainNOSPD.h.

Member Enumeration Documentation

◆ DecompMethod

Method to decompose into rotation increment and strain increment.

Enumerator
TaylorExpansion 
EigenSolution 

Definition at line 63 of file ComputeFiniteStrainNOSPD.h.

64  {
65  TaylorExpansion,
66  EigenSolution
67  };

Constructor & Destructor Documentation

◆ ComputeFiniteStrainNOSPD()

ComputeFiniteStrainNOSPD::ComputeFiniteStrainNOSPD ( const InputParameters &  parameters)

Definition at line 39 of file ComputeFiniteStrainNOSPD.C.

40  : ComputeStrainBaseNOSPD(parameters),
41  _strain_rate(declareProperty<RankTwoTensor>("strain_rate")),
42  _strain_increment(declareProperty<RankTwoTensor>("strain_increment")),
43  _rotation_increment(declareProperty<RankTwoTensor>("rotation_increment")),
44  _deformation_gradient_old(getMaterialPropertyOld<RankTwoTensor>("deformation_gradient")),
45  _mechanical_strain_old(getMaterialPropertyOld<RankTwoTensor>("mechanical_strain")),
46  _total_strain_old(getMaterialPropertyOld<RankTwoTensor>("total_strain")),
48  _Fhat(2),
49  _decomposition_method(getParam<MooseEnum>("decomposition_method").getEnum<DecompMethod>())
50 {
51  for (unsigned int i = 0; i < _eigenstrains_old.size(); ++i)
52  _eigenstrains_old[i] = &getMaterialPropertyOld<RankTwoTensor>(_eigenstrain_names[i]);
53 }

Member Function Documentation

◆ computeBondStretch()

void ComputeStrainBaseNOSPD::computeBondStretch ( )
overrideprotectedvirtualinherited

Definition at line 156 of file ComputeStrainBaseNOSPD.C.

157 {
158  _total_stretch[0] = _current_length / _origin_length - 1.0;
159  _total_stretch[1] = _total_stretch[0];
160  _mechanical_stretch[0] = _total_stretch[0];
161 
162  Real factor = 1.0;
163  if (_plane_strain)
165 
166  for (auto es : _eigenstrains)
167  _mechanical_stretch[0] -= 0.5 * factor * ((*es)[0](2, 2) + (*es)[1](2, 2));
168 
169  _mechanical_stretch[1] = _mechanical_stretch[0];
170 }

Referenced by ComputeStrainBaseNOSPD::computeProperties().

◆ computeProperties()

void ComputeStrainBaseNOSPD::computeProperties ( )
overrideprotectedvirtualinherited

Definition at line 145 of file ComputeStrainBaseNOSPD.C.

146 {
147  setupMeshRelatedData(); // function from base class
148  computeBondCurrentLength(); // current length of a bond from base class
149  computeBondStretch(); // stretch of a bond
150 
151  for (_qp = 0; _qp < _nnodes; ++_qp)
152  computeQpStrain();
153 }

◆ computeQpDeformationGradient()

void ComputeStrainBaseNOSPD::computeQpDeformationGradient ( )
protectedvirtualinherited

Function to compute bond-associated deformation gradient.

Reimplemented in ComputeForceStabilizedSmallStrainNOSPD.

Definition at line 71 of file ComputeStrainBaseNOSPD.C.

72 {
73  _shape2[_qp].zero();
74  _deformation_gradient[_qp].zero();
75  _ddgraddu[_qp].zero();
76  _ddgraddv[_qp].zero();
77  _ddgraddw[_qp].zero();
78  _multi[_qp] = 0.0;
79 
80  // for cases when current bond was broken, assign shape tensor and deformation gradient to unity
81  if (_bond_status_var->getElementalValue(_current_elem) < 0.5)
82  {
83  _shape2[_qp] = RankTwoTensor::initIdentity;
84  _deformation_gradient[_qp] = _shape2[_qp];
85  }
86  else
87  {
88  if (_dim == 2)
89  _shape2[_qp](2, 2) = _deformation_gradient[_qp](2, 2) = 1.0;
90 
91  const Node * cur_nd = _current_elem->node_ptr(_qp);
92  const Node * end_nd = _current_elem->node_ptr(1 - _qp); // two nodes for edge2 element
93  std::vector<dof_id_type> neighbors = _pdmesh.getNeighbors(cur_nd->id());
94  std::vector<dof_id_type> bonds = _pdmesh.getBonds(cur_nd->id());
95 
96  unsigned int nb =
97  std::find(neighbors.begin(), neighbors.end(), end_nd->id()) - neighbors.begin();
98  std::vector<dof_id_type> dg_neighbors = _pdmesh.getDefGradNeighbors(cur_nd->id(), nb);
99 
100  // calculate the shape tensor and prepare the deformation gradient tensor
101  Real dgnodes_vsum = 0.0;
102  RealGradient ori_vec(_dim), cur_vec(_dim);
103 
104  for (unsigned int j = 0; j < dg_neighbors.size(); ++j)
105  if (_bond_status_var->getElementalValue(_pdmesh.elemPtr(bonds[dg_neighbors[j]])) > 0.5)
106  {
107  Node * node_j = _pdmesh.nodePtr(neighbors[dg_neighbors[j]]);
108  Real vol_j = _pdmesh.getPDNodeVolume(neighbors[dg_neighbors[j]]);
109  dgnodes_vsum += vol_j;
110  ori_vec = *node_j - *_pdmesh.nodePtr(cur_nd->id());
111 
112  for (unsigned int k = 0; k < _dim; ++k)
113  cur_vec(k) = ori_vec(k) + _disp_var[k]->getNodalValue(*node_j) -
114  _disp_var[k]->getNodalValue(*cur_nd);
115 
116  Real ori_len = ori_vec.norm();
117  for (unsigned int k = 0; k < _dim; ++k)
118  {
119  for (unsigned int l = 0; l < _dim; ++l)
120  {
121  _shape2[_qp](k, l) += _horiz_rad[_qp] / ori_len * ori_vec(k) * ori_vec(l) * vol_j;
122  _deformation_gradient[_qp](k, l) +=
123  _horiz_rad[_qp] / ori_len * cur_vec(k) * ori_vec(l) * vol_j;
124  }
125  // calculate derivatives of deformation_gradient w.r.t displacements of node i
126  _ddgraddu[_qp](0, k) += -_horiz_rad[_qp] / ori_len * ori_vec(k) * vol_j;
127  _ddgraddv[_qp](1, k) += -_horiz_rad[_qp] / ori_len * ori_vec(k) * vol_j;
128  if (_dim == 3)
129  _ddgraddw[_qp](2, k) += -_horiz_rad[_qp] / ori_len * ori_vec(k) * vol_j;
130  }
131  }
132  // finalize the deformation gradient and its derivatives
133  _deformation_gradient[_qp] *= _shape2[_qp].inverse();
134  _ddgraddu[_qp] *= _shape2[_qp].inverse();
135  _ddgraddv[_qp] *= _shape2[_qp].inverse();
136  _ddgraddw[_qp] *= _shape2[_qp].inverse();
137 
138  // force state multiplier
139  _multi[_qp] = _horiz_rad[_qp] / _origin_length * _node_vol[0] * _node_vol[1] * dgnodes_vsum /
140  _horiz_vol[_qp];
141  }
142 }

Referenced by ComputeSmallStrainNOSPD::computeQpStrain(), and computeQpStrain().

◆ computeQpFhat()

void ComputeFiniteStrainNOSPD::computeQpFhat ( )
protectedvirtual

Reimplemented in ComputePlaneFiniteStrainNOSPD.

Definition at line 89 of file ComputeFiniteStrainNOSPD.C.

90 {
91  // Incremental deformation gradient of current step w.r.t previous step:
92  // _Fhat = deformation_gradient * inv(deformation_gradient_old)
93  _Fhat[_qp] = _deformation_gradient[_qp] * _deformation_gradient_old[_qp].inverse();
94 }

Referenced by computeQpStrain().

◆ computeQpStrain()

void ComputeFiniteStrainNOSPD::computeQpStrain ( )
overrideprotectedvirtual

Function to compute strain tensors.

Implements ComputeStrainBaseNOSPD.

Definition at line 56 of file ComputeFiniteStrainNOSPD.C.

57 {
59 
60  computeQpFhat();
61 
62  RankTwoTensor total_strain_increment;
63 
64  // Two ways to calculate these increments: TaylorExpansion(default) or EigenSolution
65  computeQpStrainRotationIncrements(total_strain_increment, _rotation_increment[_qp]);
66 
67  _strain_increment[_qp] = total_strain_increment;
68 
69  // Remove the eigenstrain increment
71 
72  if (_dt > 0)
73  _strain_rate[_qp] = _strain_increment[_qp] / _dt;
74  else
75  _strain_rate[_qp].zero();
76 
77  // Update strain in intermediate configuration
79  _total_strain[_qp] = _total_strain_old[_qp] + total_strain_increment;
80 
81  // Rotate strain to current configuration
82  _mechanical_strain[_qp] =
83  _rotation_increment[_qp] * _mechanical_strain[_qp] * _rotation_increment[_qp].transpose();
84  _total_strain[_qp] =
85  _rotation_increment[_qp] * _total_strain[_qp] * _rotation_increment[_qp].transpose();
86 }

◆ computeQpStrainRotationIncrements()

void ComputeFiniteStrainNOSPD::computeQpStrainRotationIncrements ( RankTwoTensor e,
RankTwoTensor r 
)
protected

Function to compute strain and rotational increments.

Definition at line 97 of file ComputeFiniteStrainNOSPD.C.

99 {
100  switch (_decomposition_method)
101  {
103  {
104  // inverse of _Fhat
105  RankTwoTensor invFhat(_Fhat[_qp].inverse());
106 
107  // A = I - _Fhat^-1
108  RankTwoTensor A(RankTwoTensor::initIdentity);
109  A -= invFhat;
110 
111  // Cinv - I = A A^T - A - A^T;
112  RankTwoTensor Cinv_I = A * A.transpose() - A - A.transpose();
113 
114  // strain rate D from Taylor expansion, Chat = (-1/2(Chat^-1 - I) + 1/4*(Chat^-1 - I)^2 + ...
115  total_strain_increment = -Cinv_I * 0.5 + Cinv_I * Cinv_I * 0.25;
116 
117  const Real a[3] = {invFhat(1, 2) - invFhat(2, 1),
118  invFhat(2, 0) - invFhat(0, 2),
119  invFhat(0, 1) - invFhat(1, 0)};
120 
121  Real q = (a[0] * a[0] + a[1] * a[1] + a[2] * a[2]) / 4.0;
122  Real trFhatinv_1 = invFhat.trace() - 1.0;
123  const Real p = trFhatinv_1 * trFhatinv_1 / 4.0;
124 
125  // cos theta_a
126  const Real C1 =
127  std::sqrt(p + 3.0 * Utility::pow<2>(p) * (1.0 - (p + q)) / Utility::pow<2>(p + q) -
128  2.0 * Utility::pow<3>(p) * (1.0 - (p + q)) / Utility::pow<3>(p + q));
129 
130  Real C2;
131  if (q > 0.01)
132  // (1-cos theta_a)/4q
133  C2 = (1.0 - C1) / (4.0 * q);
134  else
135  // alternate form for small q
136  C2 = 0.125 + q * 0.03125 * (Utility::pow<2>(p) - 12.0 * (p - 1.0)) / Utility::pow<2>(p) +
137  Utility::pow<2>(q) * (p - 2.0) * (Utility::pow<2>(p) - 10.0 * p + 32.0) /
138  Utility::pow<3>(p) +
139  Utility::pow<3>(q) *
140  (1104.0 - 992.0 * p + 376.0 * Utility::pow<2>(p) - 72.0 * Utility::pow<3>(p) +
141  5.0 * Utility::pow<4>(p)) /
142  (512.0 * Utility::pow<4>(p));
143  const Real C3 =
144  0.5 * std::sqrt((p * q * (3.0 - q) + Utility::pow<3>(p) + Utility::pow<2>(q)) /
145  Utility::pow<3>(p + q)); // sin theta_a/(2 sqrt(q))
146 
147  // Calculate incremental rotation. Note that this value is the transpose of that from Rashid,
148  // 93, so we transpose it before storing
149  RankTwoTensor R_incr;
150  R_incr.addIa(C1);
151  for (unsigned int i = 0; i < 3; ++i)
152  for (unsigned int j = 0; j < 3; ++j)
153  R_incr(i, j) += C2 * a[i] * a[j];
154 
155  R_incr(0, 1) += C3 * a[2];
156  R_incr(0, 2) -= C3 * a[1];
157  R_incr(1, 0) -= C3 * a[2];
158  R_incr(1, 2) += C3 * a[0];
159  R_incr(2, 0) += C3 * a[1];
160  R_incr(2, 1) -= C3 * a[0];
161 
162  rotation_increment = R_incr.transpose();
163  break;
164  }
166  {
167  std::vector<Real> e_value(3);
168  RankTwoTensor e_vector, N1, N2, N3;
169 
170  RankTwoTensor Chat = _Fhat[_qp].transpose() * _Fhat[_qp];
171  Chat.symmetricEigenvaluesEigenvectors(e_value, e_vector);
172 
173  const Real lambda1 = std::sqrt(e_value[0]);
174  const Real lambda2 = std::sqrt(e_value[1]);
175  const Real lambda3 = std::sqrt(e_value[2]);
176 
177  N1.vectorOuterProduct(e_vector.column(0), e_vector.column(0));
178  N2.vectorOuterProduct(e_vector.column(1), e_vector.column(1));
179  N3.vectorOuterProduct(e_vector.column(2), e_vector.column(2));
180 
181  RankTwoTensor Uhat = N1 * lambda1 + N2 * lambda2 + N3 * lambda3;
182  RankTwoTensor invUhat(Uhat.inverse());
183 
184  rotation_increment = _Fhat[_qp] * invUhat;
185 
186  total_strain_increment =
187  N1 * std::log(lambda1) + N2 * std::log(lambda2) + N3 * std::log(lambda3);
188  break;
189  }
190 
191  default:
192  mooseError("ComputeFiniteStrainNOSPD Error: Invalid decomposition type! Please specify : "
193  "TaylorExpansion or "
194  "EigenSolution.");
195  }
196 }

Referenced by computeQpStrain().

◆ decompositionType()

MooseEnum ComputeFiniteStrainNOSPD::decompositionType ( )
static

Definition at line 17 of file ComputeFiniteStrainNOSPD.C.

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

Referenced by validParams< ComputeFiniteStrainNOSPD >().

◆ initQpStatefulProperties()

void ComputeStrainBaseNOSPD::initQpStatefulProperties ( )
overridevirtualinherited

Definition at line 57 of file ComputeStrainBaseNOSPD.C.

58 {
59  _mechanical_strain[_qp].zero();
60  _total_strain[_qp].zero();
61  _deformation_gradient[_qp].zero();
62  _deformation_gradient[_qp].addIa(1.0);
63 
64  if (_qrule->n_points() < 2) // Gauss_Lobatto: order = 2n-3 (n is number of qp)
65  mooseError(
66  "NOSPD models require Gauss_Lobatto rule and a minimum of 2 quadrature points, i.e., "
67  "order of FIRST");
68 }

◆ subtractEigenstrainIncrementFromStrain()

void ComputeFiniteStrainNOSPD::subtractEigenstrainIncrementFromStrain ( RankTwoTensor strain)
protected

Function to compute the mechanical strain tensor by subtracting thermal strain from the total strain.

Definition at line 199 of file ComputeFiniteStrainNOSPD.C.

200 {
201  for (unsigned int i = 0; i < _eigenstrains.size(); ++i)
202  {
203  strain -= (*_eigenstrains[i])[_qp];
204  strain += (*_eigenstrains_old[i])[_qp];
205  }
206 }

Referenced by computeQpStrain().

Member Data Documentation

◆ _Cijkl

const MaterialProperty<RankFourTensor>& ComputeStrainBaseNOSPD::_Cijkl
protectedinherited

Material properties to fetch.

Definition at line 48 of file ComputeStrainBaseNOSPD.h.

Referenced by ComputeStrainBaseNOSPD::computeBondStretch().

◆ _ddgraddu

MaterialProperty<RankTwoTensor>& ComputeStrainBaseNOSPD::_ddgraddu
protectedinherited

◆ _ddgraddv

MaterialProperty<RankTwoTensor>& ComputeStrainBaseNOSPD::_ddgraddv
protectedinherited

◆ _ddgraddw

MaterialProperty<RankTwoTensor>& ComputeStrainBaseNOSPD::_ddgraddw
protectedinherited

◆ _decomposition_method

const DecompMethod ComputeFiniteStrainNOSPD::_decomposition_method
private

Definition at line 68 of file ComputeFiniteStrainNOSPD.h.

Referenced by computeQpStrainRotationIncrements().

◆ _deformation_gradient

MaterialProperty<RankTwoTensor>& ComputeStrainBaseNOSPD::_deformation_gradient
protectedinherited

◆ _deformation_gradient_old

const MaterialProperty<RankTwoTensor>& ComputeFiniteStrainNOSPD::_deformation_gradient_old
protected

Material properties to fetch.

Definition at line 52 of file ComputeFiniteStrainNOSPD.h.

Referenced by ComputePlaneFiniteStrainNOSPD::computeQpFhat(), and computeQpFhat().

◆ _eigenstrain_names

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

◆ _eigenstrains

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

◆ _eigenstrains_old

std::vector<const MaterialProperty<RankTwoTensor> *> ComputeFiniteStrainNOSPD::_eigenstrains_old
protected

◆ _Fhat

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

'Incremental' deformation gradient

Definition at line 59 of file ComputeFiniteStrainNOSPD.h.

Referenced by ComputePlaneFiniteStrainNOSPD::computeQpFhat(), computeQpFhat(), and computeQpStrainRotationIncrements().

◆ _mechanical_strain

MaterialProperty<RankTwoTensor>& ComputeStrainBaseNOSPD::_mechanical_strain
protectedinherited

◆ _mechanical_strain_old

const MaterialProperty<RankTwoTensor>& ComputeFiniteStrainNOSPD::_mechanical_strain_old
protected

Definition at line 53 of file ComputeFiniteStrainNOSPD.h.

Referenced by computeQpStrain().

◆ _multi

MaterialProperty<Real>& ComputeStrainBaseNOSPD::_multi
protectedinherited

◆ _plane_strain

const bool ComputeStrainBaseNOSPD::_plane_strain
protectedinherited

Plane strain problem or not, this is only used for mechanical stretch calculation.

Definition at line 45 of file ComputeStrainBaseNOSPD.h.

Referenced by ComputeStrainBaseNOSPD::computeBondStretch().

◆ _rotation_increment

MaterialProperty<RankTwoTensor>& ComputeFiniteStrainNOSPD::_rotation_increment
protected

Definition at line 48 of file ComputeFiniteStrainNOSPD.h.

Referenced by computeQpStrain().

◆ _shape2

MaterialProperty<RankTwoTensor>& ComputeStrainBaseNOSPD::_shape2
protectedinherited

◆ _strain_increment

MaterialProperty<RankTwoTensor>& ComputeFiniteStrainNOSPD::_strain_increment
protected

Definition at line 47 of file ComputeFiniteStrainNOSPD.h.

Referenced by computeQpStrain().

◆ _strain_rate

MaterialProperty<RankTwoTensor>& ComputeFiniteStrainNOSPD::_strain_rate
protected

Material properties to store.

Definition at line 46 of file ComputeFiniteStrainNOSPD.h.

Referenced by computeQpStrain().

◆ _total_strain

MaterialProperty<RankTwoTensor>& ComputeStrainBaseNOSPD::_total_strain
protectedinherited

◆ _total_strain_old

const MaterialProperty<RankTwoTensor>& ComputeFiniteStrainNOSPD::_total_strain_old
protected

Definition at line 54 of file ComputeFiniteStrainNOSPD.h.

Referenced by computeQpStrain().


The documentation for this class was generated from the following files:
ComputeStrainBaseNOSPD::_mechanical_strain
MaterialProperty< RankTwoTensor > & _mechanical_strain
Definition: ComputeStrainBaseNOSPD.h:62
ComputeFiniteStrainNOSPD::_mechanical_strain_old
const MaterialProperty< RankTwoTensor > & _mechanical_strain_old
Definition: ComputeFiniteStrainNOSPD.h:53
ComputeStrainBaseNOSPD::_ddgraddv
MaterialProperty< RankTwoTensor > & _ddgraddv
Definition: ComputeStrainBaseNOSPD.h:58
ComputeStrainBaseNOSPD::_shape2
MaterialProperty< RankTwoTensor > & _shape2
Material properties to store.
Definition: ComputeStrainBaseNOSPD.h:54
libMesh::RealGradient
VectorValue< Real > RealGradient
Definition: GrainForceAndTorqueInterface.h:17
ComputeFiniteStrainNOSPD::DecompMethod::EigenSolution
ComputeStrainBaseNOSPD::_plane_strain
const bool _plane_strain
Plane strain problem or not, this is only used for mechanical stretch calculation.
Definition: ComputeStrainBaseNOSPD.h:45
ComputeFiniteStrainNOSPD::_decomposition_method
const DecompMethod _decomposition_method
Definition: ComputeFiniteStrainNOSPD.h:68
ComputeFiniteStrainNOSPD::_strain_rate
MaterialProperty< RankTwoTensor > & _strain_rate
Material properties to store.
Definition: ComputeFiniteStrainNOSPD.h:46
ComputeFiniteStrainNOSPD::subtractEigenstrainIncrementFromStrain
void subtractEigenstrainIncrementFromStrain(RankTwoTensor &strain)
Function to compute the mechanical strain tensor by subtracting thermal strain from the total strain.
Definition: ComputeFiniteStrainNOSPD.C:199
ComputeFiniteStrainNOSPD::DecompMethod::TaylorExpansion
ComputeFiniteStrainNOSPD::_rotation_increment
MaterialProperty< RankTwoTensor > & _rotation_increment
Definition: ComputeFiniteStrainNOSPD.h:48
ComputeStrainBaseNOSPD::_deformation_gradient
MaterialProperty< RankTwoTensor > & _deformation_gradient
Definition: ComputeStrainBaseNOSPD.h:55
ElasticityTensorTools::getIsotropicPoissonsRatio
T getIsotropicPoissonsRatio(const RankFourTensorTempl< T > &elasticity_tensor)
Get the Poisson's modulus for an isotropic elasticity tensor param elasticity_tensor the tensor (must...
Definition: ElasticityTensorTools.h:115
ComputeStrainBaseNOSPD::_Cijkl
const MaterialProperty< RankFourTensor > & _Cijkl
Material properties to fetch.
Definition: ComputeStrainBaseNOSPD.h:48
ComputeStrainBaseNOSPD::computeQpDeformationGradient
virtual void computeQpDeformationGradient()
Function to compute bond-associated deformation gradient.
Definition: ComputeStrainBaseNOSPD.C:71
ComputeFiniteStrainNOSPD::_total_strain_old
const MaterialProperty< RankTwoTensor > & _total_strain_old
Definition: ComputeFiniteStrainNOSPD.h:54
ComputeStrainBaseNOSPD::_eigenstrain_names
std::vector< MaterialPropertyName > _eigenstrain_names
Definition: ComputeStrainBaseNOSPD.h:49
ComputeStrainBaseNOSPD::_multi
MaterialProperty< Real > & _multi
Definition: ComputeStrainBaseNOSPD.h:64
ComputeFiniteStrainNOSPD::_eigenstrains_old
std::vector< const MaterialProperty< RankTwoTensor > * > _eigenstrains_old
Definition: ComputeFiniteStrainNOSPD.h:55
ComputeStrainBaseNOSPD::ComputeStrainBaseNOSPD
ComputeStrainBaseNOSPD(const InputParameters &parameters)
Definition: ComputeStrainBaseNOSPD.C:34
ComputeStrainBaseNOSPD::_ddgraddw
MaterialProperty< RankTwoTensor > & _ddgraddw
Definition: ComputeStrainBaseNOSPD.h:59
ComputeStrainBaseNOSPD::_eigenstrains
std::vector< const MaterialProperty< RankTwoTensor > * > _eigenstrains
Definition: ComputeStrainBaseNOSPD.h:50
ComputeFiniteStrainNOSPD::_deformation_gradient_old
const MaterialProperty< RankTwoTensor > & _deformation_gradient_old
Material properties to fetch.
Definition: ComputeFiniteStrainNOSPD.h:52
ComputeFiniteStrainNOSPD::computeQpFhat
virtual void computeQpFhat()
Definition: ComputeFiniteStrainNOSPD.C:89
ComputeStrainBaseNOSPD::computeQpStrain
virtual void computeQpStrain()=0
Function to compute strain tensors.
ComputeFiniteStrainNOSPD::_strain_increment
MaterialProperty< RankTwoTensor > & _strain_increment
Definition: ComputeFiniteStrainNOSPD.h:47
ComputeStrainBaseNOSPD::computeBondStretch
virtual void computeBondStretch() override
Definition: ComputeStrainBaseNOSPD.C:156
RankTwoTensorTempl< Real >
ComputeStrainBaseNOSPD::_ddgraddu
MaterialProperty< RankTwoTensor > & _ddgraddu
Definition: ComputeStrainBaseNOSPD.h:57
ComputeStrainBaseNOSPD::_total_strain
MaterialProperty< RankTwoTensor > & _total_strain
Definition: ComputeStrainBaseNOSPD.h:61
ComputeFiniteStrainNOSPD::computeQpStrainRotationIncrements
void computeQpStrainRotationIncrements(RankTwoTensor &e, RankTwoTensor &r)
Function to compute strain and rotational increments.
Definition: ComputeFiniteStrainNOSPD.C:97
ComputeFiniteStrainNOSPD::_Fhat
std::vector< RankTwoTensor > _Fhat
'Incremental' deformation gradient
Definition: ComputeFiniteStrainNOSPD.h:59