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

Kernel class for ordinary state based peridynamic solid mechanics models for small strain. More...

#include <MechanicsOSPD.h>

Inheritance diagram for MechanicsOSPD:
[legend]

Public Member Functions

 MechanicsOSPD (const InputParameters &parameters)
 
virtual void computeOffDiagJacobian (MooseVariableFEBase &jvar) override
 
virtual void initialSetup () override
 
virtual void prepare () override
 

Protected Member Functions

virtual void computeLocalResidual () override
 
virtual void computeNonlocalResidual () override
 
virtual void computeLocalJacobian () override
 
virtual void computeNonlocalJacobian () override
 
void computeLocalOffDiagJacobian (unsigned int coupled_component) override
 Function to compute local contribution to the off-diagonal Jacobian at the current nodes. More...
 
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...
 

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_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 > & _bond_force_ij
 Bond based material properties. More...
 
const MaterialProperty< Real > & _bond_force_i_j
 
const MaterialProperty< Real > & _bond_dfdU_ij
 
const MaterialProperty< Real > & _bond_dfdU_i_j
 
const MaterialProperty< Real > & _bond_dfdT_ij
 
const MaterialProperty< Real > & _bond_dfdT_i_j
 
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 ordinary state based peridynamic solid mechanics models for small strain.

Definition at line 22 of file MechanicsOSPD.h.

Constructor & Destructor Documentation

◆ MechanicsOSPD()

MechanicsOSPD::MechanicsOSPD ( const InputParameters &  parameters)

Definition at line 30 of file MechanicsOSPD.C.

31  : MechanicsBasePD(parameters),
32  _bond_force_ij(getMaterialProperty<Real>("bond_force_ij")),
33  _bond_force_i_j(getMaterialProperty<Real>("bond_force_i_j")),
34  _bond_dfdU_ij(getMaterialProperty<Real>("bond_dfdU_ij")),
35  _bond_dfdU_i_j(getMaterialProperty<Real>("bond_dfdU_i_j")),
36  _bond_dfdT_ij(getMaterialProperty<Real>("bond_dfdT_ij")),
37  _bond_dfdT_i_j(getMaterialProperty<Real>("bond_dfdT_i_j")),
38  _component(getParam<unsigned int>("component"))
39 {
40 }

Member Function Documentation

◆ computeLocalJacobian()

void MechanicsOSPD::computeLocalJacobian ( )
overrideprotectedvirtual

Definition at line 108 of file MechanicsOSPD.C.

109 {
110  const Real val =
113 
114  for (_i = 0; _i < _test.size(); ++_i)
115  for (_j = 0; _j < _phi.size(); ++_j)
116  _local_ke(_i, _j) += (_i == _j ? 1 : -1) * val * _bond_status_ij;
117 }

◆ computeLocalOffDiagJacobian()

void MechanicsOSPD::computeLocalOffDiagJacobian ( 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 120 of file MechanicsOSPD.C.

121 {
122  if (coupled_component == 3)
123  {
124  for (_i = 0; _i < _test.size(); ++_i)
125  for (_j = 0; _j < _phi.size(); ++_j)
126  _local_ke(_i, _j) +=
127  (_i == 1 ? 1 : -1) * _cur_ori_ij(_component) * _bond_dfdT_ij[0] * _bond_status_ij;
128  }
129  else
130  {
131  const Real val =
132  _cur_ori_ij(_component) * _cur_ori_ij(coupled_component) * _bond_dfdU_ij[0] -
133  _bond_force_ij[0] * _cur_ori_ij(_component) * _cur_ori_ij(coupled_component) / _cur_len_ij;
134  for (_i = 0; _i < _test.size(); ++_i)
135  for (_j = 0; _j < _phi.size(); ++_j)
136  _local_ke(_i, _j) += (_i == _j ? 1 : -1) * val * _bond_status_ij;
137  }
138 }

◆ computeLocalResidual()

void MechanicsOSPD::computeLocalResidual ( )
overrideprotectedvirtual

Definition at line 43 of file MechanicsOSPD.C.

44 {
45  // H term
46  _local_re(0) = -_bond_force_ij[0] * _cur_ori_ij(_component) * _bond_status_ij;
47  _local_re(1) = -_local_re(0);
48 }

◆ computeNonlocalJacobian()

void MechanicsOSPD::computeNonlocalJacobian ( )
overrideprotectedvirtual

Definition at line 141 of file MechanicsOSPD.C.

142 {
143  for (unsigned int cur_nd = 0; cur_nd < _nnodes; ++cur_nd)
144  {
145  // calculation of jacobian contribution to current_node's neighbors
146  std::vector<dof_id_type> ivardofs(_nnodes);
147  ivardofs[cur_nd] = _ivardofs_ij[cur_nd];
148  std::vector<dof_id_type> neighbors = _pdmesh.getNeighbors(_current_elem->node_id(cur_nd));
149  std::vector<dof_id_type> bonds = _pdmesh.getBonds(_current_elem->node_id(cur_nd));
150  for (unsigned int k = 0; k < neighbors.size(); ++k)
151  {
152  const Node * node_k = _pdmesh.nodePtr(neighbors[k]);
153  ivardofs[1 - cur_nd] = node_k->dof_number(_sys.number(), _var.number(), 0);
154  const Real vol_k = _pdmesh.getPDNodeVolume(neighbors[k]);
155 
156  // obtain bond ik's origin length and current orientation
157  RealGradient origin_ori_ijk = *node_k - *_pdmesh.nodePtr(_current_elem->node_id(cur_nd));
158 
159  RealGradient cur_ori_ijk(_dim);
160  for (unsigned int j = 0; j < _dim; ++j)
161  cur_ori_ijk(j) = origin_ori_ijk(j) + _disp_var[j]->getNodalValue(*node_k) -
162  _disp_var[j]->getNodalValue(*_current_elem->node_ptr(cur_nd));
163 
164  cur_ori_ijk /= cur_ori_ijk.norm();
165 
166  // bond status for bond k
167  const Real bond_status_ijk = _bond_status_var->getElementalValue(_pdmesh.elemPtr(bonds[k]));
168 
169  const Real val = (cur_nd == 0 ? 1 : -1) * cur_ori_ijk(_component) * _cur_ori_ij(_component) *
170  _bond_dfdU_i_j[cur_nd];
171 
172  _local_ke.zero();
173  for (_i = 0; _i < _test.size(); ++_i)
174  for (_j = 0; _j < _phi.size(); ++_j)
175  _local_ke(_i, _j) += (_i == _j ? 1 : -1) * val / origin_ori_ijk.norm() * vol_k *
176  _bond_status_ij * bond_status_ijk;
177 
178  _assembly.cacheJacobianBlock(_local_ke, ivardofs, _ivardofs_ij, _var.scalingFactor());
179 
180  if (_has_diag_save_in)
181  {
182  unsigned int rows = _test.size();
183  DenseVector<Real> diag(rows);
184  for (unsigned int i = 0; i < rows; ++i)
185  diag(i) = _local_ke(i, i);
186 
187  diag(1 - cur_nd) = 0;
188 
189  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
190  for (unsigned int i = 0; i < _diag_save_in.size(); ++i)
191  {
192  std::vector<dof_id_type> diag_save_in_dofs(2);
193  diag_save_in_dofs[cur_nd] = _current_elem->node_ptr(cur_nd)->dof_number(
194  _diag_save_in[i]->sys().number(), _diag_save_in[i]->number(), 0);
195  diag_save_in_dofs[1 - cur_nd] =
196  node_k->dof_number(_diag_save_in[i]->sys().number(), _diag_save_in[i]->number(), 0);
197 
198  _diag_save_in[i]->sys().solution().add_vector(diag, diag_save_in_dofs);
199  }
200  }
201 
202  const Real val2 = _bond_force_i_j[cur_nd] *
203  (1.0 - cur_ori_ijk(_component) * cur_ori_ijk(_component)) /
204  cur_ori_ijk.norm();
205 
206  _local_ke.zero();
207  for (_i = 0; _i < _test.size(); ++_i)
208  for (_j = 0; _j < _phi.size(); ++_j)
209  _local_ke(_i, _j) += (_i == _j ? 1 : -1) * val2 / origin_ori_ijk.norm() * vol_k *
210  _bond_status_ij * bond_status_ijk;
211 
212  _assembly.cacheJacobianBlock(_local_ke, ivardofs, ivardofs, _var.scalingFactor());
213 
214  if (_has_diag_save_in)
215  {
216  unsigned int rows = _test.size();
217  DenseVector<Real> diag(rows);
218  for (unsigned int i = 0; i < rows; ++i)
219  diag(i) = _local_ke(i, i);
220 
221  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
222  for (unsigned int i = 0; i < _diag_save_in.size(); ++i)
223  {
224  std::vector<dof_id_type> diag_save_in_dofs(2);
225  diag_save_in_dofs[cur_nd] = _current_elem->node_ptr(cur_nd)->dof_number(
226  _diag_save_in[i]->sys().number(), _diag_save_in[i]->number(), 0);
227  diag_save_in_dofs[1 - cur_nd] =
228  node_k->dof_number(_diag_save_in[i]->sys().number(), _diag_save_in[i]->number(), 0);
229 
230  _diag_save_in[i]->sys().solution().add_vector(diag, diag_save_in_dofs);
231  }
232  }
233  }
234  }
235 }

◆ computeNonlocalResidual()

void MechanicsOSPD::computeNonlocalResidual ( )
overrideprotectedvirtual

Definition at line 51 of file MechanicsOSPD.C.

52 {
53  // P and Q terms
54  for (unsigned int cur_nd = 0; cur_nd < 2; ++cur_nd)
55  {
56  // calculation of residual contribution to current_node's neighbors
57  std::vector<dof_id_type> ivardofs(2);
58  ivardofs[cur_nd] = _ivardofs_ij[cur_nd];
59  std::vector<dof_id_type> neighbors = _pdmesh.getNeighbors(_current_elem->node_id(cur_nd));
60  std::vector<dof_id_type> bonds = _pdmesh.getBonds(_current_elem->node_id(cur_nd));
61  for (unsigned int k = 0; k < neighbors.size(); ++k)
62  {
63  const Node * node_k = _pdmesh.nodePtr(neighbors[k]);
64  ivardofs[1 - cur_nd] = node_k->dof_number(_sys.number(), _var.number(), 0);
65  const Real vol_k = _pdmesh.getPDNodeVolume(neighbors[k]);
66 
67  // obtain bond ik's origin length and current orientation
68  RealGradient origin_ori_ijk = *node_k - *_pdmesh.nodePtr(_current_elem->node_id(cur_nd));
69 
70  RealGradient cur_ori_ijk(_dim);
71  for (unsigned int j = 0; j < _dim; ++j)
72  cur_ori_ijk(j) = origin_ori_ijk(j) + _disp_var[j]->getNodalValue(*node_k) -
73  _disp_var[j]->getNodalValue(*_current_elem->node_ptr(cur_nd));
74 
75  cur_ori_ijk /= cur_ori_ijk.norm();
76 
77  // bond status for bond k
78  const Real bond_status_ijk = _bond_status_var->getElementalValue(_pdmesh.elemPtr(bonds[k]));
79 
80  _local_re(0) = (cur_nd == 0 ? -1 : 1) * _bond_force_i_j[cur_nd] * vol_k /
81  origin_ori_ijk.norm() * cur_ori_ijk(_component) * _bond_status_ij *
82  bond_status_ijk;
83  _local_re(1) = -_local_re(0);
84 
85  // cache the residual contribution to node_i and its neighbor k using their global dof indices
86  _assembly.cacheResidualNodes(_local_re, ivardofs);
87 
88  // save in the displacement residuals
89  if (_has_save_in)
90  {
91  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
92  for (unsigned int i = 0; i < _save_in.size(); ++i)
93  {
94  std::vector<dof_id_type> save_in_dofs(2);
95  save_in_dofs[cur_nd] = _current_elem->node_ptr(cur_nd)->dof_number(
96  _save_in[i]->sys().number(), _save_in[i]->number(), 0);
97  save_in_dofs[1 - cur_nd] =
98  node_k->dof_number(_save_in[i]->sys().number(), _save_in[i]->number(), 0);
99 
100  _save_in[i]->sys().solution().add_vector(_local_re, save_in_dofs);
101  }
102  }
103  }
104  }
105 }

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

void MechanicsOSPD::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 238 of file MechanicsOSPD.C.

240 {
241  std::vector<dof_id_type> jvardofs_ij(_nnodes);
242  for (unsigned int nd = 0; nd < _nnodes; ++nd)
243  jvardofs_ij[nd] = _current_elem->node_ptr(nd)->dof_number(_sys.number(), jvar_num, 0);
244 
245  for (unsigned int cur_nd = 0; cur_nd < _nnodes; ++cur_nd)
246  {
247  // calculation of jacobian contribution to current_node's neighbors
248  std::vector<dof_id_type> ivardofs(_nnodes), jvardofs(_nnodes);
249  ivardofs[cur_nd] = _ivardofs_ij[cur_nd];
250  jvardofs[cur_nd] = jvardofs_ij[cur_nd];
251  std::vector<dof_id_type> neighbors = _pdmesh.getNeighbors(_current_elem->node_id(cur_nd));
252  std::vector<dof_id_type> bonds = _pdmesh.getBonds(_current_elem->node_id(cur_nd));
253  for (unsigned int k = 0; k < neighbors.size(); ++k)
254  {
255  const Node * node_k = _pdmesh.nodePtr(neighbors[k]);
256  ivardofs[1 - cur_nd] = node_k->dof_number(_sys.number(), _var.number(), 0);
257  jvardofs[1 - cur_nd] = node_k->dof_number(_sys.number(), jvar_num, 0);
258  const Real vol_k = _pdmesh.getPDNodeVolume(neighbors[k]);
259 
260  // obtain bond k's origin length and current orientation
261  RealGradient origin_ori_ijk = *node_k - *_pdmesh.nodePtr(_current_elem->node_id(cur_nd));
262 
263  RealGradient cur_ori_ijk;
264  for (unsigned int j = 0; j < _dim; ++j)
265  cur_ori_ijk(j) = origin_ori_ijk(j) + _disp_var[j]->getNodalValue(*node_k) -
266  _disp_var[j]->getNodalValue(*_current_elem->node_ptr(cur_nd));
267 
268  cur_ori_ijk /= cur_ori_ijk.norm();
269 
270  // bond status for bond k
271  const Real bond_status_ijk = _bond_status_var->getElementalValue(_pdmesh.elemPtr(bonds[k]));
272 
273  _local_ke.zero();
274  if (coupled_component == 3)
275  {
276  const Real val =
277  cur_ori_ijk(_component) * _bond_dfdT_i_j[cur_nd] / origin_ori_ijk.norm() * vol_k;
278  _local_ke(0, cur_nd) += (cur_nd == 0 ? -1 : 1) * val * _bond_status_ij * bond_status_ijk;
279  _local_ke(1, cur_nd) += (cur_nd == 0 ? 1 : -1) * val * _bond_status_ij * bond_status_ijk;
280  }
281  else
282  {
283  const Real val = (cur_nd == 0 ? 1 : -1) * cur_ori_ijk(_component) *
284  _cur_ori_ij(coupled_component) * _bond_dfdU_i_j[cur_nd] /
285  origin_ori_ijk.norm() * vol_k;
286  for (_i = 0; _i < _test.size(); ++_i)
287  for (_j = 0; _j < _phi.size(); ++_j)
288  _local_ke(_i, _j) += (_i == _j ? 1 : -1) * val * _bond_status_ij * bond_status_ijk;
289  }
290 
291  _assembly.cacheJacobianBlock(_local_ke, ivardofs, jvardofs_ij, _var.scalingFactor());
292  }
293  }
294 }

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

Member Data Documentation

◆ _bond_dfdT_i_j

const MaterialProperty<Real>& MechanicsOSPD::_bond_dfdT_i_j
protected

Definition at line 44 of file MechanicsOSPD.h.

Referenced by computePDNonlocalOffDiagJacobian().

◆ _bond_dfdT_ij

const MaterialProperty<Real>& MechanicsOSPD::_bond_dfdT_ij
protected

Definition at line 43 of file MechanicsOSPD.h.

Referenced by computeLocalOffDiagJacobian().

◆ _bond_dfdU_i_j

const MaterialProperty<Real>& MechanicsOSPD::_bond_dfdU_i_j
protected

Definition at line 42 of file MechanicsOSPD.h.

Referenced by computeNonlocalJacobian(), and computePDNonlocalOffDiagJacobian().

◆ _bond_dfdU_ij

const MaterialProperty<Real>& MechanicsOSPD::_bond_dfdU_ij
protected

Definition at line 41 of file MechanicsOSPD.h.

Referenced by computeLocalJacobian(), and computeLocalOffDiagJacobian().

◆ _bond_force_i_j

const MaterialProperty<Real>& MechanicsOSPD::_bond_force_i_j
protected

Definition at line 40 of file MechanicsOSPD.h.

Referenced by computeNonlocalJacobian(), and computeNonlocalResidual().

◆ _bond_force_ij

const MaterialProperty<Real>& MechanicsOSPD::_bond_force_ij
protected

Bond based material properties.

Definition at line 39 of file MechanicsOSPD.h.

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

◆ _component

const unsigned int MechanicsOSPD::_component
protected

◆ _cur_len_ij

Real MechanicsBasePD::_cur_len_ij
protectedinherited

◆ _cur_ori_ij

RealGradient MechanicsBasePD::_cur_ori_ij
protectedinherited

◆ _disp_var

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

◆ _ivardofs_ij

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

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

◆ _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
MechanicsOSPD::_component
const unsigned int _component
The index of displacement component.
Definition: MechanicsOSPD.h:48
libMesh::RealGradient
VectorValue< Real > RealGradient
Definition: GrainForceAndTorqueInterface.h:17
MechanicsOSPD::_bond_force_ij
const MaterialProperty< Real > & _bond_force_ij
Bond based material properties.
Definition: MechanicsOSPD.h:39
MechanicsOSPD::_bond_dfdU_i_j
const MaterialProperty< Real > & _bond_dfdU_i_j
Definition: MechanicsOSPD.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
MechanicsOSPD::_bond_force_i_j
const MaterialProperty< Real > & _bond_force_i_j
Definition: MechanicsOSPD.h:40
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
MechanicsOSPD::_bond_dfdT_ij
const MaterialProperty< Real > & _bond_dfdT_ij
Definition: MechanicsOSPD.h:43
MechanicsOSPD::_bond_dfdT_i_j
const MaterialProperty< Real > & _bond_dfdT_i_j
Definition: MechanicsOSPD.h:44
MechanicsOSPD::_bond_dfdU_ij
const MaterialProperty< Real > & _bond_dfdU_ij
Definition: MechanicsOSPD.h:41
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
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