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[], PetscInt, 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: #define LANDAU_MAX_GRIDS 3
24: #else
25: #define LANDAU_MAX_SPECIES 3
26: #define LANDAU_MAX_GRIDS 2
27: #endif
28: #else
29: #define LANDAU_MAX_GRIDS 3
30: #endif
32: #if !defined(LANDAU_MAX_Q)
33: #if defined(LANDAU_MAX_NQ)
34: #error"LANDAU_MAX_NQ but not LANDAU_MAX_Q. Use -DLANDAU_MAX_Q=4 for Q3 elements"
35: #endif
36: #if LANDAU_DIM==2
37: #define LANDAU_MAX_Q 5
38: #else
39: #define LANDAU_MAX_Q 3
40: #endif
41: #else
42: #undef LANDAU_MAX_NQ
43: #endif
45: #if LANDAU_DIM==2
46: #define LANDAU_MAX_Q_FACE LANDAU_MAX_Q
47: #define LANDAU_MAX_NQ (LANDAU_MAX_Q*LANDAU_MAX_Q)
48: #else
49: #define LANDAU_MAX_Q_FACE (LANDAU_MAX_Q*LANDAU_MAX_Q)
50: #define LANDAU_MAX_NQ (LANDAU_MAX_Q*LANDAU_MAX_Q*LANDAU_MAX_Q)
51: #endif
53: typedef enum {LANDAU_CUDA, LANDAU_KOKKOS, LANDAU_CPU} LandauDeviceType;
55: // static data
56: typedef struct {
57: void *invJ; // nip*dim*dim
58: void *D; // nq*nb*dim
59: void *B; // nq*nb
60: void *alpha; // ns
61: void *beta; // ns
62: void *invMass; // ns
63: void *w; // nip
64: void *x; // nip
65: void *y; // nip
66: void *z; // nip
67: void *Eq_m; // ns - dynamic
68: void *f; // nip*Nf - dynamic (IP)
69: void *dfdx; // nip*Nf - dynamic (IP)
70: void *dfdy; // nip*Nf - dynamic (IP)
71: void *dfdz; // nip*Nf - dynamic (IP)
72: int dim_,ns_,nip_,nq_,nb_;
73: void *NCells;
74: void *species_offset;
75: void *mat_offset;
76: void *elem_offset;
77: void *ip_offset;
78: void *ipf_offset;
79: void *ipfdf_data;
80: void *maps;
81: } LandauStaticData;
83: typedef struct {
84: PetscBool interpolate; /* Generate intermediate mesh elements */
85: PetscBool gpu_assembly;
86: MPI_Comm comm; /* global communicator to use for errors and diagnostics */
87: double times[1];
88: PetscBool initialized;
89: PetscBool use_matrix_mass;
90: /* FE */
91: PetscFE fe[LANDAU_MAX_SPECIES];
92: /* geometry */
93: PetscReal radius[LANDAU_MAX_GRIDS];
94: PetscReal re_radius; /* RE: radius of refinement along v_perp=0, z>0 */
95: PetscReal vperp0_radius1; /* RE: radius of refinement along v_perp=0 */
96: PetscReal vperp0_radius2; /* RE: radius of refinement along v_perp=0 after origin AMR refinement */
97: PetscBool sphere; // not used
98: PetscBool inflate; // not used
99: PetscReal i_radius[LANDAU_MAX_GRIDS]; // not used
100: PetscReal e_radius; // not used
101: PetscInt num_sections; // not used
102: /* AMR */
103: PetscBool use_p4est;
104: PetscInt numRERefine; /* RE: refinement along v_perp=0, z > 0 */
105: PetscInt nZRefine1; /* RE: origin refinement after v_perp=0 refinement */
106: PetscInt nZRefine2; /* RE: origin refinement after origin AMR refinement */
107: PetscInt numAMRRefine[LANDAU_MAX_GRIDS]; /* normal AMR - refine from origin */
108: PetscInt postAMRRefine[LANDAU_MAX_GRIDS]; /* uniform refinement of AMR */
109: /* relativistic */
110: PetscBool use_energy_tensor_trick;
111: PetscBool use_relativistic_corrections;
112: /* physics */
113: PetscReal thermal_temps[LANDAU_MAX_SPECIES];
114: PetscReal masses[LANDAU_MAX_SPECIES]; /* mass of each species */
115: PetscReal charges[LANDAU_MAX_SPECIES]; /* charge of each species */
116: PetscReal n[LANDAU_MAX_SPECIES]; /* number density of each species */
117: PetscReal m_0; /* reference mass */
118: PetscReal v_0; /* reference velocity */
119: PetscReal n_0; /* reference number density */
120: PetscReal t_0; /* reference time */
121: PetscReal Ez;
122: PetscReal epsilon0;
123: PetscReal k;
124: PetscReal lnLam;
125: PetscReal electronShift;
126: PetscInt num_species;
127: PetscInt num_grids;
128: PetscInt species_offset[LANDAU_MAX_GRIDS+1];
129: PetscInt mat_offset[LANDAU_MAX_GRIDS+1];
130: /* cache */
131: Mat J;
132: Mat M;
133: Vec X;
134: /* derived type */
135: void *data;
136: PetscBool aux_bool; /* helper */
137: /* computing */
138: LandauDeviceType deviceType;
139: PetscInt subThreadBlockSize;
140: PetscInt numConcurrency; /* number of SMs in Cuda to use */
141: DM pack;
142: DM plex[LANDAU_MAX_GRIDS];
143: LandauStaticData SData_d; /* static geometric data on device, this holds host pointers */
144: /* diagnostics */
145: PetscInt verbose;
146: PetscLogEvent events[20];
147: } LandauCtx;
149: typedef int LandauIdx;
150: typedef struct {
151: PetscReal scale;
152: LandauIdx gid; // Lanadu matrix index (<10,000)
153: } pointInterpolationP4est;
154: typedef struct _lP4estVertexMaps {
155: LandauIdx (*gIdx)[LANDAU_MAX_SPECIES][LANDAU_MAX_NQ]; // #elems * LANDAU_MAX_NQ (spoof for max , Nb) on device,
156: LandauIdx num_elements;
157: LandauIdx num_reduced;
158: LandauIdx num_face; // (Q or Q^2 for 3D)
159: LandauDeviceType deviceType;
160: PetscInt Nf;
161: PetscInt Nq;
162: pointInterpolationP4est (*c_maps)[LANDAU_MAX_Q_FACE];
163: struct _lP4estVertexMaps*d_self;
164: void *vp1,*vp2,*vp3;
165: PetscInt numgrids;
166: } P4estVertexMaps;
168: PETSC_EXTERN PetscErrorCode LandauCreateColoring(Mat, DM, PetscContainer *);
169: #if defined(PETSC_HAVE_CUDA)
170: PETSC_EXTERN PetscErrorCode LandauCUDAJacobian(DM[], const PetscInt, const PetscInt, const PetscInt[], PetscReal[], PetscScalar[], const PetscInt, const PetscScalar[],
171: const LandauStaticData *, const PetscInt, const PetscReal, const PetscLogEvent[],const PetscInt[], const PetscInt[], Mat[], Mat);
172: PETSC_EXTERN PetscErrorCode LandauCUDACreateMatMaps(P4estVertexMaps *, pointInterpolationP4est (*)[LANDAU_MAX_Q_FACE], PetscInt[], PetscInt, PetscInt);
173: PETSC_EXTERN PetscErrorCode LandauCUDADestroyMatMaps(P4estVertexMaps *, PetscInt);
174: PETSC_EXTERN PetscErrorCode LandauCUDAStaticDataSet(DM, const PetscInt, const PetscInt, PetscInt[], PetscInt[], PetscInt[], PetscReal [], PetscReal [], PetscReal[], PetscReal[], PetscReal[], PetscReal[], PetscReal[],
175: PetscReal[], LandauStaticData *);
176: PETSC_EXTERN PetscErrorCode LandauCUDAStaticDataClear(LandauStaticData *);
177: #endif
178: #if defined(PETSC_HAVE_KOKKOS)
179: PETSC_EXTERN PetscErrorCode LandauKokkosJacobian(DM[], const PetscInt, const PetscInt, const PetscInt[], PetscReal[], PetscScalar[], const PetscInt, const PetscScalar[],
180: const LandauStaticData *, const PetscInt, const PetscReal, const PetscLogEvent[], const PetscInt[], const PetscInt[], Mat[], Mat);
181: PETSC_EXTERN PetscErrorCode LandauKokkosCreateMatMaps(P4estVertexMaps *, pointInterpolationP4est (*)[LANDAU_MAX_Q_FACE], PetscInt[], PetscInt, PetscInt);
182: PETSC_EXTERN PetscErrorCode LandauKokkosDestroyMatMaps(P4estVertexMaps *, PetscInt);
183: PETSC_EXTERN PetscErrorCode LandauKokkosStaticDataSet(DM, const PetscInt, const PetscInt, PetscInt[], PetscInt[], PetscInt[], PetscReal [], PetscReal [], PetscReal[], PetscReal[], PetscReal[], PetscReal[], PetscReal[],
184: PetscReal[], LandauStaticData *);
185: PETSC_EXTERN PetscErrorCode LandauKokkosStaticDataClear(LandauStaticData*);
186: #endif
188: #endif /* PETSCLANDAU_H */