XGCa
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
electric_field.hpp
Go to the documentation of this file.
1 #ifndef ELECTRIC_FIELD_HPP
2 #define ELECTRIC_FIELD_HPP
3 #include "space_settings.hpp"
4 #include "NamelistReader.hpp"
5 #include "grid.hpp"
6 #include "grid_field_pack.hpp"
7 #include "timer_macro.hpp"
8 #include "get_potential_grad.hpp"
10 #ifdef SONIC_GK
11 #include "get_sonic_fields.hpp"
12 #endif
13 #ifdef EXPLICIT_EM
14 #include "globals.hpp"
15 #endif
16 
17 // Fortran calls
18 #ifdef XGC1
20 # ifdef EXPLICIT_EM
22 # endif
23  );
24 #else
26 # ifdef SONIC_GK
30 # endif
31  );
32 #endif
33 
34 // Electric field class
35 template<class Device>
36 struct ElectricField {
37 
38  // Dimensions
39  int nnode;
40  int nphi;
41 
42  // Gyromatrix info (to be moved soon)
43  int nrho;
44  double inv_drho;
45 
46  // Optimization options
48 
49  // Field variables
50  bool turb_efield;
51  bool E00_efield;
52 
53  // Host Views
56 
57  //For EM
66 
67  //Kokkos::View<double*,Kokkos::LayoutRight,HostType> dpot_n0_h;
68 
71 #ifdef XGC1
73 #else
74 # ifdef SONIC_GK
78 # endif
79 #endif
82 
83  // Default constructor
85 
86  // Makes sense to use the Grid in the constructor here, but the Grid still copies from fortran so until that changes,
87  // it's simpler to not use it here so that ElectricField can be defined before the restart, which uses psn%As
88  ElectricField(NLReader::NamelistReader& nlr, int nnode_in, int nplanes_in, int nrho_in, double inv_drho_in) : //Grid<HostTypeType>& grid) :
89  /*nnode(grid.nnode),
90  nphi(grid.nplanes),
91  nrho(grid.nrho),*/
92  nnode(nnode_in),
93  nphi(nplanes_in),
94  nrho(nrho_in),
95  inv_drho(inv_drho_in),
96  // Host Views
97 #ifdef XGC1
98  pot("pot",nphi,nrho,nnode),
99  E("E",nphi,nrho,nnode),
100 # ifdef DELTAF_CONV
101  E00_ff_h("E00_ff_h",nnode),
102  ddpotdt("ddpotdt",nphi,nrho,nnode),
103 # endif
104 # ifdef EXPLICIT_EM
105  dAh("dAh",nphi,nrho,nnode),
106  Ah("Ah",nphi,nrho,nnode),
107  Ah_h("Ah_h",4,nnode), // -1:2 in fortran
108  As_h("As_h",4,nnode), // -1:2 in fortran
109 # endif
110  dpot_h("dpot_h",4,nnode), // -1:2 in fortran
111  dpot_ff_h("dpot_ff_h",nnode),
112 #else
113  E("E",nrho, nnode),
114  dpot_h("dpot_h",1,nnode),
115 # ifdef SONIC_GK
116  dEr_B2("dEr_B2",nrho, nnode),
117  dEz_B2("dEz_B2",nrho, nnode),
118  du2_E("du2_E",nrho, nnode),
119 # endif
120 #endif
121  pot0_h("pot0_h",nnode)
122  {
123  nlr.use_namelist("sml_param");
124  turb_efield = nlr.get<bool>("sml_turb_efield",true);
125  E00_efield = nlr.get<bool>("sml_00_efield",true);
126 
127 #ifdef EXPLICIT_EM
128  // Pass these (or sml) as arguments
129  bool em_mixed_variable = nlr.get<bool>("sml_em_mixed_variable",true);
130  bool em_control_variate = nlr.get<bool>("sml_em_control_variate",false);
131  int em_pullback_mode = nlr.get<int>("sml_em_pullback_mode",3);
132 
133  if (!em_mixed_variable)
134  exit_XGC("\nError: sml_em_mixed_variable is explicitly set to .false. which is incompatible with EXPLICIT_EM\n");
135 #endif
136 
137  if(nlr.namelist_present("performance_param")){
138  nlr.use_namelist("performance_param");
139  calculate_phi_ff_on_device = nlr.get<bool>("calculate_phi_ff_on_device", false);
140  }else{
142  }
143 
144 #ifdef EXPLICIT_EM
145  if (em_control_variate){
147  }
148 
149  if (em_mixed_variable){
152  if (em_pullback_mode==4){
154  }
155  }
156 #endif
157 
158  // Set pointers in Fortran to rho arrays
159 #ifdef XGC1
160  set_rho_ff_pointers(nrho, nnode, pot0_h.data(), dpot_h.data(), dpot_ff_h.data(), pot.data(), E.data()
161 # ifdef EXPLICIT_EM
162  , Ah_h.data(), Ah.data(), dAh.data(), As_h.data(), As.data(), dAs.data(), Epar_em.data()
163 # endif
164  );
165 #else
166  set_rho_pointers(nrho, nnode, pot0_h.data(), dpot_h.data(), E.data()
167 # ifdef SONIC_GK
168  , dEr_B2.data()
169  , dEz_B2.data()
170  , du2_E.data()
171 # endif
172  );
173 #endif
174  }
175 
176 
177  /**************************************************/
178  /************ Copy to GridFieldPack *************/
179  /**************************************************/
180 
181  // Should these functions be in GridFieldPack instead? - ALS
182 
183 
184  // Shallow copy if DeviceType and HostType are the same
185  template<typename T>
186  void copy_to_pack(T& device_grid_field, const T& host_grid_field) const{
187  device_grid_field = host_grid_field;
188  }
189 
190  // Allocate device grid_field and copy data from host
191  template<typename T1, typename T2>
192  void copy_to_pack(T1& device_grid_field, const T2& host_grid_field) const{
193  // Allocate device array
194  device_grid_field = T1("grid_field",host_grid_field.nphi(), host_grid_field.nrhop1()-1, host_grid_field.nnode());
195  // Deep copy from host to device
196  Kokkos::deep_copy(device_grid_field.f, host_grid_field.f);
197  }
198 
199  // Resizes allocations for large arrays that have a phi dimension (i.e used for electron push)
200  template<PhiInterpType PIT>
202  // Re-initialize new gridpack - deallocates whatever Views were already allocated
203  GridFieldPack<Device, PIT> gfpack(kintype_in, turb_efield, nrho, inv_drho);
204 #ifdef XGC1
205  // Copy field data to device (or shallow copy if no device)
206  if(kintype_in==DriftKin){
207  GPTLstart("GAT_FIELD_INFO");
208  if(pol_decomp.decompose_fields){
209  if(near_field){
210  // Gathering the "near field", defined as the field decomposition region that contains this rank's
211  // global domain decomposition
212  int local_pid = pol_decomp.field_decomp.find_domain_owner(pol_decomp.plane_index, grid.nplanes, pol_decomp.node_offset, grid.nnode);
213  gfpack.phi_offset = pol_decomp.field_decomp.all_first_plane(local_pid);
214  gfpack.node_offset = pol_decomp.field_decomp.all_first_node(local_pid);
215  }else{
216  gfpack.phi_offset = pol_decomp.field_decomp.first_plane;
217  gfpack.node_offset = pol_decomp.field_decomp.first_node;
218  }
219  }
220  PlaneFieldGatherer pfg(pol_decomp, grid, near_field);
221 
222 # ifdef DELTAF_CONV
223  if(ddpotdt.f.size()>0) pfg.gather_phi_ff_on_device(ddpotdt.f, gfpack.ddpotdt_phi_ff.f);
224  copy_to_pack(gfpack.E00_ff, E00_ff_h); // No gather_phi needed due to axisymmetry
225 # endif
226 
228  // Whether to ignore dpot
229  bool ignore_poloidal_dpot = !(turb_efield);
230  bool use_field00 = E00_efield;
231  bool potential_is_requested = false; // Not used in push
233  pfg.calculate_phi_ff_on_device(sml,grid, magnetic_field, my_subview(dpot_h.f, 2), pot0_h, gfpack.E_phi_ff, unused_pot, potential_is_requested, use_field00, ignore_poloidal_dpot);
234 # ifdef EXPLICIT_EM
236  pfg.calculate_phi_ff_on_device(sml,grid, magnetic_field, my_subview(Ah_h.f, 2), pot0_em_h, gfpack.dAh_phi_ff, gfpack.Ah_phi_ff);
237  if(As.f.size()>0) pfg.calculate_phi_ff_on_device(sml,grid, magnetic_field, my_subview(As_h.f, 2), pot0_em_h, gfpack.dAs_phi_ff, gfpack.As_phi_ff);
238 # endif
239  }else{
240  pfg.gather_phi_ff_on_device(E.f, gfpack.E_phi_ff.f);
241 # ifdef EXPLICIT_EM
242  pfg.gather_phi_ff_on_device(dAh.f, gfpack.dAh_phi_ff.f);
243  pfg.gather_phi_ff_on_device(Ah.f, gfpack.Ah_phi_ff.f);
244  if(dAs.f.size()>0) pfg.gather_phi_ff_on_device(dAs.f, gfpack.dAs_phi_ff.f);
245  if(As.f.size()>0) pfg.gather_phi_ff_on_device(As.f, gfpack.As_phi_ff.f);
246 # endif
247  }
248 # ifdef EXPLICIT_EM
249  if(Epar_em.f.size()>0) pfg.gather_phi_ff_on_device(Epar_em.f, gfpack.Epar_em_phi_ff.f);
250 # endif
251  GPTLstop("GAT_FIELD_INFO");
252  }else{
253  GPTLstart("Copy_rho_ff_fields_to_device");
254  copy_to_pack(gfpack.E_rho_ff, E);
255 # ifdef EXPLICIT_EM
256  copy_to_pack(gfpack.dAs_rho_ff, dAs);
257  copy_to_pack(gfpack.dAh_rho_ff, dAh);
258  copy_to_pack(gfpack.As_rho_ff, As);
259  copy_to_pack(gfpack.Ah_rho_ff, Ah);
260  copy_to_pack(gfpack.Epar_em_rho_ff, Epar_em);
261 # endif
262  GPTLstop("Copy_rho_ff_fields_to_device");
263  }
264 #else
265  copy_to_pack(gfpack.E_rho, E);
266 # ifdef SONIC_GK
267  copy_to_pack(gfpack.dEr_B2_rho, dEr_B2);
268  copy_to_pack(gfpack.dEz_B2_rho, dEz_B2);
269  copy_to_pack(gfpack.du2_E_rho, du2_E);
270 # endif
271 #endif
272  return gfpack;
273  }
274 
275  /* Set of functions used to allocate/deallocate gpu memory for Ah during particle scatters */
276  // Resizes allocation of Ah_cv
277  template<PhiInterpType PIT>
279  GridFieldPack<Device, PIT> gfpack(kintype_in, turb_efield, nrho, inv_drho);
280 #ifdef EXPLICIT_EM
281  PlaneFieldGatherer pfg(pol_decomp, grid);
282  pfg.gather_phi_ff_on_device(Ah_cv.f, gfpack.Ah_phi_ff.f);
283 #endif
284  return gfpack;
285  }
286 
287  template<PhiInterpType PIT>
289  GridFieldPack<Device, PIT> gfpack(kintype_in, turb_efield, nrho, inv_drho);
290 #ifdef EXPLICIT_EM
291  if(kintype_in==DriftKin){
292  PlaneFieldGatherer pfg(pol_decomp, grid);
293  pfg.gather_phi_ff_on_device(Ah.f, gfpack.Ah_phi_ff.f);
294  }else{
295  copy_to_pack(gfpack.Ah_rho_ff, Ah);
296  }
297 #endif
298  return gfpack;
299  }
300 
301 
302  // Resizes allocation of dpot_ff and copies to device
303  template<PhiInterpType PIT>
305  GridFieldPack<Device, PIT> gfpack(kintype_in, turb_efield, nrho, inv_drho);
306 #ifdef XGC1
307  copy_to_pack(gfpack.dpot_ff, dpot_ff_h);
308 #else
309  const int PLANE_ZERO=0;
310  copy_to_pack(gfpack.dpot, dpot_h.subfield(PLANE_ZERO));
311 #endif
312  return gfpack;
313  }
314 
326  GPTLstart("GET_POT_GRAD_EXCL_DESTR"); // Timer excludes the destructors which are non-negligible
327 
328  // Allocate temporary views
329  GPTLstart("GET_POT_GRAD_ALLOC_VIEWS");
330  int n_input_field_planes = dpot_h.nphi();
331  bool gyroaverage_requested = true;
332  GetPotentialGradTemp<DeviceType, HostType> tmp(sml, grid, pol_decomp, magnetic_field, n_input_field_planes, gyroaverage_requested, species.gyro_avg_matrices);
333  GPTLstop("GET_POT_GRAD_ALLOC_VIEWS");
334 
335  /****************************************/
336  // Calculation of the electric potential E
337  /****************************************/
338  GPTLstart("GET_POT_GRAD_PHI");
339 
341  (dpot_h.unmanaged(), !(turb_efield));
342 
343  // Set up E00
344  if(E00_efield){
345  gpg_field_args.field00 = Field00<DeviceType>(pot0_h, sml.grad_psitheta);
346  gpg_field_args.field00.calculate_gradient(grid, tmp.grad_matrices);
347  }
348  // Write results of E00
349  if(sml.reduced_deltaf && sml.electron_on) gpg_field_args.field00.set_field(grid, tmp.ff, E00_ff_h.f);
350 
351  // Set up EMPar
352  if(sml.explicit_electromagnetic && (sml.em_mixed_variable && sml.em_pullback_mode==4)) gpg_field_args.request_para_em(Epar_em);
353 
354  // Request potential and/or gradient outputs
355  // XGC1 algorithm uses potential but XGCa doesnt
356  if(!sml.is_XGCa) gpg_field_args.request_potential(pot);
357  gpg_field_args.request_gradient(E);
358 
359  // Get potentials
360  get_field_grad(grid, pol_decomp, magnetic_field, gpg_field_args, tmp);
361 
362  GPTLstop("GET_POT_GRAD_PHI");
363 
364  if(sml.explicit_electromagnetic){
365  GPTLstart("GET_POT_GRAD_APAR");
366 
367  /****************************************/
368  // Calculation of the Hamiltonian vector potential A_h
369  /****************************************/
371  (Ah_h.unmanaged());
372  gpg_field_args.request_potential(Ah);
373  gpg_field_args.request_gradient(dAh);
374  get_field_grad(grid, pol_decomp, magnetic_field, gpg_field_args, tmp);
375 
376  /****************************************/
377  // Calculation of the symplectic vector potential A_s
378  /****************************************/
379  if (sml.em_mixed_variable){
381  (As_h.unmanaged());
382  gpg_field_args.request_potential(As);
383  gpg_field_args.request_gradient(dAs);
384  get_field_grad(grid, pol_decomp, magnetic_field, gpg_field_args, tmp);
385  }
386 
387  /****************************************/
388  // Calculation of the Hamiltonian vector potential A_h for the control-variate method
389  /****************************************/
390  if (sml.em_control_variate){
391  TIMER("GET_POT_APAR_CV",
392  get_field_Ah_cv_ff(Ah_cv.data()) );
393  }
394  GPTLstop("GET_POT_GRAD_APAR");
395  }
396 
397 #ifdef SONIC_GK
398  GPTLstart("GET_POT_GRAD_SONIC");
399  View<double*,CLayout,HostType> Er_B2_h("Er_B2_h",grid.nnode);
400  View<double*,CLayout,HostType> Ez_B2_h("Ez_B2_h",grid.nnode);
401  View<double*,CLayout,HostType> u2_E_h("u2_E_h",grid.nnode);
402 
403  TIMER("GET_POT_SONIC",
404  get_sonic_fields(grid, E, Er_B2_h, Ez_B2_h, u2_E_h) );
405 
407  gpg_field_args.request_gradient(dEr_B2);
408  get_field_grad(grid, pol_decomp, magnetic_field, gpg_field_args, tmp);
409 
411  gpg_field_args.request_gradient(dEz_B2);
412  get_field_grad(grid, pol_decomp, magnetic_field, gpg_field_args, tmp);
413 
415  gpg_field_args.request_gradient(du2_E);
416  get_field_grad(grid, pol_decomp, magnetic_field, gpg_field_args, tmp);
417 
418  GPTLstop("GET_POT_GRAD_SONIC");
419 #endif
420 
421  GPTLstop("GET_POT_GRAD_EXCL_DESTR");
422  }
423 
424 };
425 
426 #endif
bool decompose_fields
Whether to decompose fields.
Definition: domain_decomposition.hpp:35
GridField< HostType, VarType::Scalar, PIT_GLOBAL, TorType::OnePlane, KinType::GyroKin > pot
Definition: electric_field.hpp:69
static int GPTLstart(const char *name)
Definition: timer_macro.hpp:9
int first_plane
First plane belonging to this rank, including ghost planes.
Definition: field_decomposition.hpp:35
void copy_to_pack(T1 &device_grid_field, const T2 &host_grid_field) const
Definition: electric_field.hpp:192
void get_potential_grad(const Simulation< DeviceType > &sml, const Grid< DeviceType > &grid, const DomainDecomposition< DeviceType > &pol_decomp, const MagneticField< DeviceType > &magnetic_field, const Species< DeviceType > &species)
Definition: electric_field.hpp:325
const bool explicit_electromagnetic
Equivalent to the preprocessor flag for now.
Definition: sml.hpp:29
GridField< HostType, vec2d_if_axisym< PIT_GLOBAL >), PIT_GLOBAL, TorType::OnePlane, KinType::GyroKin > dAh
Definition: electric_field.hpp:59
T get(const string &param, const T default_val, int val_ind=0)
Definition: NamelistReader.hpp:372
GridField< HostType, VarType::Scalar, PIT_GLOBAL, TorType::OnePlane, KinType::GyroKin > Ah
Definition: electric_field.hpp:60
bool em_mixed_variable
Switch for use of mixed-variable formulation.
Definition: sml.hpp:68
GridField< HostType, VarType::Vector2D, PIT_GLOBAL, TorType::OnePlane, KinType::DriftKin > E00_ff_h
Radial electric field from &lt;phi&gt; in field-following format.
Definition: electric_field.hpp:54
double inv_drho
Definition: electric_field.hpp:44
Definition: sml.hpp:8
GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::MultiplePlanes, KinType::DriftKin > As_h
Definition: electric_field.hpp:61
int plane_index
Offset of local plane.
Definition: domain_decomposition.hpp:33
Definition: globals.hpp:82
Definition: plane_field_gatherer.hpp:9
FieldFollowingCoordinates ff
Definition: get_potential_grad.hpp:76
KOKKOS_INLINE_FUNCTION int find_domain_owner(int global_plane_index, int nplanes_total, int global_node_index, int nnodes_total) const
Definition: field_decomposition.hpp:136
int nrho
Definition: electric_field.hpp:43
Definition: NamelistReader.hpp:193
ElectricField(NLReader::NamelistReader &nlr, int nnode_in, int nplanes_in, int nrho_in, double inv_drho_in)
Definition: electric_field.hpp:88
Definition: magnetic_field.hpp:12
Definition: get_potential_grad.hpp:248
int phi_offset
Offset for phi_ff field decomposition.
Definition: grid_field_pack.hpp:29
Definition: electric_field.hpp:36
GridField< Device, vec2d_if_axisym< PIT >), PIT, TorType::OnePlane, KinType::GyroKin > E_rho
Definition: grid_field_pack.hpp:58
GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::OnePlane, KinType::DriftKin > pot0_h
Definition: electric_field.hpp:81
int nplanes
Number of planes.
Definition: grid.hpp:217
GridFieldPack< Device, PIT > copy_Ah_ff_to_device(KinType kintype_in, const DomainDecomposition< DeviceType > &pol_decomp, const Grid< DeviceType > &grid) const
Definition: electric_field.hpp:288
Definition: grid_field_pack.hpp:20
int node_offset
Offset of first mesh node belonging to this MPI rank.
Definition: domain_decomposition.hpp:39
GradientMatrices< DeviceType > grad_matrices
Definition: get_potential_grad.hpp:75
void calculate_phi_ff_on_device(const Simulation< DeviceType > &sml, const Grid< DeviceType > &grid, const MagneticField< DeviceType > &magnetic_field, const View< Field< VarType::Scalar, PhiInterpType::None > *, CLayout, HostType > &dpot_h, const GT0 &pot0_h, GT &phi_ff, GT2 &pot_phi_ff, bool potential_is_requested=true, bool use_field00=false, bool ignore_poloidal_dpot=false)
Definition: plane_field_gatherer.hpp:117
void get_field_grad(const Grid< DeviceType > &grid, const DomainDecomposition< DeviceType > &pol_decomp, const MagneticField< DeviceType > &magnetic_field, GetPotGradFieldArgs< DeviceIn, DeviceOut, VT, PIT, TT, KT > &args, GetPotentialGradTemp< DeviceType, DeviceOut > &tmp)
Definition: get_potential_grad.cpp:428
GridField< HostType, VarType::Scalar, PIT_GLOBAL, TorType::OnePlane, KinType::GyroKin > Epar_em
Definition: electric_field.hpp:65
int em_pullback_mode
Definition: sml.hpp:70
bool grad_psitheta
Definition: sml.hpp:84
void get_field_Ah_cv_ff(Field< VarType::Scalar, PIT_GLOBAL > *psn_Ah_cv_ff)
ElectricField()
Definition: electric_field.hpp:84
#define TIMER(N, F)
Definition: timer_macro.hpp:24
void copy_to_pack(T &device_grid_field, const T &host_grid_field) const
Definition: electric_field.hpp:186
GridField< HostType, VarType::Scalar, PIT_GLOBAL, TorType::OnePlane, KinType::DriftKin > ddpotdt
Time derivative of phi-&lt;phi&gt; - Not field-following?
Definition: electric_field.hpp:55
void use_namelist(const string &namelist)
Definition: NamelistReader.hpp:354
void set_rho_pointers(int nrho, int n, Field< VarType::Scalar, PhiInterpType::None > *pot0, Field< VarType::Scalar, PhiInterpType::None > *dpot, Field< VarType::Vector2D, PhiInterpType::None > *E_rho)
constexpr PhiInterpType PIT_GLOBAL
Definition: globals.hpp:98
GridField< HostType, VarType::Scalar, PIT_GLOBAL, TorType::OnePlane, KinType::DriftKin > Ah_cv
Definition: electric_field.hpp:64
void gather_phi_ff_on_device(View< FT **, CLayout, HostType > &tmp, View< FT *, CLayout, HostType > &tmp_full, const T_h &rho_ff_h, View< FT **, CLayout, DeviceType > &phi_ff)
Definition: plane_field_gatherer.hpp:203
Definition: globals.hpp:83
GridFieldPack< Device, PIT > copy_push_fields_to_device(KinType kintype_in, const Simulation< DeviceType > &sml, const DomainDecomposition< DeviceType > &pol_decomp, const MagneticField< DeviceType > &magnetic_field, const Grid< DeviceType > &grid, bool near_field=false) const
Definition: electric_field.hpp:201
GridField< HostType, vec2d_if_axisym< PIT_GLOBAL >), PIT_GLOBAL, TorType::OnePlane, KinType::GyroKin > dAs
Definition: electric_field.hpp:62
GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::MultiplePlanes, KinType::DriftKin > Ah_h
Definition: electric_field.hpp:58
bool E00_efield
Flux-surface averaged potential not used for calculating the electric field if .false.
Definition: electric_field.hpp:51
int node_offset
Offset for phi_ff field decomposition.
Definition: grid_field_pack.hpp:30
KinType
Definition: globals.hpp:81
GridField< Device, VarType::Scalar, PIT, TorType::OnePlane, KinType::DriftKin > dpot
Definition: grid_field_pack.hpp:59
View< int *, CLayout, HostType > all_first_plane
First plane of each rank.
Definition: field_decomposition.hpp:43
bool namelist_present(const string &namelist)
Definition: NamelistReader.hpp:350
bool electron_on
Use kinetic electrons.
Definition: sml.hpp:54
const bool reduced_deltaf
Equivalent to the preprocessor flag for now.
Definition: sml.hpp:23
void get_sonic_fields(const Grid< DeviceType > &grid, const GridField< HostType, vec2d_if_axisym< PIT_GLOBAL >(), PIT_GLOBAL, TorType::OnePlane, KinType::GyroKin > &E, View< double *, CLayout, HostType > &Er_B2_h, View< double *, CLayout, HostType > &Ez_B2_h, View< double *, CLayout, HostType > &u2_E_h)
bool calculate_phi_ff_on_device
Definition: electric_field.hpp:47
void exit_XGC(std::string msg)
Definition: globals.hpp:36
Kokkos::View< T *, Kokkos::LayoutRight, Device > my_subview(const Kokkos::View< T ****, Kokkos::LayoutRight, Device > &view, int i, int j, int k)
Definition: my_subview.hpp:8
Definition: magnetic_field.F90:1
bool turb_efield
E-field calculated only with if .false., psndpot will still contain all (n=0,|m|&gt;0) and (|n|&gt;0...
Definition: electric_field.hpp:50
int nphi
Definition: electric_field.hpp:40
const bool is_XGCa
Equivalent to the preprocessor flag for now.
Definition: sml.hpp:15
GridField< HostType, VarType::Scalar, PIT_GLOBAL, TorType::OnePlane, KinType::GyroKin > As
Definition: electric_field.hpp:63
int nnode
Definition: electric_field.hpp:39
FieldDecomposition< Device > field_decomp
Definition: domain_decomposition.hpp:36
GyroAverageMatrices< HostType > gyro_avg_matrices
Definition: species.hpp:126
Definition: species.hpp:72
GridFieldPack< Device, PIT > copy_dpot_to_device(KinType kintype_in) const
Definition: electric_field.hpp:304
GridFieldPack< Device, PIT > copy_Ah_cv_to_device_as_Ah(KinType kintype_in, const DomainDecomposition< DeviceType > &pol_decomp, const Grid< DeviceType > &grid) const
Definition: electric_field.hpp:278
void calculate_gradient(const Grid< DeviceType > &grid, const GradientMatrices< DeviceType > &grad_matrices)
Definition: get_potential_grad.hpp:143
GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::MultiplePlanes, KinType::DriftKin > dpot_h
Definition: electric_field.hpp:80
int nnode
Number of grid nodes.
Definition: grid.hpp:216
View< int *, CLayout, HostType > all_first_node
First node of each rank.
Definition: field_decomposition.hpp:41
bool em_control_variate
Switch for use of control variate method.
Definition: sml.hpp:69
int first_node
First mesh node belonging to this rank, including ghost nodes.
Definition: field_decomposition.hpp:32
GridField< HostType, vec2d_if_axisym< PIT_GLOBAL >), PIT_GLOBAL, TorType::OnePlane, KinType::GyroKin > E
Definition: electric_field.hpp:70
static int GPTLstop(const char *name)
Definition: timer_macro.hpp:10
Definition: get_potential_grad.hpp:54