Kernel class for fictitious force stabilized conventional correspondence material model for small strain. More...

#include <ForceStabilizedSmallStrainMechanicsNOSPD.h>

Inheritance diagram for ForceStabilizedSmallStrainMechanicsNOSPD:

Public Member Functions
	ForceStabilizedSmallStrainMechanicsNOSPD (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	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...

virtual RankTwoTensor	computeDSDU (unsigned int component, unsigned int nd)
	Function to compute derivative of stress with respect to displacements. More...

Protected Attributes
const MaterialProperty< Real > &	_sf_coeff
	Bond based material property for fictitious stabilization force. More...

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 > &	_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 fictitious force stabilized conventional correspondence material model for small strain.

Definition at line 23 of file ForceStabilizedSmallStrainMechanicsNOSPD.h.

Constructor & Destructor Documentation

◆ ForceStabilizedSmallStrainMechanicsNOSPD()

ForceStabilizedSmallStrainMechanicsNOSPD::ForceStabilizedSmallStrainMechanicsNOSPD ( const InputParameters & parameters )

Definition at line 33 of file ForceStabilizedSmallStrainMechanicsNOSPD.C.

   : MechanicsBaseNOSPD(parameters),
     _sf_coeff(getMaterialProperty<Real>("stabilization_force_coeff")),
     _component(getParam<unsigned int>("component"))
 {
 }

Member Function Documentation

◆ computeDSDU()

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

protectedvirtualinherited

Function to compute derivative of stress with respect to displacements.

Parameters

component	The index of displacement component
nd	The 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.

 {
   // compute the derivative of stress w.r.t the solution components for small strain
   RankTwoTensor dSdU;
   if (component == 0)
     dSdU = _Jacobian_mult[nd] * 0.5 * (_ddgraddu[nd].transpose() + _ddgraddu[nd]);
   else if (component == 1)
     dSdU = _Jacobian_mult[nd] * 0.5 * (_ddgraddv[nd].transpose() + _ddgraddv[nd]);
   else if (component == 2)
     dSdU = _Jacobian_mult[nd] * 0.5 * (_ddgraddw[nd].transpose() + _ddgraddw[nd]);
  
   return dSdU;
 }

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

◆ computeLocalJacobian()

void ForceStabilizedSmallStrainMechanicsNOSPD::computeLocalJacobian ( )

overrideprotectedvirtual

Definition at line 64 of file ForceStabilizedSmallStrainMechanicsNOSPD.C.

 {
   // excludes dTi/dUj and dTj/dUi contribution which was considered as nonlocal contribution
   for (_i = 0; _i < _test.size(); ++_i)
     for (_j = 0; _j < _phi.size(); ++_j)
       _local_ke(_i, _j) += (_i == 0 ? -1 : 1) * _multi[_j] *
                            (computeDSDU(_component, _j) * _shape2[_j].inverse()).row(_component) *
                            _origin_vec_ij * _bond_status_ij;
  
   // compute local stabilization force jacobian
   Real diag = 0.5 * (_sf_coeff[0] + _sf_coeff[1]);
   for (_i = 0; _i < _test.size(); ++_i)
     for (_j = 0; _j < _phi.size(); ++_j)
       _local_ke(_i, _j) += (_i == _j ? 1 : -1) * diag * _bond_status_ij;
 }

◆ computeLocalOffDiagJacobian()

void ForceStabilizedSmallStrainMechanicsNOSPD::computeLocalOffDiagJacobian ( unsigned int )

overrideprotectedvirtual

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

Parameters

coupled_component The coupled variable number

Reimplemented from MechanicsBasePD.

Definition at line 138 of file ForceStabilizedSmallStrainMechanicsNOSPD.C.

 {
   _local_ke.zero();
   if (coupled_component == 3)
   {
     std::vector<RankTwoTensor> dSdT(_nnodes);
     for (unsigned int nd = 0; nd < _nnodes; ++nd)
       for (unsigned int es = 0; es < _deigenstrain_dT.size(); ++es)
         dSdT[nd] = -(_Jacobian_mult[nd] * (*_deigenstrain_dT[es])[nd]);
  
     for (_i = 0; _i < _test.size(); ++_i)
       for (_j = 0; _j < _phi.size(); ++_j)
         _local_ke(_i, _j) += (_i == 0 ? -1 : 1) * _multi[_j] *
                              (dSdT[_j] * _shape2[_j].inverse()).row(_component) * _origin_vec_ij *
                              _bond_status_ij;
   }
   else
   {
     for (_i = 0; _i < _test.size(); ++_i)
       for (_j = 0; _j < _phi.size(); ++_j)
         _local_ke(_i, _j) +=
             (_i == 0 ? -1 : 1) * _multi[_j] *
             (computeDSDU(coupled_component, _j) * _shape2[_j].inverse()).row(_component) *
             _origin_vec_ij * _bond_status_ij;
   }
 }

◆ computeLocalResidual()

void ForceStabilizedSmallStrainMechanicsNOSPD::computeLocalResidual ( )

overrideprotectedvirtual

Definition at line 42 of file ForceStabilizedSmallStrainMechanicsNOSPD.C.

 {
   Real sforce = 0.5 * (_sf_coeff[0] + _sf_coeff[1]) *
                 (_disp_var[_component]->getNodalValue(*_current_elem->node_ptr(1)) -
                  _disp_var[_component]->getNodalValue(*_current_elem->node_ptr(0)));
  
   // For small strain assumptions, stress measures, i.e., Cauchy stress (Sigma), the first
   // Piola-Kirchhoff stress (P), and the second Piola-Kirchhoff stress (S) are approximately the
   // same. Thus, the nodal force state tensors are calculated using the Cauchy stresses,
   // i.e., T = Sigma * inv(Shape) * xi * multi.
   // Cauchy stress is calculated as Sigma = C * E.
  
   std::vector<Real> nodal_force_comp(_nnodes);
   for (unsigned int nd = 0; nd < _nnodes; ++nd)
     nodal_force_comp[nd] = _multi[nd] * (_stress[nd] * _shape2[nd].inverse()).row(_component) *
                            (nd == 0 ? 1 : -1) * _origin_vec_ij;
  
   _local_re(0) = -(nodal_force_comp[0] - nodal_force_comp[1] + sforce) * _bond_status_ij;
   _local_re(1) = -_local_re(0);
 }

◆ computeNonlocalJacobian()

void ForceStabilizedSmallStrainMechanicsNOSPD::computeNonlocalJacobian ( )

overrideprotectedvirtual

Definition at line 81 of file ForceStabilizedSmallStrainMechanicsNOSPD.C.

 {
   // includes dTi/dUj and dTj/dUi contributions
   for (unsigned int cur_nd = 0; cur_nd < _nnodes; ++cur_nd)
   {
     // calculation of jacobian contribution to current_node's neighbors
     std::vector<dof_id_type> dof(_nnodes);
     dof[0] = _current_elem->node_ptr(cur_nd)->dof_number(_sys.number(), _var.number(), 0);
     std::vector<dof_id_type> neighbors = _pdmesh.getNeighbors(_current_elem->node_id(cur_nd));
     std::vector<dof_id_type> bonds = _pdmesh.getBonds(_current_elem->node_id(cur_nd));
     for (unsigned int k = 0; k < neighbors.size(); ++k)
     {
       Node * node_k = _pdmesh.nodePtr(neighbors[k]);
       dof[1] = node_k->dof_number(_sys.number(), _var.number(), 0);
       Real vol_k = _pdmesh.getPDNodeVolume(neighbors[k]);
  
       // obtain bond ik's origin vector
       RealGradient origin_vec_ijk = *node_k - *_pdmesh.nodePtr(_current_elem->node_id(cur_nd));
  
       RankTwoTensor dFdUk;
       dFdUk.zero();
       for (unsigned int j = 0; j < _dim; ++j)
         dFdUk(_component, j) =
             vol_k * _horiz_rad[cur_nd] / origin_vec_ijk.norm() * origin_vec_ijk(j);
  
       dFdUk *= _shape2[cur_nd].inverse();
       RankTwoTensor dSxdUkx = _Jacobian_mult[cur_nd] * 0.5 * (dFdUk.transpose() + dFdUk);
  
       // bond status for bond k
       Real bond_status_ijk = _bond_status_var->getElementalValue(_pdmesh.elemPtr(bonds[k]));
  
       _local_ke.resize(_test.size(), _phi.size());
       _local_ke.zero();
       for (_i = 0; _i < _test.size(); ++_i)
         for (_j = 0; _j < _phi.size(); ++_j)
           _local_ke(_i, _j) = (_i == 0 ? -1 : 1) * (_j == 0 ? 0 : 1) * _multi[cur_nd] *
                               (dSxdUkx * _shape2[cur_nd].inverse()).row(_component) *
                               _origin_vec_ij * _bond_status_ij * bond_status_ijk;
  
       _assembly.cacheJacobianBlock(_local_ke, _ivardofs_ij, dof, _var.scalingFactor());
  
       if (_has_diag_save_in)
       {
         unsigned int rows = _test.size();
         DenseVector<Real> diag(rows);
         for (unsigned int i = 0; i < rows; ++i)
           diag(i) = _local_ke(i, i);
  
         Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
         for (unsigned int i = 0; i < _diag_save_in.size(); ++i)
           _diag_save_in[i]->sys().solution().add_vector(diag, _diag_save_in[i]->dofIndices());
       }
     }
   }
 }

◆ computeOffDiagJacobian()

void MechanicsBasePD::computeOffDiagJacobian ( MooseVariableFEBase & jvar )

overridevirtualinherited

Definition at line 71 of file MechanicsBasePD.C.

 {
   prepare();
  
   if (jvar.number() == _var.number())
     computeJacobian();
   else
   {
     unsigned int coupled_component = 0;
     bool active = false;
  
     for (unsigned int i = 0; i < _dim; ++i)
       if (jvar.number() == _disp_var[i]->number())
       {
         coupled_component = i;
         active = true;
       }
  
     if (_temp_coupled && jvar.number() == _temp_var->number())
     {
       coupled_component = 3;
       active = true;
     }
  
     if (_out_of_plane_strain_coupled && jvar.number() == _out_of_plane_strain_var->number())
     {
       coupled_component = 4;
       active = true;
     }
  
     if (active)
     {
       DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), jvar.number());
       _local_ke.resize(ke.m(), ke.n());
       _local_ke.zero();
  
       computeLocalOffDiagJacobian(coupled_component);
  
       ke += _local_ke;
  
       if (_use_full_jacobian)
         computePDNonlocalOffDiagJacobian(jvar.number(), coupled_component);
     }
   }
 }

◆ computePDNonlocalOffDiagJacobian()

void ForceStabilizedSmallStrainMechanicsNOSPD::computePDNonlocalOffDiagJacobian	(	unsigned int	,
		unsigned int
	)

overrideprotectedvirtual

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

Parameters

jvar_num	The number of the first coupled variable
coupled_component	The component number of the second coupled variable

Reimplemented from MechanicsBasePD.

Definition at line 167 of file ForceStabilizedSmallStrainMechanicsNOSPD.C.

 {
   if (coupled_component != 3)
   {
     for (unsigned int cur_nd = 0; cur_nd < _nnodes; ++cur_nd)
     {
       // calculation of jacobian contribution to current_node's neighbors
       std::vector<dof_id_type> jvardofs_ijk(_nnodes);
       jvardofs_ijk[0] = _current_elem->node_ptr(cur_nd)->dof_number(_sys.number(), jvar_num, 0);
       std::vector<dof_id_type> neighbors = _pdmesh.getNeighbors(_current_elem->node_id(cur_nd));
       std::vector<dof_id_type> bonds = _pdmesh.getBonds(_current_elem->node_id(cur_nd));
       for (unsigned int k = 0; k < neighbors.size(); ++k)
       {
         Node * node_k = _pdmesh.nodePtr(neighbors[k]);
         jvardofs_ijk[1] = node_k->dof_number(_sys.number(), jvar_num, 0);
         Real vol_k = _pdmesh.getPDNodeVolume(neighbors[k]);
  
         // obtain bond k's origin vector
         RealGradient origin_vec_ijk = *node_k - *_pdmesh.nodePtr(_current_elem->node_id(cur_nd));
  
         RankTwoTensor dFdUk;
         dFdUk.zero();
         for (unsigned int j = 0; j < _dim; ++j)
           dFdUk(coupled_component, j) =
               _horiz_rad[cur_nd] / origin_vec_ijk.norm() * origin_vec_ijk(j) * vol_k;
  
         dFdUk *= _shape2[cur_nd].inverse();
         RankTwoTensor dSxdUky = _Jacobian_mult[cur_nd] * 0.5 * (dFdUk.transpose() + dFdUk);
  
         // bond status for bond k
         Real bond_status_ijk = _bond_status_var->getElementalValue(_pdmesh.elemPtr(bonds[k]));
  
         _local_ke.zero();
         for (_i = 0; _i < _test.size(); ++_i)
           for (_j = 0; _j < _phi.size(); ++_j)
             _local_ke(_i, _j) = (_i == 0 ? -1 : 1) * (_j == 0 ? 0 : 1) * _multi[cur_nd] *
                                 (dSxdUky * _shape2[cur_nd].inverse()).row(_component) *
                                 _origin_vec_ij * _bond_status_ij * bond_status_ijk;
  
         _assembly.cacheJacobianBlock(_local_ke, _ivardofs_ij, jvardofs_ijk, _var.scalingFactor());
       }
     }
   }
 }

◆ initialSetup()

void MechanicsBasePD::initialSetup ( )

overridevirtualinherited

Definition at line 47 of file MechanicsBasePD.C.

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

◆ prepare()

void MechanicsBasePD::prepare ( )

overridevirtualinherited

Definition at line 53 of file MechanicsBasePD.C.

 {
   PeridynamicsKernelBase::prepare();
  
   _ivardofs_ij.resize(_nnodes);
  
   for (unsigned int i = 0; i < _nnodes; ++i)
     _ivardofs_ij[i] = _current_elem->node_ptr(i)->dof_number(_sys.number(), _var.number(), 0);
  
   for (unsigned int i = 0; i < _dim; ++i)
     _cur_ori_ij(i) = _origin_vec_ij(i) + _disp_var[i]->getNodalValue(*_current_elem->node_ptr(1)) -
                      _disp_var[i]->getNodalValue(*_current_elem->node_ptr(0));
  
   _cur_len_ij = _cur_ori_ij.norm();
   _cur_ori_ij /= _cur_len_ij;
 }