23from CCZ4Solver
import CCZ4Solver_FV_GlobalAdaptiveTimeStepWithEnclaveTasking
24from CCZ4Solver
import CCZ4Solver_FD4_GlobalAdaptiveTimeStepWithEnclaveTasking
25from CCZ4Solver
import (
26 CCZ4Solver_FD4_SecondOrderFormulation_GlobalAdaptiveTimeStepWithEnclaveTasking,
29from CCZ4Solver
import CCZ4Solver_FV_GlobalAdaptiveTimeStep
30from CCZ4Solver
import CCZ4Solver_FD4_GlobalAdaptiveTimeStep
40 NONE =
"::peano4::utils::LoopPlacement::Serial"
41 PARALLEL_FOR =
"::peano4::utils::LoopPlacement::Nested"
42 SUBTASKS =
"::peano4::utils::LoopPlacement::SpreadOut"
48 Construct the Finite Volume (limiter) scheme
50 We assume that the underlying Finite Differences scheme has a patch
51 size of 6x6x6. To make the Finite Volume scheme's time stepping (and
52 accuracy) match this patch size, we have to employ a 16 times finer
55 It is interesting to see that the limiter does not really have a min
56 and max mesh size. The point is that the higher order solver dictates
57 the mesh structure, and we then follow this structure with the
66 amend_priorities: bool,
67 parallelisation_of_kernels: KernelParallelisation,
74 Pass in the patch size of the FD4 scheme or, if you are using RKDG,
75 hand in the number 1. The Finite Volume patch then will be 16 times
79 SubPatchSize =
int(patch_size * 4 * 4)
83 patch_size=SubPatchSize,
84 min_volume_h=1.0 / 65536,
86 pde_terms_without_state=
True,
91#include "toolbox/blockstructured/Restriction.h"
92#include "toolbox/blockstructured/Interpolation.h"
106 initial_conditions=exahype2.solvers.PDETerms.User_Defined_Implementation,
107 refinement_criterion=exahype2.solvers.PDETerms.Empty_Implementation,
108 boundary_conditions=exahype2.solvers.PDETerms.Empty_Implementation,
111 if parallelisation_of_kernels == KernelParallelisation.SUBTASKS:
123 compute_max_eigenvalue_of_next_time_step=
True,
124 solver_variant=exahype2.solvers.fv.rusanov.kernels.SolverVariant.Multicore,
125 kernel_variant=exahype2.solvers.fv.rusanov.kernels.KernelVariant.PatchWiseAoSHeap,
134 Mask out exterior cells
140 +
" and repositories::instanceOf"
142 +
".isCellOverlappingWithBHImpactArea(marker) )"
149 +
" and repositories::instanceOf"
151 +
".isCellOverlappingWithBHImpactArea(marker) )"
158 +
" and repositories::instanceOf"
160 +
".isCellOverlappingWithBHImpactArea(marker) )"
164 return "{} and repositories::instanceOf{}.isOneAdjacentCellOverlappingWithBHImpactArea(marker)".format(
170 return "({} and repositories::instanceOf{}.areBothAdjacentCellsOverlappingWithBHImpactArea(marker))".format(
175 return "({} and repositories::instanceOf{}.areBothAdjacentCellsOverlappingWithBHImpactArea(marker))".format(
182 Not really a lot of things to do here. The only exception that is
183 really important is that we have to ensure that we only solve stuff
184 inside the local domain of the FV. By default, ExaHyPE 2 solves the
185 PDE everywhere. If data is not stored persistently or loaded from
186 the persistent stacks, it still solves, as it then would assume that
187 such data arises from dynamic AMR. In this particular case, we have
188 to really mask out certain subdomains.
190 It is not just a nice optimisation to do so. It is absolutely key,
191 as the application of the compute kernel on garbage would mean that
192 we end up with invalid eigenvalues.
208 4th order Finite Differences solver with a limiter
210 The FD scheme is our primary solver, i.e. the one we wanna use (almost)
211 everywhere. Hence, we have to be sure that we use the right boundary
212 conditions. In our case, that should be Sommerfeld ones. As we work
213 with a limited FD4 solver, we have to ensure that we get the
214 @ref benchmarks_exahype2_ccz4_single_black_hole "solver coupling" right.
216 The primary solver is a plain CCZ4 FD4 solver with enclave tasking. There
217 are some modification though:
219 1. We ensure that we use Sommerfeld boundary conditions. This means that
220 we have to replace the whole generic boundary treatment with a bespoke
222 2. The coupling has to be injected. We model the coupling as separate
223 (enclave) tasks to avoid that we totally serialise the solver steps
224 whenever we encounter a coupling.
234 parallelisation_of_interpolation: KernelParallelisation,
235 parallelisation_of_kernels: KernelParallelisation,
236 interpolation_method,
242 CCZ4Solver_FD4_GlobalAdaptiveTimeStepWithEnclaveTasking.__init__(
245 patch_size=patch_size,
247 min_meshcell_h=min_cell_size / patch_size,
248 max_meshcell_h=max_cell_size / patch_size,
249 pde_terms_without_state=
False,
253#include "toolbox/blockstructured/Projection.h"
254#include "toolbox/blockstructured/Restriction.h"
255#include "toolbox/blockstructured/Interpolation.h"
256#include "toolbox/blockstructured/Copy.h"
261#include "matrixdata/Interpolator_"""
264 + str(patch_size * 16)
271 initial_conditions=exahype2.solvers.PDETerms.User_Defined_Implementation,
272 refinement_criterion=exahype2.solvers.PDETerms.User_Defined_Implementation,
273 boundary_conditions=exahype2.solvers.PDETerms.None_Implementation,
277 if parallelisation_of_kernels == KernelParallelisation.SUBTASKS:
279 "WARNING: Unfortunately, we have not yet written parallelised kernels for FD4"
287 patch_size, patch_size * 16, 3, 1
289 solver_matrix_interpolator.generateData()
294 Tailor action set behaviour
296 We first make a few additional cells skeleton cells. The rationale
297 behind additional skeletons is given in the @ref benchmarks_exahype2_ccz4_single_black_hole "generic overview".
298 Given the first remark there on FD4-FV coupling, one would be tempted
301 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
302 self._action_set_update_cell.additional_skeleton_guard = " " "(
303 repositories::instanceOf" " " + self._name_without_FD4_extension + " " "_FV.isCellOverlappingWithBHImpactArea(marker)
305 not repositories::instanceOf" " " + self._name_without_FD4_extension + " " "_FV.areAllFaceConnectedCellsOverlappingWithBHImpactArea(marker)
308 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
310 Once we study the other items (notably the fourth), we see that it is
311 reasonable to make all the overlap region a skeleton within the FD4
319 double Qinf[59]={1.0, 0.0, 0.0, 1.0, 0.0, 1.0, //q0-5
320 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, //q6-11
321 0.0, 0.0, 0.0, 0.0, //q12-15
322 1.0, 0.0, 0.0, 0.0, //q16-19
323 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, //q20-25
324 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, //q26-34
325 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,
326 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, //q35-52
327 0.0, 1.0, 0.0, 0.0, 0.0, 0.0 //q53-58
328 }; //approximate background solution at infinity
329 ::exahype2::fd::applySommerfeldConditions(
331 const double * __restrict__ Q,
332 const tarch::la::Vector<Dimensions,double>& faceCentre,
333 const tarch::la::Vector<Dimensions,double>& gridCellH,
338 return repositories::{{SOLVER_INSTANCE}}.maxEigenvalue( Q, faceCentre, gridCellH, t, dt, normal );
341 double * __restrict__ Q,
342 const tarch::la::Vector<Dimensions,double>& faceCentre,
343 const tarch::la::Vector<Dimensions,double>& gridCellH
350 Q[0] = 1.0; Q[3] = 1.0; Q[5] = 1.0;
351 Q[16] = 0.95; Q[54] = 0.95;
355 {{FACE_METADATA_ACCESSOR}}.getOldTimeStamp(marker.getSelectedFaceNumber()<Dimensions ? 1 : 0),
356 repositories::{{SOLVER_INSTANCE}}.getMinTimeStepSize(),
357 {{NUMBER_OF_GRID_CELLS_PER_PATCH_PER_AXIS}},
359 {{NUMBER_OF_UNKNOWNS}},
360 {{NUMBER_OF_AUXILIARY_VARIABLES}},
361 marker.getSelectedFaceNumber(),
363 fineGridFace{{UNKNOWN_IDENTIFIER}}Old.value,
364 fineGridFace{{UNKNOWN_IDENTIFIER}}New.value
371 enclave_task_cell_label=
"fineGridCell"
374 compute_kernel_implementation=
"""
375const int sizeOfPatch = (repositories::instanceOf"""
377 +
"""_FV.NumberOfFiniteVolumesPerAxisPerPatch+2)
378 * (repositories::instanceOf"""
380 +
"""_FV.NumberOfFiniteVolumesPerAxisPerPatch+2)
381 * (repositories::instanceOf"""
383 +
"""_FV.NumberOfFiniteVolumesPerAxisPerPatch+2)
384 * (repositories::instanceOf"""
386 +
"""_FV.NumberOfUnknowns + repositories::instanceOf"""
388 +
"""_FV.NumberOfAuxiliaryVariables);
389const int sizeOfFace = 2
390 * (repositories::instanceOf"""
392 +
"""_FV.NumberOfFiniteVolumesPerAxisPerPatch+2)
393 * (repositories::instanceOf"""
395 +
"""_FV.NumberOfFiniteVolumesPerAxisPerPatch+2)
396 * (repositories::instanceOf"""
398 +
"""_FV.NumberOfUnknowns + repositories::instanceOf"""
400 +
"""_FV.NumberOfAuxiliaryVariables);
402double* interpolatedFVDataWithHalo = ::tarch::allocateMemory<double>(sizeOfPatch, ::tarch::MemoryLocation::Heap);
404bool faceIsReal0 = repositories::instanceOf"""
406 +
"""_FV.areBothAdjacentCellsOverlappingWithBHImpactArea(marker,0);
407bool faceIsReal1 = repositories::instanceOf"""
409 +
"""_FV.areBothAdjacentCellsOverlappingWithBHImpactArea(marker,1);
410bool faceIsReal2 = repositories::instanceOf"""
412 +
"""_FV.areBothAdjacentCellsOverlappingWithBHImpactArea(marker,2);
413bool faceIsReal3 = repositories::instanceOf"""
415 +
"""_FV.areBothAdjacentCellsOverlappingWithBHImpactArea(marker,3);
416bool faceIsReal4 = repositories::instanceOf"""
418 +
"""_FV.areBothAdjacentCellsOverlappingWithBHImpactArea(marker,4);
419bool faceIsReal5 = repositories::instanceOf"""
421 +
"""_FV.areBothAdjacentCellsOverlappingWithBHImpactArea(marker,5);
423double* realFace0 = faceIsReal0 ? fineGridFaces"""
425 +
"""_FVQNew(0).value : nullptr;
426double* realFace1 = faceIsReal1 ? fineGridFaces"""
428 +
"""_FVQNew(1).value : nullptr;
429double* realFace2 = faceIsReal2 ? fineGridFaces"""
431 +
"""_FVQNew(2).value : nullptr;
432double* realFace3 = faceIsReal3 ? fineGridFaces"""
434 +
"""_FVQNew(3).value : nullptr;
435double* realFace4 = faceIsReal4 ? fineGridFaces"""
437 +
"""_FVQNew(4).value : nullptr;
438double* realFace5 = faceIsReal5 ? fineGridFaces"""
440 +
"""_FVQNew(5).value : nullptr;
442double* dummyFace0 = faceIsReal0 ? nullptr : ::tarch::allocateMemory<double>(sizeOfFace, ::tarch::MemoryLocation::Heap);
443double* dummyFace1 = faceIsReal1 ? nullptr : ::tarch::allocateMemory<double>(sizeOfFace, ::tarch::MemoryLocation::Heap);
444double* dummyFace2 = faceIsReal2 ? nullptr : ::tarch::allocateMemory<double>(sizeOfFace, ::tarch::MemoryLocation::Heap);
445double* dummyFace3 = faceIsReal3 ? nullptr : ::tarch::allocateMemory<double>(sizeOfFace, ::tarch::MemoryLocation::Heap);
446double* dummyFace4 = faceIsReal4 ? nullptr : ::tarch::allocateMemory<double>(sizeOfFace, ::tarch::MemoryLocation::Heap);
447double* dummyFace5 = faceIsReal5 ? nullptr : ::tarch::allocateMemory<double>(sizeOfFace, ::tarch::MemoryLocation::Heap);
449::toolbox::blockstructured::interpolateCellDataAssociatedToVolumesIntoOverlappingCell_"""
452 repositories::instanceOf"""
454 +
"""_FD4.NumberOfGridCellsPerPatchPerAxis,
455 repositories::instanceOf"""
457 +
"""_FV.NumberOfFiniteVolumesPerAxisPerPatch,
460 repositories::instanceOf"""
462 +
"""_FV.NumberOfUnknowns + repositories::instanceOf"""
464 +
"""_FV.NumberOfAuxiliaryVariables,
467 """InterpolationMatrixData::Data,
468 InterpolationMatrixData::Indices,
469 InterpolationMatrixData::Indptr,"""
475 interpolatedFVDataWithHalo,
481::toolbox::blockstructured::projectPatchHaloOntoFaces(
482 repositories::instanceOf"""
484 +
"""_FV.NumberOfFiniteVolumesPerAxisPerPatch,
486 repositories::instanceOf"""
488 +
"""_FV.NumberOfUnknowns,
489 repositories::instanceOf"""
491 +
"""_FV.NumberOfAuxiliaryVariables,
492 interpolatedFVDataWithHalo,
493 faceIsReal0 ? realFace0 : dummyFace0,
494 faceIsReal1 ? realFace1 : dummyFace1,
495 faceIsReal2 ? realFace2 : dummyFace2,
496 faceIsReal3 ? realFace3 : dummyFace3,
497 faceIsReal4 ? realFace4 : dummyFace4,
498 faceIsReal5 ? realFace5 : dummyFace5
501::tarch::freeMemory(interpolatedFVDataWithHalo, tarch::MemoryLocation::Heap );
502if (dummyFace0!=nullptr) ::tarch::freeMemory(dummyFace0, tarch::MemoryLocation::Heap );
503if (dummyFace1!=nullptr) ::tarch::freeMemory(dummyFace1, tarch::MemoryLocation::Heap );
504if (dummyFace2!=nullptr) ::tarch::freeMemory(dummyFace2, tarch::MemoryLocation::Heap );
505if (dummyFace3!=nullptr) ::tarch::freeMemory(dummyFace3, tarch::MemoryLocation::Heap );
506if (dummyFace4!=nullptr) ::tarch::freeMemory(dummyFace4, tarch::MemoryLocation::Heap );
507if (dummyFace5!=nullptr) ::tarch::freeMemory(dummyFace5, tarch::MemoryLocation::Heap );
516 repositories::instanceOf"""
518 +
"""_FV.isCellOverlappingWithBHImpactArea(marker)
520 not repositories::instanceOf"""
522 +
"""_FV.areAllFaceConnectedCellsOverlappingWithBHImpactArea(marker)
524 not marker.willBeRefined()
532 repositories::instanceOf{0}_FV.isCellOverlappingWithBHImpactArea(marker)
534 not marker.willBeRefined()
536 logTraceIn( "touchCellFirstTime(...)-inject" );
538 repositories::instanceOf{0}_FV.incNumberOfPatches();
540 //if ( repositories::instanceOf{0}_FV.areAllFaceConnectedCellsOverlappingWithBHImpactArea(marker) ) {{
541 ::toolbox::blockstructured::restrictCellIntoOverlappingCell_inject(
542 repositories::instanceOf{0}_FV.NumberOfFiniteVolumesPerAxisPerPatch,
543 repositories::instanceOf{0}_FD4.NumberOfGridCellsPerPatchPerAxis,
544 repositories::instanceOf{0}_FV.NumberOfUnknowns + repositories::instanceOf{0}_FV.NumberOfAuxiliaryVariables,
545 fineGridCell{0}_FVQ.value,
546 fineGridCell{0}_FD4Q.value
550 ::toolbox::blockstructured::restrictCellIntoOverlappingCell_inject_and_average(
551 repositories::instanceOf{0}_FV.NumberOfFiniteVolumesPerAxisPerPatch,
552 repositories::instanceOf{0}_FD4.NumberOfGridCellsPerPatchPerAxis,
553 repositories::instanceOf{0}_FV.NumberOfUnknowns + repositories::instanceOf{0}_FV.NumberOfAuxiliaryVariables,
554 fineGridCell{0}_FVQ.value,
555 fineGridCell{0}_FD4Q.value
559 logTraceOut( "touchCellFirstTime(...)-inject" );
570 Construct 4th order Finite Differences solver without a limiter
583 CCZ4Solver_FD4_GlobalAdaptiveTimeStepWithEnclaveTasking.__init__(
586 patch_size=patch_size,
588 min_meshcell_h=min_cell_size / patch_size,
589 max_meshcell_h=max_cell_size / patch_size,
590 pde_terms_without_state=
False,
594#include "toolbox/blockstructured/Projection.h"
595#include "toolbox/blockstructured/Restriction.h"
596#include "toolbox/blockstructured/Interpolation.h"
597#include "toolbox/blockstructured/Copy.h"
604 initial_conditions=exahype2.solvers.PDETerms.User_Defined_Implementation,
605 refinement_criterion=exahype2.solvers.PDETerms.User_Defined_Implementation,
606 boundary_conditions=exahype2.solvers.PDETerms.None_Implementation,
615 Tailor action set behaviour
617 We first make a few additional cells skeleton cells. The rationale
618 behind additional skeletons is given in the @ref benchmarks_exahype2_ccz4_single_black_hole "generic overview".
619 Given the first remark there on FD4-FV coupling, one would be tempted
622 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
623 self._action_set_update_cell.additional_skeleton_guard = " " "(
624 repositories::instanceOf" " " + self._name_without_FD4_extension + " " "_FV.isCellOverlappingWithBHImpactArea(marker)
626 not repositories::instanceOf" " " + self._name_without_FD4_extension + " " "_FV.areAllFaceConnectedCellsOverlappingWithBHImpactArea(marker)
629 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
631 Once we study the other items (notably the fourth), we see that it is
632 reasonable to make all the overlap region a skeleton within the FD4
640 double Qinf[59]={1.0, 0.0, 0.0, 1.0, 0.0, 1.0, //q0-5
641 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, //q6-11
642 0.0, 0.0, 0.0, 0.0, //q12-15
643 1.0, 0.0, 0.0, 0.0, //q16-19
644 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, //q20-25
645 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, //q26-34
646 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,
647 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, //q35-52
648 0.0, 1.0, 0.0, 0.0, 0.0, 0.0 //q53-58
649 }; //approximate background solution at infinity
650 ::exahype2::fd::applySommerfeldConditions(
652 const double * __restrict__ Q,
653 const tarch::la::Vector<Dimensions,double>& faceCentre,
654 const tarch::la::Vector<Dimensions,double>& gridCellH,
659 return repositories::{{SOLVER_INSTANCE}}.maxEigenvalue( Q, faceCentre, gridCellH, t, dt, normal );
662 double * __restrict__ Q,
663 const tarch::la::Vector<Dimensions,double>& faceCentre,
664 const tarch::la::Vector<Dimensions,double>& gridCellH
671 Q[0] = 1.0; Q[3] = 1.0; Q[5] = 1.0;
672 Q[16] = 0.95; Q[54] = 0.95;
676 {{FACE_METADATA_ACCESSOR}}.getOldTimeStamp(marker.getSelectedFaceNumber()<Dimensions ? 1 : 0),
677 repositories::{{SOLVER_INSTANCE}}.getMinTimeStepSize(),
678 {{NUMBER_OF_GRID_CELLS_PER_PATCH_PER_AXIS}},
680 {{NUMBER_OF_UNKNOWNS}},
681 {{NUMBER_OF_AUXILIARY_VARIABLES}},
682 marker.getSelectedFaceNumber(),
684 fineGridFace{{UNKNOWN_IDENTIFIER}}Old.value,
685 fineGridFace{{UNKNOWN_IDENTIFIER}}New.value
694 A finite volume solver
696 This solver is not appropriate to simulate black holes as a stand-alone
697 solver, as it is way too diffusive. If you use it without another scheme,
698 you typically see the black hole disappear after a brief period. So we
699 have it in here merely for performance tests.
712 Construct the Finite Volume solver
714 @param patch_size: Integer
715 Defines how big the individual patches are. If you pass in 10, each
716 Finite Volume patch will have the dimensions 10x10x10.
717 @param min_cell_size: Float
718 This parameter refers to the cell size, i.e. the size of a whole
719 patch. We use this one here, to make the signature the same as for
720 the FD and DG solver variants. The superclass constructor argues
721 over finite volume sizes, and we hence have to recalibrate this
722 parameter with patch_size.
728 patch_size=patch_size,
729 min_volume_h=min_cell_size / patch_size,
730 max_volume_h=max_cell_size / patch_size,
731 pde_terms_without_state=
True,
735#include "toolbox/blockstructured/Projection.h"
736#include "toolbox/blockstructured/Restriction.h"
737#include "toolbox/blockstructured/Interpolation.h"
738#include "toolbox/blockstructured/Copy.h"
745 initial_conditions=exahype2.solvers.PDETerms.User_Defined_Implementation,
746 refinement_criterion=exahype2.solvers.PDETerms.Empty_Implementation,
747 boundary_conditions=exahype2.solvers.PDETerms.Empty_Implementation,
756 Tailor action set behaviour
758 We first make a few additional cells skeleton cells. The rationale
759 behind additional skeletons is given in the @ref benchmarks_exahype2_ccz4_single_black_hole "generic overview".
760 Given the first remark there on FD4-FV coupling, one would be tempted
763 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
764 self._action_set_update_cell.additional_skeleton_guard = " " "(
765 repositories::instanceOf" " " + self._name_without_FD4_extension + " " "_FV.isCellOverlappingWithBHImpactArea(marker)
767 not repositories::instanceOf" " " + self._name_without_FD4_extension + " " "_FV.areAllFaceConnectedCellsOverlappingWithBHImpactArea(marker)
770 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
772 Once we study the other items (notably the fourth), we see that it is
773 reasonable to make all the overlap region a skeleton within the FD4
783 Update preconfigured solver parameters
785 The default parameters of CCZ4 are tailored towards gauge waves or similar.
786 For the single black hole, two parameters have to be changed. That's bs and
787 sk which both have to be 1.0.
790 solver.double_constants[
"CCZ4bs"] = 1.0
791 solver.double_constants[
"CCZ4sk"] = 1.0
add_all_solver_constants(self)
Add domain-specific constants.
CCZ4 solver using fourth-order finite differences and global adaptive time stepping incl enclave task...
CCZ4 solver using finite volumes and global adaptive time stepping incl enclave tasking.
4th order Finite Differences solver with a limiter
_action_set_postprocess_solution
__init__(self, name, patch_size, min_cell_size, max_cell_size, KernelParallelisation parallelisation_of_interpolation, KernelParallelisation parallelisation_of_kernels, interpolation_method)
Constructor.
create_action_sets(self)
Tailor action set behaviour.
_name_without_FD4_extension
_parallelisation_of_interpolation
_action_set_preprocess_solution
Construct 4th order Finite Differences solver without a limiter.
_name_without_FD4_extension
create_action_sets(self)
Tailor action set behaviour.
__init__(self, name, patch_size, min_cell_size, max_cell_size)
Constructor.
__init__(self, name, patch_size, min_cell_size, max_cell_size)
Construct the Finite Volume solver.
create_action_sets(self)
Tailor action set behaviour.
Construct the Finite Volume (limiter) scheme.
create_action_sets(self)
Not really a lot of things to do here.
_source_term_implementation
_fused_compute_kernel_call_cpu
_provide_face_data_to_compute_kernels_default_guard(self)
_store_cell_data_default_guard(self)
Mask out exterior cells.
_provide_cell_data_to_compute_kernels_default_guard(self)
Default logic when to create cell data or not.
_store_face_data_default_guard(self)
Extend the guard via ands only.
_load_cell_data_default_guard(self)
Extend the guard via ands only.
__init__(self, name, patch_size, bool amend_priorities, KernelParallelisation parallelisation_of_kernels)
Construct the limiter.
_load_face_data_default_guard(self)
Extend the guard via ands only.
_action_set_merge_enclave_task_outcome
set_implementation(self, boundary_conditions, refinement_criterion, initial_conditions, memory_location, use_split_loop, additional_action_set_includes, additional_user_includes)
If you pass in User_Defined, then the generator will create C++ stubs that you have to befill manuall...
_user_action_set_includes
_source_term_implementation
_source_term_implementation
set_implementation(self, flux=None, ncp=None, eigenvalues=None, boundary_conditions=None, refinement_criterion=None, initial_conditions=None, source_term=None, memory_location=None, use_split_loop=False, additional_action_set_includes="", additional_user_includes="")
If you pass in User_Defined, then the generator will create C++ stubs that you have to befill manuall...
_user_action_set_includes
_action_set_postprocess_solution
_action_set_preprocess_solution
_action_set_handle_boundary
set_implementation(self, flux, ncp, source_term, eigenvalues, boundary_conditions, refinement_criterion, initial_conditions, memory_location, additional_action_set_includes, additional_user_includes)
If you pass in User_Defined, then the generator will create C++ stubs that you have to befill manuall...
Preprocess the solution and reconstruct halo for this.
set_implementation(self, flux=None, ncp=None, source_term=None, eigenvalues=None, boundary_conditions=None, refinement_criterion=None, initial_conditions=None, memory_location=None, additional_action_set_includes="", additional_user_includes="", KOSigma=None)
If you pass in User_Defined, then the generator will create C++ stubs that you have to befill manuall...
update_solver_parameters_for_single_black_hole(solver)
Update preconfigured solver parameters.
create_compute_Riemann_kernel_for_Rusanov(flux_implementation, ncp_implementation, source_implementation, compute_max_eigenvalue_of_next_time_step, SolverVariant solver_variant, KernelVariant kernel_variant)
Return only the unqualified function call, i.e., without any namespaces.