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

ComputeCrackTipEnrichmentSmallStrain calculates the sum of standard strain and enrichement strain. More...

#include <ComputeCrackTipEnrichmentSmallStrain.h>

Inheritance diagram for ComputeCrackTipEnrichmentSmallStrain:
[legend]

Public Member Functions

 ComputeCrackTipEnrichmentSmallStrain (const InputParameters &parameters)
 
virtual ~ComputeCrackTipEnrichmentSmallStrain ()
 
virtual unsigned int crackTipEnrichementFunctionAtPoint (const Point &point, std::vector< Real > &B)
 calculate the enrichment function values at point More...
 
virtual unsigned int crackTipEnrichementFunctionDerivativeAtPoint (const Point &point, std::vector< RealVectorValue > &dB)
 calculate the enrichment function derivatives at point More...
 
void rotateFromCrackFrontCoordsToGlobal (const RealVectorValue &vector, RealVectorValue &rotated_vector, const unsigned int point_index)
 rotate a vector from crack front coordinate to global cooridate More...
 

Protected Member Functions

virtual void computeProperties () override
 
virtual void computeQpProperties () override
 
void initialSetup () override
 
virtual void initQpStatefulProperties () override
 
virtual void displacementIntegrityCheck ()
 

Protected Attributes

std::vector< Real > _enrich_disp
 enrichment displacement More...
 
std::vector< RealVectorValue > _grad_enrich_disp
 gradient of enrichment displacement More...
 
std::vector< std::vector< MooseVariableFEBase * > > _enrich_variable
 enrichment displacement variables More...
 
const VariablePhiValue & _phi
 the current shape functions More...
 
const VariablePhiGradient & _grad_phi
 gradient of the shape function More...
 
unsigned int _ndisp
 Coupled displacement variables. More...
 
std::vector< const VariableValue * > _disp
 
std::vector< const VariableGradient * > _grad_disp
 
std::string _base_name
 
MaterialProperty< RankTwoTensor > & _mechanical_strain
 
MaterialProperty< RankTwoTensor > & _total_strain
 
std::vector< MaterialPropertyName > _eigenstrain_names
 
std::vector< const MaterialProperty< RankTwoTensor > * > _eigenstrains
 
const MaterialProperty< RankTwoTensor > * _global_strain
 
bool _volumetric_locking_correction
 
const Real & _current_elem_volume
 

Private Attributes

std::vector< Real > _B
 enrichment function value More...
 
std::vector< RealVectorValue > _dBX
 derivatives of enrichment function respect to global cooridnate More...
 
std::vector< RealVectorValue > _dBx
 derivatives of enrichment function respect to crack front cooridnate More...
 
std::vector< std::vector< Real > > _BI
 enrichment function at node I More...
 
const std::vector< std::vector< Real > > * _fe_phi
 shape function More...
 
const std::vector< std::vector< RealGradient > > * _fe_dphi
 gradient of shape function More...
 
NonlinearSystem * _nl
 
const NumericVector< Number > * _sln
 

Detailed Description

ComputeCrackTipEnrichmentSmallStrain calculates the sum of standard strain and enrichement strain.

Definition at line 30 of file ComputeCrackTipEnrichmentSmallStrain.h.

Constructor & Destructor Documentation

◆ ComputeCrackTipEnrichmentSmallStrain()

ComputeCrackTipEnrichmentSmallStrain::ComputeCrackTipEnrichmentSmallStrain ( const InputParameters &  parameters)

Definition at line 31 of file ComputeCrackTipEnrichmentSmallStrain.C.

33  : ComputeStrainBase(parameters),
34  EnrichmentFunctionCalculation(&getUserObject<CrackFrontDefinition>("crack_front_definition")),
35  _enrich_disp(3),
38  _phi(_assembly.phi()),
39  _grad_phi(_assembly.gradPhi()),
40  _B(4),
41  _dBX(4),
42  _dBx(4)
43 {
44  for (unsigned int i = 0; i < _enrich_variable.size(); ++i)
45  _enrich_variable[i].resize(_ndisp);
46 
47  const std::vector<NonlinearVariableName> & nl_vnames =
48  getParam<std::vector<NonlinearVariableName>>("enrichment_displacements");
49 
50  if (_ndisp == 2 && nl_vnames.size() != 8)
51  mooseError("The number of enrichment displacements should be total 8 for 2D.");
52  else if (_ndisp == 3 && nl_vnames.size() != 12)
53  mooseError("The number of enrichment displacements should be total 12 for 3D.");
54 
55  NonlinearSystem & nl = _fe_problem.getNonlinearSystem();
56  _nl = &(_fe_problem.getNonlinearSystem());
57 
58  for (unsigned int j = 0; j < _ndisp; ++j)
59  for (unsigned int i = 0; i < 4; ++i)
60  _enrich_variable[i][j] = &(nl.getVariable(0, nl_vnames[j * 4 + i]));
61 
62  if (_ndisp == 2)
63  _BI.resize(4); // QUAD4
64  else if (_ndisp == 3)
65  _BI.resize(8); // HEX8
66 
67  for (unsigned int i = 0; i < _BI.size(); ++i)
68  _BI[i].resize(4);
69 }
const VariablePhiValue & _phi
the current shape functions
std::vector< std::vector< MooseVariableFEBase * > > _enrich_variable
enrichment displacement variables
std::vector< Real > _B
enrichment function value
std::vector< RealVectorValue > _grad_enrich_disp
gradient of enrichment displacement
EnrichmentFunctionCalculation(const CrackFrontDefinition *crack_front_definition)
ComputeStrainBase(const InputParameters &parameters)
std::vector< RealVectorValue > _dBX
derivatives of enrichment function respect to global cooridnate
const VariablePhiGradient & _grad_phi
gradient of the shape function
std::vector< std::vector< Real > > _BI
enrichment function at node I
std::vector< Real > _enrich_disp
enrichment displacement
unsigned int _ndisp
Coupled displacement variables.
std::vector< RealVectorValue > _dBx
derivatives of enrichment function respect to crack front cooridnate

◆ ~ComputeCrackTipEnrichmentSmallStrain()

virtual ComputeCrackTipEnrichmentSmallStrain::~ComputeCrackTipEnrichmentSmallStrain ( )
inlinevirtual

Definition at line 35 of file ComputeCrackTipEnrichmentSmallStrain.h.

35 {}

Member Function Documentation

◆ computeProperties()

void ComputeCrackTipEnrichmentSmallStrain::computeProperties ( )
overrideprotectedvirtual

Definition at line 121 of file ComputeCrackTipEnrichmentSmallStrain.C.

122 {
123  FEType fe_type(Utility::string_to_enum<Order>("first"),
124  Utility::string_to_enum<FEFamily>("lagrange"));
125  const unsigned int dim = _current_elem->dim();
126  std::unique_ptr<FEBase> fe(FEBase::build(dim, fe_type));
127  fe->attach_quadrature_rule(_qrule);
128  _fe_phi = &(fe->get_phi());
129  _fe_dphi = &(fe->get_dphi());
130 
131  if (isBoundaryMaterial())
132  fe->reinit(_current_elem, _current_side);
133  else
134  fe->reinit(_current_elem);
135 
136  for (unsigned int i = 0; i < _BI.size(); ++i)
137  crackTipEnrichementFunctionAtPoint(*(_current_elem->get_node(i)), _BI[i]);
138 
139  ComputeStrainBase::computeProperties();
140 }
virtual unsigned int crackTipEnrichementFunctionAtPoint(const Point &point, std::vector< Real > &B)
calculate the enrichment function values at point
const std::vector< std::vector< Real > > * _fe_phi
shape function
std::vector< std::vector< Real > > _BI
enrichment function at node I
const std::vector< std::vector< RealGradient > > * _fe_dphi
gradient of shape function

◆ computeQpProperties()

void ComputeCrackTipEnrichmentSmallStrain::computeQpProperties ( )
overrideprotectedvirtual

Definition at line 72 of file ComputeCrackTipEnrichmentSmallStrain.C.

73 {
74  crackTipEnrichementFunctionAtPoint(_q_point[_qp], _B);
75  unsigned int crack_front_point_index =
77 
78  for (unsigned int i = 0; i < 4; ++i)
79  rotateFromCrackFrontCoordsToGlobal(_dBx[i], _dBX[i], crack_front_point_index);
80 
81  _sln = _nl->currentSolution();
82 
83  for (unsigned int m = 0; m < _ndisp; ++m)
84  {
85  _enrich_disp[m] = 0.0;
86  _grad_enrich_disp[m].zero();
87  for (unsigned int i = 0; i < _current_elem->n_nodes(); ++i)
88  {
89  Node * node_i = _current_elem->get_node(i);
90  for (unsigned int j = 0; j < 4; ++j)
91  {
92  dof_id_type dof = node_i->dof_number(_nl->number(), _enrich_variable[j][m]->number(), 0);
93  Real soln = (*_sln)(dof);
94  _enrich_disp[m] += (*_fe_phi)[i][_qp] * (_B[j] - _BI[i][j]) * soln;
95  RealVectorValue grad_B(_dBX[j]);
96  _grad_enrich_disp[m] +=
97  ((*_fe_dphi)[i][_qp] * (_B[j] - _BI[i][j]) + (*_fe_phi)[i][_qp] * grad_B) * soln;
98  }
99  }
100  }
101 
102  RankTwoTensor grad_tensor_enrich(
104 
105  RankTwoTensor enrich_strain = (grad_tensor_enrich + grad_tensor_enrich.transpose()) / 2.0;
106 
107  RankTwoTensor grad_tensor((*_grad_disp[0])[_qp], (*_grad_disp[1])[_qp], (*_grad_disp[2])[_qp]);
108 
109  _total_strain[_qp] = (grad_tensor + grad_tensor.transpose()) / 2.0;
110 
111  _total_strain[_qp] += enrich_strain;
112 
113  _mechanical_strain[_qp] = _total_strain[_qp];
114 
115  // Remove the Eigen strain
116  for (auto es : _eigenstrains)
117  _mechanical_strain[_qp] -= (*es)[_qp];
118 }
std::vector< std::vector< MooseVariableFEBase * > > _enrich_variable
enrichment displacement variables
std::vector< Real > _B
enrichment function value
std::vector< RealVectorValue > _grad_enrich_disp
gradient of enrichment displacement
virtual unsigned int crackTipEnrichementFunctionAtPoint(const Point &point, std::vector< Real > &B)
calculate the enrichment function values at point
void rotateFromCrackFrontCoordsToGlobal(const RealVectorValue &vector, RealVectorValue &rotated_vector, const unsigned int point_index)
rotate a vector from crack front coordinate to global cooridate
std::vector< RealVectorValue > _dBX
derivatives of enrichment function respect to global cooridnate
const std::vector< std::vector< Real > > * _fe_phi
shape function
std::vector< std::vector< Real > > _BI
enrichment function at node I
MaterialProperty< RankTwoTensor > & _mechanical_strain
std::vector< Real > _enrich_disp
enrichment displacement
unsigned int _ndisp
Coupled displacement variables.
std::vector< RealVectorValue > _dBx
derivatives of enrichment function respect to crack front cooridnate
virtual unsigned int crackTipEnrichementFunctionDerivativeAtPoint(const Point &point, std::vector< RealVectorValue > &dB)
calculate the enrichment function derivatives at point
MaterialProperty< RankTwoTensor > & _total_strain
std::vector< const MaterialProperty< RankTwoTensor > * > _eigenstrains
std::vector< const VariableGradient * > _grad_disp

◆ crackTipEnrichementFunctionAtPoint()

unsigned int EnrichmentFunctionCalculation::crackTipEnrichementFunctionAtPoint ( const Point &  point,
std::vector< Real > &  B 
)
virtualinherited

calculate the enrichment function values at point

Returns
the closest crack front index

Definition at line 19 of file EnrichmentFunctionCalculation.C.

Referenced by computeProperties(), CrackTipEnrichmentStressDivergenceTensors::computeQpJacobian(), CrackTipEnrichmentStressDivergenceTensors::computeQpOffDiagJacobian(), computeQpProperties(), and CrackTipEnrichmentStressDivergenceTensors::computeQpResidual().

21 {
22  unsigned int crack_front_point_index =
24 
25  if (MooseUtils::absoluteFuzzyEqual(_r, 0.0))
26  mooseError("EnrichmentFunctionCalculation: the distance between a point and the crack "
27  "tip/front is zero.");
28 
29  Real st = std::sin(_theta);
30  Real st2 = std::sin(_theta / 2.0);
31  Real ct2 = std::cos(_theta / 2.0);
32  Real sr = std::sqrt(_r);
33 
34  B[0] = sr * st2;
35  B[1] = sr * ct2;
36  B[2] = sr * st2 * st;
37  B[3] = sr * ct2 * st;
38 
39  return crack_front_point_index;
40 }
void calculateRThetaToCrackFront(const Point qp, const unsigned int point_index, Real &r, Real &theta) const
calculate r and theta in the crack front polar cooridnate
const CrackFrontDefinition * _crack_front_definition

◆ crackTipEnrichementFunctionDerivativeAtPoint()

unsigned int EnrichmentFunctionCalculation::crackTipEnrichementFunctionDerivativeAtPoint ( const Point &  point,
std::vector< RealVectorValue > &  dB 
)
virtualinherited

calculate the enrichment function derivatives at point

Returns
the closest crack front index

Definition at line 43 of file EnrichmentFunctionCalculation.C.

Referenced by CrackTipEnrichmentStressDivergenceTensors::computeQpJacobian(), CrackTipEnrichmentStressDivergenceTensors::computeQpOffDiagJacobian(), computeQpProperties(), and CrackTipEnrichmentStressDivergenceTensors::computeQpResidual().

45 {
46  unsigned int crack_front_point_index =
48 
49  if (MooseUtils::absoluteFuzzyEqual(_r, 0.0))
50  mooseError("EnrichmentFunctionCalculation: the distance between a point and the crack "
51  "tip/front is zero.");
52 
53  Real st = std::sin(_theta);
54  Real ct = std::cos(_theta);
55  Real st2 = std::sin(_theta / 2.0);
56  Real ct2 = std::cos(_theta / 2.0);
57  Real st15 = std::sin(1.5 * _theta);
58  Real ct15 = std::cos(1.5 * _theta);
59  Real sr = std::sqrt(_r);
60 
61  dB[0](0) = -0.5 / sr * st2;
62  dB[0](1) = 0.5 / sr * ct2;
63  dB[0](2) = 0.0;
64  dB[1](0) = 0.5 / sr * ct2;
65  dB[1](1) = 0.5 / sr * st2;
66  dB[1](2) = 0.0;
67  dB[2](0) = -0.5 / sr * st15 * st;
68  dB[2](1) = 0.5 / sr * (st2 + st15 * ct);
69  dB[2](2) = 0.0;
70  dB[3](0) = -0.5 / sr * ct15 * st;
71  dB[3](1) = 0.5 / sr * (ct2 + ct15 * ct);
72  dB[3](2) = 0.0;
73 
74  return crack_front_point_index;
75 }
void calculateRThetaToCrackFront(const Point qp, const unsigned int point_index, Real &r, Real &theta) const
calculate r and theta in the crack front polar cooridnate
const CrackFrontDefinition * _crack_front_definition

◆ displacementIntegrityCheck()

void ComputeStrainBase::displacementIntegrityCheck ( )
protectedvirtualinherited

Reimplemented in Compute2DFiniteStrain, Compute2DIncrementalStrain, and Compute2DSmallStrain.

Definition at line 86 of file ComputeStrainBase.C.

Referenced by ComputeStrainBase::initialSetup().

87 {
88  // Checking for consistency between mesh size and length of the provided displacements vector
89  if (_ndisp != _mesh.dimension())
90  paramError(
91  "displacements",
92  "The number of variables supplied in 'displacements' must match the mesh dimension.");
93 }
unsigned int _ndisp
Coupled displacement variables.

◆ initialSetup()

void ComputeStrainBase::initialSetup ( )
overrideprotectedinherited

Definition at line 67 of file ComputeStrainBase.C.

Referenced by ComputeIncrementalStrainBase::initialSetup(), and ComputeAxisymmetric1DSmallStrain::initialSetup().

68 {
70  // fetch coupled variables and gradients (as stateful properties if necessary)
71  for (unsigned int i = 0; i < _ndisp; ++i)
72  {
73  _disp[i] = &coupledValue("displacements", i);
74  _grad_disp[i] = &coupledGradient("displacements", i);
75  }
76 
77  // set unused dimensions to zero
78  for (unsigned i = _ndisp; i < 3; ++i)
79  {
80  _disp[i] = &_zero;
81  _grad_disp[i] = &_grad_zero;
82  }
83 }
virtual void displacementIntegrityCheck()
std::vector< const VariableValue * > _disp
unsigned int _ndisp
Coupled displacement variables.
std::vector< const VariableGradient * > _grad_disp

◆ initQpStatefulProperties()

void ComputeStrainBase::initQpStatefulProperties ( )
overrideprotectedvirtualinherited

Reimplemented in ComputeCosseratIncrementalSmallStrain, and ComputeIncrementalStrainBase.

Definition at line 96 of file ComputeStrainBase.C.

97 {
98  _mechanical_strain[_qp].zero();
99  _total_strain[_qp].zero();
100 }
MaterialProperty< RankTwoTensor > & _mechanical_strain
MaterialProperty< RankTwoTensor > & _total_strain

◆ rotateFromCrackFrontCoordsToGlobal()

void EnrichmentFunctionCalculation::rotateFromCrackFrontCoordsToGlobal ( const RealVectorValue &  vector,
RealVectorValue &  rotated_vector,
const unsigned int  point_index 
)
inherited

rotate a vector from crack front coordinate to global cooridate

Parameters
rotated_vectorrotated vector

Definition at line 78 of file EnrichmentFunctionCalculation.C.

Referenced by CrackTipEnrichmentStressDivergenceTensors::computeQpJacobian(), CrackTipEnrichmentStressDivergenceTensors::computeQpOffDiagJacobian(), computeQpProperties(), and CrackTipEnrichmentStressDivergenceTensors::computeQpResidual().

81 {
82  rotated_vector = _crack_front_definition->rotateFromCrackFrontCoordsToGlobal(vector, point_index);
83 }
const CrackFrontDefinition * _crack_front_definition
RealVectorValue rotateFromCrackFrontCoordsToGlobal(const RealVectorValue vector, const unsigned int point_index) const
rotate a vector from crack front cartesian coordinate to global cartesian coordinate ...

Member Data Documentation

◆ _B

std::vector<Real> ComputeCrackTipEnrichmentSmallStrain::_B
private

enrichment function value

Definition at line 59 of file ComputeCrackTipEnrichmentSmallStrain.h.

Referenced by computeQpProperties().

◆ _base_name

std::string ComputeStrainBase::_base_name
protectedinherited

Definition at line 43 of file ComputeStrainBase.h.

Referenced by ComputeStrainBase::ComputeStrainBase().

◆ _BI

std::vector<std::vector<Real> > ComputeCrackTipEnrichmentSmallStrain::_BI
private

enrichment function at node I

Definition at line 65 of file ComputeCrackTipEnrichmentSmallStrain.h.

Referenced by ComputeCrackTipEnrichmentSmallStrain(), computeProperties(), and computeQpProperties().

◆ _current_elem_volume

const Real& ComputeStrainBase::_current_elem_volume
protectedinherited

◆ _dBX

std::vector<RealVectorValue> ComputeCrackTipEnrichmentSmallStrain::_dBX
private

derivatives of enrichment function respect to global cooridnate

Definition at line 61 of file ComputeCrackTipEnrichmentSmallStrain.h.

Referenced by computeQpProperties().

◆ _dBx

std::vector<RealVectorValue> ComputeCrackTipEnrichmentSmallStrain::_dBx
private

derivatives of enrichment function respect to crack front cooridnate

Definition at line 63 of file ComputeCrackTipEnrichmentSmallStrain.h.

Referenced by computeQpProperties().

◆ _disp

std::vector<const VariableValue *> ComputeStrainBase::_disp
protectedinherited

◆ _eigenstrain_names

std::vector<MaterialPropertyName> ComputeStrainBase::_eigenstrain_names
protectedinherited

◆ _eigenstrains

std::vector<const MaterialProperty<RankTwoTensor> *> ComputeStrainBase::_eigenstrains
protectedinherited

◆ _enrich_disp

std::vector<Real> ComputeCrackTipEnrichmentSmallStrain::_enrich_disp
protected

enrichment displacement

Definition at line 43 of file ComputeCrackTipEnrichmentSmallStrain.h.

Referenced by computeQpProperties().

◆ _enrich_variable

std::vector<std::vector<MooseVariableFEBase *> > ComputeCrackTipEnrichmentSmallStrain::_enrich_variable
protected

enrichment displacement variables

Definition at line 49 of file ComputeCrackTipEnrichmentSmallStrain.h.

Referenced by ComputeCrackTipEnrichmentSmallStrain(), and computeQpProperties().

◆ _fe_dphi

const std::vector<std::vector<RealGradient> >* ComputeCrackTipEnrichmentSmallStrain::_fe_dphi
private

gradient of shape function

Definition at line 69 of file ComputeCrackTipEnrichmentSmallStrain.h.

Referenced by computeProperties().

◆ _fe_phi

const std::vector<std::vector<Real> >* ComputeCrackTipEnrichmentSmallStrain::_fe_phi
private

shape function

Definition at line 67 of file ComputeCrackTipEnrichmentSmallStrain.h.

Referenced by computeProperties(), and computeQpProperties().

◆ _global_strain

const MaterialProperty<RankTwoTensor>* ComputeStrainBase::_global_strain
protectedinherited

◆ _grad_disp

std::vector<const VariableGradient *> ComputeStrainBase::_grad_disp
protectedinherited

◆ _grad_enrich_disp

std::vector<RealVectorValue> ComputeCrackTipEnrichmentSmallStrain::_grad_enrich_disp
protected

gradient of enrichment displacement

Definition at line 46 of file ComputeCrackTipEnrichmentSmallStrain.h.

Referenced by computeQpProperties().

◆ _grad_phi

const VariablePhiGradient& ComputeCrackTipEnrichmentSmallStrain::_grad_phi
protected

gradient of the shape function

Definition at line 55 of file ComputeCrackTipEnrichmentSmallStrain.h.

◆ _mechanical_strain

MaterialProperty<RankTwoTensor>& ComputeStrainBase::_mechanical_strain
protectedinherited

◆ _ndisp

unsigned int ComputeStrainBase::_ndisp
protectedinherited

◆ _nl

NonlinearSystem* ComputeCrackTipEnrichmentSmallStrain::_nl
private

◆ _phi

const VariablePhiValue& ComputeCrackTipEnrichmentSmallStrain::_phi
protected

the current shape functions

Definition at line 52 of file ComputeCrackTipEnrichmentSmallStrain.h.

◆ _sln

const NumericVector<Number>* ComputeCrackTipEnrichmentSmallStrain::_sln
private

Definition at line 71 of file ComputeCrackTipEnrichmentSmallStrain.h.

Referenced by computeQpProperties().

◆ _total_strain

MaterialProperty<RankTwoTensor>& ComputeStrainBase::_total_strain
protectedinherited

◆ _volumetric_locking_correction

bool ComputeStrainBase::_volumetric_locking_correction
protectedinherited

The documentation for this class was generated from the following files: