XGC1
 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 "smoothing.hpp"
9 #include "get_potential_grad.hpp"
10 #include "plane_field_gatherer.hpp"
12 #ifdef SONIC_GK
13 #include "get_sonic_fields.hpp"
14 #endif
15 #ifdef EXPLICIT_EM
16 #include "globals.hpp"
17 #endif
18 
19 // Fortran calls
20 #ifdef XGC1
22 # ifdef EXPLICIT_EM
24  double* dpot_es
25 # endif
26  );
27 #else
29 # ifdef SONIC_GK
33 # endif
34  );
35 #endif
36 
37 // Electric field class
38 template<class Device>
39 struct ElectricField {
40 
41  // Dimensions
42  int nnode;
43  int nphi;
44 
45  // Gyromatrix info (to be moved soon)
46  int nrho;
47 
48  // Optimization options
50 
51  // Field variables
52  bool turb_efield;
53  bool E00_efield;
54 
55  // Host Views
58 
59  //For EM
70 
73 //#ifdef XGC1
75 //#else
76 //# ifdef SONIC_GK
80 //# endif
81 //#endif
85  View<double***,CLayout, HostType> save_dpot0; // Used for adiabatic response
86  View<double***,CLayout, HostType> save_dpot; // Used for adiabatic response
87  View<double***,CLayout, HostType> dpotsave; // Used for split-weight scheme
88  View<double*,CLayout, HostType> dpot_es; // Used for a diagnostic
89 
90  // Host View for the adjustable loop voltage current drive
94 
95  // Default constructor
97 
98  // Makes sense to use the Grid in the constructor here, but the Grid still copies from fortran so until that changes,
99  // it's simpler to not use it here so that ElectricField can be defined before the restart, which uses psn%As
100  ElectricField(NLReader::NamelistReader& nlr, int nnode_in, int nplanes_in, int nrho_in) : //Grid<HostTypeType>& grid) :
101  /*nnode(grid.nnode),
102  nphi(grid.nplanes),
103  nrho(grid.nrho),*/
104  nnode(nnode_in),
105  nphi(nplanes_in),
106  nrho(nrho_in),
107  // Host Views
108 #ifdef XGC1
109  pot("pot",nphi,nrho,nnode),
110  E("E",nphi,nrho,nnode),
111 # if defined(DELTAF_CONV) && !defined(EXPLICIT_EM)
112  E00_ff_h("E00_ff_h",nnode),
113  ddpotdt("ddpotdt",nphi,nrho,nnode),
114  dpotsave("dpotsave", 2, 2, nnode), // View. dims: RK2, ff, nnode
115 # endif
116  save_dpot("save_dpot", 2, 2, nnode), // View. dims: RK2, ff, nnode
117  save_dpot0("save_dpot0", 2, 2, nnode), // View. dims: RK2, ff, nnode
118 # ifdef EXPLICIT_EM
119  dpot_es("dpot_es", nnode),
120  dAh("dAh",nphi,nrho,nnode),
121  Ah("Ah",nphi,nrho,nnode),
122  Ah_h("Ah_h",4,nnode), // -1:2 in fortran
123  As_h("As_h",4,nnode), // -1:2 in fortran
124  As_h_backup("As_h_backup",4,nnode), // -1:2 in fortran
125 # endif
126  dpot_h("dpot_h",4,nnode), // -1:2 in fortran
127  dpot_ff_h("dpot_ff_h",nnode),
128  dpot_n0_h("dpot_n0_h",nnode),
129 #else
130  E("E",nrho, nnode),
131  dpot_h("dpot_h",1,nnode),
132 # ifdef SONIC_GK
133  dEr_B2("dEr_B2",nrho, nnode),
134  dEz_B2("dEz_B2",nrho, nnode),
135  du2_E("du2_E",nrho, nnode),
136 # endif
137  dpot_n0_h("dpot_n0_h",nnode),
138 #endif
139  pot0_h("pot0_h",nnode),
140  loop_voltage_h("loop_voltage_h",nnode),
141  jpara_diff_previous_h("jpara_diff_previous_h",nnode),
142  jpara_diff_integral_h("jpara_diff_integral_h",nnode)
143  {
144  nlr.use_namelist("sml_param");
145  turb_efield = nlr.get<bool>("sml_turb_efield",true);
146  E00_efield = nlr.get<bool>("sml_00_efield",true);
147 
148 #ifdef EXPLICIT_EM
149  // Pass these (or sml) as arguments
150  bool em_mixed_variable = nlr.get<bool>("sml_em_mixed_variable",true);
151  bool em_control_variate = nlr.get<bool>("sml_em_control_variate",false);
152  int em_pullback_mode = nlr.get<int>("sml_em_pullback_mode",3);
153 
154  if (!em_mixed_variable)
155  exit_XGC("\nError: sml_em_mixed_variable is explicitly set to .false. which is incompatible with EXPLICIT_EM\n");
156 #endif
157 
158  if(nlr.namelist_present("performance_param")){
159  nlr.use_namelist("performance_param");
160  calculate_phi_ff_on_device = nlr.get<bool>("calculate_phi_ff_on_device", false);
161  }else{
163  }
164 
165 #ifdef EXPLICIT_EM
166  if (em_control_variate){
169  }
170 
171  if (em_mixed_variable){
174  if (em_pullback_mode==4){
176  }
177  }
178 #endif
179 
180  // Set pointers in Fortran to rho arrays
181 #ifdef XGC1
182  set_rho_ff_pointers(nrho, nnode, pot0_h.data(), dpot_n0_h.data(), dpot_h.data(), dpot_ff_h.data(), pot.data(), E.data()
183 # ifdef EXPLICIT_EM
184  , Ah_h.data(), Ah.data(), dAh.data(), As_h.data(), As_h_backup.data(), As.data(), dAs.data(), Epar_em.data(), dpot_es.data()
185 # endif
186  );
187 #else
188  set_rho_pointers(nrho, nnode, pot0_h.data(), dpot_h.data(), E.data()
189 # ifdef SONIC_GK
190  , dEr_B2.data()
191  , dEz_B2.data()
192  , du2_E.data()
193 # endif
194  );
195 #endif
196  // Initialize loop voltage
197  nlr.use_namelist("sml_param");
198  const double initial_loop_vol = nlr.get<double>("sml_loop_voltage",0.0);
199  //
200  for (int inode=0; inode<nnode; inode++){
201  loop_voltage_h(inode) = initial_loop_vol;
202  }
203  }
204 
205 
206  /**************************************************/
207  /************ Copy to GridFieldPack *************/
208  /**************************************************/
209 
210  // Resizes allocations for large arrays that have a phi dimension (i.e used for electron push)
211  template<PhiInterpType PIT>
213  // Re-initialize new gridpack - deallocates whatever Views were already allocated
214  GridFieldPack<Device, PIT> gfpack(kintype_in, turb_efield);
215  //
216  // Copy axisymmetric loop voltage to temporary device-view
217  //
218  ScalarGridField loop_voltage("loop_voltage",loop_voltage_h.nnode());
219  Kokkos::deep_copy(loop_voltage.f,loop_voltage_h.f);
220  //
221 #ifdef XGC1
222  // Copy field data to device (or shallow copy if no device)
223  if(kintype_in==DriftKin){
224  GPTLstart("GAT_FIELD_INFO");
225  if(pol_decomp.decompose_fields){
226  if(near_field){
227  // Gathering the "near field", defined as the field decomposition region that contains this rank's
228  // global domain decomposition
229  int local_pid = pol_decomp.field_decomp.find_domain_owner(pol_decomp.plane_index, grid.nplanes, pol_decomp.node_offset, grid.nnode);
230  gfpack.phi_offset = pol_decomp.field_decomp.all_first_plane(local_pid);
231  gfpack.node_offset = pol_decomp.field_decomp.all_first_node(local_pid);
232  }else{
233  gfpack.phi_offset = pol_decomp.field_decomp.first_plane;
234  gfpack.node_offset = pol_decomp.field_decomp.first_node;
235  }
236  }
237  PlaneFieldGatherer pfg(pol_decomp, grid, near_field);
238 
239 # if defined(DELTAF_CONV) && !defined(EXPLICIT_EM)
240  if(ddpotdt.f.size()>0) pfg.gather_phi_ff_on_device(ddpotdt.f, gfpack.ddpotdt_phi_ff.f);
241  grid_field_copy(gfpack.E00_ff, E00_ff_h); // No gather_phi needed due to axisymmetry
242 # endif
243 
245  // Allocate temporary views
246  GPTLstart("GET_POT_GRAD_ALLOC_VIEWS");
247  int n_input_field_planes = pot0_h.nphi();
248  bool gradient_requested = true;
249  bool gyroaverage_requested = false;
250  GetPotentialGradTemp<DeviceType, DeviceType> gpg_tmp(sml, grid, pol_decomp, magnetic_field, n_input_field_planes, gradient_requested, gyroaverage_requested);
251  GPTLstop("GET_POT_GRAD_ALLOC_VIEWS");
252 
253  // Whether to ignore dpot
254  bool ignore_poloidal_dpot = !(turb_efield);
255  bool use_field00 = E00_efield;
256  bool potential_is_requested = false; // Not used in push
258  pfg.calculate_phi_ff_on_device(sml,grid, magnetic_field, smoothing, gpg_tmp, my_subview(dpot_h.f, 2), pot0_h, gfpack.E_phi_ff, unused_pot, potential_is_requested, use_field00, ignore_poloidal_dpot);
259 # ifdef EXPLICIT_EM
261  pfg.calculate_phi_ff_on_device(sml,grid, magnetic_field, smoothing, gpg_tmp, my_subview(Ah_h.f, 2), pot0_em_h, gfpack.dAh_phi_ff, gfpack.Ah_phi_ff);
262  if(As.f.size()>0) pfg.calculate_phi_ff_on_device(sml,grid, magnetic_field, smoothing, gpg_tmp, my_subview(As_h.f, 2), pot0_em_h, gfpack.dAs_phi_ff, gfpack.As_phi_ff);
263 # endif
264  }else{
265  pfg.gather_phi_ff_on_device(E.f, gfpack.E_phi_ff.f);
266 # ifdef EXPLICIT_EM
267  pfg.gather_phi_ff_on_device(dAh.f, gfpack.dAh_phi_ff.f);
268  pfg.gather_phi_ff_on_device(Ah.f, gfpack.Ah_phi_ff.f);
269  if(dAs.f.size()>0) pfg.gather_phi_ff_on_device(dAs.f, gfpack.dAs_phi_ff.f);
270  if(As.f.size()>0) pfg.gather_phi_ff_on_device(As.f, gfpack.As_phi_ff.f);
271 # endif
272  }
273 # ifdef EXPLICIT_EM
274  if(Epar_em.f.size()>0) pfg.gather_phi_ff_on_device(Epar_em.f, gfpack.Epar_em_phi_ff.f);
275 # endif
276  GPTLstop("GAT_FIELD_INFO");
277  }else{
278  GPTLstart("Copy_rho_ff_fields_to_device");
279  grid_field_copy(gfpack.E_rho_ff, E);
280 # ifdef EXPLICIT_EM
281  grid_field_copy(gfpack.dAs_rho_ff, dAs);
282  grid_field_copy(gfpack.dAh_rho_ff, dAh);
283  grid_field_copy(gfpack.As_rho_ff, As);
284  grid_field_copy(gfpack.Ah_rho_ff, Ah);
285  grid_field_copy(gfpack.Epar_em_rho_ff, Epar_em);
286 # endif
287  GPTLstop("Copy_rho_ff_fields_to_device");
288  }
289  //
290  // Loop voltage is identical for gyrokinetic and drift-kinetic species,
291  // but different between XGC1 and XGCa
292  //
294  gfpack.loop_voltage = GridField<Device, VarType::Scalar,PIT, TorType::OnePlane, KinType::DriftKin>("gfpack_loop_voltage",loop_voltage.nnode());
295  ff.cnvt_grid_real2ff(grid,loop_voltage,gfpack.loop_voltage);
296 #else
297  grid_field_copy(gfpack.E_rho, E);
298 # ifdef SONIC_GK
299  grid_field_copy(gfpack.dEr_B2_rho, dEr_B2);
300  grid_field_copy(gfpack.dEz_B2_rho, dEz_B2);
301  grid_field_copy(gfpack.du2_E_rho, du2_E);
302 # endif
303  grid_field_copy(gfpack.loop_voltage,loop_voltage);
304 #endif
305  //grid_field_copy(gfpack.loop_voltage,loop_voltage_h);
306  return gfpack;
307  }
308 
309  /* Set of functions used to allocate/deallocate gpu memory for Ah during particle scatters */
310  // Resizes allocation of Ah_cv
311  template<PhiInterpType PIT>
313  GridFieldPack<Device, PIT> gfpack(kintype_in, turb_efield);
314 #ifdef EXPLICIT_EM
315  PlaneFieldGatherer pfg(pol_decomp, grid);
316  pfg.gather_phi_ff_on_device(Ah_cv.f, gfpack.Ah_phi_ff.f);
317 #endif
318  return gfpack;
319  }
320 
321  template<PhiInterpType PIT>
323  GridFieldPack<Device, PIT> gfpack(kintype_in, turb_efield);
324 #ifdef EXPLICIT_EM
325  if(kintype_in==DriftKin){
326  PlaneFieldGatherer pfg(pol_decomp, grid);
327  pfg.gather_phi_ff_on_device(Ah.f, gfpack.Ah_phi_ff.f);
328  }else{
329  grid_field_copy(gfpack.Ah_rho_ff, Ah);
330  }
331 #endif
332  return gfpack;
333  }
334 
335 
336  // Resizes allocation of dpot_ff and copies to device
337  template<PhiInterpType PIT>
339  GridFieldPack<Device, PIT> gfpack(kintype_in, turb_efield);
340 #ifdef XGC1
342 #else
343  const int PLANE_ZERO=0;
344  grid_field_copy(gfpack.dpot, dpot_h.subfield(PLANE_ZERO));
345 #endif
346  return gfpack;
347  }
348 
359  GPTLstart("GET_AH_RHO_FF_EXCL_DESTR"); // Timer excludes the destructors which are non-negligible
360 
361  // Allocate temporary views
362  GPTLstart("GET_AH_RHO_FF_ALLOC_VIEWS");
363  int n_input_field_planes = Ah_h.nphi();
364  bool gradient_requested = false;
365  bool gyroaverage_requested = true;
366  GetPotentialGradTemp<DeviceType, HostType> tmp(sml, grid, pol_decomp, magnetic_field, n_input_field_planes, gradient_requested, gyroaverage_requested, species.gyro_avg_matrices);
367  GPTLstop("GET_AH_RHO_FF_ALLOC_VIEWS");
368 
369  /****************************************/
370  // Calculation of the Hamiltonian vector potential A_h
371  /****************************************/
373  FieldArgs gpg_field_args(Ah_h.unmanaged());
374  gpg_field_args.request_potential(Ah);
375  get_field_grad(sml, grid, pol_decomp, magnetic_field, smoothing, gpg_field_args, tmp);
376 
377  GPTLstop("GET_AH_RHO_FF_EXCL_DESTR");
378  }
379 
391  GPTLstart("GET_POT_GRAD_EXCL_DESTR"); // Timer excludes the destructors which are non-negligible
392 
393  // Allocate temporary views
394  GPTLstart("GET_POT_GRAD_ALLOC_VIEWS");
395  int n_input_field_planes = dpot_h.nphi();
396  bool gradient_requested = true;
397  bool gyroaverage_requested = true;
398 #ifdef MULTI_RATE
399  gyroaverage_requested = false; // This is a hack. This variable controls the size of the tmp arrays for get_pot_grad,
400  // but the rho_ff arrays will still exist and the rho=0 values will be correct
401 #endif
402  GetPotentialGradTemp<DeviceType, HostType> tmp(sml, grid, pol_decomp, magnetic_field, n_input_field_planes, gradient_requested, gyroaverage_requested, species.gyro_avg_matrices);
403  GPTLstop("GET_POT_GRAD_ALLOC_VIEWS");
404 
405  /****************************************/
406  // Calculation of the electric potential E
407  /****************************************/
408  GPTLstart("GET_POT_GRAD_PHI");
409 
411  FieldArgs gpg_field_args(dpot_h.unmanaged(), !(turb_efield));
412 
413  // Set up E00
414  if(E00_efield){
415  gpg_field_args.field00 = Field00<DeviceType>(pot0_h, sml.grad_psitheta);
416  gpg_field_args.field00.calculate_gradient(grid, tmp.grad_matrices);
417  }
418  // Write results of E00
419  if(sml.reduced_deltaf && sml.electron_on) gpg_field_args.field00.set_field(grid, tmp.ff, E00_ff_h.f);
420 
421  // Set up EMPar
422  if(sml.explicit_electromagnetic && (sml.em_mixed_variable && sml.em_pullback_mode==4)) gpg_field_args.request_para_em(Epar_em, sml.em_pullback_dampfac, gyroaverage_requested);
423 
424  // Request potential and/or gradient outputs
425  // XGC1 algorithm uses potential but XGCa doesnt
426  if(!sml.is_XGCa) gpg_field_args.request_potential(pot);
427  gpg_field_args.request_gradient(E);
428 
429  // Get potentials
430  get_field_grad(sml, grid, pol_decomp, magnetic_field, smoothing, gpg_field_args, tmp);
431 
432  GPTLstop("GET_POT_GRAD_PHI");
433 
434  if(sml.explicit_electromagnetic){
435  GPTLstart("GET_POT_GRAD_APAR");
436 
437  /****************************************/
438  // Calculation of the Hamiltonian vector potential A_h
439  /****************************************/
440  gpg_field_args = FieldArgs(Ah_h.unmanaged());
441  gpg_field_args.request_potential(Ah);
442 #ifndef MULTI_RATE
443  gpg_field_args.request_gradient(dAh);
444 #endif
445  get_field_grad(sml, grid, pol_decomp, magnetic_field, smoothing, gpg_field_args, tmp);
446 
447  /****************************************/
448  // Calculation of the symplectic vector potential A_s
449  /****************************************/
450  if (sml.em_mixed_variable){
451  gpg_field_args = FieldArgs(As_h.unmanaged());
452  gpg_field_args.request_potential(As);
453 #ifndef MULTI_RATE
454  gpg_field_args.request_gradient(dAs);
455 #endif
456  get_field_grad(sml, grid, pol_decomp, magnetic_field, smoothing, gpg_field_args, tmp);
457  }
458 
459  /****************************************/
460  // Calculation of the Hamiltonian vector potential A_h for the control-variate method
461  /****************************************/
462  if (sml.em_control_variate){
463  TIMER("GET_POT_APAR_CV",
464  get_field_Ah_cv_ff(grid, pol_decomp, Ah_cv_h, Ah_cv) );
465  }
466  GPTLstop("GET_POT_GRAD_APAR");
467  }
468 
469 #ifdef SONIC_GK
470  GPTLstart("GET_POT_GRAD_SONIC");
471  View<double*,CLayout,HostType> Er_B2_h("Er_B2_h",grid.nnode);
472  View<double*,CLayout,HostType> Ez_B2_h("Ez_B2_h",grid.nnode);
473  View<double*,CLayout,HostType> u2_E_h("u2_E_h",grid.nnode);
474 
475  TIMER("GET_POT_SONIC",
476  get_sonic_fields(grid, E, Er_B2_h, Ez_B2_h, u2_E_h) );
477 
478  gpg_field_args = FieldArgs(Er_B2_h);
479  gpg_field_args.request_gradient(dEr_B2);
480  get_field_grad(sml, grid, pol_decomp, magnetic_field, smoothing, gpg_field_args, tmp);
481 
482  gpg_field_args = FieldArgs(Ez_B2_h);
483  gpg_field_args.request_gradient(dEz_B2);
484  get_field_grad(sml, grid, pol_decomp, magnetic_field, smoothing, gpg_field_args, tmp);
485 
486  gpg_field_args = FieldArgs(u2_E_h);
487  gpg_field_args.request_gradient(du2_E);
488  get_field_grad(sml, grid, pol_decomp, magnetic_field, smoothing, gpg_field_args, tmp);
489 
490  GPTLstop("GET_POT_GRAD_SONIC");
491 #endif
492 
493  GPTLstop("GET_POT_GRAD_EXCL_DESTR");
494  }
495 
497  nlr.use_namelist("sml_param");
498  int field_n_rho = nlr.get<int>("sml_grid_nrho", 6);
499 
500 #ifdef EXPLICIT_EM
501  int n_fields = 3 + 1; // E, dAs, dAh, + scalars As,Ah
502 #else
503  int n_fields = 1;
504 #endif
505  double host_field_GB_per_vertex = (53 + 8*grid.nplanes + 11*(field_n_rho + 1))*8*BYTES_TO_GB;
506  double cpu_memory_usage = grid.nnode*host_field_GB_per_vertex;
507 
508  double device_field_GB_per_vertex = std::max(n_fields*8*grid.nplanes, n_fields*8*(field_n_rho+1))*8*BYTES_TO_GB;
509  double gpu_memory_usage = grid.nnode*device_field_GB_per_vertex;
510 
511  MemoryPrediction memory_prediction{"Collisions", cpu_memory_usage, gpu_memory_usage};
512 
513  return memory_prediction;
514  }
515 };
516 
517 #endif
GridField< HostType, VarType::Scalar, PhiInterpType::Planes, TorType::OnePlane, KinType::DriftKin > ddpotdt
Time derivative of - Not field-following?
Definition: electric_field.hpp:57
void get_potential_grad(const Simulation< DeviceType > &sml, const Grid< DeviceType > &grid, const DomainDecomposition< DeviceType > &pol_decomp, const MagneticField< DeviceType > &magnetic_field, Smoothing &smoothing, const Species< DeviceType > &species)
Definition: electric_field.hpp:390
bool decompose_fields
Whether to decompose fields.
Definition: domain_decomposition.hpp:51
GridField< HostType, VarType::Scalar, PIT_GLOBAL, TorType::OnePlane, KinType::GyroKin > pot
Definition: electric_field.hpp:71
Definition: field_following_coordinates.hpp:9
static int GPTLstart(const char *name)
Definition: timer_macro.hpp:9
GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::OnePlane, KinType::DriftKin > jpara_diff_integral_h
Definition: electric_field.hpp:93
GridField< Device, VarType::Scalar, PIT, TorType::OnePlane, KinType::DriftKin > loop_voltage
Definition: grid_field_pack.hpp:64
int first_plane
First plane belonging to this rank, including ghost planes.
Definition: field_decomposition.hpp:35
GridField< HostType, vec2d_if_axisym< PIT_GLOBAL >), PIT_GLOBAL, TorType::OnePlane, KinType::GyroKin > dAh
Definition: electric_field.hpp:61
ElectricField(NLReader::NamelistReader &nlr, int nnode_in, int nplanes_in, int nrho_in)
Definition: electric_field.hpp:100
void set_rho_ff_pointers(int nrho, int n, Field< VarType::Scalar, PhiInterpType::None > *pot0, Field< VarType::Scalar, PhiInterpType::None > *dpot_n0, Field< VarType::Scalar, PhiInterpType::None > *dpot, Field< VarType::Scalar, PhiInterpType::Planes > *dpot_ff, Field< VarType::Scalar, PhiInterpType::Planes > *pot_rho_ff, Field< VarType::Vector, PhiInterpType::Planes > *E_rho_ff)
T get(const string &param, const T default_val, int val_ind=0)
Definition: NamelistReader.hpp:373
constexpr double BYTES_TO_GB
Definition: constants.hpp:12
GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::OnePlane, KinType::DriftKin > loop_voltage_h
Definition: electric_field.hpp:91
GridField< HostType, VarType::Scalar, PIT_GLOBAL, TorType::OnePlane, KinType::GyroKin > Ah
Definition: electric_field.hpp:62
bool em_mixed_variable
Switch for use of mixed-variable formulation.
Definition: sml.hpp:78
GridField< HostType, VarType::Vector2D, PIT_GLOBAL, TorType::OnePlane, KinType::DriftKin > E00_ff_h
Radial electric field from in field-following format.
Definition: electric_field.hpp:56
View< double *, CLayout, HostType > dpot_es
Definition: electric_field.hpp:88
GridField< HostType, VarType::Scalar, PhiInterpType::Planes, TorType::OnePlane, KinType::DriftKin > dpot_ff_h
Definition: electric_field.hpp:74
Definition: sml.hpp:8
static constexpr bool is_XGCa
Equivalent to the preprocessor flag for now.
Definition: sml.hpp:19
static constexpr bool explicit_electromagnetic
Equivalent to the preprocessor flag for now.
Definition: sml.hpp:37
GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::MultiplePlanes, KinType::DriftKin > As_h
Definition: electric_field.hpp:63
int plane_index
Offset of local plane.
Definition: domain_decomposition.hpp:48
Definition: globals.hpp:89
Definition: plane_field_gatherer.hpp:9
FieldFollowingCoordinates ff
Definition: get_potential_grad.hpp:75
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:46
Definition: NamelistReader.hpp:193
Definition: magnetic_field.hpp:12
static MemoryPrediction estimate_memory_usage(NLReader::NamelistReader &nlr, const Grid< DeviceType > &grid)
Definition: electric_field.hpp:496
Definition: get_potential_grad.hpp:283
int phi_offset
Offset for phi_ff field decomposition.
Definition: grid_field_pack.hpp:28
void cnvt_grid_real2ff(const Grid< DeviceType > &grid, const View< double **, CLayout, DeviceType > &input, const View< double **, Kokkos::LayoutRight, DeviceType > &output) const
Definition: field_following_coordinates.hpp:121
Definition: electric_field.hpp:39
void get_field_grad(const Simulation< DeviceType > &sml, const Grid< DeviceType > &grid, const DomainDecomposition< DeviceType > &pol_decomp, const MagneticField< DeviceType > &magnetic_field, Smoothing &smoothing, GetPotGradFieldArgs< DeviceIn, DeviceOut, VT, PIT, TT, KT > &args, GetPotentialGradTemp< DeviceType, DeviceOut > &tmp)
Definition: get_potential_grad.cpp:453
GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::OnePlane, KinType::DriftKin > pot0_h
Definition: electric_field.hpp:83
int nplanes
Number of planes.
Definition: grid.hpp:274
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:322
static constexpr bool reduced_deltaf
Equivalent to the preprocessor flag for now.
Definition: sml.hpp:25
Definition: grid_field_pack.hpp:19
int node_offset
Offset of first mesh node belonging to this MPI rank.
Definition: domain_decomposition.hpp:55
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, Smoothing &smoothing, bool near_field=false) const
Definition: electric_field.hpp:212
GridField< Device, vec2d_if_axisym< PIT >), PIT, TorType::MultiplePlanes, KinType::DriftKin > E_phi_ff
Definition: grid_field_pack.hpp:33
GradientMatrices< DeviceType > grad_matrices
Definition: get_potential_grad.hpp:74
Projection< HostType > half_plane_ff
Definition: grid.hpp:318
void request_para_em(const GridField< DeviceOut, VarType::Scalar, PhiInterpType::Planes, TorType::OnePlane, KinType::GyroKin > &output_field, double em_pullback_dampfac_in, bool gyroaverage_requested)
Definition: get_potential_grad.hpp:319
View< double ***, CLayout, HostType > save_dpot
Definition: electric_field.hpp:86
int em_pullback_mode
Definition: sml.hpp:80
bool grad_psitheta
Definition: sml.hpp:101
GridField< HostType, VarType::Scalar, PhiInterpType::Planes, TorType::OnePlane, KinType::DriftKin > Ah_cv
Definition: electric_field.hpp:68
ElectricField()
Definition: electric_field.hpp:96
#define TIMER(N, F)
Definition: timer_macro.hpp:24
GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::MultiplePlanes, KinType::DriftKin > Ah_cv_h
Definition: electric_field.hpp:67
void use_namelist(const string &namelist)
Definition: NamelistReader.hpp:355
constexpr PhiInterpType PIT_GLOBAL
Definition: globals.hpp:103
GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::MultiplePlanes, KinType::DriftKin > As_h_backup
Definition: electric_field.hpp:64
Definition: memory_prediction.hpp:4
void grid_field_copy(T &dest, const T &src)
Definition: grid_field.hpp:440
Definition: smoothing.hpp:8
GridField< HostType, VarType::Vector2D, PIT_GLOBAL, TorType::OnePlane, KinType::GyroKin > dEz_B2
Gradient of Ez/B^2.
Definition: electric_field.hpp:78
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:197
Field00< DeviceType > field00
Definition: get_potential_grad.hpp:290
Definition: globals.hpp:90
GridField< HostType, vec2d_if_axisym< PIT_GLOBAL >), PIT_GLOBAL, TorType::OnePlane, KinType::GyroKin > dAs
Definition: electric_field.hpp:65
void request_potential(const GridField< DeviceOut, VarType::Scalar, PIT, TT, KT > &potential_in)
Definition: get_potential_grad.hpp:309
GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::MultiplePlanes, KinType::DriftKin > Ah_h
Definition: electric_field.hpp:60
GridField< Device, vec2d_if_axisym< PIT >), PIT, TorType::OnePlane, KinType::GyroKin > E_rho_ff
Definition: grid_field_pack.hpp:34
bool E00_efield
Flux-surface averaged potential not used for calculating the electric field if .false.
Definition: electric_field.hpp:53
GridField< HostType, VarType::Vector2D, PIT_GLOBAL, TorType::OnePlane, KinType::GyroKin > dEr_B2
Gradient of Er/B^2.
Definition: electric_field.hpp:77
int node_offset
Offset for phi_ff field decomposition.
Definition: grid_field_pack.hpp:29
GridField< HostType, VarType::Vector2D, PIT_GLOBAL, TorType::OnePlane, KinType::GyroKin > du2_E
Gradient of the square (dot product) of the E x B drift.
Definition: electric_field.hpp:79
KinType
Definition: globals.hpp:88
GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::OnePlane, KinType::DriftKin > jpara_diff_previous_h
Definition: electric_field.hpp:92
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:351
bool electron_on
Use kinetic electrons.
Definition: sml.hpp:64
GridField< HostType, VarType::Scalar, PhiInterpType::Planes, TorType::OnePlane, KinType::GyroKin > Epar_em
Definition: electric_field.hpp:69
GridField< Device, VarType::Scalar, PIT, TorType::OnePlane, KinType::DriftKin > dpot_ff
Definition: grid_field_pack.hpp:35
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:49
void exit_XGC(std::string msg)
Definition: globals.hpp:37
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:52
View< double ***, CLayout, HostType > save_dpot0
Definition: electric_field.hpp:85
int nphi
Definition: electric_field.hpp:43
GridField< HostType, VarType::Scalar, PIT_GLOBAL, TorType::OnePlane, KinType::GyroKin > As
Definition: electric_field.hpp:66
int nnode
Definition: electric_field.hpp:42
FieldDecomposition< Device > field_decomp
Definition: domain_decomposition.hpp:52
GyroAverageMatrices< HostType > gyro_avg_matrices
Definition: species.hpp:129
Definition: species.hpp:75
void request_gradient(const GridField< DeviceOut, VT, PIT, TT, KT > &gradient_in)
Definition: get_potential_grad.hpp:314
void get_field_Ah_cv_ff(const Grid< DeviceType > &grid, const DomainDecomposition< DeviceType > &pol_decomp, const GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::MultiplePlanes, KinType::DriftKin > &Ah_cv_h, const GridField< HostType, VarType::Scalar, PhiInterpType::Planes, TorType::OnePlane, KinType::DriftKin > &Ah_cv_ff)
Definition: get_potential_grad.cpp:9
GridFieldPack< Device, PIT > copy_dpot_to_device(KinType kintype_in) const
Definition: electric_field.hpp:338
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:312
void calculate_gradient(const Grid< DeviceType > &grid, const GradientMatrices< DeviceType > &grad_matrices)
Definition: get_potential_grad.hpp:153
GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::MultiplePlanes, KinType::DriftKin > dpot_h
Definition: electric_field.hpp:82
int nnode
Number of grid nodes.
Definition: grid.hpp:273
double em_pullback_dampfac
Damping term gamma on -b.grad(phi) in pullback mode 4.
Definition: sml.hpp:86
void get_Ah_rho_ff(const Simulation< DeviceType > &sml, const Grid< DeviceType > &grid, const DomainDecomposition< DeviceType > &pol_decomp, const MagneticField< DeviceType > &magnetic_field, Smoothing &smoothing, const Species< DeviceType > &species)
Definition: electric_field.hpp:358
GridField< HostType, VarType::Scalar, PhiInterpType::None, TorType::OnePlane, KinType::DriftKin > dpot_n0_h
Definition: electric_field.hpp:84
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:79
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:72
void set_field(const Grid< DeviceType > &grid, const FieldFollowingCoordinates &ff, View< Field< VarType::Vector2D, PIT > *, CLayout, DeviceOut > &field00_ff_h)
Definition: get_potential_grad.hpp:169
View< double ***, CLayout, HostType > dpotsave
Definition: electric_field.hpp:87
void calculate_phi_ff_on_device(const Simulation< DeviceType > &sml, const Grid< DeviceType > &grid, const MagneticField< DeviceType > &magnetic_field, Smoothing &smoothing, GetPotentialGradTemp< DeviceType, DeviceType > &tmp, 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
static int GPTLstop(const char *name)
Definition: timer_macro.hpp:10
Definition: get_potential_grad.hpp:53