4 #define VELOCITY_SPACE_DIMENSION 2
13 #include <Kokkos_Core.hpp>
16 #include <petscversion.h>
18 #include <petscdmplex.h>
19 #if PETSC_VERSION_GE(3,14,0)
21 #include "petscdmswarm.h"
36 Kokkos::View<DMObjects*, Kokkos::LayoutRight, HostType>
obj;
41 Kokkos::resize(
obj, 0);
46 for (
int id = 0;
id <
dm_len;
id++){
47 DMDestroy(&
obj(
id).dm);
48 MatDestroy(&
obj(
id).MM);
57 for (
int id = 0;
id <
dm_len;
id++){
58 obj(
id).dm = original.
obj(
id).dm;
59 obj(
id).MM = original.
obj(
id).MM;
64 if (
this == &original){
70 for (
int id = 0;
id <
dm_len;
id++){
71 obj(
id).dm = original.
obj(
id).dm;
72 obj(
id).MM = original.
obj(
id).MM;
95 PetscInt order,
const PetscInt Nspecies,
96 const PetscReal vPara_min,
const PetscReal vPara_max,
97 const PetscInt num_vPara,
const PetscReal vPerp_min,
98 const PetscReal vPerp_max,
const PetscInt num_vPerp)
100 #if PETSC_VERSION_LT(3,14,0)
102 printf(
"\n\n\nERROR: initialize_pseudo_inverse_interpolation (should never reach here)\n");
112 PetscReal lo[2], hi[2];
113 cells[0] = num_vPara - 1;
114 cells[1] = num_vPerp - 1;
120 PetscFunctionBeginUser;
122 std::string order_str = std::to_string(order);
123 ierr = PetscOptionsSetValue(NULL,
"-petscspace_degree", order_str.c_str());CHKERRQ(ierr);
128 DMBoundaryType periodicity[2] = {DM_BOUNDARY_NONE, DM_BOUNDARY_NONE};
129 PetscBool interpolate = PETSC_TRUE;
131 for (
int id = 0;
id <
dm_len;
id++){
132 ierr = DMPlexCreateBoxMesh(PETSC_COMM_SELF, dim, PETSC_FALSE, cells, lo, hi, periodicity, interpolate, &
obj(
id).dm);CHKERRQ(ierr);
133 ierr = PetscObjectSetName((PetscObject)
obj(
id).dm,
"Potential Grid");CHKERRQ(ierr);
134 if (
id == dm_len - 1) {ierr = DMViewFromOptions(
obj(
id).dm, NULL,
"-dm_view");CHKERRQ(ierr);}
135 ierr = DMSetFromOptions(
obj(
id).dm);CHKERRQ(ierr);
136 ierr = PetscFECreateDefault(PETSC_COMM_SELF, dim, Nspecies, PETSC_FALSE,
"", PETSC_DECIDE, &fe);CHKERRQ(ierr);
137 ierr = PetscFESetFromOptions(fe);CHKERRQ(ierr);
138 ierr = PetscObjectSetName((PetscObject)fe,
"fe");CHKERRQ(ierr);
139 ierr = DMSetField(
obj(
id).dm, field, NULL, (PetscObject)fe);CHKERRQ(ierr);
140 ierr = DMCreateDS(
obj(
id).dm);CHKERRQ(ierr);
141 ierr = PetscFEDestroy(&fe);CHKERRQ(ierr);
143 ierr = DMCreateMassMatrix(
obj(
id).dm,
obj(
id).dm, &
obj(
id).MM);CHKERRQ(ierr);
146 PetscFunctionReturn(0);
159 PetscFunctionBeginUser;
161 for (
int id = 0;
id <
dm_len;
id++){
162 if(
obj(
id).dm != NULL){
163 ierr = DMDestroy(&
obj(
id).dm);CHKERRQ(ierr);
166 if(
obj(
id).MM != NULL){
167 ierr = MatDestroy(&
obj(
id).MM);CHKERRQ(ierr);
172 PetscFunctionReturn(0);
PetscErrorCode initialize_pseudo_inverse_interpolation(int nthreads, PetscInt order, const PetscInt Nspecies, const PetscReal vPara_min, const PetscReal vPara_max, const PetscInt num_vPara, const PetscReal vPerp_min, const PetscReal vPerp_max, const PetscInt num_vPerp)
Definition: DM_wrapper.hpp:94
DMWrapper()
Definition: DM_wrapper.hpp:38
DM dm
PETSc DM mesh object representing the velocity grid when pseudo-inverse is used.
Definition: DM_wrapper.hpp:28
#define VELOCITY_SPACE_DIMENSION
Definition: DM_wrapper.hpp:4
bool is_owner
Definition: DM_wrapper.hpp:34
~DMWrapper()
Definition: DM_wrapper.hpp:44
Definition: DM_wrapper.hpp:33
DMWrapper(const DMWrapper &original)
Definition: DM_wrapper.hpp:53
DMWrapper & operator=(const DMWrapper &original)
Definition: DM_wrapper.hpp:63
Kokkos::View< DMObjects *, Kokkos::LayoutRight, HostType > obj
Stores #threads DMObjects.
Definition: DM_wrapper.hpp:36
PetscErrorCode destroy_DM()
Definition: DM_wrapper.hpp:157
int dm_len
Number of PETSc mesh objects to store.
Definition: DM_wrapper.hpp:35
Definition: DM_wrapper.hpp:27
Mat MM
PETSc mass matrix object for the velocity grid when pseudo-inverse is used.
Definition: DM_wrapper.hpp:29