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