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