xref: /petsc/include/petsclandau.h (revision 608e5a7a8a35b40c4e75a0aae313cbdeb3df9fee)
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 
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
27 
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
40 
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
48 
49 typedef enum {LANDAU_CUDA, LANDAU_KOKKOS, LANDAU_CPU} LandauDeviceType;
50 
51 /* typedef PetscReal LandauIPReal; */
52 /* typedef struct { */
53 /*   LandauIPReal  *coefs; */
54 /*   int           dim_,ns_,nip_; */
55 /* } LandauIPFdF; */
56 
57 typedef struct {
58   void  *invJ;  // nip*dim*dim
59   void  *D;     // nq*nb*dim
60   void  *B;     // nq*nb
61   void  *alpha; // ns
62   void  *beta;  // ns
63   void  *invMass; // ns
64   void  *mass_w;  // nip
65   void  *w; // nip
66   void  *x; // nip
67   void  *y; // nip
68   void  *z; // nip
69   void  *Eq_m; // ns - dynamic
70   void  *f; //  nip*Nf - dynamic (IP)
71   void  *dfdx; // nip*Nf - dynamic (IP)
72   void  *dfdy; // nip*Nf - dynamic (IP)
73   void  *dfdz; // nip*Nf - dynamic (IP)
74   void  *IPf;  // Ncells*Nb*Nf - dynamic (vertex in cells)
75   int   dim_,ns_,nip_,nq_,nb_;
76 } LandauGeomData;
77 
78 typedef struct {
79   PetscBool      interpolate;                  /* Generate intermediate mesh elements */
80   PetscBool      gpu_assembly;
81   PetscFE        fe[LANDAU_MAX_SPECIES];
82   /* geometry  */
83   PetscReal      i_radius;
84   PetscReal      e_radius;
85   PetscInt       num_sections;
86   PetscReal      radius;
87   PetscReal      re_radius;           /* radius of refinement along v_perp=0, z>0 */
88   PetscReal      vperp0_radius1;      /* radius of refinement along v_perp=0 */
89   PetscReal      vperp0_radius2;      /* radius of refinement along v_perp=0 after origin AMR refinement */
90   PetscBool      sphere;
91   PetscBool      inflate;
92   PetscInt       numRERefine;       /* refinement along v_perp=0, z > 0 */
93   PetscInt       nZRefine1;          /* origin refinement after v_perp=0 refinement */
94   PetscInt       nZRefine2;          /* origin refinement after origin AMR refinement */
95   PetscInt       maxRefIts;         /* normal AMR - refine from origin */
96   PetscInt       postAMRRefine;     /* uniform refinement of AMR */
97   /* discretization - AMR */
98   PetscErrorCode (*errorIndicator)(PetscInt, PetscReal, PetscReal [], PetscInt, const PetscInt[], const PetscScalar[], const PetscScalar[], PetscReal *, void *);
99   PetscReal      refineTol[LANDAU_MAX_SPECIES];
100   PetscReal      coarsenTol[LANDAU_MAX_SPECIES];
101   /* physics */
102   PetscReal      thermal_temps[LANDAU_MAX_SPECIES];
103   PetscReal      masses[LANDAU_MAX_SPECIES];  /* mass of each species  */
104   PetscReal      charges[LANDAU_MAX_SPECIES]; /* charge of each species  */
105   PetscReal      n[LANDAU_MAX_SPECIES];       /* number density of each species  */
106   PetscReal      m_0;      /* reference mass */
107   PetscReal      v_0;      /* reference velocity */
108   PetscReal      n_0;      /* reference number density */
109   PetscReal      t_0;      /* reference time */
110   PetscReal      Ez;
111   PetscReal      epsilon0;
112   PetscReal      k;
113   PetscReal      lnLam;
114   PetscReal      electronShift; /* for tests */
115   PetscInt       num_species;
116   /* cache */
117   Mat            J;
118   Mat            M;
119   Vec            X;
120   /* derived type */
121   void          *data;
122   PetscBool      aux_bool;  /* helper */
123   /* computing */
124   LandauDeviceType deviceType;
125   PetscInt       subThreadBlockSize;
126   PetscInt       numConcurrency; /* number of SMs in Cuda to use */
127   MPI_Comm       comm; /* global communicator to use for errors and diagnostics */
128   LandauGeomData *SData_d; /* static geometric data on device, but this pointer is a host pointer */
129   double         times[1];
130   PetscBool      init;
131   PetscBool      use_matrix_mass;
132   DM             dmv;
133   DM             plex;
134   /* diagnostics */
135   PetscInt       verbose;
136   PetscLogEvent  events[20];
137 } LandauCtx;
138 
139 typedef int LandauIdx;
140 typedef struct {
141   PetscReal scale;
142   LandauIdx gid;   // Lanadu matrix index (<10,000)
143 } pointInterpolationP4est;
144 typedef struct _lP4estVertexMaps {
145   LandauIdx                (*gIdx)[LANDAU_MAX_SPECIES][LANDAU_MAX_NQ]; // #elems *  LANDAU_MAX_NQ (spoof for max , Nb) on device,
146   LandauIdx                num_elements;
147   LandauIdx                num_reduced;
148   LandauIdx                num_face;  // (Q or Q^2 for 3D)
149   LandauDeviceType         deviceType;
150   PetscInt                 Nf;
151   PetscInt                 Nq;
152   pointInterpolationP4est (*c_maps)[LANDAU_MAX_Q_FACE];
153   struct _lP4estVertexMaps*data;
154   void                    *vp1,*vp2,*vp3;
155 } P4estVertexMaps;
156 
157 PETSC_EXTERN PetscErrorCode LandauCreateColoring(Mat, DM, PetscContainer *);
158 #if defined(PETSC_HAVE_CUDA)
159 PETSC_EXTERN PetscErrorCode LandauCUDAJacobian(DM, const PetscInt, PetscReal[], PetscScalar[], const PetscInt, const PetscScalar[], LandauGeomData *, const PetscInt, PetscReal, const PetscLogEvent[], Mat);
160 PETSC_EXTERN PetscErrorCode LandauCUDACreateMatMaps(P4estVertexMaps *, pointInterpolationP4est (*)[LANDAU_MAX_Q_FACE], PetscInt, PetscInt);
161 PETSC_EXTERN PetscErrorCode LandauCUDADestroyMatMaps(P4estVertexMaps *);
162 PETSC_EXTERN PetscErrorCode LandauCUDAStaticDataSet(DM, const PetscInt, PetscReal [], PetscReal [], PetscReal[], PetscReal[], PetscReal[], PetscReal[], PetscReal[], PetscReal[], PetscReal[], LandauGeomData *);
163 PETSC_EXTERN PetscErrorCode LandauCUDAStaticDataClear(LandauGeomData *);
164 #endif
165 #if defined(PETSC_HAVE_KOKKOS)
166 PETSC_EXTERN PetscErrorCode LandauKokkosJacobian(DM, const PetscInt, PetscReal[], PetscScalar[],  const PetscInt, const PetscScalar[], LandauGeomData *, const PetscInt, PetscReal, const PetscLogEvent[], Mat);
167 PETSC_EXTERN PetscErrorCode LandauKokkosCreateMatMaps(P4estVertexMaps *, pointInterpolationP4est (*)[LANDAU_MAX_Q_FACE], PetscInt, PetscInt);
168 PETSC_EXTERN PetscErrorCode LandauKokkosDestroyMatMaps(P4estVertexMaps *);
169 PETSC_EXTERN PetscErrorCode LandauKokkosStaticDataSet(DM, const PetscInt, PetscReal [], PetscReal [], PetscReal[], PetscReal[], PetscReal[], PetscReal[], PetscReal[], PetscReal[], PetscReal[], LandauGeomData *);
170 PETSC_EXTERN PetscErrorCode LandauKokkosStaticDataClear(LandauGeomData *);
171 #endif
172 
173 #endif /* PETSCLANDAU_H */
174