Merge pull request #25179 from sourabhkadambi/columnar_3D_hex

Adding 3D columnar option in PolycrystalHex userobject

modules/phase_field/doc/content/source/userobjects/PolycrystalHex.md

@@ -1,12 +1,11 @@
 # PolycrystalHex
 
-The PolycrystalHex object generates a set of points in 2D or 3D space such that hexagonal patterns are generated from a resulting Voronoi Tesselation. The number of grains must be choosen such that the $\sqrt[dim]{n}$ root is an integer. Furthermore, a suitable number must be chosen that is supported by the mesh resolution and grain boundary width. Once the grain structure has been generated, a coloring algorithm is used to assign order parameters to grains so that a reduced number of order parameters may be used. Typical numbers are given here:
+The PolycrystalHex object generates a set of points in 2D space such that hexagonal patterns are generated from a resulting Voronoi Tesselation. The number of grains must be chosen such that the $\sqrt[dim=2]{n}$ root is an even integer. Furthermore, a suitable number must be chosen that is supported by the mesh resolution and grain boundary width. At present, the object is supported for generating hexagons in 2D and hexagonal prisms in 3D with the `columnar_3D` option. It works with periodic boundary conditions. Once the grain structure has been generated, a coloring algorithm is used to assign order parameters to grains so that a reduced number of order parameters may be used. Typical numbers are given here:
 
 !table style=border:4px solid black;width:350px;
 | Dimension | Recommended OPs |
 |-----------|-----------------|
 | 2D        | 8               |
-| 3D        | 25              |
 
 See [Polycrystal Initial Conditions](ICs/PolycrystalICs.md) for more information.
 

modules/phase_field/src/userobjects/PolycrystalHex.C

@@ -33,21 +33,41 @@ PolycrystalHex::PolycrystalHex(const InputParameters & parameters)
     _x_offset(getParam<Real>("x_offset")),
     _perturbation_percent(getParam<Real>("perturbation_percent"))
 {
+  if (_columnar_3D == false && _dim == 3)
+    mooseError(
+        "PolycrystalHex is supported on 2D domains or 3D domains with the columnar_3D option");
+  if (_grain_num % 2 != 0)
+    mooseError("PolycrystalHex requires an even square number for 2D or columnar 3D");
+
   _random.seed(_tid, getParam<unsigned int>("rand_seed"));
 }
 
 void
 PolycrystalHex::precomputeGrainStructure()
 {
-  const unsigned int root = MathUtils::round(std::pow(_grain_num, 1.0 / _dim));
+  unsigned int d = _dim;
+
+  if (_columnar_3D && _dim == 3)
+    d -= 1;
+
+  // check if periodic boundary condition is set
+  for (unsigned int j = 0; j < d; ++j)
+  {
+    for (unsigned int i = 0; i < _vars.size(); ++i)
+      if (!_mesh.isTranslatedPeriodic(_vars[i]->number(), j))
+        mooseError("PolycrystalHex works only with periodic BCs");
+  }
+
+  const unsigned int root = MathUtils::round(std::pow(_grain_num, 1.0 / d));
 
   // integer power the rounded root and check if we recover the grain number
   unsigned int grain_pow = root;
-  for (unsigned int i = 1; i < _dim; ++i)
+  for (unsigned int i = 1; i < d; ++i)
     grain_pow *= root;
 
   if (_grain_num != grain_pow)
-    mooseError("PolycrystalHex requires a square or cubic number depending on the mesh dimension");
+    mooseError(
+        "PolycrystalHex requires a square number for 2D or columnar 3D and a cubic number for 3D");
 
   // Set up domain bounds with mesh tools
   for (const auto i : make_range(Moose::dim))
@@ -67,8 +87,9 @@ PolycrystalHex::precomputeGrainStructure()
   // Assign the relative center points positions, defining the grains according to a hexagonal
   // pattern
   unsigned int count = 0;
-  for (unsigned int k = 0; k < (_dim == 3 ? root : 1); ++k)
-    for (unsigned int j = 0; j < (_dim >= 2 ? root : 1); ++j)
+
+  for (unsigned int k = 0; k < (d == 3 ? root : 1); ++k)
+    for (unsigned int j = 0; j < (d >= 2 ? root : 1); ++j)
       for (unsigned int i = 0; i < root; ++i)
       {
         // set x-coordinate
@@ -99,6 +120,9 @@ PolycrystalHex::precomputeGrainStructure()
         _centerpoints[grain](i) = _top_right(i);
       if (_centerpoints[grain](i) < _bottom_left(i))
         _centerpoints[grain](i) = _bottom_left(i);
+
+      if (_columnar_3D && _dim == 3)
+        _centerpoints[grain](2) = _bottom_left(2) + _range(2) * 0.5;
     }
 
   buildSearchTree();

modules/phase_field/test/tests/initial_conditions/BlockRestriction.i

@@ -25,10 +25,6 @@
   []
   [tricrystal_2]
   []
-  [PolycrystalVariables]
-    var_name_base = 'coloring_'
-    op_num = 8
-  []
   [random_0]
   []
   [random_1]
@@ -75,12 +71,6 @@
       var_name_base = 'tricrystal_'
       op_num = 3
     []
-    [PolycrystalColoringIC]
-      block = '0'
-      polycrystal_ic_uo = hex_ic
-      var_name_base = 'coloring_'
-      op_num = 8
-    []
     [PolycrystalRandomIC]
       block = '0'
       var_name_base = 'random_'
@@ -110,13 +100,6 @@
     op_num = 8
     grain_num = 4
   []
-  [hex_ic]
-    type = PolycrystalHex
-    coloring_algorithm = bt
-    var_name_base = 'coloring_'
-    op_num = 8
-    grain_num = 4
-  []
 []
 
 [Problem]

modules/phase_field/test/tests/initial_conditions/HexPolycrystalIC.i

@@ -34,7 +34,7 @@
 [BCs]
   [./Periodic]
     [./all]
-      auto_direction = 'x'
+      auto_direction = 'x y'
     [../]
   [../]
 []

modules/phase_field/test/tests/initial_conditions/HexPolycrystalIC_3D_columnar.i

@@ -0,0 +1,71 @@
+[Mesh]
+  type = GeneratedMesh
+  dim = 3
+  nx = 20
+  ny = 20
+  nz = 1
+  xmax = 1
+  ymax = 1
+  zmax = 0.1
+[]
+
+[GlobalParams]
+  op_num = 4
+  grain_num = 4
+  var_name_base = gr
+  int_width = 0.05
+[]
+
+[Variables]
+  [./PolycrystalVariables]
+  [../]
+[]
+
+[ICs]
+  [./PolycrystalICs]
+    [./PolycrystalColoringIC]
+      polycrystal_ic_uo = hex_ic
+    [../]
+  [../]
+  [./bnds]
+    type = BndsCalcIC
+    variable = bnds
+  [../]
+[]
+
+[AuxVariables]
+  [./bnds]
+    order = FIRST
+    family = LAGRANGE
+  [../]
+[]
+
+[UserObjects]
+  [./hex_ic]
+    type = PolycrystalHex
+    coloring_algorithm = bt
+    columnar_3D = true
+  [../]
+[]
+
+[BCs]
+  [./Periodic]
+    [./all]
+      auto_direction = 'x y'
+    [../]
+  [../]
+[]
+
+[Problem]
+  type = FEProblem
+  solve = false
+[]
+
+[Executioner]
+  type = Transient
+  num_steps = 0
+[]
+
+[Outputs]
+  exodus = true
+[]

modules/phase_field/test/tests/initial_conditions/tests

@@ -179,7 +179,17 @@
       # Testing optimal hex coloring with a backtracking algorithm
       method = '!DBG'
 
-      detail = 'hexagonal structure, and'
+      detail = 'hexagonal structure in 2D,'
+    []
+
+    [HexPolycrystalIC_3D]
+      type = Exodiff
+      input = HexPolycrystalIC_3D_columnar.i
+      exodiff = HexPolycrystalIC_3D_columnar_out.e
+      # Testing optimal hex coloring with a backtracking algorithm
+      method = '!DBG'
+
+      detail = 'hexagonal structure in 3D, and'
     []
 
     [TricrystalTripleJunctionIC]