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 "CardinalEnums.h"
22 : #include "MooseTypes.h"
23 : #include "NekBoundaryCoupling.h"
24 : #include "NekVolumeCoupling.h"
25 :
26 : #include "inipp.hpp"
27 : #include "nekrs.hpp"
28 : #include "bcMap.hpp"
29 : #include "udf.hpp"
30 : #include "inipp.hpp"
31 : #include "mesh.h"
32 :
33 : #include "libmesh/point.h"
34 :
35 : #include <string>
36 : #include <vector>
37 :
38 : /**
39 : * \brief Cardinal-specific nekRS API
40 : *
41 : * nekRS ships with a rudimentary API in their nekrs namespace, but we need additional
42 : * functionality from within Cardinal. Many of these functions are quite basic and could
43 : * eventually be ported back into nekRS itself.
44 : */
45 : namespace nekrs
46 : {
47 :
48 : static int build_only;
49 :
50 : /// Allocate memory for the host mesh parameters
51 : void initializeHostMeshParameters();
52 :
53 : /// Update the mesh parameters on host
54 : void updateHostMeshParameters();
55 :
56 : dfloat * getSgeo();
57 : dfloat * getVgeo();
58 :
59 : /**
60 : * Check that the field specified can be accessed, e.g., if a user is requesting
61 : * to access temperature, the problem must have a temperature variable
62 : * @param[in] field field to check
63 : */
64 : void checkFieldValidity(const field::NekFieldEnum & field);
65 : void checkFieldValidity(const field::NekWriteEnum & field);
66 :
67 : /**
68 : * Compute y+ on the NekRS mesh
69 : * @param[in] boundary_id boundary(s) on which to compute y+
70 : * @param[in] index index into nek::scPtr where the wall distance is stored
71 : * @return max, min, average y+
72 : */
73 : std::vector<dfloat> yPlus(const std::vector<int> & boundary_id, const unsigned int & index);
74 :
75 : /**
76 : * Compute the three components of the viscous drag along a given boundary
77 : * @param[in] boundary boundary IDs for drag computation
78 : * @return components of drag (force fluid exerts on walls)
79 : */
80 : std::vector<dfloat> viscousDrag(const std::vector<int> & boundary);
81 :
82 : /**
83 : * Set the absolute tolerance for checking energy conservation in data transfers to Nek
84 : * @param[in] tol tolerance
85 : */
86 : void setAbsoluteTol(double tol);
87 :
88 : /**
89 : * Return the reference units for a usrwrk slot
90 : * @param[in] slot usrwrk slot
91 : * @return value by which to multiply the usrwrk slot to go from non-dimensional form into
92 : * dimensional form
93 : */
94 : Real scratchUnits(const int slot);
95 :
96 : /**
97 : * Inform backend if dimensionalization should be performed
98 : * @param[in] n if dimensionalize should be performed
99 : */
100 : void nondimensional(const bool n);
101 :
102 : /**
103 : * Set the relative tolerance for checking energy conservation in data transfers to Nek
104 : * @param[in] tol tolerance
105 : */
106 : void setRelativeTol(double tol);
107 :
108 : /**
109 : * Nek's runtime statistics are formed by collecting a timer of both the initialization
110 : * and accumulated run time. We unfortunately have to split this across multiple classes,
111 : * so if we want correct times we need to have NekInitAction save the value of the time
112 : * spent on initialization.
113 : * @param[in] time time spent on initialization
114 : */
115 : void setNekSetupTime(const double & time);
116 :
117 : /**
118 : * Get time spent on initialization
119 : * @return time spent on initialization
120 : */
121 : double getNekSetupTime();
122 :
123 : /**
124 : * Set the start time used by NekRS
125 : * @param[in] start start time
126 : */
127 : void setStartTime(const double & start);
128 :
129 : /**
130 : * Whether NekRS itself has been initialized yet
131 : * @return whether NekRS is initialized
132 : */
133 : bool isInitialized();
134 :
135 : /**
136 : * Write a field file containing a specific slot of the nrs->usrwrk scratch space;
137 : * this will write the field to the 'temperature' slot in a field file.
138 : * @param[in] slot index in the nrs->usrwrk array to write
139 : * @param[in] prefix prefix for file name
140 : * @param[in] time simulation time to write file for
141 : * @param[in] step time step index
142 : * @param[in] write_coords whether to write the mesh coordinates
143 : */
144 : void write_usrwrk_field_file(const int & slot, const std::string & prefix, const dfloat & time, const int & step, const bool & write_coords);
145 :
146 : /**
147 : * Write a field file containing pressure, velocity, and scalars with given prefix
148 : * @param[in] prefix three-character prefix
149 : * @param[in] time time
150 : * @param[in] step time step index
151 : */
152 : void write_field_file(const std::string & prefix, const dfloat time, const int & step);
153 :
154 : /**
155 : * Indicate whether NekRS was run in build-only mode (this doesn't actually
156 : * cause NekRS to run in build-only mode, but only provides an interface to
157 : * this information elsewhere).
158 : * @param[in] buildOnly whether NekRS is to be run in build-only mode
159 : */
160 : void buildOnly(int buildOnly);
161 :
162 : /**
163 : * Whether NekRS was run in JIT build-only mode
164 : * @return whether NekRS was run in build-only mode
165 : */
166 : int buildOnly();
167 :
168 : /**
169 : * Interpolate a volume between NekRS's GLL points and a given-order receiving/sending mesh
170 : */
171 : void interpolateVolumeHex3D(const double * I, double * x, int N, double * Ix, int M);
172 :
173 : /**
174 : * Whether nekRS's input file has CHT
175 : * @return whether nekRS input files model CHT
176 : */
177 : bool hasCHT();
178 :
179 : /**
180 : * Whether nekRS's input file indicates a moving mesh
181 : * @return whether nekRS's input file indicates a moving mesh
182 : */
183 : bool hasMovingMesh();
184 :
185 : /**
186 : * Whether nekRS's input file indicates a variable time stepping scheme
187 : * @return whether nekRS's input file indicates a variable time stepping
188 : */
189 : bool hasVariableDt();
190 :
191 : /**
192 : * Whether nekRS's input file has the blending mesh solver
193 : * @return whether nekRS's input file has a non-user [MESH] solver
194 : */
195 : bool hasBlendingSolver();
196 :
197 : /**
198 : * Whether nekRS's input file has the user mesh solver
199 : * @return whether nekRS's input file has [MESH] solver = user
200 : */
201 : bool hasUserMeshSolver();
202 :
203 : /**
204 : * Whether nekRS's input file intends to terminate the simulation based on a wall time
205 : * @return whether a wall time is used in nekRS to end the simulation
206 : */
207 : bool endControlElapsedTime();
208 :
209 : /**
210 : * Whether nekRS's input file intends to terminate the simulation based on an end time
211 : * @return whether an end time is used in nekRS to end the simulation
212 : */
213 : bool endControlTime();
214 :
215 : /**
216 : * Whether nekRS's input file intends to terminate the simulation based on a number of steps
217 : * @return whether a time step interval is used in nekRS to end the simulation
218 : */
219 : bool endControlNumSteps();
220 :
221 : /**
222 : * Offset increment for indexing into multi-volume arrays for the scalar fields.
223 : * This assumes that all scalars are the same length as the temperature scalar.
224 : * TODO: evaluate whether this works if nekRS uses CHT
225 : * @return scalar field offset
226 : */
227 : int scalarFieldOffset();
228 :
229 : /**
230 : * Offset increment for indexing into the velocity array
231 : * @return velocity field offset
232 : */
233 : int velocityFieldOffset();
234 :
235 : /**
236 : * Offset increment to use for generic slice indexing
237 : * @return field offset
238 : */
239 : int fieldOffset();
240 :
241 : /**
242 : * Get the "entire" NekRS mesh. For cases with a temperature scalar, this returns
243 : * nrs->meshT, which will cover both the fluid and solid regions if CHT is present.
244 : * For flow-only cases, this will return the flow mesh.
245 : * @return entire NekRS mesh
246 : */
247 : mesh_t * entireMesh();
248 :
249 : /**
250 : * Get the mesh for the flow solve
251 : * @return flow mesh
252 : */
253 : mesh_t * flowMesh();
254 :
255 : /**
256 : * Get the mesh for the temperature scalar
257 : * @return temperature mesh
258 : */
259 : mesh_t * temperatureMesh();
260 :
261 : /**
262 : * Get the mesh to act on
263 : * @param[in] pp_mesh which NekRS mesh to operate on
264 : * @return mesh to act on
265 : */
266 : mesh_t * getMesh(const nek_mesh::NekMeshEnum pp_mesh);
267 :
268 : /**
269 : * Get the process rank
270 : * @return process rank
271 : */
272 : int commRank();
273 :
274 : /**
275 : * Get the communicator size
276 : * @return communicator size
277 : */
278 : int commSize();
279 :
280 : /**
281 : * Whether nekRS's input file indicates that the problem has a temperature variable
282 : * @return whether the nekRS problem includes a temperature variable
283 : */
284 : bool hasTemperatureVariable();
285 :
286 : /**
287 : * Whether nekRS actually solves for temperature (as opposed to setting its solver to 'none')
288 : * @return whether nekRS will solve for temperature
289 : */
290 : bool hasTemperatureSolve();
291 :
292 : /**
293 : * Whether nekRS's input file indicates that the problem has a scalar0(scalarId) variable
294 : * @param[in] scalarId scalar number, i.e. for scalar03 scalarId=3
295 : * @return whether the nekRS problem includes the scalar0(scalarId) variable
296 : */
297 : bool hasScalarVariable(int scalarId);
298 :
299 : /**
300 : * Whether nekRS contains an OCCA kernel to apply a source to the passive scalar equations
301 : * @return whether nekRS has an OCCA kernel for apply a passive scalar source
302 : */
303 : bool hasHeatSourceKernel();
304 :
305 : /**
306 : * Whether the scratch space has already been allocated by the user
307 : * @return whether scratch space is already allocated
308 : */
309 : bool scratchAvailable();
310 :
311 : /**
312 : * Initialize scratch space for data to get sent into NekRS
313 : * @param[in] n_slots number of slots (for volume arrays) to allocate
314 : */
315 : void initializeScratch(const unsigned int & n_slots);
316 :
317 : /// Free the scratch space
318 : void freeScratch();
319 :
320 : /**
321 : * Get the viscosity used in the definition of the Reynolds number; note that
322 : * for dimensional cases, this is only guaranteed to be correct if the viscosity is constant.
323 : * @return constant dynamic viscosity
324 : */
325 : double viscosity();
326 :
327 : /**
328 : * Get the Prandtl number; note that for dimensional cases, this is only guaranteed
329 : * to be correct if the density, viscosity, heat capacity, and conductivity are constant.
330 : * @return constant Prandtl number
331 : */
332 : double Pr();
333 :
334 : /// Copy the deformation from host to device
335 : void copyDeformationToDevice();
336 :
337 : template <typename T>
338 : void allgatherv(const std::vector<int> & base_counts,
339 : const T * input,
340 : T * output,
341 : const int multiplier = 1);
342 :
343 : /**
344 : * Determine the receiving counts and displacements for all gather routines
345 : * @param[in] base_counts unit-wise receiving counts for each process
346 : * @param[out] counts receiving counts from each process
347 : * @param[out] displacement displacement for each process's counts
348 : * @param[in] multiplier optional multiplier on the face-based data
349 : */
350 : void displacementAndCounts(const std::vector<int> & base_counts,
351 : int * counts,
352 : int * displacement,
353 : const int multiplier);
354 :
355 : /**
356 : * Form the 2-D interpolation matrix from a starting GLL quadrature rule to an ending
357 : * GLL quadrature rule.
358 : * @param[out] I interpolation matrix
359 : * @param[in] starting_points number of points in the source quadrature rule
360 : * @param[in] ending_points number of points in the end quadrature rule
361 : */
362 : void interpolationMatrix(double * I, int starting_points, int ending_points);
363 :
364 : /**
365 : * Interpolate face data onto a new set of points
366 : * @param[in] scratch available scratch space for the calculation
367 : * @param[in] I interpolation matrix
368 : * @param[in] x face data to be interpolated
369 : * @param[in] N number of points in 1-D to be interpolated
370 : * @param[out] Ix interpolated data
371 : * @param[in] M resulting number of interpolated points in 1-D
372 : */
373 : void interpolateSurfaceFaceHex3D(
374 : double * scratch, const double * I, double * x, int N, double * Ix, int M);
375 :
376 : /**
377 : * Compute the face centroid given a local element ID and face ID (NOTE: returns in dimensional
378 : * form)
379 : * @param[in] local_elem_id local element ID on this rank
380 : * @param[in] local_face_id local face ID on the element
381 : * @return centroid
382 : */
383 : Point centroidFace(int local_elem_id, int local_face_id);
384 :
385 : /**
386 : * Compute the centroid given a local element ID (NOTE: returns in dimensional form)
387 : * @param[in] local_elem_id local element ID on this rank
388 : * @return centroid
389 : */
390 : Point centroid(int local_elem_id);
391 :
392 : /**
393 : * Get the coordinate given a local element ID and local node ID (NOTE: returns in dimensional form)
394 : * @param[in] local_elem_id local element ID on this rank
395 : * @param[in] local_node_id local node ID on this element
396 : * @return point
397 : */
398 : Point gllPoint(int local_elem_id, int local_node_id);
399 :
400 : /**
401 : * Get the coordinate given a local element ID, a local face ID, and local node ID (NOTE: returns in
402 : * dimensional form)
403 : * @param[in] local_elem_id local element ID on this rank
404 : * @param[in] local_face_id local face ID on this element
405 : * @param[in] local_node_id local node ID on this element
406 : * @return point
407 : */
408 : Point gllPointFace(int local_elem_id, int local_face_id, int local_node_id);
409 :
410 : /**
411 : * Integrate the scratch space over boundaries
412 : * @param[in] slot slot in scratch space
413 : * @param[in] boundary boundaries over which to integrate the scratch space
414 : * @param[in] pp_mesh portion of NekRS mesh to integrate over
415 : * @return boundary integrated scratch space, with one value per sideset
416 : */
417 : std::vector<double> usrwrkSideIntegral(const unsigned int & slot,
418 : const std::vector<int> & boundary,
419 : const nek_mesh::NekMeshEnum pp_mesh);
420 :
421 : /**
422 : * Volume integrate the scratch space
423 : * @param[in] slot slot in scratch space to i ntegrat
424 : * @param[in] pp_mesh NekRS mesh to integrate over
425 : * @return volume integrated scratch space
426 : */
427 : double usrwrkVolumeIntegral(const unsigned int & slot, const nek_mesh::NekMeshEnum pp_mesh);
428 :
429 : /**
430 : * Scale a slot in the usrwrk by a fixed value (multiplication)
431 : * @param[in] slot slot in usrwrk to modify
432 : * @param[in] value value to multiply on scratch slot
433 : */
434 : void scaleUsrwrk(const unsigned int & slot, const dfloat & value);
435 :
436 : /**
437 : * Compute the area of a set of boundary IDs
438 : * @param[in] boundary_id nekRS boundary IDs for which to perform the integral
439 : * @param[in] pp_mesh which NekRS mesh to operate on
440 : * @return area integral
441 : */
442 : double area(const std::vector<int> & boundary_id, const nek_mesh::NekMeshEnum pp_mesh);
443 :
444 : /**
445 : * Compute the area integral of a given integrand over a set of boundary IDs
446 : * @param[in] boundary_id nekRS boundary IDs for which to perform the integral
447 : * @param[in] integrand field to integrate
448 : * @param[in] pp_mesh which NekRS mesh to operate on
449 : * @return area integral of a field
450 : */
451 : double sideIntegral(const std::vector<int> & boundary_id, const field::NekFieldEnum & integrand,
452 : const nek_mesh::NekMeshEnum pp_mesh);
453 :
454 : /**
455 : * Compute the volume over the entire scalar mesh
456 : * @param[in] pp_mesh which NekRS mesh to operate on
457 : * @return volume integral
458 : */
459 : double volume(const nek_mesh::NekMeshEnum pp_mesh);
460 :
461 : /**
462 : * Dimensionalize a volume
463 : * @param[in] integral integral to dimensionalize
464 : */
465 : void dimensionalizeVolume(double & integral);
466 :
467 : /**
468 : * Dimensionalize an area
469 : * @param[in] integral integral to dimensionalize
470 : */
471 : void dimensionalizeArea(double & integral);
472 :
473 : /**
474 : * Dimensionalize a given integral of f over volume, i.e. fdV
475 : * @param[in] integrand field to dimensionalize
476 : * @param[in] volume volume of the domain (only used for dimensionalizing temperature)
477 : * @param[in] integral integral to dimensionalize
478 : */
479 : void dimensionalizeVolumeIntegral(const field::NekFieldEnum & integrand,
480 : const Real & volume,
481 : double & integral);
482 :
483 : /**
484 : * Dimensionalize a given integral of f over a side, i.e. fdS
485 : * @param[in] integrand field to dimensionalize
486 : * @param[in] area area of the boundary
487 : * @param[in] integral integral to dimensionalize
488 : */
489 : void dimensionalizeSideIntegral(const field::NekFieldEnum & integrand,
490 : const Real & area,
491 : double & integral);
492 :
493 : /**
494 : * Dimensionalize a given integral of f over a side, i.e. fdS
495 : * @param[in] integrand field to dimensionalize
496 : * @param[in] boundary_id boundary IDs for the integral
497 : * @param[in] integral integral to dimensionalize
498 : * @param[in] pp_mesh which NekRS mesh to operate on
499 : */
500 : void dimensionalizeSideIntegral(const field::NekFieldEnum & integrand,
501 : const std::vector<int> & boundary_id,
502 : double & integral,
503 : const nek_mesh::NekMeshEnum pp_mesh);
504 :
505 : /**
506 : * Compute the volume integral of a given integrand over the entire scalar mesh
507 : * @param[in] integrand field to integrate
508 : * @param[in] volume volume of the domain (only used for dimensionalizing temperature)
509 : * @param[in] pp_mesh which NekRS mesh to operate on
510 : * @return volume integral of a field
511 : */
512 : double volumeIntegral(const field::NekFieldEnum & integrand,
513 : const double & volume,
514 : const nek_mesh::NekMeshEnum pp_mesh);
515 :
516 : /**
517 : * Compute the mass flowrate over a set of boundary IDs
518 : * @param[in] boundary_id nekRS boundary IDs for which to compute the mass flowrate
519 : * @param[in] pp_mesh which NekRS mesh to operate on
520 : * @return mass flowrate
521 : */
522 : double massFlowrate(const std::vector<int> & boundary_id,
523 : const nek_mesh::NekMeshEnum pp_mesh);
524 :
525 : /**
526 : * Compute the mass flux weighted integral of a given integrand over a set of boundary IDs
527 : * @param[in] boundary_id nekRS boundary IDs for which to perform the integral
528 : * @param[in] integrand field to integrate and weight by mass flux
529 : * @param[in] pp_mesh which NekRS mesh to operate on
530 : * @return mass flux weighted area average of a field
531 : */
532 : double sideMassFluxWeightedIntegral(const std::vector<int> & boundary_id,
533 : const field::NekFieldEnum & integrand,
534 : const nek_mesh::NekMeshEnum pp_mesh);
535 :
536 : /**
537 : * Compute the integral of pressure on a surface, multiplied by the unit normal
538 : * of the surface with a specified direction vector. This represents the force
539 : * that the fluid exerts ON the boundary.
540 : * @param[in] boundary_id NekRS boundary IDs for which to perform the integral
541 : * @param[in] direction unit vector to dot with the boundary surface normal
542 : * @param[in] pp_mesh which NekRS mesh to operate on
543 : * @return pressure surface force, along a particular direction
544 : */
545 : double pressureSurfaceForce(const std::vector<int> & boundary_id, const Point & direction, const nek_mesh::NekMeshEnum pp_mesh);
546 :
547 : /**
548 : * Compute the heat flux over a set of boundary IDs
549 : * @param[in] boundary_id nekRS boundary IDs for which to perform the integral
550 : * @param[in] pp_mesh which NekRS mesh to operate on
551 : * @return heat flux area integral
552 : */
553 : double heatFluxIntegral(const std::vector<int> & boundary_id,
554 : const nek_mesh::NekMeshEnum pp_mesh);
555 :
556 : /**
557 : * Limit the temperature in nekRS to within the range of [min_T, max_T]
558 : * @param[in] min_T minimum temperature allowable in nekRS
559 : * @param[in] max_T maximum temperature allowable in nekRS
560 : */
561 : void limitTemperature(const double * min_T, const double * max_T);
562 :
563 : /**
564 : * Compute the gradient of a volume field
565 : * @param[in] offset in the gradient field for each component (grad_x, grad_y, or grad_z)
566 : * @param[in] e element ID to compute gradient
567 : * @param[in] f field to compute the gradient of
568 : * @param[in] pp_mesh which NekRS mesh to operate on
569 : * @param[out] grad_f gradient of field
570 : */
571 : void gradient(const int offset,
572 : const int e,
573 : const double * f,
574 : double * grad_f,
575 : const nek_mesh::NekMeshEnum pp_mesh);
576 :
577 : /**
578 : * Find the extreme value of a given field over the entire nekRS domain
579 : * @param[in] field field to find the minimum value of
580 : * @param[in] pp_mesh which NekRS mesh to operate on
581 : * @param[in] max whether to take the maximum (or if false, the minimum)
582 : * @return max or min value of field in volume
583 : */
584 : double volumeExtremeValue(const field::NekFieldEnum & field,
585 : const nek_mesh::NekMeshEnum pp_mesh,
586 : const bool max);
587 :
588 : /**
589 : * Find the extreme of a given field over a set of boundary IDs
590 : * @param[in] boundary_id nekRS boundary IDs for which to find the extreme value
591 : * @param[in] field field to find the maximum value of
592 : * @param[in] pp_mesh which NekRS mesh to operate on
593 : * @param[in] max whether to take the maximum (or if false, the minimum)
594 : * @return max or min value of field on boundary
595 : */
596 : double sideExtremeValue(const std::vector<int> & boundary_id, const field::NekFieldEnum & field,
597 : const nek_mesh::NekMeshEnum pp_mesh, const bool max);
598 :
599 : /**
600 : * Number of faces per element; because NekRS only supports HEX20, this should be 6
601 : * @return number of faces per mesh element
602 : */
603 : int Nfaces();
604 :
605 : /**
606 : * Whether the specific boundary is a flux boundary
607 : * @param[in] boundary boundary ID
608 : * @return whether boundary is a flux boundary
609 : */
610 : bool isHeatFluxBoundary(const int boundary);
611 :
612 : /**
613 : * Whether the specific boundary is a moving mesh boundary
614 : * @param[in] boundary boundary ID
615 : * @return whether boundary is a moving mesh boundary
616 : */
617 : bool isMovingMeshBoundary(const int boundary);
618 :
619 : /**
620 : * Whether the specific boundary is a specified temperature boundary
621 : * @param[in] boundary boundary ID
622 : * @return whether boundary is a temperature boundary
623 : */
624 : bool isTemperatureBoundary(const int boundary);
625 :
626 : /**
627 : * String name indicating the temperature boundary condition type on a given boundary
628 : * @param[in] boundary boundary ID
629 : * @return string name of boundary condition type
630 : */
631 : const std::string temperatureBoundaryType(const int boundary);
632 :
633 : /**
634 : * Polynomial order used in nekRS solution
635 : * @return polynomial order
636 : */
637 : int polynomialOrder();
638 :
639 : /**
640 : * Total number of volume elements in nekRS mesh summed over all processes
641 : * @return number of volume elements
642 : */
643 : int Nelements();
644 :
645 : /**
646 : * Mesh dimension
647 : * @return mesh dimension
648 : */
649 : int dim();
650 :
651 : /**
652 : * \brief Number of vertices required to define an element face
653 : * Vertices refer to the points required to place the "corners" of an element face,
654 : * and _not_ the quadrature points. For instance, for hexahedral elements, the number of vertices
655 : * per face is 4 regardless of the polynomial order.
656 : * @return Number of vertices per element face
657 : */
658 : int NfaceVertices();
659 :
660 : /**
661 : * Total number of element faces on a boundary of the nekRS mesh summed over all processes
662 : * @return number of boundary element faces
663 : */
664 : int NboundaryFaces();
665 :
666 : /**
667 : * Number of boundary IDs in the nekRS mesh
668 : * @return number of boundary IDs
669 : */
670 : int NboundaryID();
671 :
672 : /**
673 : * Whether the provided boundary IDs are all valid in the nekRS mesh
674 : * @param[in] boundary_id vector of boundary IDs to check
675 : * @param[out] first_invalid_id first invalid ID encountered for printing an error on the MOOSE side
676 : * @param[out] n_boundaries maximum valid boundary ID for printing an error on the MOOSE side
677 : * @return whether all boundaries are valid
678 : */
679 : bool
680 : validBoundaryIDs(const std::vector<int> & boundary_id, int & first_invalid_id, int & n_boundaries);
681 :
682 : /**
683 : * Store the rank-local element, element-local face, and rank ownership for boundary coupling
684 : * @param[in] boundary_id boundaries through which nekRS will be coupled
685 : * @param[out] N total number of surface elements
686 : */
687 : void storeBoundaryCoupling(const std::vector<int> & boundary_id, int & N);
688 :
689 : /**
690 : * Characteristic scales assumed in nekRS if using a non-dimensional solution; initial values
691 : * are applied, which will be overridden by the DimensionalizeAction in Cardinal.
692 : */
693 : struct characteristicScales
694 : {
695 : double U_ref = 1;
696 : double T_ref = 0;
697 : double dT_ref = 1;
698 : double P_ref = 1;
699 : double L_ref = 1;
700 : double A_ref = 1;
701 : double V_ref = 1;
702 : double rho_ref = 1;
703 : double Cp_ref = 1;
704 : double flux_ref = 1;
705 : double source_ref = 1;
706 : double t_ref = 1;
707 : double s01_ref = 0;
708 : double ds01_ref = 1;
709 : double s02_ref = 0;
710 : double ds02_ref = 1;
711 : double s03_ref = 0;
712 : double ds03_ref = 1;
713 : };
714 :
715 : /**
716 : * Get pointer to various solution functions (for reading only) based on enumeration
717 : * @param[in] field field to return a pointer to
718 : * @return function pointer to method that returns said field as a function of GLL index
719 : */
720 : double (*solutionPointer(const field::NekFieldEnum & field))(int, int);
721 : double (*solutionPointer(const field::NekWriteEnum & field))(int, int);
722 :
723 : /**
724 : * Get the scalar01 solution at given GLL index
725 : * @param[in] id GLL index
726 : * @return scalar01 value at index
727 : */
728 : double get_scalar01(const int id, const int surf_offset);
729 :
730 : /**
731 : * Get the scalar02 solution at given GLL index
732 : * @param[in] id GLL index
733 : * @return scalar02 value at index
734 : */
735 : double get_scalar02(const int id, const int surf_offset);
736 :
737 : /**
738 : * Get the scalar03 solution at given GLL index
739 : * @param[in] id GLL index
740 : * @return scalar03 value at index
741 : */
742 : double get_scalar03(const int id, const int surf_offset);
743 :
744 : /**
745 : * Get the usrwrk zeroth slice at given GLL index
746 : * @param[in] id GLL index
747 : * @return zeroth slice of usrwrk value at index
748 : */
749 : double get_usrwrk00(const int id, const int surf_offset);
750 :
751 : /**
752 : * Get the usrwrk first slice at given GLL index
753 : * @param[in] id GLL index
754 : * @return first slice of usrwrk value at index
755 : */
756 : double get_usrwrk01(const int id, const int surf_offset);
757 :
758 : /**
759 : * Get the usrwrk second slice at given GLL index
760 : * @param[in] id GLL index
761 : * @return second slice of usrwrk value at index
762 : */
763 : double get_usrwrk02(const int id, const int surf_offset);
764 :
765 : /**
766 : * Get the temperature solution at given GLL index
767 : * @param[in] id GLL index
768 : * @return temperature value at index
769 : */
770 : double get_temperature(const int id, const int surf_offset);
771 :
772 : /**
773 : * Get the pressure solution at given GLL index
774 : * @param[in] id GLL index
775 : * @return pressure value at index
776 : */
777 : double get_pressure(const int id, const int surf_offset);
778 :
779 : /**
780 : * Return unity, for cases where the integrand or operator we are generalizing acts on 1
781 : * @param[in] id GLL index
782 : * @return unity
783 : */
784 : double get_unity(const int id, const int surf_offset);
785 :
786 : /**
787 : * Get the x-velocity at given GLL index
788 : * @param[in] id GLL index
789 : * @return x-velocity at index
790 : */
791 : double get_velocity_x(const int id, const int surf_offset);
792 :
793 : /**
794 : * Get the y-velocity at given GLL index
795 : * @param[in] id GLL index
796 : * @return y-velocity at index
797 : */
798 : double get_velocity_y(const int id, const int surf_offset);
799 :
800 : /**
801 : * Get the z-velocity at given GLL index
802 : * @param[in] id GLL index
803 : * @return z-velocity at index
804 : */
805 : double get_velocity_z(const int id, const int surf_offset);
806 :
807 : /**
808 : * Get the magnitude of the velocity solution at given GLL index
809 : * @param[in] id GLL index
810 : * @return velocity magnitude at index
811 : */
812 : double get_velocity(const int id, const int surf_offset);
813 :
814 : /**
815 : * Get the x-velocity squared at given GLL index
816 : * @param[in] id GLL index
817 : * @return square of x-velocity at index
818 : */
819 : double get_velocity_x_squared(const int id, const int surf_offset);
820 :
821 : /**
822 : * Get the y-velocity squared at given GLL index
823 : * @param[in] id GLL index
824 : * @return square of y-velocity at index
825 : */
826 : double get_velocity_y_squared(const int id, const int surf_offset);
827 :
828 : /**
829 : * Get the z-velocity squared at given GLL index
830 : * @param[in] id GLL index
831 : * @return square of z-velocity at index
832 : */
833 : double get_velocity_z_squared(const int id, const int surf_offset);
834 :
835 : /**
836 : * Initialize the characteristic scales for a nondimesional solution
837 : * @param[in] U reference velocity
838 : * @param[in] T reference temperature
839 : * @param[in] dT reference temperature range
840 : * @param[in] L reference length scale
841 : * @param[in] rho reference density
842 : * @param[in] Cp reference heat capacity
843 : * @param[in] s01 reference scalar01
844 : * @param[in] ds01 reference s01 range
845 : * @param[in] s02 reference scalar02
846 : * @param[in] ds02 reference s02 range
847 : * @param[in] s03 reference scalar03
848 : * @param[in] ds03 reference s03 range
849 : */
850 : void initializeDimensionalScales(const double U,
851 : const double T,
852 : const double dT,
853 : const double L,
854 : const double rho,
855 : const double Cp,
856 : const double s01,
857 : const double ds01,
858 : const double s02,
859 : const double ds02,
860 : const double s03,
861 : const double ds03);
862 :
863 : /**
864 : * \brief Return the reference divisor scale that defines the non-dimensional field
865 : *
866 : * All fields in NekRS are assumed non-dimensionalized according to the general form
867 : * f (non-dimensional) = (f - f_ref) / df
868 : *
869 : * so that to define a nondimensionalization requires two scales: the divisor scale
870 : * (df) and the additive scale (f_ref).
871 : *
872 : * @param[in] field physical interpretation of value
873 : * @param[out] value nondimensional divisor scale (df)
874 : */
875 : Real nondimensionalDivisor(const field::NekFieldEnum & field);
876 : Real nondimensionalDivisor(const field::NekWriteEnum & field);
877 :
878 : /**
879 : * All fields in NekRS are assumed non-dimensionalized according to the general form
880 : * f (non-dimensional) = (f - f_ref) / df
881 : *
882 : * so that to define a nondimensionalization requires two scales: the divisor scale
883 : * (df) and the additive scale (f_ref).
884 : *
885 : * @param[in] field physical interpretation of value
886 : * @param[out] value nondimensional additive scale (f_ref)
887 : *
888 : */
889 : Real nondimensionalAdditive(const field::NekFieldEnum & field);
890 : Real nondimensionalAdditive(const field::NekWriteEnum & field);
891 :
892 : /**
893 : * Get the reference length scale
894 : * @return reference length scale
895 : */
896 : double referenceLength();
897 :
898 : /**
899 : * Get the reference time scale
900 : * @return reference time scale
901 : */
902 : double referenceTime();
903 :
904 : /**
905 : * Get the reference area scale
906 : * @return reference area scale
907 : */
908 : double referenceArea();
909 :
910 : /**
911 : * Get the reference volume scale
912 : * @return reference volume scale
913 : */
914 : double referenceVolume();
915 :
916 : // useful concept from Stack Overflow for templating MPI calls
917 : template <typename T>
918 : MPI_Datatype resolveType();
919 :
920 : /**
921 : * Helper function for MPI_Allgatherv of results in NekRS
922 : * @param[in] base_counts once multiplied by 'multiplier', the number of counts on each rank
923 : * @param[in] input rank-local data
924 : * @param[out] output collected result
925 : * @param[in] multiplier constant multiplier to set on each count indicator
926 : */
927 : template <typename T>
928 : void
929 93115 : allgatherv(const std::vector<int> & base_counts, const T * input, T * output, const int multiplier)
930 : {
931 93115 : int * recvCounts = (int *)calloc(commSize(), sizeof(int));
932 93115 : int * displacement = (int *)calloc(commSize(), sizeof(int));
933 93115 : displacementAndCounts(base_counts, recvCounts, displacement, multiplier);
934 :
935 93115 : MPI_Allgatherv(input,
936 : recvCounts[commRank()],
937 : resolveType<T>(),
938 : output,
939 : (const int *)recvCounts,
940 : (const int *)displacement,
941 : resolveType<T>(),
942 : platform->comm.mpiComm);
943 :
944 93115 : free(recvCounts);
945 93115 : free(displacement);
946 93115 : }
947 :
948 : } // end namespace nekrs
|
