Actual source code: petsclandau.h
1: #if !defined(PETSCLANDAU_H)
2: #define PETSCLANDAU_H
4: #include <petscdmplex.h>
5: #include <petscts.h>
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 *);
15: /* the Fokker-Planck-Landau context */
16: #if !defined(LANDAU_DIM)
17: #define LANDAU_DIM 2
18: #endif
20: #if !defined(LANDAU_MAX_SPECIES)
21: #if LANDAU_DIM==2
22: #define LANDAU_MAX_SPECIES 10
23: #else
24: #define LANDAU_MAX_SPECIES 3
25: #endif
26: #endif
28: #if !defined(LANDAU_MAX_Q)
29: #if defined(LANDAU_MAX_NQ)
30: #error"LANDAU_MAX_NQ but not LANDAU_MAX_Q. Use -DLANDAU_MAX_Q=4 for Q3 elements"
31: #endif
32: #if LANDAU_DIM==2
33: #define LANDAU_MAX_Q 5
34: #else
35: #define LANDAU_MAX_Q 3
36: #endif
37: #else
38: #undef LANDAU_MAX_NQ
39: #endif
41: #if LANDAU_DIM==2
42: #define LANDAU_MAX_Q_FACE LANDAU_MAX_Q
43: #define LANDAU_MAX_NQ (LANDAU_MAX_Q*LANDAU_MAX_Q)
44: #else
45: #define LANDAU_MAX_Q_FACE (LANDAU_MAX_Q*LANDAU_MAX_Q)
46: #define LANDAU_MAX_NQ (LANDAU_MAX_Q*LANDAU_MAX_Q*LANDAU_MAX_Q)
47: #endif
49: typedef enum {LANDAU_CUDA, LANDAU_KOKKOS, LANDAU_CPU} LandauDeviceType;
50: typedef struct {
51: PetscBool interpolate; /* Generate intermediate mesh elements */
52: PetscBool gpu_assembly;
53: PetscFE fe[LANDAU_MAX_SPECIES];
54: /* geometry */
55: PetscReal i_radius;
56: PetscReal e_radius;
57: PetscInt num_sections;
58: PetscReal radius;
59: PetscReal re_radius; /* radius of refinement along v_perp=0, z>0 */
60: PetscReal vperp0_radius1; /* radius of refinement along v_perp=0 */
61: PetscReal vperp0_radius2; /* radius of refinement along v_perp=0 after origin AMR refinement */
62: PetscBool sphere;
63: PetscBool inflate;
64: PetscInt numRERefine; /* refinement along v_perp=0, z > 0 */
65: PetscInt nZRefine1; /* origin refinement after v_perp=0 refinement */
66: PetscInt nZRefine2; /* origin refinement after origin AMR refinement */
67: PetscInt maxRefIts; /* normal AMR - refine from origin */
68: PetscInt postAMRRefine; /* uniform refinement of AMR */
69: /* discretization - AMR */
70: PetscErrorCode (*errorIndicator)(PetscInt, PetscReal, PetscReal [], PetscInt, const PetscInt[], const PetscScalar[], const PetscScalar[], PetscReal *, void *);
71: PetscReal refineTol[LANDAU_MAX_SPECIES];
72: PetscReal coarsenTol[LANDAU_MAX_SPECIES];
73: /* physics */
74: PetscReal thermal_temps[LANDAU_MAX_SPECIES];
75: PetscReal masses[LANDAU_MAX_SPECIES]; /* mass of each species */
76: PetscReal charges[LANDAU_MAX_SPECIES]; /* charge of each species */
77: PetscReal n[LANDAU_MAX_SPECIES]; /* number density of each species */
78: PetscReal m_0; /* reference mass */
79: PetscReal v_0; /* reference velocity */
80: PetscReal n_0; /* reference number density */
81: PetscReal t_0; /* reference time */
82: PetscReal Ez;
83: PetscReal epsilon0;
84: PetscReal k;
85: PetscReal lnLam;
86: PetscReal electronShift; /* for tests */
87: PetscInt num_species;
88: /* diagnostics */
89: PetscInt verbose;
90: PetscLogEvent events[20];
91: DM dmv;
92: /* cache */
93: Mat J;
94: Mat M;
95: Vec X;
96: PetscReal normJ; /* used to see if function changed */
97: /* derived type */
98: void *data;
99: PetscBool aux_bool; /* helper */
100: /* computing */
101: LandauDeviceType deviceType;
102: PetscInt subThreadBlockSize;
103: PetscInt numConcurrency; /* number of SMs in Cuda to use */
104: MPI_Comm comm; /* global communicator to use for errors and diagnostics */
105: } LandauCtx;
107: typedef int LandauIdx;
108: typedef struct {
109: PetscReal scale;
110: LandauIdx gid; // Lanadu matrix index (<10,000)
111: } pointInterpolationP4est;
112: typedef struct _lP4estVertexMaps {
113: LandauIdx (*gIdx)[LANDAU_MAX_SPECIES][LANDAU_MAX_NQ]; // #elems * LANDAU_MAX_NQ (spoof for max , Nb) on device,
114: LandauIdx num_elements;
115: LandauIdx num_reduced;
116: LandauIdx num_face; // (Q or Q^2 for 3D)
117: LandauDeviceType deviceType;
118: PetscInt Nf;
119: PetscInt Nq;
120: pointInterpolationP4est (*c_maps)[LANDAU_MAX_Q_FACE];
121: struct _lP4estVertexMaps*data;
122: void *vp1,*vp2,*vp3;
123: } P4estVertexMaps;
125: typedef PetscReal LandauIPReal;
126: typedef struct {
127: LandauIPReal *w;
128: LandauIPReal *x;
129: LandauIPReal *y;
130: LandauIPReal *z;
131: LandauIPReal *coefs;
132: int dim_,ns_,nip_;
133: } LandauIPData;
135: PETSC_EXTERN PetscErrorCode LandauCreateColoring(Mat, DM, PetscContainer *);
136: PETSC_EXTERN int LandauGetIPDataSize(const LandauIPData *const);
137: #if defined(PETSC_HAVE_CUDA)
138: PETSC_EXTERN PetscErrorCode LandauCUDAJacobian(DM, const PetscInt, const PetscReal [], const PetscReal [], const PetscReal[], const PetscReal[],
139: const LandauIPData *const, const PetscReal [], PetscReal *, PetscReal, const PetscLogEvent[], Mat);
140: PETSC_EXTERN PetscErrorCode LandauCUDACreateMatMaps(P4estVertexMaps *, pointInterpolationP4est (*)[LANDAU_MAX_Q_FACE], PetscInt, PetscInt);
141: PETSC_EXTERN PetscErrorCode LandauCUDADestroyMatMaps(P4estVertexMaps *);
143: #endif
144: #if defined(PETSC_HAVE_KOKKOS)
145: PETSC_EXTERN PetscErrorCode LandauKokkosJacobian(DM, const PetscInt, PetscReal*, PetscReal*, PetscReal*, PetscReal*,
146: const LandauIPData* const, PetscReal*, const PetscInt, PetscReal*, PetscReal,
147: const PetscLogEvent*, Mat);
148: PETSC_EXTERN PetscErrorCode LandauKokkosCreateMatMaps(P4estVertexMaps *, pointInterpolationP4est (*)[LANDAU_MAX_Q_FACE], PetscInt, PetscInt);
149: PETSC_EXTERN PetscErrorCode LandauKokkosDestroyMatMaps(P4estVertexMaps *);
150: #endif
152: #endif /* PETSCLANDAU_H */