XGCa
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
species.hpp
Go to the documentation of this file.
1 #ifndef SPECIES_HPP
2 #define SPECIES_HPP
3 #include <Cabana_AoSoA.hpp>
4 #include <Cabana_DeepCopy.hpp>
5 #include <Kokkos_Core.hpp>
6 #include "NamelistReader.hpp"
7 #include "timer_macro.hpp"
8 #include "magnetic_field.hpp"
9 #include "grid.hpp"
10 #include "domain_decomposition.hpp"
11 #include "particles.hpp"
12 #include "space_settings.hpp"
13 #include "distribution.hpp"
14 #include "profile.hpp"
16 
17 extern "C" void set_spall_num_and_ptr(int idx, int n_ptl, int n_vecs, VecParticles* ptl);
18 extern "C" void set_min_max_num(int isp, int n_ptl);
19 extern "C" void adjust_n_ptl_for_core_ptl(int* n_ptl);
20 
25  return false;
26 }
27 
32  return false;
33 }
34 
35 // Used for Cabana slices (getting one particle property at a time)
36 namespace PtlSlice{
37 #ifdef ESC_PTL
38 enum{Ph=0,Ct,Gid,Flag};
39 #else
40 enum{Ph=0,Ct,Gid};
41 #endif
42 }
43 
44 struct PtlMvmt{
45  // Options for particle location management when looping over particles
46  enum SendOpt{
47  NoSend=0,
50  };
51 
52  enum ReturnOpt{
56  };
57 
60 
61  PtlMvmt(SendOpt send_opt, ReturnOpt return_opt) : send_opt(send_opt), return_opt(return_opt){}
62 };
63 
67 };
68 
69 // Species class
70 template<class Device>
71 class Species{
72  public:
73 
74  int idx;
75  bool is_electron;
76  bool is_adiabatic;
79  double mass;
80  double charge;
81  double charge_eu;
82  double c_m;
83  double c2_2m;
84 
85  bool is_deltaf;
86 
87  int ncycles;
89 
91  int n_ptl;
92  Cabana::AoSoA<ParticleDataTypes,HostType,VEC_LEN> particles;
93 
94  // Device particles
95  Cabana::AoSoA<ParticleDataTypes,Device,VEC_LEN> particles_d;
96 
98 
101 
102  /*** Could be its own class inside species? ***/
106 
107  // phase0 (for ion restoration)
108  Cabana::AoSoA<PhaseDataTypes,HostType,VEC_LEN> phase0;
109  Cabana::AoSoA<PhaseDataTypes,Device,VEC_LEN> phase0_d;
110 
111  // For electron restoration
112  Cabana::AoSoA<ParticleDataTypes,HostType,VEC_LEN> backup_particles;
113  int n_backup_particles; // Number of particles stored in backup_particles (can't be deduced from its size due to buffer)
114  /****/
115 
117 
119 
120  Eq::Profile<Device> eq_temp; // Equilibrium temperature
121  Eq::Profile<Device> eq_den; // Equilibrium density
122  Eq::Profile<Device> eq_flow; // Equilibrium flow
123 
124  Species(int idx_in, int nonadiabatic_idx_in, bool is_electron_in, bool is_adiabatic_in, KinType kintype_in, double mass_in, double charge_in, double charge_eu_in, bool is_deltaf_in,
125  int ncycles_in);
126 
127  template<class GridDevice>
128  Species(NLReader::NamelistReader& nlr, const Grid<GridDevice> &grid, const MagneticField<DeviceType> &magnetic_field, const DomainDecomposition<DeviceType>& pol_decomp, int idx_in, int nonadiabatic_idx_in);
129 
130  // Electron or ion default constructor
131  Species(SpeciesType sp_type, int n_ptl)
132  : idx(sp_type==ELECTRON ? 0 : 1),
133  is_electron(sp_type==ELECTRON),
134  mass(is_electron ? 3.344e-30 : PROTON_MASS),
136  charge_eu(is_electron ? -1.0 : 1.0),
137  is_deltaf(true),
138  is_adiabatic(false),
139  nonadiabatic_idx(idx), // Since is_adiabatic is false above
141  ncycles(is_electron ? 70 : 1),
142  c_m(charge/mass),
143  c2_2m(0.5*charge*charge/mass),
145  n_ptl(n_ptl),
146  backup_particles("backup_particles", 0),
147  particles("particles", add_vec_buffer(n_ptl)),
150  eq_temp(1.0e3,-0.1),
151  eq_den(1.0e19,-0.1),
152  owns_particles_d(false),
156 
157  // Special constructor for tests that involve tracking particles in memory
158  // The idea is to use these particles to test functions that reorder particles,
159  // e.g. sort, shift, and cleaning
160  Species(int n_ptl_in)
161  : n_ptl(n_ptl_in),
162  is_electron(true),
163  is_adiabatic(false),
164  particles("particles", add_vec_buffer(n_ptl_in)),
166  owns_particles_d(false),
171  {
172 
173  // Slice particle properties
174  auto ph = Cabana::slice<PtlSlice::Ph>(particles);
175  auto ct = Cabana::slice<PtlSlice::Ct>(particles);
176  auto gid = Cabana::slice<PtlSlice::Gid>(particles);
177 #ifdef ESC_PTL
178  auto flag = Cabana::slice<PtlSlice::Flag>(particles);
179 #endif
180 
181  // Offset gid if using MPI
182 #ifdef USE_MPI
183  long long int gid_offset = n_ptl*SML_COMM_RANK;
184 #else
185  long long int gid_offset = 0;
186 #endif
187 
188  // Assign trackable values
189  for (int i=0;i<n_ptl;i++){
190  // Set GID in order
191  gid(i) = gid_offset + i+1; // 1-indexed
192 
193  // Value of properties is gid + 0.1*(property index)
194  // First particle: (1.0, 1.1, ... 1.8)
195  // Second particle: (2.0, 2.1, ... 2.8)
196  for (int j=0;j<6;j++) ph(i, j) = gid(i) + (j)*0.1;
197  for (int j=0;j<3;j++) ct(i, j) = gid(i) + (j+6)*0.1;
198  }
199 
200  // Buffer particles: same but with gid = -1
201  if(n_ptl>0){
202  for (int i=n_ptl;i<add_vec_buffer(n_ptl);i++){
203  gid(i) = -1;
204  for (int j=0;j<6;j++) ph(i, j) = gid(i) + (j)*0.1;
205  for (int j=0;j<3;j++) ct(i, j) = gid(i) + (j+6)*0.1;
206  }
207  }
208  }
209 
210  void resize_particles(int new_n_ptl){
211  n_ptl = new_n_ptl;
212 
213 #ifndef USE_GPU
214  // If CPU-only, particles_d points to the same location as particles. If particles is resized, then Cabana will not deallocate the first allocation
215  // since it is still used by particles_d. So, reset particles_d before resize, and point it back to particles only afterwards
216  particles_d = Cabana::AoSoA<ParticleDataTypes,Device,VEC_LEN>();
217 #endif
218 
219  particles.reserve(minimum_ptl_reservation); // Can only raise reservation (no-op if AoSoA is already larger)
220  particles.resize(add_vec_buffer(n_ptl));
221 
222 #ifndef USE_GPU
223  // Point particles_d back to particles after resize
225 #endif
226 
228  }
229 
231 #ifdef USE_GPU
232  particles.reserve(minimum_ptl_reservation); // Can only raise reservation (no-op if AoSoA is already larger)
233  particles.resize(particles_d.size());
234 #else
235  // Point particles back to particles_d after resize
237 #endif
238 
240  }
241 
242  /* If using CPU-only, then "device" particles are a shallow copy of host particles so that
243  * there is no unnecessary duplication. When "device" particles are resized, Cabana will keep the
244  * original allocation if there is a second reference (i.e. host particles).
245  * To resolve this, we free the host particles here so that there is no second reference
246  * */
248  particles = Cabana::AoSoA<ParticleDataTypes,HostType,VEC_LEN>();
249  }
250 
251  /* Resizes device particles if on GPU, or just creates a shallow copy if CPU only
252  * */
254  if(!owns_particles_d) exit_XGC("\nSpecies tried to resize device particles, but doesn't own the device array.");
255 
256 #ifdef USE_GPU
257  particles_d.reserve(minimum_ptl_reservation); // Can only raise reservation (no-op if AoSoA is already larger)
258  // Resize device particles to match host particles
259  particles_d.resize(particles.size());
260 #else
261  // If kernels are on CPU, do shallow copy
263 #endif
264  }
265 
266  /* Resizes device particles
267  * */
268  void resize_device_particles(int new_n_ptl){
269  if(!owns_particles_d) exit_XGC("\nSpecies tried to resize device particles, but doesn't own the device array.");
270 
271  n_ptl = new_n_ptl;
272 
273  particles_d.reserve(minimum_ptl_reservation); // Can only raise reservation (no-op if AoSoA is already larger)
274 
275  // Resize device particles to match host particles
277  }
278 
279  /* Copies particles to device - deep copy if using GPU, otherwise shallow copy
280  * Also takes the opportunity to set the buffer particles to realistic values
281  * */
283  if(!owns_particles_d) exit_XGC("\nSpecies tried to copy particles to device, but doesn't own the device array.");
284 
285 #ifdef PRINT_COPA_ARTIFACTS
286  if(is_rank_zero()) printf("\nCopying Cabana AoSoA particle data to GPU with Cabana::deep_copy\n");
287 #endif
288 
289 #ifdef USE_GPU
290  // Copy to device
291  Cabana::deep_copy(particles_d, particles);
292 #else
293  // No operation required if CPU-only
294 #endif
295 
296  // Copy last particle to fill remainder of trailing vector in AoSoA
298  }
299 
300  /* Copies particles from device - deep copy if using GPU, otherwise no copy is necessary
301  * */
303  if(!owns_particles_d) exit_XGC("\nSpecies tried to copy particles from device, but doesn't own the device array.");
304 
305 #ifdef PRINT_COPA_ARTIFACTS
306  if(is_rank_zero()) printf("\nCopying Cabana AoSoA particle data from GPU with Cabana::deep_copy\n");
307 #endif
308 
309 #ifdef USE_GPU
310  // Copy particles to host
311  Cabana::deep_copy(particles, particles_d);
312 #else
313  // No operation required if CPU-only
314 #endif
315  }
316 
317  /* Copies particles to device if they are resident on the device
318  * */
321  }
322 
323  /* Copies particles from device if they are resident on the device
324  * */
327  }
328 
329  /* Copies particles to device if they are NOT resident on the device
330  * */
333  }
334 
335  /* Copies particles from device if they are NOT resident on the device
336  * */
339  }
340 
349  if (n_ptl>0){
350  int last_ptl_index = n_ptl - 1;
351  auto ph = Cabana::slice<PtlSlice::Ph>(particles_d);
352  auto ct = Cabana::slice<PtlSlice::Ct>(particles_d);
353  auto gid = Cabana::slice<PtlSlice::Gid>(particles_d);
354 #ifdef ESC_PTL
355  auto flag = Cabana::slice<PtlSlice::Flag>(particles_d);
356 #endif
357 
358  Kokkos::parallel_for("set_buffer_particles_d", Kokkos::RangePolicy<ExSpace>( n_ptl, add_vec_buffer(n_ptl) ), KOKKOS_LAMBDA( const int i ){
359  // Buffer particles: same as last particle, gid = -1
360  for (int j=0;j<6;j++) ph(i, j) = ph(last_ptl_index, j);
361  for (int j=0;j<3;j++) ct(i, j) = ct(last_ptl_index, j);
362  gid(i) = -1;
363 #ifdef ESC_PTL
364  flag(i) = flag(last_ptl_index);
365 #endif
366  });
367  }
368  }
369 
378  if (n_ptl>0){
379  int last_ptl_index = n_ptl - 1;
380  auto ph = Cabana::slice<PtlSlice::Ph>(phase0_d);
381 
382  Kokkos::parallel_for("set_buffer_phase0", Kokkos::RangePolicy<ExSpace>( n_ptl, add_vec_buffer(n_ptl) ), KOKKOS_LAMBDA( const int i ){
383  // copy final real particle
384  for (int j=0;j<6;j++) ph(i, j) = ph(last_ptl_index, j);
385  });
386  }
387  }
388 
389  // Options for custom launch bounds since kokkos defaults are suboptimal for electron push kernel
390  enum class LaunchBounds{
391  Default,
392  Custom
393  };
394 
400  template<typename F>
401  inline void for_all_particles(const std::string label, F lambda_func) const {
402 #ifdef PRINT_COPA_ARTIFACTS
403  if(is_rank_zero()) printf("\nLaunching GPU kernel '%s' on Cabana AoSoA of particles\n", label.c_str());
404 #endif
405  Kokkos::RangePolicy<ExSpace> particle_range_policy( 0, p_range<DeviceType>(n_ptl) );
406  Kokkos::parallel_for(label, Opt::require(particle_range_policy, Async), lambda_func);
407  }
408 
416  template<typename F>
417  inline void for_all_particles(const std::string label, F lambda_func,
418  const PtlMvmt mvmt, LaunchBounds launch_bounds=LaunchBounds::Default) {
419  if(!owns_particles_d) exit_XGC("\nSpecies tried to loop over particles on device, but doesn't own the device array.");
420 
421 #ifdef PRINT_COPA_ARTIFACTS
422  if(is_rank_zero()) printf("\nLaunching GPU kernel '%s' on Cabana AoSoA of particles\n", label.c_str());
423 #endif
424 
425  bool use_streaming = stream_particles;
426 #ifndef USE_STREAMS
427  use_streaming = false; // Just to be safe, turn streams off here
428 #endif
429 
430  bool send_ptl = ( mvmt.send_opt==PtlMvmt::Send ||
432  bool return_ptl = ( mvmt.return_opt==PtlMvmt::Return ||
434 
435  // Don't need to stream if particles are already present and don't need to be returned
436  if((!send_ptl) && (!return_ptl)) use_streaming = false;
437 
438  if(use_streaming){
439 #ifdef USE_STREAMS
440  Streamed::Option stream_option = Streamed::Normal; // Send to device and back
441  if(!send_ptl) stream_option = Streamed::NoSend;
442  if(!return_ptl) stream_option = Streamed::NoReturn;
443 
444  // Execute streaming parallel_for
445  Streamed::parallel_for(label, n_ptl, lambda_func, stream_option, particles, particles_d);
446 #endif
447  }else{
448  if(send_ptl) TIMER("copy_ptl_to_device",copy_particles_to_device() );
449 
450  if (launch_bounds==LaunchBounds::Custom) {
451 #ifdef USE_EPUSH_LAUNCH_BOUNDS
452 # if !defined(PUSH_MAX_THREADS_PER_BLOCK) || !defined(PUSH_MIN_WARPS_PER_EU)
453 # error "USE_EPUSH_LAUNCH_BOUNDS requires PUSH_MAX_THREADS_PER_BLOCK and PUSH_MIN_WARPS_PER_EU to be defined"
454 # endif
455  Kokkos::RangePolicy<ExSpace, Kokkos::LaunchBounds<PUSH_MAX_THREADS_PER_BLOCK, PUSH_MIN_WARPS_PER_EU>>
456  particle_range_policy( 0, p_range<DeviceType>(n_ptl) );
457  Kokkos::parallel_for(label, Opt::require(particle_range_policy, Async), lambda_func);
458 #else
459  exit_XGC("\nERROR: LaunchBounds::Custom specified, but USE_EPUSH_LAUNCH_BOUNDS is not defined\n");
460 #endif
461  } else {
462  Kokkos::RangePolicy<ExSpace>
463  particle_range_policy( 0, p_range<DeviceType>(n_ptl) );
464  Kokkos::parallel_for(label, Opt::require(particle_range_policy, Async), lambda_func);
465  }
466 
467  if(return_ptl) TIMER("copy_ptl_from_device", copy_particles_from_device() );
468  }
469  }
470 
471  KOKKOS_INLINE_FUNCTION VecParticles* ptl() const{
472  return (VecParticles*)(&particles_d.access(0));
473  }
474 
475  KOKKOS_INLINE_FUNCTION VecPhase* ph0() const{
476  return (VecPhase*)(&phase0_d.access(0));
477  }
478 
480  for_all_particles("copy_to_phase0", KOKKOS_LAMBDA( const int idx ){
481  AoSoAIndices<DeviceType> inds(idx);
482  VecParticles* ptl_loc = species.ptl();
483  VecPhase* ph0_loc = species.ph0();
484  for (int i_simd = 0; i_simd<SIMD_SIZE; i_simd++){
485  int p_vec = inds.a + i_simd;
486  ph0_loc[inds.s].r[p_vec] = ptl_loc[inds.s].ph.r[p_vec];
487  ph0_loc[inds.s].z[p_vec] = ptl_loc[inds.s].ph.z[p_vec];
488  ph0_loc[inds.s].phi[p_vec] = ptl_loc[inds.s].ph.phi[p_vec];
489  ph0_loc[inds.s].rho[p_vec] = ptl_loc[inds.s].ph.rho[p_vec];
490  ph0_loc[inds.s].w1[p_vec] = ptl_loc[inds.s].ph.w1[p_vec];
491  ph0_loc[inds.s].w2[p_vec] = ptl_loc[inds.s].ph.w2[p_vec];
492  }
493  });
494  }
495 
498  // If particles are resident on the device, then use the host particle allocation as the backup
499  // Otherwise, resize the backup_particle array
502  }else{
503  backup_particles.resize(particles_d.size());
504  }
505 
506  // Copy particle data to backup
507  Cabana::deep_copy(backup_particles, particles_d);
509  } else {
510  // For ions, save phase to phase0
511  phase0_d.resize(particles_d.size());
512  copy_to_phase0(*this); // Separate function to avoid implicit copy into lambda
513  }
515  }
516 
519  // If particles are resident on device, don't need to resize host particles since they are already used as the backup.
520  // If not, resize host particle array to fit backup particles
524  }else{
525  // If particles are not resident on device, resize host particle array to fit backup particles
527  }
528 
529  // Resize device particle array to fit backed up particles
531 
532  // Restore particle data from backup
533  Cabana::deep_copy(particles_d, backup_particles);
534 
536  // If the backup particles are simply pointing to the host particles, then
537  // point them to their own 0-sized allocation when finished using them
538  // so that they don't make a copy when host particles get resized
539  backup_particles = Cabana::AoSoA<ParticleDataTypes,HostType,VEC_LEN>("backup_particles", 0);
540  }
541  } else {
542  // When ipc==2, copy phase0 to device
543  // First, resize device phase0 and reset pointer
544  phase0_d.resize(phase0.size());
545 
546  // Next copy data to phase0 on device
547  Cabana::deep_copy(phase0_d, phase0);
548 
549  // Fill buffer with realistic ptl data
551  }
552  particles_are_backed_up = false;
553  }
554 
555  KOKKOS_INLINE_FUNCTION void restore_phase_from_phase0(const AoSoAIndices<Device>& inds, SimdParticles& part_one) const {
556  VecPhase* ph0_loc = ph0();
557  for (int i_simd = 0; i_simd<SIMD_SIZE; i_simd++){
558  int p_vec = inds.a + i_simd;
559  part_one.ph.r[i_simd] = ph0_loc[inds.s].r[p_vec];
560  part_one.ph.z[i_simd] = ph0_loc[inds.s].z[p_vec];
561  part_one.ph.phi[i_simd] = ph0_loc[inds.s].phi[p_vec];
562  part_one.ph.rho[i_simd] = ph0_loc[inds.s].rho[p_vec];
563  part_one.ph.w1[i_simd] = ph0_loc[inds.s].w1[p_vec];
564  part_one.ph.w2[i_simd] = ph0_loc[inds.s].w2[p_vec];
565  }
566  }
567 
568  long long int get_total_n_ptl(){
569 #ifdef USE_MPI
570  long long int tmp_n_ptl = n_ptl;
571  long long int out_n_ptl = 0;
572  MPI_Allreduce(&tmp_n_ptl, &out_n_ptl, 1, MPI_LONG_LONG_INT, MPI_SUM, SML_COMM_WORLD);
573  return out_n_ptl;
574 #else
575  return (long long int)(n_ptl);
576 #endif
577  }
578 
580 #ifdef USE_MPI
581  int tmp_n_ptl = n_ptl;
582  int out_n_ptl = 0;
583  MPI_Allreduce(&tmp_n_ptl, &out_n_ptl, 1, MPI_INT, MPI_MAX, SML_COMM_WORLD);
584  return out_n_ptl;
585 #else
586  return n_ptl;
587 #endif
588  }
589 };
590 
591 #include "species.tpp"
592 #endif
Cabana::AoSoA< PhaseDataTypes, HostType, VEC_LEN > phase0
Definition: species.hpp:108
bool stream_particles
Whether to stream particles between host and device if possible.
Definition: species.hpp:100
Definition: globals.hpp:77
KOKKOS_INLINE_FUNCTION VecPhase * ph0() const
Definition: species.hpp:475
KOKKOS_INLINE_FUNCTION int divide_and_round_up(int a, int b)
Definition: globals.hpp:126
bool owns_particles_d
Whether the species owns the device particle allocation right now.
Definition: species.hpp:97
subroutine adjust_n_ptl_for_core_ptl(n_ptl)
Definition: load.F90:445
void set_spall_num_and_ptr(int idx, int n_ptl, int n_vecs, VecParticles *ptl)
bool is_rank_zero()
Definition: globals.hpp:26
Distribution< Device > f0
Species distribution in velocity space on local mesh nodes.
Definition: species.hpp:116
MPI_Comm SML_COMM_WORLD
Definition: my_mpi.cpp:4
Cabana::AoSoA< ParticleDataTypes, HostType, VEC_LEN > backup_particles
Copy of particles to be restored for RK2.
Definition: species.hpp:112
bool is_electron
Whether this species is the electrons.
Definition: species.hpp:75
void for_all_particles(const std::string label, F lambda_func, const PtlMvmt mvmt, LaunchBounds launch_bounds=LaunchBounds::Default)
Definition: species.hpp:417
void save_backup_particles()
Definition: species.hpp:496
double c2_2m
c2/2m
Definition: species.hpp:83
double rho[VEC_LEN]
Definition: particles.hpp:72
void copy_to_phase0(Species< Device > &species)
Definition: species.hpp:479
Definition: species.hpp:54
Simd< double > w1
Definition: particles.hpp:22
double c_m
c/m
Definition: species.hpp:82
Definition: species.hpp:53
bool default_streaming_option()
Definition: species.hpp:24
Eq::Profile< Device > eq_den
Definition: species.hpp:121
Definition: globals.hpp:82
KOKKOS_INLINE_FUNCTION VecParticles * ptl() const
Definition: species.hpp:471
Definition: NamelistReader.hpp:163
KinType kintype
Whether the species is gyrokinetic or drift kinetic.
Definition: species.hpp:78
Definition: magnetic_field.hpp:12
int add_vec_buffer(int n_ptl)
Definition: particles.hpp:164
int idx
Index in all_species.
Definition: species.hpp:74
Definition: particles.hpp:68
int a
The index in the inner array of the AoSoA.
Definition: particles.hpp:120
Definition: grid.hpp:10
Definition: particles.hpp:85
bool particles_are_backed_up
Whether particles are currently backed up.
Definition: species.hpp:105
int nonadiabatic_idx
Index of species skipping adiabatic species (for compatibility with fortran arrays) ...
Definition: species.hpp:77
bool default_residence_option()
Definition: species.hpp:31
int n_ptl
Number of particles.
Definition: species.hpp:91
Definition: streamed_parallel_for.hpp:16
Definition: streamed_parallel_for.hpp:14
void set_buffer_phase0_d()
Definition: species.hpp:377
long long int get_total_n_ptl()
Definition: species.hpp:568
void set_buffer_particles_d()
Definition: species.hpp:348
Simd< double > rho
Definition: particles.hpp:21
Definition: species.hpp:55
int p_range< DeviceType >(int num_particle)
Definition: particles.hpp:157
double charge_eu
Particle charge in eu.
Definition: species.hpp:81
Definition: species.hpp:47
void resize_particles(int new_n_ptl)
Definition: species.hpp:210
double mass
Particle mass.
Definition: species.hpp:79
Species(int idx_in, int nonadiabatic_idx_in, bool is_electron_in, bool is_adiabatic_in, KinType kintype_in, double mass_in, double charge_in, double charge_eu_in, bool is_deltaf_in, int ncycles_in)
Definition: species.tpp:22
void for_all_particles(const std::string label, F lambda_func) const
Definition: species.hpp:401
Cabana::AoSoA< ParticleDataTypes, Device, VEC_LEN > particles_d
Particles on device.
Definition: species.hpp:95
Definition: species.hpp:44
double w2[VEC_LEN]
Definition: particles.hpp:74
#define TIMER(N, F)
Definition: timer_macro.hpp:29
RKRestorationMethod
Definition: species.hpp:64
idx
Definition: diag_f0_df_port1.hpp:32
void copy_particles_to_device_if_not_resident()
Definition: species.hpp:331
RKRestorationMethod RK_restoration_method
Currently, electrons must use first method and ions must use second.
Definition: species.hpp:103
Simd< double > r
Definition: particles.hpp:18
void resize_host_particles_to_match_device()
Definition: species.hpp:230
Definition: species.hpp:65
ReturnOpt return_opt
Definition: species.hpp:59
ReturnOpt
Definition: species.hpp:52
Option
Definition: streamed_parallel_for.hpp:13
void restore_particles_from_backup()
Definition: species.hpp:517
Definition: globals.hpp:83
SendOpt send_opt
Definition: species.hpp:58
double charge
Particle charge.
Definition: species.hpp:80
SimdPhase ph
Definition: particles.hpp:59
void copy_particles_from_device()
Definition: species.hpp:302
void copy_particles_from_device_if_not_resident()
Definition: species.hpp:337
void unassign_host_particles()
Definition: species.hpp:247
int ncycles_between_sorts
Number of subcycles between sorts.
Definition: species.hpp:88
Definition: particles.hpp:58
Cabana::AoSoA< PhaseDataTypes, Device, VEC_LEN > phase0_d
Definition: species.hpp:109
int SML_COMM_RANK
Definition: my_mpi.cpp:5
KinType
Definition: globals.hpp:81
Species(SpeciesType sp_type, int n_ptl)
Definition: species.hpp:131
bool is_deltaf
Whether this species is deltaf.
Definition: species.hpp:85
VecPhase ph
Definition: particles.hpp:86
Definition: species.hpp:40
Definition: species.hpp:66
constexpr double PROTON_MASS
Definition: globals.hpp:121
void set_min_max_num(int isp, int n_ptl)
int minimum_ptl_reservation
The minimum reservation size for particles.
Definition: species.hpp:90
int s
The index in the outer array of the AoSoA.
Definition: particles.hpp:119
Simd< double > z
Definition: particles.hpp:19
void copy_particles_to_device_if_resident()
Definition: species.hpp:319
Definition: species.hpp:48
void resize_device_particles(int new_n_ptl)
Definition: species.hpp:268
Definition: species.hpp:49
void exit_XGC(std::string msg)
Definition: globals.hpp:36
void copy_particles_from_device_if_resident()
Definition: species.hpp:325
bool is_adiabatic
Whether this species is adiabatic.
Definition: species.hpp:76
Simd< double > phi
Definition: particles.hpp:20
Definition: magnetic_field.F90:1
static constexpr const Kokkos::Experimental::WorkItemProperty::HintLightWeight_t Async
Definition: space_settings.hpp:74
int n_backup_particles
Definition: species.hpp:113
Eq::Profile< Device > eq_flow
Definition: species.hpp:122
Definition: streamed_parallel_for.hpp:15
SendOpt
Definition: species.hpp:46
void copy_particles_to_device()
Definition: species.hpp:282
KOKKOS_INLINE_FUNCTION void restore_phase_from_phase0(const AoSoAIndices< Device > &inds, SimdParticles &part_one) const
Definition: species.hpp:555
Species(int n_ptl_in)
Definition: species.hpp:160
double phi[VEC_LEN]
Definition: particles.hpp:71
Simd< double > w2
Definition: particles.hpp:23
double r[VEC_LEN]
Definition: particles.hpp:69
Definition: species.hpp:71
void parallel_for(const std::string name, int n_ptl, Function func, Option option, HostAoSoA aosoa_h, DeviceAoSoA aosoa_d)
Definition: streamed_parallel_for.hpp:252
Definition: species.hpp:40
Eq::Profile< Device > eq_temp
Definition: species.hpp:120
bool particles_resident_on_device
Whether the particles can reside on device.
Definition: species.hpp:99
PtlMvmt(SendOpt send_opt, ReturnOpt return_opt)
Definition: species.hpp:61
int ncycles
Number of subcycles.
Definition: species.hpp:87
Definition: profile.hpp:65
int collision_grid_index
Which collision grid to use.
Definition: species.hpp:118
Definition: particles.hpp:118
SpeciesType
Definition: globals.hpp:76
double z[VEC_LEN]
Definition: particles.hpp:70
constexpr double UNIT_CHARGE
Charge of an electron (C)
Definition: globals.hpp:118
Definition: distribution.hpp:10
void resize_device_particles()
Definition: species.hpp:253
int get_max_n_ptl()
Definition: species.hpp:579
LaunchBounds
Definition: species.hpp:390
Cabana::AoSoA< ParticleDataTypes, HostType, VEC_LEN > particles
Particles.
Definition: species.hpp:92
double w1[VEC_LEN]
Definition: particles.hpp:73
Definition: species.hpp:40