Potential off-by-one error in HexagonalGridPositions.C

[lewisgross1296]

Bug Description

Recently, #26935 added the capability to generate various positions within a Hexagonal Lattice (which I am very grateful for). I went to try it out today and encountered

*** ERROR ***
/home/lgross/gcmr/mwes/small_inf_assembly/coolant_channel_generator.i:21: (Positions/hex_grid/lattice_flat_to_flat):
    Bundle pitch is too small to fit this many rings at that pin pitch

but I think either I am misunderstanding nr or there is an off-by-one error.

When counting rings, does the central ring count or not? For example, I would say this hexagaonal lattice has 6 rings

I would say that the flat-to-flat here is $5pitch\sqrt{3}$. Currently, https://github.com/idaholab/moose/blob/659e12791ceaa627c75956d97d009b3f9f8963c0/modules/reactor/src/positions/HexagonalGridPositions.C#L61 would throw an error because $6pitch\sqrt{3} > 5pitch\sqrt{3}$. I believe the correct condition to throw an error should be $(nr-1)p\sqrt{3} > lattice$ $flat$ $to$ $flat$

I've attached some diagrams for a hexagon of arbitrary number of rings nr. Here a hexagon is shown where the edge length is $(nr-1)p$ Here is the 30-60-90 triangle created that relates the side length to the flat to flat distance. Note the base of this triangle is half the flat to flat distance.

Here the math tells us $$\big( \frac{1}{2} dftf \big)^2 + \big( \frac{1}{2} (nr-1)p \big)^2 = \big((nr-1)p \big)^2$$ $$\big(\frac{1}{2} dftf \big)^2 = \big((nr-1)p \big)^2 - \big( \frac{1}{2} (nr-1)p \big)^2$$ $$\frac{1}{2} dftf = (nr-1)p \frac{\sqrt{3}}{2}$$ Thus, the minimum flat to flat distance is given by $$dftf = (nr-1)p \sqrt{3}$$

Steps to Reproduce

Run this file coolant_channel_generator.i

pin_lattice_pitch = 0.02
[GlobalParams]
  quad_center_elements = true
[]

[Mesh]
  [cmg]
    type = CartesianMeshGenerator
    dx = 1
    dim = 1
  []
[]

[Positions]
  active = 'hex_grid'
  [hex_grid]
    type = HexagonalGridPositions
    center = '0.0 0.0 0.0'
    nr = 3
    pin_pitch = ${pin_lattice_pitch}
    lattice_flat_to_flat = ${fparse pin_lattice_pitch*sqrt(3)*2}
    pattern =  '0 1 0;
               1 0 0 1;
              0 0 1 0 2;
               1 0 0 1;
                0 1 0'
    include_in_pattern = '1'
    outputs = 'out'
  []
[]

# The following content is adding postprocessor(s) to check sideset areas.
# The reactor module is unfortunately quite brittle in its assignment of sideset
# IDs, so we want to be extra sure that any changes to sideset numbering are detected
# in our test suite.
[Problem]
  solve = false
[]

[Executioner]
  type = Steady
[]

[Outputs]
  [out]
    type = JSON
    execute_on = 'FINAL'
    file_base = 'coolant_channel_positions'
  []
[]

Impact

This would be a very small PR and I could likely get around it by using a different lattice_flat_to_flat, but I think it should be resolved.

[Optional] Diagnostics

No response

[lewisgross1296]

Tagging @GiudGiud

[lewisgross1296]

I just tried the MWE with $nr=2$ and I think that confirms my understanding of ring counting

*** ERROR ***
The following error occurred in the object "hex_grid", of type "HexagonalGridPositions".

Number of rows in pattern 5 should be equal to twice the number of hexagonal rings minus one

AKA the above code snippet should be considered 3 rings and the error checking condition is off-by-one.

[GiudGiud]

First picture is 5.5 rings Yes it possible this check is off by 1

[lewisgross1296]

Ah gotcha, 5.5 also makes sense. I only call it 6 rings because I made the image using a HexLattice in OpenMC and there it is defined as 6. The hexagonal prism I filled it into clips the last ring to "a half." I think this kind of reflective boundary (only really done in simulations only) is the limit of what should be allowed for a minimum flat to flat.

[GiudGiud]

Happy to review a patch if you feel comfortable making it

[lewisgross1296]

Sure! This seems very doable, so I'll give it a try.

[lewisgross1296]

Note: all this math above is for a half-ring situation. I should redo this when the number of rings is only an integer and there isn't a smooth boundary