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 "QuadSubChannelMesh.h"
      11             : #include <cmath>
      12             : #include <limits>
      13             : #include "libmesh/edge_edge2.h"
      14             : #include "libmesh/unstructured_mesh.h"
      15             : 
      16             : registerMooseObject("SubChannelApp", QuadSubChannelMesh);
      17             : 
      18             : InputParameters
      19         412 : QuadSubChannelMesh::validParams()
      20             : {
      21         412 :   InputParameters params = SubChannelMesh::validParams();
      22         412 :   params.addClassDescription("Creates an subchannel mesh container for a square "
      23             :                              "lattice arrangement");
      24         412 :   return params;
      25           0 : }
      26             : 
      27         206 : QuadSubChannelMesh::QuadSubChannelMesh(const InputParameters & params) : SubChannelMesh(params) {}
      28             : 
      29           0 : QuadSubChannelMesh::QuadSubChannelMesh(const QuadSubChannelMesh & other_mesh)
      30             :   : SubChannelMesh(other_mesh),
      31           0 :     _nx(other_mesh._nx),
      32           0 :     _ny(other_mesh._ny),
      33           0 :     _n_channels(other_mesh._n_channels),
      34           0 :     _n_gaps(other_mesh._n_gaps),
      35           0 :     _n_pins(other_mesh._n_pins),
      36           0 :     _side_gap(other_mesh._side_gap),
      37           0 :     _nodes(other_mesh._nodes),
      38           0 :     _pin_nodes(other_mesh._pin_nodes),
      39           0 :     _gapnodes(other_mesh._gapnodes),
      40           0 :     _gap_to_chan_map(other_mesh._gap_to_chan_map),
      41           0 :     _gap_to_pin_map(other_mesh._gap_to_pin_map),
      42           0 :     _chan_to_gap_map(other_mesh._chan_to_gap_map),
      43           0 :     _chan_to_pin_map(other_mesh._chan_to_pin_map),
      44           0 :     _pin_to_chan_map(other_mesh._pin_to_chan_map),
      45           0 :     _sign_id_crossflow_map(other_mesh._sign_id_crossflow_map),
      46           0 :     _gij_map(other_mesh._gij_map),
      47           0 :     _subch_type(other_mesh._subch_type)
      48             : {
      49           0 :   if (_nx < 2 && _ny < 2)
      50           0 :     mooseError(name(),
      51             :                ": The number of subchannels cannot be less than 2 in both directions (x and y). "
      52             :                "Smallest assembly allowed is either 2X1 or 1X2. ");
      53           0 : }
      54             : 
      55             : std::unique_ptr<MooseMesh>
      56           0 : QuadSubChannelMesh::safeClone() const
      57             : {
      58           0 :   return _app.getFactory().copyConstruct(*this);
      59             : }
      60             : 
      61             : void
      62         206 : QuadSubChannelMesh::buildMesh()
      63             : {
      64         206 : }
      65             : 
      66             : void
      67         204 : QuadSubChannelMesh::computeAssemblyHydraulicParameters()
      68             : {
      69         204 :   _assembly_flow_area = 0.0;
      70         204 :   _assembly_wetted_perimeter = 0.0;
      71         204 :   _assembly_hydraulic_diameter = 0.0;
      72             : 
      73         204 :   const Real z = _z_grid.empty() ? 0.0 : _z_grid.front();
      74             : 
      75        3747 :   for (unsigned int i_ch = 0; i_ch < _n_channels; i_ch++)
      76             :   {
      77        3543 :     _assembly_flow_area += getSubchannelFlowArea(i_ch, z);
      78        3543 :     _assembly_wetted_perimeter += getSubchannelWettedPerimeter(i_ch);
      79             :   }
      80             : 
      81         204 :   if (_assembly_wetted_perimeter == 0.0)
      82           0 :     mooseError(name(), ": Assembly wetted perimeter is zero; cannot compute hydraulic diameter.");
      83             : 
      84         204 :   _assembly_hydraulic_diameter = 4.0 * _assembly_flow_area / _assembly_wetted_perimeter;
      85         204 : }
      86             : 
      87             : Real
      88       52497 : QuadSubChannelMesh::getSubchannelFlowArea(unsigned int i_chan, Real z) const
      89             : {
      90             :   Real standard_area = 0.0;
      91             :   Real rod_area = 0.0;
      92             :   Real additional_area = 0.0;
      93             : 
      94       52497 :   const auto subch_type = getSubchannelType(i_chan);
      95       52497 :   if (subch_type == EChannelType::CORNER)
      96             :   {
      97       12596 :     standard_area = 0.25 * _pitch * _pitch;
      98       12596 :     rod_area = 0.25 * 0.25 * libMesh::pi * _pin_diameter * _pin_diameter;
      99       12596 :     additional_area = _pitch * _side_gap + _side_gap * _side_gap;
     100             :   }
     101       39901 :   else if (subch_type == EChannelType::EDGE)
     102             :   {
     103       24446 :     standard_area = 0.5 * _pitch * _pitch;
     104       24446 :     rod_area = 0.5 * 0.25 * libMesh::pi * _pin_diameter * _pin_diameter;
     105       24446 :     additional_area = _pitch * _side_gap;
     106             :   }
     107             :   else
     108             :   {
     109       15455 :     standard_area = _pitch * _pitch;
     110       15455 :     rod_area = 0.25 * libMesh::pi * _pin_diameter * _pin_diameter;
     111             :   }
     112             : 
     113       52497 :   Real flow_area = standard_area + additional_area - rod_area;
     114             : 
     115             :   unsigned int blockage_index = 0;
     116      104994 :   for (const auto & i_blockage : _index_blockage)
     117             :   {
     118       52497 :     if (i_chan == i_blockage && z >= _z_blockage.front() && z <= _z_blockage.back())
     119         292 :       flow_area *= _reduction_blockage[blockage_index];
     120       52497 :     blockage_index++;
     121             :   }
     122             : 
     123       52497 :   return flow_area;
     124             : }
     125             : 
     126             : Real
     127       51457 : QuadSubChannelMesh::getSubchannelWettedPerimeter(unsigned int i_chan) const
     128             : {
     129       51457 :   const Real rod_circumference = libMesh::pi * _pin_diameter;
     130       51457 :   const auto subch_type = getSubchannelType(i_chan);
     131             : 
     132       51457 :   if (subch_type == EChannelType::CORNER)
     133       12196 :     return 0.25 * rod_circumference + _pitch + 2.0 * _side_gap;
     134       39261 :   else if (subch_type == EChannelType::EDGE)
     135       23966 :     return 0.5 * rod_circumference + _pitch;
     136             :   else
     137             :     return rod_circumference;
     138             : }
     139             : 
     140             : unsigned int
     141       98088 : QuadSubChannelMesh::getSubchannelIndexFromPoint(const Point & p) const
     142             : {
     143       98088 :   Real offset_x = (_nx - 1) * _pitch / 2.0;
     144       98088 :   Real offset_y = (_ny - 1) * _pitch / 2.0;
     145       98088 :   unsigned int i = (p(0) + offset_x + 0.5 * _pitch) / _pitch;
     146       98088 :   unsigned int j = (p(1) + offset_y + 0.5 * _pitch) / _pitch;
     147       98088 :   return j * _nx + i;
     148             : }
     149             : 
     150             : unsigned int
     151      875842 : QuadSubChannelMesh::channelIndex(const Point & pt) const
     152             : {
     153             :   // this is identical to getSubchannelIndexFromPoint, but when it is given a point "outside" the
     154             :   // normal subchannel geometry (i.e. a point that lies in a gap around the lattice) we still report
     155             :   // a valid subchannel index this is needed for transferring the solution onto a visualization mesh
     156             : 
     157      875842 :   Real offset_x = (_nx - 1) * _pitch / 2.0;
     158      875842 :   Real offset_y = (_ny - 1) * _pitch / 2.0;
     159      875842 :   int i = (pt(0) + offset_x + 0.5 * _pitch) / _pitch;
     160      875842 :   int j = (pt(1) + offset_y + 0.5 * _pitch) / _pitch;
     161             : 
     162      875842 :   i = std::max(0, i);
     163      875842 :   i = std::min(i, (int)(_nx - 1));
     164             : 
     165      875842 :   j = std::max(0, j);
     166      875842 :   j = std::min(j, (int)(_ny - 1));
     167             : 
     168      875842 :   return j * _nx + i;
     169             : }
     170             : 
     171             : unsigned int
     172       25888 : QuadSubChannelMesh::getPinIndexFromPoint(const Point & p) const
     173             : {
     174       25888 :   if (_n_pins == 0)
     175           0 :     mooseError(name(), ": Cannot compute a pin index because this mesh has no pins.");
     176             : 
     177       25888 :   Real offset_x = (_nx - 2) * _pitch / 2.0;
     178       25888 :   Real offset_y = (_ny - 2) * _pitch / 2.0;
     179       25888 :   unsigned int i = (p(0) + offset_x) / _pitch;
     180       25888 :   unsigned int j = (p(1) + offset_y) / _pitch;
     181       25888 :   return j * (_nx - 1) + i;
     182             : }
     183             : 
     184             : unsigned int
     185      860086 : QuadSubChannelMesh::pinIndex(const Point & p) const
     186             : {
     187      860086 :   if (_n_pins == 0)
     188           0 :     mooseError(name(), ": Cannot compute a pin index because this mesh has no pins.");
     189             : 
     190      860086 :   Real offset_x = (_nx - 2) * _pitch / 2.0;
     191      860086 :   Real offset_y = (_ny - 2) * _pitch / 2.0;
     192      860086 :   unsigned int i = (p(0) + offset_x) / _pitch;
     193      860086 :   unsigned int j = (p(1) + offset_y) / _pitch;
     194      860086 :   return j * (_nx - 1) + i;
     195             : }
     196             : 
     197             : void
     198          44 : QuadSubChannelMesh::generatePinCenters(
     199             :     unsigned int nx, unsigned int ny, Real pitch, Real elev, std::vector<Point> & pin_centers)
     200             : {
     201             :   mooseAssert(nx >= 2, "Number of channels in x-direction must be 2 or more.");
     202             :   mooseAssert(ny >= 2, "Number of channels in y-direction must be 2 or more.");
     203             : 
     204          44 :   Real offset_x = (nx - 2) * pitch / 2.0;
     205          44 :   Real offset_y = (ny - 2) * pitch / 2.0;
     206         167 :   for (unsigned int iy = 0; iy < ny - 1; iy++)
     207         656 :     for (unsigned int ix = 0; ix < nx - 1; ix++)
     208         533 :       pin_centers.push_back(Point(pitch * ix - offset_x, pitch * iy - offset_y, elev));
     209          44 : }