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
 
virtual void initQpStatefulProperties () override
 

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< MooseVariable * > > _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
 
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 22 of file ComputeCrackTipEnrichmentSmallStrain.h.

Constructor & Destructor Documentation

ComputeCrackTipEnrichmentSmallStrain::ComputeCrackTipEnrichmentSmallStrain ( const InputParameters &  parameters)

Definition at line 27 of file ComputeCrackTipEnrichmentSmallStrain.C.

29  : ComputeStrainBase(parameters),
30  EnrichmentFunctionCalculation(&getUserObject<CrackFrontDefinition>("crack_front_definition")),
31  _enrich_disp(3),
34  _phi(_assembly.phi()),
35  _grad_phi(_assembly.gradPhi()),
36  _B(4),
37  _dBX(4),
38  _dBx(4)
39 {
40  for (unsigned int i = 0; i < _enrich_variable.size(); ++i)
41  _enrich_variable[i].resize(_ndisp);
42 
43  const std::vector<NonlinearVariableName> & nl_vnames =
44  getParam<std::vector<NonlinearVariableName>>("enrichment_displacements");
45 
46  if (_ndisp == 2 && nl_vnames.size() != 8)
47  mooseError("The number of enrichment displacements should be total 8 for 2D.");
48  else if (_ndisp == 3 && nl_vnames.size() != 12)
49  mooseError("The number of enrichment displacements should be total 12 for 3D.");
50 
51  NonlinearSystem & nl = _fe_problem.getNonlinearSystem();
52  _nl = &(_fe_problem.getNonlinearSystem());
53 
54  for (unsigned int j = 0; j < _ndisp; ++j)
55  for (unsigned int i = 0; i < 4; ++i)
56  _enrich_variable[i][j] = &(nl.getVariable(0, nl_vnames[j * 4 + i]));
57 
58  if (_ndisp == 2)
59  _BI.resize(4); // QUAD4
60  else if (_ndisp == 3)
61  _BI.resize(8); // HEX8
62 
63  for (unsigned int i = 0; i < _BI.size(); ++i)
64  _BI[i].resize(4);
65 }
const VariablePhiValue & _phi
the current shape functions
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
std::vector< std::vector< MooseVariable * > > _enrich_variable
enrichment displacement variables
virtual ComputeCrackTipEnrichmentSmallStrain::~ComputeCrackTipEnrichmentSmallStrain ( )
inlinevirtual

Definition at line 27 of file ComputeCrackTipEnrichmentSmallStrain.h.

27 {}

Member Function Documentation

void ComputeCrackTipEnrichmentSmallStrain::computeProperties ( )
overrideprotectedvirtual

Definition at line 117 of file ComputeCrackTipEnrichmentSmallStrain.C.

Referenced by ~ComputeCrackTipEnrichmentSmallStrain().

118 {
119  FEType fe_type(Utility::string_to_enum<Order>("first"),
120  Utility::string_to_enum<FEFamily>("lagrange"));
121  const unsigned int dim = _current_elem->dim();
122  UniquePtr<FEBase> fe(FEBase::build(dim, fe_type));
123  fe->attach_quadrature_rule(_qrule);
124  _fe_phi = &(fe->get_phi());
125  _fe_dphi = &(fe->get_dphi());
126 
127  if (isBoundaryMaterial())
128  fe->reinit(_current_elem, _current_side);
129  else
130  fe->reinit(_current_elem);
131 
132  for (unsigned int i = 0; i < _BI.size(); ++i)
133  crackTipEnrichementFunctionAtPoint(*(_current_elem->get_node(i)), _BI[i]);
134 
135  ComputeStrainBase::computeProperties();
136 }
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
void ComputeCrackTipEnrichmentSmallStrain::computeQpProperties ( )
overrideprotectedvirtual

Definition at line 68 of file ComputeCrackTipEnrichmentSmallStrain.C.

Referenced by ~ComputeCrackTipEnrichmentSmallStrain().

69 {
70  crackTipEnrichementFunctionAtPoint(_q_point[_qp], _B);
71  unsigned int crack_front_point_index =
73 
74  for (unsigned int i = 0; i < 4; ++i)
75  rotateFromCrackFrontCoordsToGlobal(_dBx[i], _dBX[i], crack_front_point_index);
76 
77  _sln = _nl->currentSolution();
78 
79  for (unsigned int m = 0; m < _ndisp; ++m)
80  {
81  _enrich_disp[m] = 0.0;
82  _grad_enrich_disp[m].zero();
83  for (unsigned int i = 0; i < _current_elem->n_nodes(); ++i)
84  {
85  Node * node_i = _current_elem->get_node(i);
86  for (unsigned int j = 0; j < 4; ++j)
87  {
88  dof_id_type dof = node_i->dof_number(_nl->number(), _enrich_variable[j][m]->number(), 0);
89  Real soln = (*_sln)(dof);
90  _enrich_disp[m] += (*_fe_phi)[i][_qp] * (_B[j] - _BI[i][j]) * soln;
91  RealVectorValue grad_B(_dBX[j]);
92  _grad_enrich_disp[m] +=
93  ((*_fe_dphi)[i][_qp] * (_B[j] - _BI[i][j]) + (*_fe_phi)[i][_qp] * grad_B) * soln;
94  }
95  }
96  }
97 
98  RankTwoTensor grad_tensor_enrich(
100 
101  RankTwoTensor enrich_strain = (grad_tensor_enrich + grad_tensor_enrich.transpose()) / 2.0;
102 
103  RankTwoTensor grad_tensor((*_grad_disp[0])[_qp], (*_grad_disp[1])[_qp], (*_grad_disp[2])[_qp]);
104 
105  _total_strain[_qp] = (grad_tensor + grad_tensor.transpose()) / 2.0;
106 
107  _total_strain[_qp] += enrich_strain;
108 
109  _mechanical_strain[_qp] = _total_strain[_qp];
110 
111  // Remove the Eigen strain
112  for (auto es : _eigenstrains)
113  _mechanical_strain[_qp] -= (*es)[_qp];
114 }
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
std::vector< std::vector< MooseVariable * > > _enrich_variable
enrichment displacement variables
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 17 of file EnrichmentFunctionCalculation.C.

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

19 {
20  unsigned int crack_front_point_index =
22 
23  if (MooseUtils::absoluteFuzzyEqual(_r, 0.0))
24  mooseError("EnrichmentFunctionCalculation: the distance between a point and the crack "
25  "tip/front is zero.");
26 
27  Real st = std::sin(_theta);
28  Real st2 = std::sin(_theta / 2.0);
29  Real ct2 = std::cos(_theta / 2.0);
30  Real sr = std::sqrt(_r);
31 
32  B[0] = sr * st2;
33  B[1] = sr * ct2;
34  B[2] = sr * st2 * st;
35  B[3] = sr * ct2 * st;
36 
37  return crack_front_point_index;
38 }
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
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 41 of file EnrichmentFunctionCalculation.C.

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

43 {
44  unsigned int crack_front_point_index =
46 
47  if (MooseUtils::absoluteFuzzyEqual(_r, 0.0))
48  mooseError("EnrichmentFunctionCalculation: the distance between a point and the crack "
49  "tip/front is zero.");
50 
51  Real st = std::sin(_theta);
52  Real ct = std::cos(_theta);
53  Real st2 = std::sin(_theta / 2.0);
54  Real ct2 = std::cos(_theta / 2.0);
55  Real st15 = std::sin(1.5 * _theta);
56  Real ct15 = std::cos(1.5 * _theta);
57  Real sr = std::sqrt(_r);
58 
59  dB[0](0) = -0.5 / sr * st2;
60  dB[0](1) = 0.5 / sr * ct2;
61  dB[0](2) = 0.0;
62  dB[1](0) = 0.5 / sr * ct2;
63  dB[1](1) = 0.5 / sr * st2;
64  dB[1](2) = 0.0;
65  dB[2](0) = -0.5 / sr * st15 * st;
66  dB[2](1) = 0.5 / sr * (st2 + st15 * ct);
67  dB[2](2) = 0.0;
68  dB[3](0) = -0.5 / sr * ct15 * st;
69  dB[3](1) = 0.5 / sr * (ct2 + ct15 * ct);
70  dB[3](2) = 0.0;
71 
72  return crack_front_point_index;
73 }
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
void ComputeStrainBase::initQpStatefulProperties ( )
overrideprotectedvirtualinherited

Reimplemented in ComputeCosseratIncrementalSmallStrain, and ComputeIncrementalStrainBase.

Definition at line 78 of file ComputeStrainBase.C.

Referenced by ComputeStrainBase::~ComputeStrainBase().

79 {
80  _mechanical_strain[_qp].zero();
81  _total_strain[_qp].zero();
82 }
MaterialProperty< RankTwoTensor > & _mechanical_strain
MaterialProperty< RankTwoTensor > & _total_strain
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 76 of file EnrichmentFunctionCalculation.C.

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

79 {
80  rotated_vector = _crack_front_definition->rotateFromCrackFrontCoordsToGlobal(vector, point_index);
81 }
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

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

enrichment function value

Definition at line 51 of file ComputeCrackTipEnrichmentSmallStrain.h.

Referenced by computeQpProperties().

std::string ComputeStrainBase::_base_name
protectedinherited

Definition at line 33 of file ComputeStrainBase.h.

Referenced by ComputeStrainBase::ComputeStrainBase().

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

enrichment function at node I

Definition at line 57 of file ComputeCrackTipEnrichmentSmallStrain.h.

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

const Real& ComputeStrainBase::_current_elem_volume
protectedinherited
std::vector<RealVectorValue> ComputeCrackTipEnrichmentSmallStrain::_dBX
private

derivatives of enrichment function respect to global cooridnate

Definition at line 53 of file ComputeCrackTipEnrichmentSmallStrain.h.

Referenced by computeQpProperties().

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

derivatives of enrichment function respect to crack front cooridnate

Definition at line 55 of file ComputeCrackTipEnrichmentSmallStrain.h.

Referenced by computeQpProperties().

std::vector<const VariableValue *> ComputeStrainBase::_disp
protectedinherited
std::vector<MaterialPropertyName> ComputeStrainBase::_eigenstrain_names
protectedinherited
std::vector<const MaterialProperty<RankTwoTensor> *> ComputeStrainBase::_eigenstrains
protectedinherited
std::vector<Real> ComputeCrackTipEnrichmentSmallStrain::_enrich_disp
protected

enrichment displacement

Definition at line 35 of file ComputeCrackTipEnrichmentSmallStrain.h.

Referenced by computeQpProperties().

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

enrichment displacement variables

Definition at line 41 of file ComputeCrackTipEnrichmentSmallStrain.h.

Referenced by ComputeCrackTipEnrichmentSmallStrain(), and computeQpProperties().

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

gradient of shape function

Definition at line 61 of file ComputeCrackTipEnrichmentSmallStrain.h.

Referenced by computeProperties().

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

shape function

Definition at line 59 of file ComputeCrackTipEnrichmentSmallStrain.h.

Referenced by computeProperties(), and computeQpProperties().

std::vector<const VariableGradient *> ComputeStrainBase::_grad_disp
protectedinherited
std::vector<RealVectorValue> ComputeCrackTipEnrichmentSmallStrain::_grad_enrich_disp
protected

gradient of enrichment displacement

Definition at line 38 of file ComputeCrackTipEnrichmentSmallStrain.h.

Referenced by computeQpProperties().

const VariablePhiGradient& ComputeCrackTipEnrichmentSmallStrain::_grad_phi
protected

gradient of the shape function

Definition at line 47 of file ComputeCrackTipEnrichmentSmallStrain.h.

MaterialProperty<RankTwoTensor>& ComputeStrainBase::_mechanical_strain
protectedinherited
unsigned int ComputeStrainBase::_ndisp
protectedinherited
NonlinearSystem* ComputeCrackTipEnrichmentSmallStrain::_nl
private
const VariablePhiValue& ComputeCrackTipEnrichmentSmallStrain::_phi
protected

the current shape functions

Definition at line 44 of file ComputeCrackTipEnrichmentSmallStrain.h.

const NumericVector<Number>* ComputeCrackTipEnrichmentSmallStrain::_sln
private

Definition at line 63 of file ComputeCrackTipEnrichmentSmallStrain.h.

Referenced by computeQpProperties().

MaterialProperty<RankTwoTensor>& ComputeStrainBase::_total_strain
protectedinherited
bool ComputeStrainBase::_volumetric_locking_correction
protectedinherited

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