www.mooseframework.org
ComputeCrackTipEnrichmentSmallStrain.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
11 #include "MooseMesh.h"
12 #include "libmesh/fe_interface.h"
13 #include "libmesh/string_to_enum.h"
14 
15 #include "libmesh/quadrature_gauss.h"
16 
18 
21 {
23  params.addClassDescription(
24  "Computes the crack tip enrichment at a point within a small strain formulation.");
25  params.addRequiredParam<std::vector<NonlinearVariableName>>("enrichment_displacements",
26  "The enrichment displacement");
27  params.addRequiredParam<UserObjectName>("crack_front_definition",
28  "The CrackFrontDefinition user object name");
29  return params;
30 }
31 
33  const InputParameters & parameters)
34  : ComputeStrainBase(parameters),
35  EnrichmentFunctionCalculation(&getUserObject<CrackFrontDefinition>("crack_front_definition")),
36  _enrich_disp(3),
37  _grad_enrich_disp(3),
38  _enrich_variable(4),
39  _phi(_assembly.phi()),
40  _grad_phi(_assembly.gradPhi()),
41  _B(4),
42  _dBX(4),
43  _dBx(4)
44 {
45  for (unsigned int i = 0; i < _enrich_variable.size(); ++i)
46  _enrich_variable[i].resize(_ndisp);
47 
48  const std::vector<NonlinearVariableName> & nl_vnames =
49  getParam<std::vector<NonlinearVariableName>>("enrichment_displacements");
50 
51  if (_ndisp == 2 && nl_vnames.size() != 8)
52  mooseError("The number of enrichment displacements should be total 8 for 2D.");
53  else if (_ndisp == 3 && nl_vnames.size() != 12)
54  mooseError("The number of enrichment displacements should be total 12 for 3D.");
55 
56  _nl = &(_fe_problem.getNonlinearSystem(/*nl_sys_num=*/0));
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 }
70 
71 void
73 {
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 
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  const Node * node_i = _current_elem->node_ptr(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  auto grad_tensor_enrich = RankTwoTensor::initializeFromRows(
104 
105  RankTwoTensor enrich_strain = (grad_tensor_enrich + grad_tensor_enrich.transpose()) / 2.0;
106 
107  auto grad_tensor = RankTwoTensor::initializeFromRows(
108  (*_grad_disp[0])[_qp], (*_grad_disp[1])[_qp], (*_grad_disp[2])[_qp]);
109 
110  _total_strain[_qp] = (grad_tensor + grad_tensor.transpose()) / 2.0;
111 
112  _total_strain[_qp] += enrich_strain;
113 
115 
116  // Remove the Eigen strain
117  for (auto es : _eigenstrains)
118  _mechanical_strain[_qp] -= (*es)[_qp];
119 }
120 
121 void
123 {
124  FEType fe_type(Utility::string_to_enum<Order>("first"),
125  Utility::string_to_enum<FEFamily>("lagrange"));
126  const unsigned int dim = _current_elem->dim();
127  std::unique_ptr<FEBase> fe(FEBase::build(dim, fe_type));
128  fe->attach_quadrature_rule(const_cast<QBase *>(_qrule));
129  _fe_phi = &(fe->get_phi());
130  _fe_dphi = &(fe->get_dphi());
131 
132  if (isBoundaryMaterial())
133  fe->reinit(_current_elem, _current_side);
134  else
135  fe->reinit(_current_elem);
136 
137  for (unsigned int i = 0; i < _BI.size(); ++i)
139 
141 }
const MooseArray< Point > & _q_point
std::vector< std::vector< MooseVariableFEBase * > > _enrich_variable
enrichment displacement variables
std::vector< Real > _B
enrichment function value
unsigned int dim
std::vector< RealVectorValue > _grad_enrich_disp
gradient of enrichment displacement
virtual void computeProperties() override
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
virtual bool isBoundaryMaterial() const override
const std::vector< std::vector< Real > > * _fe_phi
shape function
void addRequiredParam(const std::string &name, const std::string &doc_string)
static RankTwoTensorTempl initializeFromRows(const libMesh::TypeVector< Real > &row0, const libMesh::TypeVector< Real > &row1, const libMesh::TypeVector< Real > &row2)
registerMooseObject("XFEMApp", ComputeCrackTipEnrichmentSmallStrain)
Class used in fracture integrals to define geometric characteristics of the crack front...
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.
unsigned int number() const
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)
const NumericVector< Number > *const & currentSolution() const override
virtual NonlinearSystem & getNonlinearSystem(const unsigned int sys_num)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static InputParameters validParams()
void mooseError(Args &&... args) const
virtual unsigned int crackTipEnrichementFunctionDerivativeAtPoint(const Point &point, std::vector< RealVectorValue > &dB)
calculate the enrichment function derivatives at point
void addClassDescription(const std::string &doc_string)
ComputeStrainBase is the base class for strain tensors.
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
ComputeCrackTipEnrichmentSmallStrain calculates the sum of standard strain and enrichement strain...
MooseVariableFieldBase & getVariable(THREAD_ID tid, const std::string &var_name) const
MaterialProperty< RankTwoTensor > & _total_strain
std::vector< const MaterialProperty< RankTwoTensor > * > _eigenstrains
uint8_t dof_id_type
const std::vector< std::vector< RealGradient > > * _fe_dphi
gradient of shape function
std::vector< const VariableGradient * > _grad_disp
Gradient of displacements.