XGC1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends 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 #include "diffusion.hpp"
8 #include "radiation.hpp"
9 #include "pellet.hpp"
10 #include "thermal_bath.hpp"
11 
12 enum class CollisionMode{
13  None = 0,
14  MonteCarloNonConserving, // monte-carlo non-conserving
15  NonLinearFokkerPlanckLandau, // fully non-linear multi-species Fokker-Planck-Landau [Hager JCP 315, 644-660 (2016)]
16  PETScLandauAMR // PETSc multi-species Landau with adaptive mesh refinement (ongoing implementation)
17 };
18 
19 struct Sources{
20  private:
21 
22  StepTrigger step_trigger; // Triggers when f_source is called
23 
24  public:
25 
27 
28  // List of sources (other collisions-related members will be folded into a child class of Source)
37 
41 
42  Sources(NLReader::NamelistReader& nlr, const DomainDecomposition<DeviceType>& pol_decomp, const MagneticField<DeviceType>& magnetic_field, const Grid<DeviceType>& grid, const VelocityGrid& vgrid, Plasma& plasma, bool overwrite_existing_files) {
43  nlr.use_namelist("sml_param");
44  f_source_period = nlr.get<int>("sml_f_source_period", 1);
45 
47 
48  // Set up for each source:
49  const int CALL_AT_F_SOURCE = 1; // Set these triggers so that these sources are called every time f_source is called
50  bool sml_radiation = nlr.get<bool>("sml_radiation",false);
51  // data base (must be configured with the rad_param namelist)
52  bool sml_neutral = nlr.get<bool>("sml_neutral",false);
53  bool sml_source = nlr.get<bool>("sml_source",false);
54  bool sml_current_drive = nlr.get<bool>("sml_current_drive_on",false);
55  bool sml_thermal_bath_on = nlr.get<bool>("sml_thermal_bath_on",false);
56  if (sml_radiation) sml_neutral = true; // radiation overrides neutrals (should throw error instead)
57 
58  if(sml_radiation){
59  radiation = Radiation(nlr, grid, magnetic_field);
60  }
61 
62  if(sml_neutral){
63  nlr.use_namelist("neu_param");
64  int start_time = nlr.get<int>("neu_start_time",2);
65  neutrals.step_trigger = StepTrigger(CALL_AT_F_SOURCE, start_time);
66  }
67 
68  if(sml_source){
69  nlr.use_namelist("src_param");
70  int start_time = nlr.get<int>("src_start_time",1);
71  heat_and_torque.step_trigger = StepTrigger(CALL_AT_F_SOURCE, start_time);
72  }
73 
74  if(sml_current_drive){
75  nlr.use_namelist("src_param");
76  int start_time = nlr.get<int>("src_current_drive_start_time",1);
77  current_drive.step_trigger = StepTrigger(CALL_AT_F_SOURCE, start_time);
78  }
79 
80  if(sml_thermal_bath_on){
82  }
83 
84  nlr.use_namelist("diff_param");
85  bool diff_on = nlr.get<bool>("diff_on", false);
86  if(diff_on){
87  // If diffusion is on, call diffusion constructor
88  diffusion = Diffusion(nlr, magnetic_field);
89 
90  int start_time = nlr.get<int>("diff_start_time", 1);
91  diffusion.step_trigger = StepTrigger(CALL_AT_F_SOURCE, start_time);
92  }
93 
94  nlr.use_namelist("src_param");
95  bool pellet_on = nlr.get<bool>("src_pellet_on", false);
96  if(pellet_on){
97  // If pellet is on, call pellet constructor
98  pellet = Pellet(nlr, grid, plasma, magnetic_field);
99  }
100 
101  // Collisions
102  nlr.use_namelist("col_param");
103  int col_mode_in = nlr.get<int>("col_mode", 0);
104  if(col_mode_in==0) col_mode = CollisionMode::None;
105  else if(col_mode_in==1) col_mode = CollisionMode::MonteCarloNonConserving;
106  else if(col_mode_in==4) col_mode = CollisionMode::NonLinearFokkerPlanckLandau;
107  else if(col_mode_in==5) col_mode = CollisionMode::PETScLandauAMR;
108  else{
109  if(col_mode_in==2){
110  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");
111  }
112  exit_XGC("\nInvalid col_mode specified. Acceptable options are 0, 1, 4, or 5.\n");
113  }
114 
116  col_grid = CollisionGrid<DeviceType>(nlr, magnetic_field, grid, overwrite_existing_files);
118  col_mc = MonteCarloCollider<DeviceType>(nlr, pol_decomp, magnetic_field, grid, plasma.n_nonadiabatic_species);
119  }
120 
122  collisions.step_trigger = StepTrigger(CALL_AT_F_SOURCE, 0); // start time is zero
123  }
124  }
125 
126  // Close streams for diagnostics associated with a source. This could be part of a destructor.
129  col_grid.io_stream->Close();
130  }
131  }
132 
133  /* Returns whether the sources should be executed on this step
134  * */
135  bool is_triggered(int step){
136  return step_trigger.is_triggered(step);
137  }
138 };
139 
140 #endif
Source neutrals
Definition: sources.hpp:29
Sources(NLReader::NamelistReader &nlr, const DomainDecomposition< DeviceType > &pol_decomp, const MagneticField< DeviceType > &magnetic_field, const Grid< DeviceType > &grid, const VelocityGrid &vgrid, Plasma &plasma, bool overwrite_existing_files)
Definition: sources.hpp:42
Definition: sources.hpp:19
bool is_rank_zero()
Definition: globals.hpp:27
T get(const string &param, const T default_val, int val_ind=0)
Definition: NamelistReader.hpp:373
Source collisions
Definition: sources.hpp:33
int f_source_period
Frequency of f_sources.
Definition: sources.hpp:26
Definition: velocity_grid.hpp:8
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:1224
Definition: thermal_bath.hpp:8
Definition: NamelistReader.hpp:193
Definition: magnetic_field.hpp:12
Definition: pellet.hpp:12
MonteCarloCollider< DeviceType > col_mc
Monte carlo collisions, only used by MonteCarloNonConserving.
Definition: sources.hpp:39
Source current_drive
Definition: sources.hpp:34
Definition: source.hpp:7
Definition: diffusion.hpp:15
CollisionGrid< DeviceType > col_grid
Collision grid, only used by NonLinearFokkerPlanckLandau.
Definition: sources.hpp:40
void use_namelist(const string &namelist)
Definition: NamelistReader.hpp:355
Source heat_and_torque
Definition: sources.hpp:32
bool is_triggered(int step) const
Definition: step_trigger.hpp:25
CollisionMode
Definition: sources.hpp:12
int n_nonadiabatic_species
Number of nonadiabatic species.
Definition: plasma.hpp:95
Diffusion diffusion
Definition: sources.hpp:30
StepTrigger step_trigger
Definition: source.hpp:8
CollisionMode col_mode
Collision operator model.
Definition: sources.hpp:38
void exit_XGC(std::string msg)
Definition: globals.hpp:37
Pellet pellet
Definition: sources.hpp:35
Definition: magnetic_field.F90:1
Definition: plasma.hpp:14
Radiation radiation
Definition: sources.hpp:31
ThermalBath thermal_bath
Definition: sources.hpp:36
std::shared_ptr< XGC_IO_Stream > io_stream
Definition: col_grid.hpp:214
void close_streams()
Definition: sources.hpp:127
StepTrigger step_trigger
Definition: sources.hpp:22
bool is_triggered(int step)
Definition: sources.hpp:135
Definition: step_trigger.hpp:4
bool diag_on
Definition: col_grid.hpp:213
Definition: radiation.hpp:13