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 "MooseMesh.h"
22 : #include "MooseApp.h"
23 : #include "CardinalEnums.h"
24 : #include "NekBoundaryCoupling.h"
25 : #include "NekVolumeCoupling.h"
26 : #include "NekInterface.h"
27 : #include "NekUtility.h"
28 :
29 : /**
30 : * Representation of a nekRS surface mesh as a native MooseMesh. This is
31 : * constructed by interpolating from the surface Gauss-Lobatto-Legendre points
32 : * in nekRS to either a first-order (Quad4) or second-order (Quad9) mesh. This
33 : * mesh is only constructed for a user-specified set of boundaries in the nekRS
34 : * mesh with the 'boundary' parameter. Therefore, this class contains a mixture of
35 : * information related to the nekRS mesh (that nekRS solves its equations on)
36 : * versus the surface mesh constructed for data transfer with MOOSE (which is
37 : * only used by nekRS for the purpose of transferring its solution). All information
38 : * specific to the mesh nekRS actually uses for its solution are prefaced with either
39 : * '_nek' or 'nek' to help with this distinction.
40 : *
41 : * The nekRS mesh is currently implemented as a replicated mesh. On the nekRS side,
42 : * an Allgather is used to get the surface geometry information on each
43 : * nekRS process such that access from MOOSE can be performed on each process.
44 : *
45 : * TODO: The extension to higher than a second-order representation requires
46 : * some modifications to the formation of the mesh, as well as the interpolation
47 : * matrices used in NekRSProblem, because for 3rd order of higher, the equispaced
48 : * libMesh nodes no longer are a subset of the GLL ndoes.
49 : */
50 : class NekRSMesh : public MooseMesh
51 : {
52 : public:
53 : NekRSMesh(const InputParameters & parameters);
54 : static InputParameters validParams();
55 :
56 44 : NekRSMesh(const NekRSMesh & /* other_mesh */) = default;
57 :
58 : NekRSMesh & operator=(const NekRSMesh & other_mesh) = delete;
59 : virtual std::unique_ptr<MooseMesh> safeClone() const override;
60 :
61 : /// Save the initial volumetric mesh for volume mirror-based moving mesh problems
62 : void saveInitialVolMesh();
63 :
64 : /// Initialize previous displacement values to zero for boundary mirror-based moving mesh problems
65 : void initializePreviousDisplacements();
66 :
67 : /**
68 : * NekRS mesh polynomial order
69 : * @return NekRS polynomial order
70 : */
71 : int nekPolynomialOrder() const { return _nek_polynomial_order; }
72 :
73 : /**
74 : * Number of elements to build for each NekRS volume element
75 : * @return number of MOOSE elements per NekRS volume element
76 : */
77 7281291 : int nBuildPerVolumeElem() const { return _n_build_per_volume_elem; }
78 :
79 : /**
80 : * Number of elements to build for each NekRS surface element
81 : * @return number of MOOSE elements per NekRS surface element
82 : */
83 1154294 : int nBuildPerSurfaceElem() const { return _n_build_per_surface_elem; }
84 :
85 : /**
86 : * Get the initial mesh x coordinates
87 : * @return initial mesh x coordinates
88 : */
89 1024 : const std::vector<double> & nek_initial_x() const { return _initial_x; }
90 :
91 : /**
92 : * Get the initial mesh y coordinates
93 : * @return initial mesh y coordinates
94 : */
95 1024 : const std::vector<double> & nek_initial_y() const { return _initial_y; }
96 :
97 : /**
98 : * Get the initial mesh z coordinates
99 : * @return initial mesh z coordinates
100 : */
101 1024 : const std::vector<double> & nek_initial_z() const { return _initial_z; }
102 :
103 : /**
104 : * Get the previous x displacement
105 : * @return previous x displacement values
106 : */
107 : std::vector<double> & prev_disp_x() { return _prev_disp_x; }
108 :
109 : /**
110 : * Get the previous y displacement
111 : * @return previous y displacement values
112 : */
113 : std::vector<double> & prev_disp_y() { return _prev_disp_y; }
114 :
115 : /**
116 : * Get the previous z displacement
117 : * @return previous z displacement values
118 : */
119 : std::vector<double> & prev_disp_z() { return _prev_disp_z; }
120 :
121 : /**
122 : * Get the boundary coupling data structure
123 : * @return boundary coupling data structure
124 : */
125 33723 : const NekBoundaryCoupling & boundaryCoupling() const { return _boundary_coupling; }
126 :
127 : /**
128 : * Get the volume coupling data structure
129 : * @return volume coupling data structure
130 : */
131 1507378 : const NekVolumeCoupling & volumeCoupling() const { return _volume_coupling; }
132 :
133 : /// Add all the elements in the mesh to the MOOSE data structures
134 : virtual void addElems();
135 :
136 : /**
137 : * Get the order of the surface mesh; note that this is zero-indexed, so 0 = first, 1 = second
138 : * \return order
139 : */
140 : const order::NekOrderEnum & order() const { return _order; }
141 :
142 : /**
143 : * Get the number of quadrature points per coordinate direction in MOOSE's representation of
144 : * nekRS's mesh \return number of quadrature points per coordinate direction
145 : */
146 : int numQuadraturePoints1D() const;
147 :
148 : /**
149 : * Get the number of quadrature points per coordiante direction in nekRS's mesh
150 : * \return number of quadrature points per coordinate direction
151 : */
152 : int nekNumQuadraturePoints1D() const;
153 :
154 : /**
155 : * \brief Get the number of NekRS elements we rebuild in the MOOSE mesh
156 : *
157 : * This function is used to perform the data transfer routines in NekRSProblem
158 : * agnostic of whether we have surface or volume coupling.
159 : * return number of elements
160 : */
161 : const int & numElems() const { return _n_elems; }
162 :
163 : /**
164 : * \brief Get the number of vertices per element in MOOSE's representation of nekRS's mesh
165 : *
166 : * This function is used to perform the data transfer routines in NekRSProblem
167 : * agnostic of whether we have surface or volume coupling.
168 : * return number of vertices per element
169 : */
170 : const int & numVerticesPerElem() const { return _n_vertices_per_elem; }
171 :
172 : /**
173 : * \brief Get the libMesh node index from nekRS's GLL index ordering
174 : *
175 : * This function is used to perform the data transfer routines in NekRSProblem
176 : * agnostic of whether we have surface or volume coupling.
177 : * @param[in] gll_index nekRS GLL index
178 : * @return node index
179 : */
180 66611488 : int nodeIndex(const int gll_index) const { return (*_node_index)[gll_index]; }
181 :
182 : /**
183 : * Get the number of surface elements in MOOSE's representation of nekRS's mesh
184 : * \return number of surface elements
185 : */
186 : const int & numSurfaceElems() const { return _n_surface_elems; }
187 :
188 : /**
189 : * Get the total number of surface elements in nekRS's mesh
190 : * \return number of surface elements
191 : */
192 : const int & nekNumSurfaceElems() const { return _nek_n_surface_elems; }
193 :
194 : /**
195 : * Get the number of vertices per surface element in MOOSE's representation of nekRS's mesh
196 : * \return number of vertices per surface
197 : */
198 : const int & numVerticesPerSurface() const { return _n_vertices_per_surface; }
199 :
200 : /**
201 : * Get the number of volume elements in MOOSE's representation of nekRS's mesh
202 : * \return number of volume elements
203 : */
204 : const int & numVolumeElems() const { return _n_volume_elems; }
205 :
206 : /**
207 : * Get the total number of volume elements in nekRS's mesh
208 : * \return number of volume elements
209 : */
210 : const int & nekNumVolumeElems() const { return _nek_n_volume_elems; }
211 :
212 : /**
213 : * Get the number of vertices per volume element in MOOSE's representation of nekRS's mesh
214 : * \return number of vertices per volume
215 : */
216 : const int & numVerticesPerVolume() const { return _n_vertices_per_volume; }
217 :
218 : /**
219 : * Get the boundary ID for which nekRS and MOOSE are coupled
220 : * \return boundary ID
221 : */
222 290 : const std::vector<int> * boundary() const { return _boundary; }
223 :
224 : /**
225 : * Get whether the mesh permits volume-based coupling
226 : * \return whether mesh contains volume elements
227 : */
228 641129 : const bool & volume() const { return _volume; }
229 :
230 : /// Create a new element for a boundary mesh
231 : Elem * boundaryElem() const;
232 :
233 : /// Create a new element for a volume mesh
234 : Elem * volumeElem() const;
235 :
236 : virtual void buildMesh() override;
237 :
238 : /**
239 : * If running NekRS in JIT mode, we still need to make a mesh based on requirements
240 : * in MOOSE, so we just make a dummy mesh of a single Quad4 element
241 : */
242 : virtual void buildDummyMesh();
243 :
244 : /**
245 : * For the case of surface coupling only (i.e. no volume coupling), we create a surface
246 : * mesh for the elements on the specified boundary IDs
247 : */
248 : virtual void extractSurfaceMesh();
249 :
250 : /**
251 : * For the case of volume coupling only (i.e. no surface coupling), we create a volume
252 : * mesh for all volume elements
253 : */
254 : virtual void extractVolumeMesh();
255 :
256 : /**
257 : * Get the libMesh node index from nekRS's GLL index ordering
258 : * @param[in] gll_index nekRS GLL index
259 : * @return node index
260 : */
261 3837312 : int boundaryNodeIndex(const int gll_index) const { return _bnd_node_index[gll_index]; }
262 :
263 : /**
264 : * Get the libMesh node index from nekRS's GLL index ordering
265 : * @param[in] gll_index nekRS GLL index
266 : * @return node index
267 : */
268 14927532 : int volumeNodeIndex(const int gll_index) const { return _vol_node_index[gll_index]; }
269 :
270 : /**
271 : * Get the scaling factor applied to the nekRS mesh
272 : * @return scaling factor
273 : */
274 27914 : const Real & scaling() const { return _scaling; }
275 :
276 : /// Print diagnostic information related to the mesh
277 : virtual void printMeshInfo() const;
278 :
279 : /**
280 : * Processor id (rank) owning the given boundary element
281 : * @return processor id
282 : */
283 : int boundaryElemProcessorID(const int elem_id);
284 :
285 : /**
286 : * Processor id (rank) owning the given volume element
287 : * @return processor id
288 : */
289 : int volumeElemProcessorID(const int elem_id);
290 :
291 : /**
292 : * Get the number of faces of this global element that are on a coupling boundary
293 : * @param[in] elem_id global element ID
294 : * @return number of faces on a coupling boundary
295 : */
296 : int facesOnBoundary(const int elem_id) const;
297 :
298 : /**
299 : * Whether the mesh mirror is an exact representation of the NekRS mesh
300 : * @param return whether mesh mirror is exact
301 : */
302 20629786 : bool exactMirror() const { return _exact; }
303 :
304 : /**
305 : * Get the corner indices for the GLL points to be used in the mesh mirror
306 : * @return mapping of mesh mirror nodes to GLL points
307 : */
308 1950956 : std::vector<std::vector<int>> cornerIndices() const { return _corner_indices; }
309 :
310 : /**
311 : * Get the number of MOOSE elements we build for each NekRS element
312 : * @return MOOSE elements built for each NekRS element
313 : */
314 112344120 : int nMoosePerNek() const { return _n_moose_per_nek; }
315 :
316 : /**
317 : * Copy a new boundary displacement value for a given element's face
318 : * @param[in] src the source displacement at element e and face f
319 : * @param[in] e the element to which src belongs
320 : * @param[in] field the displacement field we are updating
321 : */
322 : void updateDisplacement(const int e, const double * src, const field::NekWriteEnum field);
323 :
324 : protected:
325 : /// Store the rank-local element and rank ownership for volume coupling
326 : void storeVolumeCoupling();
327 :
328 : /**
329 : * Store the rank-local element and rank ownership for boundary coupling;
330 : * this loops over the NekRS mesh and fetches relevant information on the boundaries
331 : */
332 : void storeBoundaryCoupling();
333 :
334 : /**
335 : * Sideset ID corresponding to a given volume element with give local face ID
336 : * @param[in] elem_id element local rank ID
337 : * @param[in] face_id element-local face ID
338 : * @return sideset ID (-1 means not one a boundary)
339 : */
340 : int boundary_id(const int elem_id, const int face_id);
341 :
342 : /**
343 : * Get the vertices defining the surface mesh interpolation from the
344 : * stored coupling information and store in _x, _y, and _z
345 : */
346 : void faceVertices();
347 :
348 : /**
349 : * Get the vertices defining the volume mesh interpolation from the
350 : * stored coupling information and store in _x, _y, and _z
351 : */
352 : void volumeVertices();
353 :
354 : /// Initialize members for the mesh and determine the GLL-to-node mapping
355 : void initializeMeshParams();
356 :
357 : /**
358 : * \brief Whether nekRS is coupled through volumes to MOOSE
359 : *
360 : * Unlike the case with _boundary, nekRS has no concept of volume/block IDs,
361 : * so we cannot have the user provide a vector of volumes that they want to
362 : * construct, so the best we can do is use a boolean here to turn on/off the
363 : * volume-based coupling for the entire mesh.
364 : */
365 : const bool & _volume;
366 :
367 : /// Boundary ID(s) through which to couple Nek to MOOSE
368 : const std::vector<int> * _boundary;
369 :
370 : /**
371 : * \brief Order of the surface interpolation between nekRS and MOOSE
372 : *
373 : * Order of the interpolation to be performed between nekRS and MOOSE;
374 : * options = FIRST, SECOND. For a first-order interpolation, nekRS's
375 : * solution is interpolated onto a first-order surface mesh (i.e. Quad4),
376 : * while for a second-order interpolation, nekRS's solution is interpolated
377 : * onto a second-order surface mesh (i.e. Quad9). Note that this is zero-indexed
378 : * so that an integer value of 0 = first-order, 1 = second-order, etc.
379 : **/
380 : const order::NekOrderEnum _order;
381 :
382 : /**
383 : * Whether the NekRS mesh mirror is an exact replica of the NekRS mesh.
384 : * If false (the default), then we build one MOOSE element for each NekRS element,
385 : * and the order of the MOOSE element is selected with 'order'. If true,
386 : * then we instead build one MOOSE element for each "first-order element"
387 : * within each NekRS high-order spectral element. In other words, if the NekRS
388 : * mesh is polynomial order 7, then we would build 7^2 MOOSE surface elements
389 : * for each NekRS surface element, and 7^3 MOOSE volume elements for each NekRS
390 : * volume element. The order of these elements will be first order.
391 : */
392 : const bool & _exact;
393 :
394 : /// Number of vertices per surface
395 : int _n_vertices_per_surface;
396 :
397 : /// Number of vertices per volume element
398 : int _n_vertices_per_volume;
399 :
400 : /**
401 : * \brief Spatial scaling factor to apply to the mesh
402 : *
403 : * nekRS is dimension agnostic - depending on the values used for the material
404 : * properties, the units of the mesh are arbitrary. Other apps that nekRS might
405 : * be coupled to could be in different units - to allow each app to use the
406 : * units that it wants, we can simply scale the NekRSMesh by a constant factor.
407 : * This will also adjust the heat flux coming in to nekRS by an appropriate factor.
408 : * For instance, if nekRS solves a problem in units of meters, but a BISON solution
409 : * is done on a mesh in units of centimeters, this scaling factor should be set to
410 : * 100. Note that other postprocessors will still be calculated on the nekRS mesh,
411 : * which will be in whatever units nekRS is internally using.
412 : */
413 : const Real & _scaling;
414 :
415 : /// Block ID for the fluid portion of the mesh mirror
416 : const unsigned int & _fluid_block_id;
417 :
418 : /// Block ID for the solid portion of the mesh mirror
419 : const unsigned int & _solid_block_id;
420 :
421 : /// Order of the nekRS solution
422 : int _nek_polynomial_order;
423 :
424 : /**
425 : * Number of NekRS surface elements in MooseMesh. The total number of surface
426 : * elements in the mesh mirror is _n_surface_elems * _n_build_per_surface_elem.
427 : */
428 : int _n_surface_elems;
429 :
430 : /// Number of MOOSE surface elements to build per NekRS surface element
431 : int _n_build_per_surface_elem;
432 :
433 : /**
434 : * Number of NekRS volume elements in MooseMesh. The total number of volume
435 : * elements in the mesh mirror is _n_volume_elems * _n_build_per_volume_elem.
436 : */
437 : int _n_volume_elems;
438 :
439 : /// Number of MOOSE volume elements to build per NekRS volume element
440 : int _n_build_per_volume_elem;
441 :
442 : /**
443 : * Number of NekRS elements we build in the MooseMesh. The total number of
444 : * elements in the mesh mirror is _n_elems * _n_moose_per_nek.
445 : */
446 : int _n_elems;
447 :
448 : /**
449 : * Number of MOOSE elements corresponding to each NekRS element, which depends on whether
450 : * building a boundary/volume mesh
451 : */
452 : int _n_moose_per_nek;
453 :
454 : /// Function returning the processor id which should own each element
455 : int (NekRSMesh::*_elem_processor_id)(const int elem_id);
456 :
457 : /// Number of vertices per element, which depends on whether building a boundary/volume mesh
458 : int _n_vertices_per_elem;
459 :
460 : /// Mapping of GLL nodes to libMesh nodes, which depends on whether building a boundary/volume mesh
461 : std::vector<int> * _node_index;
462 :
463 : /// Total number of surface elements in the nekRS problem
464 : int _nek_n_surface_elems;
465 :
466 : /// Total number of volume elements in the nekRS problem
467 : int _nek_n_volume_elems;
468 :
469 : /**
470 : * \brief "Phase" for each element (fluid = 0, solid = 1)
471 : *
472 : * TODO: This could be improved so that we don't need to collect this information,
473 : * but would require rewriting of a lot of the code used to assemble the members in
474 : * NekVolumeCoupling.
475 : */
476 : std::vector<int> _phase;
477 :
478 : /**
479 : * \brief \f$x\f$ coordinates of the current GLL points (which can move in time), for this rank
480 : *
481 : * This is ordered according to nekRS's internal geometry layout, and is indexed
482 : * first by the element and then by the node.
483 : */
484 : std::vector<double> _x;
485 :
486 : /**
487 : * \brief \f$y\f$ coordinates of the current GLL points (which can move in time), for this rank
488 : *
489 : * This is ordered according to nekRS's internal geometry layout, and is indexed
490 : * first by the element and then by the node.
491 : */
492 : std::vector<double> _y;
493 :
494 : /**
495 : * \brief \f$z\f$ coordinates of the current GLL points (which can move in time), for this rank
496 : *
497 : * This is ordered according to nekRS's internal geometry layout, and is indexed
498 : * first by the element and then by the node.
499 : */
500 : std::vector<double> _z;
501 :
502 : /**
503 : * \brief \f$x\f$ coordinates of the initial GLL points, for this rank
504 : *
505 : * This is ordered according to nekRS's internal geometry layout, and is indexed
506 : * first by the element and then by the node.
507 : */
508 : std::vector<double> _initial_x;
509 :
510 : /**
511 : * \brief \f$y\f$ coordinates of the initial GLL points, for this rank
512 : *
513 : * This is ordered according to nekRS's internal geometry layout, and is indexed
514 : * first by the element and then by the node.
515 : */
516 : std::vector<double> _initial_y;
517 :
518 : /**
519 : * \brief \f$z\f$ coordinates of the initial GLL points, for this rank
520 : *
521 : * This is ordered according to nekRS's internal geometry layout, and is indexed
522 : * first by the element and then by the node.
523 : */
524 : std::vector<double> _initial_z;
525 :
526 : ///@{
527 : /**
528 : * \f$x\f$, \f$y\f$, \f$z\f$ displacements of the boundary mesh mirror
529 : * for calculating displacement, in this rank
530 : **/
531 : std::vector<double> _prev_disp_x;
532 : std::vector<double> _prev_disp_y;
533 : std::vector<double> _prev_disp_z;
534 : ///@}
535 : /**
536 : * \brief Mapping of boundary GLL indices to MooseMesh node indices
537 : *
538 : * In nekRS, the GLL points are ordered by \f$x\f$, \f$y\f$, and \f$z\f$ coordinates,
539 : * but in order to construct sensible elements in Moose, we need to reorder these
540 : * points so that they match a libMesh-friendly node ordering. Without such a mapping,
541 : * we would construct triangles with zero/negative Jacobians instead of quad elements.
542 : * By indexing in the GLL index, this returns the node index.
543 : **/
544 : std::vector<int> _bnd_node_index;
545 :
546 : /**
547 : * \brief Mapping of volume GLL indices to MooseMesh node indices
548 : *
549 : * In nekRS, the GLL points are ordered by \f$x\f$, \f$y\f$, and \f$z\f$ coordinates,
550 : * but in order to construct sensible elements in Moose, we need to reorder these
551 : * points so that they match a libMesh-friendly node ordering. Without such a mapping,
552 : * we would construct triangles with zero/negative Jacobians instead of hex elements.
553 : * By indexing in the GLL index, this returns the node index.
554 : **/
555 : std::vector<int> _vol_node_index;
556 :
557 : /**
558 : * \brief Mapping of side indices to libMesh side indices
559 : *
560 : * nekRS uses its own side mapping that differs from that assumed in libMesh. In order
561 : * to assign the correct sideset IDs to the MooseMesh, we need to know the mapping between
562 : * these different conventions. By indexing in the nekRS side index, this returns the
563 : * libMesh side index.
564 : */
565 : std::vector<int> _side_index;
566 :
567 : /// Function pointer to the type of new element to add
568 : Elem * (NekRSMesh::*_new_elem)() const;
569 :
570 : /// Data structure holding mapping information for boundary coupling
571 : NekBoundaryCoupling _boundary_coupling;
572 :
573 : /// Data structure holding mapping information for volume coupling
574 : NekVolumeCoupling _volume_coupling;
575 :
576 : /// Pointer to NekRS's internal mesh data structure
577 : mesh_t * _nek_internal_mesh = nullptr;
578 :
579 : /// Corner indices for GLL points of mesh mirror elements
580 : std::vector<std::vector<int>> _corner_indices;
581 : };
|
