WarpX
HybridPICModel.H
Go to the documentation of this file.
1 /* Copyright 2023 The WarpX Community
2  *
3  * This file is part of WarpX.
4  *
5  * Authors: Roelof Groenewald (TAE Technologies)
6  *
7  * License: BSD-3-Clause-LBNL
8  */
9 
10 #ifndef WARPX_HYBRIDPICMODEL_H_
11 #define WARPX_HYBRIDPICMODEL_H_
12 
13 #include "HybridPICModel_fwd.H"
14 
16 
19 #include "Utils/WarpXConst.H"
21 
22 #include <AMReX_Array.H>
23 #include <AMReX_REAL.H>
24 
30 {
31 public:
32  HybridPICModel (int nlevs_max); // constructor
33 
35  void ReadParameters ();
36 
38  void AllocateMFs (int nlevs_max);
39  void AllocateLevelMFs (int lev, const amrex::BoxArray& ba, const amrex::DistributionMapping& dm,
40  int ncomps, const amrex::IntVect& ngJ, const amrex::IntVect& ngRho,
41  const amrex::IntVect& jx_nodal_flag, const amrex::IntVect& jy_nodal_flag,
42  const amrex::IntVect& jz_nodal_flag, const amrex::IntVect& rho_nodal_flag);
43 
45  void ClearLevel (int lev);
46 
47  void InitData ();
48 
55  void GetCurrentExternal (
56  amrex::Vector<std::array< std::unique_ptr<amrex::MultiFab>, 3>> const& edge_lengths
57  );
58  void GetCurrentExternal (
59  std::array< std::unique_ptr<amrex::MultiFab>, 3> const& edge_lengths,
60  int lev
61  );
62 
73  amrex::Vector<std::array< std::unique_ptr<amrex::MultiFab>, 3>> const& Bfield,
74  amrex::Vector<std::array< std::unique_ptr<amrex::MultiFab>, 3>> const& edge_lengths
75  );
77  std::array< std::unique_ptr<amrex::MultiFab>, 3> const& Bfield,
78  std::array< std::unique_ptr<amrex::MultiFab>, 3> const& edge_lengths,
79  int lev
80  );
81 
86  void HybridPICSolveE (
87  amrex::Vector<std::array< std::unique_ptr<amrex::MultiFab>, 3>>& Efield,
88  amrex::Vector<std::array< std::unique_ptr<amrex::MultiFab>, 3>> const& Jfield,
89  amrex::Vector<std::array< std::unique_ptr<amrex::MultiFab>, 3>> const& Bfield,
90  amrex::Vector<std::unique_ptr<amrex::MultiFab>> const& rhofield,
91  amrex::Vector<std::array< std::unique_ptr<amrex::MultiFab>, 3>> const& edge_lengths,
92  bool include_resistivity_term);
93 
94  void HybridPICSolveE (
95  std::array< std::unique_ptr<amrex::MultiFab>, 3>& Efield,
96  std::array< std::unique_ptr<amrex::MultiFab>, 3> const& Jfield,
97  std::array< std::unique_ptr<amrex::MultiFab>, 3> const& Bfield,
98  std::unique_ptr<amrex::MultiFab> const& rhofield,
99  std::array< std::unique_ptr<amrex::MultiFab>, 3> const& edge_lengths,
100  int lev, bool include_resistivity_term);
101 
102  void HybridPICSolveE (
103  std::array< std::unique_ptr<amrex::MultiFab>, 3>& Efield,
104  std::array< std::unique_ptr<amrex::MultiFab>, 3> const& Jfield,
105  std::array< std::unique_ptr<amrex::MultiFab>, 3> const& Bfield,
106  std::unique_ptr<amrex::MultiFab> const& rhofield,
107  std::array< std::unique_ptr<amrex::MultiFab>, 3> const& edge_lengths,
108  int lev, PatchType patch_type, bool include_resistivity_term);
109 
110  void BfieldEvolveRK (
111  amrex::Vector<std::array< std::unique_ptr<amrex::MultiFab>, 3>>& Bfield,
112  amrex::Vector<std::array< std::unique_ptr<amrex::MultiFab>, 3>>& Efield,
113  amrex::Vector<std::array< std::unique_ptr<amrex::MultiFab>, 3>> const& Jfield,
114  amrex::Vector<std::unique_ptr<amrex::MultiFab>> const& rhofield,
115  amrex::Vector<std::array< std::unique_ptr<amrex::MultiFab>, 3>> const& edge_lengths,
116  amrex::Real dt, DtType a_dt_type,
117  amrex::IntVect ng, std::optional<bool> nodal_sync);
118 
119  void BfieldEvolveRK (
120  amrex::Vector<std::array< std::unique_ptr<amrex::MultiFab>, 3>>& Bfield,
121  amrex::Vector<std::array< std::unique_ptr<amrex::MultiFab>, 3>>& Efield,
122  amrex::Vector<std::array< std::unique_ptr<amrex::MultiFab>, 3>> const& Jfield,
123  amrex::Vector<std::unique_ptr<amrex::MultiFab>> const& rhofield,
124  amrex::Vector<std::array< std::unique_ptr<amrex::MultiFab>, 3>> const& edge_lengths,
125  amrex::Real dt, int lev, DtType dt_type,
126  amrex::IntVect ng, std::optional<bool> nodal_sync);
127 
128  void FieldPush (
129  amrex::Vector<std::array< std::unique_ptr<amrex::MultiFab>, 3>>& Bfield,
130  amrex::Vector<std::array< std::unique_ptr<amrex::MultiFab>, 3>>& Efield,
131  amrex::Vector<std::array< std::unique_ptr<amrex::MultiFab>, 3>> const& Jfield,
132  amrex::Vector<std::unique_ptr<amrex::MultiFab>> const& rhofield,
133  amrex::Vector<std::array< std::unique_ptr<amrex::MultiFab>, 3>> const& edge_lengths,
134  amrex::Real dt, DtType dt_type,
135  amrex::IntVect ng, std::optional<bool> nodal_sync);
136 
143  void CalculateElectronPressure ( DtType a_dt_type);
144  void CalculateElectronPressure (int lev, DtType a_dt_type);
145 
155  std::unique_ptr<amrex::MultiFab> const& Pe_field,
156  amrex::MultiFab* const& rho_field ) const;
157 
158  // Declare variables to hold hybrid-PIC model parameters
160  int m_substeps = 10;
161 
163  amrex::Real m_elec_temp;
165  amrex::Real m_n0_ref = 1.0;
167  amrex::Real m_gamma = 5.0/3.0;
168 
170  amrex::Real m_n_floor = 1.0;
171 
173  std::string m_eta_expression = "0.0";
174  std::unique_ptr<amrex::Parser> m_resistivity_parser;
177 
179  amrex::Real m_eta_h = 0.0;
180 
182  std::string m_Jx_ext_grid_function = "0.0";
183  std::string m_Jy_ext_grid_function = "0.0";
184  std::string m_Jz_ext_grid_function = "0.0";
185  std::array< std::unique_ptr<amrex::Parser>, 3> m_J_external_parser;
186  std::array< amrex::ParserExecutor<4>, 3> m_J_external;
188 
189  // Declare multifabs specifically needed for the hybrid-PIC model
195 
196  // Helper functions to retrieve hybrid-PIC multifabs
197  [[nodiscard]] amrex::MultiFab*
198  get_pointer_current_fp_ampere (int lev, int direction) const
199  {
200  return current_fp_ampere[lev][direction].get();
201  }
202 
203  [[nodiscard]] amrex::MultiFab*
204  get_pointer_current_fp_external (int lev, int direction) const
205  {
206  return current_fp_external[lev][direction].get();
207  }
208 
209  [[nodiscard]] amrex::MultiFab*
211  {
212  return electron_pressure_fp[lev].get();
213  }
214 
233 };
234 
242 
244  static amrex::Real get_pressure (amrex::Real const n0,
245  amrex::Real const T0,
246  amrex::Real const gamma,
247  amrex::Real const rho) {
248  return n0 * T0 * std::pow((rho/PhysConst::q_e)/n0, gamma);
249  }
250 };
251 
252 #endif // WARPX_HYBRIDPICMODEL_H_
#define AMREX_FORCE_INLINE
#define AMREX_GPU_HOST_DEVICE
PatchType
Definition: WarpXAlgorithmSelection.H:59
DtType
Definition: WarpXDtType.H:11
This class contains the parameters needed to evaluate hybrid field solutions (kinetic ions with fluid...
Definition: HybridPICModel.H:30
void ReadParameters()
Definition: HybridPICModel.cpp:23
amrex::GpuArray< int, 3 > Bz_IndexType
Definition: HybridPICModel.H:226
amrex::GpuArray< int, 3 > Bx_IndexType
Definition: HybridPICModel.H:222
void FillElectronPressureMF(std::unique_ptr< amrex::MultiFab > const &Pe_field, amrex::MultiFab *const &rho_field) const
Fill the electron pressure multifab given the kinetic particle charge density (and assumption of quas...
Definition: HybridPICModel.cpp:517
void HybridPICSolveE(amrex::Vector< std::array< std::unique_ptr< amrex::MultiFab >, 3 >> &Efield, amrex::Vector< std::array< std::unique_ptr< amrex::MultiFab >, 3 >> const &Jfield, amrex::Vector< std::array< std::unique_ptr< amrex::MultiFab >, 3 >> const &Bfield, amrex::Vector< std::unique_ptr< amrex::MultiFab >> const &rhofield, amrex::Vector< std::array< std::unique_ptr< amrex::MultiFab >, 3 >> const &edge_lengths, bool include_resistivity_term)
Function to update the E-field using Ohm's law (hybrid-PIC model).
Definition: HybridPICModel.cpp:428
void FieldPush(amrex::Vector< std::array< std::unique_ptr< amrex::MultiFab >, 3 >> &Bfield, amrex::Vector< std::array< std::unique_ptr< amrex::MultiFab >, 3 >> &Efield, amrex::Vector< std::array< std::unique_ptr< amrex::MultiFab >, 3 >> const &Jfield, amrex::Vector< std::unique_ptr< amrex::MultiFab >> const &rhofield, amrex::Vector< std::array< std::unique_ptr< amrex::MultiFab >, 3 >> const &edge_lengths, amrex::Real dt, DtType dt_type, amrex::IntVect ng, std::optional< bool > nodal_sync)
Definition: HybridPICModel.cpp:678
amrex::Real m_n0_ref
Definition: HybridPICModel.H:165
void CalculateCurrentAmpere(amrex::Vector< std::array< std::unique_ptr< amrex::MultiFab >, 3 >> const &Bfield, amrex::Vector< std::array< std::unique_ptr< amrex::MultiFab >, 3 >> const &edge_lengths)
Function to calculate the total current based on Ampere's law while neglecting displacement current (...
Definition: HybridPICModel.cpp:398
std::string m_Jz_ext_grid_function
Definition: HybridPICModel.H:184
amrex::GpuArray< int, 3 > Ez_IndexType
Definition: HybridPICModel.H:232
void BfieldEvolveRK(amrex::Vector< std::array< std::unique_ptr< amrex::MultiFab >, 3 >> &Bfield, amrex::Vector< std::array< std::unique_ptr< amrex::MultiFab >, 3 >> &Efield, amrex::Vector< std::array< std::unique_ptr< amrex::MultiFab >, 3 >> const &Jfield, amrex::Vector< std::unique_ptr< amrex::MultiFab >> const &rhofield, amrex::Vector< std::array< std::unique_ptr< amrex::MultiFab >, 3 >> const &edge_lengths, amrex::Real dt, DtType a_dt_type, amrex::IntVect ng, std::optional< bool > nodal_sync)
Definition: HybridPICModel.cpp:546
amrex::MultiFab * get_pointer_electron_pressure_fp(int lev) const
Definition: HybridPICModel.H:210
amrex::MultiFab * get_pointer_current_fp_external(int lev, int direction) const
Definition: HybridPICModel.H:204
std::array< amrex::ParserExecutor< 4 >, 3 > m_J_external
Definition: HybridPICModel.H:186
amrex::Real m_elec_temp
Definition: HybridPICModel.H:163
amrex::Vector< std::array< std::unique_ptr< amrex::MultiFab >, 3 > > current_fp_external
Definition: HybridPICModel.H:193
amrex::GpuArray< int, 3 > Jx_IndexType
Definition: HybridPICModel.H:216
void ClearLevel(int lev)
Definition: HybridPICModel.cpp:118
amrex::GpuArray< int, 3 > Jz_IndexType
Definition: HybridPICModel.H:220
std::string m_Jy_ext_grid_function
Definition: HybridPICModel.H:183
std::unique_ptr< amrex::Parser > m_resistivity_parser
Definition: HybridPICModel.H:174
HybridPICModel(int nlevs_max)
Definition: HybridPICModel.cpp:17
std::string m_eta_expression
Definition: HybridPICModel.H:173
void AllocateMFs(int nlevs_max)
Definition: HybridPICModel.cpp:57
void AllocateLevelMFs(int lev, const amrex::BoxArray &ba, const amrex::DistributionMapping &dm, int ncomps, const amrex::IntVect &ngJ, const amrex::IntVect &ngRho, const amrex::IntVect &jx_nodal_flag, const amrex::IntVect &jy_nodal_flag, const amrex::IntVect &jz_nodal_flag, const amrex::IntVect &rho_nodal_flag)
Definition: HybridPICModel.cpp:66
amrex::Vector< std::unique_ptr< amrex::MultiFab > > electron_pressure_fp
Definition: HybridPICModel.H:194
amrex::Vector< std::array< std::unique_ptr< amrex::MultiFab >, 3 > > current_fp_ampere
Definition: HybridPICModel.H:192
amrex::Real m_gamma
Definition: HybridPICModel.H:167
std::array< std::unique_ptr< amrex::Parser >, 3 > m_J_external_parser
Definition: HybridPICModel.H:185
amrex::Vector< std::unique_ptr< amrex::MultiFab > > rho_fp_temp
Definition: HybridPICModel.H:190
std::string m_Jx_ext_grid_function
Definition: HybridPICModel.H:182
int m_substeps
Definition: HybridPICModel.H:160
bool m_resistivity_has_J_dependence
Definition: HybridPICModel.H:176
amrex::GpuArray< int, 3 > Ex_IndexType
Definition: HybridPICModel.H:228
amrex::Real m_n_floor
Definition: HybridPICModel.H:170
bool m_external_field_has_time_dependence
Definition: HybridPICModel.H:187
void GetCurrentExternal(amrex::Vector< std::array< std::unique_ptr< amrex::MultiFab >, 3 >> const &edge_lengths)
Function to evaluate the external current expressions and populate the external current multifab....
Definition: HybridPICModel.cpp:250
void InitData()
Definition: HybridPICModel.cpp:129
amrex::ParserExecutor< 2 > m_eta
Definition: HybridPICModel.H:175
amrex::MultiFab * get_pointer_current_fp_ampere(int lev, int direction) const
Definition: HybridPICModel.H:198
amrex::Vector< std::array< std::unique_ptr< amrex::MultiFab >, 3 > > current_fp_temp
Definition: HybridPICModel.H:191
void CalculateElectronPressure(DtType a_dt_type)
Function to calculate the electron pressure at a given timestep type using the simulation charge dens...
Definition: HybridPICModel.cpp:488
amrex::GpuArray< int, 3 > Jy_IndexType
Definition: HybridPICModel.H:218
amrex::Real m_eta_h
Definition: HybridPICModel.H:179
amrex::GpuArray< int, 3 > By_IndexType
Definition: HybridPICModel.H:224
amrex::GpuArray< int, 3 > Ey_IndexType
Definition: HybridPICModel.H:230
static constexpr auto q_e
elementary charge [C]
Definition: constant.H:50
float dt
Definition: stencil.py:442
int gamma
boosted frame
Definition: stencil.py:431
This struct contains only static functions to compute the electron pressure using the particle densit...
Definition: HybridPICModel.H:241
AMREX_GPU_HOST_DEVICE static AMREX_FORCE_INLINE amrex::Real get_pressure(amrex::Real const n0, amrex::Real const T0, amrex::Real const gamma, amrex::Real const rho)
Definition: HybridPICModel.H:244