1 #if !defined(PETSCLANDAU_H) 2 #define PETSCLANDAU_H 3 4 #include <petscdmplex.h> /*I "petscdmplex.h" I*/ 5 #include <petscts.h> 6 7 PETSC_EXTERN PetscErrorCode LandauPrintNorms(Vec, PetscInt); 8 PETSC_EXTERN PetscErrorCode LandauCreateVelocitySpace(MPI_Comm,PetscInt,const char[],Vec*,Mat*,DM*); 9 PETSC_EXTERN PetscErrorCode LandauDestroyVelocitySpace(DM*); 10 PETSC_EXTERN PetscErrorCode LandauAddMaxwellians(DM, Vec, PetscReal, PetscReal[], PetscReal[], void *); 11 PETSC_EXTERN PetscErrorCode LandauCreateMassMatrix(DM dm, Mat *Amat); 12 PETSC_EXTERN PetscErrorCode LandauIFunction(TS, PetscReal,Vec,Vec,Vec,void *); 13 PETSC_EXTERN PetscErrorCode LandauIJacobian(TS, PetscReal,Vec,Vec,PetscReal,Mat,Mat,void *); 14 15 /* the Fokker-Planck-Landau context */ 16 #if !defined(LANDAU_DIM) 17 #define LANDAU_DIM 2 18 #endif 19 #if !defined(LANDAU_MAX_SPECIES) 20 #if LANDAU_DIM==2 21 #define LANDAU_MAX_SPECIES 10 22 #else 23 #define LANDAU_MAX_SPECIES 3 24 #endif 25 #endif 26 #if !defined(LANDAU_MAX_NQ) 27 #if LANDAU_DIM==2 28 #define LANDAU_MAX_NQ 25 29 #else 30 #define LANDAU_MAX_NQ 27 31 #endif 32 #endif 33 typedef enum {LANDAU_CUDA, LANDAU_KOKKOS, LANDAU_CPU} LandauDeviceType; 34 typedef struct { 35 PetscBool interpolate; /* Generate intermediate mesh elements */ 36 PetscBool simplex; 37 PetscFE fe[LANDAU_MAX_SPECIES]; 38 /* geometry */ 39 PetscReal i_radius; 40 PetscReal e_radius; 41 PetscInt num_sections; 42 PetscReal radius; 43 PetscReal re_radius; /* radius of refinement along v_perp=0, z>0 */ 44 PetscReal vperp0_radius1; /* radius of refinement along v_perp=0 */ 45 PetscReal vperp0_radius2; /* radius of refinement along v_perp=0 after origin AMR refinement */ 46 PetscBool sphere; 47 PetscBool inflate; 48 PetscInt numRERefine; /* refinement along v_perp=0, z > 0 */ 49 PetscInt nZRefine1; /* origin refinement after v_perp=0 refinement */ 50 PetscInt nZRefine2; /* origin refinement after origin AMR refinement */ 51 PetscInt maxRefIts; /* normal AMR - refine from origin */ 52 PetscInt postAMRRefine; /* uniform refinement of AMR */ 53 /* discretization - AMR */ 54 PetscErrorCode (*errorIndicator)(PetscInt, PetscReal, PetscReal [], PetscInt, const PetscInt[], const PetscScalar[], const PetscScalar[], PetscReal *, void *); 55 PetscReal refineTol[LANDAU_MAX_SPECIES]; 56 PetscReal coarsenTol[LANDAU_MAX_SPECIES]; 57 /* physics */ 58 PetscReal thermal_temps[LANDAU_MAX_SPECIES]; 59 PetscReal masses[LANDAU_MAX_SPECIES]; /* mass of each species */ 60 PetscReal charges[LANDAU_MAX_SPECIES]; /* charge of each species */ 61 PetscReal n[LANDAU_MAX_SPECIES]; /* number density of each species */ 62 PetscReal m_0; /* reference mass */ 63 PetscReal v_0; /* reference velocity */ 64 PetscReal n_0; /* reference number density */ 65 PetscReal t_0; /* reference time */ 66 PetscReal Ez; 67 PetscReal epsilon0; 68 PetscReal k; 69 PetscReal lnLam; 70 PetscReal electronShift; /* for tests */ 71 PetscInt num_species; 72 /* diagnostics */ 73 PetscInt verbose; 74 PetscLogEvent events[20]; 75 DM dmv; 76 /* cache */ 77 Mat J; 78 Mat M; 79 Vec X; 80 PetscReal normJ; /* used to see if function changed */ 81 /* derived type */ 82 void *data; 83 PetscBool aux_bool; /* helper */ 84 /* computing */ 85 LandauDeviceType deviceType; 86 PetscInt subThreadBlockSize; 87 } LandauCtx; 88 89 typedef PetscReal LandauIPReal; 90 typedef struct { 91 LandauIPReal *w_data; 92 LandauIPReal *x; 93 LandauIPReal *y; 94 LandauIPReal *z; 95 LandauIPReal *f; 96 LandauIPReal *dfx; 97 LandauIPReal *dfy; 98 LandauIPReal *dfz; 99 int dim_,ns_,nip_; 100 } LandauIPData; 101 102 PETSC_EXTERN PetscErrorCode LandauAssembleOpenMP(PetscInt cStart, PetscInt cEnd, PetscInt totDim, DM plex, PetscSection section, PetscSection globalSection, Mat JacP, PetscScalar elemMats[], PetscContainer container); 103 PETSC_EXTERN PetscErrorCode LandauCreateColoring(Mat, DM, PetscContainer *); 104 PETSC_EXTERN PetscErrorCode LandauFormJacobian_Internal(Vec, Mat, const PetscInt, void *); 105 PETSC_EXTERN int LandauGetIPDataSize(const LandauIPData *const); 106 #if defined(PETSC_HAVE_CUDA) 107 PETSC_EXTERN PetscErrorCode LandauCUDAJacobian(DM, const PetscInt, const PetscReal [], const PetscReal [], const PetscReal[], const PetscReal[], 108 const LandauIPData *const, const PetscReal [],const PetscInt, const PetscLogEvent[], Mat); 109 #endif 110 #if defined(PETSC_HAVE_KOKKOS) 111 /* TODO: this won't work if PETSc is built with C++ */ 112 #if !defined(__cplusplus) 113 PETSC_EXTERN PetscErrorCode LandauKokkosJacobian(DM, const PetscInt, const PetscReal [], const PetscReal [], const PetscReal[], const PetscReal[], 114 const LandauIPData *const, const PetscReal [],const PetscInt, const PetscLogEvent[], Mat); 115 #endif 116 #endif 117 118 #endif /* PETSCLANDAU_H */ 119