Line data Source code
1 : #pragma once
2 :
3 : #include "GeneralUserObject.h"
4 : #include "MaterialBase.h"
5 : #include "MooseMesh.h"
6 :
7 : #include "moab/Core.hpp"
8 : #include "moab/Skinner.hpp"
9 : #include "moab/GeomTopoTool.hpp"
10 : #include "MBTagConventions.hpp"
11 :
12 : /**
13 : * \brief Skins the [Mesh] according to individual bins for temperature, density, and subdomain ID
14 : *
15 : * Skins a [Mesh] according to temperature, density, and subdomain. The MOAB surfaces bounding
16 : * those grouped elements are then generated, providing geometry information needed for DAGMC
17 : * to then track particles on this new geometry.
18 : */
19 : class MoabSkinner : public GeneralUserObject
20 : {
21 : public:
22 : MoabSkinner(const InputParameters & parameters);
23 :
24 : static InputParameters validParams();
25 :
26 : virtual void execute() override;
27 :
28 : virtual void initialize() override;
29 :
30 : virtual void finalize() override;
31 :
32 0 : virtual void threadJoin(const UserObject & /* uo */) override {}
33 :
34 : /**
35 : * Wrap the error handling in MOAB to print errors to user
36 : * @param[in] input MOAB error code
37 : * @return error mode
38 : */
39 : virtual moab::ErrorCode check(const moab::ErrorCode input) const;
40 :
41 : std::string materialName(const unsigned int & block,
42 : const unsigned int & density,
43 : const unsigned int & temp) const;
44 :
45 : /// Perform the skinning operation
46 : virtual void update();
47 :
48 : /**
49 : * Set the names to be used for naming the subdomains in the skinned mesh.
50 : * This is used when material names are supplied externally rather than through
51 : * the 'material_blocks' and 'material_names' input parameters.
52 : * @param[in] names names for subdomains
53 : */
54 38 : virtual void setMaterialNames(std::vector<std::string> names) { _material_names = names; }
55 :
56 : /**
57 : * Get the total number of bins
58 : * @return total number of bins
59 : */
60 : unsigned int nBins() const;
61 :
62 : /**
63 : * Get the bin index for the temperature
64 : * @param[in] elem element
65 : * @return temperature bin index
66 : */
67 : virtual unsigned int getTemperatureBin(const Elem * const elem) const;
68 :
69 : /**
70 : * Get the bin index for the density
71 : * @param[in] elem element
72 : * @return density bin index
73 : */
74 : virtual unsigned int getDensityBin(const Elem * const elem) const;
75 :
76 : /**
77 : * Get the bin index for the subdomain
78 : * @param[in] elem element
79 : * @return subdomain bin index
80 : */
81 92470 : virtual unsigned int getSubdomainBin(const Elem * const elem) const
82 : {
83 92470 : return _blocks.at(elem->subdomain_id());
84 : }
85 :
86 : /**
87 : * Override the user parameter for use_displaced
88 : * @param[in] use whether to use the displaced mesh
89 : */
90 : void setUseDisplacedMesh(const bool & use);
91 :
92 : /**
93 : * Set the length multiplier to get from [Mesh] units into centimeters
94 : * @param[in] scale multiplier
95 : */
96 39 : virtual void setScaling(const Real & scale) { _scaling = scale; }
97 :
98 : /**
99 : * Set the verbosity level
100 : * @param[in] verbose whether to print diagnostic information
101 : */
102 39 : virtual void setVerbosity(const bool & verbose) { _verbose = verbose; }
103 :
104 : /**
105 : * Indicate whether this userobject is run by itself (for testing purposes)
106 : * or controlled by some other class.
107 : */
108 39 : virtual void makeDependentOnExternalAction() { _standalone = false; }
109 :
110 : /**
111 : * Get variable number in the auxiliary system
112 : * @param[in] name variable name
113 : * @param[in] param_name parameter name, for printing a helpful error message
114 : * @return variable number
115 : */
116 : unsigned int getAuxiliaryVariableNumber(const std::string & name,
117 : const std::string & param_name) const;
118 :
119 : /// Clear mesh data
120 : void reset();
121 :
122 : /**
123 : * Get total bin index given individual indices for the temperature, density, and subdomain bins
124 : * @param[in] temp_bin temperature bin
125 : * @param[in] density_bin density bin
126 : * @param[in] subdomain_bin subdomain ID bin
127 : * @return total bin index
128 : */
129 : virtual unsigned int getBin(const unsigned int & temp_bin,
130 : const unsigned int & density_bin,
131 : const unsigned int & subdomain_bin) const;
132 :
133 : /**
134 : * Whether the skinner builds a graveyard
135 : * @return whether a graveyard is built
136 : */
137 0 : virtual const bool & hasGraveyard() const { return _build_graveyard; }
138 :
139 : /**
140 : * Set the graveyard setting
141 : * @param[in] build whether to build a graveyard
142 : */
143 : void setGraveyard(bool build);
144 :
145 : /**
146 : * Number of density bins; if greater than 1, this means we must be re-generating
147 : * OpenMC materials during the course of the simulation.
148 : * @return number of density bins
149 : */
150 0 : virtual unsigned int nDensityBins() const { return _n_density_bins; }
151 :
152 : /**
153 : * Whether density skinning is applied
154 : * @return using density skinning
155 : */
156 42 : virtual bool hasDensitySkinning() const { return _bin_by_density; }
157 :
158 : /**
159 : * Get pointer to underlying moab interface
160 : * @return pointer to moab interface
161 : */
162 : const std::shared_ptr<moab::Interface> & moabPtr() const { return _moab; }
163 :
164 : protected:
165 : /// Boundary condition types that can be assigned to DAGMC surfaces
166 : enum class BoundaryConditionType
167 : {
168 : Transmission,
169 : Vacuum,
170 : Reflective
171 : };
172 :
173 : std::unique_ptr<NumericVector<Number>> _serialized_solution;
174 :
175 : /// MOAB interface
176 : std::shared_ptr<moab::Interface> _moab;
177 :
178 : /// Whether to print diagnostic information
179 : bool _verbose;
180 :
181 : /// Name of the temperature variable
182 : const std::string & _temperature_name;
183 :
184 : /// Lower bound of temperature bins
185 : const Real & _temperature_min;
186 :
187 : /// Upper bound of temperature bins
188 : const Real & _temperature_max;
189 :
190 : /// Number of temperature bins
191 : const unsigned int & _n_temperature_bins;
192 :
193 : /// Temperature bin width
194 : const Real _temperature_bin_width;
195 :
196 : /// Whether elements are binned by density (in addition to temperature and block)
197 : const bool _bin_by_density;
198 :
199 : /// Material names corresponding to each subdomain. These are used to name the
200 : /// new skinned volumes in MOAB
201 : std::vector<std::string> _material_names;
202 :
203 : /// Faceting tolerence needed by DAGMC
204 : const Real & _faceting_tol;
205 :
206 : /// Geometry tolerence needed by DAGMC
207 : const Real & _geom_tol;
208 :
209 : /// Multiplier on bounding box for inner surface of graveyard
210 : const Real & _graveyard_scale_inner;
211 :
212 : /// Multiplier on bounding box for outer surface of graveyard
213 : const Real & _graveyard_scale_outer;
214 :
215 : /// Whether to output the MOAB mesh skins to a .h5m file
216 : const bool & _output_skins;
217 :
218 : /// Whether to output the MOAB mesh to a .h5m file
219 : const bool & _output_full;
220 :
221 : /**
222 : * Whether to build a graveyard as two additional cube surfaces surrounding the mesh.
223 : * This is only needed if the skinned geometry is fed into a Monte Carlo code.
224 : */
225 : bool _build_graveyard;
226 :
227 : /// Whether to assign a material to the implicit complement region
228 : bool _set_implicit_complement_material = false;
229 :
230 : /// OpenMC material name or ID which will be assigned to the implicit complement
231 : std::string _implicit_complement_group_name;
232 :
233 : /// Whether the skinned mesh should be generated from a displaced mesh
234 : bool _use_displaced;
235 :
236 : /// Length multiplier to get from [Mesh] units into OpenMC's centimeters
237 : Real _scaling;
238 :
239 : /// Count number of times output files have been written
240 : unsigned int _n_write;
241 :
242 : /// Whether this class runs by itself, or is controlled by an external class
243 : bool _standalone;
244 :
245 : /// True when buildTetMesh() is called
246 : bool _tet_mesh_built;
247 :
248 : /// TET4 clone of the MOOSE mesh. Present only when the source mesh contains non-tetrahedral elements. Rebuilt each update() cycle.
249 : std::unique_ptr<MeshBase> _tet_mesh;
250 :
251 : /// Encode the whether the surface normal faces into or out of the volume
252 : enum Sense
253 : {
254 : BACKWARDS = -1,
255 : FORWARDS = 1
256 : };
257 :
258 : /// Encode MOAB information about volumes needed when creating surfaces
259 : struct VolData
260 : {
261 : moab::EntityHandle vol;
262 : Sense sense;
263 : };
264 :
265 : /// Get the MooseMesh (displaced or not, depending on _use_displaced)
266 : MooseMesh & getMooseMesh();
267 :
268 : /**
269 : * Return the mesh used for DAGMC geometry construction.
270 : * @return *_tet_mesh when the source mesh contained non-tetrahedral elements
271 : * and an internal conversion was performed; otherwise returns getMooseMesh().getMesh().
272 : */
273 : MeshBase & getDAGMCGeometryMesh();
274 :
275 : /// If the mesh has non-tetrahedral elements, clone it into _tet_mesh and convert to all-TET4
276 : void buildTetMesh();
277 :
278 : /**
279 : * Copy the libMesh [Mesh] into a MOAB mesh. This first loops through all of the
280 : * nodes, and rebuilds each as a MOAB vertex. Then, we loop over all of the elements
281 : * and rebuild each as a TET4 (if the libMesh mesh has TET10 elements, they are each
282 : * rebuilt into 8 TET4 elements). When the source mesh contained non-tetrahedral elements,
283 : * the geometry mesh is the internally-converted all-TET4 clone rather than the MOOSE mesh
284 : */
285 : void createMOABElems();
286 :
287 : /// Helper method to create MOAB tags
288 : virtual void createTags();
289 :
290 : /**
291 : * Helper method to create MOAB group entity set
292 : * @param[in] id ID for the group
293 : * @param[in] name name for the group
294 : * @param[in] group_set group of entities
295 : */
296 : void
297 : createGroup(const unsigned int & id, const std::string & name, moab::EntityHandle & group_set);
298 :
299 : /// Helper method to create MOAB volume entity set
300 : void
301 : createVol(const unsigned int & id, moab::EntityHandle & volume_set, moab::EntityHandle group_set);
302 :
303 : /// Helper method to create MOAB surface entity set
304 : void createSurf(const unsigned int & id,
305 : moab::EntityHandle & surface_set,
306 : moab::Range & faces,
307 : const std::vector<VolData> & voldata,
308 : BoundaryConditionType bc_type = BoundaryConditionType::Transmission);
309 :
310 : /// Helper method to create MOAB surfaces with no overlaps; bc_type is recorded on
311 : /// every surface created or matched by this call
312 : void createSurfaces(moab::Range & faces,
313 : VolData & voldata,
314 : unsigned int & surf_id,
315 : BoundaryConditionType bc_type = BoundaryConditionType::Transmission);
316 :
317 : /**
318 : * Create a MOAB surface from a bounding box
319 : */
320 : void createSurfaceFromBox(const BoundingBox & box,
321 : const VolData & voldata,
322 : unsigned int & surf_id,
323 : bool normalout,
324 : const Real & factor);
325 :
326 : /**
327 : * Create MOAB nodes from a bounding box
328 : * @param[in] box bounding box
329 : * @param[in] factor multiplicative factor to resize the bounding box sides
330 : * @return nodes
331 : */
332 : std::vector<moab::EntityHandle> createNodesFromBox(const BoundingBox & box,
333 : const Real & factor) const;
334 :
335 : /// Create 3 tri faces stemming from one corner of a cude (an open tetrahedron)
336 : void createCornerTris(const std::vector<moab::EntityHandle> & verts,
337 : unsigned int corner,
338 : unsigned int v1,
339 : unsigned int v2,
340 : unsigned int v3,
341 : bool normalout,
342 : moab::Range & surface_tris);
343 :
344 : /// Create MOAB tri surface element
345 : moab::EntityHandle createTri(const std::vector<moab::EntityHandle> & vertices,
346 : unsigned int v1,
347 : unsigned int v2,
348 : unsigned int v3);
349 :
350 : /// Add parent-child metadata relating a surface to its volume
351 : void updateSurfData(moab::EntityHandle surface_set, const VolData & data);
352 :
353 : /// Generic method to set the tags that DAGMC requires
354 : void
355 : setTags(moab::EntityHandle ent, std::string name, std::string category, unsigned int id, int dim);
356 :
357 : /// Helper function to wrap moab::tag_set_data for a string
358 : void setTagData(moab::Tag tag, moab::EntityHandle ent, std::string data, unsigned int SIZE);
359 :
360 : /// Helper function to wrap moab::tag_set_data for a generic pointer
361 : void setTagData(moab::Tag tag, moab::EntityHandle ent, void * data);
362 :
363 : /**
364 : * Get the node numberings for the MOAB TET4 elements to build for each [Mesh] element
365 : * @param[in] type element type
366 : */
367 : const std::vector<std::vector<unsigned int>> & getTetSets(ElemType type) const;
368 :
369 : /**
370 : * \brief Build a graveyard volume around the domain
371 : *
372 : * The graveyard is a containing volume which bounds the volume of interest. This is
373 : * only needed if the skinned geometry is going to be input into a Monte Carlo solver. For
374 : * performance reasons, a cubic shell is optimal. So, here we build two cubic surfaces,
375 : * both larger than the bounding box of the "actual" geometry. We name this region
376 : * "mat:Graveyard", so that when OpenMC parses the geometry it knows to assign "void"
377 : * to this region, and set vacuum BCs on the outer surfaces of the cubic shell. The
378 : * remaining space between the "actual" geometry and the inner graveyard surface is
379 : * treated as the implicit complement of the rest of the geometry (e.g. a transmissive region).
380 : */
381 : void buildGraveyard(unsigned int & vol_id, unsigned int & surf_id);
382 :
383 : /// Store a mapping from [Mesh] subdomain IDs to an index, to be used for binning by block ID
384 : virtual void findBlocks();
385 :
386 : /// Sort all the elements in the [Mesh] into bins for temperature, density, and subdomain.
387 : virtual void sortElemsByResults();
388 :
389 : /// Group the binned elems into local temperature regions and find their surfaces
390 : void findSurfaces();
391 :
392 : /**
393 : * Convert sideset names or numeric IDs to mesh BoundaryIDs and validate that each
394 : * is a sideset. Both string names and integer IDs are accepted.
395 : * @param[in] names sideset names or IDs
396 : * @param[in] param_name input parameter name, for error messages
397 : * @return resolved boundary IDs
398 : */
399 : std::set<BoundaryID> boundaryNamesToIDs(const std::vector<BoundaryName> & names,
400 : const std::string & param_name);
401 :
402 : /// Error if the same sideset ID appears in both vacuum and reflective BC sets.
403 : void checkBoundaryConditionOverlap() const;
404 :
405 : /**
406 : * Determine which BC type, if any, applies to a side of an element based on the
407 : * sidesets that side belongs to.
408 : * @param[in] elem element (in the DAGMC geometry mesh)
409 : * @param[in] side side index on the element
410 : * @param[in] boundary_info boundary info of the DAGMC geometry mesh
411 : * @return boundary condition type to assign to the side
412 : */
413 : BoundaryConditionType boundaryConditionType(const Elem * const elem,
414 : const unsigned int side,
415 : const libMesh::BoundaryInfo & boundary_info) const;
416 :
417 : /**
418 : * Classify the skinned triangles for one region into transmission, vacuum, and
419 : * reflective sets by checking each triangle against the BC sidesets.
420 : * @param[in] region MOAB tets forming the current local region
421 : * @param[in] skin skinned triangles bounding the region
422 : * @param[out] transmission_tris triangles with no boundary condition (transmission)
423 : * @param[out] vacuum_tris triangles assigned vacuum boundary conditions
424 : * @param[out] reflective_tris triangles assigned reflective boundary conditions
425 : */
426 : void splitSkinByBoundaryCondition(const moab::Range & region,
427 : const moab::Range & skin,
428 : moab::Range & transmission_tris,
429 : moab::Range & vacuum_tris,
430 : moab::Range & reflective_tris);
431 :
432 : /**
433 : * Classify and create DAGMC surfaces from find_skin() result.
434 : * Splits the skin into transmission, vacuum, and reflective subsets and calls
435 : * createSurfaces() for each type.
436 : * @param[in] region MOAB tets forming the current local region
437 : * @param[in] skin skinned triangles bounding the region
438 : * @param[in] voldata volume to associate with the created surfaces
439 : * @param[in,out] surf_id running counter of surface IDs
440 : */
441 : void createSurfacesFromSkin(const moab::Range & region,
442 : moab::Range & skin,
443 : VolData & voldata,
444 : unsigned int & surf_id);
445 :
446 : /**
447 : * Record that a DAGMC surface meshset has been assigned a BC type.
448 : * Called from createSurf() and from the overlap path in createSurfaces().
449 : * @param[in] surface_set surface meshset
450 : * @param[in] bc_type boundary condition type assigned to the surface
451 : */
452 : void recordBoundaryConditionSurface(moab::EntityHandle surface_set,
453 : BoundaryConditionType bc_type);
454 :
455 : /**
456 : * Get the boundary condition type previously recorded for a surface, if any
457 : * @param[in] surface_set surface meshset
458 : * @return recorded boundary condition type (Transmission if the surface has no explicit BC)
459 : */
460 : BoundaryConditionType recordedBoundaryCondition(moab::EntityHandle surface_set) const;
461 :
462 : /** Return the first group ID available after material, graveyard, and implicit-complement groups.
463 : * @return first available boundary condition group ID
464 : */
465 : unsigned int firstBoundaryConditionGroupID() const;
466 :
467 : /// Create the DAGMC BC group entity sets from the surfaces recorded
468 : /// during skinning and add them to the MOAB geometry
469 : void createBoundaryConditionGroups();
470 :
471 : /**
472 : * Get the DAGMC group name used to assign a boundary condition type
473 : * @param[in] bc_type boundary condition type
474 : * @return DAGMC group name
475 : */
476 : std::string boundaryConditionGroupName(BoundaryConditionType bc_type) const;
477 :
478 : /// Group a given bin into local regions
479 : /// NB elems in param is a copy, localElems is a reference
480 : void groupLocalElems(std::set<dof_id_type> elems, std::vector<moab::Range> & localElems);
481 :
482 : /// Clear MOAB entity sets
483 : bool resetMOAB();
484 :
485 : /// Find the surfaces for the provided range and add to group
486 : void findSurface(const moab::Range & region,
487 : moab::EntityHandle group,
488 : unsigned int & vol_id,
489 : unsigned int & surf_id,
490 : moab::EntityHandle & volume_set);
491 :
492 : /// Write MOAB volume and/or skin meshes to file
493 : virtual void write();
494 :
495 : /// Moab skinner for finding temperature surfaces
496 : std::unique_ptr<moab::Skinner> skinner;
497 :
498 : /// Topology tool for setting surface sense
499 : std::unique_ptr<moab::GeomTopoTool> gtt;
500 :
501 : /// Map from libmesh id to MOAB element entity handles
502 : std::map<dof_id_type, std::vector<moab::EntityHandle>> _id_to_elem_handles;
503 :
504 : /// Map from libMesh id to MOAB vertex handles
505 : std::unordered_map<dof_id_type, moab::EntityHandle> _node_id_to_handle;
506 :
507 : /// Save the first tet entity handle
508 : moab::EntityHandle offset;
509 :
510 : /// Name of the MOOSE variable containing the density
511 : std::string _density_name;
512 :
513 : /// Resolved vacuum BC sideset IDs (from 'vacuum_bcs_surfaces' input parameter)
514 : std::set<BoundaryID> _vacuum_bcs_surface_ids;
515 :
516 : /// Resolved reflective BC sideset IDs (from 'reflective_bcs_surfaces' input parameter)
517 : std::set<BoundaryID> _reflective_bcs_surface_ids;
518 :
519 : /// Map from DAGMC surface meshset to its assigned BC type, recorded during skinning.
520 : /// Transmission surfaces are not recorded (they need no DAGMC group).
521 : std::map<moab::EntityHandle, BoundaryConditionType> _surface_bc_types;
522 :
523 : /// Reverse map from MOAB tet entity handle to libMesh element ID. Populated in createMOABElems()
524 : std::unordered_map<moab::EntityHandle, dof_id_type> _elem_handle_to_id;
525 :
526 : /// Whether to assign boundary conditions to surfaces
527 : bool _set_bcs;
528 :
529 : /// Lower bound of density bins
530 : Real _density_min;
531 :
532 : /// Upper bound of density bins
533 : Real _density_max;
534 :
535 : /// Density bin width
536 : Real _density_bin_width;
537 :
538 : /// Number of density bins
539 : unsigned int _n_density_bins;
540 :
541 : /// Number of block bins
542 : unsigned int _n_block_bins;
543 :
544 : /// Mapping from total bin ID to a set of elements sorted into that bin
545 : std::vector<std::set<dof_id_type>> _elem_bins;
546 :
547 : /// Blocks in the [Mesh]
548 : std::map<SubdomainID, unsigned int> _blocks;
549 :
550 : /// Entity handle to represent the set of all tets
551 : moab::EntityHandle _all_tets;
552 :
553 : /// Save some topological data: map from surface handle to vol handle and sense
554 : std::map<moab::EntityHandle, std::vector<VolData>> surfsToVols;
555 :
556 : /// Tag for dimension for geometry
557 : moab::Tag geometry_dimension_tag;
558 :
559 : /// Tag for entitiy set ID
560 : moab::Tag id_tag;
561 :
562 : /// Tag for faceting tolerance
563 : moab::Tag faceting_tol_tag;
564 :
565 : /// Tag needed by DAGMC
566 : moab::Tag geometry_resabs_tag;
567 :
568 : /// Tag for type of entity set
569 : moab::Tag category_tag;
570 :
571 : /// Tag for name of entity set
572 : moab::Tag name_tag;
573 :
574 : /// Bounds of the temperature bins
575 : std::vector<Real> _temperature_bin_bounds;
576 :
577 : /// Bounds of the density bins
578 : std::vector<Real> _density_bin_bounds;
579 :
580 : /// Node ordering for a TET4 MOAB element, based on libMesh node numberings
581 : std::vector<std::vector<unsigned int>> _tet4_nodes;
582 :
583 : /**
584 : * Node ordering for eight TET4 MOAB elements, based on libMesh node numberings
585 : * for a TET10 element. We re-build the libMesh element into first-order MOAB elements.
586 : */
587 : std::vector<std::vector<unsigned int>> _tet10_nodes;
588 :
589 : /// Auxiliary variable number for temperature
590 : unsigned int _temperature_var_num;
591 :
592 : /// Auxiliary variable number for density
593 : unsigned int _density_var_num;
594 :
595 : /// Number of nodes per MOAB tet (which are first order, so TET4)
596 : const unsigned int NODES_PER_MOAB_TET = 4;
597 :
598 : /// Tolerance to use for comparing values to bin bounds
599 : const Real BIN_TOLERANCE = 1e-6;
600 :
601 : /// Map from mesh SubdomainID to OpenMC material name
602 : std::map<SubdomainID, std::string> _block_id_to_material_name;
603 : };
|