www.mooseframework.org
InteractionIntegral.h
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 
10 #pragma once
11 
13 
14 // Forward Declarations
16 
22 template <bool is_ad>
24 {
25 public:
27 
29 
30  virtual void initialSetup() override;
31  virtual void initialize() override;
32  virtual void execute() override;
33  virtual void finalize() override;
34  virtual void threadJoin(const UserObject & y) override;
35 
40  static MooseEnum qFunctionType();
41 
46  static MooseEnum sifModeType();
47 
48 protected:
58  Real computeQpIntegral(const std::size_t crack_front_point_index,
59  const Real scalar_q,
60  const RealVectorValue & grad_of_scalar_q);
61 
73  void computeAuxFields(RankTwoTensor & aux_stress,
74  RankTwoTensor & grad_disp,
75  RankTwoTensor & aux_strain,
76  RankTwoTensor & aux_disp);
84  void computeTFields(RankTwoTensor & aux_stress, RankTwoTensor & grad_disp);
85 
87  std::size_t _ndisp;
99  std::vector<MooseVariableFEBase *> _fe_vars;
101  const FEType & _fe_type;
103  std::vector<const VariableValue *> _disp;
105  std::vector<const VariableGradient *> _grad_disp;
107  const bool _has_temp;
123  const bool _fgm_crack;
125  std::size_t _ring_index;
143  std::vector<Real> _q_curr_elem;
145  const std::vector<std::vector<Real>> * _phi_curr_elem;
147  const std::vector<std::vector<RealGradient>> * _dphi_curr_elem;
157  unsigned int _qp;
158 
161  const enum class QMethod { Geometry, Topology } _q_function_type;
162 
165  const enum class PositionType { Angle, Distance } _position_type;
166 
168  const enum class SifMethod { KI, KII, KIII, T } _sif_mode;
169 
180 
186 };
187 
const VariableGradient * _additional_eigenstrain_gradient_02
Gradient used to add contribution to the interaction integral (XZ)
unsigned int _qp
Current quadrature point index.
const std::vector< std::vector< Real > > * _phi_curr_elem
Vector of shape function values for the current element.
OutputTools< Real >::VariableGradient VariableGradient
const GenericMaterialProperty< RankTwoTensor, is_ad > *const _total_deigenstrain_dT
Derivative of the total eigenstrain with respect to temperature.
const std::vector< std::vector< RealGradient > > * _dphi_curr_elem
Vector of gradients of shape function values for the current element.
VectorPostprocessorValue & _x
Vectors computed by this VectorPostprocessor: x,y,z coordinates, position of nodes along crack front...
enum InteractionIntegralTempl::SifMethod _sif_mode
bool _has_additional_eigenstrain
Whether the user chooses to add other eigenstrain influence (e.g. irradiation-induced) ...
std::vector< const VariableValue * > _disp
Displacement variables.
const CrackFrontDefinition *const _crack_front_definition
Pointer to the crack front definition object.
const MaterialProperty< Real > * _functionally_graded_youngs_modulus
Spatial elasticity modulus variable for FGM.
const MaterialProperty< Real > * _functionally_graded_youngs_modulus_crack_dir_gradient
Spatial derivative of the youngs modulus in the crack direction.
void computeAuxFields(RankTwoTensor &aux_stress, RankTwoTensor &grad_disp, RankTwoTensor &aux_strain, RankTwoTensor &aux_disp)
Compute the auxiliary fields, including the auxiliary stress and the gradient of the auxiliary disp...
virtual void initialize() override
VectorPostprocessorValue & _y
Real _kappa
Kappa Lame constant.
std::vector< const VariableGradient * > _grad_disp
Gradient of displacements.
const VariableGradient * _additional_eigenstrain_gradient_01
Gradient used to add contribution to the interaction integral (XY)
Real _r
Radial distance from the current point to the crack front.
Real _poissons_ratio
Poisson&#39;s ratio of the material.
const VariableGradient * _additional_eigenstrain_gradient_11
Gradient used to add contribution to the interaction integral (YY)
const std::vector< double > y
virtual void execute() override
const VariableGradient & _grad_temp
Gradient of temperature.
static InputParameters validParams()
QMethod
Enum used to select the method used to compute the q function used in the fracture integrals...
bool _treat_as_2d
Whether to treat a model as 2D for computation of fracture integrals.
InteractionIntegralTempl< false > InteractionIntegral
Class used in fracture integrals to define geometric characteristics of the crack front...
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 a fracture integral for a point on the cra...
bool _using_mesh_cutter
Is the crack defined by an XFEM cutter mesh.
static MooseEnum qFunctionType()
Get the MooseEnum defining the options for how the q function is defined.
InteractionIntegralTempl< true > ADInteractionIntegral
const GenericMaterialProperty< RankTwoTensor, is_ad > & _stress
Reference to the stress tensor computed by the material models.
const MaterialProperty< RealVectorValue > * _body_force
VectorPostprocessorValue & _z
typename GenericMaterialPropertyStruct< T, is_ad >::type GenericMaterialProperty
std::vector< MooseVariableFEBase * > _fe_vars
Vector of all coupled variables.
virtual void finalize() override
Real _theta
Angle of current point relative to the crack front.
bool _has_symmetry_plane
Whether the crack plane is also a symmetry plane in the model.
This vectorpostprocessor computes the Interaction Integral, which is used to compute various fracture...
InteractionIntegralTempl(const InputParameters &parameters)
Real _K_factor
Conversion factor applied to convert interaction integral to stress intensity factor K...
const bool _has_temp
Whether the temperature variable is coupled.
const bool _fgm_crack
Whether to consider interaction integral and material properties for a crack in functionally graded m...
std::vector< Real > VectorPostprocessorValue
static MooseEnum sifModeType()
Get the MooseEnum defining the options for the integral to be computed.
const VariableGradient * _additional_eigenstrain_gradient_12
Gradient used to add contribution to the interaction integral (YZ)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual void initialSetup() override
const InputParameters & parameters() const
VectorPostprocessorValue & _position
SifMethod
Enum used to select the type of integral to be performed.
const GenericMaterialProperty< RankThreeTensor, is_ad > * _eigenstrain_gradient
Pointers to optionally-used eigenstrain gradient and body force.
const FEType & _fe_type
FEType object defining order and family of displacement variables.
std::vector< Real > _q_curr_elem
Vector of q function values for the nodes in the current element.
void computeTFields(RankTwoTensor &aux_stress, RankTwoTensor &grad_disp)
Compute the auxiliary fields, including the auxiliary stress and the gradient of the auxiliary displa...
enum InteractionIntegralTempl::QMethod _q_function_type
Real _shear_modulus
Shear modulus.
virtual void threadJoin(const UserObject &y) override
enum InteractionIntegralTempl::PositionType _position_type
const GenericMaterialProperty< RankTwoTensor, is_ad > & _strain
Reference to the strain tensor computed by the material models.
std::size_t _ndisp
Number of displacement components.
const VariableGradient * _additional_eigenstrain_gradient_22
Gradient used to add contribution to the interaction integral (ZZ)
std::size_t _ring_index
Index of the ring for the integral computed by this object.
PositionType
Enum used to define how the distance along the crack front is measured (angle or distance) ...
VectorPostprocessorValue & _interaction_integral
Real _youngs_modulus
Young&#39;s modulus of the material.
const VariableGradient * _additional_eigenstrain_gradient_00
Gradient used to add contribution to the interaction integral (XX)