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

Kernel class for Form II of the horizon-stabilized peridynamic correspondence model for finite strain. More...

#include <HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD.h>

Inheritance diagram for HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD:
[legend]

Public Types

typedef DerivativeMaterialPropertyNameInterface::SymbolName SymbolName
 

Public Member Functions

 HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD (const InputParameters &parameters)
 
virtual void computeOffDiagJacobian (unsigned int jvar) override
 
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 void computeLocalResidual () override
 
virtual void computeNonlocalResidual () override
 
virtual void computeLocalJacobian () override
 
virtual void computeNonlocalJacobian () override
 
virtual void computeLocalOffDiagJacobian (unsigned int jvar_num, unsigned int coupled_component) override
 Function to compute local contribution to the off-diagonal Jacobian at the current nodes. More...
 
virtual void computePDNonlocalOffDiagJacobian (unsigned int jvar_num, unsigned int coupled_component) override
 Function to compute nonlocal contribution to the off-diagonal Jacobian at the current nodes. More...
 
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

const unsigned int _component
 The index of displacement component. More...
 
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

Kernel class for Form II of the horizon-stabilized peridynamic correspondence model for finite strain.

Definition at line 20 of file HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD.h.

Constructor & Destructor Documentation

◆ HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD()

HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD ( const InputParameters parameters)

Definition at line 31 of file HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD.C.

32  : MechanicsFiniteStrainBaseNOSPD(parameters), _component(getParam<unsigned int>("component"))
33 {
34 }
MechanicsFiniteStrainBaseNOSPD(const InputParameters &parameters)
const unsigned int _component
The index of displacement component.

Member Function Documentation

◆ computeDinvFTDU()

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

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(), computeLocalJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), and 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 
)
protectedinherited

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(), computeLocalJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), and 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)
protectedinherited

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 MechanicsFiniteStrainBaseNOSPD::computeDSDU(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeNonlocalJacobian(), computeNonlocalJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computePDNonlocalOffDiagJacobian(), and 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 
)
overrideprotectedvirtualinherited

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(), computeLocalJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), and 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

◆ computeLocalJacobian()

void HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeLocalJacobian ( )
overrideprotectedvirtual

Definition at line 116 of file HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD.C.

117 {
118  for (unsigned int nd = 0; nd < _nnodes; ++nd)
119  {
120  std::vector<dof_id_type> ivardofs(_nnodes);
121  ivardofs[0] = _current_elem->node_ptr(nd)->dof_number(_sys.number(), _var.number(), 0);
122  std::vector<dof_id_type> neighbors = _pdmesh.getNeighbors(_current_elem->node_id(nd));
123  std::vector<dof_id_type> bonds = _pdmesh.getBonds(_current_elem->node_id(nd));
124 
125  dof_id_type nb_index =
126  std::find(neighbors.begin(), neighbors.end(), _current_elem->node_id(1 - nd)) -
127  neighbors.begin();
128  std::vector<dof_id_type> dg_neighbors =
129  _pdmesh.getBondDeformationGradientNeighbors(_current_elem->node_id(nd), nb_index);
130 
131  RankTwoTensor dPxdUx =
132  computeDJDU(_component, nd) * _stress[nd] * _dgrad[nd].inverse().transpose() +
133  _dgrad[nd].det() * computeDSDU(_component, nd) * _dgrad[nd].inverse().transpose() +
134  _dgrad[nd].det() * _stress[nd] * computeDinvFTDU(_component, nd);
135 
136  RealGradient origin_vec_nb;
137  Real node_vol_nb;
138 
139  for (unsigned int nb = 0; nb < dg_neighbors.size(); ++nb)
140  if (_bond_status_var->getElementalValue(_pdmesh.elemPtr(bonds[dg_neighbors[nb]])) > 0.5)
141  {
142  ivardofs[1] = _pdmesh.nodePtr(neighbors[dg_neighbors[nb]])
143  ->dof_number(_sys.number(), _var.number(), 0);
144  origin_vec_nb = _pdmesh.getNodeCoord(neighbors[dg_neighbors[nb]]) -
145  _pdmesh.getNodeCoord(_current_elem->node_id(nd));
146  node_vol_nb = _pdmesh.getNodeVolume(neighbors[dg_neighbors[nb]]);
147 
148  for (unsigned int i = 0; i < _nnodes; ++i)
149  for (unsigned int j = 0; j < _nnodes; ++j)
150  _local_ke(i, j) = (i == 0 ? -1 : 1) * (j == 0 ? 1 : 0) * _multi[nd] *
151  _horizon_radius[nd] / origin_vec_nb.norm() *
152  (dPxdUx * _shape2[nd].inverse()).row(_component) * origin_vec_nb *
153  node_vol_nb * _bond_status;
154 
155  addJacobian(_assembly, _local_ke, ivardofs, ivardofs, _var.scalingFactor());
156  }
157  }
158  _local_ke.zero();
159 }
auto norm() const -> decltype(std::norm(Real()))
RankTwoTensor computeDinvFTDU(unsigned int component, unsigned int nd)
Function to compute derivative of deformation gradient inverse with respect to displacements.
const MaterialProperty< RankTwoTensor > & _stress
virtual RankTwoTensor computeDSDU(unsigned int component, unsigned int nd) override
Function to compute derivative of stress with respect to displacements for small strain problems...
const unsigned int _component
The index of displacement component.
Real computeDJDU(unsigned int component, unsigned int nd)
Function to compute derivative of determinant of deformation gradient with respect to displacements...
const MaterialProperty< RankTwoTensor > & _shape2
const MaterialProperty< Real > & _multi
Material point based material properties.
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")
uint8_t dof_id_type
const MaterialProperty< RankTwoTensor > & _dgrad

◆ computeLocalOffDiagJacobian()

void HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeLocalOffDiagJacobian ( unsigned int  ,
unsigned int   
)
overrideprotectedvirtual

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

Parameters
coupled_componentThe coupled variable number

Reimplemented from MechanicsBasePD.

Definition at line 273 of file HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD.C.

275 {
276  if (_temp_coupled && jvar_num == _temp_var->number()) // temperature is coupled
277  {
278  std::vector<RankTwoTensor> dSdT(_nnodes);
279  for (unsigned int nd = 0; nd < _nnodes; ++nd)
280  for (unsigned int es = 0; es < _deigenstrain_dT.size(); ++es)
281  dSdT[nd] = -_dgrad[nd].det() * _Jacobian_mult[nd] * (*_deigenstrain_dT[es])[nd] *
282  _dgrad[nd].inverse().transpose();
283 
284  for (unsigned int nd = 0; nd < _nnodes; ++nd)
285  {
286  std::vector<dof_id_type> ivardofs(_nnodes), jvardofs(_nnodes);
287  ivardofs[0] = _current_elem->node_ptr(nd)->dof_number(_sys.number(), _var.number(), 0);
288  jvardofs[0] = _current_elem->node_ptr(nd)->dof_number(_sys.number(), jvar_num, 0);
289  std::vector<dof_id_type> neighbors = _pdmesh.getNeighbors(_current_elem->node_id(nd));
290  std::vector<dof_id_type> bonds = _pdmesh.getBonds(_current_elem->node_id(nd));
291 
292  dof_id_type nb_index =
293  std::find(neighbors.begin(), neighbors.end(), _current_elem->node_id(1 - nd)) -
294  neighbors.begin();
295  std::vector<dof_id_type> dg_neighbors =
296  _pdmesh.getBondDeformationGradientNeighbors(_current_elem->node_id(nd), nb_index);
297 
298  RealGradient origin_vec_nb;
299  Real node_vol_nb;
300 
301  for (unsigned int nb = 0; nb < dg_neighbors.size(); ++nb)
302  if (_bond_status_var->getElementalValue(_pdmesh.elemPtr(bonds[dg_neighbors[nb]])) > 0.5)
303  {
304  ivardofs[1] = _pdmesh.nodePtr(neighbors[dg_neighbors[nb]])
305  ->dof_number(_sys.number(), _var.number(), 0);
306  jvardofs[1] =
307  _pdmesh.nodePtr(neighbors[dg_neighbors[nb]])->dof_number(_sys.number(), jvar_num, 0);
308  origin_vec_nb = _pdmesh.getNodeCoord(neighbors[dg_neighbors[nb]]) -
309  _pdmesh.getNodeCoord(_current_elem->node_id(nd));
310  node_vol_nb = _pdmesh.getNodeVolume(neighbors[dg_neighbors[nb]]);
311 
312  for (unsigned int i = 0; i < _nnodes; ++i)
313  for (unsigned int j = 0; j < _nnodes; ++j)
314  _local_ke(i, j) = (i == 0 ? -1 : 1) * (j == 0 ? 1 : 0) * _multi[nd] *
315  _horizon_radius[nd] / origin_vec_nb.norm() *
316  (dSdT[nd] * _shape2[nd].inverse()).row(_component) * origin_vec_nb *
317  node_vol_nb * _bond_status;
318 
319  addJacobian(_assembly, _local_ke, ivardofs, jvardofs, _var.scalingFactor());
320  }
321  }
322  _local_ke.zero();
323  }
324  else if (_out_of_plane_strain_coupled &&
325  jvar_num == _out_of_plane_strain_var
326  ->number()) // weak plane stress case, out_of_plane_strain is coupled
327  {
328  std::vector<RankTwoTensor> dSdE33(_nnodes);
329  for (unsigned int nd = 0; nd < _nnodes; ++nd)
330  {
331  for (unsigned int i = 0; i < 3; ++i)
332  for (unsigned int j = 0; j < 3; ++j)
333  dSdE33[nd](i, j) = _Jacobian_mult[nd](i, j, 2, 2);
334 
335  dSdE33[nd] = _dgrad[nd].det() * dSdE33[nd] * _dgrad[nd].inverse().transpose();
336  }
337 
338  for (unsigned int nd = 0; nd < _nnodes; ++nd)
339  {
340  std::vector<dof_id_type> ivardofs(_nnodes), jvardofs(_nnodes);
341  ivardofs[0] = _current_elem->node_ptr(nd)->dof_number(_sys.number(), _var.number(), 0);
342  jvardofs[0] = _current_elem->node_ptr(nd)->dof_number(_sys.number(), jvar_num, 0);
343  std::vector<dof_id_type> neighbors = _pdmesh.getNeighbors(_current_elem->node_id(nd));
344  std::vector<dof_id_type> bonds = _pdmesh.getBonds(_current_elem->node_id(nd));
345 
346  dof_id_type nb_index =
347  std::find(neighbors.begin(), neighbors.end(), _current_elem->node_id(1 - nd)) -
348  neighbors.begin();
349  std::vector<dof_id_type> dg_neighbors =
350  _pdmesh.getBondDeformationGradientNeighbors(_current_elem->node_id(nd), nb_index);
351 
352  RealGradient origin_vec_nb;
353  Real node_vol_nb;
354 
355  for (unsigned int nb = 0; nb < dg_neighbors.size(); ++nb)
356  if (_bond_status_var->getElementalValue(_pdmesh.elemPtr(bonds[dg_neighbors[nb]])) > 0.5)
357  {
358  ivardofs[1] = _pdmesh.nodePtr(neighbors[dg_neighbors[nb]])
359  ->dof_number(_sys.number(), _var.number(), 0);
360  jvardofs[1] =
361  _pdmesh.nodePtr(neighbors[dg_neighbors[nb]])->dof_number(_sys.number(), jvar_num, 0);
362  origin_vec_nb = _pdmesh.getNodeCoord(neighbors[dg_neighbors[nb]]) -
363  _pdmesh.getNodeCoord(_current_elem->node_id(nd));
364  node_vol_nb = _pdmesh.getNodeVolume(neighbors[dg_neighbors[nb]]);
365 
366  for (unsigned int i = 0; i < _nnodes; ++i)
367  for (unsigned int j = 0; j < _nnodes; ++j)
368  _local_ke(i, j) = (i == 0 ? -1 : 1) * (j == 0 ? 1 : 0) * _multi[nd] *
369  _horizon_radius[nd] / origin_vec_nb.norm() *
370  (dSdE33[nd] * _shape2[nd].inverse()).row(_component) *
371  origin_vec_nb * node_vol_nb * _bond_status;
372 
373  addJacobian(_assembly, _local_ke, ivardofs, jvardofs, _var.scalingFactor());
374  }
375  }
376  _local_ke.zero();
377  }
378  else // off-diagonal Jacobian with respect to other displacement variables
379  {
380  for (unsigned int nd = 0; nd < _nnodes; ++nd)
381  {
382  std::vector<dof_id_type> ivardofs(_nnodes), jvardofs(_nnodes);
383  ivardofs[0] = _current_elem->node_ptr(nd)->dof_number(_sys.number(), _var.number(), 0);
384  jvardofs[0] = _current_elem->node_ptr(nd)->dof_number(_sys.number(), jvar_num, 0);
385  std::vector<dof_id_type> neighbors = _pdmesh.getNeighbors(_current_elem->node_id(nd));
386  std::vector<dof_id_type> bonds = _pdmesh.getBonds(_current_elem->node_id(nd));
387 
388  dof_id_type nb_index =
389  std::find(neighbors.begin(), neighbors.end(), _current_elem->node_id(1 - nd)) -
390  neighbors.begin();
391  std::vector<dof_id_type> dg_neighbors =
392  _pdmesh.getBondDeformationGradientNeighbors(_current_elem->node_id(nd), nb_index);
393 
394  RankTwoTensor dPxdUy =
395  computeDJDU(coupled_component, nd) * _stress[nd] * _dgrad[nd].inverse().transpose() +
396  _dgrad[nd].det() * computeDSDU(coupled_component, nd) * _dgrad[nd].inverse().transpose() +
397  _dgrad[nd].det() * _stress[nd] * computeDinvFTDU(coupled_component, nd);
398 
399  RealGradient origin_vec_nb;
400  Real node_vol_nb;
401 
402  for (unsigned int nb = 0; nb < dg_neighbors.size(); ++nb)
403  if (_bond_status_var->getElementalValue(_pdmesh.elemPtr(bonds[dg_neighbors[nb]])) > 0.5)
404  {
405  ivardofs[1] = _pdmesh.nodePtr(neighbors[dg_neighbors[nb]])
406  ->dof_number(_sys.number(), _var.number(), 0);
407  jvardofs[1] =
408  _pdmesh.nodePtr(neighbors[dg_neighbors[nb]])->dof_number(_sys.number(), jvar_num, 0);
409  origin_vec_nb = _pdmesh.getNodeCoord(neighbors[dg_neighbors[nb]]) -
410  _pdmesh.getNodeCoord(_current_elem->node_id(nd));
411  node_vol_nb = _pdmesh.getNodeVolume(neighbors[dg_neighbors[nb]]);
412 
413  for (unsigned int i = 0; i < _nnodes; ++i)
414  for (unsigned int j = 0; j < _nnodes; ++j)
415  _local_ke(i, j) = (i == 0 ? -1 : 1) * (j == 0 ? 1 : 0) * _multi[nd] *
416  _horizon_radius[nd] / origin_vec_nb.norm() *
417  (dPxdUy * _shape2[nd].inverse()).row(_component) * origin_vec_nb *
418  node_vol_nb * _bond_status;
419 
420  addJacobian(_assembly, _local_ke, ivardofs, jvardofs, _var.scalingFactor());
421  }
422  }
423  _local_ke.zero();
424  }
425 }
const bool _temp_coupled
Temperature variable.
auto norm() const -> decltype(std::norm(Real()))
RankTwoTensor computeDinvFTDU(unsigned int component, unsigned int nd)
Function to compute derivative of deformation gradient inverse with respect to displacements.
const MaterialProperty< RankFourTensor > & _Jacobian_mult
unsigned int number() const
const MaterialProperty< RankTwoTensor > & _stress
void inverse(const std::vector< std::vector< Real >> &m, std::vector< std::vector< Real >> &m_inv)
const bool _out_of_plane_strain_coupled
Parameters for out-of-plane strain in weak plane stress formulation.
std::vector< const MaterialProperty< RankTwoTensor > * > _deigenstrain_dT
virtual RankTwoTensor computeDSDU(unsigned int component, unsigned int nd) override
Function to compute derivative of stress with respect to displacements for small strain problems...
const unsigned int _component
The index of displacement component.
Real computeDJDU(unsigned int component, unsigned int nd)
Function to compute derivative of determinant of deformation gradient with respect to displacements...
MooseVariable * _out_of_plane_strain_var
const MaterialProperty< RankTwoTensor > & _shape2
const MaterialProperty< Real > & _multi
Material point based material properties.
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")
MooseVariable * _temp_var
uint8_t dof_id_type
const MaterialProperty< RankTwoTensor > & _dgrad

◆ computeLocalResidual()

void HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeLocalResidual ( )
overrideprotectedvirtual

Definition at line 37 of file HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD.C.

38 {
39  // For finite strain formulation, the _stress tensor gotten from material class is the
40  // Cauchy stress (Sigma). the first Piola-Kirchhoff stress (P) is then obtained as
41  // P = J * Sigma * inv(F)^T.
42  // Nodal force states are based on the first Piola-Kirchhoff stress tensors (P).
43  // i.e., T = (J * Sigma * inv(F)^T) * inv(Shape) * xi * multi.
44  // Cauchy stress is calculated as Sigma_n+1 = Sigma_n + R * (C * dt * D) * R^T
45 
46  for (unsigned int nd = 0; nd < _nnodes; ++nd)
47  for (unsigned int i = 0; i < _nnodes; ++i)
48  _local_re(i) += (i == 0 ? -1 : 1) * _multi[nd] * _horizon_radius[nd] / _origin_vec.norm() *
49  ((_dgrad[nd].det() * _stress[nd] * _dgrad[nd].inverse().transpose()) *
50  _shape2[nd].inverse())
51  .row(_component) *
52  _origin_vec * _node_vol[1 - nd] * _bond_status;
53 }
const MaterialProperty< RankTwoTensor > & _stress
void inverse(const std::vector< std::vector< Real >> &m, std::vector< std::vector< Real >> &m_inv)
const unsigned int _component
The index of displacement component.
const MaterialProperty< RankTwoTensor > & _shape2
const MaterialProperty< Real > & _multi
Material point based material properties.
const MaterialProperty< RankTwoTensor > & _dgrad

◆ computeNonlocalJacobian()

void HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeNonlocalJacobian ( )
overrideprotectedvirtual

Definition at line 162 of file HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD.C.

163 {
164  // includes dTi/dUj and dTj/dUi contributions
165  for (unsigned int nd = 0; nd < _nnodes; ++nd)
166  {
167  RankFourTensor dSdFhat = computeDSDFhat(nd);
168  RankTwoTensor invF = _dgrad[nd].inverse();
169  Real detF = _dgrad[nd].det();
170  // calculation of jacobian contribution to current_node's neighbors
171  std::vector<dof_id_type> ivardofs(_nnodes), jvardofs(_nnodes);
172  ivardofs[0] = _current_elem->node_ptr(nd)->dof_number(_sys.number(), _var.number(), 0);
173  jvardofs[0] = _current_elem->node_ptr(nd)->dof_number(_sys.number(), _var.number(), 0);
174  std::vector<dof_id_type> neighbors = _pdmesh.getNeighbors(_current_elem->node_id(nd));
175  std::vector<dof_id_type> bonds = _pdmesh.getBonds(_current_elem->node_id(nd));
176 
177  dof_id_type nb_index =
178  std::find(neighbors.begin(), neighbors.end(), _current_elem->node_id(1 - nd)) -
179  neighbors.begin();
180  std::vector<dof_id_type> dg_neighbors =
181  _pdmesh.getBondDeformationGradientNeighbors(_current_elem->node_id(nd), nb_index);
182 
183  RealGradient origin_vec_nb1;
184  Real node_vol_nb1;
185 
186  for (unsigned int nb1 = 0; nb1 < dg_neighbors.size(); ++nb1)
187  if (_bond_status_var->getElementalValue(_pdmesh.elemPtr(bonds[dg_neighbors[nb1]])) > 0.5)
188  {
189  ivardofs[1] = _pdmesh.nodePtr(neighbors[dg_neighbors[nb1]])
190  ->dof_number(_sys.number(), _var.number(), 0);
191  origin_vec_nb1 = _pdmesh.getNodeCoord(neighbors[dg_neighbors[nb1]]) -
192  _pdmesh.getNodeCoord(_current_elem->node_id(nd));
193  node_vol_nb1 = _pdmesh.getNodeVolume(neighbors[dg_neighbors[nb1]]);
194 
195  Real vol_nb2, dJdU;
196  RealGradient origin_vec_nb2;
197  RankTwoTensor dFdUk, dPxdUkx, dSdU, dinvFTdU;
198 
199  for (unsigned int nb2 = 0; nb2 < dg_neighbors.size(); ++nb2)
200  if (_bond_status_var->getElementalValue(_pdmesh.elemPtr(bonds[dg_neighbors[nb2]])) > 0.5)
201  {
202  ivardofs[1] = _pdmesh.nodePtr(neighbors[dg_neighbors[nb2]])
203  ->dof_number(_sys.number(), _var.number(), 0);
204  vol_nb2 = _pdmesh.getNodeVolume(neighbors[dg_neighbors[nb2]]);
205 
206  origin_vec_nb2 = _pdmesh.getNodeCoord(neighbors[dg_neighbors[nb2]]) -
207  _pdmesh.getNodeCoord(_current_elem->node_id(nd));
208 
209  dFdUk.zero();
210  for (unsigned int i = 0; i < _dim; ++i)
211  dFdUk(_component, i) =
212  _horizon_radius[nd] / origin_vec_nb2.norm() * origin_vec_nb2(i) * vol_nb2;
213 
214  dFdUk *= _shape2[nd].inverse();
215 
216  // calculate dJ/du
217  dJdU = 0.0;
218  for (unsigned int i = 0; i < 3; ++i)
219  for (unsigned int j = 0; j < 3; ++j)
220  dJdU += detF * invF(j, i) * dFdUk(i, j);
221 
222  // calculate dS/du
223  dSdU = dSdFhat * dFdUk * _dgrad_old[nd].inverse();
224 
225  // calculate dinv(F)Tdu
226  dinvFTdU.zero();
227  for (unsigned int i = 0; i < 3; ++i)
228  for (unsigned int J = 0; J < 3; ++J)
229  for (unsigned int k = 0; k < 3; ++k)
230  for (unsigned int L = 0; L < 3; ++L)
231  dinvFTdU(i, J) += -invF(J, k) * invF(L, i) * dFdUk(k, L);
232 
233  // calculate the derivative of first Piola-Kirchhoff stress w.r.t displacements
234  dPxdUkx = dJdU * _stress[nd] * invF.transpose() + detF * dSdU * invF.transpose() +
235  detF * _stress[nd] * dinvFTdU;
236 
237  for (unsigned int i = 0; i < _nnodes; ++i)
238  for (unsigned int j = 0; j < _nnodes; ++j)
239  _local_ke(i, j) = (i == 0 ? -1 : 1) * (j == 0 ? 0 : 1) * _multi[nd] *
240  _horizon_radius[nd] / origin_vec_nb1.norm() *
241  (dPxdUkx * _shape2[nd].inverse()).row(_component) *
242  origin_vec_nb1 * node_vol_nb1 * _bond_status;
243 
244  addJacobian(_assembly, _local_ke, ivardofs, jvardofs, _var.scalingFactor());
245 
246  if (_has_diag_save_in)
247  {
248  unsigned int rows = _nnodes;
249  DenseVector<Real> diag(rows);
250  for (unsigned int i = 0; i < rows; ++i)
251  diag(i) = _local_ke(i, i);
252 
253  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
254  for (unsigned int i = 0; i < _diag_save_in.size(); ++i)
255  {
256  std::vector<dof_id_type> diag_save_in_dofs(2);
257  diag_save_in_dofs[0] = _current_elem->node_ptr(nd)->dof_number(
258  _diag_save_in[i]->sys().number(), _diag_save_in[i]->number(), 0);
259  diag_save_in_dofs[1] = _pdmesh.nodePtr(neighbors[dg_neighbors[nb2]])
260  ->dof_number(_diag_save_in[i]->sys().number(),
261  _diag_save_in[i]->number(),
262  0);
263 
264  _diag_save_in[i]->sys().solution().add_vector(diag, diag_save_in_dofs);
265  }
266  }
267  }
268  }
269  }
270 }
RankFourTensor computeDSDFhat(unsigned int nd)
Function to compute derivative of stress with respect to derived deformation gradient.
auto norm() const -> decltype(std::norm(Real()))
const MaterialProperty< RankTwoTensor > & _stress
const MaterialProperty< RankTwoTensor > & _dgrad_old
Material point based material property.
const unsigned int _component
The index of displacement component.
const MaterialProperty< RankTwoTensor > & _shape2
const MaterialProperty< Real > & _multi
Material point based material properties.
RankTwoTensorTempl< Real > transpose() const
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")
static const std::string k
Definition: NS.h:130
uint8_t dof_id_type
const MaterialProperty< RankTwoTensor > & _dgrad

◆ computeNonlocalResidual()

void HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computeNonlocalResidual ( )
overrideprotectedvirtual

Definition at line 56 of file HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD.C.

57 {
58  for (unsigned int nd = 0; nd < _nnodes; ++nd)
59  {
60  // calculation of residual contribution to current_node's neighbors
61  std::vector<dof_id_type> ivardofs(_nnodes);
62  ivardofs[0] = _ivardofs[nd];
63  std::vector<dof_id_type> neighbors = _pdmesh.getNeighbors(_current_elem->node_id(nd));
64  std::vector<dof_id_type> bonds = _pdmesh.getBonds(_current_elem->node_id(nd));
65 
66  dof_id_type nb_index =
67  std::find(neighbors.begin(), neighbors.end(), _current_elem->node_id(1 - nd)) -
68  neighbors.begin();
69  std::vector<dof_id_type> dg_neighbors =
70  _pdmesh.getBondDeformationGradientNeighbors(_current_elem->node_id(nd), nb_index);
71 
72  RealGradient origin_vec_nb;
73  Real node_vol_nb;
74 
75  for (unsigned int nb = 0; nb < dg_neighbors.size(); ++nb)
76  if (neighbors[dg_neighbors[nb]] != _current_elem->node_id(1 - nd) &&
77  _bond_status_var->getElementalValue(_pdmesh.elemPtr(bonds[dg_neighbors[nb]])) > 0.5)
78  {
79  ivardofs[1] = _pdmesh.nodePtr(neighbors[dg_neighbors[nb]])
80  ->dof_number(_sys.number(), _var.number(), 0);
81  origin_vec_nb = _pdmesh.getNodeCoord(neighbors[dg_neighbors[nb]]) -
82  _pdmesh.getNodeCoord(_current_elem->node_id(nd));
83  node_vol_nb = _pdmesh.getNodeVolume(neighbors[dg_neighbors[nb]]);
84 
85  for (unsigned int i = 0; i < _nnodes; ++i)
86  _local_re(i) = (i == 0 ? -1 : 1) * _multi[nd] * _horizon_radius[nd] /
87  origin_vec_nb.norm() *
88  ((_dgrad[nd].det() * _stress[nd] * _dgrad[nd].inverse().transpose()) *
89  _shape2[nd].inverse())
90  .row(_component) *
91  origin_vec_nb * node_vol_nb * _bond_status;
92 
93  // cache the residual contribution
94  addResiduals(_assembly, _local_re, ivardofs, _var.scalingFactor());
95 
96  if (_has_save_in)
97  {
98  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
99  for (unsigned int i = 0; i < _save_in.size(); ++i)
100  {
101  std::vector<dof_id_type> save_in_dofs(_nnodes);
102  save_in_dofs[0] = _current_elem->node_ptr(nd)->dof_number(
103  _save_in[i]->sys().number(), _save_in[i]->number(), 0);
104  save_in_dofs[1] =
105  _pdmesh.nodePtr(neighbors[dg_neighbors[nb]])
106  ->dof_number(_save_in[i]->sys().number(), _save_in[i]->number(), 0);
107 
108  _save_in[i]->sys().solution().add_vector(_local_re, save_in_dofs);
109  }
110  }
111  }
112  }
113 }
auto norm() const -> decltype(std::norm(Real()))
const MaterialProperty< RankTwoTensor > & _stress
void inverse(const std::vector< std::vector< Real >> &m, std::vector< std::vector< Real >> &m_inv)
const unsigned int _component
The index of displacement component.
std::vector< dof_id_type > _ivardofs
Current variable dof numbers for nodes i and j.
const MaterialProperty< RankTwoTensor > & _shape2
const MaterialProperty< Real > & _multi
Material point based material properties.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
uint8_t dof_id_type
const MaterialProperty< RankTwoTensor > & _dgrad

◆ 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()

void HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::computePDNonlocalOffDiagJacobian ( unsigned int  ,
unsigned int   
)
overrideprotectedvirtual

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 from MechanicsBasePD.

Definition at line 428 of file HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD.C.

430 {
431  if (_temp_coupled && jvar_num == _temp_var->number())
432  {
433  // no nonlocal contribution from temperature
434  }
436  {
437  // no nonlocal contribution from out of plane strain
438  }
439  else
440  {
441  for (unsigned int nd = 0; nd < _nnodes; ++nd)
442  {
443  RankFourTensor dSdFhat = computeDSDFhat(nd);
444  RankTwoTensor invF = _dgrad[nd].inverse();
445  Real detF = _dgrad[nd].det();
446  // calculation of jacobian contribution to current_node's neighbors
447  // NOT including the contribution to nodes i and j, which is considered as local off-diagonal
448  std::vector<dof_id_type> ivardofs(_nnodes), jvardofs(_nnodes);
449  ivardofs[0] = _current_elem->node_ptr(nd)->dof_number(_sys.number(), _var.number(), 0);
450  jvardofs[0] = _current_elem->node_ptr(nd)->dof_number(_sys.number(), jvar_num, 0);
451  std::vector<dof_id_type> neighbors = _pdmesh.getNeighbors(_current_elem->node_id(nd));
452  std::vector<dof_id_type> bonds = _pdmesh.getBonds(_current_elem->node_id(nd));
453 
454  dof_id_type nb_index =
455  std::find(neighbors.begin(), neighbors.end(), _current_elem->node_id(1 - nd)) -
456  neighbors.begin();
457  std::vector<dof_id_type> dg_neighbors =
458  _pdmesh.getBondDeformationGradientNeighbors(_current_elem->node_id(nd), nb_index);
459 
460  RealGradient origin_vec_nb1;
461  Real node_vol_nb1;
462 
463  for (unsigned int nb1 = 0; nb1 < dg_neighbors.size(); ++nb1)
464  if (_bond_status_var->getElementalValue(_pdmesh.elemPtr(bonds[dg_neighbors[nb1]])) > 0.5)
465  {
466  ivardofs[1] = _pdmesh.nodePtr(neighbors[dg_neighbors[nb1]])
467  ->dof_number(_sys.number(), _var.number(), 0);
468  origin_vec_nb1 = _pdmesh.getNodeCoord(neighbors[dg_neighbors[nb1]]) -
469  _pdmesh.getNodeCoord(_current_elem->node_id(nd));
470  node_vol_nb1 = _pdmesh.getNodeVolume(neighbors[dg_neighbors[nb1]]);
471 
472  Real vol_nb2, dJdU;
473  RealGradient origin_vec_nb2;
474  RankTwoTensor dFdUk, dPxdUky, dSdU, dinvFTdU;
475 
476  for (unsigned int nb2 = 0; nb2 < dg_neighbors.size(); ++nb2)
477  if (_bond_status_var->getElementalValue(_pdmesh.elemPtr(bonds[dg_neighbors[nb2]])) >
478  0.5)
479  {
480  jvardofs[1] = _pdmesh.nodePtr(neighbors[dg_neighbors[nb2]])
481  ->dof_number(_sys.number(), jvar_num, 0);
482  vol_nb2 = _pdmesh.getNodeVolume(neighbors[dg_neighbors[nb2]]);
483 
484  origin_vec_nb2 = _pdmesh.getNodeCoord(neighbors[dg_neighbors[nb2]]) -
485  _pdmesh.getNodeCoord(_current_elem->node_id(nd));
486 
487  dFdUk.zero();
488  for (unsigned int i = 0; i < _dim; ++i)
489  dFdUk(coupled_component, i) =
490  _horizon_radius[nd] / origin_vec_nb2.norm() * origin_vec_nb2(i) * vol_nb2;
491 
492  dFdUk *= _shape2[nd].inverse();
493 
494  // calculate dJ/du
495  dJdU = 0.0;
496  for (unsigned int i = 0; i < 3; ++i)
497  for (unsigned int j = 0; j < 3; ++j)
498  dJdU += detF * invF(j, i) * dFdUk(i, j);
499 
500  // calculate dS/du
501  dSdU = dSdFhat * dFdUk * _dgrad_old[nd].inverse();
502 
503  // calculate dinv(F)Tdu
504  dinvFTdU.zero();
505  for (unsigned int i = 0; i < 3; ++i)
506  for (unsigned int J = 0; J < 3; ++J)
507  for (unsigned int k = 0; k < 3; ++k)
508  for (unsigned int L = 0; L < 3; ++L)
509  dinvFTdU(i, J) += -invF(J, k) * invF(L, i) * dFdUk(k, L);
510 
511  // calculate the derivative of first Piola-Kirchhoff stress w.r.t displacements
512  dPxdUky = dJdU * _stress[nd] * invF.transpose() + detF * dSdU * invF.transpose() +
513  detF * _stress[nd] * dinvFTdU;
514 
515  for (unsigned int i = 0; i < _nnodes; ++i)
516  for (unsigned int j = 0; j < _nnodes; ++j)
517  _local_ke(i, j) = (i == 0 ? -1 : 1) * (j == 0 ? 0 : 1) * _multi[nd] *
518  _horizon_radius[nd] / origin_vec_nb1.norm() *
519  (dPxdUky * _shape2[nd].inverse()).row(_component) *
520  origin_vec_nb1 * node_vol_nb1 * _bond_status;
521 
522  addJacobian(_assembly, _local_ke, ivardofs, jvardofs, _var.scalingFactor());
523  }
524  }
525  }
526  }
527 }
RankFourTensor computeDSDFhat(unsigned int nd)
Function to compute derivative of stress with respect to derived deformation gradient.
const bool _temp_coupled
Temperature variable.
auto norm() const -> decltype(std::norm(Real()))
unsigned int number() const
const MaterialProperty< RankTwoTensor > & _stress
const bool _out_of_plane_strain_coupled
Parameters for out-of-plane strain in weak plane stress formulation.
const MaterialProperty< RankTwoTensor > & _dgrad_old
Material point based material property.
const unsigned int _component
The index of displacement component.
MooseVariable * _out_of_plane_strain_var
const MaterialProperty< RankTwoTensor > & _shape2
const MaterialProperty< Real > & _multi
Material point based material properties.
RankTwoTensorTempl< Real > transpose() const
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")
MooseVariable * _temp_var
static const std::string k
Definition: NS.h:130
uint8_t dof_id_type
const MaterialProperty< RankTwoTensor > & _dgrad

◆ 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 HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::validParams ( )
static

Definition at line 16 of file HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD.C.

17 {
19  params.addClassDescription(
20  "Class for calculating the residual and the Jacobian for Form II "
21  "of the horizon-stabilized peridynamic correspondence model under finite strain assumptions");
22 
23  params.addRequiredParam<unsigned int>(
24  "component",
25  "An integer corresponding to the variable this kernel acts on (0 for x, 1 for y, 2 for z)");
26 
27  return params;
28 }
void addRequiredParam(const std::string &name, const std::string &doc_string)
void addClassDescription(const std::string &doc_string)

Member Data Documentation

◆ _component

const unsigned int HorizonStabilizedFormIIFiniteStrainMechanicsNOSPD::_component
protected

◆ _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
protectedinherited

◆ _disp_var

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

◆ _E_inc

const MaterialProperty<RankTwoTensor>& MechanicsFiniteStrainBaseNOSPD::_E_inc
protectedinherited

◆ _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(), computeLocalJacobian(), ForceStabilizedSmallStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), HorizonStabilizedFormISmallStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), computeLocalOffDiagJacobian(), HorizonStabilizedFormIISmallStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), ForceStabilizedSmallStrainMechanicsNOSPD::computeLocalResidual(), HorizonStabilizedFormISmallStrainMechanicsNOSPD::computeLocalResidual(), computeLocalResidual(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeLocalResidual(), HorizonStabilizedFormIISmallStrainMechanicsNOSPD::computeLocalResidual(), ForceStabilizedSmallStrainMechanicsNOSPD::computeNonlocalJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeNonlocalJacobian(), HorizonStabilizedFormISmallStrainMechanicsNOSPD::computeNonlocalJacobian(), HorizonStabilizedFormIISmallStrainMechanicsNOSPD::computeNonlocalJacobian(), computeNonlocalJacobian(), HorizonStabilizedFormIISmallStrainMechanicsNOSPD::computeNonlocalResidual(), computeNonlocalResidual(), ForceStabilizedSmallStrainMechanicsNOSPD::computePDNonlocalOffDiagJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computePDNonlocalOffDiagJacobian(), HorizonStabilizedFormISmallStrainMechanicsNOSPD::computePDNonlocalOffDiagJacobian(), 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
protectedinherited

◆ _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(), computeLocalJacobian(), ForceStabilizedSmallStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), HorizonStabilizedFormISmallStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), computeLocalOffDiagJacobian(), HorizonStabilizedFormIISmallStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), ForceStabilizedSmallStrainMechanicsNOSPD::computeLocalResidual(), HorizonStabilizedFormISmallStrainMechanicsNOSPD::computeLocalResidual(), computeLocalResidual(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeLocalResidual(), HorizonStabilizedFormIISmallStrainMechanicsNOSPD::computeLocalResidual(), ForceStabilizedSmallStrainMechanicsNOSPD::computeNonlocalJacobian(), HorizonStabilizedFormISmallStrainMechanicsNOSPD::computeNonlocalJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computeNonlocalJacobian(), computeNonlocalJacobian(), HorizonStabilizedFormIISmallStrainMechanicsNOSPD::computeNonlocalJacobian(), HorizonStabilizedFormIISmallStrainMechanicsNOSPD::computeNonlocalResidual(), computeNonlocalResidual(), WeakPlaneStressNOSPD::computePDNonlocalOffDiagJacobian(), ForceStabilizedSmallStrainMechanicsNOSPD::computePDNonlocalOffDiagJacobian(), HorizonStabilizedFormISmallStrainMechanicsNOSPD::computePDNonlocalOffDiagJacobian(), HorizonStabilizedFormIFiniteStrainMechanicsNOSPD::computePDNonlocalOffDiagJacobian(), 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: