Line data    Source code

       1             : //* This file is part of the MOOSE framework
       2             : //* https://mooseframework.inl.gov
       3             : //*
       4             : //* All rights reserved, see COPYRIGHT for full restrictions
       5             : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
       6             : //*
       7             : //* Licensed under LGPL 2.1, please see LICENSE for details
       8             : //* https://www.gnu.org/licenses/lgpl-2.1.html
       9             : 
      10             : #include "SCMTriDuctMeshGenerator.h"
      11             : #include "TriSubChannelMesh.h"
      12             : #include <cmath>
      13             : #include "libmesh/face_quad4.h"
      14             : #include "libmesh/unstructured_mesh.h"
      15             : 
      16             : registerMooseObject("SubChannelApp", SCMTriDuctMeshGenerator);
      17             : registerMooseObjectRenamed("SubChannelApp",
      18             :                            TriDuctMeshGenerator,
      19             :                            "06/30/2025 24:00",
      20             :                            SCMTriDuctMeshGenerator);
      21             : 
      22             : InputParameters
      23         120 : SCMTriDuctMeshGenerator::validParams()
      24             : {
      25         120 :   InputParameters params = MeshGenerator::validParams();
      26         120 :   params.addClassDescription(
      27             :       "Creates a mesh of 1D duct cells around a triangular lattice subassembly");
      28         240 :   params.addRequiredParam<MeshGeneratorName>("input", "The corresponding subchannel mesh");
      29         240 :   params.addParam<unsigned int>("block_id", 2, "Domain Index");
      30         240 :   params.addRequiredParam<unsigned int>("n_cells", "The number of cells in the axial direction");
      31         240 :   params.addParam<Real>("unheated_length_entry", 0.0, "Unheated length at entry [m]");
      32         240 :   params.addRequiredParam<Real>("heated_length", "Heated length [m]");
      33         240 :   params.addParam<Real>("unheated_length_exit", 0.0, "Unheated length at exit [m]");
      34         240 :   params.addRequiredParam<Real>("pitch", "Pitch [m]");
      35         240 :   params.addRequiredParam<unsigned int>("nrings", "Number of fuel Pin rings per assembly [-]");
      36         240 :   params.addRequiredParam<Real>("flat_to_flat",
      37             :                                 "Flat to flat distance for the hexagonal assembly [m]");
      38         120 :   return params;
      39           0 : }
      40             : 
      41          60 : SCMTriDuctMeshGenerator::SCMTriDuctMeshGenerator(const InputParameters & params)
      42             :   : MeshGenerator(params),
      43          60 :     _input(getMesh("input")),
      44         120 :     _n_cells(getParam<unsigned int>("n_cells")),
      45         120 :     _unheated_length_entry(getParam<Real>("unheated_length_entry")),
      46         120 :     _heated_length(getParam<Real>("heated_length")),
      47         120 :     _unheated_length_exit(getParam<Real>("unheated_length_exit")),
      48         120 :     _block_id(getParam<unsigned int>("block_id")),
      49         120 :     _pitch(getParam<Real>("pitch")),
      50         120 :     _n_rings(getParam<unsigned int>("nrings")),
      51         180 :     _flat_to_flat(getParam<Real>("flat_to_flat"))
      52             : {
      53          60 :   SubChannelMesh::generateZGrid(
      54          60 :       _unheated_length_entry, _heated_length, _unheated_length_exit, _n_cells, _z_grid);
      55          60 : }
      56             : 
      57             : std::unique_ptr<MeshBase>
      58          60 : SCMTriDuctMeshGenerator::generate()
      59             : {
      60          60 :   std::unique_ptr<MeshBase> mesh_base = std::move(_input);
      61          60 :   if (!mesh_base)
      62           0 :     mooseError("SCMTriDuctMeshGenerator has to be connected to a sub channel mesh generator.");
      63             : 
      64          60 :   mesh_base->set_mesh_dimension(3);
      65             : 
      66             :   std::vector<Point> corners;
      67          60 :   ductCorners(corners, _flat_to_flat, Point(0, 0, 0));
      68             :   std::vector<Point> xsec;
      69          60 :   ductXsec(xsec, corners, _n_rings, _pitch, _flat_to_flat);
      70             :   std::vector<Point> points;
      71          60 :   ductPoints(points, xsec, _z_grid);
      72             :   std::vector<std::vector<size_t>> elem_point_indices;
      73          60 :   ductElems(elem_point_indices, _z_grid.size(), xsec.size());
      74             :   std::vector<Node *> duct_nodes;
      75          60 :   buildDuct(mesh_base, duct_nodes, points, elem_point_indices, _block_id);
      76          60 :   mesh_base->subdomain_name(_block_id) = name();
      77             : 
      78          60 :   mesh_base->prepare_for_use();
      79             : 
      80          60 :   auto & sch_mesh = static_cast<TriSubChannelMesh &>(*_mesh);
      81          60 :   sch_mesh.setChannelToDuctMaps(duct_nodes);
      82          60 :   sch_mesh._duct_mesh_exist = true;
      83             : 
      84          60 :   return mesh_base;
      85          60 : }
      86             : 
      87             : void
      88          60 : SCMTriDuctMeshGenerator::ductCorners(std::vector<Point> & corners,
      89             :                                      Real flat_to_flat,
      90             :                                      const Point & center) const
      91             : {
      92          60 :   corners.resize(TriSubChannelMesh::N_CORNERS);
      93          60 :   Real r_corner = flat_to_flat / 2. / std::cos(libMesh::pi / 6.);
      94         420 :   for (size_t i = 0; i < TriSubChannelMesh::N_CORNERS; i++)
      95             :   {
      96         360 :     Real theta = i * libMesh::pi / 3.;
      97         360 :     Point corner = {r_corner * std::cos(theta), r_corner * std::sin(theta)};
      98         360 :     corners[i] = center + corner;
      99             :   }
     100          60 : }
     101             : 
     102             : void
     103          60 : SCMTriDuctMeshGenerator::ductXsec(std::vector<Point> & xsec,
     104             :                                   const std::vector<Point> & corners,
     105             :                                   unsigned int nrings,
     106             :                                   Real pitch,
     107             :                                   Real flat_to_flat) const
     108             : {
     109          60 :   xsec.clear();
     110             : 
     111          60 :   Real r_corner = flat_to_flat / 2. / std::cos(libMesh::pi / 6.);
     112          60 :   Real start_offset = (r_corner - (nrings - 2) * pitch) * std::sin(libMesh::pi / 6.);
     113          60 :   Real side_length = (corners[0] - corners[1]).norm();
     114             : 
     115         420 :   for (size_t i = 0; i < corners.size(); i++)
     116             :   {
     117         360 :     auto left = corners[i];
     118         360 :     auto right = corners[(i + 1) % corners.size()];
     119         360 :     xsec.push_back(left);
     120         360 :     auto direc = (right - left).unit();
     121        1674 :     for (Real offset_from_corner = start_offset; offset_from_corner < side_length;
     122        1314 :          offset_from_corner += pitch)
     123        1314 :       xsec.push_back(left + direc * offset_from_corner);
     124             :   }
     125          60 : }
     126             : 
     127             : size_t
     128      120474 : SCMTriDuctMeshGenerator::ductPointIndex(unsigned int points_per_layer,
     129             :                                         unsigned int layer,
     130             :                                         unsigned int point) const
     131             : {
     132      120474 :   return layer * points_per_layer + point;
     133             : }
     134             : 
     135             : void
     136          60 : SCMTriDuctMeshGenerator::ductPoints(std::vector<Point> & points,
     137             :                                     const std::vector<Point> & xsec,
     138             :                                     const std::vector<Real> & z_layers) const
     139             : {
     140          60 :   points.resize(xsec.size() * z_layers.size());
     141         750 :   for (size_t i = 0; i < z_layers.size(); i++)
     142       26124 :     for (size_t j = 0; j < xsec.size(); j++)
     143       25434 :       points[ductPointIndex(xsec.size(), i, j)] = Point(xsec[j](0), xsec[j](1), z_layers[i]);
     144          60 : }
     145             : 
     146             : void
     147          60 : SCMTriDuctMeshGenerator::ductElems(std::vector<std::vector<size_t>> & elem_point_indices,
     148             :                                    unsigned int n_layers,
     149             :                                    unsigned int points_per_layer) const
     150             : {
     151          60 :   elem_point_indices.clear();
     152         690 :   for (unsigned int i = 0; i < n_layers - 1; i++)
     153             :   {
     154             :     unsigned int bottom = i;
     155         630 :     unsigned int top = i + 1;
     156       24390 :     for (unsigned int j = 0; j < points_per_layer; j++)
     157             :     {
     158             :       unsigned int left = j;
     159       23760 :       unsigned int right = (j + 1) % points_per_layer;
     160       71280 :       elem_point_indices.push_back({ductPointIndex(points_per_layer, bottom, left),
     161       23760 :                                     ductPointIndex(points_per_layer, bottom, right),
     162       23760 :                                     ductPointIndex(points_per_layer, top, right),
     163       23760 :                                     ductPointIndex(points_per_layer, top, left)});
     164             :     }
     165             :   }
     166          60 : }
     167             : 
     168             : void
     169          60 : SCMTriDuctMeshGenerator::buildDuct(std::unique_ptr<MeshBase> & mesh,
     170             :                                    std::vector<Node *> & duct_nodes,
     171             :                                    const std::vector<Point> & points,
     172             :                                    const std::vector<std::vector<size_t>> & elem_point_indices,
     173             :                                    SubdomainID block) const
     174             : {
     175       25494 :   for (size_t i = 0; i < points.size(); i++)
     176       25434 :     duct_nodes.push_back(mesh->add_point(points[i]));
     177             : 
     178       23820 :   for (auto & elem_indices : elem_point_indices)
     179             :   {
     180       47520 :     auto elem = mesh->add_elem(new Quad4());
     181       23760 :     elem->subdomain_id() = block;
     182      118800 :     for (size_t i = 0; i < elem_indices.size(); i++)
     183       95040 :       elem->set_node(i, duct_nodes[elem_indices[i]]);
     184             :   }
     185          60 : }