Line data Source code
1 : /********************************************************************/
2 : /* SOFTWARE COPYRIGHT NOTIFICATION */
3 : /* Cardinal */
4 : /* */
5 : /* (c) 2021 UChicago Argonne, LLC */
6 : /* ALL RIGHTS RESERVED */
7 : /* */
8 : /* Prepared by UChicago Argonne, LLC */
9 : /* Under Contract No. DE-AC02-06CH11357 */
10 : /* With the U. S. Department of Energy */
11 : /* */
12 : /* Prepared by Battelle Energy Alliance, LLC */
13 : /* Under Contract No. DE-AC07-05ID14517 */
14 : /* With the U. S. Department of Energy */
15 : /* */
16 : /* See LICENSE for full restrictions */
17 : /********************************************************************/
18 :
19 : #pragma once
20 :
21 : #include "OpenMCProblemBase.h"
22 : #include "SymmetryPointGenerator.h"
23 :
24 : /// Tally/filter includes.
25 : #include "TallyBase.h"
26 : #include "FilterBase.h"
27 :
28 : #ifdef ENABLE_DAGMC
29 : #include "MoabSkinner.h"
30 : #include "DagMC.hpp"
31 : #endif
32 :
33 : /// Forward declarations to avoid cyclic dependencies.
34 : class OpenMCCellMaterialFill;
35 : class OpenMCVolumeCalculation;
36 :
37 : /**
38 : * Mapping of OpenMC to a collection of MOOSE elements, with temperature and/or
39 : * density feedback. The mappind is established authomatically by looping over
40 : * all the MOOSE elements and finding the OpenMC cell at each element's centroid.
41 : *
42 : * All feedback into OpenMC is performed via element averages. The
43 : * 'temperature_blocks' parameter is used to indicate which MOOSE blocks should
44 : * provide temperature feedback, while the 'density_blocks' parameter is used to
45 : * indicate which MOOSE blocks should provide density feedback. Tallies are
46 : * automatically added to OpenMC using either cell or mesh tallies.
47 : *
48 : * Each OpenMC cell shall not have ambiguous data transfers. That is, a cell
49 : * should map to a set of elements that are ALL/NONE providing temperature
50 : * feedback, ALL/NONE providing density feedback, ALL/NONE providing cell
51 : * tallies, and ALL/NONE being uncoupled altogether.
52 : *
53 : * TODO: If this is too restrictive in the future, we could implement some type
54 : * of weighted averaging process. Also, if a cell maps to a phase and an
55 : * unmapped region, perhaps we want to allow that.
56 : *
57 : * Other considerations you should be aware of:
58 : * - The density being transferred into OpenMC from MOOSE is in units of kg/m3; this
59 : * is the unit employed by the MOOSE fluid properties module.
60 : * - The temperature being transferred into OpenMC from MOOSE is in units of K; this
61 : * is the unit employed by the MOOSE fluid and solid properties modules.
62 : * - If your geometry is highly fine-scale (such as TRISOs), you might be able to get a
63 : * huge speedup in your runtime if you (i) build your OpenMC model by repeating the same
64 : * TRISO universe in each of your repeatable-units (e.g. pebbles, compacts, plates)
65 : * AND (ii) leverage the 'identical_cell_fills' option.
66 : */
67 : class OpenMCCellAverageProblem : public OpenMCProblemBase
68 : {
69 : public:
70 : OpenMCCellAverageProblem(const InputParameters & params);
71 : static InputParameters validParams();
72 :
73 : virtual void initialSetup() override;
74 : virtual void externalSolve() override;
75 : virtual void syncSolutions(ExternalProblem::Direction direction) override;
76 2733 : virtual bool converged(unsigned int) override { return true; }
77 :
78 : /**
79 : * Read a 2d vector of subdomain names, and check that there are no duplications
80 : * and that all provided values exist on the mesh.
81 : * @param[in] name string name for the 2d vector parameter
82 : * @param[out] names subdomain names
83 : * @param[out] flattened_ids flattened 1d vector of subdomain IDs
84 : */
85 : void read2DBlockParameters(const std::string name,
86 : std::vector<std::vector<SubdomainName>> & names,
87 : std::vector<SubdomainID> & flattened_ids);
88 :
89 : /**
90 : * Check that the specified blocks are in the mesh
91 : * @param[in] name name for throwing an error
92 : * @param[in] ids block IDs to check
93 : * @param[in] names block subdomain names for throwing an error
94 : */
95 : void checkBlocksInMesh(const std::string name,
96 : const std::vector<SubdomainID> & ids,
97 : const std::vector<SubdomainName> & names) const;
98 :
99 : /// Initialize the mapping of OpenMC to the MooseMesh and perform additional setup actions
100 : void setupProblem();
101 :
102 : /**
103 : * Add the tally variable(s) (to receive OpenMC tally values), temperature variable(s)
104 : * (to write into OpenMC cells), and density variable(s) (to write into OpenMC materials)
105 : */
106 : virtual void addExternalVariables() override;
107 :
108 : /**
109 : * Get the cell volume from a stochastic calculation
110 : * @param[in] cell_info cell index, instance pair
111 : * @return stochastically-computed OpenMC cell volume
112 : */
113 : virtual Real cellVolume(const cellInfo & cell_info) const;
114 :
115 : /**
116 : * Reference to stochastic volume calculation
117 : * @return reference to stochastic volume calculation
118 : */
119 609314 : virtual const OpenMCVolumeCalculation * volumeCalculation() const { return _volume_calc; }
120 :
121 : /**
122 : * Get the mapping of cells to MOOSE elements
123 : * @return mapping of cells to MOOSE elements
124 : */
125 12106 : virtual const std::map<cellInfo, std::vector<unsigned int>> & cellToElem() const
126 : {
127 12106 : return _cell_to_elem;
128 : }
129 :
130 : /**
131 : * Get the MOOSE subdomains associated with an OpenMC cell
132 : * @param info the cell info
133 : * @return MOOSE subdomains associated with an OpenMC cell
134 : */
135 14482 : virtual std::unordered_set<SubdomainID> getCellToElementSub(const cellInfo & info)
136 : {
137 14482 : return _cell_to_elem_subdomain.at(info);
138 : }
139 :
140 : /**
141 : * Whether transformations are applied to the [Mesh] points when mapping to OpenMC
142 : * @return whether transformations are applied
143 : */
144 2312467 : virtual bool hasPointTransformations() const { return _symmetry != nullptr; }
145 :
146 : /**
147 : * Get all the scores added to the tally
148 : * @return scores
149 : */
150 108 : virtual const std::vector<std::string> & getTallyScores() const { return _all_tally_scores; }
151 :
152 : /**
153 : * Get the number of tallies scoring a particular score.
154 : * @param[in] score the score to check
155 : * @return the number of tallies accumulating 'score'
156 : */
157 528 : unsigned int getNumScoringTallies(const std::string & score) const
158 : {
159 528 : return _score_count.count(score) > 0 ? _score_count.at(score) : 0;
160 : }
161 :
162 : /**
163 : * Check to see if this problem contains a specific tally score.
164 : * @param[in] score the tally score
165 : * @return whether this problem contains the tally score in a tally object
166 : */
167 356 : bool hasScore(const std::string & score)
168 : {
169 356 : return std::find(_all_tally_scores.begin(), _all_tally_scores.end(), score) !=
170 356 : _all_tally_scores.end();
171 : }
172 :
173 : /**
174 : * Get a tally object by its name.
175 : * @param[in] name the name of the TallyBase
176 : * @return a Cardinal wrapped tally
177 : */
178 : const TallyBase * getTally(const std::string & name);
179 :
180 : /**
181 : * Get the variable(s) associated with an OpenMC tally score.
182 : * @param[in] score the OpenMC score
183 : * @param[in] tally_name the name of the tally object to fetch score variables from
184 : * @param[in] tid the thread ID associated with the current MOOSE object
185 : * @param[in] output the output variable (relative error, standard deviation, etc.) to fetch
186 : * @param[in] skip_func_exp whether functional expansion filter bins should be skipped or not when
187 : * fetching variable values
188 : * @return a vector of variable values associated with score
189 : */
190 : std::vector<const MooseVariableFE<Real> *> getTallyScoreVariables(const std::string & score,
191 : const std::string & tally_name,
192 : THREAD_ID tid,
193 : const std::string & output = "",
194 : bool skip_func_exp = false);
195 :
196 : /**
197 : * Get the variable value(s) associated with an OpenMC tally score.
198 : * @param[in] score the OpenMC score
199 : * @param[in] tid the thread ID associated with the current MOOSE object
200 : * @param[in] tally_name the name of the tally object to fetch score variable values from
201 : * @param[in] output the output variable (relative error, standard deviation, etc.) to fetch
202 : * @param[in] skip_func_exp whether functional expansion filter bins should be skipped or not when
203 : * fetching variable values
204 : * @return a vector of variable values associated with score
205 : */
206 : std::vector<const VariableValue *> getTallyScoreVariableValues(const std::string & score,
207 : const std::string & tally_name,
208 : THREAD_ID tid,
209 : const std::string & output = "",
210 : bool skip_func_exp = false);
211 :
212 : /**
213 : * Get the variable value(s) associated with an OpenMC tally score.
214 : * @param[in] score the OpenMC score
215 : * @param[in] tally_name the name of the tally object to fetch score variable neighbor values from
216 : * @param[in] tid the thread ID associated with the current MOOSE object
217 : * @param[in] output the output variable (relative error, standard deviation, etc.) to fetch
218 : * @param[in] skip_func_exp whether functional expansion filter bins should be skipped or not when
219 : * fetching variable values
220 : * @return a vector of variable values associated with score
221 : */
222 : std::vector<const VariableValue *>
223 : getTallyScoreNeighborVariableValues(const std::string & score,
224 : const std::string & tally_name,
225 : THREAD_ID tid,
226 : const std::string & output = "",
227 : bool skip_func_exp = false);
228 :
229 : /**
230 : * Whether a tally contains a specified output or not.
231 : * @param[in] score the tally score to check
232 : * @param[in] output the additional output (unrelaxed standard deviation, relative error, or
233 : * tally)
234 : * @return whether an added tally has the output or not
235 : */
236 : bool hasOutput(const std::string & score, const std::string & output) const;
237 :
238 : /**
239 : * Apply transformations to point
240 : * @param[in] pt point
241 : * @return transformed point
242 : */
243 2263313 : virtual Point transformPoint(const Point & pt) const
244 : {
245 2263313 : return this->hasPointTransformations() ? _symmetry->transformPoint(pt) : pt;
246 : }
247 :
248 : /**
249 : * This class uses elem->volume() in order to normalize the tally values. However,
250 : * elem->volume() is expensive, so whenever MOOSE does integration, they set
251 : * _current_elem_volume to the volume as set by the sum of the quadrature weights.
252 : * The quadrature rule that MOOSE provides when you only have CONSTANT MONOMIALS is
253 : * insufficient for exactly integrating the element Jacobian mapping type (which
254 : * is FIRST LAGRANGE for a first order element), so you get an error relative to
255 : * the libmesh volume computation.
256 : *
257 : * So, we need to make sure that a minimum order quadrature rule is used
258 : * so that the total tally as computed by an
259 : * ElementIntegralVariablePostprocessor actually matches the specified total
260 : * (for low quadrature orders, there can be an error up to about 5% or so in total
261 : * power). This override simply forces the volume quadrature order to be 2 or higher
262 : * when using Gauss (default), monomial, or Gauss-Lobatto quadrature.
263 : *
264 : * For other quadrature rules, the approximations made in elem->volume() are never
265 : * going to match the volume integrations in MOOSE (unless the quadrature order is
266 : * very very high). For these orders, we print an error message informing the user
267 : * that they should switch to a different order.
268 : */
269 : virtual void createQRules(QuadratureType type,
270 : Order order,
271 : Order volume_order,
272 : Order face_order,
273 : SubdomainID block,
274 : bool allow_negative_weights = true) override;
275 :
276 : /**
277 : * Type definition for cells contained within a parent cell; the first value
278 : * is the cell index, while the second is the set of cell instances
279 : */
280 : typedef std::unordered_map<int32_t, std::vector<int32_t>> containedCells;
281 :
282 : /**
283 : * Get the cell index from the element ID; will return UNMAPPED for unmapped elements
284 : * @param[in] elem_id element ID
285 : * @return cell index
286 : */
287 : int32_t elemToCellIndex(const int & elem_id) const { return elemToCellInfo(elem_id).first; }
288 :
289 : /**
290 : * Get the cell ID from the element ID. Note that this function requires that the elem_id
291 : * maps to an OpenMC cell, or else an error will be raised from OpenMC in cellID.
292 : * @param[in] elem_id element ID
293 : * @return cell ID
294 : */
295 25909658 : int32_t elemToCellID(const int & elem_id) const { return cellID(elemToCellIndex(elem_id)); }
296 :
297 : /**
298 : * Get the cell instance from the element ID; will return UNMAPPED for unmapped elements
299 : * @param[in] elem_id element ID
300 : * @return cell instance
301 : */
302 : int32_t elemToCellInstance(const int & elem_id) const { return elemToCellInfo(elem_id).second; }
303 :
304 : /**
305 : * Get the cell index, instance pair from element ID; if the element doesn't map to an OpenMC
306 : * cell, the index and instance are both set to UNMAPPED
307 : * @param[in] elem_id element ID
308 : * @return cell index, instance pair
309 : */
310 37571068 : cellInfo elemToCellInfo(const int & elem_id) const { return _elem_to_cell[elem_id]; }
311 :
312 : /**
313 : * Get the fields coupled for each cell; because we require that each cell maps to a consistent
314 : * set, we simply look up the coupled fields of the first element that this cell maps to. Note
315 : * that this function requires a valid instance, index pair for cellInfo - you cannot pass in an
316 : * unmapped cell, i.e. (UNMAPPED, UNMAPPED)
317 : * @param[in] cell_info cell index, instance pair
318 : * @return coupling fields
319 : */
320 : coupling::CouplingFields cellFeedback(const cellInfo & cell_info) const;
321 :
322 : /**
323 : * Whether a cell has density feedback
324 : * @param[in] cell_info cell index,instance pair
325 : * @return if cell has density feedback
326 : */
327 5300198 : bool hasDensityFeedback(const cellInfo & cell_info) const
328 : {
329 : std::vector<coupling::CouplingFields> phase = {coupling::density,
330 5300198 : coupling::density_and_temperature};
331 5300198 : return std::find(phase.begin(), phase.end(), cellFeedback(cell_info)) != phase.end();
332 5300198 : }
333 :
334 : /**
335 : * Whether a cell has temperature feedback
336 : * @param[in] cell_info cell index,instance pair
337 : * @return if cell has temperature feedback
338 : */
339 4096360 : bool hasTemperatureFeedback(const cellInfo & cell_info) const
340 : {
341 : std::vector<coupling::CouplingFields> phase = {coupling::temperature,
342 4096360 : coupling::density_and_temperature};
343 4096360 : return std::find(phase.begin(), phase.end(), cellFeedback(cell_info)) != phase.end();
344 4096360 : }
345 :
346 : /**
347 : * Checks if the [Problem/Filters] block contains a specific filter.
348 : * @param[in] filter_name the MOOSE object name of the filter
349 : * @return whether the problem contains the specified filter
350 : */
351 : bool hasFilter(const std::string & filter_name) const { return _filters.count(filter_name) > 0; }
352 :
353 : /**
354 : * Get a filter added by the [Problem/Filters] block by it's MOOSE object name.
355 : * @param[in] filter_name the MOOSE object name of the filter
356 : * @return the filter object
357 : */
358 : std::shared_ptr<FilterBase> & getFilter(const std::string & filter_name)
359 : {
360 1054 : return _filters.at(filter_name);
361 : }
362 :
363 : /**
364 : * Get the local tally
365 : * @return local tally
366 : */
367 : const std::vector<std::shared_ptr<TallyBase>> & getLocalTallies() const { return _local_tallies; }
368 :
369 : /**
370 : * Get the temperature of a cell; for cells not filled with materials, this will return
371 : * the temperature of the first material-type cell
372 : * @param[in] cell_info cell index, instance pair
373 : */
374 : double cellTemperature(const cellInfo & cell_info) const;
375 :
376 : /**
377 : * Get the density of a cell; for cells not filled with materials, this will return
378 : * the density of the first material-type cell
379 : * @param[in] cell_info cell index, instance pair
380 : * @param[in] elem element to fetch multigroup reference densities
381 : */
382 : double cellDensity(const cellInfo & cell_info, const Elem * elem) const;
383 :
384 : /**
385 : * Get the volume that each OpenMC cell mapped to
386 : * @param[in] cell_info cell index, instance pair
387 : */
388 : double cellMappedVolume(const cellInfo & cell_info) const;
389 :
390 : /// Reconstruct the DAGMC geometry after skinning
391 : void reloadDAGMC();
392 :
393 : /**
394 : * Add a Filter object using the filter system.
395 : * @param[in] type the new tally type
396 : * @param[in] name the name of the new tally
397 : * @param[in] moose_object_pars the input parameters of the new tally
398 : */
399 : void addFilter(const std::string & type,
400 : const std::string & name,
401 : InputParameters & moose_object_pars);
402 :
403 : /**
404 : * Add a Tally object using the tally system.
405 : * @param[in] type the new tally type
406 : * @param[in] name the name of the new tally
407 : * @param[in] moose_object_pars the input parameters of the new tally
408 : * @return a shared pointer to the TallyBase
409 : */
410 : std::shared_ptr<TallyBase>
411 : addTally(const std::string & type, const std::string & name, InputParameters & moose_object_pars);
412 :
413 : /**
414 : * Multiplier on the normalized tally results; for fixed source runs,
415 : * we multiply the tally (which has units of eV/source)
416 : * by the source strength and the eV to joule conversion, while for k-eigenvalue runs, we
417 : * multiply the normalized tally (which is unitless and has an integral
418 : * value of 1.0) by the power.
419 : * @param[in] score_name name of the score
420 : * @param[in] local_mean_tally the mean tally associated with score_name
421 : */
422 : Real tallyMultiplier(const std::string & score_name, const Real & local_mean_tally) const;
423 :
424 : /**
425 : * Get the reference density of an element when running in multi-group mode.
426 : * @param[in] elem the element
427 : * @return the reference density (kg/m3) or unity (if not running in multi-group mode)
428 : */
429 : const Real getReferenceDensity(const Elem * elem) const;
430 :
431 : /**
432 : * Check whether a vector extracted with getParam is empty
433 : * @param[in] vector vector
434 : * @param[in] name name to use for printing error if empty
435 : */
436 : template <typename T>
437 4766 : void checkEmptyVector(const std::vector<T> & vector, const std::string & name) const
438 : {
439 4660 : if (vector.empty())
440 8 : mooseError(name + " cannot be empty!");
441 4758 : }
442 :
443 4058 : int fixedPointIteration() const { return _fixed_point_iteration; }
444 :
445 : /**
446 : * Checks if the problem uses adaptivity or not.
447 : * @return if the problem uses adaptivity.
448 : */
449 2838 : bool hasAdaptivity() const { return _has_adaptivity; }
450 :
451 : /**
452 : * Checks if the problem is using a MoabSkinner or not.
453 : * @return if the problem uses skinning
454 : */
455 157 : bool hasSkinner() const { return _using_skinner; }
456 :
457 : /// Constant flag to indicate that a cell/element was unmapped
458 : static constexpr int32_t UNMAPPED{-1};
459 :
460 : /// Spatial dimension of the Monte Carlo problem
461 : static constexpr int DIMENSION{3};
462 :
463 : /// Get a modifyable non-const reference to the Moose mesh
464 : virtual MooseMesh & getMooseMesh();
465 :
466 : /// Get a modifyable const reference to the Moose mesh
467 : virtual const MooseMesh & getMooseMesh() const;
468 :
469 : /**
470 : * Whether a moving mesh is used
471 : * @return whether the [Mesh] is moving
472 : */
473 : const bool & useDisplaced() const { return _use_displaced; }
474 :
475 : /**
476 : * Get the number of material-fill cells contained within the given cell
477 : * @param[in] cell_info cell index, instance pair
478 : * @return number of contained cells filled by a material
479 : */
480 : int numContainedMaterialCells(const cellInfo & cell_info) const;
481 :
482 : /**
483 : * Get the first material cell contained in the given cell
484 : * @param[in] cell_info cell index, instance pair
485 : * @return material cell index, instance pair
486 : */
487 : cellInfo firstContainedMaterialCell(const cellInfo & cell_info) const;
488 :
489 : protected:
490 : /**
491 : * A function to re-initialize coupling and apply feedback to the OpenMC problem.
492 : * Applied before OpenMC is executed in either: i) normal Picard iterations, or
493 : * ii) criticality searches. This function performs these operations in the
494 : * following order:
495 : * 1. Updates the OpenMC geometry using the skinner;
496 : * 2. Resets tallies and the cell->element maps to take into account
497 : * mesh/geometry changes;
498 : * 3. Updates the nuclide composition of OpenMC materials;
499 : * 4. Sets cell temperatures and densities;
500 : * 5. Exports OpenMC properties;
501 : * 6. Reinitializes multi-group cross sections to take into account changing
502 : * temperatures and densities.
503 : */
504 : void reinitCouplingAndApplyFeedback();
505 :
506 : /**
507 : * Implement critSearchStep() to re-generate the cell-to-element (and dual)
508 : * mapping. This sends new temperatures and densities to OpenMC from the
509 : * re-mapped elements to ensure the state remains critical under changes to the
510 : * model with feedbacks.
511 : */
512 : virtual void critSearchStep() override;
513 :
514 : /**
515 : * Get the cell level in OpenMC to use for coupling
516 : * @param[in] c point
517 : * @return cell level
518 : */
519 : unsigned int getCellLevel(const Point & c) const;
520 :
521 : /**
522 : * Read the names of the MOOSE variables used for sending feedback into OpenMC
523 : * @param[in] param feedback term to read
524 : * @param[in] default_name default name to use for MOOSE variables holding this field
525 : * @param[out] vars_to_specified_blocks map from MOOSE variable names to the blocks on which they
526 : * are defined
527 : * @param[out] specified_blocks user-specified blocks for feedback
528 : */
529 : void
530 : readBlockVariables(const std::string & param,
531 : const std::string & default_name,
532 : std::map<std::string, std::vector<SubdomainName>> & vars_to_specified_blocks,
533 : std::vector<SubdomainID> & specified_blocks);
534 :
535 : /**
536 : * Whether this cell has an identical fill
537 : * @param[in] cell_info cell index, instance pair
538 : * @return whether this cell has an identical fill
539 : */
540 : bool cellHasIdenticalFill(const cellInfo & cell_info) const;
541 :
542 : /**
543 : * When using the 'identical_cell_fills' feature, this is used to determine the
544 : * new instance of a cell contained in 'cell_info'
545 : * @param[in] cell_info cell index, instance pair of the containing cell
546 : * @param[in] cc_idx index in the openmc::model::cells array for a cell contained in 'cell_info'
547 : * @param[in] cc_instance_idx_to_shift the index in the containedCells instance array for the
548 : * instance we want to shift
549 : * @return a shifted instance
550 : */
551 : int containedCellInstanceShift(const cellInfo & cell_info,
552 : int32_t cc_idx,
553 : int32_t cc_instance_idx_to_shift) const;
554 :
555 : /**
556 : * When using the 'identical_cell_fills' feature, this is use to get a list of
557 : * contained cells and their instances for a given 'cell_info'.
558 : * @param[in] cell_info the containing cell
559 : * @return the cells contained in 'cell_info'. Instances are not shifted.
560 : */
561 : const containedCells & unshiftedContainedCells(const cellInfo & cell_info) const;
562 :
563 : /**
564 : * Whether this cell overlaps with ANY value in the given subdomain set
565 : * @param[in] cell_info cell index, instance pair
566 : * @param[in] id subdomain IDs
567 : * @return whether the cell overlaps with the subdomain
568 : */
569 : bool cellMapsToSubdomain(const cellInfo & cell_info,
570 : const std::unordered_set<SubdomainID> & id) const;
571 :
572 : /**
573 : * Delete the OpenMC DAGMC geometry and re-generate the CSG geometry data structures in-place.
574 : */
575 : void updateOpenMCGeometry();
576 :
577 : /**
578 : * Get a list of each material in the problem, sorted by subdomain. This function also checks
579 : * that there is just one OpenMC material in each subdomain, necessary for the DAGMC skinning.
580 : * @return material in each subdomain
581 : */
582 : std::vector<std::string> getMaterialInEachSubdomain() const;
583 :
584 : /**
585 : * Apply transformations and scale point from MOOSE into the OpenMC domain
586 : * @param[in] pt point
587 : * @return transformed point
588 : */
589 : Point transformPointToOpenMC(const Point & pt) const;
590 :
591 : /**
592 : * For geometries with fine-scale details (e.g. TRISO), Cardinal's default settings can
593 : * take a very long time to initialize the problem (but we can't change those defaults
594 : * because they are not 100% applicable all the time). So, we print out a message to
595 : * the user to point them in the right direction if their initialization is taking a
596 : * long time.
597 : * @param[in] start time to use for evaluating whether we've exceeded our limit for printing the
598 : * message
599 : */
600 : void
601 : printTrisoHelp(const std::chrono::time_point<std::chrono::high_resolution_clock> & start) const;
602 :
603 : /**
604 : * Print to the console the names of the auxvariables used for I/O with OpenMC.
605 : * We only print these tables once, upon initialization, because this data does
606 : * not change if re-initializing the spatial mapping for moving-mesh problems,
607 : * adaptive refinement, skinning, etc.
608 : */
609 : void printAuxVariableIO();
610 :
611 : /**
612 : * Get all the material indices within the set of cells
613 : * @param[in] contained_cells set of cells
614 : * @return contained materials
615 : */
616 : std::vector<int32_t> materialsInCells(const containedCells & contained_cells) const;
617 :
618 : /// Loop over the mapped cells, and build a map between subdomains to OpenMC materials
619 : void subdomainsToMaterials();
620 :
621 : /**
622 : * Get a set of all subdomains that have at least 1 element coupled to an OpenMC cell
623 : * @return subdomains with at least 1 element coupled to OpenMC
624 : */
625 : std::set<SubdomainID> coupledSubdomains() const;
626 :
627 : /**
628 : * Gather a vector of values to be summed for each cell
629 : * @param[in] local local values to be summed for the cells
630 : * @param[out] global global mapping of the summed values to the cells
631 : */
632 : template <typename T>
633 : void gatherCellSum(std::vector<T> & local, std::map<cellInfo, T> & global) const;
634 :
635 : /**
636 : * Gather a vector of values to be pushed back to for each cell
637 : * @param[in] local local values to be pushed back for the cells
638 : * @param[in] n_local number of local values contributed to each cell
639 : * @param[out] global global mapping of the pushed back values to the cells
640 : */
641 : template <typename T>
642 : void gatherCellVector(std::vector<T> & local,
643 : std::vector<unsigned int> & n_local,
644 : std::map<cellInfo, std::vector<T>> & global);
645 :
646 : /**
647 : * Get the feedback which this element provides to OpenMC
648 : * @param[in] elem
649 : * @return coupling phase
650 : */
651 : coupling::CouplingFields elemFeedback(const Elem * elem) const;
652 :
653 : /**
654 : * Read the parameters needed for triggers
655 : * @param[in] params input parameters
656 : */
657 : void getTallyTriggerParameters(const InputParameters & params);
658 :
659 : /**
660 : * Read the block parameters based on user settings
661 : * @param[in] name name of input parameter representing a vector of subdomain names
662 : * @param[in] blocks list of block ids to write
663 : */
664 : void readBlockParameters(const std::string name, std::unordered_set<SubdomainID> & blocks);
665 :
666 : /**
667 : * Cache the material cells contained within each coupling cell;
668 : * depending on user settings, this may attempt to take shortcuts
669 : * by assuming each cell has the same fills
670 : */
671 : void cacheContainedCells();
672 :
673 : /**
674 : * Fill the cached contained cells data structure for a given cell
675 : * @param[in] cell_info cell index, instance pair
676 : * @param[in] hint location hint used to accelerate the search
677 : * @param[out] map contained cell map
678 : */
679 : void setContainedCells(const cellInfo & cell_info,
680 : const Point & hint,
681 : std::map<cellInfo, containedCells> & map);
682 :
683 : /**
684 : * Check that the structure of the contained material cells for two cell matches;
685 : * i.e. this checks that the keys are the same and that the *number* of instances
686 : * of each filling material cell match.
687 : * @param[in] cell_info cell index, instance pair
688 : * @param[in] reference map we want to check against
689 : * @param[in] compare map we want to check
690 : */
691 : void checkContainedCellsStructure(const cellInfo & cell_info,
692 : containedCells & reference,
693 : containedCells & compare) const;
694 :
695 : /**
696 : * Set a minimum order for a volume quadrature rule
697 : * @param[in] volume_order order of the volume quadrature rule
698 : * @param[in] type string type of quadrature rule for printing a console message
699 : */
700 : void setMinimumVolumeQRules(Order & volume_order, const std::string & type);
701 :
702 : /// For keeping the output neat when using verbose
703 : std::string printNewline() const
704 : {
705 : if (_verbose)
706 : return "\n";
707 : else
708 : return "";
709 : }
710 :
711 : /// Loop over the elements in the MOOSE mesh and store the type of feedback applied by each.
712 : void storeElementPhase();
713 :
714 : /**
715 : * Loop over all the OpenMC cells and count the number of MOOSE elements to which the cell
716 : * is mapped based on phase.
717 : */
718 : void getCellMappedPhase();
719 :
720 : /// This function is used to ensure that each OpenMC cell only maps to a single phase
721 : void checkCellMappedPhase();
722 :
723 : /// Loop over all the OpenMC cells and get the element subdomain IDs that map to each cell
724 : void getCellMappedSubdomains();
725 :
726 : /**
727 : * Loop over all the OpenMC cells and compute the volume of the MOOSE elements that each
728 : * cell maps to
729 : */
730 : void computeCellMappedVolumes();
731 :
732 : /// Set up the mapping from MOOSE elements to OpenMC cells
733 : void initializeElementToCellMapping();
734 :
735 : /// Populate maps of MOOSE elements to OpenMC cells
736 : void mapElemsToCells();
737 :
738 : /**
739 : * A function which validates local tallies. This is done to ensure that at least one of the
740 : * tallies contains a heating score when running in eigenvalue mode. This must be done outside
741 : * of the constructor as tallies are added from an external system.
742 : */
743 : void validateLocalTallies();
744 :
745 : /// Add OpenMC tallies to facilitate the coupling
746 : void initializeTallies();
747 :
748 : /**
749 : * Reset any tallies previously added by Cardinal, by deleting them from OpenMC.
750 : * Also delete any mesh filters and meshes added to OpenMC for mesh filters.
751 : */
752 : void resetTallies();
753 :
754 : /**
755 : * Get one point inside each cell, for accelerating the particle search routine.
756 : * This function will get the centroid of the first global element in the lowest
757 : * rank in the cell.
758 : */
759 : void getPointInCell();
760 :
761 : /**
762 : * Compute the product of volume with a field across ranks and sum into a global map
763 : * @param[in] var_num variable to weight with volume, mapped by subdomain ID
764 : * @param[in] phase phases to compute the operation for
765 : * @param[in] scaling a scaling factor to apply, mapped by subdomain ID
766 : * @return volume-weighted field for each cell, in a global sense
767 : */
768 : std::map<cellInfo, Real> computeVolumeWeightedCellInput(
769 : const std::map<SubdomainID, std::pair<unsigned int, std::string>> & var_num,
770 : const std::vector<coupling::CouplingFields> * phase = nullptr,
771 : const std::map<SubdomainID, Real> * scaling = nullptr) const;
772 :
773 : /**
774 : * Send temperature from MOOSE to OpenMC by computing a volume average
775 : * and applying a single temperature per OpenMC cell
776 : */
777 : void sendTemperatureToOpenMC() const;
778 :
779 : /**
780 : * Send density from MOOSE to OpenMC by computing a volume average
781 : * and applying a single density per OpenMC cell.
782 : */
783 : void sendDensityToOpenMC() const;
784 :
785 : /**
786 : * Check if a mapped location is in the outer universe of a lattice
787 : * @param[in] level lattice level
788 : * @return whether the location is in the outer universe
789 : */
790 : void latticeOuterCheck(const Point & c, int level) const;
791 :
792 : /**
793 : * Report an error for a mapped location in an outer universe of a lattice
794 : * @param[in] c Mapped location
795 : * @param[in] level level of the mapped cell
796 : */
797 : void latticeOuterError(const Point & c, int level) const;
798 :
799 : /**
800 : * Find the OpenMC cell at a given point in space
801 : * @param[in] point point
802 : * @return whether OpenMC reported an error
803 : */
804 : bool findCell(const Point & point);
805 :
806 : /**
807 : * Checks that the contained material cells exactly match between a reference obtained
808 : * by calling openmc::Cell::get_contained_cells for each cell and a shortcut
809 : * approach that assumes all identical cells (which aren't simply just material fills)
810 : * has exactly the same contained material cells.
811 : * @param[in] reference reference map to compare against
812 : * @param[in] compare shortcut map to compare
813 : */
814 : void compareContainedCells(std::map<cellInfo, containedCells> & reference,
815 : std::map<cellInfo, containedCells> & compare) const;
816 :
817 : /**
818 : * Return all IDs of all Cardinal-mapped Tallies
819 : * @return all Cardinal-mapped Tally IDs
820 : */
821 : virtual std::vector<int32_t> getMappedTallyIDs() const override;
822 :
823 : /// A reference to the serialized auxvariable solution.
824 : NumericVector<Number> & _serialized_solution;
825 :
826 : /**
827 : * Whether to automatically compute the mapping of OpenMC cell IDs and
828 : * instances to the [Mesh].
829 : */
830 : const bool & _output_cell_mapping;
831 :
832 : /**
833 : * Where to get the initial OpenMC temperatures and densities from;
834 : * can be either hdf5 (from a properties.h5 file), xml (whatever is already
835 : * set in the XML files), or moose (meaning whatever ICs are set on the 'temperature_variables'
836 : * and 'density_variables'
837 : */
838 : const coupling::OpenMCInitialCondition _initial_condition;
839 :
840 : /// Type of relaxation to apply to the OpenMC tallies
841 : const relaxation::RelaxationEnum _relaxation;
842 :
843 : /**
844 : * Type of trigger to apply to k eigenvalue to indicate when
845 : * the simulation is complete. These can be used to on-the-fly adjust the number
846 : * of active batches in order to reach some desired criteria (which is specified
847 : * by this parameter).
848 : */
849 : const trigger::TallyTriggerTypeEnum _k_trigger;
850 :
851 : /**
852 : * Coordinate level in the OpenMC domain to use for mapping cells to mesh.
853 : * When using 'lowest_cell_level', this parameter indicates that the lowest
854 : * cell level is used, up until _cell_level.
855 : */
856 : unsigned int _cell_level;
857 :
858 : /**
859 : * Whether OpenMC properties (temperature and density) should be exported
860 : * after being updated in syncSolutions.
861 : */
862 : const bool & _export_properties;
863 :
864 : /// Whether or not the problem uses a skinner to regenerate the OpenMC geometry.
865 : const bool _using_skinner;
866 :
867 : /**
868 : * When the mesh changes during the simulation (either from adaptive mesh refinement
869 : * or deformation), the mapping from OpenMC cells to the [Mesh] must be re-established
870 : * after each OpenMC run.
871 : */
872 : bool _need_to_reinit_coupling;
873 :
874 : /**
875 : * If known a priori by the user, whether the tally cells (which are not simply material
876 : * fills) have EXACTLY the same contained material cells. This is a big optimization for
877 : * TRISO problems in setting up homogenized temperature/density feedback to OpenMC.
878 : *
879 : * The concept can best be explained with a pebble bed reactor.
880 : * If every pebble is filled with an identical TRISO universe, then the material fills
881 : * in each pebble are identical to one another except for a constant offset. This idea
882 : * can be used to then skip all but the first two openmc::Cell::get_contained_cells
883 : * calls (which are required in order to figure out the pattern by which pebble N is
884 : * incremented relative to pebble 1).
885 : *
886 : * When using this parameter, we HIGHLY recommend setting 'check_identical_cell_fills =
887 : * true' the first time you run your model. This will figure out the material cell fills using a
888 : * method that calls openmc::Cell::get_contained_cells for every tally cell, i.e. without assuming
889 : * anything about repeated structure in your OpenMC model. Setting 'identical_cell_fills'
890 : * without also setting 'check_identical_cell_fills = true' may result in SILENT
891 : * errors!!! So it is essential to be sure you've removed any error sources before you turn the
892 : * error check off to actually leverage the speedup.
893 : *
894 : * Note: for any tally cells that are just filled with a material, we use the approach
895 : * where openmc::Cell::get_contained_cells is called in full.
896 : *
897 : * This optimization will not work (and 'check_identical_cell_fills = true' *will*
898 : * catch these) for:
899 : * - any situation where tallied, non-material-fill pebbles have different fills
900 : * (such as if you have different TRISO lattices in each pebble)
901 : * - any situation where there is a "gap" in the incrementing of the material fill
902 : * instances (such as if pebble 89 does not map to 'tally_blocks', then the instance
903 : * shift for pebble 90 relative to pebble 1 is 89, when it should have been 90).
904 : */
905 : const bool _has_identical_cell_fills;
906 :
907 : /**
908 : * Whether we should rigorously check that each tally cell has identical fills;
909 : * this is SLOW for large TRISO problems, but is essential to ensure the accuracy of
910 : * 'identical_cell_fills'. Please set 'check_identical_cell_fills' to 'true' at least
911 : * once before running production cases to be sure the optimization can be applied.
912 : */
913 : const bool & _check_identical_cell_fills;
914 :
915 : /**
916 : * Whether it can be assumed that all of the tallies (both those set by the user
917 : * in the XML file, as well as those created automatically by Cardinal) are
918 : * spatially separate. This means that once a particle scores to one tally bin, it wouldn't
919 : * score to ANY other tally bins. This can dramatically increase tracking rates
920 : * for problems with many tallies.
921 : */
922 : const bool & _assume_separate_tallies;
923 :
924 : /**
925 : * Whether to map density according to each individual OpenMC cell (in which case an
926 : * error is thrown if you don't have a unique material in each cell) or by material.
927 : */
928 : bool _map_density_by_cell;
929 :
930 : /**
931 : * Whether the problem has density feedback blocks specified; note that this is NOT necessarily
932 : * indicative that the mapping was successful in finding any cells corresponding to those blocks
933 : */
934 : const bool _specified_density_feedback;
935 :
936 : /**
937 : * Whether the problem has temperature feedback blocks specified; note that this is NOT
938 : * necessarily indicative that the mapping was successful in finding any cells corresponding to
939 : * those blocks
940 : */
941 : const bool _specified_temperature_feedback;
942 :
943 : /// Whether any cell tallies exist.
944 : bool _has_cell_tallies = false;
945 :
946 : /// Whether any spatial mapping from OpenMC's cells to the mesh is needed
947 : bool _needs_to_map_cells;
948 :
949 : /**
950 : * A map of the filter objects created by the [Problem/Filters] block. The key for each filter is
951 : * it's corresponding MOOSE name to allow tallies to look up filters.
952 : */
953 : std::map<std::string, std::shared_ptr<FilterBase>> _filters;
954 :
955 : /// A vector of the tally objects created by the [Problem/Tallies] block.
956 : std::vector<std::shared_ptr<TallyBase>> _local_tallies;
957 :
958 : /// A list of all of the scores contained by the local tallies added in the [Tallies] block.
959 : std::vector<std::string> _all_tally_scores;
960 :
961 : /// Number of tallies scoring a particular score.
962 : std::map<std::string, unsigned int> _score_count;
963 :
964 : /// A vector of auxvariable ids added by the [Tallies] block.
965 : std::vector<std::vector<unsigned int>> _tally_var_ids;
966 :
967 : /// A vector of external (output-based) auxvariable ids added by the [Tallies] block.
968 : std::vector<std::vector<std::vector<unsigned int>>> _tally_ext_var_ids;
969 :
970 : /// Blocks in MOOSE mesh that provide density feedback
971 : std::vector<SubdomainID> _density_blocks;
972 :
973 : /// Blocks in MOOSE mesh that provide temperature feedback
974 : std::vector<SubdomainID> _temp_blocks;
975 :
976 : /// Blocks for which the cell fills are identical
977 : std::unordered_set<SubdomainID> _identical_cell_fill_blocks;
978 :
979 : /// Mapping of MOOSE elements to the OpenMC cell they map to (if any)
980 : std::vector<cellInfo> _elem_to_cell{};
981 :
982 : /// Phase of each cell
983 : std::map<cellInfo, coupling::CouplingFields> _cell_phase;
984 :
985 : /// Number of elements in the MOOSE mesh that exclusively provide density feedback
986 : int _n_moose_density_elems;
987 :
988 : /// Number of elements in the MOOSE mesh that exclusively provide temperature feedback
989 : int _n_moose_temp_elems;
990 :
991 : /// Number of elements in the MOOSE mesh which provide temperature+density feedback
992 : int _n_moose_temp_density_elems;
993 :
994 : /// Number of no-coupling elements in the MOOSE mesh
995 : int _n_moose_none_elems;
996 :
997 : /**
998 : * Number of MOOSE elements that exclusively provide temperature feedback,
999 : * and which successfully mapped to OpenMC cells
1000 : */
1001 : int _n_mapped_temp_elems;
1002 :
1003 : /**
1004 : * Number of MOOSE elements that exclusively provide density feedback,
1005 : * and which successfully mapped to OpenMC cells
1006 : */
1007 : int _n_mapped_density_elems;
1008 :
1009 : /**
1010 : * Number of MOOSE elements that provide temperature+density feedback,
1011 : * and which successfully mapped to OpenMC cells
1012 : */
1013 : int _n_mapped_temp_density_elems;
1014 :
1015 : /// Number of no-coupling elements mapped to OpenMC cells
1016 : int _n_mapped_none_elems;
1017 :
1018 : /// Total volume of uncoupled MOOSE mesh elements
1019 : Real _uncoupled_volume;
1020 :
1021 : /// Whether non-material cells are mapped
1022 : bool _material_cells_only{true};
1023 :
1024 : /// Mapping of OpenMC cell indices to a vector of MOOSE element IDs
1025 : std::map<cellInfo, std::vector<unsigned int>> _cell_to_elem;
1026 :
1027 : /// Mapping of OpenMC cell indices to a vector of MOOSE element IDs, on each local rank
1028 : std::map<cellInfo, std::vector<unsigned int>> _local_cell_to_elem;
1029 :
1030 : /// Mapping of OpenMC cell indices to the subdomain IDs each maps to
1031 : std::map<cellInfo, std::unordered_set<SubdomainID>> _cell_to_elem_subdomain;
1032 :
1033 : /// Mapping of elem subdomains to materials
1034 : std::map<SubdomainID, std::set<int32_t>> _subdomain_to_material;
1035 :
1036 : /**
1037 : * A point inside the cell, taken simply as the centroid of the first global
1038 : * element inside the cell. This is stored to accelerate the particle search.
1039 : */
1040 : std::map<cellInfo, Point> _cell_to_point;
1041 :
1042 : /**
1043 : * Volume associated with the mapped element space for each OpenMC cell; the unit
1044 : * for this volume is whatever is used in the [Mesh] block
1045 : */
1046 : std::map<cellInfo, Real> _cell_to_elem_volume;
1047 :
1048 : /**
1049 : * Volume associated with the actual OpenMC cell, computed by an optional
1050 : * OpenMCVolumeCalculation user object
1051 : */
1052 : std::map<cellInfo, Real> _cell_volume;
1053 :
1054 : /**
1055 : * Material-type cells contained within a cell; this is only populated if a cell
1056 : * is NOT indicated as having an identical fill
1057 : */
1058 : std::map<cellInfo, containedCells> _cell_to_contained_material_cells;
1059 :
1060 : /// Number of material-type cells contained within a cell
1061 : std::map<cellInfo, int32_t> _cell_to_n_contained;
1062 :
1063 : /// Whether the present transfer is the first transfer
1064 : static bool _first_transfer;
1065 :
1066 : /// Whether the diagnostic tables on initialization have already been printed
1067 : static bool _printed_initial;
1068 :
1069 : /// Whether a warning has already been printed about very long setup times (for TRISOs)
1070 : static bool _printed_triso_warning;
1071 :
1072 : /// Dummy particle to reduce number of allocations of particles for cell lookup routines
1073 : openmc::Particle _particle;
1074 :
1075 : /// Number of particles simulated in the first iteration in Dufek-Gudowski relaxation
1076 : unsigned int _n_particles_1;
1077 :
1078 : /// Mapping from temperature variable name to the subdomains on which to read it from
1079 : std::map<std::string, std::vector<SubdomainName>> _temp_vars_to_blocks;
1080 :
1081 : /// Mapping from density variable name to the subdomains on which to read it from
1082 : std::map<std::string, std::vector<SubdomainName>> _density_vars_to_blocks;
1083 :
1084 : /// Optional volume calculation for cells which map to MOOSE
1085 : OpenMCVolumeCalculation * _volume_calc;
1086 :
1087 : /// Userobject that maps from a partial-symmetry OpenMC model to a whole-domain [Mesh]
1088 : const SymmetryPointGenerator * _symmetry;
1089 :
1090 : /// Number of temperature-only feedback elements in each mapped OpenMC cell (global)
1091 : std::map<cellInfo, int> _n_temp;
1092 :
1093 : /// Number of density-only feedback elements in each mapped OpenMC cell (global)
1094 : std::map<cellInfo, int> _n_rho;
1095 :
1096 : /// Number of temperature+density feedback elements in each mapped OpenMC cell (global)
1097 : std::map<cellInfo, int> _n_temp_rho;
1098 :
1099 : /// Number of none elements in each mapped OpenMC cell (global)
1100 : std::map<cellInfo, int> _n_none;
1101 :
1102 : /// The tally to be used for normalizing all other tallies when running an eigenvalue calculation.
1103 : std::shared_ptr<TallyBase> _source_rate_norm_tally;
1104 :
1105 : /// The score index into "_source_rate_norm_tally".
1106 : int _source_rate_score = -1;
1107 :
1108 : #ifdef ENABLE_DAGMC
1109 : /// Optional skinner to re-generate the OpenMC geometry on-the-fly for DAGMC models
1110 : MoabSkinner * _skinner = nullptr;
1111 :
1112 : /// Pointer to DAGMC
1113 : std::shared_ptr<moab::DagMC> _dagmc = nullptr;
1114 : #endif
1115 :
1116 : /// Total number of unique OpenMC cell IDs + instances combinations
1117 : long unsigned int _n_openmc_cells;
1118 :
1119 : /// ID of the OpenMC universe corresponding to the DAGMC universe
1120 : int32_t _dagmc_universe_id;
1121 :
1122 : /// Whether the DAGMC universe is the root universe or not.
1123 : bool _dagmc_root_universe = true;
1124 :
1125 : /// ID of the OpenMC cell corresponding to the cell which uses the DAGMC universe as a fill.
1126 : int32_t _cell_using_dagmc_universe_id;
1127 :
1128 : /// The number of OpenMC surfaces before skinning occurs. This is required to properly reinitialize
1129 : /// the CSG geometry contained in the OpenMC model.
1130 : const int32_t _initial_num_openmc_surfaces;
1131 :
1132 : /// Conversion rate from eV to Joule
1133 : static constexpr Real EV_TO_JOULE = 1.6022e-19;
1134 :
1135 : /// Tolerance for setting zero tally
1136 : static constexpr Real ZERO_TALLY_THRESHOLD = 1e-12;
1137 :
1138 : private:
1139 : /**
1140 : * Update the number of particles according to the Dufek-Gudowski relaxation scheme
1141 : */
1142 : void dufekGudowskiParticleUpdate();
1143 :
1144 : /// Flattened cell IDs collected after parallel communication
1145 : std::vector<int32_t> _flattened_ids;
1146 :
1147 : /// Flattened cell instancess collected after parallel communication
1148 : std::vector<int32_t> _flattened_instances;
1149 :
1150 : /// Offsets for each cell instance in an identically-repeated universe
1151 : containedCells _instance_offsets;
1152 :
1153 : /// Offset for each cell relative to the first identical-fill cell
1154 : std::map<cellInfo, int32_t> _n_offset;
1155 :
1156 : /// First identical-fill cell
1157 : cellInfo _first_identical_cell;
1158 :
1159 : /// Materials in the first identical-fill cell
1160 : std::vector<int32_t> _first_identical_cell_materials;
1161 :
1162 : /// Whether OpenMCCellAverageProblem should use the displaced mesh
1163 : bool _use_displaced;
1164 :
1165 : /// Mapping from subdomain IDs to which aux variable to read temperature (K) from
1166 : std::map<SubdomainID, std::pair<unsigned int, std::string>> _subdomain_to_temp_vars;
1167 :
1168 : /// Mapping from subdomain IDs to which aux variable to read density (kg/m3) from
1169 : std::map<SubdomainID, std::pair<unsigned int, std::string>> _subdomain_to_density_vars;
1170 :
1171 : /// Mapping from subdomain IDs to the reference density (kg/m3).
1172 : std::map<SubdomainID, Real> _subdomain_to_ref_density;
1173 :
1174 : /// Mapping from cell index to the OpenMC Cell Material Modifier that contains its material list
1175 : std::map<int32_t, OpenMCCellMaterialFill *> _cell_material_modifiers;
1176 : };
|
