neams-th-coe / cardinal Goto Github PK

Currently, CellTemperatureAux reads a cell's temperature by calling openmc_cell_get_temperature. However, temperature can only be read for cells that have a material fill. For a TRISO pebble example I am writing, because we use the set_contained feature to set the temperature of all the included cells in the pebble, when we try to call openmc_cell_get_temperature, we encounter an error due to the cell not being filled with a material.

We should make CellTemperatureAux more robust by reading the temperature of the first filled material-type cell (which is guaranteed to be the same as all the other filled cells due to the nature with which we use set_contained = true).

For coupling to SAM, we will need postprocessors in Cardinal that integrate nekRS's pressure along boundaries. This will be implemented by simply adding pressure as part of the NekFieldEnum enumeration, so that we can also get other postprocessors to act on pressure, such as volume/side max/min, etc.

Do we want to update the build instructions so it doesn't try to clone the SAM repo by default? I would expect many users not to have access to this.

For fluids applications, it is often convenient to run in non-dimensional form. However, this significantly complicates any data transfers to/from nekRS, as well as the interpretation of the various postprocessors we might use to query the nekRS solution. We should add a systematic method to inform Cardinal as to the reference scales selected for the nekRS non-dimensionalization, and then be sure to dimensionalize/non-dimensionalize all data going in or out of nekRS, respectively.

This will require obtaining the reference scales from the user (since nekRS does not save this information explicitly - it's up to the user to define correctly in the .par and .udf files), and then using those in Cardinal to transform dimensional -> non-dimensional quantities and vice versa.

Cardinal does not have very much documentation - as the number of users grows, it will be helpful to have some tutorials and a high-level description of how Cardinal works.

Right now these tolerances are hard-coded values or rely on values from OpenMC's constants. It would perhaps be convenient to allow users to specify what tolerances are applied in these checks.

To facilitate seamless coupling of Cardinal to Sockeye, we can add Sockeye as an optional submodule dependency in the same manner as was done with SAM.

For the OpenMC inputs in the Cardinal repository, we should update the OpenMC input files to use the new model generation script to take advantage of lattices. This isn't strictly necessary, but is a good idea to make sure comments like

This model was generated with ./scripts/openmc_pebble_bed_model.py -e 1.0

do actually apply to the inputs they're referring to when we change the generation script.

For some examples I'm running, if I use too large a time step or too high a thermal conductivity, the nekRS-MOOSE thermal solution diverges - the flux and temperature on the shared surface between these two applications oscillates with growing amplitude. I think this could probably be fixed by incorporating relaxation, which I believe is already supported in the MOOSE framework.

So we should make a test that demonstrates that relaxation is effective, and if not, incorporate some type of relaxation specific to Cardinal.

... copying over from Gitlab

Somewhere in the gap between CI on the gitlab repo and CI on this github repo, the boundary_and_volume/prism test is failing with about 1e-3 relative differences in temperature. We should resolve this discrepancy, possibly by re-golding if appropriate, to ensure that the CI tests pass for new PRs.

While NEKRS_HOME is read from the environment for Cardinal's makefiles, the Summit submission scripts currently just set it internally, which can be hard to keep up-to-date if the case directories change positions in the tree or you copy a script for a new case in a different location.

... copying over from Gitlab

Currently, both NekMesh and NekRSMesh use a replicated mesh approach, where all the mesh information is replicated on every parallel process. However, the external mesh-type objects should know about the parallel decomposition in the external application, so we can take advantage of this and move to a distributed mesh.

... copying over from Gitlab

In #22 , @roystgnr fixed the way that standard flags are passed from libmesh_CXX to CMAKE_CXX_COMPILER:

the autoconf macros libMesh uses will choose a $libmesh_CXX like
"mpicxx -std=c++17" (with spaces, so it can be directly included in
command lines), but cmake wants to see "mpicxx;-std=c++17"

I mentioned a few other possiblities for further enhancements:

Ideally, we should split a token like mpicxx -std=c++17 into two components such that:

mpicxx gets assigned to -DCMAKE_CXX_COMPILER
-std=c++17 is parsed and assigned to -DCXX_STANDARD
If you like, I can put together a PR that does this

And Roy had some feedback on that:

This alone probably wouldn't be safe: the autoconf C++ standard detection module also looks for +std=c++17 (HP's aCC) and/or "-h std=c++17" (Cray's crayCC), and who knows what else might be possible. Ironically, the flag that selects a C++ standard is not at all standardized.

But your later suggestions, where you pull out any standard you can parse, but then make sure to save any other flags, in a more appropriate variable, seems like a pure improvement over this PR. (assuming CMAKE_CXX_FLAGS_INIT is the right variable; I'm not familiar enough to be certain)

This issue is my reminder to follow-up and implement our consensus decision.

In order to superimpose a heat flux due to radiation on the CHT heat flux from nekRS, it would be helpful to change the current data transfers (temperature into nekRS, flux into MOOSE) to the other direction - flux from neKRS->MOOSE, and temperature from MOOSE->nekRS.

If you try to uniformly run all the tests in Cardinal with

./run_tests --parallel 64

a lot of the nekRS tests will error because there's too many processors for how many elements exist in the problem. For those tests, we can set max_parallel to a reasonable value so that the tests don't fail depending on the size of the nekRS mesh.

From the perspective of mesh generation, it may be very convenient if we could run multiple, different, nekRS applications within the same Cardinal run. For instance, suppose we have a mesh for a single fuel bundle and a mesh for the inter-wrapper flow space between a group of 3 bundles. Rather than repeat the mesh generation process just to get all these meshes in the same file, we could achieve coupling by running multiple nekRS simulations as sub-apps.

However, if the nekRS simulations need to use different meshes, I don't see a current way for our passing of --nekrs-setup on the command line to distinguish between multiple MultiApp blocks in the MOOSE input file.

It would be nice if we could instead set the value of --nekrs-setup in the nek.i input files.

The power normalization from OpenMC currently normalizes the power by the volumes of the OpenMC cells, when the power should instead be normalized by the volume of the MOOSE receiving app's cells. For instance, consider this single-pebble illustration:

• OpenMC model of a single sphere of radius 1.5. The volume of this cell is 14.137 cm^3
• MOOSE mesh of a single sphere of radius 1.5. Because the mesh is coarse, the volume of the MOOSE version of the pebble is 13.22 cm^3
• If the OpenMC power is 10 W, OpenMCProblem normalizes by 14.137 cm^3 to get a power of 0.707 W/cm^3
• But then when this power gets sent to MOOSE, the total power seen by the MOOSE mesh is instead 9.34 W (when it should be 10).

There's nothing technically wrong with OpenMCProblem, but all of the inputs in large-problems use sphere volumes of 14.137, when sphere.e actually has a volume of 13.22 cm^3. We should update all inputs with this error to be sure that we're setting the desired power in MOOSE.

One feature that might be nice to provide is the specification of tally triggers in the OpenMCCellAverageProblem input. These triggers will halt an OpenMC run once a tally's standard deviation or relative error goes below the specified value for all bins in the tally. This would be especially useful when applying mesh tallies where the number of particles/batches needed to reduce the variance of the heating tally is difficult to predict.

OpenMC uses units of centimeters, while MOOSE and nekRS are both dimension-agnostic. However, from convention, it's more convenient to solve MOOSE and nekRS in SI units (meters for length). If we add a scaling factor to OpenMCProblem to correct the

1. Mesh used for unstructured tallies
2. Cell identification for cell-based tallies
3. Power normalization by volume

then in theory we could run OpenMC in cm and MOOSE and nekRS in some other system of units.

... copying over from Gitlab

Additional to-do items

• Add capability to OpenMCProblem for both cell and mesh tallies
• Add capability to OpenMCCellAverageProblem for mesh tallies

People may want to use Cardinal's Nek wrapping for purposes other than CHT coupling - such as extracting a turbulent viscosity and projecting it onto a Pronghorn mesh for a lower-order model. We should add an input parameter that writes various Nek fields into the Exodus output file used in the wrapper. The input file could look like:

[Problems]
  type = NekRSProblem
  ...
  output = 'mu_t pressure'
[]

where output is an enumeration of various fields we want to write onto the Exodus mesh.

After running into the bug fixed by #24 the hard way, I decided to let ./run_tests --valgrind tell me whether there was anything else I ought to be looking out for. And the answer was, e.g.:

test:nek_errors/no_temp_var.no_temperature_variable: ==1465809== ERROR SUMMARY: 498 errors from 344 contexts (suppressed: 17 from 17)

And I don't have time to dig through that minefield. But although most of the errors are probably red herrings (like the 17 OpenMPI reports that were suppressed by MOOSE settings), any one of them could be another segfault-under-different-circumstances waiting to happen. Ideally someone will look through and fix or suppress all these some day, but enough of them are upstream (I don't think the explicit output I got even made it past nekrs::setup) that that's a daunting task.

You can run MOOSE in --mesh-only mode to just create the Mesh to be used for the simulation. This could be very useful as a preprocessing activity when using NekRS, since there's no way to visualize .re2 files at the moment.

When you try to run a non-Cardinal NekRS input with --mesh-only, however, you get some errors related to boundary conditions not being set up in the correctly multiphysics sense. Normally, these errors are helpful because they help identify errors in your coupling setup. But if you just want to create a mesh, they mean that you need to modify your nekRS input files, even if you really only intend to look at the mesh in Paraview.

We can move these error checks to other parts of the setup that don't get called when running in --mesh-only mode.

For OpenMCCellAverageProblem, we currently check that any cell filled with a material is fissionable before trying to add tallies to it. For cells filled with universes or lattices, we can extend this check to check that the cell contains fissile material by recursing down through all the various fills.

Once NekRS allows symmetry conditions on non-aligned planes (i.e. planes aside from the x, y, and z axes), we will need to update the cht1.md tutorial input files and remove some warnings about this simplification in the model.

For coupling with OpenMC, we will need to extract the density from nekRS. We can do this in two approaches:

1. Solve nekRS with variable density properties, and then simply re-use the same transfer routines that we use for temperatures, replacing nrs->cds->S by the density copied back from the device.

2. Run nekRS in incompressible form, and apply a correlation for rho(T) based on the computed temperature.

For most of our applications, we'll probably prefer option (2) because density changes are often fairly small. For option (2), we have two options:

2a: Apply a density correlation as an auxkernel in the OpenMC input file.
2b: Apply a correlation within nekRS itself (we could add a function in NekInterface for properties of interest).

When compiling Cardinal, nekRS can throw a fatal error (the one I most often see is Unsupported compiler!), but then the rest of the compilation keeps going like nothing bad happened. At the end, no executable is produced, but the actual error message didn't produce a fatal error as it should have.

For normalizing the boundary heat flux or a volumetric heat source, right now Cardinal requires the user to add a flux_integral and a source_integral Receiver-type postprocessors. Like the receiving heat flux and temperature variables, we can just add these automatically to simplify input file setup a bit.

nekRS now uses three CMake variables - NEKRS_VERSION, NEKRS_SUBVERSION, and GITCOMMITHASH to print a header. Within Cardinal's build system, we need to manually specify these three variables. Right now, this is done in cardinal/Makefile by writing directly to CXXFLAGS:

CXXFLAGS += -DNEKRS_VERSION=21
CXXFLAGS += -DNEKRS_SUBVERSION=1
CXXFLAGS += -DGITCOMMITHASH=\"73b4f2ee2ab1e7dc9f154da33b8f0368854d50e2\"

However, it would make more sense and be much clearer if we set these variables within NEKRS_CXXFLAGS in cardinal/config/nekrs.mk. However, when I naively just move the above to nekrs.mk (replacing CXXFLAGS by NEKRS_CXXFLAGS, I get errors related to those variables still not being found, like this:

/home/anovak/cardinal/contrib/nekRS/src/lib/nekrs.cpp:34:18: error: ‘NEKRS_VERSION’ was not declared in this scope; did you mean ‘OKL_VERSION’?
   34 |        << "v" << NEKRS_VERSION << "." << NEKRS_SUBVERSION
      |                  ^~~~~~~~~~~~~
      |                  OKL_VERSION

And, if GITCOMMITHASH is actually a mix of letters and numbers, I also get errors related to a missing operator, possibly due to some weirdness with the escaped quotes:

/home/anovak/cardinal/contrib/nekRS/src/lib/nekrs.cpp: In function ‘void printHeader()’:
<command-line>: error: unable to find numeric literal operator ‘operator""b4f2ee2ab1e7dc9f154da33b8f0368854d50e2’
/home/anovak/temp/cardinal_nekrs/cardinal/contrib/nekRS/src/lib/nekrs.cpp:35:19: note: in expansion of macro ‘GITCOMMITHASH’
   35 |        << " (" << GITCOMMITHASH << ")" << endl
      |                   ^~~~~~~~~~~~~

where b4f2ee2ab1e7dc9f154da33b8f0368854d50e2 is the part of the actual git commit hash after any leading numbers (if any).

This is weird behavior that would be nice to fix.

Currently, we ensure a conservative transfer of heat source (i.e. heat source computed by OpenMC exactly matches the power when integrated on the receiving MOOSE mesh, to conserve power) by

1. Manually setting the quadrature rule in the OpenMCProblem/OpenMCCellAverageProblem to be the same as in the receiving MOOSE application
2. Using exactly the same mesh (and variable basis/order) in the OpenMC wrapping and the receiving MOOSE application

This can be tedious and error prone to set up, when instead we could use the from_postprocessors_to_be_preserved and to_postprocessors_to_be_preserved options in MOOSE. We should update the tests and examples in Cardinal to use this simpler approach.

It would be helpful to view the velocity as one of the fields that the Nek postprocessors can operate on.

In order to utilize Picard iteration (as opposed to the loose coupling that we currently have), nekRS must be able to reset it's internal state back to the start of a time step, multiple times. Having Picard iteration would be advantageous when using subcycling, so that MOOSE might only run every 100 time steps, as it might be more accurate than having loose coupling with an interpolated heat flux (from MOOSE -> nekRS) if the true solution is changing quickly such that the O(dt) error arising from loose coupling is not small.

... copying over from Gitlab

The current GitLab pipelines tests don't check whether the simulation results have changed, but only if the run finishes without an error. To protect us from accidentally making changes that cause the simulation to change, we should be checking that the simulation results are independent of changes we make.

So, I think we should change the spherical_heat_conduction tests to Exodiff-type CIVET tests.

... copying over from Gitlab

The mesh tests in Cardinal currently check the node positions based on the element node IDs. These tests were golded with a particular number of processes, and if you run with a different amount, the partitioning information from nekRS is different (so elements are added in a different order). Therefore, these tests will only work for the processor count for which they were golded.

In line with #45, our tests should be agnostic of the number of ranks.

After we finish making the new problems repo, we should update the download script to point to the new location. It would also be convenient to have this script download the auxiliary large files (i.e. meshes stored in Box somewhere).

In order to force data transfers into NekApp only when the master application has new information, we can use the minimize_transfers_in = true option. This requires the manual addition of a Receiver postprocessor.

For mesh deformation problems, we will always need to use this minimize_transfers_in feature if NekApp is run with a smaller time step than the master application because we add displacement to the NekRS mesh. In other words, if a tensor mechanics solve computes a displacement of 1 meter and runs with a time step of 1 second, and we run a NekApp with a time step of 0.1 seconds, then if we didn't restrict the incoming mesh transfer we'd apply the 1 meter displacement 10 times (1 second / 0.1 seconds), giving the incorrect displacement.

Because we're always going to be using this minimize_transfers_in feature for mesh movement problems, this feature is no longer a "fringe" concept, so the setup should therefore be simplified by automation.

After NekRS adds the ability to take variable time step sizes, the NekTimeStepper will need to be updated to also utilize this feature. We will also need to update the documentation, specifically a warning printed in the cht1.md tutorial.

We don't currently have any CIVET tests for OpenMC-MOOSE coupling, and for the distributed transfers it would be useful to have baseline models that we expect to match that can be quickly run with ./run_tests and exodiffed.

It would be good to have the capacity to pass 0D information between codes at boundaries. For example: Imagine a flow loop being simulated in SAM, with only a portion of a pipe being simulated in NekRS. We would want:

• to pass avg. temperature, avg. pressure and total mass flow rate (or avg. velocity, area and density) from NekRS to SAM on both pipe surfaces (inlet and outlet); Similar information back from SAM to NekRS;

• OR to pass forces like the average pressure gradient between NekRS and SAM. temperature, velocity and pressure from SAM to NekRS.

For a general feature I think we would to have an API that allows a number of N properties to be exchanged at the boundaries. Postprocessors can then be added in SAM or other system analysis codes to pass to NekRS.

To keep the various coupling types distinct, I propose adding a new problem, like NekScalarProblem or Nek0DProblem that would be used to couple Nek to systems-level feedback to differentiation from NekRSProblem, that couples for 2-D boundary and 3-D volume feedback.

Basic checklist:

• Build input models for the 3 coupling options (NekRS b/w two 1-D domains, NekRS at exit of 1-D domain, NekRS at inlet of 1-D domain)
• Add robust error checking
• Add online documentation
• Add online tutorial
• Ensure that the correct transformations exist for non-dimensional Nek models
• Ensure compatibility with distributed mesh implementation

Possible more advanced features:

• Allow a multi-inlet, multi-outlet NekRS domain to be coupled, by allowing different inlet and outlet BCs to be set based on boundary

Allow nekRS mesh to deform based on displacements provided by MOOSE. This will involve adding new receiver variables for displacement, computing a mesh velocity based on a finite difference approximation of displacement, and then interpolating this mesh velocity onto the nekRS GLL points.

This issue can be closed when a PR that addresses the following has been submitted:

Recent additions for #92 use 3 mesh deformation functions for x,y,z coordinates that need to be condensed into one using enums. will be tackled in a separate issue

• A call to nekrs::outfld needs to be added to the NekRSProblem constructor to write the initial mesh to a field file to match the exodus output.

We're planning on having a discussion around options for a more flexible tally specification in the OpenMCProblem. It's pretty rigid and fairly opaque right now, requiring users to have all tally cells on the same level of the geometry and to know what level that is. As mentioned in #33, we also expect that the tally/tallies in the OpenMC problem enclose all of the fissionable material in the problem and are non-overlapping spatially.

Until we can discuss this further, @aprilnovak made the good suggestion that we perform a check for unique bins in the cell tally to protect against a situation in which the incorrect pebble_cell_level has been specified in the model, resulting in a cell tally either too high or too low in the geometry hierarchy.

If we want to use a lower-resolution T/H model (such as a porous media or subchannel application) to get a good initial condition for nekRS, we need to be able to set an initial temperature, velocity, and pressure field in nekRS. From the perspective of just getting an initial condition, one possible workflow could be:

1. Run Pronghorn as a master application to estimate both the fluid and solid temperatures.
2. Run BISON-nekRS as a sub-application after Pronghorn finishes a pseudo-steady calculation.

This requires the ability to transfer initial conditions into nekRS, which is different from setting either a heat source or heat flux boundary condition because we would directly be modifying nrs->S, nrs->o_S, etc.

Using METHOD=oprof (to build with performance monitoring instrumentation), devel (to build with assertions and debugging symbols enabled), and/or dbg (to use more expensive libMesh and libstdc++ assertions, and to disable optimization to make stepping through a debugger easier) are extremely useful, but we can't get a usable Cardinal build with them right now.

The culprit seems to be contrib/openmc/cmake/Modules/FindLIBMESH.cmake - it's looking for libmesh-opt.pc, specifically, with no way to override that at the moment, but that ought to be easy enough to fix and send upstream. If I need more debugging support (which I probably will, shortly) then I'll do so.

Building with more than one method at a time would be harder - use multiple out-of-source builds for dependency libraries, maybe - so that's not going to be a near-term goal.

The nekRS interface right now currently assumes that nekRS is solving a fluid-only problem. But, it's possible that someone might want to have nekRS solve the heat equation in one part of the geometry, while MOOSE handles the solid phase elsewhere. While theoretically the NekInterface should work regardless of whether nekRS's problem includes CHT or not, there might be some small hiccups in the setup. So, we should add a test covering this option to make sure this capability is actually supported as-is (and otherwise support it if not).

This will help streamline the coupling of nekRS and SAM via Cardinal.

tag @vcleite2019 @rui-hu

Many of the nekRS .par files used in the tests have some old Nek5000 syntax that doesn't do anything in nekRS. We should remove it to avoid confusing people and to be clear what in the input files actually affects the nekRS solution.

It is expected that the cell/mesh tally in the OpenMCProblem encompasses all fissionable material in the problem in a non-overlapping manner. To check that this is true, we can add a global kappa-fission tally to the OpenMC problem and check that the sum of the tallies used in transfers equals this global tally value. This provides an extra check for users to ensure that the tallies have been setup correctly.

The current wrapping of OpenMC, in OpenMCProblem, is very dependent on the system being a pebble bed reactor. The temperature that gets sent into OpenMC is computed based on a provided point cloud in pebble_centers by doing a nearest-point average. This is great for feedback to sphere domains, but doesn't translate well to many other geometries. Consider these two cases:

• Suppose we wanted to transfer the temperature for the pellet and cladding in a typical LWR pincell geometry to OpenMC. While we can provide a point cloud that will pick up the pellet, there's no straightforward choice for what "point" to associate with the nearest-point average in the cladding. Even if we could get the centroids of all the OpenMC cells that are in the cladding, there will be some elements towards the edges of the fuel pellet that will "pick up" the cladding point cloud. So, the point cloud approach won't allow anything with an annular geometry (like cladding).
• More generally any geometry where the maximum distance of the cell surface from its centroid is greater than the distance between that cell's centroid and a neighboring cell, is not going to work with the nearest point idea. Consider this 3-element setup. The centroid of each element is marked with a *.

|-------|-------|
|  el 1 |  el 2 |
|       |   *   |
|   *   |       |
|       |-------|
|       | el 3  |
|       |   *   |
|-------|-------|

There will be some points in the lower right corner of element 1 that will "pick up" the centroid of element 3, and contribute to the average there. This could end up having elements with completely different materials contributing temperatures to different areas of the problem that they should not be.

In addition, OpenMCProblem has a few limitations that should be addressed:

• OpenMCProblem does not have any density feedback
• OpenMCProblem does not have any capability for fluid temperature feedback (because of the general limitations of the nearest point transfer)
• OpenMCProblem requires the entire fissile part of the OpenMC geometry to be coupled to MOOSE; while this will cover most people's use cases, we may at some point only want to do multiphysics feedback for a single bundle with OpenMC (and the rest of the core with Griffin). This manifests in how we normalize the tallies and the errors we throw if the OpenMC mapped domain does not match the global tally.
• OpenMCProblem required you to automatically provide the volumes of the MOOSE elements that each OpenMC cell maps to. This can be tedious to set up if OpenMC isn't coupling to 300k identical pebbles.

We should add a new OpenMCCellAverageProblem that address all of the above limitations by allowing general density and temperature feedback to MOOSE. Volume normalization will be performed by automatically computing volumes by querying the MOOSE mesh elements.

A great way to verify this implementation will be a comparison against Enrico's cell averaging approach, since the idea is quite similar.

Remaining to-do items:

• Implement distributed mesh version
• More extensive testing in parallel
• Add mesh-based tallies

When running OpenMC, it would be helpful to have a number of postprocessors for both solve diagnostics as well as to control termination of the solve by MOOSE.

• When running OpenMC, it would be helpful to be able to view the maximum tally relative error as the solve progresses. Of course we can use the OpenMC python API to view quantities like this, but it would also be helpful during the course of the simulation (say, if the power distribution is changing once feedback is turned on).
• Convergence of multiphysics calculations is often based on k-effective. We should make a postprocessor that just reads the k-value from OpenMC.

Right now the default value of the pebble_cell_level for the OpenMCProblem is 0. Because it is not a required parameter (it is not needed if a mesh tally is being used), the cell tally will happily proceed with this value and potentially incorrectly tally heating values. We've had this happen a few times now.

If a cell tally is being used, we should have a check that the pebble_cell_level is explicitly set in the input file by the user, indicating that they understand at what level of the geometry the tally should occur.