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

Base kernel class for bond-associated correspondence material models. More...

#include <MechanicsBaseNOSPD.h>

Inheritance diagram for MechanicsBaseNOSPD:
[legend]

Public Member Functions

 MechanicsBaseNOSPD (const InputParameters &parameters)
 
virtual void computeOffDiagJacobian (MooseVariableFEBase &jvar) override
 
virtual void computeLocalOffDiagJacobian (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
 

Protected Member Functions

virtual RankTwoTensor computeDSDU (unsigned int component, unsigned int nd)
 Function to compute derivative of stress 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_ij
 Current variable dof numbers for nodes i and j. More...
 
RealGradient _cur_ori_ij
 Vector of bond in current configuration. More...
 
Real _cur_len_ij
 Current bond length. More...
 
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 bond-associated correspondence material models.

Definition at line 22 of file MechanicsBaseNOSPD.h.

Constructor & Destructor Documentation

◆ MechanicsBaseNOSPD()

MechanicsBaseNOSPD::MechanicsBaseNOSPD ( const InputParameters &  parameters)

Definition at line 29 of file MechanicsBaseNOSPD.C.

30  : MechanicsBasePD(parameters),
31  _multi(getMaterialProperty<Real>("multi")),
32  _stress(getMaterialProperty<RankTwoTensor>("stress")),
33  _shape2(getMaterialProperty<RankTwoTensor>("rank_two_shape_tensor")),
34  _dgrad(getMaterialProperty<RankTwoTensor>("deformation_gradient")),
35  _ddgraddu(getMaterialProperty<RankTwoTensor>("ddeformation_gradient_du")),
36  _ddgraddv(getMaterialProperty<RankTwoTensor>("ddeformation_gradient_dv")),
37  _ddgraddw(getMaterialProperty<RankTwoTensor>("ddeformation_gradient_dw")),
38  _Jacobian_mult(getMaterialProperty<RankFourTensor>("Jacobian_mult")),
39  _eigenstrain_names(getParam<std::vector<MaterialPropertyName>>("eigenstrain_names")),
41 {
42  if (_temp_coupled)
43  for (unsigned int i = 0; i < _deigenstrain_dT.size(); ++i)
44  _deigenstrain_dT[i] =
45  &getMaterialPropertyDerivative<RankTwoTensor>(_eigenstrain_names[i], _temp_var->name());
46 }

Member Function Documentation

◆ computeDSDU()

RankTwoTensor MechanicsBaseNOSPD::computeDSDU ( unsigned int  component,
unsigned int  nd 
)
protectedvirtual

Function to compute derivative of stress 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

Reimplemented in FiniteStrainMechanicsNOSPD.

Definition at line 49 of file MechanicsBaseNOSPD.C.

50 {
51  // compute the derivative of stress w.r.t the solution components for small strain
52  RankTwoTensor dSdU;
53  if (component == 0)
54  dSdU = _Jacobian_mult[nd] * 0.5 * (_ddgraddu[nd].transpose() + _ddgraddu[nd]);
55  else if (component == 1)
56  dSdU = _Jacobian_mult[nd] * 0.5 * (_ddgraddv[nd].transpose() + _ddgraddv[nd]);
57  else if (component == 2)
58  dSdU = _Jacobian_mult[nd] * 0.5 * (_ddgraddw[nd].transpose() + _ddgraddw[nd]);
59 
60  return dSdU;
61 }

Referenced by GeneralizedPlaneStrainOffDiagNOSPD::computeDispFullOffDiagJacobianScalar(), GeneralizedPlaneStrainOffDiagNOSPD::computeDispPartialOffDiagJacobianScalar(), SmallStrainMechanicsNOSPD::computeLocalJacobian(), ForceStabilizedSmallStrainMechanicsNOSPD::computeLocalJacobian(), WeakPlaneStressNOSPD::computeLocalOffDiagJacobian(), SmallStrainMechanicsNOSPD::computeLocalOffDiagJacobian(), and ForceStabilizedSmallStrainMechanicsNOSPD::computeLocalOffDiagJacobian().

◆ computeLocalOffDiagJacobian()

virtual void MechanicsBasePD::computeLocalOffDiagJacobian ( 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 MechanicsBPD, MechanicsOSPD, FiniteStrainMechanicsNOSPD, ForceStabilizedSmallStrainMechanicsNOSPD, SmallStrainMechanicsNOSPD, and WeakPlaneStressNOSPD.

Definition at line 35 of file MechanicsBasePD.h.

35 {};

Referenced by MechanicsBasePD::computeOffDiagJacobian().

◆ computeOffDiagJacobian()

void MechanicsBasePD::computeOffDiagJacobian ( MooseVariableFEBase &  jvar)
overridevirtualinherited

Definition at line 71 of file MechanicsBasePD.C.

72 {
73  prepare();
74 
75  if (jvar.number() == _var.number())
76  computeJacobian();
77  else
78  {
79  unsigned int coupled_component = 0;
80  bool active = false;
81 
82  for (unsigned int i = 0; i < _dim; ++i)
83  if (jvar.number() == _disp_var[i]->number())
84  {
85  coupled_component = i;
86  active = true;
87  }
88 
89  if (_temp_coupled && jvar.number() == _temp_var->number())
90  {
91  coupled_component = 3;
92  active = true;
93  }
94 
95  if (_out_of_plane_strain_coupled && jvar.number() == _out_of_plane_strain_var->number())
96  {
97  coupled_component = 4;
98  active = true;
99  }
100 
101  if (active)
102  {
103  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), jvar.number());
104  _local_ke.resize(ke.m(), ke.n());
105  _local_ke.zero();
106 
107  computeLocalOffDiagJacobian(coupled_component);
108 
109  ke += _local_ke;
110 
111  if (_use_full_jacobian)
112  computePDNonlocalOffDiagJacobian(jvar.number(), coupled_component);
113  }
114  }
115 }

◆ 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 MechanicsOSPD, FiniteStrainMechanicsNOSPD, ForceStabilizedSmallStrainMechanicsNOSPD, SmallStrainMechanicsNOSPD, and WeakPlaneStressNOSPD.

Definition at line 42 of file MechanicsBasePD.h.

43  {};

Referenced by MechanicsBasePD::computeOffDiagJacobian().

◆ initialSetup()

void MechanicsBasePD::initialSetup ( )
overridevirtualinherited

Definition at line 47 of file MechanicsBasePD.C.

48 {
49  _orientation = &_assembly.getFE(FEType(), 1)->get_dxyzdxi();
50 }

◆ prepare()

void MechanicsBasePD::prepare ( )
overridevirtualinherited

Definition at line 53 of file MechanicsBasePD.C.

54 {
56 
57  _ivardofs_ij.resize(_nnodes);
58 
59  for (unsigned int i = 0; i < _nnodes; ++i)
60  _ivardofs_ij[i] = _current_elem->node_ptr(i)->dof_number(_sys.number(), _var.number(), 0);
61 
62  for (unsigned int i = 0; i < _dim; ++i)
63  _cur_ori_ij(i) = _origin_vec_ij(i) + _disp_var[i]->getNodalValue(*_current_elem->node_ptr(1)) -
64  _disp_var[i]->getNodalValue(*_current_elem->node_ptr(0));
65 
66  _cur_len_ij = _cur_ori_ij.norm();
68 }

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

Member Data Documentation

◆ _cur_len_ij

Real MechanicsBasePD::_cur_len_ij
protectedinherited

◆ _cur_ori_ij

RealGradient MechanicsBasePD::_cur_ori_ij
protectedinherited

◆ _ddgraddu

const MaterialProperty<RankTwoTensor>& MechanicsBaseNOSPD::_ddgraddu
protected

◆ _ddgraddv

const MaterialProperty<RankTwoTensor>& MechanicsBaseNOSPD::_ddgraddv
protected

◆ _ddgraddw

const MaterialProperty<RankTwoTensor>& MechanicsBaseNOSPD::_ddgraddw
protected

◆ _deigenstrain_dT

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

◆ _dgrad

const MaterialProperty<RankTwoTensor>& MechanicsBaseNOSPD::_dgrad
protected

◆ _disp_var

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

◆ _eigenstrain_names

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

Definition at line 45 of file MechanicsBaseNOSPD.h.

Referenced by MechanicsBaseNOSPD().

◆ _ivardofs_ij

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

◆ _Jacobian_mult

const MaterialProperty<RankFourTensor>& MechanicsBaseNOSPD::_Jacobian_mult
protected

◆ _multi

const MaterialProperty<Real>& MechanicsBaseNOSPD::_multi
protected

◆ _ndisp

unsigned int MechanicsBasePD::_ndisp
protectedinherited

number of displacement components

Definition at line 58 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 66 of file MechanicsBasePD.h.

Referenced by MechanicsBasePD::initialSetup().

◆ _out_of_plane_strain_coupled

const bool MechanicsBasePD::_out_of_plane_strain_coupled
protectedinherited

Parameters for out-of-plane strain in weak plane stress formulation.

Definition at line 61 of file MechanicsBasePD.h.

Referenced by MechanicsBasePD::computeOffDiagJacobian().

◆ _out_of_plane_strain_var

MooseVariable* MechanicsBasePD::_out_of_plane_strain_var
protectedinherited

Definition at line 62 of file MechanicsBasePD.h.

Referenced by MechanicsBasePD::computeOffDiagJacobian().

◆ _shape2

const MaterialProperty<RankTwoTensor>& MechanicsBaseNOSPD::_shape2
protected

◆ _stress

const MaterialProperty<RankTwoTensor>& MechanicsBaseNOSPD::_stress
protected

◆ _temp_coupled

const bool MechanicsBasePD::_temp_coupled
protectedinherited

◆ _temp_var

MooseVariable* MechanicsBasePD::_temp_var
protectedinherited

The documentation for this class was generated from the following files:
MechanicsBasePD::MechanicsBasePD
MechanicsBasePD(const InputParameters &parameters)
Definition: MechanicsBasePD.C:29
MechanicsBasePD::_cur_ori_ij
RealGradient _cur_ori_ij
Vector of bond in current configuration.
Definition: MechanicsBasePD.h:72
MechanicsBasePD::_temp_var
MooseVariable * _temp_var
Definition: MechanicsBasePD.h:54
MechanicsBasePD::_temp_coupled
const bool _temp_coupled
Temperature variable.
Definition: MechanicsBasePD.h:53
MechanicsBaseNOSPD::_dgrad
const MaterialProperty< RankTwoTensor > & _dgrad
Definition: MechanicsBaseNOSPD.h:40
MechanicsBaseNOSPD::_shape2
const MaterialProperty< RankTwoTensor > & _shape2
Definition: MechanicsBaseNOSPD.h:39
MechanicsBaseNOSPD::_multi
const MaterialProperty< Real > & _multi
Material point based material properties.
Definition: MechanicsBaseNOSPD.h:37
MechanicsBaseNOSPD::_Jacobian_mult
const MaterialProperty< RankFourTensor > & _Jacobian_mult
Definition: MechanicsBaseNOSPD.h:44
MechanicsBaseNOSPD::_stress
const MaterialProperty< RankTwoTensor > & _stress
Definition: MechanicsBaseNOSPD.h:38
MechanicsBaseNOSPD::_ddgraddw
const MaterialProperty< RankTwoTensor > & _ddgraddw
Definition: MechanicsBaseNOSPD.h:43
MechanicsBaseNOSPD::_ddgraddu
const MaterialProperty< RankTwoTensor > & _ddgraddu
Definition: MechanicsBaseNOSPD.h:41
MechanicsBaseNOSPD::_ddgraddv
const MaterialProperty< RankTwoTensor > & _ddgraddv
Definition: MechanicsBaseNOSPD.h:42
MechanicsBasePD::_orientation
const std::vector< RealGradient > * _orientation
Vector of bond in current configuration.
Definition: MechanicsBasePD.h:66
MechanicsBasePD::computeLocalOffDiagJacobian
virtual void computeLocalOffDiagJacobian(unsigned int)
Function to compute local contribution to the off-diagonal Jacobian at the current nodes.
Definition: MechanicsBasePD.h:35
MechanicsBasePD::_out_of_plane_strain_var
MooseVariable * _out_of_plane_strain_var
Definition: MechanicsBasePD.h:62
PeridynamicsKernelBase::prepare
virtual void prepare()
Function to precalculate data which will be used in the derived classes.
Definition: PeridynamicsKernelBase.C:43
MechanicsBasePD::_out_of_plane_strain_coupled
const bool _out_of_plane_strain_coupled
Parameters for out-of-plane strain in weak plane stress formulation.
Definition: MechanicsBasePD.h:61
MechanicsBasePD::_ivardofs_ij
std::vector< dof_id_type > _ivardofs_ij
Current variable dof numbers for nodes i and j.
Definition: MechanicsBasePD.h:69
MechanicsBaseNOSPD::_eigenstrain_names
const std::vector< MaterialPropertyName > _eigenstrain_names
Definition: MechanicsBaseNOSPD.h:45
MaterialTensorCalculatorTools::component
Real component(const SymmTensor &symm_tensor, unsigned int index)
Definition: MaterialTensorCalculatorTools.C:16
MechanicsBasePD::computePDNonlocalOffDiagJacobian
virtual void computePDNonlocalOffDiagJacobian(unsigned int, unsigned int)
Function to compute nonlocal contribution to the off-diagonal Jacobian at the current nodes.
Definition: MechanicsBasePD.h:42
MechanicsBasePD::prepare
virtual void prepare() override
Definition: MechanicsBasePD.C:53
RankTwoTensorTempl< Real >
MechanicsBaseNOSPD::_deigenstrain_dT
std::vector< const MaterialProperty< RankTwoTensor > * > _deigenstrain_dT
Definition: MechanicsBaseNOSPD.h:46
MechanicsBasePD::_cur_len_ij
Real _cur_len_ij
Current bond length.
Definition: MechanicsBasePD.h:75
MechanicsBasePD::_disp_var
std::vector< MooseVariable * > _disp_var
displacement variables
Definition: MechanicsBasePD.h:50