Hipace
Hipace.H
Go to the documentation of this file.
1 /* Copyright 2020-2022
2  *
3  * This file is part of HiPACE++.
4  *
5  * Authors: AlexanderSinn, Andrew Myers, Axel Huebl, MaxThevenet
6  * Remi Lehe, Severin Diederichs, WeiqunZhang
7  * License: BSD-3-Clause-LBNL
8  */
9 #ifndef HIPACE_H_
10 #define HIPACE_H_
11 
12 #include "fields/Fields.H"
17 #include "utils/AdaptiveTimeStep.H"
18 #include "utils/GridCurrent.H"
19 #include "laser/MultiLaser.H"
20 #include "utils/Constants.H"
21 #include "utils/Parser.H"
22 #include "utils/MultiBuffer.H"
23 #include "diagnostics/Diagnostic.H"
25 
26 #include <AMReX_AmrCore.H>
27 #ifdef AMREX_USE_LINEAR_SOLVERS
28 # include <AMReX_MLALaplacian.H>
29 # include <AMReX_MLMG.H>
30 #endif
31 
32 #ifdef HIPACE_USE_OPENPMD
33 # include <openPMD/openPMD.hpp>
34 #endif
35 
36 #include <memory>
37 
38 namespace hpmg { class MultiGrid; }
39 
42 {
45  Hipace_early_init (Hipace* instance);
46 
50 
53  inline static int m_depos_order_xy = 2;
56  inline static int m_depos_order_z = 0;
59  inline static int m_depos_derivative_type = 2;
60  /* if the loop over depos_order is included in the loop over particles
61  */
62  inline static bool m_outer_depos_loop = false;
63 
64  /* Number of mesh refinement levels */
65  int m_N_level = 1;
66 };
67 
69 class Hipace final : public Hipace_early_init
70 {
71 public:
74  Hipace ();
75 
77  static Hipace& GetInstance ();
78 
80  void InitData ();
81 
83  void Evolve ();
84 
86  void MakeGeometry ();
87 
93  void WriteDiagnostics (int output_step, const OpenPMDWriterCallType call_type);
94 
97 
103  void SolveOneSlice (int islice, int step);
104 
110  void InitializeSxSyWithBeam (const int lev);
111 
118  void ExplicitMGSolveBxBy (const int lev, const int which_slice);
119 
124  void ResetAllQuantities ();
125 
127  void ResetLaser ();
128 
132  void doCoulombCollision ();
133 
137  static bool HeadRank ()
138  {
140  }
141 
146  static std::string Version ();
147 
161  inline static Hipace* m_instance = nullptr;
163  inline static bool m_normalized_units = false;
170 
173  inline static int m_numprocs = 0;
175  inline static int m_max_step = 0;
177  inline static amrex::Real m_max_time = std::numeric_limits<amrex::Real>::infinity();
180  inline static amrex::Real m_physical_time = 0.0;
182  inline static amrex::Real m_initial_time = 0.0;
183 
184  bool m_has_last_step = false;
186  inline static int m_verbose = 0;
189  inline static amrex::Real m_predcorr_B_error_tolerance = 4e-2;
192  inline static int m_predcorr_max_iterations = 30;
195  amrex::Real m_predcorr_avg_iterations = 0.;
198  amrex::Real m_predcorr_avg_B_error = 0.;
201  inline static amrex::Real m_predcorr_B_mixing_factor = 0.05;
203  inline static bool m_do_beam_jx_jy_deposition = true;
205  inline static bool m_do_beam_jz_minus_rho = false;
207  inline static bool m_deposit_rho = false;
209  inline static bool m_deposit_rho_individual = false;
211  inline static bool m_interpolate_neutralizing_background = false;
213 #ifdef AMREX_USE_GPU
214  inline static bool m_do_tiling = false;
215 #else
216  inline static bool m_do_tiling = true;
217 #endif
219  inline static bool m_explicit = true;
221  inline static amrex::Real m_MG_tolerance_rel = 1.e-4;
223  inline static amrex::Real m_MG_tolerance_abs = std::numeric_limits<amrex::Real>::min();
225  inline static int m_MG_verbose = 0;
227  inline static bool m_use_amrex_mlmg = false;
229  inline static bool m_use_laser = false;
232  inline static amrex::Real m_background_density_SI = 0.;
234  amrex::Real m_dt = 0.0;
236  inline static int m_ncollisions = 0;
243 #ifdef HIPACE_USE_OPENPMD
245  OpenPMDWriter m_openpmd_writer;
246 #endif
250  bool m_salame_do_advance = true;
254  bool m_salame_overloaded = false;
256  amrex::Real m_salame_zeta_initial = 0;
262  bool m_comms_buffer_on_gpu = false;
264  int m_comms_buffer_max_leading_slices = std::numeric_limits<int>::max();
266  int m_comms_buffer_max_trailing_slices = std::numeric_limits<int>::max();
267 
268 private:
269 
270 #ifdef AMREX_USE_LINEAR_SOLVERS
275 #endif
281  std::vector<std::string> m_collision_names;
284 
285  void InitDiagnostics (const int step);
286  void FillFieldDiagnostics (const int lev, int islice);
287  void FillBeamDiagnostics (const int step);
288  void WriteDiagnostics (const int step);
289  void FlushDiagnostics ();
290 
293 
309  void PredictorCorrectorLoopToSolveBxBy (const int islice, const int current_N_level,
310  const int step);
311 };
312 
313 #endif
OpenPMDWriterCallType
Whether the beam, the field data is written, or if it is just flushing the stored data.
Definition: OpenPMDWriter.H:30
class handling the adaptive time step
Definition: AdaptiveTimeStep.H:18
This class holds data for all diagnostics.
Definition: Diagnostic.H:55
amrex::Vector< std::string > & getBeamNames()
return names of the beams to output
Definition: Diagnostic.H:66
Main class handling all field data structures and operations.
Definition: Fields.H:85
class handling a current directly written to the grid
Definition: GridCurrent.H:16
Singleton class, that intialize, runs and finalizes the simulation.
Definition: Hipace.H:70
int m_salame_last_slice
Definition: Hipace.H:252
void ExplicitMGSolveBxBy(const int lev, const int which_slice)
Knowing the sources Sx, Sy and chi, apply MG to solve for Bx, By.
Definition: Hipace.cpp:694
amrex::Vector< CoulombCollision > m_all_collisions
Definition: Hipace.H:283
void InitData()
Definition: Hipace.cpp:172
Hipace()
Definition: Hipace.cpp:64
static int m_ncollisions
Definition: Hipace.H:236
void FlushDiagnostics()
Definition: Hipace.cpp:1017
static amrex::Real m_predcorr_B_error_tolerance
Definition: Hipace.H:189
MultiLaser m_multi_laser
Definition: Hipace.H:240
Diagnostic m_diags
Definition: Hipace.H:279
void InitDiagnostics(const int step)
Definition: Hipace.cpp:956
static amrex::Real m_initial_time
Definition: Hipace.H:182
PhysConst get_phys_const()
Return a copy of member struct for physical constants.
Definition: Hipace.H:96
amrex::Vector< std::string > & getDiagBeamNames()
get diagnostics Component names of Fields to output
Definition: Hipace.H:292
MultiBeam m_multi_beam
Definition: Hipace.H:167
static int m_verbose
Definition: Hipace.H:186
static amrex::Real m_MG_tolerance_abs
Definition: Hipace.H:223
GridCurrent m_grid_current
Definition: Hipace.H:242
void WriteDiagnostics(int output_step, const OpenPMDWriterCallType call_type)
Dump simulation data to file.
int m_comms_buffer_max_leading_slices
Definition: Hipace.H:264
bool m_salame_overloaded
Definition: Hipace.H:254
void FillBeamDiagnostics(const int step)
Definition: Hipace.cpp:984
amrex::Vector< amrex::Geometry > m_slice_geom
Definition: Hipace.H:155
static int m_MG_verbose
Definition: Hipace.H:225
amrex::Vector< amrex::Geometry > m_3D_geom
Definition: Hipace.H:149
static bool m_use_amrex_mlmg
Definition: Hipace.H:227
int m_salame_n_iter
Definition: Hipace.H:248
amrex::Vector< amrex::BoxArray > m_slice_ba
Definition: Hipace.H:159
amrex::Vector< std::unique_ptr< hpmg::MultiGrid > > m_hpmg
Definition: Hipace.H:277
std::vector< std::string > m_collision_names
Definition: Hipace.H:281
amrex::Real m_dt
Definition: Hipace.H:234
static amrex::Real m_background_density_SI
Definition: Hipace.H:232
void Evolve()
Definition: Hipace.cpp:324
amrex::Real m_predcorr_avg_iterations
Definition: Hipace.H:195
static bool HeadRank()
Returns true on the head rank, otherwise false.
Definition: Hipace.H:137
static amrex::Real m_MG_tolerance_rel
Definition: Hipace.H:221
amrex::ParserExecutor< 3 > m_salame_target_func
Definition: Hipace.H:260
static std::string Version()
Definition: HipaceVersion.cpp:16
bool m_has_last_step
Definition: Hipace.H:184
amrex::Real m_predcorr_avg_B_error
Definition: Hipace.H:198
void doCoulombCollision()
does Coulomb collisions between plasmas and beams
Definition: Hipace.cpp:923
amrex::Vector< amrex::DistributionMapping > m_slice_dm
Definition: Hipace.H:157
static amrex::Real m_max_time
Definition: Hipace.H:177
void FillFieldDiagnostics(const int lev, int islice)
Definition: Hipace.cpp:974
static bool m_interpolate_neutralizing_background
Definition: Hipace.H:211
void PredictorCorrectorLoopToSolveBxBy(const int islice, const int current_N_level, const int step)
Predictor-corrector loop to calculate Bx and By.
Definition: Hipace.cpp:819
amrex::Vector< amrex::DistributionMapping > m_3D_dm
Definition: Hipace.H:151
static amrex::Real m_predcorr_B_mixing_factor
Definition: Hipace.H:201
static int m_predcorr_max_iterations
Definition: Hipace.H:192
bool m_comms_buffer_on_gpu
Definition: Hipace.H:262
amrex::Real m_salame_zeta_initial
Definition: Hipace.H:256
void SolveOneSlice(int islice, int step)
Full evolve on 1 slice.
Definition: Hipace.cpp:442
static bool m_deposit_rho
Definition: Hipace.H:207
bool m_salame_do_advance
Definition: Hipace.H:250
static int m_numprocs
Definition: Hipace.H:173
static bool m_do_beam_jz_minus_rho
Definition: Hipace.H:205
Fields m_fields
Definition: Hipace.H:165
static bool m_deposit_rho_individual
Definition: Hipace.H:209
static bool m_explicit
Definition: Hipace.H:219
void ResetAllQuantities()
Reset plasma and field slice quantities to initial value.
Definition: Hipace.cpp:614
void ResetLaser()
reset all laser slices to 0
Definition: Hipace.cpp:628
static bool m_do_tiling
Definition: Hipace.H:216
static int m_max_step
Definition: Hipace.H:175
void MakeGeometry()
Definition: Hipace.cpp:233
MultiPlasma m_multi_plasma
Definition: Hipace.H:169
static bool m_use_laser
Definition: Hipace.H:229
amrex::Vector< amrex::BoxArray > m_3D_ba
Definition: Hipace.H:153
static amrex::Real m_physical_time
Definition: Hipace.H:180
AdaptiveTimeStep m_adaptive_time_step
Definition: Hipace.H:238
static bool m_normalized_units
Definition: Hipace.H:163
void InitializeSxSyWithBeam(const int lev)
Initialize Sx and Sy with the beam contributions.
Definition: Hipace.cpp:636
MultiBuffer m_multi_buffer
Definition: Hipace.H:171
static Hipace & GetInstance()
Definition: Hipace.cpp:58
static Hipace * m_instance
Definition: Hipace.H:161
amrex::Parser m_salame_parser
Definition: Hipace.H:258
static bool m_do_beam_jx_jy_deposition
Definition: Hipace.H:203
int m_comms_buffer_max_trailing_slices
Definition: Hipace.H:266
Definition: MultiBeam.H:15
Definition: MultiBuffer.H:16
Definition: MultiLaser.H:99
Definition: MultiPlasma.H:18
Multigrid solver.
Definition: HpMultiGrid.H:34
Definition: Hipace.H:38
Helper struct to initialize m_phys_const and Parser before amrex::AmrCore.
Definition: Hipace.H:42
static int m_depos_order_z
Definition: Hipace.H:56
static int m_depos_derivative_type
Definition: Hipace.H:59
static int m_depos_order_xy
Definition: Hipace.H:53
static bool m_outer_depos_loop
Definition: Hipace.H:62
int m_N_level
Definition: Hipace.H:65
Hipace_early_init(Hipace *instance)
Definition: Hipace.cpp:31
PhysConst m_phys_const
Definition: Hipace.H:49
Struct containing physical constants, our main strategy to handle both SI and normalized units.
Definition: Constants.H:40