SNAP 8 Experimental Reactor (S8ER) Multiphysics Model

Contact: Isaac Naupa, iaguirre6.at.gatech.edu Contact: Stefano Terlizzi, sbt5572.at.psu.edu

Model link: S8ER Model

commentnote

For citing purposes, please cite (Garcia et al., 2022) and (Naupa et al., 2022).

Introduction

The SNAP8ER is a good stress test for the MOOSE-based computational tools, specifically GRIFFIN due to its small size, high leakage, sensitivity to few group parameters, and material composition.

Some aims of this model was to reduce the use of resources external to MOOSE, i.e. Cubit and to simplify the workflow from cross section generation to full griffin model. For this reason, the mesh is fully generated inside MOOSE using a variety of Mesh Generator Objects. Also for this reason, the DFEM-SN Transport Scheme, rather than the CFEM-Diffusion Scheme is used in the neutronic model. This allows us to have more a direct workflow from generation of the few group parameters and our griffin model, such as avoiding the generation of SPH correction factors, often needed for good agreement between griffin and reference solutions when using CFEM-Diffusion. DFEM-SN is also heavilyy optimized for parallel computation and comes in handy for larger problems such as this.

commentnote

To request access to Bison or Griffin, please submit a request on the INL Modeling and Software Website.

Generation of Few Group Parameters

For generation of the reference solution and neutronic data the monte carlo code Serpent was used. It must be noted that this model is quite sensitive to the energy structure and spatial resolution used in the generation of the few group parameters. Various sensitivity studies were conducted to find the ideal parameters that were conducive to good agreement between Griffin and Serpent.

Serpent Model

For the serpent model, only the the active core region including core barrel and internal reflectors were modeled. External reflector control drums are not included in this model, but will be included in future updates.

S8ER Core Snapshot, from (AID, 1973)

S8ER Serpent Model Core Snapshot, from (Garcia et al., 2022)

% --- SNAP 8 Core Model ---------

/*
General comments:
Core Modeling completed to allow for preliminary dry critical configuration
comparison with NAA-SR-9642. Work total excess reactivity is 5.7%, this was done
via implementing SS316 layering in core cylinder, Sm2O3 poisoning in fuel, and
removal of hexagonal vertices in reflector material.
*/

% --- Problem title:

set title "SNAP 8"

% --- Cross section library file path:

set acelib "/hpc-common/data/serpent/xsdata/s2v0_endfb71/s2v0_endfb71.xsdata"

% --- Materials:

include "SNAP8ER_f300c300r300.mat"
% ------------------------------------------------------------

/************************
 * Geometry definitions *
 ************************/
% ---------------------
% Fuel Pin Definitions
% ---------------------
% --- Fuel Pin (0.56in OD, 0.01in Hastelloy N clad, 0.0022in Ceramic Coating, 0.0016 He gap)
% Ceramic thickness can be found here: SNAP and Al Fuel Summary Report, pg. 4
% Other dimensions: NAA-SR-9642, pg. 19
pin 100
UZrH     0.67564
ceramic  0.681228
Sm2O3    0.6812587
intatm   0.6858
hasteN   0.7112
nak

pin 200
UZrH     0.67564
ceramic  0.681228
Sm2O3    0.6812587
intatm   0.6858
hasteN   0.7112
nak

pin 400
UZrH     0.67564
ceramic  0.681228
Sm2O3    0.6812587
intatm   0.6858
hasteN   0.7112
nak

pin 500
UZrH     0.67564
ceramic  0.681228
Sm2O3    0.6812587
intatm   0.6858
hasteN   0.7112
nak

pin 600
UZrH     0.67564
ceramic  0.681228
Sm2O3    0.6812587
intatm   0.6858
hasteN   0.7112
nak

pin 700
UZrH     0.67564
ceramic  0.681228
Sm2O3    0.6812587
intatm   0.6858
hasteN   0.7112
nak

pin 800
UZrH     0.67564
ceramic  0.681228
Sm2O3    0.6812587
intatm   0.6858
hasteN   0.7112
nak

pin 900
UZrH     0.67564
ceramic  0.681228
Sm2O3    0.6812587
intatm   0.6858
hasteN   0.7112
nak

% --- Void Pin
pin 300
nak

% ---------------------
% Surface Definitions
% ---------------------

% --- Cylindrical Surface for reactor core boundary, 9.352in OD (NAA-SR-9642, pg. 13)

surf S5 cyl 0.0 0.0 11.87704 % -18.3769 18.3769

surf SUG cyl 0.0 0.0 11.87704 % 18.3769  19.2500 % for upper grid plate
surf SLG cyl 0.0 0.0 11.87704 %-19.1707 -18.3769 % for lower grid plate
surf S12 cyl 0.0 0.0 11.718036 %-18.3769 18.3769
surf S13 cyl 0.0 0.0 11.93 %-18.3769 18.3769
surf S14 cyl 0.0 0.0 11.6926 %-18.3769 18.3769
surf SCube sqc 0.0 0.0 11.87704

% --- surfaces for drums
surf sDrum1 cyl  23.972012 0.0 11.9126% -18.3769 18.3769
surf sDrum4 cyl -23.972012 0.0 11.9126% -18.3769 18.3769
surf sDrum2 cyl  11.9860  20.7604 11.9126% -18.3769 18.3769
surf sDrum3 cyl -11.9860  20.7604 11.9126% -18.3769 18.3769
surf sDrum5 cyl -11.9860 -20.7604 11.9126% -18.3769 18.3769
surf sDrum6 cyl  11.9860 -20.7604 11.9126% -18.3769 18.3769

% --- surfaces for void near drums
surf sVDrum1 cyl  23.972012 0.0    11.95% -18.3769 18.3769
surf sVDrum4 cyl -23.972012 0.0    11.95% -18.3769 18.3769
surf sVDrum2 cyl  11.9860  20.7604 11.95% -18.3769 18.3769
surf sVDrum3 cyl -11.9860  20.7604 11.95% -18.3769 18.3769
surf sVDrum5 cyl -11.9860 -20.7604 11.95% -18.3769 18.3769
surf sVDrum6 cyl  11.9860 -20.7604 11.95% -18.3769 18.3769

% --- Cutoff at the end of hexagonal vertex for drums
surf sCut1 plane     0       20.9     0 436.81
surf sCut2 plane   -18.0999  10.450   0 436.81
surf sCut3 plane   -18.0999 -10.450   0 436.81
surf sCut4 plane     0      -20.9     0 436.81
surf sCut5 plane    18.0999 -10.45    0 436.81
surf sCut6 plane    18.0999  10.45    0 436.81

% --- Cutoff for stationary reflectors
surf sStatCut1 plane     0       18.9     0 357.21
surf sStatCut2 plane   -16.3679   9.450   0 357.21
surf sStatCut3 plane   -16.3679  -9.450   0 357.21
surf sStatCut4 plane     0      -18.9     0 357.21
surf sStatCut5 plane    16.3679  -9.45    0 357.21
surf sStatCut6 plane    16.3679   9.45    0 357.21

% --- surfaces for empty shims
surf sShimE1 plane  17.8206  0       0 317.5752
surf sShimE2 plane   8.9103  15.4331 0 317.5752
surf sShimE3 plane  -8.9103  15.4331 0 317.5752
surf sShimE4 plane -17.8206  0       0 317.5752
surf sShimE5 plane  -8.9103 -15.4331 0 317.5752
surf sShimE6 plane   8.9103 -15.4331 0 317.5752

% --- surfaces for B shims
surf sShimB1 plane  21.9354   0      0 481.1618
surf sShimB2 plane  10.9677  18.9966 0 481.1618
surf sShimB3 plane -10.9677  18.9966 0 481.1618
surf sShimB4 plane -21.9354   0      0 481.1618
surf sShimB5 plane -10.9677 -18.9966 0 481.1618
surf sShimB6 plane  10.9677 -18.9966 0 481.1618

%--- BEGIN CUBOID DEFINITIONS
surf sCuboid1 rect 19.7002 21.9354 -7.712 7.712 %-15.24 15.24
surf sCuboid2 rect 19.7002 21.9354 -7.712 7.712 %-15.24 15.24
surf sCuboid3 rect 19.7002 21.9354 -7.712 7.712 %-15.24 15.24
surf sCuboid4 rect 19.7002 21.9354 -7.712 7.712 %-15.24 15.24
surf sCuboid5 rect 19.7002 21.9354 -7.712 7.712 %-15.24 15.24
surf sCuboid6 rect 19.7002 21.9354 -7.712 7.712 %-15.24 15.24

% --- surfaces for internal reflectors
surf srefl1 plane  0      11.0668 0 122.4736
surf srefl2 plane -9.5841  5.5334 0 122.4736
surf srefl3 plane -9.5841 -5.5334 0 122.4736
surf srefl4 plane  0     -11.0668 0 122.4736
surf srefl5 plane  9.5841 -5.5334 0 122.4736
surf srefl6 plane  9.5841  5.5334 0 122.4736

% --- Begin Surface Definitions for Housing --- %
surf sHouseD1 cyl  23.972012 0.0    11.75385% -18.3769 18.3769
surf sHouseD4 cyl -23.972012 0.0    11.75385% -18.3769 18.3769
surf sHouseD2 cyl  11.9860  20.7604 11.75385% -18.3769 18.3769
surf sHouseD3 cyl -11.9860  20.7604 11.75385% -18.3769 18.3769
surf sHouseD5 cyl -11.9860 -20.7604 11.75385% -18.3769 18.3769
surf sHouseD6 cyl  11.9860 -20.7604 11.75385 -18.3769 18.3769
surf sHouseE1 plane  17.6619  0       0 311.9409
surf sHouseE2 plane   8.8309  15.2956 0 311.9409
surf sHouseE3 plane  -8.8309  15.2956 0 311.9409
surf sHouseE4 plane -17.6619  0       0 311.9409
surf sHouseE5 plane  -8.8309 -15.2956 0 311.9409
surf sHouseE6 plane   8.8309 -15.2956 0 311.9409
surf sHCut1 plane     0       20.7413  0 430.1995
surf sHCut2 plane   -17.9624  10.3706  0 430.1995
surf sHCut3 plane   -17.9624 -10.3706  0 430.1995
surf sHCut4 plane     0      -20.741   0 430.1995
surf sHCut5 plane    17.9624 -10.3706  0 430.1995
surf sHCut6 plane    17.9624  10.3706  0 430.1995
surf S8House hexxc 0.0 0.0 19.54145% -18.3769 18.3769
surf sHrefl1 plane  0      11.0414 0 121.9121
surf sHrefl2 plane -9.5621  5.5207 0 121.9121
surf sHrefl3 plane -9.5621 -5.5207 0 121.9121
surf sHrefl4 plane  0     -11.0414 0 121.9121
surf sHrefl5 plane  9.5621 -5.5207 0 121.9121
surf sHrefl6 plane  9.5621  5.5207 0 121.9121

% --- Hexagonal surfrace for reflector core boundaries, 9.352in OD + 0.0818 reflector radial thickness at thinnest point
%     + 4.68 drum radius =23.972012 cm radial distance flat to flat x-hexagonal (NAA-SR-9642, pg. 13)
%     Note that the actual thickest portion of the drum is noted as 3 inches which makes the half distance from flat point
%     to flat point 9.352 OD + 0.0818 + 3 in drum thickness = 19.704812 cm distance flat to flat x-hexagonal (AI-AEC-13070
%     Table 2 Sheet 2 of 4)

surf S8 hexxc 0.0 0.0 19.7002% -18.3769 18.3769

% ---------------------
% Cell Definitions
% ---------------------

% --- Latice x-type hexagonal, pitch = 0.57in (NAA-SR-9642, pg. 19)
%     Lattice universe is part of universe "core"
% lat Uni Type x_o y_o UNI
lat C3critload 2 0.0 0.0 21 21 1.4478
300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300
  300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300
    300 300 300 300 300 300 300 300 300 300 300 800 800 800 800 800 800 800 300 300 300 % top row, 7 inside, ending in position 4
      300 300 300 300 300 300 300 300 300 800 700 700 700 700 700 700 700 700 800 300 300 % 10 inside, ending in position 300
        300 300 300 300 300 300 300 300 800 700 600 600 600 600 600 600 600 700 800 300 300 % 11 inside, ending in position 300
          300 300 300 300 300 300 300 800 700 600 500 500 500 500 500 500 600 700 800 300 300 % 12 inside, ending in position 300
            300 300 300 300 300 300 800 700 600 500 400 400 400 400 400 500 600 700 800 300 300 % 1300 inside, ending in position 300
              300 300 300 300 300 800 700 600 500 400 900 900 900 900 400 500 600 700 800 300 300 % 14 inside, ending in position 300
                300 300 300 300 800 700 600 500 400 900 200 200 200 900 400 500 600 700 800 300 300 % 15 inside, ending in position 300
                  300 300 300 800 700 600 500 400 900 200 100 100 200 900 400 500 600 700 800 300 300 % 16 inside, ending in position 300
                    300 300 300 700 600 500 400 900 200 100 100 100 200 900 400 500 600 700 300 300 300 % middle row, 15 inside, starting in position 4
                      300 300 800 700 600 500 400 900 200 100 100 200 900 400 500 600 700 800 300 300 300 % 15 inside, starting in position 300
                        300 300 800 700 600 500 400 900 200 200 200 900 400 500 600 700 800 300 300 300 300 % 14 inside, starting in position 300
                          300 300 800 700 600 500 400 900 900 900 900 400 500 600 700 800 300 300 300 300 300 % 1300 inside, starting in position 300
                            300 300 800 700 600 500 400 400 400 400 400 500 600 700 800 300 300 300 300 300 300 % 12 inside, starting in position 300
                              300 300 800 700 600 500 500 500 500 500 500 600 700 800 300 300 300 300 300 300 300 % 11 inside, starting in position 300
                                300 300 800 700 600 600 600 600 600 600 600 700 800 300 300 300 300 300 300 300 300 % 10 inside, starting in position 300
                                  300 300 800 700 700 700 700 700 700 700 700 800 300 300 300 300 300 300 300 300 300 % 9 inside, starting in position 300
                                    300 300 300 800 800 800 800 800 800 800 300 300 300 300 300 300 300 300 300 300 300 % bottom row, 7 are inside, starting in position 4
                                      300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300
                                        300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300 300

% --- These cells define the reactor i.e. cutting off the "core"
%     universe with cylindrical boundaries

% --- fill definitions

cell cInternRefl1   1200  reflMix  srefl1 -S14% -sHouseZ1 sHouseZ2
cell cInternRefl2   1200  reflMix  srefl2 -S14% -sHouseZ1 sHouseZ2
cell cInternRefl3   1200  reflMix  srefl3 -S14% -sHouseZ1 sHouseZ2
cell cInternRefl4   1200  reflMix  srefl4 -S14% -sHouseZ1 sHouseZ2
cell cInternRefl5   1200  reflMix  srefl5 -S14% -sHouseZ1 sHouseZ2
cell cInternRefl6   1200  reflMix  srefl6 -S14% -sHouseZ1 sHouseZ2
cell cHouseRefl1    1200 reflMix (sHrefl1 -srefl1 -S12):(S14 -S12 srefl1)
cell cHouseRefl2    1200 reflMix (sHrefl2 -srefl2 -S12):(S14 -S12 srefl2)
cell cHouseRefl3    1200 reflMix (sHrefl3 -srefl3 -S12):(S14 -S12 srefl3)
cell cHouseRefl4    1200 reflMix (sHrefl4 -srefl4 -S12):(S14 -S12 srefl4)
cell cHouseRefl5    1200 reflMix (sHrefl5 -srefl5 -S12):(S14 -S12 sHrefl5)
cell cHouseRefl6    1200 reflMix (sHrefl6 -srefl6 -S12):(S14 -S12 srefl6)

cell cIntRefCore reactor fill 1200 (-S12 sHrefl1):(-S12 sHrefl2):(-S12 sHrefl3):(-S12 sHrefl4):(-S12 sHrefl5):(-S12 sHrefl6)
cell cCore reactor fill C3critload -S12 -sHrefl1 -sHrefl2 -sHrefl3 -sHrefl4 -sHrefl5 -sHrefl6
cell cCoreWall 1200  ss316 S12 -S5
cell cBarrel reactor fill 1200 S12 -S5
cell voidCell reactor void S5 -SCube

%surf cuty px 0
%surf cuthex plane 0 0 0 25.98076211350 15 0 0 0 1
% --- Cell cIN  is filled with universe "core", also its important to keep in mind that
%     the "0" universe is the universe for which outside needs to be defined.
%     Serpent gives the warning that the  '0' universe should be the only one defining outside
%     although this is not strictly true based on serpent-2 documentation.
%cell cIN 0 fill reactor cuty -cuthex -SCube
cell cIN 0 fill reactor -SCube

% --- Cell cOUT  is defined as everything outside the cubic cell
cell cOuT 0 outside SCube

% ------------------------------------------------------------

/******************
 * Run parameters *
 ******************/

% --- Boundary condition (1 = black, 2 = reflective, 3 = periodic)

set bc 1 1 2

% --- Neutron population: 100000 neutrons per cycle, 60 active / 20 inactive cycles

set pop 1000000 100 40

% --- XY-plot (3)

plot 31 1000 1000  %-19.0
plot 31 1000 1000  18.2
plot 31 1000 1000 -17.8
plot 21 1000 1000
%plot 11 1000 1000
% --- XY-meshplot (3), which is 700 by 700 pixels and covers the whole geometry

mesh 3 900 900
mesh 2 900 900
%branch Fhi  stp UZrH -6.0968 600
%branch dens repm UZrH UZrH_dens
%set power 1000
%dep daystep 10 20 30
set gcu 300 100 200 400 500 600 700 800 900 1200
ene sxtngroup 1 2.53E-08  1.00E-07  4.00E-07  1.00E-06  3.00E-06  1.00E-05  3.00E-05  1.00E-04  5.50E-04  3.00E-03  1.70E-02  1.00E-01  4.00E-01  9.00E-01  1.40E+00  3.00E+00  1.00E+01
set nfg sxtngroup

det ring1 de eqLet du 100
det ring2 de eqLet du 200
det ring3 de eqLet du 900
det ring4 de eqLet du 400
det ring5 de eqLet du 500
det ring6 de eqLet du 600
det ring7 de eqLet du 700
det ring8 de eqLet du 800
det ring9 de eqLet du 300
det ref de eqLet du 1200
ene eqLet 3 1500 1e-9 2e1

%det ringpowers dr -8 void dx -10.71 10.71 17 dy -10.71 10.71 17
(microreactors/s8er/serpent/SNAP8ER_C3_2D_NODRUMS_WET_1_1.i)

A radial ring spatial resolution was used for generation of the few group parameters 8 ring fuel layers, 1 coolant ring layer, and 1 internal reflector and core barrel ring layer for a total of 10 radial ring layers. Serpent uses universes to define spatial domains for generation of few group parameters, below is a table defining each universe. This information will be used in the griffin model to map neutronic data to the mesh.

Table 1: Radial Ring Universe Mapping

Radial LayerUniverse IDUniverse Material/Region
1100Fuel Channel
2200Fuel Channel
3900Fuel Channel
4400Fuel Channel
5500Fuel Channel
6600Fuel Channel
7700Fuel Channel
8800Fuel Channel
9300Coolant Channel
101200Internal Reflector + Core Barrel

A 16 energy few group structure documented in (Swenson, 1969) was used for the generation of few group parameters.

Table 2: Few Group Energy Structure

GroupLower Limit (eV)
010E+06
13E+06
21.4E+06
30.9E+06
40.4E+06
50.1E+06
617E+03
73E+03
80.55E+03
9100
1030
1110
123
131
140.4
150.1
160.025

State Points

Reference data at each state point are generated for each tabulation parameter: fuel and coolant temperature, and each tabulation point: 300K, 600K, 900K, 1200K. This leads to a total of state points. The reflector temperature and the coolant density are maintained constant for each statepoint. Future updates to the model will include these effects.

Table 3: Serpent Model State Points

CaseGrid IndexFuel Temperature (K)Coolant Temperature (K)
11 1300300
21 2300600
31 3300900
41 43001200
52 1600300
62 2600600
72 3600900
82 46001200
93 1900300
103 2900600
113 3900900
123 49001200
134 11200300
144 21200600
154 31200900
164 412001200

Preparation of Neutronic Data for Griffin

Now having all the reference neutronic data generated by serpent, the data must be prepared for use in our griffin model. This is done through ISOXML, a general toolkit aimed at handling large multigroup libraries. This is also comes in handy for facilitating the coupling process and handling feedback effects by using the data generated for our various state points.

<ISOXML Format="SERPENT2" Name="core_2D_nodrums_16G_WET_XS" NGroup="16" NDLGroup="6" Description="core_2D_nodrums_16G_WET">
<Controls>
   <ScatMatrixTol>1E-32</ScatMatrixTol>
   <ForceConsistency>false</ForceConsistency>
   <Debug>false</Debug>
   <EnergyUnit>MeV</EnergyUnit>
</Controls>

<Librarywise DataType="macro" L="1" LType="1">
  <Tabulation>Tfuel Tcool</Tabulation>
  <Tfuel>300 600 900 1200</Tfuel>
  <Tcool>300 600 900 1200</Tcool>
  <AddData>Total Fission Absorption NGamma</AddData>
  <Directory>./</Directory>
  <BaseFile>SNAP8ER_C3_2D_NODRUMS_WET</BaseFile>
  <BaseExtension>i</BaseExtension>
</Librarywise>

</ISOXML>
(microreactors/s8er/serpent/core_2D_nodrums_16G_WET.xml)

Full documentation on the serpent to ISOXML workflow can be found here. ISOXML will use the grid indexes for each state point and parse the serpent results files and construct the multigroup library for use in our griffin model.

Griffin Model

In serpent, the fuel, ceramic, absorber, gap, cladding, coolant, internal reflectors, and barrel are all modeled explicitly. However, the homogenized universe data is generated at the channel level, meaning in griffin the mesh will also be done at the channel level to correspond with the homogenized universe data.

Geometry

The main tools used to generate this mesh were MOOSE based meshgenerators, primarily the PolygonConcentricCircleMeshGenerator for a fuel channel, HexIDPatternMeshGenerator for a lattice of fuel channels, the PeripheralRingMeshGenerator for the internal reflector and barrel, and the PlaneDeletionGenerator for adjustments to the geometry.

[Mesh]
  [fuel_elem1]
    type = PolygonConcentricCircleMeshGenerator
    num_sides = '6'
    num_sectors_per_side = '2 2 2 2 2 2'
    polygon_size = '${hex_apothem}'
    preserve_volumes = true
    background_intervals = '2'
    background_block_ids = '1 1'
    quad_center_elements = true
    smoothing_max_it = 3
  []
  [fuel_elem2]
    type = PolygonConcentricCircleMeshGenerator
    num_sides = '6'
    num_sectors_per_side = '2 2 2 2 2 2'
    polygon_size = '${hex_apothem}'
    preserve_volumes = true
    background_intervals = '2'
    background_block_ids = '2 2'
    quad_center_elements = true
    smoothing_max_it = 3
  []
  [fuel_elem3]
    type = PolygonConcentricCircleMeshGenerator
    num_sides = '6'
    num_sectors_per_side = '2 2 2 2 2 2'
    polygon_size = '${hex_apothem}'
    preserve_volumes = true
    background_intervals = '2'
    background_block_ids = '3 3'
    quad_center_elements = true
    smoothing_max_it = 3
  []
  [fuel_elem4]
    type = PolygonConcentricCircleMeshGenerator
    num_sides = '6'
    num_sectors_per_side = '2 2 2 2 2 2'
    polygon_size = '${hex_apothem}'
    preserve_volumes = true
    background_intervals = '2'
    background_block_ids = '4 4'
    quad_center_elements = true
    smoothing_max_it = 3
  []
  [fuel_elem5]
    type = PolygonConcentricCircleMeshGenerator
    num_sides = '6'
    num_sectors_per_side = '2 2 2 2 2 2'
    polygon_size = '${hex_apothem}'
    preserve_volumes = true
    background_intervals = '2'
    background_block_ids = '5 5'
    quad_center_elements = true
    smoothing_max_it = 3
  []
  [fuel_elem6]
    type = PolygonConcentricCircleMeshGenerator
    num_sides = '6'
    num_sectors_per_side = '2 2 2 2 2 2'
    polygon_size = '${hex_apothem}'
    preserve_volumes = true
    background_intervals = '2'
    background_block_ids = '6 6'
    quad_center_elements = true
    smoothing_max_it = 3
  []
  [fuel_elem7]
    type = PolygonConcentricCircleMeshGenerator
    num_sides = '6'
    num_sectors_per_side = '2 2 2 2 2 2'
    polygon_size = '${hex_apothem}'
    preserve_volumes = true
    background_intervals = '2'
    background_block_ids = '7 7'
    quad_center_elements = true
    smoothing_max_it = 3
  []
  [fuel_elem8]
    type = PolygonConcentricCircleMeshGenerator
    num_sides = '6'
    num_sectors_per_side = '2 2 2 2 2 2'
    polygon_size = '${hex_apothem}'
    preserve_volumes = true
    background_intervals = '2'
    background_block_ids = '8 8'
    quad_center_elements = true
    smoothing_max_it = 3
  []
  [air_elem1]
    type = PolygonConcentricCircleMeshGenerator
    num_sides = '6'
    num_sectors_per_side = '2 2 2 2 2 2'
    polygon_size = '${hex_apothem}'
    preserve_volumes = true
    background_intervals = '2'
    background_block_ids = '9 9'
    quad_center_elements = true
    smoothing_max_it = 3
  []
  [fuel_assem]
    type = PatternedHexMeshGenerator
    inputs = 'fuel_elem1 fuel_elem2 fuel_elem3 fuel_elem4 fuel_elem5 fuel_elem6 fuel_elem7 fuel_elem8 air_elem1'
    # Pattern ID  0 1 2 3 4 5 6 7 8
    background_intervals = 1
    background_block_id = '10'
    background_block_name = 'AIR'
    hexagon_size = '${active_apothem}'
    hexagon_size_style = 'apothem'
    pattern_boundary = 'hexagon'
    pattern = '8 7 7 7 7 7 7 7 8;
            7 6 6 6 6 6 6 6 6 7 ;
            7 6 5 5 5 5 5 5 5 6 7 ;
            7 6 5 4 4 4 4 4 4 5 6 7 ;
            7 6 5 4 3 3 3 3 3 4 5 6 7 ;
            7 6 5 4 3 2 2 2 2 3 4 5 6 7 ;
            7 6 5 4 3 2 1 1 1 2 3 4 5 6 7 ;
            7 6 5 4 3 2 1 0 0 1 2 3 4 5 6 7 ;
            8 6 5 4 3 2 1 0 0 0 1 2 3 4 5 6 8 ;
            7 6 5 4 3 2 1 0 0 1 2 3 4 5 6 7 ;
            7 6 5 4 3 2 1 1 1 2 3 4 5 6 7 ;
            7 6 5 4 3 2 2 2 2 3 4 5 6 7 ;
            7 6 5 4 3 3 3 3 3 4 5 6 7 ;
            7 6 5 4 4 4 4 4 4 5 6 7 ;
            7 6 5 5 5 5 5 5 5 6 7 ;
            7 6 6 6 6 6 6 6 6 7 ;
            8 7 7 7 7 7 7 7 8'
    external_boundary_id = '10'
    external_boundary_name = 'active_outer'
  []
  [int_ref]
    type = PeripheralRingMeshGenerator
    input = fuel_assem
    peripheral_ring_radius = '${nonactive_corner_half}'
    input_mesh_external_boundary = 'active_outer'
    external_boundary_name = 'core_outer'
    external_boundary_id = '11'
    peripheral_ring_block_id = 11
    peripheral_ring_block_name = int_ref
    peripheral_layer_num = 10
  []
  [cut_1]
    type = PlaneDeletionGenerator
    point = '0 ${cut_apothem} 0'
    normal = '0 ${cut_apothem} 0'
    input = int_ref
    new_boundary = 'core_outer'
  []
  [cut_2]
    type = PlaneDeletionGenerator
    point = '${cut_x} ${cut_y} 0'
    normal = '${cut_x} ${cut_y} 0'
    input = cut_1
    new_boundary = 'core_outer'
  []
  [cut_3]
    type = PlaneDeletionGenerator
    point = '${cut_x} -${cut_y} 0'
    normal = '${cut_x} -${cut_y} 0'
    input = cut_2
    new_boundary = 'core_outer'
  []
  [cut_4]
    type = PlaneDeletionGenerator
    point = '0 -${cut_apothem} 0'
    normal = '0 -${cut_apothem} 0'
    input = cut_3
    new_boundary = 'core_outer'
  []
  [cut_5]
    type = PlaneDeletionGenerator
    point = '-${cut_x} -${cut_y} 0'
    normal = '-${cut_x} -${cut_y} 0'
    input = cut_4
    new_boundary = 'core_outer'
  []
  [cut_6]
    type = PlaneDeletionGenerator
    point = '-${cut_x} ${cut_y} 0'
    normal = '-${cut_x} ${cut_y} 0'
    input = cut_5
    new_boundary = 'core_outer'
  []
  [scube_1]
    type = PlaneDeletionGenerator
    point = '${scube_apothem} 0  0'
    normal = '${scube_apothem} 0  0'
    input = cut_6
    new_boundary = 'core_outer'
  []
  [scube_2]
    type = PlaneDeletionGenerator
    point = '${scube_x} ${scube_y} 0'
    normal = '${scube_x} ${scube_y} 0'
    input = scube_1
    new_boundary = 'core_outer'
  []
  [scube_3]
    type = PlaneDeletionGenerator
    point = '${scube_x} -${scube_y} 0'
    normal = '${scube_x} -${scube_y} 0'
    input = scube_2
    new_boundary = 'core_outer'
  []
  [scube_4]
    type = PlaneDeletionGenerator
    point = '-${scube_apothem} 0 0'
    normal = '-${scube_apothem} 0 0'
    input = scube_3
    new_boundary = 'core_outer'
  []
  [scube_5]
    type = PlaneDeletionGenerator
    point = '-${scube_x} -${scube_y} 0'
    normal = '-${scube_x} -${scube_y} 0'
    input = scube_4
    new_boundary = 'core_outer'
  []
  [scube_6]
    type = PlaneDeletionGenerator
    point = '-${scube_x} ${scube_y} 0'
    normal = '-${scube_x} ${scube_y} 0'
    input = scube_5
    new_boundary = 'core_outer'
  []
  [set_mat_id]
    type = SubdomainExtraElementIDGenerator
    input = 'scube_6'
    subdomains = '1 2 3 4 5 6 7 8 9 10 11'
    extra_element_id_names = 'material_id'
    extra_element_ids = '100 200 900 400 500 600 700 800 300 300 1200'
  []
[]
(microreactors/s8er/core_2D_griffin_coupled.i)

  1. 9 individual fuel channels are modeled, 1 for each corresponding radial ring

  2. The fuel channels are mapped into the lattice in the appropiate arrangement from the homogenized universe data.

  3. A peripheral ring is enclosed around the lattice to represent the internal reflector + core barrel (last radial layer)

  4. The mesh is cut along the ring at various points to preserve the geometry and dimension of the reference serpent model.

  5. The blocks are mapped to their appropriate material ID.

A table with corresponding mapping between the Serpent Universe ID and the Griffin Block ID is shown below.

Table 4: Homogenized Universe Mapping

Serpent Universe IDMaterial IDBlock IDMaterial/Region
1001001Fuel Channel
2002002Fuel Channel
9009003Fuel Channel
4004004Fuel Channel
5005005Fuel Channel
6006006Fuel Channel
7007007Fuel Channel
8008008Fuel Channel
3003009, 10Coolant Channel
1200120011Internal Reflector + Core Barrel

S8ER Griffin Mesh, from (Naupa et al., 2022)

TransportSystems

As stated previously, the Griffin model uses the newly optimized DFEM-SN transport scheme to solve for the eigenvalue. Here a 16 group energy structure is used with vacuum boundary conditions on the outer edge of the core barrel. The solver uses 9 azimuthal and 3 polar angles, typically the default for these are 6 and 2 respectively but because this system is highly anisotropic, a higher order of angles is required for solution fidelity. The benefit with using the DFEM-SN scheme rather than the stand CFEM-Diffusion scheme, is the removal of several layers of calculation for the same fidelity, however at the expense of a more expensive computation. The DFEM-SN scheme is optimized and scales well with problem size.

[TransportSystems]
  # Eigenvalue Problem
  particle = neutron
  equation_type = eigenvalue

  # Sixteen Group XS Structure
  G = 16
  VacuumBoundary = 'core_outer'

  # DFEM-Transport
  [dsn]
    scheme = DFEM-SN
    #n_delay_groups = 6
    family = MONOMIAL
    order = FIRST
    AQtype = Gauss-Chebyshev
    NPolar = 3 # use >=2 for final runs (4 sawtooth nodes sufficient)
    NAzmthl = 9 # use >=6 for final runs (4 sawtooth nodes sufficient)
    NA = 1
    sweep_type = asynchronous_parallel_sweeper
    using_array_variable = true
    collapse_scattering = true
    hide_angular_flux = true
    # verbose = 3
  []
[]
(microreactors/s8er/core_2D_griffin_coupled.i)

Materials

Through the use of the CoupledFeedbackMatIDNeutronicsMaterial, by the mapping of materialIDs to the mesh done earlier, this object conveniently handles the direct mapping of neutronic data into the mesh. The object also has a feature to correct the neutronic data used in a material by adjusting the cross sections for a correction volume. This is done for the internal reflector and core barrel region to more appropriately match the reference model.

[Materials]
  [core]
    type = CoupledFeedbackMatIDNeutronicsMaterial
    block = '1 2 3 4 5 6 7 8 9 10'
    library_file = './serpent/core_2D_nodrums_16G_WET_XS.xml'
    library_name = 'core_2D_nodrums_16G_WET_XS'
    isotopes = 'pseudo'
    densities = '1.0'
    plus = 1
    is_meter = true
    grid_names = 'Tfuel Tcool'
    grid_variables = 'griffin_Tfuel griffin_Tcool'
    #volume_adjuster = core_vol_adjust  #to be fixed
  []
  [int_ref]
    type = CoupledFeedbackMatIDNeutronicsMaterial
    block = '11'
    library_file = './serpent/core_2D_nodrums_16G_WET_XS.xml'
    library_name = 'core_2D_nodrums_16G_WET_XS'
    isotopes = 'pseudo'
    densities = '1.0'
    plus = 1
    is_meter = true
    grid_names = 'Tfuel Tcool'
    grid_variables = 'griffin_Tfuel griffin_Tcool'
    volume_adjuster = int_ref_vol_adjust
  []
[]
(microreactors/s8er/core_2D_griffin_coupled.i)

Bison Model

Because bison is only taking the power distribution from the griffin model, we have the freedom to model each fuel channel component explicitly and hence we are able to capture unique heat transfer phenomena such as the gap radiative transfer and model thermal properties for each component. This model only captures heat conduction but future updates will include feedback from thermomechanical effects such as thermal expansion.

S8ER Bison Mesh, from (Naupa et al., 2022)

Geometry

The process for generating the Bison mesh is almost identical to the process described for the Griffin mesh, with a few extra steps. This time since we are modeling the explicit fuel element geometry, we use the PolygonConcentricCircleMeshGenerator to model concentric rings corresponding to the outer radii of each material in the fuel element.

[Mesh]
  [fuel_elem1]
    type = PolygonConcentricCircleMeshGenerator
    num_sides = '6'
    num_sectors_per_side = '4 4 4 4 4 4'
    polygon_size = '${hex_apothem}'
    preserve_volumes = true
    ring_radii = '${fuel_radius} ${ceramic_radius_outer} ${gap_radius_outer} ${clad_radius_outer} '
    ring_intervals = '1 1 1 1'
    ring_block_ids = '1 2 3 4'
    background_intervals = '1'
    background_block_ids = '5'
    quad_center_elements = true
    smoothing_max_it = 3
  []
  [coolant_elem1]
    type = PolygonConcentricCircleMeshGenerator
    num_sides = '6'
    num_sectors_per_side = '4 4 4 4 4 4'
    polygon_size = '${hex_apothem}'
    preserve_volumes = true
    background_intervals = '1'
    background_block_ids = '6'
    quad_center_elements = true
    smoothing_max_it = 3
    external_boundary_id = '6'
  []
  [fuel_assem]
    type = PatternedHexMeshGenerator
    inputs = 'fuel_elem1 coolant_elem1'
    # Pattern ID  0 1
    background_intervals = 1
    background_block_id = '7'
    background_block_name = 'coolant'
    hexagon_size = '${active_apothem}'
    hexagon_size_style = 'apothem'
    pattern_boundary = 'hexagon'
    pattern = '1 0 0 0 0 0 0 0 1;
            0 0 0 0 0 0 0 0 0 0 ;
            0 0 0 0 0 0 0 0 0 0 0 ;
            0 0 0 0 0 0 0 0 0 0 0 0 ;
            0 0 0 0 0 0 0 0 0 0 0 0 0 ;
            0 0 0 0 0 0 0 0 0 0 0 0 0 0 ;
            0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ;
            0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ;
            1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 ;
            0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ;
            0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ;
            0 0 0 0 0 0 0 0 0 0 0 0 0 0 ;
            0 0 0 0 0 0 0 0 0 0 0 0 0 ;
            0 0 0 0 0 0 0 0 0 0 0 0 ;
            0 0 0 0 0 0 0 0 0 0 0 ;
            0 0 0 0 0 0 0 0 0 0 ;
            1 0 0 0 0 0 0 0 1'
  []
  [block_rename1]
    type = RenameBlockGenerator
    input = fuel_assem
    old_block = '1 2 3 4 5 6 7'
    new_block = 'fuel ceramic gap clad coolant coolant coolant'
  []
  [bound_rename1]
    type = RenameBoundaryGenerator
    input = block_rename1
    old_boundary = '3 5 7'
    new_boundary = 'gap_inner gap_outer clad_outer'
  []
  [inner_gap1]
    type = SideSetsBetweenSubdomainsGenerator
    input = bound_rename1
    primary_block = ceramic
    paired_block = gap
    new_boundary = 'gap_inner'
  []
  [outer_gap1]
    type = SideSetsBetweenSubdomainsGenerator
    input = inner_gap1
    primary_block = clad
    paired_block = gap
    new_boundary = 'gap_outer'
  []
  [del_gap1]
    type = BlockDeletionGenerator
    block = gap
    input = outer_gap1
  []
  [int_ref]
    type = PeripheralRingMeshGenerator
    input = del_gap1
    peripheral_ring_radius = '${nonactive_corner_half}'
    input_mesh_external_boundary = 10000
    peripheral_ring_block_id = 8
    peripheral_ring_block_name = int_ref
    peripheral_layer_num = 10
  []
  [cut_1]
    type = PlaneDeletionGenerator
    point = '0 ${cut_apothem} 0'
    normal = '0 ${cut_apothem} 0'
    input = int_ref
    new_boundary = 10000
  []
  [cut_2]
    type = PlaneDeletionGenerator
    point = '${cut_x} ${cut_y} 0'
    normal = '${cut_x} ${cut_y} 0'
    input = cut_1
    new_boundary = 10000
  []
  [cut_3]
    type = PlaneDeletionGenerator
    point = '${cut_x} -${cut_y} 0'
    normal = '${cut_x} -${cut_y} 0'
    input = cut_2
    new_boundary = 10000
  []
  [cut_4]
    type = PlaneDeletionGenerator
    point = '0 -${cut_apothem} 0'
    normal = '0 -${cut_apothem} 0'
    input = cut_3
    new_boundary = 10000
  []
  [cut_5]
    type = PlaneDeletionGenerator
    point = '-${cut_x} -${cut_y} 0'
    normal = '-${cut_x} -${cut_y} 0'
    input = cut_4
    new_boundary = 10000
  []
  [cut_6]
    type = PlaneDeletionGenerator
    point = '-${cut_x} ${cut_y} 0'
    normal = '-${cut_x} ${cut_y} 0'
    input = cut_5
    new_boundary = 10000
  []
  [scube_1]
    type = PlaneDeletionGenerator
    point = '${scube_apothem} 0  0'
    normal = '${scube_apothem} 0  0'
    input = cut_6
    new_boundary = 10000
  []
  [scube_2]
    type = PlaneDeletionGenerator
    point = '${scube_x} ${scube_y} 0'
    normal = '${scube_x} ${scube_y} 0'
    input = scube_1
    new_boundary = 10000
  []
  [scube_3]
    type = PlaneDeletionGenerator
    point = '${scube_x} -${scube_y} 0'
    normal = '${scube_x} -${scube_y} 0'
    input = scube_2
    new_boundary = 10000
  []
  [scube_4]
    type = PlaneDeletionGenerator
    point = '-${scube_apothem} 0 0'
    normal = '-${scube_apothem} 0 0'
    input = scube_3
    new_boundary = 10000
  []
  [scube_5]
    type = PlaneDeletionGenerator
    point = '-${scube_x} -${scube_y} 0'
    normal = '-${scube_x} -${scube_y} 0'
    input = scube_4
    new_boundary = 10000
  []
  [scube_6]
    type = PlaneDeletionGenerator
    point = '-${scube_x} ${scube_y} 0'
    normal = '-${scube_x} ${scube_y} 0'
    input = scube_5
    new_boundary = 10000
  []
  [bound_rename2]
    type = RenameBoundaryGenerator
    input = scube_6
    old_boundary = '10000'
    new_boundary = 'void_outer'
  []
[]
(microreactors/s8er/core_2D_bison_coupled.i)

Table 5: S8ER Bison Dimensions

Material/RegionRadius (m)
UZrH Radius0.0067564)
Ceramic Outer Radius0.007130879
Gap Outer Radius0.006866719
Cladding Outer Radius0.007130879
Lattice Pitch0.007239

Kernels, AuxKernels, Functions

The Kernels used are the HeatConduction to solve the steady state heat conduction equation, as well as the CoupledForce to couple the power distribution from bison as the spatially dependent heat generation term.

[Kernels]
  # Equation 0 (Heat conduction equation with heat generation).
  [heat_conduction]
    type = ADHeatConduction
    variable = bison_temp
  []
  [heat_source]
    type = ADCoupledForce
    variable = bison_temp
    v = bison_norm_power_density
    block = fuel
  []
[]
(microreactors/s8er/core_2D_bison_coupled.i)

The AuxKernels used are a variety of function executors and variable normalizers for solving thermal properties that are temperature dependent, normalizing the power distribution for power conservation between Bison and Griffin, and for seperating the fuel, coolant, and reflector temperatures into individual variables to be used as the tabulation parameters.

[AuxKernels]

  [aux_tk]
    type = FunctionAux
    function = tk_f
    variable = aux_tk
    block = fuel
  []
  [aux_cp]
    type = FunctionAux
    function = cp_f
    variable = aux_cp
    block = fuel
  []
  [aux_T_inf]
    type = FunctionAux
    function = T_inf_f
    variable = aux_T_inf
    block = clad
  []
  [norm_Tfuel]
    type = NormalizationAux
    variable = bison_Tfuel
    source_variable = bison_temp
    execute_on = 'timestep_end' #check
  []
  [norm_Tcool]
    type = NormalizationAux
    variable = bison_Tcool
    source_variable = bison_temp
    execute_on = 'timestep_end' #check
  []
  [norm_Tintref]
    type = NormalizationAux
    variable = bison_Tintref
    source_variable = bison_temp
    execute_on = 'timestep_end' #check
  []
  [norm_power_density]
    type = NormalizationAux
    variable = bison_norm_power_density
    source_variable = bison_power_density
    normal_factor = 1.26359177225 #1.265805981
    execute_on = 'timestep_begin' #check
  []
[]
(microreactors/s8er/core_2D_bison_coupled.i)

The Functions used include the thermal conductivity, and specific heat of the fuel as a function of temperature from (Swenson, 1969).

[Functions]

  [tk_f]
    type = ParsedFunction
    symbol_names = 'bison_temp'
    symbol_values = 'temp_av'
    expression = '27.73992+bison_temp*0.027428444' #Models/SNAP10A_dimensions
  []
  [cp_f]
    type = ParsedFunction
    symbol_names = 'bison_temp'
    symbol_values = 'temp_av'
    expression = '472.27104+bison_temp*0.7275728' #Models/SNAP10A_dimensions
  []
  [T_inf_f]
    type = ParsedFunction
    expression = '${inlet_T_fluid}'
  []
[]
(microreactors/s8er/core_2D_bison_coupled.i)

Materials

The thermal properties of each material are modeled through using GenericConstantMaterial. This allows you to set thermal properties such as the density, specific heat, and thermal conductivity. In the case of the fuel we have documented formulations for the specific heat and thermal conductivity as a function of temperature. This temperature dependence can implemented through the use of the HeatConductionMaterial which allows you to set material properties as a function of temperature.

[Materials]
  [fuel_thermal_conduction]
    # Metallic Uranium Zirc Hydride fuel conduction params.
    type = ADHeatConductionMaterial
    temp = bison_temp
    specific_heat_temperature_function = cp_f
    thermal_conductivity_temperature_function = tk_f
    block = fuel
  []
  # Metallic Uranium Zirc Hydride fuel density
  [fuel_density]
    type = ADGenericConstantMaterial
    prop_names = 'density'
    prop_values = '${fuel_density}'
    block = fuel
  []
  # Fluid ambient temperature.
  # [T_infinity]
  #     #Models/SNAP10A_dimensions
  #     type = ADCoupledValueFunctionMaterial
  #     prop_name = 'T_infinity'
  #     function = 'if(x<180, 293.15+2.908240722*x, 816.48333)'
  #     v = 'aux_time'
  #     block = 2
  #     boundary = 'right'
  # []
  [clad_thermal_conduction]
    #Models/SNAP10A_dimensions
    type = ADGenericConstantMaterial
    prop_names = 'density specific_heat thermal_conductivity'
    prop_values = '${clad_density} ${clad_cp} ${ clad_tc}'
    block = clad
  []
  [gap_heat_transfer]
    #Models/SNAP10A_dimensions
    type = GenericConstantMaterial
    prop_names = 'gap_conductance gap_conductance_dT'
    prop_values = '${gap_tc} 0.0'
    boundary = 'gap_inner gap_outer'
  []
  # [gap_thermal_conduction]
  #     #Models/SNAP10A_dimensions
  #     type = ADGenericConstantMaterial
  #     prop_names = 'density specific_heat thermal_conductivity'
  #     prop_values = '${gap_dens} ${gap_cp} ${gap_tc}'
  #     block = gap
  # []
  [ceramic_thermal_conduction]
    #Models/SNAP10A_dimensions
    type = ADGenericConstantMaterial
    prop_names = 'density specific_heat thermal_conductivity'
    prop_values = '${ceramic_dens} ${ ceramic_cp} ${ ceramic_tc}'
    block = ceramic
  []
  [intref_thermal_conduction]
    #Models/SNAP10A_dimensions
    type = ADGenericConstantMaterial
    prop_names = 'density specific_heat thermal_conductivity'
    prop_values = '${intref_dens} ${intref_cp} ${ intref_tc}'
    block = int_ref
  []
  # [air_thermal_conduction]
  #     #Models/SNAP10A_dimensions
  #     type = ADGenericConstantMaterial
  #     prop_names = 'density specific_heat thermal_conductivity'
  #     prop_values = '${air_dens} ${air_cp} ${air_tc}'
  #     block = air
  # []
  [coolant_thermal_conduction]
    #Models/SNAP10A_dimensions
    type = ADGenericConstantMaterial
    prop_names = 'density specific_heat thermal_conductivity'
    prop_values = '${coolant_dens} ${coolant_cp} ${coolant_tc}'
    block = coolant
  []
[]
(microreactors/s8er/core_2D_bison_coupled.i)

Table 6: S8ER Thermal Parameters, from (Swenson, 1969), (AID, 1973)

Material/RegionDensity (kg/m^3)Specfic Heat (J/kg-K)Thermal Conductivity (W/m-K)
UZrH5963.13472.27 + T(0.72)27.73 + T(0.027)
Ceramic2242.584837.361.73
Gap0.01665193.160.34
Cladding8617.93418.6818.85
Internal Reflector1810.082721.42131.53

Multiphysics Coupling

The multiphysics coupling is done through MOOSE's MultiApp system, where Griffin is the main driver and executes Bison as a subapp. In this model the power distribution is solved by Griffin, where it is then transferred to Bison to solve the temperature distribution, where then the temperature distribution is transferred back to Griffin to calculate the new power distribution with updated neutronic data from the temperature feedback. This process is iterated by the [Executioner] until the system converges.

Multiphysics Flow Diagram, from (Naupa et al., 2022)

[MultiApps]
  [bison_coupled]
    type = FullSolveMultiApp
    app_type = GriffinApp
    input_files = 'core_2D_bison_coupled.i'
    positions = '0 0 0'
  []
[]
(microreactors/s8er/core_2D_griffin_coupled.i)

The data from each respective model is transfered using MOOSE's Transfers system. Since the meshes are tailored for each application, the differences in each mesh must be accounted for in the transfers. As described earlier, the Griffin mesh is generated at the channel level, while the Bison mesh is generated with explicit fuel components modeled. For this reason the MultiAppProjectionTransfer is used, which is handy when projecting field variables from one mesh to another. In this case the power generated in the Griffin fuel channel is transfered to only the fuel in Bison. For this reason, a volume normalization factor must be applied to the power density transfered so that power between models is conserved. The normalization factor is simply the ratio of the fuel channel to fuel volume.

The temperature distribution for each respective tabulation parameter i.e, fuel temperature and coolant temperature is transfered using the MultiAppInterpolationTransfer. This transfer is suitable for nodal variables such as variables in the LAGRANGE family.

[Transfers]
  [to_bison_power_density]
    type = MultiAppProjectionTransfer
    to_multi_app = bison_coupled
    source_variable = griffin_power_density
    variable = bison_power_density
  []
  [from_bison_Tfuel]
    type = MultiAppGeometricInterpolationTransfer
    from_multi_app = bison_coupled
    source_variable = bison_Tfuel
    variable = griffin_Tfuel
  []
  [from_bison_Tcool]
    type = MultiAppGeometricInterpolationTransfer
    from_multi_app = bison_coupled
    source_variable = bison_Tcool
    variable = griffin_Tcool
  []
  [from_bison_Tintref]
    type = MultiAppGeometricInterpolationTransfer
    from_multi_app = bison_coupled
    source_variable = bison_Tintref
    variable = griffin_Tintref
  []
[]
(microreactors/s8er/core_2D_griffin_coupled.i)

Referring back to the [Materials] block, this is where the coupled feedback is taken into account. The CoupledFeedbackMatIDNeutronicsMaterial takes the tabulated neutronic data and uses temperature from bison to interpolate between statepoints. For this reason, Multiple transfers for each tabulation parameter are required so that the neutronic data can be properly interpolated between statepoints.

[Materials]
  [core]
    type = CoupledFeedbackMatIDNeutronicsMaterial
    block = '1 2 3 4 5 6 7 8 9 10'
    library_file = './serpent/core_2D_nodrums_16G_WET_XS.xml'
    library_name = 'core_2D_nodrums_16G_WET_XS'
    isotopes = 'pseudo'
    densities = '1.0'
    plus = 1
    is_meter = true
    grid_names = 'Tfuel Tcool'
    grid_variables = 'griffin_Tfuel griffin_Tcool'
    #volume_adjuster = core_vol_adjust  #to be fixed
  []
  [int_ref]
    type = CoupledFeedbackMatIDNeutronicsMaterial
    block = '11'
    library_file = './serpent/core_2D_nodrums_16G_WET_XS.xml'
    library_name = 'core_2D_nodrums_16G_WET_XS'
    isotopes = 'pseudo'
    densities = '1.0'
    plus = 1
    is_meter = true
    grid_names = 'Tfuel Tcool'
    grid_variables = 'griffin_Tfuel griffin_Tcool'
    volume_adjuster = int_ref_vol_adjust
  []
[]
(microreactors/s8er/core_2D_griffin_coupled.i)

Running the Model

The Discontinuous Finite Element Discrete Ordinate (DFEM-SN) solver in conjunction with the asynchronous parallel sweeper implemented in griffin for the inversion of the streaming operator allows to efficiently perform neutron transport calculations in a parallel fashion on unstructured mesh (Wang et al., 2021). Additionally, it allows to streamline the modeling and simulation pipeline due to the enhanced fidelity that renders the generation of equivalence parameters not essential (agreement with Serpent within few hundreds of pcm for the effective multiplication factor). The neutronic model coupled with heat conduction runs in one minute on one Sawtooth node composed of 48 processors.

[Executioner]
  type = SweepUpdate
  verbose = true

  richardson_max_its = 50
  richardson_value = eigenvalue
  richardson_rel_tol = 1e-8
  richardson_abs_tol = 1e-8

  inner_solve_type = GMRes
  max_inner_its = 2

  fixed_point_max_its = 2
  custom_pp = eigenvalue
  custom_rel_tol = 1e-4
  force_fixed_point_solve = true
[]
(microreactors/s8er/core_2D_griffin_coupled.i)

The code is executed with the following lines.


mpiexec -n 48 griffin-opt -i core_2D_griffin_coupled.i
commentnote

More information about the SNAP Reactors may be found here.

References

  1. AID. Snap 8 summary report. 9 1973. URL: https://www.osti.gov/biblio/4393793, doi:10.2172/4393793.[BibTeX]
  2. Samuel Garcia, Isaac Naupa, Dan Kotlyar, and Ben Lindley. Validation of snap8 criticality configuration experiments using serpent. In Proceedings of ANS Winter Conference. 1 2022.[BibTeX]
  3. Isaac Naupa, Samuel Garcia, Stefano Terlizzi, Dan Kotlyar, and Ben Lindley. Validation of snap8 criticality configuration experiments using neams tools. In Submitted to Proceedings of M&C. 1 2022.[BibTeX]
  4. L D Swenson. Snap 8 development reactor nuclear analysis. 1 1969. URL: https://www.osti.gov/biblio/4709921, doi:10.2172/4709921.[BibTeX]
  5. Yaqi Wang, Zachary M. Prince, Joshua Thomas Hanophy, Changho Lee, Yeon Sang Jung, Hansol Park, Logan H. Harbour, and Javier Ortensi. Performance improvements for the griffin transport solvers. 9 2021. URL: https://www.osti.gov/biblio/1822446, doi:10.2172/1822446.[BibTeX]