XGCa
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
col_grid.hpp
Go to the documentation of this file.
1 #ifndef COL_GRID_HPP
2 #define COL_GRID_HPP
3 
4 #include <vector>
5 #include <stdio.h>
6 #include <Kokkos_Core.hpp>
7 #include <Kokkos_DualView.hpp>
8 
10 #include "velocity_grid.hpp"
11 
12 extern "C" int get_col_range_setup(FortranPtr grid_fptr, int inode);
13 
14 //> XGC's multi-species Fokker-Planck-Landau collision operator
15 // see R. Hager et al., Journal of Computational Physics 315 (2016), pp. 644-660
16 // E. S. Yoon et al., Physics of Plasmas 21 (2014), 032503
17 //
18 // The headers of the subroutines in this module refer to the corresponding parts of Hager et al.
19 // where applicable in order to make the source code easier to understand.
20 
21 
22 // <SOLVER-related>
23 const int vpic_inner_iter_max=20;
24 
25 // set use_superlu = .true. to use sparse direct SuperLU solver
26 // set use_superlu = .false. to use lapack band solver
27 const bool use_superlu=false;
28 
30  int gid; // Corresponding XGC species
31  bool is_electron; // Whether species is electron
32  double T_avg;
33  double mass;
34  double mass_au;
35  double charge_eu;
36  double vpar_beg;
38  double conv_factor;
39  double numeric_vth2; //[col_fm_core_init<col_fm_core]
40  double numeric_t; //[col_fm_core_init<col_fm_core]
41  double mesh_dz;
42  double mesh_dr;
43  double dens;
44  double ens;
45  double mom; // [col_fm_core_init<col_fm_core]
46 };
47 
48 template<class Device>
50 
51  Kokkos::DualView<CollisionSpeciesScalars**,Device> s;
52 
53  Kokkos::View<double****,HostType> pdf_n; //[get_local_f<setup_col_sp] Given initial PDF
54  Kokkos::DualView<double****,Kokkos::LayoutRight,Device> pdf_np1; //[get_local_f<setup_col_sp] To be used as PDF at next Picard step. BUT has a role of "df" at end
55 
56  Kokkos::View<double***,Kokkos::LayoutRight,HostType> mesh_r;
57  Kokkos::View<double***,Kokkos::LayoutRight,HostType> mesh_z;
58  Kokkos::View<double***,Kokkos::LayoutRight,Device> vol;
59  Kokkos::View<double***,Kokkos::LayoutRight,HostType> vol_h;
60  Kokkos::View<double****,Device> delta_r;
61  Kokkos::View<double****,Device> delta_z;
62 
63  // These three are sent to GPU, so use a dualview
64  Kokkos::DualView<double***,Device> mesh_r_half;
65  Kokkos::DualView<double***,Device> mesh_z_half;
66  Kokkos::DualView<double***,Device> local_center_volume;
67 
68  public:
69 
71 
72  CollisionSpecies<Device>(int col_f_nvr, int col_f_nvz, int col_f_sp_num, std::vector<Species<DeviceType>> &all_species, int n_species, int mb_n_nodes)
73  : s("s", n_species, mb_n_nodes),
74  pdf_n("pdf_n", mb_n_nodes, n_species, col_f_nvz, col_f_nvr),
75  pdf_np1("pdf_np1", mb_n_nodes, n_species, col_f_nvz, col_f_nvr),
76  mesh_r("mesh_r", mb_n_nodes, n_species, col_f_nvr),
77  mesh_z("mesh_z", mb_n_nodes, n_species, col_f_nvz),
78  mesh_r_half("mesh_r_half",mb_n_nodes,col_f_nvr-1,n_species),
79  mesh_z_half("mesh_z_half",mb_n_nodes,col_f_nvz-1,n_species),
80  local_center_volume("local_center_volume",mb_n_nodes,col_f_nvr-1,n_species),
81  vol("vol",col_f_nvr,n_species, mb_n_nodes),
82  vol_h("vol_h",mb_n_nodes,n_species, col_f_nvr), // Must have col_f_nvr fastest index for passing into fortran negative_f_correction routine
83  delta_r("delta_r",col_f_sp_num,col_f_nvr-1,n_species, mb_n_nodes),
84  delta_z("delta_z",col_f_sp_num,col_f_nvz-1,n_species, mb_n_nodes)
85  {
86  int col_species_cnt=0;
87  for (int isp = 0; isp<all_species.size(); isp++){
88  if (!all_species[isp].is_adiabatic){
89  for (int mesh_ind = 0; mesh_ind<mb_n_nodes; mesh_ind++){
90  s.h_view(col_species_cnt, mesh_ind).gid = isp;
91  s.h_view(col_species_cnt, mesh_ind).is_electron = all_species[isp].is_electron;
92  }
93  col_species_cnt++;
94  }
95  }
96  }
97 
98  void setup_one(int isp, int mesh_ind, Species<DeviceType> &all_species, int inode, int first_node, int nvr, int nvz, double dsmu, double dvp, double vp_max);
99  void setup_all(std::vector<Species<DeviceType>>& all_species, Kokkos::View<int**,HostType>& mesh_nodes, int mesh_batch_ind, int first_node, int nvr, int nvz, double dsmu, double dvp, double vp_max);
100 };
101 
102 template<class Device>
103 struct TmpColData{
104  Kokkos::View<double***,HostType> dist_col; // Raw pointer used so careful about padding - ALS
105  Kokkos::DualView<double***,Kokkos::LayoutRight,Device> LU_values;// Super LU Variables // Raw pointer used so careful about padding - ALS
106 
107  // These views are both on CPU and GPU - dual views
108  Kokkos::DualView<double***,Device> gammac_spall;
109  Kokkos::DualView<int[4][4],Device> mat_pos_rel_indx; // Shouldn't need to be a dual view
110 
111  Kokkos::View<double***,Device> fhalf;
112  Kokkos::View<double***,Device> dfdr;
113  Kokkos::View<double***,Device> dfdz;
114  Kokkos::View<double****,Device> EDs;
115  Kokkos::View<double******,Device> M_ab;
116  Kokkos::View<double*****,Device> M_s;
117 
119 
120  TmpColData<Device>(int col_f_nvr, int col_f_nvz, int col_f_sp_num, int col_f_ntotal_v, int LU_nnz, int mb_n_nodes)
121  // Note: layout right
122  : dist_col("dist_col",mb_n_nodes,col_f_sp_num,col_f_ntotal_v),
123  LU_values("LU_values",mb_n_nodes,col_f_sp_num,LU_nnz),
124 
125  // Note: layout left
126  gammac_spall("gammac_spall",col_f_sp_num,col_f_sp_num,mb_n_nodes),
127  fhalf("fhalf",mb_n_nodes,col_f_nvz-1,col_f_nvr-1),
128  dfdr("dfdr",mb_n_nodes,col_f_nvz-1,col_f_nvr-1),
129  dfdz("dfdz",mb_n_nodes,col_f_nvz-1,col_f_nvr-1),
130  EDs("EDs",mb_n_nodes,col_f_nvz-1,col_f_nvr-1,5),
131  M_s("M_s",mb_n_nodes,col_f_sp_num, (col_f_nvr-1)*(col_f_nvz-1), 5, col_f_nvr-1),
132  M_ab("M_ab",mb_n_nodes,col_f_sp_num, col_f_sp_num-1, (col_f_nvr-1)*(col_f_nvz-1), 3, (col_f_nvr-1)*(col_f_nvz-1)),
133  mat_pos_rel_indx("mat_pos_rel_indx")
134  {
135  // Set up mat_pos_rel_indx and send to GPU
136  int mat_loc[16] = {4,5,7,8,
137  1,2,4,5,
138  3,4,6,7,
139  0,1,3,4};
140 
141  for (int i=0;i<16;i++)
142  mat_pos_rel_indx.h_view(i/4,i%4)=mat_loc[i];
143 
144 #ifdef USE_GPU
145  Kokkos::deep_copy(mat_pos_rel_indx.d_view, mat_pos_rel_indx.h_view);
146 #endif
147  }
148 };
149 
150 
151 template<class Device>
153 
154  public:
155 
156  int col_f_nvr, col_f_nvz; //< The number of "GRIDS" in vr and vz, respectively
157  int col_f_nvrm1,col_f_nvzm1; //< To enhance understandings of compiler and conseq. to avoid temp. copy
158  int col_f_nvrm1_nvzm1; //< To enhance understandings of compiler and conseq. to avoid temp. copy
159  double col_f_vp_max; //< Max v_para of velocity-grid
160  double col_f_smu_max; //< Max v_perp of velocity-grid
161  double col_f_dvp; //< v_para grid spacing
162  double col_f_dsmu; //< v_perp grid spacing
163  int col_f_ntotal_v; //< set in col_f_setup
164  double col_f_dt; //< Time step for collision operation
165  Kokkos::DualView<int**,Device> index_map_LU;
166  Kokkos::View<int*,HostType> LU_cvalues;
167  Kokkos::View<int*,HostType> LU_rowindx;
168  Kokkos::View<int*,HostType> LU_colptr;
169  int LU_n;
170  int LU_nnz;
171  int LU_nrhs;
172  int LU_ldb;
173 
174  int col_f_sp_num; // [col_f_setup] the number of species to be collided ( = col_f_sp_index_e-col_f_sp_index_s+1 )
175 
176  double* node_cost; // pointer to array of timings for collisions
177 
178  Kokkos::View<bool*, Kokkos::HostSpace> in_range;
179 
181 
182  //> Set up global variables of the collision operator
183  // @param[in] nspecies Number of species
184  // @param[in] dt Time step, real(8)
185  CollisionGrid(const VelocityGrid& vgrid, int nspecies,double dt,int f_source_period,bool sml_symmetric_f, int nnode, int first_node, double* node_cost_ptr, FortranPtr grid_fptr)
186  : node_cost(node_cost_ptr), // Point timer array to fortran array
187  symmetric_f(sml_symmetric_f),
188  in_range("in_range",nnode)
189  {
190  // nodes outside [sml_inpsi,sml_outpsi] and in private flux region are skipped
191  for (int i=0;i<nnode;i++){
192  int ftr_node = i + first_node;
193  in_range(i) = (get_col_range_setup(grid_fptr, ftr_node)==1);
194  }
195 
196  //##### Multi-species >>> !! DEV_IMP
197  col_f_sp_num=nspecies;
198 
199  col_f_nvr=vgrid.nvr;
200  col_f_nvz=vgrid.nvz;
205  col_f_dt=dt*f_source_period; // delta time for collision operation
206  col_f_vp_max=vgrid.vp_max;
207  col_f_smu_max=vgrid.smu_max;
208  col_f_dvp=vgrid.dvp;
209  col_f_dsmu=vgrid.dsmu;
210 
211  // 2013-02-23 =============================================== SUPER LU
212  // index_map_LU, LU_rowindx, LU_cvalues
214  LU_nnz=4*4+6*2*(col_f_nvr-2)+6*2*(col_f_nvz-2)+9*(col_f_nvr-2)*(col_f_nvz-2);
215  LU_nrhs=1;
216  LU_ldb=LU_n;
217 
218  // Construct CollisionGrid views
219  LU_rowindx = Kokkos::View<int*,HostType>("LU_rowindx",LU_nnz);
220  LU_cvalues = Kokkos::View<int*,HostType>("LU_cvalues",LU_n);
221  LU_colptr = Kokkos::View<int*,HostType>("LU_colptr",LU_n+1);
222  index_map_LU = Kokkos::DualView<int**,Device>("index_map_LU",LU_n,9);
223 
224  // Local
225  Kokkos::View<int*,HostType> LU_colptr_num("LU_colptr_num",LU_n);
226 
227  //below is time independent (what does this comment me
228 
229  LU_colptr(0)=0;
230  for (int i=0; i<LU_n; i++){
231  int incr_LU;
232  //for colptr
233  int LU_i = i/col_f_nvr;
234  int LU_j = i%col_f_nvr;
235  if((LU_i==0)||(LU_i==col_f_nvz-1)){
236  if((LU_j==0)||(LU_j==col_f_nvr-1)){
237  incr_LU=4;
238  } else {
239  incr_LU=6;
240  }
241  } else {
242  if((LU_j==0)||(LU_j==col_f_nvr-1)){
243  incr_LU=6;
244  } else {
245  incr_LU=9;
246  }
247  }
248  LU_colptr(i+1)=LU_colptr(i)+incr_LU;
249  }
250 
251  //===============
252  // 2--5--8
253  // | | |
254  // 1--4--7
255  // | | |
256  // 0--3--6
257  //==============
258  int mat_pos_rel[9];
259  mat_pos_rel[0]=-col_f_nvr-1;
260  mat_pos_rel[1]=-col_f_nvr;
261  mat_pos_rel[2]=-col_f_nvr+1;
262  mat_pos_rel[3]=-1;
263  mat_pos_rel[4]=0;
264  mat_pos_rel[5]=1;
265  mat_pos_rel[6]=col_f_nvr-1;
266  mat_pos_rel[7]=col_f_nvr;
267  mat_pos_rel[8]=col_f_nvr+1;
268  int LU_i_arr[9];
269  std::vector<std::vector<int>> mat_pos_rel_indx;
270  mat_pos_rel_indx.push_back(std::vector<int>{4,5,7,8}); // Upper right
271  mat_pos_rel_indx.push_back(std::vector<int>{3,4,6,7}); // Lower right
272  mat_pos_rel_indx.push_back(std::vector<int>{3,4,5,6,7,8}); // Right
273  mat_pos_rel_indx.push_back(std::vector<int>{1,2,4,5}); // Upper left
274  mat_pos_rel_indx.push_back(std::vector<int>{0,1,3,4}); // Lower left
275  mat_pos_rel_indx.push_back(std::vector<int>{0,1,2,3,4,5}); // Left
276  mat_pos_rel_indx.push_back(std::vector<int>{1,2,4,5,7,8}); // Top
277  mat_pos_rel_indx.push_back(std::vector<int>{0,1,3,4,6,7}); // Bottom
278  mat_pos_rel_indx.push_back(std::vector<int>{0,1,2,3,4,5,6,7,8}); // Center
279 
280  for (int i=0; i<col_f_nvz; i++){
281  for (int j=0; j<col_f_nvr; j++){
282  int mat_pos=j+i*col_f_nvr;
283  int loc;
284  // INDEXING
285  if(i==0){
286  if(j==0){
287  //(J,I)=(1,1)
288  loc = 0;
289  } else if(j==col_f_nvr-1){
290  //(J,I)=(Nr,1)
291  loc = 1;
292  } else {
293  //(J,I)=(:,1)
294  loc = 2;
295  }
296  } else if(i==col_f_nvz-1){
297  if(j==0){
298  //(J,I)=(1,mesh_Nz)
299  loc = 3;
300  } else if(j==col_f_nvr-1){
301  //(J,I)=(Nr,mesh_Nz)
302  loc = 4;
303  } else {
304  //(J,I)=(:,mesh_Nz)
305  loc = 5;
306  }
307  } else {
308  if(j==0){
309  //(J,I) = (1,:)
310  loc = 6;
311  } else if(j==col_f_nvr-1){
312  //(J,I) = (mesh_Nr,:)
313  loc = 7;
314  } else {
315  //(J,I) = (:,:)
316  loc = 8;
317  }
318  }
319 
320  for (int i_elem=0; i_elem<mat_pos_rel_indx[loc].size(); i_elem++)
321  LU_i_arr[i_elem]=mat_pos+mat_pos_rel[mat_pos_rel_indx[loc][i_elem]];
322 
323  for (int i_elem=0; i_elem<mat_pos_rel_indx[loc].size(); i_elem++)
324  index_map_LU.h_view(mat_pos,mat_pos_rel_indx[loc][i_elem])=LU_colptr(LU_i_arr[i_elem])+LU_colptr_num(LU_i_arr[i_elem]);
325 
326  for (int i_elem=0; i_elem<mat_pos_rel_indx[loc].size(); i_elem++)
327  LU_colptr_num(LU_i_arr[i_elem])+=1;
328 
329  for (int i_elem=0; i_elem<mat_pos_rel_indx[loc].size(); i_elem++)
330  LU_rowindx(index_map_LU.h_view(mat_pos,mat_pos_rel_indx[loc][i_elem]))=mat_pos;
331 
332  LU_cvalues(mat_pos)=index_map_LU.h_view(mat_pos,4); //For implicit time marching
333  }
334  }
335 
336  // Send to GPU
337  index_map_LU.template modify<typename Kokkos::DualView<int**,Device>::host_mirror_space>();
338  index_map_LU.template sync<typename Kokkos::DualView<int**,Device>::execution_space>();
339  }
340 
341  // Computes the base collision frequency -- Eq. (6) in Hager et al.
342  double lambda_gamma_pair(const CollisionSpeciesScalars& sp_a, const CollisionSpeciesScalars& sp_b) const;// C(a->b)
343 
344  // Core init
345  void core_init(int spi, int mesh_ind, const CollisionSpecies<Device>& col_spall) const;
346 
347  // Get maxw_fac
348  static KOKKOS_INLINE_FUNCTION double get_maxw_fac(double mesh_dr, double mesh_r, double numeric_vth2) ;
349 
350  // Delta init
351  void core_delta_init(int mb_n_nodes, int spi, int spj, CollisionSpecies<Device>& col_spall) const;
352 
353  // LU_matrix_ftn
354  static KOKKOS_INLINE_FUNCTION void LU_matrix_ftn(int mesh_ind, int spi, int spj,int cell_i,int cell_j,const CollisionSpecies<Device>& col_spall, int mprl_col,int mat_pos,double coeff1,double coeff2,const TmpColData<Device>& tcd, const Kokkos::View<int**,Device>& index_map_LU_d) ;
355 
356  // LU_matrix
357  void LU_matrix(int mb_n_nodes, int spi, int spj, CollisionSpecies<Device>& col_spall, TmpColData<Device>& tcd, const Kokkos::View<int**,Device>& index_map_LU_d) const;
358 
359  // Picard step
360  void picard_step(int iter_inter,const Kokkos::DualView<double***,Kokkos::LayoutRight,Device>& LU_values,const Kokkos::View<double***,HostType>& dist_col,int spi, int mesh_ind, const Kokkos::View<double****,Kokkos::LayoutRight,HostType>& dist_iter) const;
361 
362  // dispatch for E_and_D calc
363  void E_and_D(int mb_n_nodes, int spi, int spj, CollisionSpecies<Device>& col_spall, TmpColData<Device>& tcd) const;
364 
365  // E_and_D_s
366  static KOKKOS_INLINE_FUNCTION void E_and_D_s(int idx, int mb_n_nodes, int nvrm1, int nvzm1, double inv_mass, const TmpColData<Device>& tcd, int spi) ;
367 
368  // dispatch for angle_avg calc
369  void angle_avg(int mb_n_nodes, int spi, int spj, CollisionSpecies<Device>& col_spall, TmpColData<Device>& tcd) const;
370 
371  // angle_avg_s
372  static KOKKOS_INLINE_FUNCTION void angle_avg_s(int idx, int mb_n_nodes, int nvrm1, int nvzm1, const CollisionSpecies<Device>& col_spall, const TmpColData<Device>& tcd, int spi) ;
373 
374  // E_and_D_ab
375  static KOKKOS_INLINE_FUNCTION void E_and_D_ab(int idx, int mb_n_nodes,int nvrm1, int nvzm1, const CollisionSpecies<Device>& col_spall, const TmpColData<Device>& tcd, int spi, int spj_mod) ;
376 
377  // angle_avg_ab
378  static KOKKOS_INLINE_FUNCTION void angle_avg_ab(int idx, int mb_n_nodes, int nvrm1, int nvzm1, const CollisionSpecies<Device>& col_spall, const TmpColData<Device>& tcd, int spi, int spj_mod) ;
379 
380  void f_df(int mb_n_nodes, CollisionSpecies<Device>& col_spall,int spi,int spj,TmpColData<Device>& tcd) const;
381 
382  // Main collision algorithm
383  Kokkos::View<int*,HostType> core(CollisionSpecies<Device>& col_spall, TmpColData<Device>& tcd, int mb_n_nodes) const;
384 
385 
386  // Carries out collisions
387  void collision(std::vector<Species<DeviceType>> &all_species, const DomainDecomposition<DeviceType>& pol_decomp) const;
388 };
389 
390 #include "col_grid.tpp"
391 
392 #endif
Kokkos::DualView< double ****, Kokkos::LayoutRight, Device > pdf_np1
Definition: col_grid.hpp:54
double lambda_gamma_pair(const CollisionSpeciesScalars &sp_a, const CollisionSpeciesScalars &sp_b) const
Definition: col_grid.tpp:237
double smu_max
Definition: velocity_grid.hpp:11
void core_init(int spi, int mesh_ind, const CollisionSpecies< Device > &col_spall) const
Definition: col_grid.tpp:658
const bool use_superlu
Definition: col_grid.hpp:27
Kokkos::View< double ******, Device > M_ab
Definition: col_grid.hpp:115
int col_f_nvrm1_nvzm1
Definition: col_grid.hpp:158
int * FortranPtr
Definition: globals.hpp:80
Kokkos::View< double *****, Device > M_s
Definition: col_grid.hpp:116
Definition: col_grid.hpp:152
double mass_au
Definition: col_grid.hpp:34
Definition: velocity_grid.hpp:5
bool symmetric_f
Definition: col_grid.hpp:180
double den_lambda_gamma
Definition: col_grid.hpp:37
Kokkos::View< double ***, Device > fhalf
Definition: col_grid.hpp:111
bool is_electron
Definition: col_grid.hpp:31
double vpar_beg
Definition: col_grid.hpp:36
void LU_matrix(int mb_n_nodes, int spi, int spj, CollisionSpecies< Device > &col_spall, TmpColData< Device > &tcd, const Kokkos::View< int **, Device > &index_map_LU_d) const
Definition: col_grid.tpp:1229
int LU_nrhs
Definition: col_grid.hpp:171
static KOKKOS_INLINE_FUNCTION void E_and_D_s(int idx, int mb_n_nodes, int nvrm1, int nvzm1, double inv_mass, const TmpColData< Device > &tcd, int spi)
Definition: col_grid.tpp:977
Kokkos::DualView< double ***, Device > mesh_r_half
Definition: col_grid.hpp:64
double mesh_dr
Definition: col_grid.hpp:42
double col_f_dsmu
Definition: col_grid.hpp:162
int nvr
Definition: velocity_grid.hpp:16
int gid
Definition: col_grid.hpp:30
Kokkos::View< int *, HostType > core(CollisionSpecies< Device > &col_spall, TmpColData< Device > &tcd, int mb_n_nodes) const
Definition: col_grid.tpp:416
void E_and_D(int mb_n_nodes, int spi, int spj, CollisionSpecies< Device > &col_spall, TmpColData< Device > &tcd) const
Definition: col_grid.tpp:789
Kokkos::View< double ***, Kokkos::LayoutRight, HostType > mesh_z
Definition: col_grid.hpp:57
double conv_factor
Definition: col_grid.hpp:38
static KOKKOS_INLINE_FUNCTION void angle_avg_s(int idx, int mb_n_nodes, int nvrm1, int nvzm1, const CollisionSpecies< Device > &col_spall, const TmpColData< Device > &tcd, int spi)
Definition: col_grid.tpp:870
Kokkos::View< double ***, Kokkos::LayoutRight, Device > vol
Definition: col_grid.hpp:58
Kokkos::View< bool *, Kokkos::HostSpace > in_range
Definition: col_grid.hpp:178
Definition: col_grid.hpp:29
Kokkos::View< double ****, Device > delta_r
Definition: col_grid.hpp:60
Kokkos::DualView< int **, Device > index_map_LU
Definition: col_grid.hpp:165
void setup_one(int isp, int mesh_ind, Species< DeviceType > &all_species, int inode, int first_node, int nvr, int nvz, double dsmu, double dvp, double vp_max)
Definition: col_grid.tpp:129
double col_f_dt
Definition: col_grid.hpp:164
Definition: col_grid.hpp:103
const int vpic_inner_iter_max
Definition: col_grid.hpp:23
int col_f_nvz
Definition: col_grid.hpp:156
Kokkos::View< double ***, Kokkos::LayoutRight, HostType > mesh_r
Definition: col_grid.hpp:56
double numeric_t
Definition: col_grid.hpp:40
Kokkos::DualView< CollisionSpeciesScalars **, Device > s
Definition: col_grid.hpp:51
void picard_step(int iter_inter, const Kokkos::DualView< double ***, Kokkos::LayoutRight, Device > &LU_values, const Kokkos::View< double ***, HostType > &dist_col, int spi, int mesh_ind, const Kokkos::View< double ****, Kokkos::LayoutRight, HostType > &dist_iter) const
Definition: col_grid.tpp:1023
Kokkos::DualView< double ***, Device > local_center_volume
Definition: col_grid.hpp:66
Kokkos::View< double ***, HostType > dist_col
Definition: col_grid.hpp:104
Kokkos::View< double ***, Device > dfdz
Definition: col_grid.hpp:113
void angle_avg(int mb_n_nodes, int spi, int spj, CollisionSpecies< Device > &col_spall, TmpColData< Device > &tcd) const
Definition: col_grid.tpp:703
double mesh_dz
Definition: col_grid.hpp:41
CollisionGrid(const VelocityGrid &vgrid, int nspecies, double dt, int f_source_period, bool sml_symmetric_f, int nnode, int first_node, double *node_cost_ptr, FortranPtr grid_fptr)
Definition: col_grid.hpp:185
double mass
Definition: col_grid.hpp:33
static KOKKOS_INLINE_FUNCTION void E_and_D_ab(int idx, int mb_n_nodes, int nvrm1, int nvzm1, const CollisionSpecies< Device > &col_spall, const TmpColData< Device > &tcd, int spi, int spj_mod)
Definition: col_grid.tpp:813
int get_col_range_setup(FortranPtr grid_fptr, int inode)
Kokkos::DualView< double ***, Kokkos::LayoutRight, Device > LU_values
Definition: col_grid.hpp:105
static KOKKOS_INLINE_FUNCTION void LU_matrix_ftn(int mesh_ind, int spi, int spj, int cell_i, int cell_j, const CollisionSpecies< Device > &col_spall, int mprl_col, int mat_pos, double coeff1, double coeff2, const TmpColData< Device > &tcd, const Kokkos::View< int **, Device > &index_map_LU_d)
Definition: col_grid.tpp:1184
Kokkos::View< double ****, Device > delta_z
Definition: col_grid.hpp:61
double vp_max
Definition: velocity_grid.hpp:7
void setup_all(std::vector< Species< DeviceType >> &all_species, Kokkos::View< int **, HostType > &mesh_nodes, int mesh_batch_ind, int first_node, int nvr, int nvz, double dsmu, double dvp, double vp_max)
Definition: col_grid.tpp:200
double charge_eu
Definition: col_grid.hpp:35
int LU_nnz
Definition: col_grid.hpp:170
double col_f_smu_max
Definition: col_grid.hpp:160
int col_f_ntotal_v
Definition: col_grid.hpp:163
void collision(std::vector< Species< DeviceType >> &all_species, const DomainDecomposition< DeviceType > &pol_decomp) const
Definition: col_grid.tpp:299
Kokkos::DualView< double ***, Device > mesh_z_half
Definition: col_grid.hpp:65
Kokkos::View< double ****, Device > EDs
Definition: col_grid.hpp:114
double numeric_vth2
Definition: col_grid.hpp:39
double * node_cost
Definition: col_grid.hpp:176
void core_delta_init(int mb_n_nodes, int spi, int spj, CollisionSpecies< Device > &col_spall) const
Definition: col_grid.tpp:1107
Kokkos::DualView< double ***, Device > gammac_spall
Definition: col_grid.hpp:108
Definition: domain_decomposition.hpp:7
double T_avg
Definition: col_grid.hpp:32
double col_f_dvp
Definition: col_grid.hpp:161
Kokkos::View< double ***, Device > dfdr
Definition: col_grid.hpp:112
double dsmu
Definition: velocity_grid.hpp:12
double dens
Definition: col_grid.hpp:43
int col_f_nvzm1
Definition: col_grid.hpp:157
double mom
Definition: col_grid.hpp:45
int col_f_nvr
Definition: col_grid.hpp:156
Definition: species.hpp:13
int nvz
Definition: velocity_grid.hpp:17
int LU_n
Definition: col_grid.hpp:169
double col_f_vp_max
Definition: col_grid.hpp:159
Kokkos::DualView< int[4][4], Device > mat_pos_rel_indx
Definition: col_grid.hpp:109
Definition: col_grid.hpp:49
Kokkos::View< int *, HostType > LU_colptr
Definition: col_grid.hpp:168
double dvp
Definition: velocity_grid.hpp:8
static KOKKOS_INLINE_FUNCTION double get_maxw_fac(double mesh_dr, double mesh_r, double numeric_vth2)
Definition: col_grid.tpp:1093
int LU_ldb
Definition: col_grid.hpp:172
void f_df(int mb_n_nodes, CollisionSpecies< Device > &col_spall, int spi, int spj, TmpColData< Device > &tcd) const
Definition: col_grid.tpp:934
Kokkos::View< int *, HostType > LU_cvalues
Definition: col_grid.hpp:166
int col_f_nvrm1
Definition: col_grid.hpp:157
Kokkos::View< double ****, HostType > pdf_n
Definition: col_grid.hpp:53
double ens
Definition: col_grid.hpp:44
static KOKKOS_INLINE_FUNCTION void angle_avg_ab(int idx, int mb_n_nodes, int nvrm1, int nvzm1, const CollisionSpecies< Device > &col_spall, const TmpColData< Device > &tcd, int spi, int spj_mod)
Definition: col_grid.tpp:721
Kokkos::View< int *, HostType > LU_rowindx
Definition: col_grid.hpp:167
Kokkos::View< double ***, Kokkos::LayoutRight, HostType > vol_h
Definition: col_grid.hpp:59
int col_f_sp_num
Definition: col_grid.hpp:174