www.mooseframework.org
Public Types | Public Member Functions | Static Public Member Functions | Public Attributes | Static Public Attributes | Protected Types | Protected Member Functions | Static Protected Member Functions | Protected Attributes | List of all members
DomainIntegralAction Class Reference

Action to set up all objects used in computation of fracture domain integrals. More...

#include <DomainIntegralAction.h>

Inheritance diagram for DomainIntegralAction:
[legend]

Public Types

typedef DataFileName DataFileParameterType
 

Public Member Functions

 DomainIntegralAction (const InputParameters &params)
 
 ~DomainIntegralAction ()
 
virtual void act () override
 
virtual void addRelationshipManagers (Moose::RelationshipManagerType input_rm_type) override
 
virtual void addRelationshipManagers (Moose::RelationshipManagerType when_type)
 
bool addRelationshipManagers (Moose::RelationshipManagerType when_type, const InputParameters &moose_object_pars)
 
void timedAct ()
 
MooseObjectName uniqueActionName () const
 
const std::string & specificTaskName () const
 
const std::set< std::string > & getAllTasks () const
 
void appendTask (const std::string &task)
 
MooseAppgetMooseApp () const
 
const std::string & type () const
 
virtual const std::string & name () const
 
std::string typeAndName () const
 
std::string errorPrefix (const std::string &error_type) const
 
void callMooseError (std::string msg, const bool with_prefix) const
 
MooseObjectParameterName uniqueParameterName (const std::string &parameter_name) const
 
const InputParametersparameters () const
 
MooseObjectName uniqueName () const
 
const T & getParam (const std::string &name) const
 
std::vector< std::pair< T1, T2 > > getParam (const std::string &param1, const std::string &param2) const
 
const T & getRenamedParam (const std::string &old_name, const std::string &new_name) const
 
getCheckedPointerParam (const std::string &name, const std::string &error_string="") const
 
bool isParamValid (const std::string &name) const
 
bool isParamSetByUser (const std::string &nm) const
 
void paramError (const std::string &param, Args... args) const
 
void paramWarning (const std::string &param, Args... args) const
 
void paramInfo (const std::string &param, Args... args) const
 
void connectControllableParams (const std::string &parameter, const std::string &object_type, const std::string &object_name, const std::string &object_parameter) const
 
void mooseError (Args &&... args) const
 
void mooseErrorNonPrefixed (Args &&... args) const
 
void mooseWarning (Args &&... args) const
 
void mooseWarningNonPrefixed (Args &&... args) const
 
void mooseDeprecated (Args &&... args) const
 
void mooseInfo (Args &&... args) const
 
PerfGraphperfGraph ()
 
std::string getDataFileName (const std::string &param) const
 
std::string getDataFileNameByName (const std::string &name, const std::string *param=nullptr) const
 
const Parallel::Communicator & comm () const
 
processor_id_type n_processors () const
 
processor_id_type processor_id () const
 

Static Public Member Functions

static InputParameters validParams ()
 

Public Attributes

const ConsoleStream _console
 

Static Public Attributes

static constexpr auto SYSTEM
 
static constexpr auto NAME
 

Protected Types

enum  INTEGRAL {
  J_INTEGRAL, C_INTEGRAL, K_FROM_J_INTEGRAL, INTERACTION_INTEGRAL_KI,
  INTERACTION_INTEGRAL_KII, INTERACTION_INTEGRAL_KIII, INTERACTION_INTEGRAL_T
}
 Enum used to select the type of integral to be performed. More...
 
enum  Q_FUNCTION_TYPE { GEOMETRY, TOPOLOGY }
 Enum used to select the method used to compute the Q function used in the fracture integrals. More...
 

Protected Member Functions

unsigned int calcNumCrackFrontPoints ()
 Compute the number of points on the crack front. More...
 
bool addRelationshipManagers (Moose::RelationshipManagerType when_type, const InputParameters &moose_object_pars)
 
void associateWithParameter (const std::string &param_name, InputParameters &params) const
 
void associateWithParameter (const InputParameters &from_params, const std::string &param_name, InputParameters &params) const
 
const T & getMeshProperty (const std::string &data_name, const std::string &prefix)
 
const T & getMeshProperty (const std::string &data_name)
 
bool hasMeshProperty (const std::string &data_name, const std::string &prefix) const
 
bool hasMeshProperty (const std::string &data_name, const std::string &prefix) const
 
bool hasMeshProperty (const std::string &data_name) const
 
bool hasMeshProperty (const std::string &data_name) const
 
std::string meshPropertyName (const std::string &data_name) const
 
PerfID registerTimedSection (const std::string &section_name, const unsigned int level) const
 
PerfID registerTimedSection (const std::string &section_name, const unsigned int level, const std::string &live_message, const bool print_dots=true) const
 
std::string timedSectionName (const std::string &section_name) const
 

Static Protected Member Functions

static std::string meshPropertyName (const std::string &data_name, const std::string &prefix)
 

Protected Attributes

std::set< INTEGRAL_integrals
 Container for enumerations describing the individual integrals computed. More...
 
const std::vector< BoundaryName > & _boundary_names
 Boundaries containing the crack front points. More...
 
std::vector< Point > _crack_front_points
 User-defined vector of crack front points. More...
 
bool _closed_loop
 Indicates whether the crack forms a closed loop. More...
 
UserObjectName _crack_front_points_provider
 Name of crack front points provider user object used to optionally define the crack points. More...
 
bool _use_crack_front_points_provider
 Whether to use a crack front points provider. More...
 
MooseEnum _direction_method_moose_enum
 Enum used to define the method to compute crack front direction. More...
 
MooseEnum _end_direction_method_moose_enum
 Enum used to define the method to compute crack front direction at ends of crack front. More...
 
const bool _have_crack_direction_vector
 Whether the crack direction vector has been provided. More...
 
const RealVectorValue _crack_direction_vector
 Vector optionally used to prescribe direction of crack extension. More...
 
const bool _have_crack_direction_vector_end_1
 Whether the crack direction vector at the 1st end of the crack has been provided. More...
 
const RealVectorValue _crack_direction_vector_end_1
 Vector optionally used to prescribe direction of crack extension at the 1st end of the crack. More...
 
const bool _have_crack_direction_vector_end_2
 Whether the crack direction vector at the 2nd end of the crack has been provided. More...
 
const RealVectorValue _crack_direction_vector_end_2
 Vector optionally used to prescribe direction of crack extension at the 2nd end of the crack. More...
 
std::vector< BoundaryName > _crack_mouth_boundary_names
 Names of boundaries optionally used to define the crack mouth location. More...
 
std::vector< BoundaryName > _intersecting_boundary_names
 Names of boundaries optionally used for improved computation of crack extension direction at ends of the crack where it intersects the prescribed boundaries. More...
 
bool _treat_as_2d
 Whether fracture computations for a 3D model should be treated as though it were a 2D model. More...
 
unsigned int _axis_2d
 Out-of-plane axis for 3D models treated as 2D. More...
 
unsigned int _ring_first
 Number of elements away from the crack tip to inside of inner ring with the topological q function. More...
 
unsigned int _ring_last
 Number of elements away from the crack tip to outside of outer ring with the topological q function. More...
 
std::vector< VariableName > _output_variables
 List of variables for which values are to be sampled and output at the crack front points. More...
 
Real _poissons_ratio
 Poisson's ratio of material. More...
 
Real _youngs_modulus
 Young's modulus of material. More...
 
std::vector< SubdomainName > _blocks
 Blocks for which the domain integrals are to be computed. More...
 
std::vector< VariableName > _displacements
 Vector of displacement variables. More...
 
VariableName _temp
 Temperature variable. More...
 
bool _has_symmetry_plane
 Whether the model has a symmetry plane passing through the plane of the crack. More...
 
unsigned int _symmetry_plane
 Identifier for which plane is the symmetry plane. More...
 
MooseEnum _position_type
 How the distance along the crack front is measured (angle or distance) More...
 
MooseEnum _q_function_type
 How the q function is evaluated (geometric distance from crack front or ring of elements) More...
 
bool _get_equivalent_k
 Whether to compute the equivalent K from the individual fracture integrals for mixed-mode fracture. More...
 
bool _use_displaced_mesh
 Whether to compute the fracture integrals on the displaced mesh. More...
 
bool _output_q
 Whether to ouput the q function as a set of AuxVariables. More...
 
std::vector< unsigned int_ring_vec
 Vector of ids for the individual rings on which the fracture integral is computed. More...
 
bool _incremental
 Whether the constitutive models for the mechanics calculations use an incremental form. More...
 
bool _convert_J_to_K
 Whether to convert the J-integral to a stress intensity factor (K) –deprecated. More...
 
bool _fgm_crack
 Whether the crack lives in a functionally-graded material. More...
 
MaterialPropertyName _functionally_graded_youngs_modulus_crack_dir_gradient
 Material property name for the Youngs modulus derivative for functionally graded materials. More...
 
MaterialPropertyName _functionally_graded_youngs_modulus
 Material property name for spatially-dependent Youngs modulus for functionally graded materials. More...
 
const bool _use_ad
 Whether to create automatic differentiation objects from the action. More...
 
const bool _used_by_xfem_to_grow_crack
 This determines if fracture integrals should be executed on nonlinear in order to grow the crack when num_xfem_updates in the executioner block is greater than 1. More...
 
std::string _registered_identifier
 
std::string _specific_task_name
 
std::set< std::string > _all_tasks
 
ActionWarehouse_awh
 
const std::string & _current_task
 
std::shared_ptr< MooseMesh > & _mesh
 
std::shared_ptr< MooseMesh > & _displaced_mesh
 
std::shared_ptr< FEProblemBase > & _problem
 
PerfID _act_timer
 
MooseApp_app
 
const std::string _type
 
const std::string _name
 
const InputParameters_pars
 
Factory_factory
 
ActionFactory_action_factory
 
MooseApp_pg_moose_app
 
const std::string _prefix
 
const Parallel::Communicator & _communicator
 
const std::string _order
 Order and family of the AuxVariables optionally created to output the values of q. More...
 
const std::string _family
 
std::vector< Real_radius_inner
 Sets of inner and outer radii of the rings used for the domain form of the fracture integrals. More...
 
std::vector< Real_radius_outer
 

Detailed Description

Action to set up all objects used in computation of fracture domain integrals.

Definition at line 22 of file DomainIntegralAction.h.

Member Enumeration Documentation

◆ INTEGRAL

◆ Q_FUNCTION_TYPE

Enum used to select the method used to compute the Q function used in the fracture integrals.

Enumerator
GEOMETRY 
TOPOLOGY 

Definition at line 51 of file DomainIntegralAction.h.

Constructor & Destructor Documentation

◆ DomainIntegralAction()

DomainIntegralAction::DomainIntegralAction ( const InputParameters params)

Definition at line 148 of file DomainIntegralAction.C.

149  : Action(params),
150  _boundary_names(getParam<std::vector<BoundaryName>>("boundary")),
151  _closed_loop(getParam<bool>("closed_loop")),
153  _order(getParam<std::string>("order")),
154  _family(getParam<std::string>("family")),
155  _direction_method_moose_enum(getParam<MooseEnum>("crack_direction_method")),
156  _end_direction_method_moose_enum(getParam<MooseEnum>("crack_end_direction_method")),
157  _have_crack_direction_vector(isParamValid("crack_direction_vector")),
159  _have_crack_direction_vector ? getParam<RealVectorValue>("crack_direction_vector") : 0.0),
160  _have_crack_direction_vector_end_1(isParamValid("crack_direction_vector_end_1")),
162  ? getParam<RealVectorValue>("crack_direction_vector_end_1")
163  : 0.0),
164  _have_crack_direction_vector_end_2(isParamValid("crack_direction_vector_end_2")),
166  ? getParam<RealVectorValue>("crack_direction_vector_end_2")
167  : 0.0),
168  _treat_as_2d(getParam<bool>("2d")),
169  _axis_2d(getParam<unsigned int>("axis_2d")),
170  _has_symmetry_plane(isParamValid("symmetry_plane")),
171  _symmetry_plane(_has_symmetry_plane ? getParam<unsigned int>("symmetry_plane")
172  : std::numeric_limits<unsigned int>::max()),
173  _position_type(getParam<MooseEnum>("position_type")),
174  _q_function_type(getParam<MooseEnum>("q_function_type")),
175  _get_equivalent_k(getParam<bool>("equivalent_k")),
176  _use_displaced_mesh(false),
177  _output_q(getParam<bool>("output_q")),
178  _incremental(getParam<bool>("incremental")),
179  _convert_J_to_K(isParamValid("convert_J_to_K") ? getParam<bool>("convert_J_to_K") : false),
180  _fgm_crack(false),
181  _use_ad(getParam<bool>("use_automatic_differentiation")),
182  _used_by_xfem_to_grow_crack(getParam<bool>("used_by_xfem_to_grow_crack"))
183 {
184 
185  if (isParamValid("functionally_graded_youngs_modulus_crack_dir_gradient") !=
186  isParamValid("functionally_graded_youngs_modulus"))
187  paramError("functionally_graded_youngs_modulus_crack_dir_gradient",
188  "You have selected to compute the interaction integral for a crack in FGM. That "
189  "selection requires the user to provide a spatially varying elasticity modulus that "
190  "defines the transition of material properties (i.e. "
191  "'functionally_graded_youngs_modulus') and its "
192  "spatial derivative in the crack direction (i.e. "
193  "'functionally_graded_youngs_modulus_crack_dir_gradient').");
194 
195  if (isParamValid("functionally_graded_youngs_modulus_crack_dir_gradient") &&
196  isParamValid("functionally_graded_youngs_modulus"))
197  {
198  _fgm_crack = true;
200  getParam<MaterialPropertyName>("functionally_graded_youngs_modulus_crack_dir_gradient");
202  getParam<MaterialPropertyName>("functionally_graded_youngs_modulus");
203  }
204 
205  if (_q_function_type == GEOMETRY)
206  {
207  if (isParamValid("radius_inner") && isParamValid("radius_outer"))
208  {
209  _radius_inner = getParam<std::vector<Real>>("radius_inner");
210  _radius_outer = getParam<std::vector<Real>>("radius_outer");
211  }
212  else
213  mooseError("DomainIntegral error: must set radius_inner and radius_outer.");
214  for (unsigned int i = 0; i < _radius_inner.size(); ++i)
215  _ring_vec.push_back(i + 1);
216  }
217  else if (_q_function_type == TOPOLOGY)
218  {
219  if (isParamValid("ring_first") && isParamValid("ring_last"))
220  {
221  _ring_first = getParam<unsigned int>("ring_first");
222  _ring_last = getParam<unsigned int>("ring_last");
223  }
224  else
225  mooseError(
226  "DomainIntegral error: must set ring_first and ring_last if q_function_type = Topology.");
227  for (unsigned int i = _ring_first; i <= _ring_last; ++i)
228  _ring_vec.push_back(i);
229  }
230  else
231  paramError("q_function_type", "DomainIntegral error: invalid q_function_type.");
232 
233  if (isParamValid("crack_front_points"))
234  _crack_front_points = getParam<std::vector<Point>>("crack_front_points");
235 
236  if (isParamValid("crack_front_points_provider"))
237  {
238  if (!isParamValid("number_points_from_provider"))
239  paramError("number_points_from_provider",
240  "DomainIntegral error: when crack_front_points_provider is used, "
241  "number_points_from_provider must be provided.");
243  _crack_front_points_provider = getParam<UserObjectName>("crack_front_points_provider");
244  }
245  else if (isParamValid("number_points_from_provider"))
246  paramError("crack_front_points_provider",
247  "DomainIntegral error: number_points_from_provider is provided but "
248  "crack_front_points_provider cannot be found.");
249  if (isParamValid("crack_mouth_boundary"))
250  _crack_mouth_boundary_names = getParam<std::vector<BoundaryName>>("crack_mouth_boundary");
251  if (isParamValid("intersecting_boundary"))
252  _intersecting_boundary_names = getParam<std::vector<BoundaryName>>("intersecting_boundary");
253  if (_radius_inner.size() != _radius_outer.size())
254  mooseError("Number of entries in 'radius_inner' and 'radius_outer' must match.");
255 
256  bool youngs_modulus_set(false);
257  bool poissons_ratio_set(false);
258 
259  // All domain integral types block restrict the objects created by this action.
260  _blocks = getParam<std::vector<SubdomainName>>("block");
261 
262  MultiMooseEnum integral_moose_enums = getParam<MultiMooseEnum>("integrals");
263  for (unsigned int i = 0; i < integral_moose_enums.size(); ++i)
264  {
265  _displacements = getParam<std::vector<VariableName>>("displacements");
266 
267  if (_displacements.size() < 2)
268  paramError(
269  "displacements",
270  "DomainIntegral error: The size of the displacements vector should at least be 2.");
271 
272  if (integral_moose_enums[i] != "JIntegral" && integral_moose_enums[i] != "CIntegral" &&
273  integral_moose_enums[i] != "KFromJIntegral")
274  {
275  // Check that parameters required for interaction integrals are defined
276  if (!(isParamValid("poissons_ratio")) || !(isParamValid("youngs_modulus")))
277  mooseError(
278  "DomainIntegral error: must set Poisson's ratio and Young's modulus for integral: ",
279  integral_moose_enums[i]);
280 
281  _poissons_ratio = getParam<Real>("poissons_ratio");
282  poissons_ratio_set = true;
283  _youngs_modulus = getParam<Real>("youngs_modulus");
284  youngs_modulus_set = true;
285  }
286 
287  _integrals.insert(INTEGRAL(int(integral_moose_enums.get(i))));
288  }
289 
290  if ((_integrals.count(J_INTEGRAL) != 0 && _integrals.count(C_INTEGRAL) != 0) ||
291  (_integrals.count(J_INTEGRAL) != 0 && _integrals.count(K_FROM_J_INTEGRAL) != 0) ||
292  (_integrals.count(C_INTEGRAL) != 0 && _integrals.count(K_FROM_J_INTEGRAL) != 0))
293  paramError("integrals",
294  "JIntegral, CIntegral, and KFromJIntegral options are mutually exclusive");
295 
296  // Acommodate deprecated parameter convert_J_to_K
297  if (_convert_J_to_K && _integrals.count(K_FROM_J_INTEGRAL) != 0)
298  {
300  _integrals.erase(J_INTEGRAL);
301  }
302 
303  if (isParamValid("temperature"))
304  _temp = getParam<VariableName>("temperature");
305 
306  if (_temp != "" && !isParamValid("eigenstrain_names"))
307  paramError(
308  "eigenstrain_names",
309  "DomainIntegral error: must provide `eigenstrain_names` when temperature is coupled.");
310 
312  _integrals.count(INTERACTION_INTEGRAL_KII) == 0 ||
314  paramError("integrals",
315  "DomainIntegral error: must calculate KI, KII and KIII to get equivalent K.");
316 
317  if (isParamValid("output_variable"))
318  {
319  _output_variables = getParam<std::vector<VariableName>>("output_variable");
320  if (_crack_front_points.size() > 0)
321  paramError("output_variables",
322  "'output_variables' not yet supported with 'crack_front_points'");
323  }
324 
325  if (_integrals.count(K_FROM_J_INTEGRAL) != 0)
326  {
327  if (!isParamValid("youngs_modulus") || !isParamValid("poissons_ratio"))
328  mooseError("DomainIntegral error: must set Young's modulus and Poisson's ratio "
329  "if K_FROM_J_INTEGRAL is selected.");
330  if (!youngs_modulus_set)
331  _youngs_modulus = getParam<Real>("youngs_modulus");
332  if (!poissons_ratio_set)
333  _poissons_ratio = getParam<Real>("poissons_ratio");
334  }
335 
336  if (_integrals.count(J_INTEGRAL) != 0 || _integrals.count(C_INTEGRAL) != 0 ||
337  _integrals.count(K_FROM_J_INTEGRAL) != 0)
338  {
339  if (isParamValid("eigenstrain_gradient"))
340  paramError("eigenstrain_gradient",
341  "'eigenstrain_gradient' cannot be specified when the computed integrals include "
342  "JIntegral, CIntegral, or KFromJIntegral");
343  if (isParamValid("body_force"))
344  paramError("body_force",
345  "'body_force' cannot be specified when the computed integrals include JIntegral, "
346  "CIntegral, or KFromJIntegral");
347  }
348  if (isParamValid("eigenstrain_gradient") && (_temp != "" || isParamValid("eigenstrain_names")))
349  paramError("eigenstrain_gradient",
350  "'eigenstrain_gradient' cannot be specified together with 'temperature' or "
351  "'eigenstrain_names'. These are for separate, mutually exclusive systems for "
352  "including the effect of eigenstrains");
353 }
const RealVectorValue _crack_direction_vector
Vector optionally used to prescribe direction of crack extension.
const bool _used_by_xfem_to_grow_crack
This determines if fracture integrals should be executed on nonlinear in order to grow the crack when...
unsigned int _symmetry_plane
Identifier for which plane is the symmetry plane.
VariableName _temp
Temperature variable.
std::vector< VariableName > _output_variables
List of variables for which values are to be sampled and output at the crack front points...
const std::vector< BoundaryName > & _boundary_names
Boundaries containing the crack front points.
UserObjectName _crack_front_points_provider
Name of crack front points provider user object used to optionally define the crack points...
const std::string _order
Order and family of the AuxVariables optionally created to output the values of q.
const bool _have_crack_direction_vector_end_2
Whether the crack direction vector at the 2nd end of the crack has been provided. ...
MaterialPropertyName _functionally_graded_youngs_modulus
Material property name for spatially-dependent Youngs modulus for functionally graded materials...
Action(const InputParameters &parameters)
const bool _have_crack_direction_vector
Whether the crack direction vector has been provided.
std::vector< BoundaryName > _intersecting_boundary_names
Names of boundaries optionally used for improved computation of crack extension direction at ends of ...
const bool _use_ad
Whether to create automatic differentiation objects from the action.
std::vector< Real > _radius_inner
Sets of inner and outer radii of the rings used for the domain form of the fracture integrals...
bool _has_symmetry_plane
Whether the model has a symmetry plane passing through the plane of the crack.
bool _use_crack_front_points_provider
Whether to use a crack front points provider.
bool isParamValid(const std::string &name) const
bool _treat_as_2d
Whether fracture computations for a 3D model should be treated as though it were a 2D model...
Real _poissons_ratio
Poisson&#39;s ratio of material.
const std::string _family
MooseEnum _direction_method_moose_enum
Enum used to define the method to compute crack front direction.
std::vector< VariableName > _displacements
Vector of displacement variables.
std::vector< BoundaryName > _crack_mouth_boundary_names
Names of boundaries optionally used to define the crack mouth location.
const T & getParam(const std::string &name) const
void paramError(const std::string &param, Args... args) const
std::vector< SubdomainName > _blocks
Blocks for which the domain integrals are to be computed.
MooseEnum _end_direction_method_moose_enum
Enum used to define the method to compute crack front direction at ends of crack front.
std::vector< unsigned int > _ring_vec
Vector of ids for the individual rings on which the fracture integral is computed.
bool _incremental
Whether the constitutive models for the mechanics calculations use an incremental form...
std::vector< Real > _radius_outer
MooseEnum _q_function_type
How the q function is evaluated (geometric distance from crack front or ring of elements) ...
const bool _have_crack_direction_vector_end_1
Whether the crack direction vector at the 1st end of the crack has been provided. ...
bool _convert_J_to_K
Whether to convert the J-integral to a stress intensity factor (K) –deprecated.
MaterialPropertyName _functionally_graded_youngs_modulus_crack_dir_gradient
Material property name for the Youngs modulus derivative for functionally graded materials.
bool _use_displaced_mesh
Whether to compute the fracture integrals on the displaced mesh.
const RealVectorValue _crack_direction_vector_end_1
Vector optionally used to prescribe direction of crack extension at the 1st end of the crack...
INTEGRAL
Enum used to select the type of integral to be performed.
void mooseError(Args &&... args) const
bool _output_q
Whether to ouput the q function as a set of AuxVariables.
bool _get_equivalent_k
Whether to compute the equivalent K from the individual fracture integrals for mixed-mode fracture...
MooseEnum _position_type
How the distance along the crack front is measured (angle or distance)
std::vector< Point > _crack_front_points
User-defined vector of crack front points.
unsigned int _ring_last
Number of elements away from the crack tip to outside of outer ring with the topological q function...
unsigned int _ring_first
Number of elements away from the crack tip to inside of inner ring with the topological q function...
bool _fgm_crack
Whether the crack lives in a functionally-graded material.
const RealVectorValue _crack_direction_vector_end_2
Vector optionally used to prescribe direction of crack extension at the 2nd end of the crack...
unsigned int _axis_2d
Out-of-plane axis for 3D models treated as 2D.
bool _closed_loop
Indicates whether the crack forms a closed loop.
Real _youngs_modulus
Young&#39;s modulus of material.
std::set< INTEGRAL > _integrals
Container for enumerations describing the individual integrals computed.

◆ ~DomainIntegralAction()

DomainIntegralAction::~DomainIntegralAction ( )

Definition at line 355 of file DomainIntegralAction.C.

355 {}

Member Function Documentation

◆ act()

void DomainIntegralAction::act ( )
overridevirtual

Implements Action.

Definition at line 358 of file DomainIntegralAction.C.

359 {
360  const std::string uo_name("crackFrontDefinition");
361  const std::string ak_base_name("q");
362  const std::string av_base_name("q");
363  const unsigned int num_crack_front_points = calcNumCrackFrontPoints();
364  const std::string aux_stress_base_name("aux_stress");
365  const std::string aux_grad_disp_base_name("aux_grad_disp");
366 
367  std::string ad_prepend = "";
368  if (_use_ad)
369  ad_prepend = "AD";
370 
371  if (_current_task == "add_user_object")
372  {
373  const std::string uo_type_name("CrackFrontDefinition");
374 
375  InputParameters params = _factory.getValidParams(uo_type_name);
377  params.set<ExecFlagEnum>("execute_on") = {EXEC_INITIAL, EXEC_TIMESTEP_END, EXEC_NONLINEAR};
378  else
379  params.set<ExecFlagEnum>("execute_on") = {EXEC_INITIAL, EXEC_TIMESTEP_END};
380 
381  params.set<MooseEnum>("crack_direction_method") = _direction_method_moose_enum;
382  params.set<MooseEnum>("crack_end_direction_method") = _end_direction_method_moose_enum;
384  params.set<RealVectorValue>("crack_direction_vector") = _crack_direction_vector;
386  params.set<RealVectorValue>("crack_direction_vector_end_1") = _crack_direction_vector_end_1;
388  params.set<RealVectorValue>("crack_direction_vector_end_2") = _crack_direction_vector_end_2;
389  if (_crack_mouth_boundary_names.size() != 0)
390  params.set<std::vector<BoundaryName>>("crack_mouth_boundary") = _crack_mouth_boundary_names;
391  if (_intersecting_boundary_names.size() != 0)
392  params.set<std::vector<BoundaryName>>("intersecting_boundary") = _intersecting_boundary_names;
393  params.set<bool>("2d") = _treat_as_2d;
394  params.set<unsigned int>("axis_2d") = _axis_2d;
396  params.set<unsigned int>("symmetry_plane") = _symmetry_plane;
397  if (_boundary_names.size() != 0)
398  params.set<std::vector<BoundaryName>>("boundary") = _boundary_names;
399  if (_crack_front_points.size() != 0)
400  params.set<std::vector<Point>>("crack_front_points") = _crack_front_points;
402  params.applyParameters(parameters(),
403  {"crack_front_points_provider, number_points_from_provider"});
404  if (_closed_loop)
405  params.set<bool>("closed_loop") = _closed_loop;
406  params.set<bool>("use_displaced_mesh") = _use_displaced_mesh;
407  if (_integrals.count(INTERACTION_INTEGRAL_T) != 0)
408  {
409  params.set<VariableName>("disp_x") = _displacements[0];
410  params.set<VariableName>("disp_y") = _displacements[1];
411  if (_displacements.size() == 3)
412  params.set<VariableName>("disp_z") = _displacements[2];
413  params.set<bool>("t_stress") = true;
414  }
415 
416  unsigned int nrings = 0;
417  if (_q_function_type == TOPOLOGY)
418  {
419  params.set<bool>("q_function_rings") = true;
420  params.set<unsigned int>("last_ring") = _ring_last;
421  params.set<unsigned int>("first_ring") = _ring_first;
422  nrings = _ring_last - _ring_first + 1;
423  }
424  else if (_q_function_type == GEOMETRY)
425  {
426  params.set<std::vector<Real>>("j_integral_radius_inner") = _radius_inner;
427  params.set<std::vector<Real>>("j_integral_radius_outer") = _radius_outer;
428  nrings = _ring_vec.size();
429  }
430 
431  params.set<unsigned int>("nrings") = nrings;
432  params.set<MooseEnum>("q_function_type") = _q_function_type;
433 
434  _problem->addUserObject(uo_type_name, uo_name, params);
435  }
436  else if (_current_task == "add_aux_variable" && _output_q)
437  {
438  for (unsigned int ring_index = 0; ring_index < _ring_vec.size(); ++ring_index)
439  {
440  std::string aux_var_type;
441  if (_family == "LAGRANGE")
442  aux_var_type = "MooseVariable";
443  else if (_family == "MONOMIAL")
444  aux_var_type = "MooseVariableConstMonomial";
445  else if (_family == "SCALAR")
446  aux_var_type = "MooseVariableScalar";
447  else
448  mooseError("Unsupported finite element family in, " + name() +
449  ". Please use LAGRANGE, MONOMIAL, or SCALAR");
450 
451  auto params = _factory.getValidParams(aux_var_type);
452  params.set<MooseEnum>("order") = _order;
453  params.set<MooseEnum>("family") = _family;
454 
456  {
457  std::ostringstream av_name_stream;
458  av_name_stream << av_base_name << "_" << _ring_vec[ring_index];
459  _problem->addAuxVariable(aux_var_type, av_name_stream.str(), params);
460  }
461  else
462  {
463  for (unsigned int cfp_index = 0; cfp_index < num_crack_front_points; ++cfp_index)
464  {
465  std::ostringstream av_name_stream;
466  av_name_stream << av_base_name << "_" << cfp_index + 1 << "_" << _ring_vec[ring_index];
467  _problem->addAuxVariable(aux_var_type, av_name_stream.str(), params);
468  }
469  }
470  }
471  }
472 
473  else if (_current_task == "add_aux_kernel" && _output_q)
474  {
475  std::string ak_type_name;
476  unsigned int nrings = 0;
477  if (_q_function_type == GEOMETRY)
478  {
479  ak_type_name = "DomainIntegralQFunction";
480  nrings = _ring_vec.size();
481  }
482  else if (_q_function_type == TOPOLOGY)
483  {
484  ak_type_name = "DomainIntegralTopologicalQFunction";
485  nrings = _ring_last - _ring_first + 1;
486  }
487 
488  InputParameters params = _factory.getValidParams(ak_type_name);
489  params.set<ExecFlagEnum>("execute_on") = {EXEC_INITIAL, EXEC_TIMESTEP_END};
490  params.set<UserObjectName>("crack_front_definition") = uo_name;
491  params.set<bool>("use_displaced_mesh") = _use_displaced_mesh;
492 
493  for (unsigned int ring_index = 0; ring_index < nrings; ++ring_index)
494  {
495  if (_q_function_type == GEOMETRY)
496  {
497  params.set<Real>("j_integral_radius_inner") = _radius_inner[ring_index];
498  params.set<Real>("j_integral_radius_outer") = _radius_outer[ring_index];
499  }
500  else if (_q_function_type == TOPOLOGY)
501  {
502  params.set<unsigned int>("ring_index") = _ring_first + ring_index;
503  }
504 
506  {
507  std::ostringstream ak_name_stream;
508  ak_name_stream << ak_base_name << "_" << _ring_vec[ring_index];
509  std::ostringstream av_name_stream;
510  av_name_stream << av_base_name << "_" << _ring_vec[ring_index];
511  params.set<AuxVariableName>("variable") = av_name_stream.str();
512  _problem->addAuxKernel(ak_type_name, ak_name_stream.str(), params);
513  }
514  else
515  {
516  for (unsigned int cfp_index = 0; cfp_index < num_crack_front_points; ++cfp_index)
517  {
518  std::ostringstream ak_name_stream;
519  ak_name_stream << ak_base_name << "_" << cfp_index + 1 << "_" << _ring_vec[ring_index];
520  std::ostringstream av_name_stream;
521  av_name_stream << av_base_name << "_" << cfp_index + 1 << "_" << _ring_vec[ring_index];
522  params.set<AuxVariableName>("variable") = av_name_stream.str();
523  params.set<unsigned int>("crack_front_point_index") = cfp_index;
524  _problem->addAuxKernel(ak_type_name, ak_name_stream.str(), params);
525  }
526  }
527  }
528  }
529 
530  else if (_current_task == "add_postprocessor")
531  {
532  for (std::set<INTEGRAL>::iterator sit = _integrals.begin(); sit != _integrals.end(); ++sit)
533  {
534  std::string pp_base_name;
535  switch (*sit)
536  {
537  case J_INTEGRAL:
538  pp_base_name = "J";
539  break;
540 
541  case C_INTEGRAL:
542  pp_base_name = "C";
543  break;
544 
545  case K_FROM_J_INTEGRAL:
546  pp_base_name = "K";
547  break;
548 
550  pp_base_name = "II_KI";
551  break;
552 
554  pp_base_name = "II_KII";
555  break;
556 
558  pp_base_name = "II_KIII";
559  break;
560 
562  pp_base_name = "II_T";
563  break;
564  }
565  const std::string pp_type_name("VectorPostprocessorComponent");
566  InputParameters params = _factory.getValidParams(pp_type_name);
567  for (unsigned int ring_index = 0; ring_index < _ring_vec.size(); ++ring_index)
568  {
570  {
571  params.set<VectorPostprocessorName>("vectorpostprocessor") =
572  pp_base_name + "_2DVPP_" + Moose::stringify(_ring_vec[ring_index]);
573  std::string pp_name = pp_base_name + +"_" + Moose::stringify(_ring_vec[ring_index]);
574  params.set<unsigned int>("index") = 0;
575  params.set<std::string>("vector_name") =
576  pp_base_name + "_" + Moose::stringify(_ring_vec[ring_index]);
577  _problem->addPostprocessor(pp_type_name, pp_name, params);
578  }
579  else
580  {
581  for (unsigned int cfp_index = 0; cfp_index < num_crack_front_points; ++cfp_index)
582  {
583  params.set<VectorPostprocessorName>("vectorpostprocessor") =
584  pp_base_name + "_" + Moose::stringify(_ring_vec[ring_index]);
585  std::string pp_name = pp_base_name + "_" + Moose::stringify(cfp_index + 1) + "_" +
586  Moose::stringify(_ring_vec[ring_index]);
587  params.set<unsigned int>("index") = cfp_index;
588  params.set<std::string>("vector_name") =
589  pp_base_name + "_" + Moose::stringify(_ring_vec[ring_index]);
590  _problem->addPostprocessor(pp_type_name, pp_name, params);
591  }
592  }
593  }
594  }
595 
596  if (_get_equivalent_k)
597  {
598  std::string pp_base_name("Keq");
599  const std::string pp_type_name("VectorPostprocessorComponent");
600  InputParameters params = _factory.getValidParams(pp_type_name);
601  for (unsigned int ring_index = 0; ring_index < _ring_vec.size(); ++ring_index)
602  {
604  {
605  params.set<VectorPostprocessorName>("vectorpostprocessor") =
606  pp_base_name + "_2DVPP_" + Moose::stringify(_ring_vec[ring_index]);
607  std::string pp_name = pp_base_name + +"_" + Moose::stringify(_ring_vec[ring_index]);
608  params.set<unsigned int>("index") = 0;
609  params.set<std::string>("vector_name") =
610  pp_base_name + "_" + Moose::stringify(_ring_vec[ring_index]);
611  _problem->addPostprocessor(pp_type_name, pp_name, params);
612  }
613  else
614  {
615  for (unsigned int cfp_index = 0; cfp_index < num_crack_front_points; ++cfp_index)
616  {
617  params.set<VectorPostprocessorName>("vectorpostprocessor") =
618  pp_base_name + "_" + Moose::stringify(_ring_vec[ring_index]);
619  std::string pp_name = pp_base_name + "_" + Moose::stringify(cfp_index + 1) + "_" +
620  Moose::stringify(_ring_vec[ring_index]);
621  params.set<unsigned int>("index") = cfp_index;
622  params.set<std::string>("vector_name") =
623  pp_base_name + "_" + Moose::stringify(_ring_vec[ring_index]);
624  _problem->addPostprocessor(pp_type_name, pp_name, params);
625  }
626  }
627  }
628  }
629 
630  for (unsigned int i = 0; i < _output_variables.size(); ++i)
631  {
632  const std::string ov_base_name(_output_variables[i]);
633  const std::string pp_type_name("CrackFrontData");
634  InputParameters params = _factory.getValidParams(pp_type_name);
635  params.set<ExecFlagEnum>("execute_on") = EXEC_TIMESTEP_END;
636  params.set<UserObjectName>("crack_front_definition") = uo_name;
638  {
639  std::ostringstream pp_name_stream;
640  pp_name_stream << ov_base_name << "_crack";
641  params.set<VariableName>("variable") = _output_variables[i];
642  _problem->addPostprocessor(pp_type_name, pp_name_stream.str(), params);
643  }
644  else
645  {
646  for (unsigned int cfp_index = 0; cfp_index < num_crack_front_points; ++cfp_index)
647  {
648  std::ostringstream pp_name_stream;
649  pp_name_stream << ov_base_name << "_crack_" << cfp_index + 1;
650  params.set<VariableName>("variable") = _output_variables[i];
651  params.set<unsigned int>("crack_front_point_index") = cfp_index;
652  _problem->addPostprocessor(pp_type_name, pp_name_stream.str(), params);
653  }
654  }
655  }
656  }
657 
658  else if (_current_task == "add_vector_postprocessor")
659  {
660  if (_integrals.count(J_INTEGRAL) != 0 || _integrals.count(C_INTEGRAL) != 0 ||
661  _integrals.count(K_FROM_J_INTEGRAL) != 0)
662  {
663  std::string vpp_base_name;
664  std::string jintegral_selection = "JIntegral";
665 
666  if (_integrals.count(J_INTEGRAL) != 0)
667  {
668  vpp_base_name = "J";
669  jintegral_selection = "JIntegral";
670  }
671  else if (_integrals.count(K_FROM_J_INTEGRAL) != 0)
672  {
673  vpp_base_name = "K";
674  jintegral_selection = "KFromJIntegral";
675  }
676  else if (_integrals.count(C_INTEGRAL) != 0)
677  {
678  vpp_base_name = "C";
679  jintegral_selection = "CIntegral";
680  }
681 
683  vpp_base_name += "_2DVPP";
684 
685  const std::string vpp_type_name("JIntegral");
686  InputParameters params = _factory.getValidParams(vpp_type_name);
687  if (!getParam<bool>("output_vpp"))
688  params.set<std::vector<OutputName>>("outputs") = {"none"};
689 
690  params.set<ExecFlagEnum>("execute_on") = EXEC_TIMESTEP_END;
691  params.set<UserObjectName>("crack_front_definition") = uo_name;
692  params.set<std::vector<SubdomainName>>("block") = {_blocks};
693  params.set<MooseEnum>("position_type") = _position_type;
694 
695  if (_integrals.count(K_FROM_J_INTEGRAL) != 0)
696  {
697  params.set<Real>("youngs_modulus") = _youngs_modulus;
698  params.set<Real>("poissons_ratio") = _poissons_ratio;
699  }
700 
702  params.set<unsigned int>("symmetry_plane") = _symmetry_plane;
703 
704  // Select the integral type to be computed in JIntegral vector postprocessor
705  params.set<MooseEnum>("integral") = jintegral_selection;
706 
707  params.set<std::vector<VariableName>>("displacements") = _displacements;
708  params.set<bool>("use_displaced_mesh") = _use_displaced_mesh;
709  for (unsigned int ring_index = 0; ring_index < _ring_vec.size(); ++ring_index)
710  {
711  params.set<unsigned int>("ring_index") = _ring_vec[ring_index];
712  params.set<MooseEnum>("q_function_type") = _q_function_type;
713 
714  std::string vpp_name = vpp_base_name + "_" + Moose::stringify(_ring_vec[ring_index]);
715  _problem->addVectorPostprocessor(vpp_type_name, vpp_name, params);
716  }
717  }
718 
719  if (_integrals.count(INTERACTION_INTEGRAL_KI) != 0 ||
720  _integrals.count(INTERACTION_INTEGRAL_KII) != 0 ||
723  {
726  paramError("symmetry_plane",
727  "In DomainIntegral, symmetry_plane option cannot be used with mode-II or "
728  "mode-III interaction integral");
729 
730  std::string vpp_base_name;
731  std::string vpp_type_name(ad_prepend + "InteractionIntegral");
732 
733  InputParameters params = _factory.getValidParams(vpp_type_name);
734  if (!getParam<bool>("output_vpp"))
735  params.set<std::vector<OutputName>>("outputs") = {"none"};
736 
738  params.set<ExecFlagEnum>("execute_on") = {EXEC_TIMESTEP_END, EXEC_NONLINEAR};
739  else
740  params.set<ExecFlagEnum>("execute_on") = {EXEC_TIMESTEP_END};
741 
742  if (isParamValid("additional_eigenstrain_00") && isParamValid("additional_eigenstrain_01") &&
743  isParamValid("additional_eigenstrain_11") && isParamValid("additional_eigenstrain_22"))
744  {
745  params.set<CoupledName>("additional_eigenstrain_00") =
746  getParam<CoupledName>("additional_eigenstrain_00");
747  params.set<CoupledName>("additional_eigenstrain_01") =
748  getParam<CoupledName>("additional_eigenstrain_01");
749  params.set<CoupledName>("additional_eigenstrain_11") =
750  getParam<CoupledName>("additional_eigenstrain_11");
751  params.set<CoupledName>("additional_eigenstrain_22") =
752  getParam<CoupledName>("additional_eigenstrain_22");
753  }
754 
755  params.set<UserObjectName>("crack_front_definition") = uo_name;
756  params.set<bool>("use_displaced_mesh") = _use_displaced_mesh;
757  params.set<std::vector<SubdomainName>>("block") = {_blocks};
758 
760  params.set<unsigned int>("symmetry_plane") = _symmetry_plane;
761 
762  if (_fgm_crack)
763  {
764  params.set<MaterialPropertyName>(
765  "functionally_graded_youngs_modulus_crack_dir_gradient") = {
767  params.set<MaterialPropertyName>("functionally_graded_youngs_modulus") = {
769  }
770 
771  params.set<Real>("poissons_ratio") = _poissons_ratio;
772  params.set<Real>("youngs_modulus") = _youngs_modulus;
773  params.set<std::vector<VariableName>>("displacements") = _displacements;
774  if (_temp != "")
775  params.set<std::vector<VariableName>>("temperature") = {_temp};
776 
777  if (parameters().isParamValid("eigenstrain_gradient"))
778  params.set<MaterialPropertyName>("eigenstrain_gradient") =
779  parameters().get<MaterialPropertyName>("eigenstrain_gradient");
780  if (parameters().isParamValid("body_force"))
781  params.set<MaterialPropertyName>("body_force") =
782  parameters().get<MaterialPropertyName>("body_force");
783 
784  for (std::set<INTEGRAL>::iterator sit = _integrals.begin(); sit != _integrals.end(); ++sit)
785  {
786  switch (*sit)
787  {
788  case J_INTEGRAL:
789  continue;
790 
791  case C_INTEGRAL:
792  continue;
793 
794  case K_FROM_J_INTEGRAL:
795  continue;
796 
798  vpp_base_name = "II_KI";
799  params.set<Real>("K_factor") =
800  0.5 * _youngs_modulus / (1.0 - std::pow(_poissons_ratio, 2.0));
801  params.set<MooseEnum>("sif_mode") = "KI";
802  break;
803 
805  vpp_base_name = "II_KII";
806  params.set<Real>("K_factor") =
807  0.5 * _youngs_modulus / (1.0 - std::pow(_poissons_ratio, 2.0));
808  params.set<MooseEnum>("sif_mode") = "KII";
809  break;
810 
812  vpp_base_name = "II_KIII";
813  params.set<Real>("K_factor") = 0.5 * _youngs_modulus / (1.0 + _poissons_ratio);
814  params.set<MooseEnum>("sif_mode") = "KIII";
815  break;
816 
818  vpp_base_name = "II_T";
819  params.set<Real>("K_factor") = _youngs_modulus / (1 - std::pow(_poissons_ratio, 2));
820  params.set<MooseEnum>("sif_mode") = "T";
821  break;
822  }
824  vpp_base_name += "_2DVPP";
825  for (unsigned int ring_index = 0; ring_index < _ring_vec.size(); ++ring_index)
826  {
827  params.set<unsigned int>("ring_index") = _ring_vec[ring_index];
828  params.set<MooseEnum>("q_function_type") = _q_function_type;
829 
830  std::string vpp_name = vpp_base_name + "_" + Moose::stringify(_ring_vec[ring_index]);
831  _problem->addVectorPostprocessor(vpp_type_name, vpp_name, params);
832  }
833  }
834  }
835 
836  if (_get_equivalent_k)
837  {
838  std::string vpp_base_name("Keq");
840  vpp_base_name += "_2DVPP";
841  const std::string vpp_type_name("MixedModeEquivalentK");
842  InputParameters params = _factory.getValidParams(vpp_type_name);
843  params.set<ExecFlagEnum>("execute_on") = EXEC_TIMESTEP_END;
844  params.set<Real>("poissons_ratio") = _poissons_ratio;
845 
846  for (unsigned int ring_index = 0; ring_index < _ring_vec.size(); ++ring_index)
847  {
848  std::string ki_name = "II_KI_";
849  std::string kii_name = "II_KII_";
850  std::string kiii_name = "II_KIII_";
851  params.set<unsigned int>("ring_index") = _ring_vec[ring_index];
853  {
854  params.set<VectorPostprocessorName>("KI_vectorpostprocessor") =
855  ki_name + "2DVPP_" + Moose::stringify(_ring_vec[ring_index]);
856  params.set<VectorPostprocessorName>("KII_vectorpostprocessor") =
857  kii_name + "2DVPP_" + Moose::stringify(_ring_vec[ring_index]);
858  params.set<VectorPostprocessorName>("KIII_vectorpostprocessor") =
859  kiii_name + "2DVPP_" + Moose::stringify(_ring_vec[ring_index]);
860  }
861  else
862  {
863  params.set<VectorPostprocessorName>("KI_vectorpostprocessor") =
864  ki_name + Moose::stringify(_ring_vec[ring_index]);
865  params.set<VectorPostprocessorName>("KII_vectorpostprocessor") =
866  kii_name + Moose::stringify(_ring_vec[ring_index]);
867  params.set<VectorPostprocessorName>("KIII_vectorpostprocessor") =
868  kiii_name + Moose::stringify(_ring_vec[ring_index]);
869  }
870  params.set<std::string>("KI_vector_name") =
871  ki_name + Moose::stringify(_ring_vec[ring_index]);
872  params.set<std::string>("KII_vector_name") =
873  kii_name + Moose::stringify(_ring_vec[ring_index]);
874  params.set<std::string>("KIII_vector_name") =
875  kiii_name + Moose::stringify(_ring_vec[ring_index]);
876  std::string vpp_name = vpp_base_name + "_" + Moose::stringify(_ring_vec[ring_index]);
877  _problem->addVectorPostprocessor(vpp_type_name, vpp_name, params);
878  }
879  }
880 
882  {
883  for (unsigned int i = 0; i < _output_variables.size(); ++i)
884  {
885  const std::string vpp_type_name("VectorOfPostprocessors");
886  InputParameters params = _factory.getValidParams(vpp_type_name);
887  params.set<ExecFlagEnum>("execute_on") = EXEC_TIMESTEP_END;
888  std::ostringstream vpp_name_stream;
889  vpp_name_stream << _output_variables[i] << "_crack";
890  std::vector<PostprocessorName> postprocessor_names;
891  for (unsigned int cfp_index = 0; cfp_index < num_crack_front_points; ++cfp_index)
892  {
893  std::ostringstream pp_name_stream;
894  pp_name_stream << vpp_name_stream.str() << "_" << cfp_index + 1;
895  postprocessor_names.push_back(pp_name_stream.str());
896  }
897  params.set<std::vector<PostprocessorName>>("postprocessors") = postprocessor_names;
898  _problem->addVectorPostprocessor(vpp_type_name, vpp_name_stream.str(), params);
899  }
900  }
901  }
902 
903  else if (_current_task == "add_material")
904  {
905  if (_temp != "")
906  {
907  std::string mater_name;
908  const std::string mater_type_name("ThermalFractureIntegral");
909  mater_name = "ThermalFractureIntegral";
910 
911  InputParameters params = _factory.getValidParams(mater_type_name);
912  params.set<std::vector<MaterialPropertyName>>("eigenstrain_names") =
913  getParam<std::vector<MaterialPropertyName>>("eigenstrain_names");
914  params.set<std::vector<VariableName>>("temperature") = {_temp};
915  params.set<std::vector<SubdomainName>>("block") = {_blocks};
916  _problem->addMaterial(mater_type_name, mater_name, params);
917  }
918  MultiMooseEnum integral_moose_enums = getParam<MultiMooseEnum>("integrals");
919  bool have_j_integral = false;
920  bool have_c_integral = false;
921 
922  for (auto ime : integral_moose_enums)
923  {
924  if (ime == "JIntegral" || ime == "CIntegral" || ime == "KFromJIntegral" ||
925  ime == "InteractionIntegralKI" || ime == "InteractionIntegralKII" ||
926  ime == "InteractionIntegralKIII" || ime == "InteractionIntegralT")
927  have_j_integral = true;
928 
929  if (ime == "CIntegral")
930  have_c_integral = true;
931  }
932  if (have_j_integral)
933  {
934  std::string mater_name;
935  const std::string mater_type_name(ad_prepend + "StrainEnergyDensity");
936  mater_name = ad_prepend + "StrainEnergyDensity";
937 
938  InputParameters params = _factory.getValidParams(mater_type_name);
939  _incremental = getParam<bool>("incremental");
940  params.set<bool>("incremental") = _incremental;
941  params.set<std::vector<SubdomainName>>("block") = {_blocks};
942  _problem->addMaterial(mater_type_name, mater_name, params);
943 
944  {
945  std::string mater_name;
946  const std::string mater_type_name(ad_prepend + "EshelbyTensor");
947  mater_name = ad_prepend + "EshelbyTensor";
948 
949  InputParameters params = _factory.getValidParams(mater_type_name);
950  _displacements = getParam<std::vector<VariableName>>("displacements");
951  params.set<std::vector<VariableName>>("displacements") = _displacements;
952  params.set<std::vector<SubdomainName>>("block") = {_blocks};
953 
954  if (have_c_integral)
955  params.set<bool>("compute_dissipation") = true;
956 
957  if (_temp != "")
958  params.set<std::vector<VariableName>>("temperature") = {_temp};
959 
960  _problem->addMaterial(mater_type_name, mater_name, params);
961  }
962  // Strain energy rate density needed for C(t)/C* integral
963  if (have_c_integral)
964  {
965  std::string mater_name;
966  const std::string mater_type_name(ad_prepend + "StrainEnergyRateDensity");
967  mater_name = ad_prepend + "StrainEnergyRateDensity";
968 
969  InputParameters params = _factory.getValidParams(mater_type_name);
970  params.set<std::vector<SubdomainName>>("block") = {_blocks};
971  params.set<std::vector<MaterialName>>("inelastic_models") =
972  getParam<std::vector<MaterialName>>("inelastic_models");
973 
974  _problem->addMaterial(mater_type_name, mater_name, params);
975  }
976  }
977  }
978 }
const RealVectorValue _crack_direction_vector
Vector optionally used to prescribe direction of crack extension.
const bool _used_by_xfem_to_grow_crack
This determines if fracture integrals should be executed on nonlinear in order to grow the crack when...
unsigned int _symmetry_plane
Identifier for which plane is the symmetry plane.
VariableName _temp
Temperature variable.
std::vector< VariableName > _output_variables
List of variables for which values are to be sampled and output at the crack front points...
const std::vector< BoundaryName > & _boundary_names
Boundaries containing the crack front points.
const std::string _order
Order and family of the AuxVariables optionally created to output the values of q.
const bool _have_crack_direction_vector_end_2
Whether the crack direction vector at the 2nd end of the crack has been provided. ...
MaterialPropertyName _functionally_graded_youngs_modulus
Material property name for spatially-dependent Youngs modulus for functionally graded materials...
std::vector< std::pair< R1, R2 > > get(const std::string &param1, const std::string &param2) const
T & set(const std::string &name, bool quiet_mode=false)
const bool _have_crack_direction_vector
Whether the crack direction vector has been provided.
std::vector< BoundaryName > _intersecting_boundary_names
Names of boundaries optionally used for improved computation of crack extension direction at ends of ...
InputParameters getValidParams(const std::string &name) const
unsigned int calcNumCrackFrontPoints()
Compute the number of points on the crack front.
void applyParameters(const InputParameters &common, const std::vector< std::string > &exclude={}, const bool allow_private=false)
const ExecFlagType EXEC_TIMESTEP_END
const bool _use_ad
Whether to create automatic differentiation objects from the action.
std::vector< Real > _radius_inner
Sets of inner and outer radii of the rings used for the domain form of the fracture integrals...
bool _has_symmetry_plane
Whether the model has a symmetry plane passing through the plane of the crack.
bool _use_crack_front_points_provider
Whether to use a crack front points provider.
virtual const std::string & name() const
bool isParamValid(const std::string &name) const
Factory & _factory
bool _treat_as_2d
Whether fracture computations for a 3D model should be treated as though it were a 2D model...
Real _poissons_ratio
Poisson&#39;s ratio of material.
const std::string _family
MooseEnum _direction_method_moose_enum
Enum used to define the method to compute crack front direction.
std::vector< VariableName > _displacements
Vector of displacement variables.
std::vector< BoundaryName > _crack_mouth_boundary_names
Names of boundaries optionally used to define the crack mouth location.
const T & getParam(const std::string &name) const
const std::string & _current_task
void paramError(const std::string &param, Args... args) const
std::vector< SubdomainName > _blocks
Blocks for which the domain integrals are to be computed.
std::string stringify(const T &t)
MooseEnum _end_direction_method_moose_enum
Enum used to define the method to compute crack front direction at ends of crack front.
std::vector< unsigned int > _ring_vec
Vector of ids for the individual rings on which the fracture integral is computed.
bool _incremental
Whether the constitutive models for the mechanics calculations use an incremental form...
const ExecFlagType EXEC_NONLINEAR
std::vector< Real > _radius_outer
MooseEnum _q_function_type
How the q function is evaluated (geometric distance from crack front or ring of elements) ...
const bool _have_crack_direction_vector_end_1
Whether the crack direction vector at the 1st end of the crack has been provided. ...
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
MaterialPropertyName _functionally_graded_youngs_modulus_crack_dir_gradient
Material property name for the Youngs modulus derivative for functionally graded materials.
bool _use_displaced_mesh
Whether to compute the fracture integrals on the displaced mesh.
std::vector< VariableName > CoupledName
const RealVectorValue _crack_direction_vector_end_1
Vector optionally used to prescribe direction of crack extension at the 1st end of the crack...
void mooseError(Args &&... args) const
bool _output_q
Whether to ouput the q function as a set of AuxVariables.
std::shared_ptr< FEProblemBase > & _problem
const InputParameters & parameters() const
bool _get_equivalent_k
Whether to compute the equivalent K from the individual fracture integrals for mixed-mode fracture...
MooseEnum _position_type
How the distance along the crack front is measured (angle or distance)
std::vector< Point > _crack_front_points
User-defined vector of crack front points.
unsigned int _ring_last
Number of elements away from the crack tip to outside of outer ring with the topological q function...
unsigned int _ring_first
Number of elements away from the crack tip to inside of inner ring with the topological q function...
bool _fgm_crack
Whether the crack lives in a functionally-graded material.
const RealVectorValue _crack_direction_vector_end_2
Vector optionally used to prescribe direction of crack extension at the 2nd end of the crack...
MooseUnits pow(const MooseUnits &, int)
unsigned int _axis_2d
Out-of-plane axis for 3D models treated as 2D.
bool _closed_loop
Indicates whether the crack forms a closed loop.
Real _youngs_modulus
Young&#39;s modulus of material.
std::set< INTEGRAL > _integrals
Container for enumerations describing the individual integrals computed.
bool isParamValid(const std::string &name) const
const ExecFlagType EXEC_INITIAL

◆ addRelationshipManagers() [1/3]

virtual void Action::addRelationshipManagers

◆ addRelationshipManagers() [2/3]

bool Action::addRelationshipManagers

◆ addRelationshipManagers() [3/3]

void DomainIntegralAction::addRelationshipManagers ( Moose::RelationshipManagerType  input_rm_type)
overridevirtual

Reimplemented from Action.

Definition at line 1016 of file DomainIntegralAction.C.

1017 {
1018  if (_integrals.count(INTERACTION_INTEGRAL_T) != 0)
1019  {
1020  InputParameters params = _factory.getValidParams("CrackFrontDefinition");
1021  addRelationshipManagers(input_rm_type, params);
1022  }
1023 }
InputParameters getValidParams(const std::string &name) const
Factory & _factory
virtual void addRelationshipManagers(Moose::RelationshipManagerType input_rm_type) override
std::set< INTEGRAL > _integrals
Container for enumerations describing the individual integrals computed.

◆ calcNumCrackFrontPoints()

unsigned int DomainIntegralAction::calcNumCrackFrontPoints ( )
protected

Compute the number of points on the crack front.

This is either the number of points in the crack front nodeset, or the number of points from the crack front points provider.

Returns
Number of points on the crack front

Definition at line 981 of file DomainIntegralAction.C.

Referenced by act().

982 {
983  unsigned int num_points = 0;
984  if (_boundary_names.size() != 0)
985  {
986  std::vector<BoundaryID> bids = _mesh->getBoundaryIDs(_boundary_names, true);
987  std::set<unsigned int> nodes;
988 
989  ConstBndNodeRange & bnd_nodes = *_mesh->getBoundaryNodeRange();
990  for (ConstBndNodeRange::const_iterator nd = bnd_nodes.begin(); nd != bnd_nodes.end(); ++nd)
991  {
992  const BndNode * bnode = *nd;
993  BoundaryID boundary_id = bnode->_bnd_id;
994 
995  for (unsigned int ibid = 0; ibid < bids.size(); ++ibid)
996  {
997  if (boundary_id == bids[ibid])
998  {
999  nodes.insert(bnode->_node->id());
1000  break;
1001  }
1002  }
1003  }
1004  num_points = nodes.size();
1005  }
1006  else if (_crack_front_points.size() != 0)
1007  num_points = _crack_front_points.size();
1009  num_points = getParam<unsigned int>("number_points_from_provider");
1010  else
1011  mooseError("Must define either 'boundary' or 'crack_front_points'");
1012  return num_points;
1013 }
const std::vector< BoundaryName > & _boundary_names
Boundaries containing the crack front points.
BoundaryID _bnd_id
bool _use_crack_front_points_provider
Whether to use a crack front points provider.
Node * _node
boundary_id_type BoundaryID
std::shared_ptr< MooseMesh > & _mesh
void mooseError(Args &&... args) const
std::vector< Point > _crack_front_points
User-defined vector of crack front points.
StoredRange< MooseMesh::const_bnd_node_iterator, const BndNode *> ConstBndNodeRange

◆ validParams()

InputParameters DomainIntegralAction::validParams ( )
static

Definition at line 32 of file DomainIntegralAction.C.

33 {
36  MultiMooseEnum integral_vec(
37  "JIntegral CIntegral KFromJIntegral InteractionIntegralKI InteractionIntegralKII "
38  "InteractionIntegralKIII InteractionIntegralT");
39  params.addClassDescription(
40  "Creates the MOOSE objects needed to compute fraction domain integrals");
41  params.addRequiredParam<MultiMooseEnum>("integrals",
42  integral_vec,
43  "Domain integrals to calculate. Choices are: " +
44  integral_vec.getRawNames());
45  params.addParam<std::vector<BoundaryName>>(
46  "boundary", {}, "Boundary containing the crack front points");
47  params.addParam<std::vector<Point>>("crack_front_points", "Set of points to define crack front");
48  params.addParam<std::string>(
49  "order", "FIRST", "Specifies the order of the FE shape function to use for q AuxVariables");
50  params.addParam<std::string>(
51  "family", "LAGRANGE", "Specifies the family of FE shape functions to use for q AuxVariables");
52  params.addParam<std::vector<Real>>("radius_inner", "Inner radius for volume integral domain");
53  params.addParam<std::vector<Real>>("radius_outer", "Outer radius for volume integral domain");
54  params.addParam<unsigned int>("ring_first",
55  "The first ring of elements for volume integral domain");
56  params.addParam<unsigned int>("ring_last",
57  "The last ring of elements for volume integral domain");
58  params.addParam<std::vector<VariableName>>(
59  "output_variable", "Variable values to be reported along the crack front");
60  params.addParam<Real>("poissons_ratio", "Poisson's ratio");
61  params.addParam<Real>("youngs_modulus", "Young's modulus");
62  params.addParam<std::vector<SubdomainName>>(
63  "block", {}, "The block ids where integrals are defined");
64  params.addParam<std::vector<VariableName>>(
65  "displacements",
66  {},
67  "The displacements appropriate for the simulation geometry and coordinate system");
68  params.addParam<VariableName>("temperature", "", "The temperature");
69  params.addParam<MaterialPropertyName>(
70  "functionally_graded_youngs_modulus_crack_dir_gradient",
71  "Gradient of the spatially varying Young's modulus provided in "
72  "'functionally_graded_youngs_modulus' in the direction of crack extension.");
73  params.addParam<MaterialPropertyName>(
74  "functionally_graded_youngs_modulus",
75  "Spatially varying elasticity modulus variable. This input is required when "
76  "using the functionally graded material capability.");
77  params.addCoupledVar("additional_eigenstrain_00",
78  "Optional additional eigenstrain variable that will be accounted for in the "
79  "interaction integral (component 00 or XX).");
80  params.addCoupledVar("additional_eigenstrain_01",
81  "Optional additional eigenstrain variable that will be accounted for in the "
82  "interaction integral (component 01 or XY).");
83  params.addCoupledVar("additional_eigenstrain_11",
84  "Optional additional eigenstrain variable that will be accounted for in the "
85  "interaction integral (component 11 or YY).");
86  params.addCoupledVar("additional_eigenstrain_22",
87  "Optional additional eigenstrain variable that will be accounted for in the "
88  "interaction integral (component 22 or ZZ).");
89  params.addCoupledVar("additional_eigenstrain_02",
90  "Optional additional eigenstrain variable that will be accounted for in the "
91  "interaction integral (component 02 or XZ).");
92  params.addCoupledVar("additional_eigenstrain_12",
93  "Optional additional eigenstrain variable that will be accounted for in the "
94  "interaction integral (component 12 or XZ).");
95  params.addParamNamesToGroup(
96  "additional_eigenstrain_00 additional_eigenstrain_01 additional_eigenstrain_11 "
97  "additional_eigenstrain_22 additional_eigenstrain_02 additional_eigenstrain_12",
98  "Generic eigenstrains for the computation of the interaction integral.");
99  MooseEnum position_type("Angle Distance", "Distance");
100  params.addParam<MooseEnum>(
101  "position_type",
102  position_type,
103  "The method used to calculate position along crack front. Options are: " +
104  position_type.getRawNames());
105  MooseEnum q_function_type("Geometry Topology", "Geometry");
106  params.addParam<MooseEnum>("q_function_type",
107  q_function_type,
108  "The method used to define the integration domain. Options are: " +
109  q_function_type.getRawNames());
110  params.addParam<bool>(
111  "equivalent_k",
112  false,
113  "Calculate an equivalent K from KI, KII and KIII, assuming self-similar crack growth.");
114  params.addParam<bool>("output_q", true, "Output q");
115  params.addRequiredParam<bool>(
116  "incremental", "Flag to indicate whether an incremental or total model is being used.");
117  params.addParam<std::vector<MaterialPropertyName>>(
118  "eigenstrain_names", "List of eigenstrains applied in the strain calculation");
119  params.addDeprecatedParam<bool>("convert_J_to_K",
120  false,
121  "Convert J-integral to stress intensity factor K.",
122  "This input parameter is deprecated and will be removed soon. "
123  "Use 'integrals = KFromJIntegral' to request output of the "
124  "conversion from the J-integral to stress intensity factors");
125  params.addParam<std::vector<MaterialName>>(
126  "inelastic_models",
127  "The material objects to use to calculate the strain energy rate density.");
128  params.addParam<MaterialPropertyName>("eigenstrain_gradient",
129  "Material defining gradient of eigenstrain tensor");
130  params.addParam<MaterialPropertyName>("body_force", "Material defining body force");
131  params.addParam<bool>("use_automatic_differentiation",
132  false,
133  "Flag to use automatic differentiation (AD) objects when possible");
134  params.addParam<bool>(
135  "used_by_xfem_to_grow_crack",
136  false,
137  "Flag to trigger domainIntregal vector postprocessors to be executed on nonlinear. This "
138  "updates the values in the vector postprocessor which will allow the crack to grow in XFEM "
139  "cutter objects that use the domainIntegral vector postprocssor values as a growth "
140  "criterion.");
141  params.addParam<bool>("output_vpp",
142  true,
143  "Flag to control the vector postprocessor outputs. Select false to "
144  "suppress the redundant csv files for each time step and ring");
145  return params;
146 }
void addDeprecatedParam(const std::string &name, const T &value, const std::string &doc_string, const std::string &deprecation_message)
void addCrackFrontDefinitionParams(InputParameters &params)
void addParam(const std::string &name, const std::initializer_list< typename T::value_type > &value, const std::string &doc_string)
std::string getRawNames() const
void addRequiredParam(const std::string &name, const std::string &doc_string)
static InputParameters validParams()
void addCoupledVar(const std::string &name, const std::string &doc_string)
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void addClassDescription(const std::string &doc_string)
void addParamNamesToGroup(const std::string &space_delim_names, const std::string group_name)

Member Data Documentation

◆ _axis_2d

unsigned int DomainIntegralAction::_axis_2d
protected

Out-of-plane axis for 3D models treated as 2D.

Definition at line 106 of file DomainIntegralAction.h.

Referenced by act().

◆ _blocks

std::vector<SubdomainName> DomainIntegralAction::_blocks
protected

Blocks for which the domain integrals are to be computed.

Definition at line 124 of file DomainIntegralAction.h.

Referenced by act(), and DomainIntegralAction().

◆ _boundary_names

const std::vector<BoundaryName>& DomainIntegralAction::_boundary_names
protected

Boundaries containing the crack front points.

Definition at line 68 of file DomainIntegralAction.h.

Referenced by act(), and calcNumCrackFrontPoints().

◆ _closed_loop

bool DomainIntegralAction::_closed_loop
protected

Indicates whether the crack forms a closed loop.

Definition at line 72 of file DomainIntegralAction.h.

Referenced by act().

◆ _convert_J_to_K

bool DomainIntegralAction::_convert_J_to_K
protected

Whether to convert the J-integral to a stress intensity factor (K) –deprecated.

Definition at line 148 of file DomainIntegralAction.h.

Referenced by DomainIntegralAction().

◆ _crack_direction_vector

const RealVectorValue DomainIntegralAction::_crack_direction_vector
protected

Vector optionally used to prescribe direction of crack extension.

Definition at line 89 of file DomainIntegralAction.h.

Referenced by act().

◆ _crack_direction_vector_end_1

const RealVectorValue DomainIntegralAction::_crack_direction_vector_end_1
protected

Vector optionally used to prescribe direction of crack extension at the 1st end of the crack.

Definition at line 93 of file DomainIntegralAction.h.

Referenced by act().

◆ _crack_direction_vector_end_2

const RealVectorValue DomainIntegralAction::_crack_direction_vector_end_2
protected

Vector optionally used to prescribe direction of crack extension at the 2nd end of the crack.

Definition at line 97 of file DomainIntegralAction.h.

Referenced by act().

◆ _crack_front_points

std::vector<Point> DomainIntegralAction::_crack_front_points
protected

User-defined vector of crack front points.

Definition at line 70 of file DomainIntegralAction.h.

Referenced by act(), calcNumCrackFrontPoints(), and DomainIntegralAction().

◆ _crack_front_points_provider

UserObjectName DomainIntegralAction::_crack_front_points_provider
protected

Name of crack front points provider user object used to optionally define the crack points.

Definition at line 74 of file DomainIntegralAction.h.

Referenced by DomainIntegralAction().

◆ _crack_mouth_boundary_names

std::vector<BoundaryName> DomainIntegralAction::_crack_mouth_boundary_names
protected

Names of boundaries optionally used to define the crack mouth location.

Definition at line 99 of file DomainIntegralAction.h.

Referenced by act(), and DomainIntegralAction().

◆ _direction_method_moose_enum

MooseEnum DomainIntegralAction::_direction_method_moose_enum
protected

Enum used to define the method to compute crack front direction.

Definition at line 83 of file DomainIntegralAction.h.

Referenced by act().

◆ _displacements

std::vector<VariableName> DomainIntegralAction::_displacements
protected

Vector of displacement variables.

Definition at line 126 of file DomainIntegralAction.h.

Referenced by act(), and DomainIntegralAction().

◆ _end_direction_method_moose_enum

MooseEnum DomainIntegralAction::_end_direction_method_moose_enum
protected

Enum used to define the method to compute crack front direction at ends of crack front.

Definition at line 85 of file DomainIntegralAction.h.

Referenced by act().

◆ _family

const std::string DomainIntegralAction::_family
protected

Definition at line 80 of file DomainIntegralAction.h.

Referenced by act().

◆ _fgm_crack

bool DomainIntegralAction::_fgm_crack
protected

Whether the crack lives in a functionally-graded material.

Definition at line 150 of file DomainIntegralAction.h.

Referenced by act(), and DomainIntegralAction().

◆ _functionally_graded_youngs_modulus

MaterialPropertyName DomainIntegralAction::_functionally_graded_youngs_modulus
protected

Material property name for spatially-dependent Youngs modulus for functionally graded materials.

Definition at line 154 of file DomainIntegralAction.h.

Referenced by act(), and DomainIntegralAction().

◆ _functionally_graded_youngs_modulus_crack_dir_gradient

MaterialPropertyName DomainIntegralAction::_functionally_graded_youngs_modulus_crack_dir_gradient
protected

Material property name for the Youngs modulus derivative for functionally graded materials.

Definition at line 152 of file DomainIntegralAction.h.

Referenced by act(), and DomainIntegralAction().

◆ _get_equivalent_k

bool DomainIntegralAction::_get_equivalent_k
protected

Whether to compute the equivalent K from the individual fracture integrals for mixed-mode fracture.

Definition at line 138 of file DomainIntegralAction.h.

Referenced by act(), and DomainIntegralAction().

◆ _has_symmetry_plane

bool DomainIntegralAction::_has_symmetry_plane
protected

Whether the model has a symmetry plane passing through the plane of the crack.

Definition at line 130 of file DomainIntegralAction.h.

Referenced by act().

◆ _have_crack_direction_vector

const bool DomainIntegralAction::_have_crack_direction_vector
protected

Whether the crack direction vector has been provided.

Definition at line 87 of file DomainIntegralAction.h.

Referenced by act().

◆ _have_crack_direction_vector_end_1

const bool DomainIntegralAction::_have_crack_direction_vector_end_1
protected

Whether the crack direction vector at the 1st end of the crack has been provided.

Definition at line 91 of file DomainIntegralAction.h.

Referenced by act().

◆ _have_crack_direction_vector_end_2

const bool DomainIntegralAction::_have_crack_direction_vector_end_2
protected

Whether the crack direction vector at the 2nd end of the crack has been provided.

Definition at line 95 of file DomainIntegralAction.h.

Referenced by act().

◆ _incremental

bool DomainIntegralAction::_incremental
protected

Whether the constitutive models for the mechanics calculations use an incremental form.

Definition at line 146 of file DomainIntegralAction.h.

Referenced by act().

◆ _integrals

std::set<INTEGRAL> DomainIntegralAction::_integrals
protected

Container for enumerations describing the individual integrals computed.

Definition at line 66 of file DomainIntegralAction.h.

Referenced by act(), addRelationshipManagers(), and DomainIntegralAction().

◆ _intersecting_boundary_names

std::vector<BoundaryName> DomainIntegralAction::_intersecting_boundary_names
protected

Names of boundaries optionally used for improved computation of crack extension direction at ends of the crack where it intersects the prescribed boundaries.

Definition at line 102 of file DomainIntegralAction.h.

Referenced by act(), and DomainIntegralAction().

◆ _order

const std::string DomainIntegralAction::_order
protected

Order and family of the AuxVariables optionally created to output the values of q.

Definition at line 79 of file DomainIntegralAction.h.

Referenced by act().

◆ _output_q

bool DomainIntegralAction::_output_q
protected

Whether to ouput the q function as a set of AuxVariables.

Definition at line 142 of file DomainIntegralAction.h.

Referenced by act().

◆ _output_variables

std::vector<VariableName> DomainIntegralAction::_output_variables
protected

List of variables for which values are to be sampled and output at the crack front points.

Definition at line 118 of file DomainIntegralAction.h.

Referenced by act(), and DomainIntegralAction().

◆ _poissons_ratio

Real DomainIntegralAction::_poissons_ratio
protected

Poisson's ratio of material.

Definition at line 120 of file DomainIntegralAction.h.

Referenced by act(), and DomainIntegralAction().

◆ _position_type

MooseEnum DomainIntegralAction::_position_type
protected

How the distance along the crack front is measured (angle or distance)

Definition at line 134 of file DomainIntegralAction.h.

Referenced by act().

◆ _q_function_type

MooseEnum DomainIntegralAction::_q_function_type
protected

How the q function is evaluated (geometric distance from crack front or ring of elements)

Definition at line 136 of file DomainIntegralAction.h.

Referenced by act(), and DomainIntegralAction().

◆ _radius_inner

std::vector<Real> DomainIntegralAction::_radius_inner
protected

Sets of inner and outer radii of the rings used for the domain form of the fracture integrals.

These are defined in corresponding pairs for each ring.

Definition at line 110 of file DomainIntegralAction.h.

Referenced by act(), and DomainIntegralAction().

◆ _radius_outer

std::vector<Real> DomainIntegralAction::_radius_outer
protected

Definition at line 111 of file DomainIntegralAction.h.

Referenced by act(), and DomainIntegralAction().

◆ _ring_first

unsigned int DomainIntegralAction::_ring_first
protected

Number of elements away from the crack tip to inside of inner ring with the topological q function.

Definition at line 114 of file DomainIntegralAction.h.

Referenced by act(), and DomainIntegralAction().

◆ _ring_last

unsigned int DomainIntegralAction::_ring_last
protected

Number of elements away from the crack tip to outside of outer ring with the topological q function.

Definition at line 116 of file DomainIntegralAction.h.

Referenced by act(), and DomainIntegralAction().

◆ _ring_vec

std::vector<unsigned int> DomainIntegralAction::_ring_vec
protected

Vector of ids for the individual rings on which the fracture integral is computed.

Definition at line 144 of file DomainIntegralAction.h.

Referenced by act(), and DomainIntegralAction().

◆ _symmetry_plane

unsigned int DomainIntegralAction::_symmetry_plane
protected

Identifier for which plane is the symmetry plane.

Definition at line 132 of file DomainIntegralAction.h.

Referenced by act().

◆ _temp

VariableName DomainIntegralAction::_temp
protected

Temperature variable.

Definition at line 128 of file DomainIntegralAction.h.

Referenced by act(), and DomainIntegralAction().

◆ _treat_as_2d

bool DomainIntegralAction::_treat_as_2d
protected

Whether fracture computations for a 3D model should be treated as though it were a 2D model.

Definition at line 104 of file DomainIntegralAction.h.

Referenced by act().

◆ _use_ad

const bool DomainIntegralAction::_use_ad
protected

Whether to create automatic differentiation objects from the action.

Definition at line 156 of file DomainIntegralAction.h.

Referenced by act().

◆ _use_crack_front_points_provider

bool DomainIntegralAction::_use_crack_front_points_provider
protected

Whether to use a crack front points provider.

Definition at line 76 of file DomainIntegralAction.h.

Referenced by act(), calcNumCrackFrontPoints(), and DomainIntegralAction().

◆ _use_displaced_mesh

bool DomainIntegralAction::_use_displaced_mesh
protected

Whether to compute the fracture integrals on the displaced mesh.

Definition at line 140 of file DomainIntegralAction.h.

Referenced by act().

◆ _used_by_xfem_to_grow_crack

const bool DomainIntegralAction::_used_by_xfem_to_grow_crack
protected

This determines if fracture integrals should be executed on nonlinear in order to grow the crack when num_xfem_updates in the executioner block is greater than 1.

Definition at line 162 of file DomainIntegralAction.h.

Referenced by act().

◆ _youngs_modulus

Real DomainIntegralAction::_youngs_modulus
protected

Young's modulus of material.

Definition at line 122 of file DomainIntegralAction.h.

Referenced by act(), and DomainIntegralAction().


The documentation for this class was generated from the following files: