www.mooseframework.org
ComputeCrackTipEnrichmentSmallStrain.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
3 /* */
4 /* All contents are licensed under LGPL V2.1 */
5 /* See LICENSE for full restrictions */
6 /****************************************************************/
7 
9 #include "MooseMesh.h"
10 #include "libmesh/fe_interface.h"
11 #include "libmesh/string_to_enum.h"
12 
13 #include "libmesh/quadrature_gauss.h"
14 
15 template <>
16 InputParameters
18 {
19  InputParameters params = validParams<ComputeStrainBase>();
20  params.addRequiredParam<std::vector<NonlinearVariableName>>("enrichment_displacements",
21  "The enrichment displacement");
22  params.addRequiredParam<UserObjectName>("crack_front_definition",
23  "The CrackFrontDefinition user object name");
24  return params;
25 }
26 
28  const InputParameters & parameters)
29  : ComputeStrainBase(parameters),
30  EnrichmentFunctionCalculation(&getUserObject<CrackFrontDefinition>("crack_front_definition")),
31  _enrich_disp(3),
32  _grad_enrich_disp(3),
33  _enrich_variable(4),
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 }
66 
67 void
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 }
115 
116 void
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 }
InputParameters validParams< ComputeCrackTipEnrichmentSmallStrain >()
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
Works on top of NodalNormalsPreprocessor.
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.
Perform calculation of enrichment function values and derivatives.
std::vector< RealVectorValue > _dBx
derivatives of enrichment function respect to crack front cooridnate
ComputeCrackTipEnrichmentSmallStrain(const InputParameters &parameters)
virtual unsigned int crackTipEnrichementFunctionDerivativeAtPoint(const Point &point, std::vector< RealVectorValue > &dB)
calculate the enrichment function derivatives at point
ComputeStrainBase is the base class for strain tensors.
MaterialProperty< RankTwoTensor > & _total_strain
std::vector< const MaterialProperty< RankTwoTensor > * > _eigenstrains
InputParameters validParams< ComputeStrainBase >()
const std::vector< std::vector< RealGradient > > * _fe_dphi
gradient of shape function
std::vector< const VariableGradient * > _grad_disp
std::vector< std::vector< MooseVariable * > > _enrich_variable
enrichment displacement variables