XGCa
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
sources.hpp
Go to the documentation of this file.
1 #ifndef SOURCES_HPP
2 #define SOURCES_HPP
3 
4 #include "step_trigger.hpp"
6 #include "col_grid.hpp"
7 
8 enum class CollisionMode{
9  None = 0,
10  MonteCarloNonConserving, // monte-carlo non-conserving
11  NonLinearFokkerPlanckLandau, // fully non-linear multi-species Fokker-Planck-Landau [Hager JCP 315, 644-660 (2016)]
12  PETScLandauAMR // PETSc multi-species Landau with adaptive mesh refinement (ongoing implementation)
13 };
14 
15 struct Sources{
16  private:
17 
18  StepTrigger step_trigger; // Triggers when f_source is called
19 
20  public:
21 
23 
24 
28 
30  nlr.use_namelist("sml_param");
31  f_source_period = nlr.get<int>("sml_f_source_period", 1);
32 
34 
35  // Collisions
36  nlr.use_namelist("col_param");
37  int col_mode_in = nlr.get<int>("col_mode", 0);
38  if(col_mode_in==0) col_mode = CollisionMode::None;
39  else if(col_mode_in==1) col_mode = CollisionMode::MonteCarloNonConserving;
40  else if(col_mode_in==4) col_mode = CollisionMode::NonLinearFokkerPlanckLandau;
41  else if(col_mode_in==5) col_mode = CollisionMode::PETScLandauAMR;
42  else{
43  if(col_mode_in==2){
44  if(is_rank_zero()) printf("\ncol_mode==2 is currently not supported. If you are reviving it, dont forget to address this potential issue: https://github.com/PrincetonUniversity/XGC-Devel/issues/266");
45  }
46  exit_XGC("\nInvalid col_mode specified. Acceptable options are 0, 1, 4, or 5.\n");
47  }
48 
50 #ifndef USE_FORTRAN_COLLISIONS
52 #endif
54  col_mc = MonteCarloCollider<DeviceType>(nlr, magnetic_field, plasma.n_nonadiabatic_species);
55  }
56  }
57 
58  /* Returns whether the sources should be executed on this step
59  * */
60  bool execute(int step){
61  return step_trigger.execute(step);
62  }
63 };
64 
65 #endif
Definition: sources.hpp:15
bool is_rank_zero()
Definition: globals.hpp:26
bool execute(int step)
Definition: sources.hpp:60
T get(const string &param, const T default_val, int val_ind=0)
Definition: NamelistReader.hpp:353
int f_source_period
Frequency of f_sources.
Definition: sources.hpp:22
Definition: velocity_grid.hpp:7
subroutine plasma(grid, itr, p, dene_out, deni_out, Te_out, Ti_out, Vparai_out)
Calculate the plasma density, temperature, and parallel velocity for a point in triangle itr using pl...
Definition: neutral_totalf.F90:1235
Definition: NamelistReader.hpp:163
Definition: magnetic_field.hpp:9
MonteCarloCollider< DeviceType > col_mc
Monte carlo collisions, only used by MonteCarloNonConserving.
Definition: sources.hpp:26
Sources(NLReader::NamelistReader &nlr, const MagneticField< DeviceType > &magnetic_field, const VelocityGrid &vgrid, const Plasma &plasma)
Definition: sources.hpp:29
CollisionGrid< DeviceType > col_grid
Collision grid, only used by NonLinearFokkerPlanckLandau.
Definition: sources.hpp:27
void use_namelist(const string &namelist)
Definition: NamelistReader.hpp:322
CollisionMode
Definition: sources.hpp:8
int n_nonadiabatic_species
Number of nonadiabatic species.
Definition: plasma.hpp:84
CollisionMode col_mode
Collision operator model.
Definition: sources.hpp:25
void exit_XGC(std::string msg)
Definition: globals.hpp:36
Definition: magnetic_field.F90:1
Definition: plasma.hpp:15
bool execute(int step)
Definition: step_trigger.hpp:18
StepTrigger step_trigger
Definition: sources.hpp:18
Definition: step_trigger.hpp:4