Merge pull request #29576 from lynnmunday/nonlocal_cracktip_matl

Nonlocal cracktip material

modules/solid_mechanics/doc/content/source/vectorpostprocessors/CrackFrontNonlocalScalarMaterial.md

@@ -0,0 +1,19 @@
+# CrackFrontNonlocalScalarMaterial
+
+!syntax description /VectorPostprocessors/CrackFrontNonlocalScalarMaterial
+
+## Description
+
+This object computes the average of a scalar material property in the region of the crack front points defined by  [CrackFrontDefinition.md].  The main use case for this `VectorPostprocessor` is to compute an average fracture toughness or $K_c$ at the crack front for use with the `MeshCut2DFractureUserObject` to grow cracks. This allows for spatially varying $K_c$ values defined by a `Material`.
+
+`CrackFrontNonlocalScalarMaterial` computes an average of the material property over a box-shaped domain at each crack tip point that is centered on the crack tip and extends [!param](/VectorPostprocessors/CrackFrontNonlocalScalarMaterial/box_length) in front of the crack tip.  The [!param](/VectorPostprocessors/CrackFrontNonlocalScalarMaterial/box_height) is the dimension normal to the crack face, and [!param](/VectorPostprocessors/CrackFrontNonlocalScalarMaterial/box_width) is the dimension tangential to the crack face.  [!param](/VectorPostprocessors/CrackFrontNonlocalScalarMaterial/box_width) is not used in 2D problems.
+
+In the following input file example, the mesh consists of a 3D plate with a hole in the middle. The CrackFrontDefinition defines crack points around the center line of the hole, `boundary=1001`. This `CrackFrontNonlocalScalarMaterial` averages a material property named `scalar_kcrit` over each 3D box at each crack front point.
+
+!listing crack_front_nonlocal_materials.i block=UserObjects VectorPostprocessors/CrackFrontNonlocalKcrit
+
+!syntax parameters /VectorPostprocessors/CrackFrontNonlocalScalarMaterial
+
+!syntax inputs /VectorPostprocessors/CrackFrontNonlocalScalarMaterial
+
+!syntax children /VectorPostprocessors/CrackFrontNonlocalScalarMaterial

modules/solid_mechanics/doc/content/source/vectorpostprocessors/CrackFrontNonlocalStress.md

@@ -13,9 +13,9 @@ where $\boldsymbol{n}$ is the normal direction vector and $\boldsymbol{\sigma}$
 
 Data produced by this VectorPostprocessor is used in conjunction with the [InteractionIntegral.md] in the XFEM module by the `MeshCut2DFractureUserObject` to grow cracks. The `CrackFrontNonlocalStress` is useful for extending cracks that are approaching free surfaces because the interaction integrals computing `KI` and `KII` are affected when the integration domain intersects the free surface.  `CrackFrontNonlocalStress` computes an average of the stress over a box-shaped domain at each crack tip point that is centered on the crack tip and extends [!param](/VectorPostprocessors/CrackFrontNonlocalStress/box_length) in front of the crack tip.  The [!param](/VectorPostprocessors/CrackFrontNonlocalStress/box_height) is the dimension normal to the crack face, and [!param](/VectorPostprocessors/CrackFrontNonlocalStress/box_width) is the dimension tangential to the crack face.  [!param](/VectorPostprocessors/CrackFrontNonlocalStress/box_width) is not used in 2D problems.  Unlike the other stress integrals, like the [InteractionIntegral.md], that use the [CrackFrontDefinition.md], `CrackFrontNonlocalStress` is not set-up by the [/DomainIntegralAction.md].
 
-In the following input file example, the mesh consists of a 3D plate with a hole in the middle. The CrackFrontDefinition defines crack points around the center line of the hole, `boundary=1001`. This `CrackFrontNonlocalStress` integrates a generic stress field set-up in the input file over the box dimensions shown.
+In the following input file example, the mesh consists of a 3D plate with a hole in the middle. The CrackFrontDefinition defines crack points around the center line of the hole, `boundary=1001`. The `CrackFrontNonlocalStress` averages a generic stress field over the box-shaped region at each crack front point, with dimensions defined in the input file.
 
-!listing crack_front_nonlocal_stress.i block=UserObjects VectorPostprocessors
+!listing crack_front_nonlocal_materials.i block=UserObjects VectorPostprocessors/CrackFrontNonlocalStress
 
 !syntax parameters /VectorPostprocessors/CrackFrontNonlocalStress
 

modules/solid_mechanics/include/vectorpostprocessors/CrackFrontNonlocalMaterialBase.h

@@ -0,0 +1,84 @@
+//* This file is part of the MOOSE framework
+//* https://www.mooseframework.org
+//*
+//* All rights reserved, see COPYRIGHT for full restrictions
+//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
+//*
+//* Licensed under LGPL 2.1, please see LICENSE for details
+//* https://www.gnu.org/licenses/lgpl-2.1.html
+
+#pragma once
+
+#include "ElementVectorPostprocessor.h"
+
+// Forward Declarations
+class CrackFrontDefinition;
+
+/**
+ * Computes the average material property in regions near points provided by the
+ * crack_front_definition vectorpostprocessor.
+ */
+class CrackFrontNonlocalMaterialBase : public ElementVectorPostprocessor
+{
+public:
+  static InputParameters validParams();
+
+  CrackFrontNonlocalMaterialBase(const InputParameters & parameters,
+                                 const std::string & property_name);
+
+  virtual void initialSetup() override;
+  virtual void initialize() override;
+  virtual void execute() override;
+  virtual void finalize() override;
+  virtual void threadJoin(const UserObject & y) override;
+
+protected:
+  /// Material property name from derived class
+  const std::string _property_name;
+  /** dimensions of the box in front of the crack tip that the stress is averaged over
+   * The box is centered in front of the crack tip
+   * _box_length distance box extends in front of the crack tip
+   * _box_width is tangent to the crack tip, centered on the crack point
+   * _box_height is normal to the crack tip, centered on the crack point
+   */
+  ///@{
+  Real _box_length;
+  Real _box_width;
+  Real _box_height;
+  ///@}
+
+  /// used to transform local coordinates to crack front coordinates
+  const CrackFrontDefinition * _crack_front_definition;
+
+  /// Base name of the material system
+  const std::string _base_name;
+
+  // volume being integrated over for each crack front
+  std::vector<Real> _volume;
+
+  /// Vectors computed by this VectorPostprocessor:
+  /// x,y,z coordinates, and position of nodes along crack front, and crack tip average scalar stress
+  ///@{
+  VectorPostprocessorValue & _x;
+  VectorPostprocessorValue & _y;
+  VectorPostprocessorValue & _z;
+  VectorPostprocessorValue & _position;
+  VectorPostprocessorValue & _avg_crack_tip_scalar;
+  ///@}
+
+  /**
+   * Determine whether a point is located within a specified crack front oriented box
+   * @param crack_front_point_index Index of the point on the crack front that the box is based on
+   * @param qp_coord Point object to determine whether it is inside the box
+   * @return  1 if point is within the box, 0 otherwise
+   */
+  Real BoxWeightingFunction(std::size_t crack_front_point_index, const Point & qp_coord) const;
+  /**
+   * Determine whether a point is located within a specified crack front oriented box
+   * @param qp quardature point to get material properties at
+   * @param crack_face_normal normal direction to crack face
+   * @return  Value of material property at this qp and direction
+   */
+  virtual Real getQPCrackFrontScalar(const unsigned int qp,
+                                     const Point crack_face_normal) const = 0;
+};

modules/solid_mechanics/include/vectorpostprocessors/CrackFrontNonlocalScalarMaterial.h

@@ -0,0 +1,31 @@
+//* This file is part of the MOOSE framework
+//* https://www.mooseframework.org
+//*
+//* All rights reserved, see COPYRIGHT for full restrictions
+//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
+//*
+//* Licensed under LGPL 2.1, please see LICENSE for details
+//* https://www.gnu.org/licenses/lgpl-2.1.html
+
+#pragma once
+
+#include "CrackFrontNonlocalMaterialBase.h"
+
+/**
+ * Computes the average material at points provided by the crack_front_definition
+ * vectorpostprocessor.
+ */
+class CrackFrontNonlocalScalarMaterial : public CrackFrontNonlocalMaterialBase
+{
+public:
+  static InputParameters validParams();
+
+  CrackFrontNonlocalScalarMaterial(const InputParameters & parameters);
+
+protected:
+  /// Property that is averaged over the crack front points
+  const MaterialProperty<Real> & _scalar;
+
+  Real getQPCrackFrontScalar(const unsigned int qp,
+                             const Point /*crack_face_normal*/) const override;
+};

modules/solid_mechanics/include/vectorpostprocessors/CrackFrontNonlocalStress.h

@@ -9,66 +9,23 @@
 
 #pragma once
 
-#include "ElementVectorPostprocessor.h"
+#include "CrackFrontNonlocalMaterialBase.h"
+#include "RankTwoTensor.h"
 
-// Forward Declarations
-class CrackFrontDefinition;
 /**
- * This vectorpostprocessor computes the average stress scalar normal to the crack front at points
- * along the crack front.
+ * Computes the average stress magnitude in the direction normal to the crack front at points
+ * provided by the crack_front_definition vectorpostprocessor.
  */
-class CrackFrontNonlocalStress : public ElementVectorPostprocessor
+class CrackFrontNonlocalStress : public CrackFrontNonlocalMaterialBase
 {
 public:
   static InputParameters validParams();
 
   CrackFrontNonlocalStress(const InputParameters & parameters);
 
-  virtual void initialSetup() override;
-  virtual void initialize() override;
-  virtual void execute() override;
-  virtual void finalize() override;
-  virtual void threadJoin(const UserObject & y) override;
-
 protected:
-  /** dimensions of the box in front of the crack tip that the stress is averaged over
-   * The box is centered in front of the crack tip
-   * _box_length distance box extends in front of the crack tip
-   * _box_width is tangent to the crack tip, centered on the crack point
-   * _box_height is normal to the crack tip, centered on the crack point
-   */
-  ///@{
-  Real _box_length;
-  Real _box_width;
-  Real _box_height;
-  ///@}
-
-  /// used to transform local coordinates to crack front coordinates
-  const CrackFrontDefinition * _crack_front_definition;
-
-  /// Base name of the material system
-  const std::string _base_name;
   /// The stress tensor
   const MaterialProperty<RankTwoTensor> & _stress;
 
-  // volume being integrated over for each crack front
-  std::vector<Real> _volume;
-
-  /// Vectors computed by this VectorPostprocessor:
-  /// x,y,z coordinates, and position of nodes along crack front, and crack tip average scalar stress
-  ///@{
-  VectorPostprocessorValue & _x;
-  VectorPostprocessorValue & _y;
-  VectorPostprocessorValue & _z;
-  VectorPostprocessorValue & _position;
-  VectorPostprocessorValue & _avg_crack_tip_stress;
-  ///@}
-
-  /**
-   * Determine whether a point is located within a specified crack front oriented box
-   * @param crack_front_point_index Index of the point on the crack front that the box is based on
-   * @param qp_coord Point object to determine whether it is inside the box
-   * @return  1 if point is within the box, 0 otherwise
-   */
-  Real BoxWeightingFunction(std::size_t crack_front_point_index, const Point & qp_coord) const;
+  Real getQPCrackFrontScalar(const unsigned int qp, const Point crack_face_normal) const override;
 };

modules/solid_mechanics/src/vectorpostprocessors/CrackFrontNonlocalMaterialBase.C

@@ -0,0 +1,148 @@
+//* This file is part of the MOOSE framework
+//* https://www.mooseframework.org
+//*
+//* All rights reserved, see COPYRIGHT for full restrictions
+//* https://github.com/idaholab/moose/blob/master/COPYRIGHT
+//*
+//* Licensed under LGPL 2.1, please see LICENSE for details
+//* https://www.gnu.org/licenses/lgpl-2.1.html
+
+#include "CrackFrontNonlocalMaterialBase.h"
+#include "CrackFrontDefinition.h"
+
+InputParameters
+CrackFrontNonlocalMaterialBase::validParams()
+{
+  InputParameters params = ElementVectorPostprocessor::validParams();
+  params.addRequiredParam<UserObjectName>("crack_front_definition",
+                                          "The CrackFrontDefinition user object name");
+  params.addRequiredParam<Real>(
+      "box_length", "Dimension of property-averaging box in direction of crack extension.");
+  params.addRequiredParam<Real>(
+      "box_height", "Dimension of property-averaging box in direction normal to crack.");
+  params.addParam<Real>("box_width", 1.0, "Distance tangent to front of crack front.");
+  params.addParam<std::string>("base_name",
+                               "Optional parameter that allows the user to define "
+                               "multiple mechanics material systems on the same "
+                               "block, i.e. for multiple phases");
+  params.set<bool>("use_displaced_mesh") = false;
+  params.set<ExecFlagEnum>("execute_on") = {EXEC_TIMESTEP_BEGIN};
+  params.addClassDescription("Computes the average material property at a crack front.");
+  return params;
+}
+
+CrackFrontNonlocalMaterialBase::CrackFrontNonlocalMaterialBase(const InputParameters & parameters,
+                                                               const std::string & property_name)
+  : ElementVectorPostprocessor(parameters),
+    _property_name(property_name),
+    _box_length(getParam<Real>("box_length")),
+    _box_width(getParam<Real>("box_width")),
+    _box_height(getParam<Real>("box_height")),
+    _base_name(isParamValid("base_name") ? getParam<std::string>("base_name") + "_" : ""),
+    _x(declareVector("x")),
+    _y(declareVector("y")),
+    _z(declareVector("z")),
+    _position(declareVector("id")),
+    // get the property name instead of materialname
+    _avg_crack_tip_scalar(declareVector("crack_tip_" + _base_name + _property_name))
+{
+  if (_mesh.dimension() == 3 && !isParamSetByUser("box_width"))
+    paramError("box_width", "Must define box_width in 3D problems.");
+  // add user object dependencies by name (the UOs do not need to exist yet for this)
+  _depend_uo.insert(getParam<UserObjectName>("crack_front_definition"));
+}
+
+void
+CrackFrontNonlocalMaterialBase::initialSetup()
+{
+  const auto uo_name = getParam<UserObjectName>("crack_front_definition");
+  _crack_front_definition =
+      &(getUserObjectByName<CrackFrontDefinition>(uo_name, /*is_dependency = */ false));
+}
+
+void
+CrackFrontNonlocalMaterialBase::initialize()
+{
+  std::size_t num_pts = _crack_front_definition->getNumCrackFrontPoints();
+
+  _volume.assign(num_pts, 0.0);
+  _x.assign(num_pts, 0.0);
+  _y.assign(num_pts, 0.0);
+  _z.assign(num_pts, 0.0);
+  _position.assign(num_pts, 0.0);
+  _avg_crack_tip_scalar.assign(num_pts, 0.0);
+}
+
+void
+CrackFrontNonlocalMaterialBase::execute()
+{
+  // icfp crack front point index
+  for (const auto icfp : index_range(_avg_crack_tip_scalar))
+  {
+    Point crack_front_normal = _crack_front_definition->getCrackFrontNormal(icfp);
+    for (unsigned int qp = 0; qp < _qrule->n_points(); qp++)
+    {
+      Real q = BoxWeightingFunction(icfp, _q_point[qp]);
+      if (q == 0)
+        continue;
+
+      Real scalar = getQPCrackFrontScalar(qp, crack_front_normal);
+      _avg_crack_tip_scalar[icfp] += _JxW[qp] * _coord[qp] * scalar * q;
+      _volume[icfp] += _JxW[qp] * _coord[qp] * q;
+    }
+  }
+}
+
+void
+CrackFrontNonlocalMaterialBase::finalize()
+{
+  gatherSum(_avg_crack_tip_scalar);
+  gatherSum(_volume);
+  for (const auto icfp : index_range(_avg_crack_tip_scalar))
+  {
+    if (_volume[icfp] != 0)
+      _avg_crack_tip_scalar[icfp] = _avg_crack_tip_scalar[icfp] / _volume[icfp];
+    else
+      _avg_crack_tip_scalar[icfp] = 0;
+
+    const auto cfp = _crack_front_definition->getCrackFrontPoint(icfp);
+    _x[icfp] = (*cfp)(0);
+    _y[icfp] = (*cfp)(1);
+    _z[icfp] = (*cfp)(2);
+    _position[icfp] = _crack_front_definition->getDistanceAlongFront(icfp);
+  }
+}
+
+void
+CrackFrontNonlocalMaterialBase::threadJoin(const UserObject & y)
+{
+  const auto & uo = static_cast<const CrackFrontNonlocalMaterialBase &>(y);
+  for (const auto i : index_range(_avg_crack_tip_scalar))
+  {
+    _volume[i] += uo._volume[i];
+    _avg_crack_tip_scalar[i] += uo._avg_crack_tip_scalar[i];
+  }
+}
+
+Real
+CrackFrontNonlocalMaterialBase::BoxWeightingFunction(std::size_t crack_front_point_index,
+                                                     const Point & qp_coord) const
+{
+  const Point * cf_pt = _crack_front_definition->getCrackFrontPoint(crack_front_point_index);
+  RealVectorValue crack_node_to_current_node = qp_coord - *cf_pt;
+
+  // crackfront coordinates are:
+  // crack_node_to_current_node_rot[0]= crack direction
+  // crack_node_to_current_node_rot[1]= normal to crack face
+  // crack_node_to_current_node_rot[2]= tangent to crack face along crack front (not used in 2D)
+  RealVectorValue crack_node_to_current_node_rot =
+      _crack_front_definition->rotateToCrackFrontCoords(crack_node_to_current_node,
+                                                        crack_front_point_index);
+  if ((crack_node_to_current_node_rot(0) > 0) &&
+      (crack_node_to_current_node_rot(0) <= _box_length) &&
+      (std::abs(crack_node_to_current_node_rot(1)) <= _box_height / 2) &&
+      (std::abs(crack_node_to_current_node_rot(2)) <= _box_width / 2))
+    return 1.0;
+
+  return 0.0;
+}