Compute radius correction factor for error checking and add azimuthal…

… node tolerance parameter

modules/reactor/include/meshgenerators/PolygonMeshGeneratorBase.h

@@ -109,6 +109,21 @@ class PolygonMeshGeneratorBase : public MeshGenerator
     Real bias;
   };
 
+  /**
+   * Makes radial correction to preserve ring area.
+   * @param azimuthal_list azimuthal angles (in degrees) of all the nodes on the circle
+   * @param full_circle whether the circle is a full or partial circle
+   * @param order order of mesh elements
+   * @param is_first_value_vertex whether the first value of the azimuthal_list belongs to a vertex
+   * instead of a midpoint
+   * @return a correction factor to preserve the area of the circle after polygonization during
+   * meshing
+   */
+  static Real radiusCorrectionFactor(const std::vector<Real> & azimuthal_list,
+                                     const bool full_circle = true,
+                                     const unsigned int order = 1,
+                                     const bool is_first_value_vertex = true);
+
 protected:
   /**
    * Creates a mesh of a slice that corresponds to a single side of the polygon to be generated.
@@ -503,21 +518,6 @@ class PolygonMeshGeneratorBase : public MeshGenerator
                    const bool generate_side_specific_boundaries = true,
                    const QUAD_ELEM_TYPE quad_elem_type = QUAD_ELEM_TYPE::QUAD4) const;
 
-  /**
-   * Makes radial correction to preserve ring area.
-   * @param azimuthal_list azimuthal angles (in degrees) of all the nodes on the circle
-   * @param full_circle whether the circle is a full or partial circle
-   * @param order order of mesh elements
-   * @param is_first_value_vertex whether the first value of the azimuthal_list belongs to a vertex
-   * instead of a midpoint
-   * @return a correction factor to preserve the area of the circle after polygonization during
-   * meshing
-   */
-  Real radiusCorrectionFactor(const std::vector<Real> & azimuthal_list,
-                              const bool full_circle = true,
-                              const unsigned int order = 1,
-                              const bool is_first_value_vertex = true) const;
-
   /**
    * Based on a pair of azimuthal angles, calculates the volume of a TRI6 element with one vertex at
    * the origin, the other two vertices on the unit circle. Here, the second vertex is on the
@@ -527,7 +527,7 @@ class PolygonMeshGeneratorBase : public MeshGenerator
    * @param azi_pair a pair of the input azimuthal angles
    * @return the volume of the TRI6 element
    */
-  Real dummyTRI6VolCalculator(const std::pair<Real, Real> & azi_pair) const;
+  static Real dummyTRI6VolCalculator(const std::pair<Real, Real> & azi_pair);
 
   /**
    * Creates peripheral area mesh for the patterned hexagon mesh. Note that the function create the

modules/reactor/src/meshgenerators/ControlDrumMeshGenerator.C

@@ -14,6 +14,7 @@
 #include "Factory.h"
 #include "libmesh/elem.h"
 #include "MooseMeshUtils.h"
+#include "PolygonMeshGeneratorBase.h"
 
 registerMooseObject("ReactorApp", ControlDrumMeshGenerator);
 
@@ -64,12 +65,17 @@ ControlDrumMeshGenerator::validParams()
       "block_names",
       "Block names for each radial and axial zone. "
       "Inner indexing is radial zones (drum inner/drum/drum outer), outer indexing is axial");
+  params.addParam<Real>("azimuthal_node_tolerance",
+                        0.1,
+                        "(in degrees) The absolute tolerance for which to shift an azimuthal node "
+                        "to match the pad start/end angles");
 
   params.addParamNamesToGroup("region_ids assembly_type", "ID assigment");
   params.addParamNamesToGroup("drum_inner_radius drum_outer_radius drum_inner_intervals "
                               "drum_intervals num_azimuthal_sectors",
                               "Drum specifications");
-  params.addParamNamesToGroup("pad_start_angle pad_end_angle", "Control pad specifications");
+  params.addParamNamesToGroup("pad_start_angle pad_end_angle azimuthal_node_tolerance",
+                              "Control pad specifications");
 
   params.addClassDescription(
       "This ControlDrumMeshGenerator object is designed to generate "
@@ -108,14 +114,6 @@ ControlDrumMeshGenerator::ControlDrumMeshGenerator(const InputParameters & param
   const auto num_sectors = getParam<unsigned int>("num_azimuthal_sectors");
   const auto assembly_pitch = getReactorParam<Real>(RGMB::assembly_pitch);
 
-  // Check drum radius parameters
-  if (_drum_inner_radius >= _drum_outer_radius)
-    paramError("drum_outer_radius", "Drum outer radius must be larger than the inner radius");
-  if (_drum_outer_radius >= assembly_pitch / 2.)
-    paramError("drum_outer_radius",
-               "Outer radius of drum region must be smaller than half the assembly pitch as "
-               "defined by 'ReactorMeshParams/assembly_pitch'");
-
   // Check drum pad parameters
   if (isParamSetByUser("pad_start_angle"))
   {
@@ -139,6 +137,22 @@ ControlDrumMeshGenerator::ControlDrumMeshGenerator(const InputParameters & param
     _has_pad_region = false;
   }
 
+  // Error checking for azimuthal node tolerance
+  const auto azimuthal_node_tolerance = getParam<Real>("azimuthal_node_tolerance");
+  if (!_has_pad_region && isParamSetByUser("azimuthal_node_tolerance"))
+    paramWarning("azimuthal_node_tolerance",
+                 "This parameter is relevant only when pad start and end angles are defined");
+  if (_has_pad_region && MooseUtils::absoluteFuzzyGreaterEqual(2. * azimuthal_node_tolerance,
+                                                               360. / (Real)num_sectors))
+    paramError("azimuthal_node_tolerance",
+               "Azimuthal node tolerance should be smaller than half the azimuthal interval size "
+               "as defined by 'num_azimuthal_sectors'");
+  if (_has_pad_region && MooseUtils::absoluteFuzzyGreaterEqual(2. * azimuthal_node_tolerance,
+                                                               _pad_end_angle - _pad_start_angle))
+    paramError("azimuthal_node_tolerance",
+               "Azimuthal node tolerance should be smaller than half the difference of the pad "
+               "angle range");
+
   // Check region IDs have correct size
   const unsigned int n_radial_regions = _has_pad_region ? 4 : 3;
   if (isParamValid("region_ids"))
@@ -194,40 +208,78 @@ ControlDrumMeshGenerator::ControlDrumMeshGenerator(const InputParameters & param
     mooseError("Both top_boundary_id and bottom_boundary_id must be provided in ReactorMeshParams "
                "if using extruded geometry");
 
-  // No subgenerators will be called if option to bypass mesh generators is enabled
-  if (!getReactorParam<bool>(RGMB::bypass_meshgen))
+  // Define azimuthal angles explicitly based on num_azimuthal_sectors and manually add pad start
+  // angle and end angle if they are not contained within these angle intervals
+  std::vector<Real> azimuthal_angles;
+  const auto custom_start_angle = _pad_start_angle;
+  const auto custom_end_angle = (_pad_end_angle > 360) ? _pad_end_angle - 360. : _pad_end_angle;
+  for (unsigned int i = 0; i < num_sectors; ++i)
   {
-    const std::string block_name_prefix =
-        RGMB::DRUM_BLOCK_NAME_PREFIX + std::to_string(_assembly_type);
-
-    // Define azimuthal angles explicitly based on num_azimuthal_sectors and manually add pad start
-    // angle and end angle if they are not contained within these angle intervals
-    std::vector<Real> azimuthal_angles;
-    const auto custom_start_angle = _pad_start_angle;
-    const auto custom_end_angle = (_pad_end_angle > 360) ? _pad_end_angle - 360. : _pad_end_angle;
-    for (unsigned int i = 0; i < num_sectors; ++i)
-    {
-      Real current_azim_angle = (Real)i * 360. / (Real)num_sectors;
-      Real next_azim_angle = (Real)(i + 1) * 360. / (Real)num_sectors;
+    Real current_azim_angle = (Real)i * 360. / (Real)num_sectors;
+    Real next_azim_angle = (Real)(i + 1) * 360. / (Real)num_sectors;
+    if (!_has_pad_region)
       azimuthal_angles.push_back(current_azim_angle);
-
-      if (!MooseUtils::absoluteFuzzyEqual(current_azim_angle, custom_start_angle) &&
-          !MooseUtils::absoluteFuzzyEqual(next_azim_angle, custom_start_angle) &&
-          (custom_start_angle > current_azim_angle) && (custom_start_angle < next_azim_angle))
+    else
+    {
+      // When pad regions are involved, check if the current azimuthal node angle coincides with
+      // pad start/end angle to within tolerance. If it does, then shift the azimuthal node location
+      // to match the pad angle location. If it does not and the pad angle falls within the
+      // azimuthal sector, then create an additional node for the pad location
+      if (MooseUtils::absoluteFuzzyLessEqual(std::abs(current_azim_angle - custom_start_angle),
+                                             azimuthal_node_tolerance))
+      {
+        azimuthal_angles.push_back(custom_start_angle);
+      }
+      else if (MooseUtils::absoluteFuzzyLessEqual(std::abs(current_azim_angle - custom_end_angle),
+                                                  azimuthal_node_tolerance))
+      {
+        azimuthal_angles.push_back(custom_end_angle);
+      }
+      else if (MooseUtils::absoluteFuzzyGreaterThan(custom_start_angle - current_azim_angle,
+                                                    azimuthal_node_tolerance) &&
+               MooseUtils::absoluteFuzzyGreaterThan(next_azim_angle - custom_start_angle,
+                                                    azimuthal_node_tolerance))
       {
         mooseWarning("pad_start_angle not contained within drum azimuthal discretization so "
                      "additional azimuthal nodes are created to align mesh with this angle");
+        azimuthal_angles.push_back(current_azim_angle);
         azimuthal_angles.push_back(custom_start_angle);
       }
-      if (!MooseUtils::absoluteFuzzyEqual(current_azim_angle, custom_end_angle) &&
-          !MooseUtils::absoluteFuzzyEqual(next_azim_angle, custom_end_angle) &&
-          (custom_end_angle > current_azim_angle) && (custom_end_angle < next_azim_angle))
+      else if (MooseUtils::absoluteFuzzyGreaterThan(custom_end_angle - current_azim_angle,
+                                                    azimuthal_node_tolerance) &&
+               MooseUtils::absoluteFuzzyGreaterThan(next_azim_angle - custom_end_angle,
+                                                    azimuthal_node_tolerance))
       {
         mooseWarning("pad_end_angle not contained within drum azimuthal discretization so "
                      "additional azimuthal nodes are created to align mesh with this angle");
+        azimuthal_angles.push_back(current_azim_angle);
         azimuthal_angles.push_back(custom_end_angle);
       }
+      else
+        azimuthal_angles.push_back(current_azim_angle);
     }
+  }
+
+  // Check drum radius parameters
+  if (_drum_inner_radius >= _drum_outer_radius)
+    paramError("drum_outer_radius", "Drum outer radius must be larger than the inner radius");
+  // Check if volume preserved radius of outer radius exceeds assembly halfpitch. In data driven
+  // mode, radius is assumed not to need volume preservation
+  auto radius_corrfac = getReactorParam<bool>(RGMB::bypass_meshgen)
+                            ? 1.0
+                            : PolygonMeshGeneratorBase::radiusCorrectionFactor(azimuthal_angles);
+  if (_drum_outer_radius * radius_corrfac >= assembly_pitch / 2.)
+    paramError("drum_outer_radius",
+               "Volume-corrected outer radius of drum region must be smaller than half the "
+               "assembly pitch as "
+               "defined by 'ReactorMeshParams/assembly_pitch'");
+
+  // No subgenerators will be called if option to bypass mesh generators is enabled
+  if (!getReactorParam<bool>(RGMB::bypass_meshgen))
+  {
+    const std::string block_name_prefix =
+        RGMB::DRUM_BLOCK_NAME_PREFIX + std::to_string(_assembly_type);
+
     {
       // Invoke AdvancedConcentricCircleGenerator to define drum mesh without background region
       auto params = _app.getFactory().getValidParams("AdvancedConcentricCircleGenerator");

modules/reactor/src/meshgenerators/PolygonMeshGeneratorBase.C

@@ -1140,7 +1140,7 @@ Real
 PolygonMeshGeneratorBase::radiusCorrectionFactor(const std::vector<Real> & azimuthal_list,
                                                  const bool full_circle,
                                                  const unsigned int order,
-                                                 const bool is_first_value_vertex) const
+                                                 const bool is_first_value_vertex)
 {
   Real tmp_acc = 0.0;
   Real tmp_acc_azi = 0.0;
@@ -1195,7 +1195,7 @@ PolygonMeshGeneratorBase::radiusCorrectionFactor(const std::vector<Real> & azimu
 }
 
 Real
-PolygonMeshGeneratorBase::dummyTRI6VolCalculator(const std::pair<Real, Real> & azi_pair) const
+PolygonMeshGeneratorBase::dummyTRI6VolCalculator(const std::pair<Real, Real> & azi_pair)
 {
   // The algorithm is copied from libMesh's face_tri6.C
   // Original license is LGPL so it can be used here.

modules/reactor/test/tests/meshgenerators/control_drum_mesh_generator/gold/drum_hex_pad_shift_nodes_out.csv

@@ -0,0 +1,2 @@
+time,area_reg1,area_reg2,area_reg3,area_reg4
+1,201.06192982975,14.74052266086,44.352835153164,86.254873870005

modules/reactor/test/tests/meshgenerators/control_drum_mesh_generator/tests

@@ -42,7 +42,7 @@
       mesh_mode = REPLICATED
     []
     [pad_azimuthal_sector_mismatch]
-      detail = 'additional azimuthal nodes for a control drum assembly mesh when the azimuthal discretization does not line up with the drum start and end angles'
+      detail = 'additional azimuthal nodes for a control drum assembly mesh when the azimuthal discretization does not line up with the pad start and end angles'
       type = 'CSVDiff'
       input = 'drum_pad.i'
       csvdiff = 'drum_hex_pad_out.csv'
@@ -52,6 +52,19 @@
       recover = false
       mesh_mode = REPLICATED
     []
+    [pad_azimuthal_sector_shift]
+      detail = 'shift azimuthal nodes for a control drum assembly mesh to match the drum pad locations when the azimuthal angles line up to within tolerance'
+      type = 'CSVDiff'
+      input = 'drum_pad.i'
+      csvdiff = 'drum_hex_pad_shift_nodes_out.csv'
+      cli_args = 'Mesh/drum/pad_start_angle=90.1
+                  Mesh/drum/pad_end_angle=179.9
+                  Mesh/drum/azimuthal_node_tolerance=0.2
+                  Outputs/file_base=drum_hex_pad_shift_nodes_out'
+      allow_warnings = true
+      recover = false
+      mesh_mode = REPLICATED
+    []
     [drum_depletion_id]
       detail = 'depletion IDs for each region in a 2D control drum mesh'
       type = 'CSVDiff'
@@ -81,8 +94,9 @@
       detail = 'if the outer radius of the control drum exceeds the radius of the inscribed circle of the assembly boundary'
       type = 'RunException'
       input = 'drum_nopad.i'
-      cli_args = 'Mesh/drum/drum_outer_radius=10.1'
-      expect_err = "Outer radius of drum region must be smaller"
+      cli_args = 'Mesh/drum/drum_outer_radius=9.9
+                  Mesh/drum/num_azimuthal_sectors=6'
+      expect_err = "Volume-corrected outer radius of drum region must be smaller"
     []
     [missing_pad_end_angle]
       detail = 'if the user sets the start angle of the drum pad but not the end angle'
@@ -114,6 +128,23 @@
                   Mesh/drum/pad_end_angle=360'
       expect_err = "The difference between 'pad_end_angle' and 'pad_start_angle' must be between"
     []
+    [node_tolerance_interval_size]
+      detail = 'if the azimuthal node tolerance is larger than half the azimuthal interval size'
+      type = 'RunException'
+      input = 'drum_pad.i'
+      cli_args = 'Mesh/drum/azimuthal_node_tolerance=1'
+      expect_err = "Azimuthal node tolerance should be smaller than half the azimuthal interval size"
+    []
+    [node_tolerance_pad_angle_range]
+      detail = 'if the azimuthal node tolerance is larger than half the difference of the pad angle range'
+      type = 'RunException'
+      input = 'drum_pad.i'
+      cli_args = 'Mesh/drum/azimuthal_node_tolerance=1
+                  Mesh/drum/num_azimuthal_sectors=36
+                  Mesh/drum/pad_start_angle=90
+                  Mesh/drum/pad_end_angle=92'
+      expect_err = "Azimuthal node tolerance should be smaller than half the difference of the pad angle range"
+    []
     [no_region_ids]
       detail = 'if region IDs are not provided for the control drum mesh'
       type = 'RunException'