12 #include "libmesh/fe_interface.h" 13 #include "libmesh/string_to_enum.h" 15 #include "libmesh/quadrature_gauss.h" 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");
28 "The CrackFrontDefinition user object name");
39 _phi(_assembly.phi()),
40 _grad_phi(_assembly.gradPhi()),
48 const std::vector<NonlinearVariableName> & nl_vnames =
49 getParam<std::vector<NonlinearVariableName>>(
"enrichment_displacements");
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.");
59 for (
unsigned int i = 0; i < 4; ++i)
67 for (
unsigned int i = 0; i <
_BI.size(); ++i)
75 unsigned int crack_front_point_index =
78 for (
unsigned int i = 0; i < 4; ++i)
83 for (
unsigned int m = 0; m <
_ndisp; ++m)
90 for (
unsigned int j = 0;
j < 4; ++
j)
93 Real soln = (*_sln)(dof);
105 RankTwoTensor enrich_strain = (grad_tensor_enrich + grad_tensor_enrich.transpose()) / 2.0;
124 FEType fe_type(Utility::string_to_enum<Order>(
"first"),
125 Utility::string_to_enum<FEFamily>(
"lagrange"));
127 std::unique_ptr<FEBase> fe(FEBase::build(
dim, fe_type));
128 fe->attach_quadrature_rule(const_cast<QBase *>(
_qrule));
137 for (
unsigned int i = 0; i <
_BI.size(); ++i)
const MooseArray< Point > & _q_point
std::vector< std::vector< MooseVariableFEBase * > > _enrich_variable
enrichment displacement variables
FEProblemBase & _fe_problem
const QBase *const & _qrule
std::vector< Real > _B
enrichment function value
virtual void computeQpProperties() override
virtual void computeProperties() override
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
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 ¶meters)
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()
const unsigned int & _current_side
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
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
const NumericVector< Number > * _sln
MaterialProperty< RankTwoTensor > & _total_strain
const Elem *const & _current_elem
std::vector< const MaterialProperty< RankTwoTensor > * > _eigenstrains
const std::vector< std::vector< RealGradient > > * _fe_dphi
gradient of shape function
std::vector< const VariableGradient * > _grad_disp
Gradient of displacements.