diff --git a/modules/solid_mechanics/doc/content/source/actions/DomainIntegralAction.md b/modules/solid_mechanics/doc/content/source/actions/DomainIntegralAction.md index b82e6962281f..b9b7ec7bdd11 100644 --- a/modules/solid_mechanics/doc/content/source/actions/DomainIntegralAction.md +++ b/modules/solid_mechanics/doc/content/source/actions/DomainIntegralAction.md @@ -8,7 +8,7 @@ The `DomainIntegral` action is used to set up all of the objects used in computi Meshed crack: The crack can be explicitly included by creating a mesh with a topology that conforms to the crack. The location of the crack tip is provided to the code by defining a nodeset that includes all nodes in the finite element mesh that are located along the crack front. For 2D analyses, this nodeset would only contain a single node at the crack tip. For 3D analyses, the mesh connectivity is used to construct a set of line segments that connect these nodes, and this is used to order the crack nodes. -XFEM: Rather than defining the topology of the crack through the mesh, XFEM can be used to cut the mesh. In this case, a set of points, which does not need to conform to points in the mesh, must be provided by the user, and is used to define the location of the crack for computing the fracture integrals. Fracture integrals are computed at the locations of these points, in the order provided by the user. +XFEM: Rather than defining the topology of the crack through the mesh, XFEM can be used to cut the mesh. In this case, a set of points, which does not need to conform to points in the mesh, must be provided by the user, and is used to define the location of the crack for computing the fracture integrals. Fracture integrals are computed at the locations of these points, in the order provided by the user. If the `DomainIntegral` is setting up an [InteractionIntegral.md] being used as a growth criterion for the XFEM crack, i.e. [!param](/DomainIntegral/DomainIntegralAction/used_by_xfem_to_grow_crack) = true, the [CrackFrontDefinition.md] must be executed before the `DomainIntegral` vectorpostprocessors. This is done by default in the `DomainIntegral` action which sets the the `CrackFrontDefinition` [!param](/UserObjects/CrackFrontDefinition/execution_order_group) to -1. ## Theory diff --git a/modules/solid_mechanics/src/actions/DomainIntegralAction.C b/modules/solid_mechanics/src/actions/DomainIntegralAction.C index c9a4128c73ed..9498aaf29d7d 100644 --- a/modules/solid_mechanics/src/actions/DomainIntegralAction.C +++ b/modules/solid_mechanics/src/actions/DomainIntegralAction.C @@ -374,7 +374,11 @@ DomainIntegralAction::act() InputParameters params = _factory.getValidParams(uo_type_name); if (_use_crack_front_points_provider && _used_by_xfem_to_grow_crack) + { params.set("execute_on") = {EXEC_INITIAL, EXEC_TIMESTEP_END, EXEC_NONLINEAR}; + // The CrackFrontDefinition updates the vpps and MUST execute before them + params.set("execution_order_group") = -1; + } else params.set("execute_on") = {EXEC_INITIAL, EXEC_TIMESTEP_END}; @@ -735,7 +739,12 @@ DomainIntegralAction::act() params.set>("outputs") = {"none"}; if (_use_crack_front_points_provider && _used_by_xfem_to_grow_crack) + { + // The CrackFrontDefinition updates this vpp and MUST execute before it + // This is enforced by setting the execution_order_group = -1 for the CrackFrontDefinition + // The CrackFrontDefinition must execute on nonlinear to update with xfem updates params.set("execute_on") = {EXEC_TIMESTEP_END, EXEC_NONLINEAR}; + } else params.set("execute_on") = {EXEC_TIMESTEP_END}; diff --git a/modules/xfem/doc/content/source/userobjects/MeshCut2DRankTwoTensorNucleation.md b/modules/xfem/doc/content/source/userobjects/MeshCut2DRankTwoTensorNucleation.md index 7570ca92403c..10e4548ca2ce 100644 --- a/modules/xfem/doc/content/source/userobjects/MeshCut2DRankTwoTensorNucleation.md +++ b/modules/xfem/doc/content/source/userobjects/MeshCut2DRankTwoTensorNucleation.md @@ -11,14 +11,14 @@ such as a principal stress or a component of stress. If the scalar exceeds a threshold specified by [!param](/UserObjects/MeshCut2DRankTwoTensorNucleation/nucleation_threshold), a line segment with a length specified by [!param](/UserObjects/MeshCut2DRankTwoTensorNucleation/nucleation_length) will be inserted into the [MeshCut2DFractureUserObject.md] cutter mesh. -Cracks are only allowed to initiate from elements on boundaries specified by [!param](/UserObjects/MeshCut2DRankTwoTensorNucleation/initiate_on_boundary). +Cracks will nucleate on the domain interior or boundary. Crack nucleation can be limited to elements on boundaries specified by [!param](/UserObjects/MeshCut2DRankTwoTensorNucleation/initiate_on_boundary). Once the nucleation criterion is reached, a line segment of the specified length is inserted into the cutter mesh, centered on the element centroid it nucleates from. The direction of the nucleated crack is normal to the direction returned by the `RankTwoTensor` scalar. For example, `MaxInPlanePrincipal` returns the direction of the maximum in-plane principal component and the crack direction is normal to this. -The nucleation length should be at least the length of the element it nucleates in so that the nucleated crack will completely cut the element. +It is an error for [!param](/UserObjects/MeshCut2DRankTwoTensorNucleation/nucleation_length) to be smaller than the length of the element the crack nucleates in, unless the user sets `nucleate_across_full_element=true`. A crack will only be nucleated if it is at least a distance specified by [!param](/UserObjects/MeshCut2DRankTwoTensorNucleation/nucleation_radius) away from existing or nucleated cracks. -If multiple cracks nucleate in the same xfem update and are within the specified `nucleation_radius`, the crack nucleated from the element with the lowest id will be retained and no other cracks within the `nucleation_radius` will be nucleated. +If multiple cracks nucleate in the same xfem update and are within the specified `nucleation_radius`, the crack nucleated from the element with the lowest id will be retained and no other cracks within the `nucleation_radius` will be nucleated. The `nucleation_radius` is a good way to limit the number of cracks that nucleate within a region, which can help limit problems related to excessive numbers of XFEM cuts in a single element. `MeshCut2DRankTwoTensorNucleation` copies several features available in the [XFEMRankTwoTensorMarkerUserObject.md]. These include the nucleation threshold being provided as a coupled variable and the computation of the maximum value of the scalar diff --git a/modules/xfem/include/userobjects/MeshCut2DRankTwoTensorNucleation.h b/modules/xfem/include/userobjects/MeshCut2DRankTwoTensorNucleation.h index 9e52fe25c553..dbd56685a65a 100644 --- a/modules/xfem/include/userobjects/MeshCut2DRankTwoTensorNucleation.h +++ b/modules/xfem/include/userobjects/MeshCut2DRankTwoTensorNucleation.h @@ -21,6 +21,18 @@ class MeshCut2DRankTwoTensorNucleation : public MeshCut2DNucleationBase virtual ~MeshCut2DRankTwoTensorNucleation() {} protected: + /// scaling factor to extend the nucleated crack beyond the element edges. + const Real _edge_extension_factor; + + /// should element be cut if the nucleation_length is smaller than the element length. + const bool _nucleate_across_full_element; + + /// is the nucleation length provided in the input file. + const bool _is_nucleation_length_provided; + + /// Length of the nucleated cracks + const Real _nucleation_length; + /// The tensor from which the scalar quantity used as a nucleating criterion is extracted const MaterialProperty & _tensor; @@ -34,13 +46,20 @@ class MeshCut2DRankTwoTensorNucleation : public MeshCut2DNucleationBase const Point _point1; const Point _point2; - /// Length of crack to be nucleated - const Real _nucleation_length; - /// Transformed Jacobian weights const MooseArray & _JxW; const MooseArray & _coord; virtual bool doesElementCrack(std::pair & cutterElemNodes) override; + + // FIXME Lynn Copy from TraceRayTools.C in rayTracing module. Remove once this function migrates + // to libmesh. + bool lineLineIntersect2D(const Point & start, + const Point & direction, + const Real length, + const Point & v0, + const Point & v1, + Point & intersection_point, + Real & intersection_distance); }; diff --git a/modules/xfem/src/userobjects/MeshCut2DNucleationBase.C b/modules/xfem/src/userobjects/MeshCut2DNucleationBase.C index 7d709b28bf12..9dcc23bbd911 100644 --- a/modules/xfem/src/userobjects/MeshCut2DNucleationBase.C +++ b/modules/xfem/src/userobjects/MeshCut2DNucleationBase.C @@ -20,7 +20,7 @@ InputParameters MeshCut2DNucleationBase::validParams() { InputParameters params = ElementUserObject::validParams(); - params.addRequiredParam>( + params.addParam>( "initiate_on_boundary", "Cracks can only initiate on elements adjacent to specified boundaries."); params.addParam("nucleation_radius", @@ -56,7 +56,7 @@ MeshCut2DNucleationBase::MeshCut2DNucleationBase(const InputParameters & paramet { std::vector initiation_boundary_names = getParam>("initiate_on_boundary"); - _initiation_boundary_ids = _mesh.getBoundaryIDs(initiation_boundary_names, true); + _initiation_boundary_ids = _mesh.getBoundaryIDs(initiation_boundary_names); } } @@ -74,21 +74,36 @@ MeshCut2DNucleationBase::execute() if (_current_elem->processor_id() != processor_id()) return; - bool isOnBoundary = false; unsigned int current_eid = _current_elem->id(); std::map>::iterator mit; mit = _nucleated_elems.find(current_eid); - for (unsigned int i = 0; i < _initiation_boundary_ids.size(); ++i) - if (_mesh.isBoundaryElem(current_eid, _initiation_boundary_ids[i])) - isOnBoundary = true; - // This does not currently allow for nucleation in an element that is already cut - if (!is_cut && isOnBoundary && doesElementCrack(cutterElemNodes)) + bool cut_element = false; + if (_initiation_boundary_ids.empty()) + { + // nucleating crack in bulk + // This does not currently allow for nucleation in an element that is already cut + cut_element = (!is_cut && doesElementCrack(cutterElemNodes)); + } + else + { + // nucleating crack on specified boundary + bool isOnBoundary = false; + for (unsigned int i = 0; i < _initiation_boundary_ids.size(); ++i) + if (_mesh.isBoundaryElem(current_eid, _initiation_boundary_ids[i])) + { + isOnBoundary = true; + break; + } + // This does not currently allow for nucleation in an element that is already cut + cut_element = (!is_cut && isOnBoundary && doesElementCrack(cutterElemNodes)); + } + + if (cut_element) { if (mit != _nucleated_elems.end()) - { mooseError("ERROR: element ", current_eid, " already marked for crack nucleation."); - } + _nucleated_elems[current_eid] = cutterElemNodes; } } @@ -97,7 +112,6 @@ void MeshCut2DNucleationBase::threadJoin(const UserObject & y) { const auto & xmuo = static_cast(y); - for (std::map>::const_iterator mit = xmuo._nucleated_elems.begin(); mit != xmuo._nucleated_elems.end(); diff --git a/modules/xfem/src/userobjects/MeshCut2DRankTwoTensorNucleation.C b/modules/xfem/src/userobjects/MeshCut2DRankTwoTensorNucleation.C index 17b21d3ed0ff..27af49a05387 100644 --- a/modules/xfem/src/userobjects/MeshCut2DRankTwoTensorNucleation.C +++ b/modules/xfem/src/userobjects/MeshCut2DRankTwoTensorNucleation.C @@ -9,11 +9,12 @@ #include "MeshCut2DRankTwoTensorNucleation.h" +#include "MooseError.h" +#include "MooseTypes.h" #include "libmesh/quadrature.h" #include "RankTwoTensor.h" #include "RankTwoScalarTools.h" #include "Assembly.h" -#include registerMooseObject("XFEMApp", MeshCut2DRankTwoTensorNucleation); @@ -23,6 +24,18 @@ MeshCut2DRankTwoTensorNucleation::validParams() InputParameters params = MeshCut2DNucleationBase::validParams(); params.addClassDescription( "Nucleate a crack in MeshCut2D UO based on a scalar extracted from a RankTwoTensor"); + params.addRangeCheckedParam( + "edge_extension_factor", + 1e-5, + "edge_extension_factor >= 0", + "Factor by which the length of cracks that have been adjusted to coincide with element edges " + "are increased to avoid missing an intersection with the edge due to numerical tolerance " + "issues."); + params.addRangeCheckedParam( + "nucleation_length", + "nucleation_length >= 0", + "Size of crack to nucleate. Must be larger than the length of the element in which the " + "crack is nucleated, unless 'nucleate_across_full_element' is set to 'true'."); params.addParam( "scalar_type", RankTwoScalarTools::scalarOptions(), @@ -31,7 +44,6 @@ MeshCut2DRankTwoTensorNucleation::validParams() params.addRequiredCoupledVar( "nucleation_threshold", "Threshold for the scalar quantity of the RankTwoTensor to nucleate new cracks"); - params.addRequiredParam("nucleation_length", "Nucleated crack length."); params.addParam( "point1", Point(0, 0, 0), @@ -40,18 +52,29 @@ MeshCut2DRankTwoTensorNucleation::validParams() "point2", Point(0, 1, 0), "End point for axis used to calculate some cylindrical material tensor quantities"); + params.addParam( + "nucleate_across_full_element", + false, + "Controls the behavior of nucleating cracks if 'nucleation_length' is smaller than the " + "length required to travserse an element with a nucleating crack. If this is set to 'false', " + "that condition will result in an error, but if set to 'true', a crack with the length " + "needed " + "to traverse the element will be inserted."); return params; } MeshCut2DRankTwoTensorNucleation::MeshCut2DRankTwoTensorNucleation( const InputParameters & parameters) : MeshCut2DNucleationBase(parameters), + _edge_extension_factor(getParam("edge_extension_factor")), + _nucleate_across_full_element(getParam("nucleate_across_full_element")), + _is_nucleation_length_provided(isParamValid("nucleation_length")), + _nucleation_length(_is_nucleation_length_provided ? getParam("nucleation_length") : 0), _tensor(getMaterialProperty(getParam("tensor"))), _nucleation_threshold(coupledValue("nucleation_threshold")), _scalar_type(getParam("scalar_type")), _point1(parameters.get("point1")), _point2(parameters.get("point2")), - _nucleation_length(getParam("nucleation_length")), _JxW(_assembly.JxW()), _coord(_assembly.coordTransformation()) { @@ -84,13 +107,125 @@ MeshCut2DRankTwoTensorNucleation::doesElementCrack( if (tensor_quantity > average_threshold) { - const Point elem_center(_current_elem->vertex_average()); does_it_crack = true; + + const Point elem_center(_current_elem->vertex_average()); Point crack_dir = point_dir.cross(Point(0, 0, 1)); - RealVectorValue point_0 = elem_center - _nucleation_length / 2.0 * crack_dir.unit(); - RealVectorValue point_1 = elem_center + _nucleation_length / 2.0 * crack_dir.unit(); + + // get max element length to use for ray length + std::unique_ptr edge; + Real circumference = 0; + for (const auto e : _current_elem->edge_index_range()) + { + _current_elem->build_edge_ptr(edge, e); + circumference += edge->length(0, 1); + } + + // Temproraries for doing things below + Real intersection_distance; + + Point point_0; + Point point_1; + bool is_point_0_on_external_boundary = false; + bool is_point_1_on_external_boundary = false; + for (const auto e : _current_elem->edge_index_range()) + { + _current_elem->build_edge_ptr(edge, e); + const auto & edge0 = edge->point(0); + const auto & edge1 = edge->point(1); + if (lineLineIntersect2D( + elem_center, -crack_dir, circumference, edge0, edge1, point_0, intersection_distance)) + is_point_0_on_external_boundary = (_current_elem->neighbor_ptr(e) == nullptr); + else if (lineLineIntersect2D(elem_center, + crack_dir, + circumference, + edge0, + edge1, + point_1, + intersection_distance)) + is_point_1_on_external_boundary = (_current_elem->neighbor_ptr(e) == nullptr); + } + + // traverse_length is the length of a cut that goes across a single elment + // extend_length is the amount that needs to be added to each side of the traverse_length cut to + // make its length equal to the nucleation length + // We also add a small amount to the crack length to make sure it fully cuts the element edge. + // This factor is based on the element length. + Real traverse_length = (point_0 - point_1).norm(); + Real extend_length = (_nucleation_length - traverse_length) / 2.0; + if (_nucleation_length < traverse_length) + { + if (_is_nucleation_length_provided && !_nucleate_across_full_element) + mooseError( + "Trying to nucleate crack smaller than element length, increase nucleation_length.\n " + "Location of crack being nucleated: ", + point_0, + "\n nucleation_length: ", + _nucleation_length, + "\n Length to traverse element: ", + traverse_length, + "\n This error can be suppressed by setting nucleate_across_full_element=True which " + "will nucleate a crack the size of the element."); + else + //_edge_extension_factor is used to make sure cut will cut both edges of the element. + extend_length = (traverse_length * _edge_extension_factor) / 2.0; + } + + // First create a cut assuming bulk nucleation + point_0 = point_0 - extend_length * crack_dir.unit(); + point_1 = point_1 + extend_length * crack_dir.unit(); + + // modify edges of cut should go so that they are on the edge of the domain and have a length + // equal to nucleation_length + if (is_point_0_on_external_boundary) + { + point_0 = point_0 - (traverse_length * _edge_extension_factor) / 2.0 * crack_dir.unit(); + } + else if (is_point_1_on_external_boundary) + { + point_1 = point_1 + (traverse_length * _edge_extension_factor) / 2.0 * crack_dir.unit(); + } cutterElemNodes = {point_0, point_1}; } return does_it_crack; } + +bool +MeshCut2DRankTwoTensorNucleation::lineLineIntersect2D(const Point & start, + const Point & direction, + const Real length, + const Point & v0, + const Point & v1, + Point & intersection_point, + Real & intersection_distance) +{ + /// The standard "looser" tolerance to use in ray tracing when having difficulty finding intersection + const Real TRACE_TOLERANCE = 1e-5; + + const auto r = direction * length; + const auto s = v1 - v0; + const auto rxs = r(0) * s(1) - r(1) * s(0); + + // Lines are parallel or colinear + if (std::abs(rxs) < TRACE_TOLERANCE) + return false; + + const auto v0mu0 = v0 - start; + const auto t = (v0mu0(0) * s(1) - v0mu0(1) * s(0)) / rxs; + if (0 >= t + TRACE_TOLERANCE || t - TRACE_TOLERANCE > 1.0) + { + return false; + } + + const auto u = (v0mu0(0) * r(1) - v0mu0(1) * r(0)) / rxs; + if (0 < u + TRACE_TOLERANCE && u - TRACE_TOLERANCE <= 1.0) + { + intersection_point = start + r * t; + intersection_distance = t * length; + return true; + } + + // Not parallel, but don't intersect + return false; +} diff --git a/modules/xfem/src/userobjects/MeshCut2DUserObjectBase.C b/modules/xfem/src/userobjects/MeshCut2DUserObjectBase.C index 5c8463dccd46..534539a6922a 100644 --- a/modules/xfem/src/userobjects/MeshCut2DUserObjectBase.C +++ b/modules/xfem/src/userobjects/MeshCut2DUserObjectBase.C @@ -172,6 +172,7 @@ MeshCut2DUserObjectBase::getCrackFrontPoints(unsigned int number_crack_front_poi " does not match the number of nodes given in " "_original_and_current_front_node_ids=" + Moose::stringify(_original_and_current_front_node_ids.size())); + for (unsigned int i = 0; i < number_crack_front_points; ++i) { dof_id_type id = _original_and_current_front_node_ids[i].second; @@ -213,22 +214,33 @@ MeshCut2DUserObjectBase::getCrackPlaneNormals(unsigned int number_crack_front_po bool found_it0 = (it0 != _original_and_current_front_node_ids.end()); bool found_it1 = (it1 != _original_and_current_front_node_ids.end()); - if (found_it0 || found_it1) + + // Newly nucleated crack elements can have one normal if they are on the edge OR + // two normals if they are in the bulk. + if (found_it0) { + Point end_pt, connecting_pt; + + end_pt = elem->node_ref(0); + connecting_pt = elem->node_ref(1); + id = it0->first; // sort by original crack front node ids + + Point fracture_dir = end_pt - connecting_pt; + // The crack normal is orthogonal to the crack extension direction (fracture_dir), + // and is defined in this implementation as the cross product of the direction of crack + // extension with the tangent direction, which is always (0, 0, 1) in 2D. + RealVectorValue normal_dir{fracture_dir(1), -fracture_dir(0), 0}; + normal_dir /= normal_dir.norm(); + crack_plane_normals.push_back(std::make_pair(id, normal_dir)); + } + if (found_it1) + { Point end_pt, connecting_pt; - if (found_it0) - { - end_pt = elem->node_ref(0); - connecting_pt = elem->node_ref(1); - id = it0->first; // sort by original crack front node ids - } - else - { - end_pt = elem->node_ref(1); - connecting_pt = elem->node_ref(0); - id = it1->first; // sort by original crack front node ids - } + + end_pt = elem->node_ref(1); + connecting_pt = elem->node_ref(0); + id = it1->first; // sort by original crack front node ids Point fracture_dir = end_pt - connecting_pt; // The crack normal is orthogonal to the crack extension direction (fracture_dir), @@ -361,14 +373,18 @@ MeshCut2DUserObjectBase::addNucleatedCracksToMesh() } _cutter_mesh->add_elem(new_elem); // now add the nucleated nodes to the crack id data struct + // edge nucleated cracks will add one node to _original_and_current_front_node_ids + // bulk nucleated cracks will add two nodes to _original_and_current_front_node_ids Point & point_0 = *node_0; const Elem * crack_front_elem_0 = (*pl)(point_0); if (crack_front_elem_0 != NULL) _original_and_current_front_node_ids.push_back(std::make_pair(node_id_0, node_id_0)); + Point & point_1 = *node_1; const Elem * crack_front_elem_1 = (*pl)(point_1); if (crack_front_elem_1 != NULL) _original_and_current_front_node_ids.push_back(std::make_pair(node_id_1, node_id_1)); + _is_mesh_modified = true; } _cutter_mesh->prepare_for_use(); diff --git a/modules/xfem/test/tests/nucleation_uo/gold/nucleate_AllEdgeCracks_out_II_KI_1_0003.csv b/modules/xfem/test/tests/nucleation_uo/gold/nucleate_AllEdgeCracks_out_II_KI_1_0003.csv index c7b807122fe2..ac9fbbb48400 100644 --- a/modules/xfem/test/tests/nucleation_uo/gold/nucleate_AllEdgeCracks_out_II_KI_1_0003.csv +++ b/modules/xfem/test/tests/nucleation_uo/gold/nucleate_AllEdgeCracks_out_II_KI_1_0003.csv @@ -1,3 +1,3 @@ II_KI_1,id,x,y,z --34.582669725973,0,0.14492706922106,1.9185539902516,0 -106.79227914941,0,0.14992426436198,1.5538911941249,0 +-32.53710316762,0,0.100001,1.9334364218096,0 +102.02665838426,0,0.100001,1.551947110631,0 diff --git a/modules/xfem/test/tests/nucleation_uo/gold/nucleate_AllEdgeCracks_out_II_KI_1_0005.csv b/modules/xfem/test/tests/nucleation_uo/gold/nucleate_AllEdgeCracks_out_II_KI_1_0005.csv index 1970714719a6..ad77bd29614f 100644 --- a/modules/xfem/test/tests/nucleation_uo/gold/nucleate_AllEdgeCracks_out_II_KI_1_0005.csv +++ b/modules/xfem/test/tests/nucleation_uo/gold/nucleate_AllEdgeCracks_out_II_KI_1_0005.csv @@ -1,7 +1,7 @@ II_KI_1,id,x,y,z --64.017070197101,0,0.14492706922106,1.9185539902516,0 -163.46849638768,0,0.14992426436198,1.5538911941249,0 -130.48886106958,0,0.14999859520155,0.35053005444629,0 -143.73418947688,0,0.14999988262991,0.65015321228526,0 -144.57385251999,0,0.14999971353134,0.94976063907878,0 -153.50224057472,0,0.14997193844276,1.2476311344495,0 +-59.72916137921,0,0.100001,1.9334364218096,0 +154.91870842731,0,0.100001,1.551947110631,0 +121.1812753263,0,0.100001,0.35026503624691,0 +133.23428222264,0,0.100001,0.65007660776467,0 +134.09608007163,0,0.100001,0.94988031680293,0 +144.31451648654,0,0.100001,1.2488152110658,0 diff --git a/modules/xfem/test/tests/nucleation_uo/gold/nucleate_bulkCrack_out_II_KI_1_0006.csv b/modules/xfem/test/tests/nucleation_uo/gold/nucleate_bulkCrack_out_II_KI_1_0006.csv new file mode 100644 index 000000000000..a20324f62342 --- /dev/null +++ b/modules/xfem/test/tests/nucleation_uo/gold/nucleate_bulkCrack_out_II_KI_1_0006.csv @@ -0,0 +1,7 @@ +II_KI_1,id,x,y,z +133.82705581083,0,0.095,0.44999999999998,0 +181.57609474945,0,0.205,0.45000000000002,0 +161.29728505566,0,0.695,0.24999999999994,0 +131.21463396511,0,0.805,0.25000000000006,0 +159.46251098133,0,0.895,0.64999999999997,0 +159.54460625526,0,0.105,1.15,0 diff --git a/modules/xfem/test/tests/nucleation_uo/gold/nucleate_edge_bulk_crack_2d_out_II_KI_1_0011.csv b/modules/xfem/test/tests/nucleation_uo/gold/nucleate_edge_bulk_crack_2d_out_II_KI_1_0011.csv new file mode 100644 index 000000000000..5fcdd39e4032 --- /dev/null +++ b/modules/xfem/test/tests/nucleation_uo/gold/nucleate_edge_bulk_crack_2d_out_II_KI_1_0011.csv @@ -0,0 +1,7 @@ +II_KI_1,id,x,y,z +178.3108853287,0,-1.7962636005232,0.59956349097406,0 +113.20564319875,0,-0.089367205247895,0.56725734686736,0 +116.03231568824,0,0.10299701379111,0.62802201693232,0 +124.41900686536,0,1.9091882432968,0.5966936906918,0 +219.03669652501,0,-2.5200487218268,0.51122488784312,0 +217.79168695061,0,2.519879839307,0.51043596218937,0 diff --git a/modules/xfem/test/tests/nucleation_uo/gold/nucleate_edge_crack_2d_out_II_KI_1_0004.csv b/modules/xfem/test/tests/nucleation_uo/gold/nucleate_edge_crack_2d_out_II_KI_1_0004.csv deleted file mode 100644 index f86f3efb8588..000000000000 --- a/modules/xfem/test/tests/nucleation_uo/gold/nucleate_edge_crack_2d_out_II_KI_1_0004.csv +++ /dev/null @@ -1,3 +0,0 @@ -II_KI_1,id,x,y,z -217.62685814112,0,0.25887298477258,0.47165266238128,0 -216.82358588287,0,1.741039702803,0.47077842129492,0 diff --git a/modules/xfem/test/tests/nucleation_uo/nucleate_2edge_cracks_2d.i b/modules/xfem/test/tests/nucleation_uo/nucleate_2edge_cracks_2d.i index 009d4808d00b..49cd6884bc00 100644 --- a/modules/xfem/test/tests/nucleation_uo/nucleate_2edge_cracks_2d.i +++ b/modules/xfem/test/tests/nucleation_uo/nucleate_2edge_cracks_2d.i @@ -109,7 +109,7 @@ [] [] -[Modules/TensorMechanics/Master] +[Physics/SolidMechanics/QuasiStatic] [./all] strain = FINITE planar_formulation = plane_strain diff --git a/modules/xfem/test/tests/nucleation_uo/nucleate_AllEdgeCracks.i b/modules/xfem/test/tests/nucleation_uo/nucleate_AllEdgeCracks.i index 69fa021ea4fe..1846ad842538 100644 --- a/modules/xfem/test/tests/nucleation_uo/nucleate_AllEdgeCracks.i +++ b/modules/xfem/test/tests/nucleation_uo/nucleate_AllEdgeCracks.i @@ -52,8 +52,9 @@ scalar_type = MaxPrincipal nucleation_threshold = 180 initiate_on_boundary = 'left' - nucleation_length = .2 nucleation_radius = .21 + edge_extension_factor = 2e-5 + nucleation_length = 0.11 [] [cut_mesh2] type = MeshCut2DFractureUserObject @@ -64,7 +65,7 @@ [] [] -[Modules/TensorMechanics/Master] +[Physics/SolidMechanics/QuasiStatic] [all] strain = FINITE planar_formulation = plane_strain @@ -141,6 +142,5 @@ [Outputs] csv=true - exodus=true execute_on = TIMESTEP_END [] diff --git a/modules/xfem/test/tests/nucleation_uo/nucleate_bulkCrack.i b/modules/xfem/test/tests/nucleation_uo/nucleate_bulkCrack.i new file mode 100644 index 000000000000..7c073bd475d9 --- /dev/null +++ b/modules/xfem/test/tests/nucleation_uo/nucleate_bulkCrack.i @@ -0,0 +1,160 @@ +[GlobalParams] + displacements = 'disp_x disp_y' +[] + +[XFEM] + geometric_cut_userobjects = 'cut_mesh2' + qrule = volfrac + output_cut_plane = true +[] + +[Mesh] +[gen] + type = GeneratedMeshGenerator + dim = 2 + nx = 10 + ny = 20 + xmin = 0 + xmax = 1.0 + ymin = 0.0 + ymax = 2.0 + elem_type = QUAD4 +[] + +[] + +[DomainIntegral] + integrals = 'InteractionIntegralKI InteractionIntegralKII' + displacements = 'disp_x disp_y' + crack_front_points_provider = cut_mesh2 + 2d=true + number_points_from_provider = 0 + crack_direction_method = CurvedCrackFront + radius_inner = '0.15' + radius_outer = '0.45' + poissons_ratio = 0.0 + youngs_modulus = 207000 + block = 0 + incremental = true + used_by_xfem_to_grow_crack = true +[] + +[AuxVariables] +[strength] + order = CONSTANT + family = MONOMIAL +[] +[] + +[ICs] +[strength] + type = VolumeWeightedWeibull + variable = strength + reference_volume = 1e-2 + weibull_modulus = 1 + median = 5000 +[] +[] + +[UserObjects] + [nucleate] + type = MeshCut2DRankTwoTensorNucleation + tensor = stress + scalar_type = MaxPrincipal + nucleation_threshold = strength + nucleation_radius = .21 + edge_extension_factor = .1 + [] + [cut_mesh2] + type = MeshCut2DFractureUserObject + mesh_file = make_edge_crack_in.e + k_critical=500000 #Large so that cracks will not grow + growth_increment = 0.11 + nucleate_uo = nucleate + [] +[] + +[Physics/SolidMechanics/QuasiStatic] + [all] + strain = FINITE + planar_formulation = plane_strain + add_variables = true + generate_output = 'stress_xx stress_yy vonmises_stress max_principal_stress' + [] +[] + +[Functions] + [bc_pull_top] + type = ParsedFunction + expression = 0.0005*t + [] +[] + +[BCs] + [top_edges] + type = FunctionDirichletBC + boundary = 'top' + variable = disp_y + function = bc_pull_top + [] + [bottom_x] + type = DirichletBC + boundary = bottom + variable = disp_x + value = 0.0 + [] + [bottom_y] + type = DirichletBC + boundary = bottom + variable = disp_y + value = 0.0 + [] +[] + +[Materials] + [elasticity_tensor] + type = ComputeIsotropicElasticityTensor + youngs_modulus = 207000 + poissons_ratio = 0.0 + [] + [stress] + type = ComputeFiniteStrainElasticStress + [] +[] + +[Executioner] + type = Transient + + solve_type = 'PJFNK' + petsc_options_iname = '-ksp_gmres_restart -pc_type -pc_hypre_type -pc_hypre_boomeramg_max_iter' + petsc_options_value = '201 hypre boomeramg 8' + + line_search = 'none' + + [./Predictor] + type = SimplePredictor + scale = 1.0 + [../] + + l_max_its = 100 + l_tol = 1e-2 + + nl_max_its = 15 + nl_rel_tol = 1e-8 + nl_abs_tol = 1e-9 + + start_time = 0.0 + dt = 1.0 + end_time = 5 + max_xfem_update = 1 +[] + +[Outputs] + csv=true + execute_on = final + # exodus=true + # [xfemcutter] + # type=XFEMCutMeshOutput + # xfem_cutter_uo=cut_mesh2 + # [] +[] diff --git a/modules/xfem/test/tests/nucleation_uo/nucleate_edge_crack_2d.i b/modules/xfem/test/tests/nucleation_uo/nucleate_edge_bulk_crack_2d.i similarity index 57% rename from modules/xfem/test/tests/nucleation_uo/nucleate_edge_crack_2d.i rename to modules/xfem/test/tests/nucleation_uo/nucleate_edge_bulk_crack_2d.i index 2907fa76c8ea..b176ae58ddda 100644 --- a/modules/xfem/test/tests/nucleation_uo/nucleate_edge_crack_2d.i +++ b/modules/xfem/test/tests/nucleation_uo/nucleate_edge_bulk_crack_2d.i @@ -13,10 +13,10 @@ [gen] type = GeneratedMeshGenerator dim = 2 - nx = 20 + nx = 60 ny = 10 - xmin = 0 - xmax = 2 + xmin = -3 + xmax = 3 ymin = 0.0 ymax = 1.0 elem_type = QUAD4 @@ -24,29 +24,61 @@ [top_left] type = BoundingBoxNodeSetGenerator new_boundary = pull_top_left - bottom_left = '-0.01 0.99 0' - top_right = '0.11 1.01 0' + bottom_left = '-3.01 0.99 0' + top_right = '-2.99 1.01 0' input = gen [] +[top_mid_left] + type = BoundingBoxNodeSetGenerator + new_boundary = pull_mid_left + bottom_left = '-1.01 0.99 0' + top_right = '-0.99 1.01 0' + input = top_left +[] +[top_mid_right] + type = BoundingBoxNodeSetGenerator + new_boundary = pull_mid_right + bottom_left = '0.99 0.99 0' + top_right = '1.01 1.01 0' + input = top_mid_left +[] [top_right] type = BoundingBoxNodeSetGenerator new_boundary = pull_top_right - bottom_left = '1.89 0.99 0' - top_right = '2.01 1.01 0' - input = top_left + bottom_left = '2.99 0.99 0' + top_right = '3.01 1.01 0' + input = top_mid_right [] -[top_middle_ss] + +[top_mid_left_ss] type = SideSetsFromBoundingBoxGenerator input = top_right - bottom_left = '0.79 0.89 0' - top_right = '1.21 1.01 0' - block_id = '0' - boundary_new = top_middle_ss + bottom_left = '-2.21 0.89 0' + top_right = '-1.79 1.01 0' + boundary_new = top_mid_left_ss + boundaries_old = top +[] +[top_mid_ss] + type = SideSetsFromBoundingBoxGenerator + input = top_mid_left_ss + bottom_left = '-0.21 0.89 0' + top_right = '0.21 1.01 0' + boundary_new = top_mid_ss + boundaries_old = top +[] +[top_mid_right_ss] + type = SideSetsFromBoundingBoxGenerator + input = top_mid_ss + bottom_left = '1.79 0.89 0' + top_right = '2.21 1.01 0' + boundary_new = top_mid_right_ss boundaries_old = top [] -[nucleate] + +[nucleation_strip] + # strip in middle of domain where cracks can nucleate type = ParsedSubdomainMeshGenerator - input = top_middle_ss + input = top_mid_right_ss combinatorial_geometry = 'y > 0.39 & y < 0.51' block_id = 10 [] @@ -64,7 +96,7 @@ poissons_ratio = 0.3 youngs_modulus = 207000 block = 0 - incremental = true + incremental = false used_by_xfem_to_grow_crack = true [] @@ -74,8 +106,9 @@ tensor = stress scalar_type = MaxPrincipal nucleation_threshold = nucleation_threshold - initiate_on_boundary = 'left right' - nucleation_length = .2 + # initiate_on_boundary = 'left right' + nucleation_radius = .41 + nucleation_length = 0.21 [] [cut_mesh2] type = MeshCut2DFractureUserObject @@ -109,38 +142,55 @@ [Functions] [nucleation_x] type = ParsedFunction - expression = '300+x*50' + expression = 'A*cos(pi*x)+D' + symbol_names = 'A D' + symbol_values = '100 200' [] [] -[Modules/TensorMechanics/Master] - [./all] - strain = FINITE - planar_formulation = plane_strain +[Physics/SolidMechanics/QuasiStatic] + [all] + strain = SMALL + planar_formulation = PLANE_STRAIN add_variables = true - generate_output = 'stress_xx stress_yy vonmises_stress max_principal_stress' - [../] + [] [] [Functions] - [bc_pull_top] + [bc_pull_edge] + type = ParsedFunction + expression = 0.0004*t + [] + [bc_pull_mid] type = ParsedFunction expression = 0.0005*t [] [] [BCs] - [top_edges] + [top_edge_nodes] type = FunctionDirichletBC boundary = 'pull_top_left pull_top_right' variable = disp_y - function = bc_pull_top + function = bc_pull_edge [] + [top_mid_nodes] + type = FunctionDirichletBC + boundary = 'pull_mid_left pull_mid_right' + variable = disp_y + function = bc_pull_mid + [] + # [top_middle] + # type = NeumannBC + # boundary = 'top_mid_left_ss top_mid_ss top_mid_right_ss' + # variable = disp_y + # value = -2000 + # [] [top_middle] - type = NeumannBC - boundary = top_middle_ss + type = DirichletBC + boundary = 'top_mid_left_ss top_mid_ss top_mid_right_ss' variable = disp_y - value = -2000 + value = 0 [] [bottom_x] type = DirichletBC @@ -163,16 +213,16 @@ poissons_ratio = 0.3 [../] [./stress] - type = ComputeFiniteStrainElasticStress + type = ComputeLinearElasticStress [../] [] [Executioner] type = Transient - solve_type = 'PJFNK' - petsc_options_iname = '-ksp_gmres_restart -pc_type -pc_hypre_type -pc_hypre_boomeramg_max_iter' - petsc_options_value = '201 hypre boomeramg 8' + solve_type = 'NEWTON' + petsc_options_iname = '-ksp_type -pc_type -pc_factor_mat_solver_package' + petsc_options_value = 'gmres lu superlu_dist' line_search = 'none' @@ -181,6 +231,9 @@ scale = 1.0 [../] + reuse_preconditioner=true + reuse_preconditioner_max_linear_its = 25 + # controls for linear iterations l_max_its = 100 l_tol = 1e-2 @@ -193,18 +246,18 @@ # time control start_time = 0.0 dt = 1.0 - end_time = 5 + end_time = 10 max_xfem_update = 100 [] [Outputs] csv=true - execute_on = TIMESTEP_END + execute_on = FINAL + # exodus=true # [xfemcutter] # type=XFEMCutMeshOutput # xfem_cutter_uo=cut_mesh2 # [] - # console = false [./console] type = Console output_linear = false diff --git a/modules/xfem/test/tests/nucleation_uo/tests b/modules/xfem/test/tests/nucleation_uo/tests index 53c34b45cba1..4e041a003b7d 100644 --- a/modules/xfem/test/tests/nucleation_uo/tests +++ b/modules/xfem/test/tests/nucleation_uo/tests @@ -1,21 +1,52 @@ [Tests] - issues = '#24714' design = 'MeshCut2DRankTwoTensorNucleation.md' - [nucleateCrack] + [nucleateGrowCrack] + issues = '#24714 #27054' type = CSVDiff - input = nucleate_edge_crack_2d.i - csvdiff = 'nucleate_edge_crack_2d_out_II_KI_1_0004.csv' + input = nucleate_edge_bulk_crack_2d.i + csvdiff = 'nucleate_edge_bulk_crack_2d_out_II_KI_1_0011.csv' + rel_err = 1e-3 # XFEM requires --enable-unique-ids in libmesh unique_id = true - requirement = 'The XFEM module shall nucleate and grow cracks based on a nucleation threshold given by an auxvariable.' + requirement = 'The system shall nucleate and grow edge and bulk cracks based on a nucleation threshold given by an auxvariable.' [] [filterNucleation] + issues = '#24714' type = CSVDiff input = nucleate_AllEdgeCracks.i csvdiff = 'nucleate_AllEdgeCracks_out_II_KI_1_0005.csv nucleate_AllEdgeCracks_out_II_KI_1_0003.csv' - rel_err = 0.1 + rel_err = 1e-1 # XFEM requires --enable-unique-ids in libmesh unique_id = true - requirement = 'The XFEM module shall not nucleate cracks within a nucleation_radius of other cracks.' + requirement = 'The system shall not nucleate edge cracks within a nucleation_radius of other cracks.' + [] + [nucleateEdgeBulkCrack] + issues = '#27054' + type = CSVDiff + input = nucleate_bulkCrack.i + csvdiff = 'nucleate_bulkCrack_out_II_KI_1_0006.csv' + # XFEM requires --enable-unique-ids in libmesh + unique_id = true + requirement = 'The system shall nucleate bulk and edge cracks that cut a single element if no nucleation length is given.' + [] + [errorSmallCut] + issues = '#27054' + type = RunException + input = nucleate_bulkCrack.i + cli_args = 'UserObjects/nucleate/nucleation_length=.001' + expect_err = 'Trying to nucleate crack smaller than element length' + # XFEM requires --enable-unique-ids in libmesh + unique_id = true + requirement = 'The system shall error if the nucleation length is smaller than the element being cut.' + [] + [suppressErrorSmallCut] + issues = '#27054' + type = CSVDiff + input = nucleate_bulkCrack.i + cli_args = 'UserObjects/nucleate/nucleation_length=.001 UserObjects/nucleate/nucleate_across_full_element=True' + csvdiff = 'nucleate_bulkCrack_out_II_KI_1_0006.csv' + # XFEM requires --enable-unique-ids in libmesh + unique_id = true + requirement = 'The system shall suppress the error if the nucleation length is smaller than the element being cut and nucleate a cut across the single element.' [] []