https://mooseframework.inl.gov
Public Types | Public Member Functions | Static Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
MechanicsFiniteStrainBaseNOSPD Class Reference

Base kernel class for finite strain correspondence models. More...

#include <MechanicsFiniteStrainBaseNOSPD.h>

Inheritance diagram for MechanicsFiniteStrainBaseNOSPD:
[legend]

Public Types

typedef DerivativeMaterialPropertyNameInterface::SymbolName SymbolName
 

Public Member Functions

 MechanicsFiniteStrainBaseNOSPD (const InputParameters &parameters)
 
virtual void computeOffDiagJacobian (unsigned int jvar) override
 
virtual void computeLocalOffDiagJacobian (unsigned int, unsigned int)
 Function to compute local contribution to the off-diagonal Jacobian at the current nodes. More...
 
virtual void computePDNonlocalOffDiagJacobian (unsigned int, unsigned int)
 Function to compute nonlocal contribution to the off-diagonal Jacobian at the current nodes. More...
 
virtual void initialSetup () override
 
virtual void prepare () override
 
const GenericMaterialProperty< U, is_ad > & getDefaultMaterialProperty (const std::string &name)
 
const GenericMaterialProperty< U, is_ad > & getDefaultMaterialPropertyByName (const std::string &name)
 
void validateDerivativeMaterialPropertyBase (const std::string &base)
 
const MaterialPropertyName derivativePropertyName (const MaterialPropertyName &base, const std::vector< SymbolName > &c) const
 
const MaterialPropertyName derivativePropertyNameFirst (const MaterialPropertyName &base, const SymbolName &c1) const
 
const MaterialPropertyName derivativePropertyNameSecond (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2) const
 
const MaterialPropertyName derivativePropertyNameThird (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2, const SymbolName &c3) const
 
GenericMaterialProperty< U, is_ad > & declarePropertyDerivative (const std::string &base, const std::vector< VariableName > &c)
 
GenericMaterialProperty< U, is_ad > & declarePropertyDerivative (const std::string &base, const std::vector< SymbolName > &c)
 
GenericMaterialProperty< U, is_ad > & declarePropertyDerivative (const std::string &base, const SymbolName &c1, const SymbolName &c2="", const SymbolName &c3="")
 
GenericMaterialProperty< U, is_ad > & declarePropertyDerivative (const std::string &base, const std::vector< VariableName > &c)
 
GenericMaterialProperty< U, is_ad > & declarePropertyDerivative (const std::string &base, const std::vector< SymbolName > &c)
 
GenericMaterialProperty< U, is_ad > & declarePropertyDerivative (const std::string &base, const SymbolName &c1, const SymbolName &c2="", const SymbolName &c3="")
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivative (const std::string &base, const std::vector< VariableName > &c)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivative (const std::string &base, const std::vector< SymbolName > &c)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivative (const std::string &base, const SymbolName &c1, const SymbolName &c2="", const SymbolName &c3="")
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivative (const std::string &base, const SymbolName &c1, unsigned int v2, unsigned int v3=libMesh::invalid_uint)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivative (const std::string &base, unsigned int v1, unsigned int v2=libMesh::invalid_uint, unsigned int v3=libMesh::invalid_uint)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivative (const std::string &base, const std::vector< VariableName > &c)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivative (const std::string &base, const std::vector< SymbolName > &c)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivative (const std::string &base, const SymbolName &c1, const SymbolName &c2="", const SymbolName &c3="")
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivative (const std::string &base, const SymbolName &c1, unsigned int v2, unsigned int v3=libMesh::invalid_uint)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivative (const std::string &base, unsigned int v1, unsigned int v2=libMesh::invalid_uint, unsigned int v3=libMesh::invalid_uint)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivativeByName (const MaterialPropertyName &base, const std::vector< VariableName > &c)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivativeByName (const MaterialPropertyName &base, const std::vector< SymbolName > &c)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivativeByName (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2="", const SymbolName &c3="")
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivativeByName (const MaterialPropertyName &base, const std::vector< VariableName > &c)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivativeByName (const MaterialPropertyName &base, const std::vector< SymbolName > &c)
 
const GenericMaterialProperty< U, is_ad > & getMaterialPropertyDerivativeByName (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2="", const SymbolName &c3="")
 
void validateCoupling (const MaterialPropertyName &base, const std::vector< VariableName > &c, bool validate_aux=true)
 
void validateCoupling (const MaterialPropertyName &base, const VariableName &c1="", const VariableName &c2="", const VariableName &c3="")
 
void validateCoupling (const MaterialPropertyName &base, const std::vector< VariableName > &c, bool validate_aux=true)
 
void validateCoupling (const MaterialPropertyName &base, const VariableName &c1="", const VariableName &c2="", const VariableName &c3="")
 
void validateNonlinearCoupling (const MaterialPropertyName &base, const VariableName &c1="", const VariableName &c2="", const VariableName &c3="")
 
void validateNonlinearCoupling (const MaterialPropertyName &base, const VariableName &c1="", const VariableName &c2="", const VariableName &c3="")
 
const MaterialPropertyName propertyName (const MaterialPropertyName &base, const std::vector< SymbolName > &c) const
 
const MaterialPropertyName propertyName (const MaterialPropertyName &base, const std::vector< SymbolName > &c) const
 
const MaterialPropertyName propertyNameFirst (const MaterialPropertyName &base, const SymbolName &c1) const
 
const MaterialPropertyName propertyNameFirst (const MaterialPropertyName &base, const SymbolName &c1) const
 
const MaterialPropertyName propertyNameSecond (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2) const
 
const MaterialPropertyName propertyNameSecond (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2) const
 
const MaterialPropertyName propertyNameThird (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2, const SymbolName &c3) const
 
const MaterialPropertyName propertyNameThird (const MaterialPropertyName &base, const SymbolName &c1, const SymbolName &c2, const SymbolName &c3) const
 

Static Public Member Functions

static InputParameters validParams ()
 

Protected Member Functions

virtual RankTwoTensor computeDSDU (unsigned int component, unsigned int nd) override
 Function to compute derivative of stress with respect to displacements for small strain problems. More...
 
RankFourTensor computeDSDFhat (unsigned int nd)
 Function to compute derivative of stress with respect to derived deformation gradient. More...
 
Real computeDJDU (unsigned int component, unsigned int nd)
 Function to compute derivative of determinant of deformation gradient with respect to displacements. More...
 
RankTwoTensor computeDinvFTDU (unsigned int component, unsigned int nd)
 Function to compute derivative of deformation gradient inverse with respect to displacements. More...
 

Protected Attributes

std::vector< MooseVariable * > _disp_var
 displacement variables More...
 
unsigned int _ndisp
 number of displacement components More...
 
const std::vector< RealGradient > * _orientation
 Vector of bond in current configuration. More...
 
std::vector< dof_id_type_ivardofs
 Current variable dof numbers for nodes i and j. More...
 
std::vector< Real_weights
 weights used for the current element to obtain the nodal stress More...
 
RealGradient _current_vec
 Vector of bond in current configuration. More...
 
RealGradient _current_unit_vec
 Unit vector of bond in current configuration. More...
 
const MaterialProperty< RankTwoTensor > & _dgrad_old
 Material point based material property. More...
 
const MaterialProperty< RankTwoTensor > & _E_inc
 
const MaterialProperty< RankTwoTensor > & _R_inc
 
const MaterialProperty< Real > & _multi
 Material point based material properties. More...
 
const MaterialProperty< RankTwoTensor > & _stress
 
const MaterialProperty< RankTwoTensor > & _shape2
 
const MaterialProperty< RankTwoTensor > & _dgrad
 
const MaterialProperty< RankTwoTensor > & _ddgraddu
 
const MaterialProperty< RankTwoTensor > & _ddgraddv
 
const MaterialProperty< RankTwoTensor > & _ddgraddw
 
const MaterialProperty< RankFourTensor > & _Jacobian_mult
 
const std::vector< MaterialPropertyName > _eigenstrain_names
 
std::vector< const MaterialProperty< RankTwoTensor > * > _deigenstrain_dT
 
const bool _temp_coupled
 Temperature variable. More...
 
MooseVariable_temp_var
 
const bool _out_of_plane_strain_coupled
 Parameters for out-of-plane strain in weak plane stress formulation. More...
 
MooseVariable_out_of_plane_strain_var
 

Detailed Description

Base kernel class for finite strain correspondence models.

Definition at line 17 of file MechanicsFiniteStrainBaseNOSPD.h.

Constructor & Destructor Documentation

◆ MechanicsFiniteStrainBaseNOSPD()

MechanicsFiniteStrainBaseNOSPD::MechanicsFiniteStrainBaseNOSPD ( const InputParameters parameters)

Definition at line 27 of file MechanicsFiniteStrainBaseNOSPD.C.

28  : MechanicsBaseNOSPD(parameters),
29  _dgrad_old(getMaterialPropertyOld<RankTwoTensor>("deformation_gradient")),
30  _E_inc(getMaterialProperty<RankTwoTensor>("strain_increment")),
31  _R_inc(getMaterialProperty<RankTwoTensor>("rotation_increment"))
32 {
33 }
const MaterialProperty< RankTwoTensor > & _dgrad_old
Material point based material property.
const MaterialProperty< RankTwoTensor > & _E_inc
MechanicsBaseNOSPD(const InputParameters &parameters)
const MaterialProperty< RankTwoTensor > & _R_inc

Member Function Documentation

◆ computeDinvFTDU()

RankTwoTensor MechanicsFiniteStrainBaseNOSPD::computeDinvFTDU ( unsigned int  component,
unsigned int  nd 
)
protected

Function to compute derivative of deformation gradient inverse with respect to displacements.

Parameters
componentThe index of displacement component
ndThe local index of element node (either 1 or 2 for Edge2 element)
Returns
The calculated derivative

Definition at line 134 of file MechanicsFiniteStrainBaseNOSPD.C.

Referenced by HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeLocalJacobian(), HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeLocalJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), and HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeLocalOffDiagJacobian().

135 {
136  // for finite formulation, compute the derivative of transpose of inverse of deformation gradient
137  // w.r.t the solution components
138  // d(inv(F)_ji)/du = d(inv(F)_ji)/dF_kl * dF_kl/du = - inv(F)_jk * inv(F)_li * dF_kl/du
139  // the bases are gi, gj, gk, gl, indictates the inverse transpose rather than the inverse
140 
141  RankTwoTensor dinvFTdU;
142  dinvFTdU.zero();
143  RankTwoTensor invF = _dgrad[nd].inverse();
144  if (component == 0)
145  {
146  dinvFTdU(0, 1) =
147  _ddgraddu[nd](0, 2) * _dgrad[nd](2, 1) - _ddgraddu[nd](0, 1) * _dgrad[nd](2, 2);
148  dinvFTdU(0, 2) =
149  _ddgraddu[nd](0, 1) * _dgrad[nd](1, 2) - _ddgraddu[nd](0, 2) * _dgrad[nd](1, 1);
150  dinvFTdU(1, 1) =
151  _ddgraddu[nd](0, 0) * _dgrad[nd](2, 2) - _ddgraddu[nd](0, 2) * _dgrad[nd](2, 0);
152  dinvFTdU(1, 2) =
153  _ddgraddu[nd](0, 2) * _dgrad[nd](1, 0) - _ddgraddu[nd](0, 0) * _dgrad[nd](1, 2);
154  dinvFTdU(2, 1) =
155  _ddgraddu[nd](0, 1) * _dgrad[nd](2, 0) - _ddgraddu[nd](0, 0) * _dgrad[nd](2, 1);
156  dinvFTdU(2, 2) =
157  _ddgraddu[nd](0, 0) * _dgrad[nd](1, 1) - _ddgraddu[nd](0, 1) * _dgrad[nd](1, 0);
158  }
159  else if (component == 1)
160  {
161  dinvFTdU(0, 0) =
162  _ddgraddv[nd](1, 1) * _dgrad[nd](2, 2) - _ddgraddv[nd](1, 2) * _dgrad[nd](2, 1);
163  dinvFTdU(0, 2) =
164  _ddgraddv[nd](1, 2) * _dgrad[nd](0, 1) - _ddgraddv[nd](0, 2) * _dgrad[nd](1, 1);
165  dinvFTdU(1, 0) =
166  _ddgraddv[nd](1, 2) * _dgrad[nd](2, 0) - _ddgraddv[nd](1, 0) * _dgrad[nd](2, 2);
167  dinvFTdU(1, 2) =
168  _ddgraddv[nd](1, 0) * _dgrad[nd](0, 2) - _ddgraddv[nd](1, 2) * _dgrad[nd](0, 0);
169  dinvFTdU(2, 0) =
170  _ddgraddv[nd](1, 0) * _dgrad[nd](2, 1) - _ddgraddv[nd](1, 1) * _dgrad[nd](2, 0);
171  dinvFTdU(2, 2) =
172  _ddgraddv[nd](1, 1) * _dgrad[nd](0, 0) - _ddgraddv[nd](1, 0) * _dgrad[nd](0, 1);
173  }
174  else if (component == 2)
175  {
176  dinvFTdU(0, 0) =
177  _ddgraddw[nd](2, 2) * _dgrad[nd](1, 1) - _ddgraddw[nd](2, 1) * _dgrad[nd](1, 2);
178  dinvFTdU(0, 1) =
179  _ddgraddw[nd](2, 1) * _dgrad[nd](0, 2) - _ddgraddw[nd](2, 2) * _dgrad[nd](0, 1);
180  dinvFTdU(1, 0) =
181  _ddgraddw[nd](2, 0) * _dgrad[nd](1, 2) - _ddgraddw[nd](2, 2) * _dgrad[nd](1, 0);
182  dinvFTdU(1, 1) =
183  _ddgraddw[nd](2, 2) * _dgrad[nd](0, 0) - _ddgraddw[nd](2, 0) * _dgrad[nd](0, 2);
184  dinvFTdU(2, 0) =
185  _ddgraddw[nd](2, 1) * _dgrad[nd](1, 0) - _ddgraddw[nd](2, 0) * _dgrad[nd](1, 1);
186  dinvFTdU(2, 1) =
187  _ddgraddw[nd](2, 0) * _dgrad[nd](0, 1) - _ddgraddw[nd](2, 1) * _dgrad[nd](0, 0);
188  }
189 
190  dinvFTdU /= _dgrad[nd].det();
191  for (unsigned int i = 0; i < 3; ++i)
192  for (unsigned int j = 0; j < 3; ++j)
193  for (unsigned int k = 0; k < 3; ++k)
194  for (unsigned int l = 0; l < 3; ++l)
195  {
196  if (component == 0)
197  dinvFTdU(i, j) -= invF(i, j) * invF(l, k) * _ddgraddu[nd](k, l);
198  else if (component == 1)
199  dinvFTdU(i, j) -= invF(i, j) * invF(l, k) * _ddgraddv[nd](k, l);
200  else if (component == 2)
201  dinvFTdU(i, j) -= invF(i, j) * invF(l, k) * _ddgraddw[nd](k, l);
202  }
203 
204  return dinvFTdU;
205 }
const MaterialProperty< RankTwoTensor > & _ddgraddu
static const std::string component
Definition: NS.h:153
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
const MaterialProperty< RankTwoTensor > & _ddgraddw
const MaterialProperty< RankTwoTensor > & _ddgraddv
static const std::string k
Definition: NS.h:130
const MaterialProperty< RankTwoTensor > & _dgrad

◆ computeDJDU()

Real MechanicsFiniteStrainBaseNOSPD::computeDJDU ( unsigned int  component,
unsigned int  nd 
)
protected

Function to compute derivative of determinant of deformation gradient with respect to displacements.

Parameters
componentThe index of displacement component
ndThe local index of element node (either 1 or 2 for Edge2 element)
Returns
The calculated derivative

Definition at line 110 of file MechanicsFiniteStrainBaseNOSPD.C.

Referenced by HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeLocalJacobian(), HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeLocalJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), and HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeLocalOffDiagJacobian().

111 {
112  // for finite formulation, compute the derivative of determinant of deformation gradient w.r.t the
113  // solution components
114  // dJ / du = dJ / dF_ij * dF_ij / du = J * inv(F)_ji * dF_ij / du
115 
116  Real dJdU = 0.0;
117  RankTwoTensor invF = _dgrad[nd].inverse();
118  Real detF = _dgrad[nd].det();
119  for (unsigned int i = 0; i < 3; ++i)
120  for (unsigned int j = 0; j < 3; ++j)
121  {
122  if (component == 0)
123  dJdU += detF * invF(j, i) * _ddgraddu[nd](i, j);
124  else if (component == 1)
125  dJdU += detF * invF(j, i) * _ddgraddv[nd](i, j);
126  else if (component == 2)
127  dJdU += detF * invF(j, i) * _ddgraddw[nd](i, j);
128  }
129 
130  return dJdU;
131 }
const MaterialProperty< RankTwoTensor > & _ddgraddu
static const std::string component
Definition: NS.h:153
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
const MaterialProperty< RankTwoTensor > & _ddgraddw
const MaterialProperty< RankTwoTensor > & _ddgraddv
const MaterialProperty< RankTwoTensor > & _dgrad

◆ computeDSDFhat()

RankFourTensor MechanicsFiniteStrainBaseNOSPD::computeDSDFhat ( unsigned int  nd)
protected

Function to compute derivative of stress with respect to derived deformation gradient.

Parameters
ndThe local index of element node (either 1 or 2 for Edge2 element)
Returns
The calculated derivative

Definition at line 61 of file MechanicsFiniteStrainBaseNOSPD.C.

Referenced by computeDSDU(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeNonlocalJacobian(), HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeNonlocalJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computePDNonlocalOffDiagJacobian(), and HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computePDNonlocalOffDiagJacobian().

62 {
63  // compute the derivative of stress w.r.t the Fhat for finite strain
65  RankFourTensor dSdFhat;
66  dSdFhat.zero();
67 
68  // first calculate the derivative of incremental Cauchy stress w.r.t the inverse of Fhat
69  // Reference: M. M. Rashid (1993), Incremental Kinematics for finite element applications, IJNME
70  RankTwoTensor S_inc = _Jacobian_mult[nd] * _E_inc[nd];
71  RankFourTensor Tp1;
72  Tp1.zero();
73  for (unsigned int i = 0; i < 3; ++i)
74  for (unsigned int j = 0; j < 3; ++j)
75  for (unsigned int k = 0; k < 3; ++k)
76  for (unsigned int l = 0; l < 3; ++l)
77  for (unsigned int m = 0; m < 3; ++m)
78  for (unsigned int n = 0; n < 3; ++n)
79  for (unsigned int r = 0; r < 3; ++r)
80  Tp1(i, j, k, l) +=
81  S_inc(m, n) *
82  (_R_inc[nd](j, n) * (0.5 * I(k, m) * I(i, l) - I(m, l) * _R_inc[nd](i, k) +
83  0.5 * _R_inc[nd](i, k) * _R_inc[nd](m, l)) +
84  _R_inc[nd](i, m) * (0.5 * I(k, n) * I(j, l) - I(n, l) * _R_inc[nd](j, k) +
85  0.5 * _R_inc[nd](j, k) * _R_inc[nd](n, l))) -
86  _R_inc[nd](l, m) * _R_inc[nd](i, n) * _R_inc[nd](j, r) *
87  _Jacobian_mult[nd](n, r, m, k);
88 
89  // second calculate derivative of inverse of Fhat w.r.t Fhat
90  // d(inv(Fhat)_kl)/dFhat_mn = - inv(Fhat)_km * inv(Fhat)_nl
91  // the bases are gk, gl, gm, gn, indictates the inverse rather than the inverse transpose
92 
93  RankFourTensor Tp2;
94  Tp2.zero();
95  RankTwoTensor invFhat = (_dgrad[nd] * _dgrad_old[nd].inverse()).inverse();
96  for (unsigned int k = 0; k < 3; ++k)
97  for (unsigned int l = 0; l < 3; ++l)
98  for (unsigned int m = 0; m < 3; ++m)
99  for (unsigned int n = 0; n < 3; ++n)
100  Tp2(k, l, m, n) += -invFhat(k, m) * invFhat(n, l);
101 
102  // assemble two calculated quantities to form the derivative of Cauchy stress w.r.t
103  // Fhat
104  dSdFhat = Tp1 * Tp2;
105 
106  return dSdFhat;
107 }
const MaterialProperty< RankFourTensor > & _Jacobian_mult
void inverse(const std::vector< std::vector< Real >> &m, std::vector< std::vector< Real >> &m_inv)
const MaterialProperty< RankTwoTensor > & _dgrad_old
Material point based material property.
const MaterialProperty< RankTwoTensor > & _E_inc
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
static const std::string k
Definition: NS.h:130
const MaterialProperty< RankTwoTensor > & _R_inc
const MaterialProperty< RankTwoTensor > & _dgrad

◆ computeDSDU()

RankTwoTensor MechanicsFiniteStrainBaseNOSPD::computeDSDU ( unsigned int  component,
unsigned int  nd 
)
overrideprotectedvirtual

Function to compute derivative of stress with respect to displacements for small strain problems.

Parameters
componentThe index of displacement component
ndThe local index of element node (either 1 or 2 for Edge2 element)
Returns
The calculated derivative

Reimplemented from MechanicsBaseNOSPD.

Definition at line 36 of file MechanicsFiniteStrainBaseNOSPD.C.

Referenced by HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeLocalJacobian(), HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeLocalJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), and HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeLocalOffDiagJacobian().

37 {
38  // compute the derivative of stress w.r.t the solution components for finite strain
39  RankTwoTensor dSdU;
40 
41  // fetch the derivative of stress w.r.t the Fhat
42  RankFourTensor DSDFhat = computeDSDFhat(nd);
43 
44  // third calculate derivative of Fhat w.r.t solution components
45  RankTwoTensor Tp3;
46  if (component == 0)
47  Tp3 = _dgrad_old[nd].inverse() * _ddgraddu[nd];
48  else if (component == 1)
49  Tp3 = _dgrad_old[nd].inverse() * _ddgraddv[nd];
50  else if (component == 2)
51  Tp3 = _dgrad_old[nd].inverse() * _ddgraddw[nd];
52 
53  // assemble the fetched and calculated quantities to form the derivative of Cauchy stress w.r.t
54  // solution components
55  dSdU = DSDFhat * Tp3;
56 
57  return dSdU;
58 }
RankFourTensor computeDSDFhat(unsigned int nd)
Function to compute derivative of stress with respect to derived deformation gradient.
const MaterialProperty< RankTwoTensor > & _ddgraddu
static const std::string component
Definition: NS.h:153
const MaterialProperty< RankTwoTensor > & _dgrad_old
Material point based material property.
const MaterialProperty< RankTwoTensor > & _ddgraddw
const MaterialProperty< RankTwoTensor > & _ddgraddv

◆ computeLocalOffDiagJacobian()

virtual void MechanicsBasePD::computeLocalOffDiagJacobian ( unsigned int  ,
unsigned int   
)
inlinevirtualinherited

Function to compute local contribution to the off-diagonal Jacobian at the current nodes.

Parameters
coupled_componentThe coupled variable number

Reimplemented in HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD, HorizonStabilizedFormIISmallStrainMechanicsNOSPD, HorizonStabilizedFormIFiniteStrainMechanicsNOSPD, HorizonStabilizedFormISmallStrainMechanicsNOSPD, MechanicsOSPD, ForceStabilizedSmallStrainMechanicsNOSPD, MechanicsBPD, and WeakPlaneStressNOSPD.

Definition at line 31 of file MechanicsBasePD.h.

Referenced by MechanicsBasePD::computeOffDiagJacobian().

32  {};

◆ computeOffDiagJacobian()

void MechanicsBasePD::computeOffDiagJacobian ( unsigned int  jvar)
overridevirtualinherited

Definition at line 73 of file MechanicsBasePD.C.

74 {
75  prepare();
76 
77  if (jvar_num == _var.number())
78  computeJacobian();
79  else
80  {
81  unsigned int coupled_component = 0;
82  bool active = false;
83 
84  for (unsigned int i = 0; i < _dim; ++i)
85  if (jvar_num == _disp_var[i]->number())
86  {
87  coupled_component = i;
88  active = true;
89  }
90 
91  if (_temp_coupled && jvar_num == _temp_var->number())
92  active = true;
93 
95  active = true;
96 
97  if (active)
98  {
99  prepareMatrixTag(_assembly, _var.number(), jvar_num);
100  computeLocalOffDiagJacobian(jvar_num, coupled_component);
101  accumulateTaggedLocalMatrix();
102 
103  if (_use_full_jacobian)
104  computePDNonlocalOffDiagJacobian(jvar_num, coupled_component);
105  }
106  }
107 }
const bool _temp_coupled
Temperature variable.
unsigned int number() const
virtual void computeLocalOffDiagJacobian(unsigned int, unsigned int)
Function to compute local contribution to the off-diagonal Jacobian at the current nodes...
const bool _out_of_plane_strain_coupled
Parameters for out-of-plane strain in weak plane stress formulation.
std::vector< MooseVariable * > _disp_var
displacement variables
virtual void prepare() override
MooseVariable * _out_of_plane_strain_var
virtual void computePDNonlocalOffDiagJacobian(unsigned int, unsigned int)
Function to compute nonlocal contribution to the off-diagonal Jacobian at the current nodes...
MooseVariable * _temp_var

◆ computePDNonlocalOffDiagJacobian()

virtual void MechanicsBasePD::computePDNonlocalOffDiagJacobian ( unsigned int  ,
unsigned int   
)
inlinevirtualinherited

Function to compute nonlocal contribution to the off-diagonal Jacobian at the current nodes.

Parameters
jvar_numThe number of the first coupled variable
coupled_componentThe component number of the second coupled variable

Reimplemented in HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD, HorizonStabilizedFormIISmallStrainMechanicsNOSPD, HorizonStabilizedFormIFiniteStrainMechanicsNOSPD, HorizonStabilizedFormISmallStrainMechanicsNOSPD, MechanicsOSPD, ForceStabilizedSmallStrainMechanicsNOSPD, and WeakPlaneStressNOSPD.

Definition at line 39 of file MechanicsBasePD.h.

Referenced by MechanicsBasePD::computeOffDiagJacobian().

40  {};

◆ initialSetup()

void MechanicsBasePD::initialSetup ( )
overridevirtualinherited

Definition at line 45 of file MechanicsBasePD.C.

46 {
47  _orientation = &_assembly.getFE(FEType(), 1)->get_dxyzdxi();
48 }
const std::vector< RealGradient > * _orientation
Vector of bond in current configuration.

◆ prepare()

void MechanicsBasePD::prepare ( )
overridevirtualinherited

Definition at line 51 of file MechanicsBasePD.C.

Referenced by MechanicsBasePD::computeOffDiagJacobian(), GeneralizedPlaneStrainOffDiagNOSPD::computeOffDiagJacobianScalar(), and GeneralizedPlaneStrainOffDiagOSPD::computeOffDiagJacobianScalar().

52 {
54 
55  _ivardofs.resize(_nnodes);
56  _weights.resize(_nnodes);
57  for (unsigned int nd = 0; nd < _nnodes; ++nd)
58  {
59  _ivardofs[nd] = _current_elem->node_ptr(nd)->dof_number(_sys.number(), _var.number(), 0);
60  _weights[nd] = _pdmesh.getNeighborWeight(
61  _current_elem->node_id(nd),
62  _pdmesh.getNeighborIndex(_current_elem->node_id(nd), _current_elem->node_id(1 - nd)));
63  }
64 
65  for (unsigned int i = 0; i < _dim; ++i)
66  _current_vec(i) = _origin_vec(i) + _disp_var[i]->getNodalValue(*_current_elem->node_ptr(1)) -
67  _disp_var[i]->getNodalValue(*_current_elem->node_ptr(0));
68 
70 }
auto norm() const -> decltype(std::norm(Real()))
RealGradient _current_unit_vec
Unit vector of bond in current configuration.
std::vector< Real > _weights
weights used for the current element to obtain the nodal stress
std::vector< MooseVariable * > _disp_var
displacement variables
RealGradient _current_vec
Vector of bond in current configuration.
std::vector< dof_id_type > _ivardofs
Current variable dof numbers for nodes i and j.
virtual void prepare()
Function to precalculate data which will be used in the derived classes.

◆ validParams()

InputParameters MechanicsFiniteStrainBaseNOSPD::validParams ( )
static

Definition at line 13 of file MechanicsFiniteStrainBaseNOSPD.C.

Referenced by HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::validParams(), and HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::validParams().

14 {
16  params.addClassDescription("Base class for kernels of the stabilized non-ordinary "
17  "state-based peridynamic correspondence models");
18 
19  params.addParam<std::vector<MaterialPropertyName>>(
20  "eigenstrain_names",
21  {},
22  "List of eigenstrains to be coupled in non-ordinary state-based mechanics kernels");
23 
24  return params;
25 }
void addParam(const std::string &name, const std::initializer_list< typename T::value_type > &value, const std::string &doc_string)
static InputParameters validParams()
void addClassDescription(const std::string &doc_string)

Member Data Documentation

◆ _current_unit_vec

RealGradient MechanicsBasePD::_current_unit_vec
protectedinherited

◆ _current_vec

RealGradient MechanicsBasePD::_current_vec
protectedinherited

◆ _ddgraddu

const MaterialProperty<RankTwoTensor>& MechanicsBaseNOSPD::_ddgraddu
protectedinherited

◆ _ddgraddv

const MaterialProperty<RankTwoTensor>& MechanicsBaseNOSPD::_ddgraddv
protectedinherited

◆ _ddgraddw

const MaterialProperty<RankTwoTensor>& MechanicsBaseNOSPD::_ddgraddw
protectedinherited

◆ _deigenstrain_dT

std::vector<const MaterialProperty<RankTwoTensor> *> MechanicsBaseNOSPD::_deigenstrain_dT
protectedinherited

◆ _dgrad

const MaterialProperty<RankTwoTensor>& MechanicsBaseNOSPD::_dgrad
protectedinherited

◆ _dgrad_old

const MaterialProperty<RankTwoTensor>& MechanicsFiniteStrainBaseNOSPD::_dgrad_old
protected

◆ _disp_var

std::vector<MooseVariable *> MechanicsBasePD::_disp_var
protectedinherited

◆ _E_inc

const MaterialProperty<RankTwoTensor>& MechanicsFiniteStrainBaseNOSPD::_E_inc
protected

Definition at line 53 of file MechanicsFiniteStrainBaseNOSPD.h.

Referenced by computeDSDFhat().

◆ _eigenstrain_names

const std::vector<MaterialPropertyName> MechanicsBaseNOSPD::_eigenstrain_names
protectedinherited

Definition at line 43 of file MechanicsBaseNOSPD.h.

Referenced by MechanicsBaseNOSPD::MechanicsBaseNOSPD().

◆ _ivardofs

std::vector<dof_id_type> MechanicsBasePD::_ivardofs
protectedinherited

◆ _Jacobian_mult

const MaterialProperty<RankFourTensor>& MechanicsBaseNOSPD::_Jacobian_mult
protectedinherited

◆ _multi

const MaterialProperty<Real>& MechanicsBaseNOSPD::_multi
protectedinherited

Material point based material properties.

Definition at line 35 of file MechanicsBaseNOSPD.h.

Referenced by GeneralizedPlaneStrainOffDiagNOSPD::computeDispFullOffDiagJacobianScalar(), GeneralizedPlaneStrainOffDiagNOSPD::computeDispPartialOffDiagJacobianScalar(), ForceStabilizedSmallStrainMechanicsNOSPD::computeLocalJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeLocalJacobian(), HorizonStabilizedFormISmallStrainMechanicsNOSPD::computeLocalJacobian(), HorizonStabilizedFormIISmallStrainMechanicsNOSPD::computeLocalJacobian(), HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeLocalJacobian(), ForceStabilizedSmallStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), HorizonStabilizedFormISmallStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), HorizonStabilizedFormIISmallStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), ForceStabilizedSmallStrainMechanicsNOSPD::computeLocalResidual(), HorizonStabilizedFormISmallStrainMechanicsNOSPD::computeLocalResidual(), HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeLocalResidual(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeLocalResidual(), HorizonStabilizedFormIISmallStrainMechanicsNOSPD::computeLocalResidual(), ForceStabilizedSmallStrainMechanicsNOSPD::computeNonlocalJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeNonlocalJacobian(), HorizonStabilizedFormISmallStrainMechanicsNOSPD::computeNonlocalJacobian(), HorizonStabilizedFormIISmallStrainMechanicsNOSPD::computeNonlocalJacobian(), HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeNonlocalJacobian(), HorizonStabilizedFormIISmallStrainMechanicsNOSPD::computeNonlocalResidual(), HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeNonlocalResidual(), ForceStabilizedSmallStrainMechanicsNOSPD::computePDNonlocalOffDiagJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computePDNonlocalOffDiagJacobian(), HorizonStabilizedFormISmallStrainMechanicsNOSPD::computePDNonlocalOffDiagJacobian(), HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computePDNonlocalOffDiagJacobian(), and HorizonStabilizedFormIISmallStrainMechanicsNOSPD::computePDNonlocalOffDiagJacobian().

◆ _ndisp

unsigned int MechanicsBasePD::_ndisp
protectedinherited

number of displacement components

Definition at line 55 of file MechanicsBasePD.h.

Referenced by MechanicsBasePD::MechanicsBasePD().

◆ _orientation

const std::vector<RealGradient>* MechanicsBasePD::_orientation
protectedinherited

Vector of bond in current configuration.

Definition at line 63 of file MechanicsBasePD.h.

Referenced by MechanicsBasePD::initialSetup().

◆ _out_of_plane_strain_coupled

const bool MechanicsBasePD::_out_of_plane_strain_coupled
protectedinherited

◆ _out_of_plane_strain_var

MooseVariable* MechanicsBasePD::_out_of_plane_strain_var
protectedinherited

◆ _R_inc

const MaterialProperty<RankTwoTensor>& MechanicsFiniteStrainBaseNOSPD::_R_inc
protected

Definition at line 54 of file MechanicsFiniteStrainBaseNOSPD.h.

Referenced by computeDSDFhat().

◆ _shape2

const MaterialProperty<RankTwoTensor>& MechanicsBaseNOSPD::_shape2
protectedinherited

Definition at line 37 of file MechanicsBaseNOSPD.h.

Referenced by GeneralizedPlaneStrainOffDiagNOSPD::computeDispFullOffDiagJacobianScalar(), GeneralizedPlaneStrainOffDiagNOSPD::computeDispPartialOffDiagJacobianScalar(), ForceStabilizedSmallStrainMechanicsNOSPD::computeLocalJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeLocalJacobian(), HorizonStabilizedFormISmallStrainMechanicsNOSPD::computeLocalJacobian(), HorizonStabilizedFormIISmallStrainMechanicsNOSPD::computeLocalJacobian(), HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeLocalJacobian(), ForceStabilizedSmallStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), HorizonStabilizedFormISmallStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), HorizonStabilizedFormIISmallStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), ForceStabilizedSmallStrainMechanicsNOSPD::computeLocalResidual(), HorizonStabilizedFormISmallStrainMechanicsNOSPD::computeLocalResidual(), HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeLocalResidual(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeLocalResidual(), HorizonStabilizedFormIISmallStrainMechanicsNOSPD::computeLocalResidual(), ForceStabilizedSmallStrainMechanicsNOSPD::computeNonlocalJacobian(), HorizonStabilizedFormISmallStrainMechanicsNOSPD::computeNonlocalJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeNonlocalJacobian(), HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeNonlocalJacobian(), HorizonStabilizedFormIISmallStrainMechanicsNOSPD::computeNonlocalJacobian(), HorizonStabilizedFormIISmallStrainMechanicsNOSPD::computeNonlocalResidual(), HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeNonlocalResidual(), WeakPlaneStressNOSPD::computePDNonlocalOffDiagJacobian(), ForceStabilizedSmallStrainMechanicsNOSPD::computePDNonlocalOffDiagJacobian(), HorizonStabilizedFormISmallStrainMechanicsNOSPD::computePDNonlocalOffDiagJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computePDNonlocalOffDiagJacobian(), HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computePDNonlocalOffDiagJacobian(), and HorizonStabilizedFormIISmallStrainMechanicsNOSPD::computePDNonlocalOffDiagJacobian().

◆ _stress

const MaterialProperty<RankTwoTensor>& MechanicsBaseNOSPD::_stress
protectedinherited

◆ _temp_coupled

const bool MechanicsBasePD::_temp_coupled
protectedinherited

◆ _temp_var

MooseVariable* MechanicsBasePD::_temp_var
protectedinherited

◆ _weights

std::vector<Real> MechanicsBasePD::_weights
protectedinherited

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