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

Kernel class for bond based peridynamic solid mechanics models. More...

#include <MechanicsBPD.h>

Inheritance diagram for MechanicsBPD:
[legend]

Public Types

typedef DerivativeMaterialPropertyNameInterface::SymbolName SymbolName
 

Public Member Functions

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

Protected Attributes

const unsigned int _component
 The index of displcement 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< Real > & _bond_local_force
 Bond based material properties. More...
 
const MaterialProperty< Real > & _bond_local_dfdU
 
const MaterialProperty< Real > & _bond_local_dfdT
 
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 bond based peridynamic solid mechanics models.

Definition at line 17 of file MechanicsBPD.h.

Constructor & Destructor Documentation

◆ MechanicsBPD()

MechanicsBPD::MechanicsBPD ( const InputParameters parameters)

Definition at line 28 of file MechanicsBPD.C.

29  : MechanicsBasePD(parameters),
30  _bond_local_force(getMaterialProperty<Real>("bond_local_force")),
31  _bond_local_dfdU(getMaterialProperty<Real>("bond_dfdU")),
32  _bond_local_dfdT(getMaterialProperty<Real>("bond_dfdT")),
33  _component(getParam<unsigned int>("component"))
34 {
35 }
const MaterialProperty< Real > & _bond_local_dfdT
Definition: MechanicsBPD.h:32
const unsigned int _component
The index of displcement component.
Definition: MechanicsBPD.h:36
MechanicsBasePD(const InputParameters &parameters)
const MaterialProperty< Real > & _bond_local_force
Bond based material properties.
Definition: MechanicsBPD.h:30
const MaterialProperty< Real > & _bond_local_dfdU
Definition: MechanicsBPD.h:31

Member Function Documentation

◆ computeLocalJacobian()

void MechanicsBPD::computeLocalJacobian ( )
overrideprotectedvirtual

Definition at line 45 of file MechanicsBPD.C.

46 {
51 
52  for (unsigned int i = 0; i < _nnodes; ++i)
53  for (unsigned int j = 0; j < _nnodes; ++j)
54  _local_ke(i, j) += (i == j ? 1 : -1) * diag * _bond_status;
55 }
auto norm() const -> decltype(std::norm(Real()))
const unsigned int _component
The index of displcement component.
Definition: MechanicsBPD.h:36
RealGradient _current_unit_vec
Unit vector of bond in current configuration.
const MaterialProperty< Real > & _bond_local_force
Bond based material properties.
Definition: MechanicsBPD.h:30
RealGradient _current_vec
Vector of bond in current configuration.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const MaterialProperty< Real > & _bond_local_dfdU
Definition: MechanicsBPD.h:31
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")

◆ computeLocalOffDiagJacobian()

void MechanicsBPD::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 58 of file MechanicsBPD.C.

59 {
60  if (_temp_coupled && jvar_num == _temp_var->number())
61  {
62  for (unsigned int i = 0; i < _nnodes; ++i)
63  for (unsigned int j = 0; j < _nnodes; ++j)
64  _local_ke(i, j) +=
65  (i == 1 ? 1 : -1) * _current_unit_vec(_component) * _bond_local_dfdT[j] * _bond_status;
66  }
67  else
68  {
69  Real val =
72  _current_unit_vec(coupled_component) / _current_vec.norm();
73  for (unsigned int i = 0; i < _nnodes; ++i)
74  for (unsigned int j = 0; j < _nnodes; ++j)
75  _local_ke(i, j) += (i == j ? 1 : -1) * val * _bond_status;
76  }
77 }
const bool _temp_coupled
Temperature variable.
const MaterialProperty< Real > & _bond_local_dfdT
Definition: MechanicsBPD.h:32
auto norm() const -> decltype(std::norm(Real()))
unsigned int number() const
const unsigned int _component
The index of displcement component.
Definition: MechanicsBPD.h:36
RealGradient _current_unit_vec
Unit vector of bond in current configuration.
const MaterialProperty< Real > & _bond_local_force
Bond based material properties.
Definition: MechanicsBPD.h:30
RealGradient _current_vec
Vector of bond in current configuration.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const MaterialProperty< Real > & _bond_local_dfdU
Definition: MechanicsBPD.h:31
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
MooseVariable * _temp_var

◆ computeLocalResidual()

void MechanicsBPD::computeLocalResidual ( )
overrideprotectedvirtual

Definition at line 38 of file MechanicsBPD.C.

39 {
40  _local_re(0) = -_bond_local_force[0] * _current_unit_vec(_component) * _bond_status;
41  _local_re(1) = -_local_re(0);
42 }
const unsigned int _component
The index of displcement component.
Definition: MechanicsBPD.h:36
RealGradient _current_unit_vec
Unit vector of bond in current configuration.
const MaterialProperty< Real > & _bond_local_force
Bond based material properties.
Definition: MechanicsBPD.h:30

◆ 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(), GeneralizedPlaneStrainOffDiagOSPD::computeOffDiagJacobianScalar(), and GeneralizedPlaneStrainOffDiagNOSPD::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 MechanicsBPD::validParams ( )
static

Definition at line 15 of file MechanicsBPD.C.

16 {
18  params.addClassDescription("Class for calculating the residual and Jacobian for the bond-based "
19  "peridynamic mechanics formulation");
20 
21  params.addRequiredParam<unsigned int>(
22  "component",
23  "An integer corresponding to the variable this kernel acts on (0 for x, 1 for y, 2 for z)");
24 
25  return params;
26 }
static InputParameters validParams()
void addRequiredParam(const std::string &name, const std::string &doc_string)
void addClassDescription(const std::string &doc_string)

Member Data Documentation

◆ _bond_local_dfdT

const MaterialProperty<Real>& MechanicsBPD::_bond_local_dfdT
protected

Definition at line 32 of file MechanicsBPD.h.

Referenced by computeLocalOffDiagJacobian().

◆ _bond_local_dfdU

const MaterialProperty<Real>& MechanicsBPD::_bond_local_dfdU
protected

Definition at line 31 of file MechanicsBPD.h.

Referenced by computeLocalJacobian(), and computeLocalOffDiagJacobian().

◆ _bond_local_force

const MaterialProperty<Real>& MechanicsBPD::_bond_local_force
protected

Bond based material properties.

Definition at line 30 of file MechanicsBPD.h.

Referenced by computeLocalJacobian(), computeLocalOffDiagJacobian(), and computeLocalResidual().

◆ _component

const unsigned int MechanicsBPD::_component
protected

The index of displcement component.

Definition at line 36 of file MechanicsBPD.h.

Referenced by computeLocalJacobian(), computeLocalOffDiagJacobian(), and computeLocalResidual().

◆ _current_unit_vec

RealGradient MechanicsBasePD::_current_unit_vec
protectedinherited

◆ _current_vec

RealGradient MechanicsBasePD::_current_vec
protectedinherited

◆ _disp_var

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

◆ _ivardofs

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

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

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