13 #include "libmesh/string_to_enum.h" 14 #include "libmesh/quadrature.h" 15 #include "libmesh/utility.h" 25 "The displacements appropriate for the simulation geometry and coordinate system");
27 "The CrackFrontDefinition user object name");
28 MooseEnum position_type(
"Angle Distance",
"Distance");
32 "The method used to calculate position along crack front. Options are: " +
34 MooseEnum integral_vec(
"JIntegral CIntegral KFromJIntegral",
"JIntegral");
37 "Domain integrals to calculate. Choices are: " +
39 params.
addParam<
unsigned int>(
"symmetry_plane",
40 "Account for a symmetry plane passing through " 41 "the plane of the crack, normal to the specified " 42 "axis (0=x, 1=y, 2=z)");
43 params.
addParam<
Real>(
"poissons_ratio",
"Poisson's ratio");
44 params.
addParam<
Real>(
"youngs_modulus",
"Young's modulus of the material.");
45 params.
set<
bool>(
"use_displaced_mesh") =
false;
47 MooseEnum q_function_type(
"Geometry Topology",
"Geometry");
50 "The method used to define the integration domain. Options are: " +
52 params.
addClassDescription(
"Computes the J-Integral, a measure of the strain energy release rate " 53 "at a crack tip, which can be used as a criterion for fracture " 54 "growth. It can, alternatively, compute the C(t) integral");
62 _J_thermal_term_vec(hasMaterialProperty<
RealVectorValue>(
"J_thermal_term_vec")
68 _Eshelby_tensor_dissipation(
70 ? &getMaterialProperty<
RankTwoTensor>(
"Eshelby_tensor_dissipation")
72 _fe_vars(getCoupledMooseVars()),
73 _fe_type(_fe_vars[0]->feType()),
74 _has_symmetry_plane(isParamValid(
"symmetry_plane")),
75 _poissons_ratio(isParamValid(
"poissons_ratio") ? getParam<
Real>(
"poissons_ratio") : 0),
76 _youngs_modulus(isParamValid(
"youngs_modulus") ? getParam<
Real>(
"youngs_modulus") : 0),
77 _ring_index(getParam<unsigned
int>(
"ring_index")),
80 _x(declareVector(
"x")),
81 _y(declareVector(
"y")),
82 _z(declareVector(
"z")),
83 _position(declareVector(
"id")),
84 _j_integral(declareVector((_integral ==
IntegralType::KFromJIntegral
98 mooseError(
"youngs_modulus and poissons_ratio must be specified if integrals = KFromJIntegral");
110 _x.assign(num_pts, 0.0);
111 _y.assign(num_pts, 0.0);
112 _z.assign(num_pts, 0.0);
116 for (
const auto * fe_var :
_fe_vars)
119 mooseError(
"displacements",
"All coupled displacements must have the same type");
131 grad_of_vector_q(0, 0) = crack_direction(0) * grad_of_scalar_q(0);
132 grad_of_vector_q(0, 1) = crack_direction(0) * grad_of_scalar_q(1);
133 grad_of_vector_q(0, 2) = crack_direction(0) * grad_of_scalar_q(2);
134 grad_of_vector_q(1, 0) = crack_direction(1) * grad_of_scalar_q(0);
135 grad_of_vector_q(1, 1) = crack_direction(1) * grad_of_scalar_q(1);
136 grad_of_vector_q(1, 2) = crack_direction(1) * grad_of_scalar_q(2);
137 grad_of_vector_q(2, 0) = crack_direction(2) * grad_of_scalar_q(0);
138 grad_of_vector_q(2, 1) = crack_direction(2) * grad_of_scalar_q(1);
139 grad_of_vector_q(2, 2) = crack_direction(2) * grad_of_scalar_q(2);
144 eq = (*_Eshelby_tensor)[
_qp].doubleContraction(grad_of_vector_q);
146 eq = (*_Eshelby_tensor_dissipation)[
_qp].doubleContraction(grad_of_vector_q);
149 Real eq_thermal = 0.0;
152 for (std::size_t i = 0; i < 3; i++)
153 eq_thermal += crack_direction(i) * scalar_q * (*_J_thermal_term_vec)[
_qp](i);
156 Real q_avg_seg = 1.0;
165 Real etot = -eq + eq_thermal;
167 return etot / q_avg_seg;
176 fe->attach_quadrature_rule(const_cast<QBase *>(
_qrule));
184 for (std::size_t icfp = 0; icfp <
_j_integral.size(); icfp++)
226 for (std::size_t i = 0; i <
_j_integral.size(); ++i)
252 const auto & uo =
static_cast<const JIntegral &
>(
y);
std::unique_ptr< FEGenericBase< Real > > build(const unsigned int dim, const FEType &fet)
IntegralType
Enum defining the type of integral to be computed.
const libMesh::FEType & _fe_type
FEType object defining order and family of displacement variables.
enum JIntegral::PositionType _position_type
VectorPostprocessorValue & _j_integral
const MooseArray< Real > & _coord
Real _poissons_ratio
Poisson's ratio of the material.
static InputParameters validParams()
Real getCrackFrontBackwardSegmentLength(const std::size_t point_index) const
Get the length of the line segment on the crack front behind the specified position.
Real computeQpIntegral(const std::size_t crack_front_point_index, const Real scalar_q, const RealVectorValue &grad_of_scalar_q)
Compute the contribution of a specific quadrature point to the J integral for a point on the crack fr...
std::vector< MooseVariableFEBase * > _fe_vars
Vector of all coupled variables.
virtual void threadJoin(const UserObject &y) override
JIntegral(const InputParameters ¶meters)
Real _youngs_modulus
Young's modulus of the material.
This vectorpostprocessor computes the J-Integral, which is a measure of the strain energy release rat...
const std::vector< std::vector< libMesh::RealGradient > > * _dphi_curr_elem
Vector of gradients of shape function values for the current element.
VectorPostprocessorValue & _position
const std::vector< double > y
Real DomainIntegralQFunction(std::size_t crack_front_point_index, std::size_t ring_index, const Node *const current_node) const
Compute the q function for the case where it is defined geometrically.
std::string getRawNames() const
Real DomainIntegralTopologicalQFunction(std::size_t crack_front_point_index, std::size_t ring_index, const Node *const current_node) const
Compute the q function for the case where it is defined through element connectivity.
VectorPostprocessorValue & _y
const MaterialProperty< RealVectorValue > *const _J_thermal_term_vec
bool isParamValid(const std::string &name) const
Class used in fracture integrals to define geometric characteristics of the crack front...
std::size_t getNumCrackFrontPoints() const
Get the number of points defining the crack front as a set of line segments.
const std::vector< std::vector< Real > > * _phi_curr_elem
Vector of shape function values for the current element.
QMethod
Enum used to select the method used to compute the q function used in the fracture integrals...
const Point * getCrackFrontPoint(const std::size_t point_index) const
Get a Point object for a specified point on the crack front.
const CrackFrontDefinition *const _crack_front_definition
std::size_t _ring_index
Index of the ring for the integral computed by this object.
virtual void execute() override
Real getDistanceAlongFront(const std::size_t point_index) const
Get the distance along the crack front from the beginning of the crack to the specified position...
const RealVectorValue & getCrackDirection(const std::size_t point_index) const
Get the unit vector of the crack extension direction at the specified position.
enum JIntegral::IntegralType _integral
std::string stringify(const T &t)
std::vector< Real > _q_curr_elem
Vector of q function values for the nodes in the current element.
PositionType
Enum used to define how the distance along the crack front is measured (angle or distance) ...
registerMooseObject("SolidMechanicsApp", JIntegral)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const QBase *const & _qrule
const Elem *const & _current_elem
const MooseArray< Real > & _JxW
virtual void initialize() override
Real getAngleAlongFront(const std::size_t point_index) const
Get the angle along the crack front from the beginning of the crack to the specified position...
unsigned int _qp
Current quadrature point index.
bool _has_symmetry_plane
Whether the crack plane is also a symmetry plane in the model.
void mooseError(Args &&... args) const
Real getCrackFrontForwardSegmentLength(const std::size_t point_index) const
Get the length of the line segment on the crack front ahead of the specified position.
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
static InputParameters validParams()
VectorPostprocessorValue & _x
Vectors computed by this VectorPostprocessor: x,y,z coordinates, position of nodes along crack front...
bool usingMeshCutter() const
Is the crack defined by a mesh cutter object.
virtual void finalize() override
bool _treat_as_2d
Whether to treat a 3D model as 2D for computation of fracture integrals.
bool treatAs2D() const
Whether the fracture computations are treated as 2D for the model.
void ErrorVector unsigned int
bool _using_mesh_cutter
Whether the crack is defined by an XFEM cutter mesh.
virtual void initialSetup() override
enum JIntegral::QMethod _q_function_type
VectorPostprocessorValue & _z