WarpX
FiniteDifferenceSolver.H
Go to the documentation of this file.
1 /* Copyright 2020 Remi Lehe
2  *
3  * This file is part of WarpX.
4  *
5  * License: BSD-3-Clause-LBNL
6  */
7 
8 #ifndef WARPX_FINITE_DIFFERENCE_SOLVER_H_
9 #define WARPX_FINITE_DIFFERENCE_SOLVER_H_
10 
14 
16 #include "Evolve/WarpXDtType.H"
19 
20 #include <AMReX_GpuContainers.H>
21 #include <AMReX_REAL.H>
22 
23 #include <AMReX_BaseFwd.H>
24 
25 #include <array>
26 #include <memory>
27 
35 {
36  public:
37 
38  // Constructor
48  int fdtd_algo,
49  std::array<amrex::Real,3> cell_size,
50  short grid_type );
51 
52  void EvolveB ( std::array< std::unique_ptr<amrex::MultiFab>, 3 >& Bfield,
53  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& Efield,
54  std::unique_ptr<amrex::MultiFab> const& Gfield,
55  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& face_areas,
56  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& area_mod,
57  std::array< std::unique_ptr<amrex::MultiFab>, 3 >& ECTRhofield,
58  std::array< std::unique_ptr<amrex::MultiFab>, 3 >& Venl,
59  std::array< std::unique_ptr<amrex::iMultiFab>, 3 >& flag_info_cell,
60  std::array< std::unique_ptr<amrex::LayoutData<FaceInfoBox> >, 3 >& borrowing,
61  int lev, amrex::Real dt );
62 
63  void EvolveE ( std::array< std::unique_ptr<amrex::MultiFab>, 3 >& Efield,
64  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& Bfield,
65  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& Jfield,
66  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& edge_lengths,
67  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& face_areas,
68  std::array< std::unique_ptr<amrex::MultiFab>, 3 >& ECTRhofield,
69  std::unique_ptr<amrex::MultiFab> const& Ffield,
70  int lev, amrex::Real dt );
71 
72  void EvolveF ( std::unique_ptr<amrex::MultiFab>& Ffield,
73  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& Efield,
74  std::unique_ptr<amrex::MultiFab> const& rhofield,
75  int rhocomp,
76  amrex::Real dt );
77 
78  void EvolveG (std::unique_ptr<amrex::MultiFab>& Gfield,
79  std::array<std::unique_ptr<amrex::MultiFab>,3> const& Bfield,
80  amrex::Real dt);
81 
82  void EvolveECTRho ( std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& Efield,
83  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& edge_lengths,
84  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& face_areas,
85  std::array< std::unique_ptr<amrex::MultiFab>, 3 >& ECTRhofield,
86  int lev );
87 
89  std::array< std::unique_ptr<amrex::MultiFab>, 3 >& Efield,
90  std::array< std::unique_ptr<amrex::MultiFab>, 3 >& Bfield,
91  amrex::Box domain_box,
92  amrex::Real dt,
93  amrex::Vector<FieldBoundaryType> field_boundary_lo,
94  amrex::Vector<FieldBoundaryType> field_boundary_hi);
95 
96  void ComputeDivE ( const std::array<std::unique_ptr<amrex::MultiFab>,3>& Efield,
97  amrex::MultiFab& divE );
98 
111  void MacroscopicEvolveE ( std::array< std::unique_ptr<amrex::MultiFab>, 3>& Efield,
112  std::array< std::unique_ptr<amrex::MultiFab>, 3> const& Bfield,
113  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& Jfield,
114  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& edge_lengths,
115  amrex::Real dt,
116  std::unique_ptr<MacroscopicProperties> const& macroscopic_properties);
117 
118  void EvolveBPML ( std::array< amrex::MultiFab*, 3 > Bfield,
119  std::array< amrex::MultiFab*, 3 > Efield,
120  amrex::Real dt,
121  bool dive_cleaning);
122 
123  void EvolveEPML ( std::array< amrex::MultiFab*, 3 > Efield,
124  std::array< amrex::MultiFab*, 3 > Bfield,
125  std::array< amrex::MultiFab*, 3 > Jfield,
126  std::array< amrex::MultiFab*, 3 > edge_lengths,
127  amrex::MultiFab* Ffield,
128  MultiSigmaBox const& sigba,
129  amrex::Real dt, bool pml_has_particles );
130 
131  void EvolveFPML ( amrex::MultiFab* Ffield,
132  std::array< amrex::MultiFab*, 3 > Efield,
133  amrex::Real dt );
134 
152  void HybridPICSolveE ( std::array< std::unique_ptr<amrex::MultiFab>, 3>& Efield,
153  std::array< std::unique_ptr<amrex::MultiFab>, 3>& Jfield,
154  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& Jifield,
155  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& Jextfield,
156  std::array< std::unique_ptr<amrex::MultiFab>, 3> const& Bfield,
157  std::unique_ptr<amrex::MultiFab> const& rhofield,
158  std::unique_ptr<amrex::MultiFab> const& Pefield,
159  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& edge_lengths,
160  int lev, HybridPICModel const* hybrid_model,
161  bool include_resistivity_term );
162 
173  std::array< std::unique_ptr<amrex::MultiFab>, 3>& Jfield,
174  std::array< std::unique_ptr<amrex::MultiFab>, 3> const& Bfield,
175  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& edge_lengths,
176  int lev );
177 
178  private:
179 
181  short m_grid_type;
182 
183 #ifdef WARPX_DIM_RZ
184  amrex::Real m_dr, m_rmin;
185  int m_nmodes;
186  // host-only
188  // device copy after init
191 #else
192  // host-only
193  amrex::Vector<amrex::Real> m_h_stencil_coefs_x, m_h_stencil_coefs_y, m_h_stencil_coefs_z;
194  // device copy after init
195  amrex::Gpu::DeviceVector<amrex::Real> m_stencil_coefs_x;
196  amrex::Gpu::DeviceVector<amrex::Real> m_stencil_coefs_y;
198 #endif
199 
200  public:
201  // The member functions below contain extended __device__ lambda.
202  // In order to compile with nvcc, they need to be public.
203 
204 #ifdef WARPX_DIM_RZ
205  template< typename T_Algo >
206  void EvolveBCylindrical (
207  std::array< std::unique_ptr<amrex::MultiFab>, 3 >& Bfield,
208  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& Efield,
209  int lev,
210  amrex::Real dt );
211 
212  template< typename T_Algo >
213  void EvolveECylindrical (
214  std::array< std::unique_ptr<amrex::MultiFab>, 3 >& Efield,
215  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& Bfield,
216  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& Jfield,
217  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& edge_lengths,
218  std::unique_ptr<amrex::MultiFab> const& Ffield,
219  int lev,
220  amrex::Real dt );
221 
222  template< typename T_Algo >
223  void EvolveFCylindrical (
224  std::unique_ptr<amrex::MultiFab>& Ffield,
225  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& Efield,
226  std::unique_ptr<amrex::MultiFab> const& rhofield,
227  int rhocomp,
228  amrex::Real dt );
229 
230  template< typename T_Algo >
232  const std::array<std::unique_ptr<amrex::MultiFab>,3>& Efield,
233  amrex::MultiFab& divE );
234 
235  template<typename T_Algo>
237  std::array< std::unique_ptr<amrex::MultiFab>, 3>& Efield,
238  std::array< std::unique_ptr<amrex::MultiFab>, 3> const& Jfield,
239  std::array< std::unique_ptr<amrex::MultiFab>, 3> const& Jifield,
240  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& Jextfield,
241  std::array< std::unique_ptr<amrex::MultiFab>, 3> const& Bfield,
242  std::unique_ptr<amrex::MultiFab> const& rhofield,
243  std::unique_ptr<amrex::MultiFab> const& Pefield,
244  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& edge_lengths,
245  int lev, HybridPICModel const* hybrid_model,
246  bool include_resistivity_term );
247 
248  template<typename T_Algo>
250  std::array< std::unique_ptr<amrex::MultiFab>, 3 >& Jfield,
251  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& Bfield,
252  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& edge_lengths,
253  int lev
254  );
255 
256 #else
257  template< typename T_Algo >
258  void EvolveBCartesian (
259  std::array< std::unique_ptr<amrex::MultiFab>, 3 >& Bfield,
260  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& Efield,
261  std::unique_ptr<amrex::MultiFab> const& Gfield,
262  int lev, amrex::Real dt );
263 
264  template< typename T_Algo >
265  void EvolveECartesian (
266  std::array< std::unique_ptr<amrex::MultiFab>, 3 >& Efield,
267  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& Bfield,
268  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& Jfield,
269  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& edge_lengths,
270  std::unique_ptr<amrex::MultiFab> const& Ffield,
271  int lev, amrex::Real dt );
272 
273  template< typename T_Algo >
274  void EvolveFCartesian (
275  std::unique_ptr<amrex::MultiFab>& Ffield,
276  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& Efield,
277  std::unique_ptr<amrex::MultiFab> const& rhofield,
278  int rhocomp,
279  amrex::Real dt );
280 
281  template< typename T_Algo >
282  void EvolveGCartesian (
283  std::unique_ptr<amrex::MultiFab>& Gfield,
284  std::array<std::unique_ptr<amrex::MultiFab>,3> const& Bfield,
285  amrex::Real dt);
286 
287  void EvolveRhoCartesianECT (
288  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& Efield,
289  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& edge_lengths,
290  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& face_areas,
291  std::array< std::unique_ptr<amrex::MultiFab>, 3 >& ECTRhofield, int lev);
292 
293  void EvolveBCartesianECT (
294  std::array< std::unique_ptr<amrex::MultiFab>, 3 >& Bfield,
295  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& face_areas,
296  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& area_mod,
297  std::array< std::unique_ptr<amrex::MultiFab>, 3 >& ECTRhofield,
298  std::array< std::unique_ptr<amrex::MultiFab>, 3 >& Venl,
299  std::array< std::unique_ptr<amrex::iMultiFab>, 3 >& flag_info_cell,
300  std::array< std::unique_ptr<amrex::LayoutData<FaceInfoBox> >, 3 >& borrowing,
301  int lev, amrex::Real dt
302  );
303 
304  template< typename T_Algo >
305  void ComputeDivECartesian (
306  const std::array<std::unique_ptr<amrex::MultiFab>,3>& Efield,
307  amrex::MultiFab& divE );
308 
309  template< typename T_Algo, typename T_MacroAlgo >
310  void MacroscopicEvolveECartesian (
311  std::array< std::unique_ptr< amrex::MultiFab>, 3>& Efield,
312  std::array< std::unique_ptr< amrex::MultiFab>, 3> const& Bfield,
313  std::array< std::unique_ptr< amrex::MultiFab>, 3> const& Jfield,
314  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& edge_lengths,
315  amrex::Real dt,
316  std::unique_ptr<MacroscopicProperties> const& macroscopic_properties);
317 
318  template< typename T_Algo >
319  void EvolveBPMLCartesian (
320  std::array< amrex::MultiFab*, 3 > Bfield,
321  std::array< amrex::MultiFab*, 3 > Efield,
322  amrex::Real dt,
323  bool dive_cleaning);
324 
325  template< typename T_Algo >
326  void EvolveEPMLCartesian (
327  std::array< amrex::MultiFab*, 3 > Efield,
328  std::array< amrex::MultiFab*, 3 > Bfield,
329  std::array< amrex::MultiFab*, 3 > Jfield,
330  std::array< amrex::MultiFab*, 3 > edge_lengths,
331  amrex::MultiFab* Ffield,
332  MultiSigmaBox const& sigba,
333  amrex::Real dt, bool pml_has_particles );
334 
335  template< typename T_Algo >
336  void EvolveFPMLCartesian ( amrex::MultiFab* Ffield,
337  std::array< amrex::MultiFab*, 3 > Efield,
338  amrex::Real dt );
339 
340  template<typename T_Algo>
341  void HybridPICSolveECartesian (
342  std::array< std::unique_ptr<amrex::MultiFab>, 3>& Efield,
343  std::array< std::unique_ptr<amrex::MultiFab>, 3> const& Jfield,
344  std::array< std::unique_ptr<amrex::MultiFab>, 3> const& Jifield,
345  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& Jextfield,
346  std::array< std::unique_ptr<amrex::MultiFab>, 3> const& Bfield,
347  std::unique_ptr<amrex::MultiFab> const& rhofield,
348  std::unique_ptr<amrex::MultiFab> const& Pefield,
349  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& edge_lengths,
350  int lev, HybridPICModel const* hybrid_model,
351  bool include_resistivity_term );
352 
353  template<typename T_Algo>
354  void CalculateCurrentAmpereCartesian (
355  std::array< std::unique_ptr<amrex::MultiFab>, 3 >& Jfield,
356  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& Bfield,
357  std::array< std::unique_ptr<amrex::MultiFab>, 3 > const& edge_lengths,
358  int lev
359  );
360 #endif
361 
362 };
363 
364 #endif // WARPX_FINITE_DIFFERENCE_SOLVER_H_
Top-level class for the electromagnetic finite-difference solver.
Definition: FiniteDifferenceSolver.H:35
int m_nmodes
Definition: FiniteDifferenceSolver.H:185
void HybridPICSolveECylindrical(std::array< std::unique_ptr< amrex::MultiFab >, 3 > &Efield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Jfield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Jifield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Jextfield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Bfield, std::unique_ptr< amrex::MultiFab > const &rhofield, std::unique_ptr< amrex::MultiFab > const &Pefield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &edge_lengths, int lev, HybridPICModel const *hybrid_model, bool include_resistivity_term)
Definition: HybridPICSolveE.cpp:406
void MacroscopicEvolveE(std::array< std::unique_ptr< amrex::MultiFab >, 3 > &Efield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Bfield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Jfield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &edge_lengths, amrex::Real dt, std::unique_ptr< MacroscopicProperties > const &macroscopic_properties)
Macroscopic E-update for non-vacuum medium using the user-selected finite-difference algorithm and ma...
Definition: MacroscopicEvolveE.cpp:37
void CalculateCurrentAmpere(std::array< std::unique_ptr< amrex::MultiFab >, 3 > &Jfield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Bfield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &edge_lengths, int lev)
Calculation of total current using Ampere's law (without displacement current): J = (curl x B) / mu0.
Definition: HybridPICSolveE.cpp:25
void EvolveG(std::unique_ptr< amrex::MultiFab > &Gfield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Bfield, amrex::Real dt)
Definition: EvolveG.cpp:40
void EvolveBPML(std::array< amrex::MultiFab *, 3 > Bfield, std::array< amrex::MultiFab *, 3 > Efield, amrex::Real dt, bool dive_cleaning)
Update the B field, over one timestep.
Definition: EvolveBPML.cpp:43
void EvolveFCylindrical(std::unique_ptr< amrex::MultiFab > &Ffield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Efield, std::unique_ptr< amrex::MultiFab > const &rhofield, int rhocomp, amrex::Real dt)
Definition: EvolveF.cpp:137
void EvolveECylindrical(std::array< std::unique_ptr< amrex::MultiFab >, 3 > &Efield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Bfield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Jfield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &edge_lengths, std::unique_ptr< amrex::MultiFab > const &Ffield, int lev, amrex::Real dt)
Definition: EvolveE.cpp:232
void HybridPICSolveE(std::array< std::unique_ptr< amrex::MultiFab >, 3 > &Efield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > &Jfield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Jifield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Jextfield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Bfield, std::unique_ptr< amrex::MultiFab > const &rhofield, std::unique_ptr< amrex::MultiFab > const &Pefield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &edge_lengths, int lev, HybridPICModel const *hybrid_model, bool include_resistivity_term)
E-update in the hybrid PIC algorithm as described in Winske et al. (2003) Eq. 10. https://link....
Definition: HybridPICSolveE.cpp:368
void EvolveF(std::unique_ptr< amrex::MultiFab > &Ffield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Efield, std::unique_ptr< amrex::MultiFab > const &rhofield, int rhocomp, amrex::Real dt)
Update the F field, over one timestep.
Definition: EvolveF.cpp:46
void ComputeDivECylindrical(const std::array< std::unique_ptr< amrex::MultiFab >, 3 > &Efield, amrex::MultiFab &divE)
Definition: ComputeDivE.cpp:125
amrex::Gpu::DeviceVector< amrex::Real > m_stencil_coefs_z
Definition: FiniteDifferenceSolver.H:190
int m_fdtd_algo
Definition: FiniteDifferenceSolver.H:180
void EvolveE(std::array< std::unique_ptr< amrex::MultiFab >, 3 > &Efield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Bfield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Jfield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &edge_lengths, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &face_areas, std::array< std::unique_ptr< amrex::MultiFab >, 3 > &ECTRhofield, std::unique_ptr< amrex::MultiFab > const &Ffield, int lev, amrex::Real dt)
Update the E field, over one timestep.
Definition: EvolveE.cpp:48
void ComputeDivE(const std::array< std::unique_ptr< amrex::MultiFab >, 3 > &Efield, amrex::MultiFab &divE)
Update the F field, over one timestep.
Definition: ComputeDivE.cpp:42
void EvolveECTRho(std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Efield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &edge_lengths, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &face_areas, std::array< std::unique_ptr< amrex::MultiFab >, 3 > &ECTRhofield, int lev)
Update the B field, over one timestep.
Definition: EvolveECTRho.cpp:49
void EvolveFPML(amrex::MultiFab *Ffield, std::array< amrex::MultiFab *, 3 > Efield, amrex::Real dt)
Update the E field, over one timestep.
Definition: EvolveFPML.cpp:42
void ApplySilverMuellerBoundary(std::array< std::unique_ptr< amrex::MultiFab >, 3 > &Efield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > &Bfield, amrex::Box domain_box, amrex::Real dt, amrex::Vector< FieldBoundaryType > field_boundary_lo, amrex::Vector< FieldBoundaryType > field_boundary_hi)
Update the B field at the boundary, using the Silver-Mueller condition.
Definition: ApplySilverMuellerBoundary.cpp:37
void EvolveEPML(std::array< amrex::MultiFab *, 3 > Efield, std::array< amrex::MultiFab *, 3 > Bfield, std::array< amrex::MultiFab *, 3 > Jfield, std::array< amrex::MultiFab *, 3 > edge_lengths, amrex::MultiFab *Ffield, MultiSigmaBox const &sigba, amrex::Real dt, bool pml_has_particles)
Update the E field, over one timestep.
Definition: EvolveEPML.cpp:46
amrex::Vector< amrex::Real > m_h_stencil_coefs_z
Definition: FiniteDifferenceSolver.H:187
short m_grid_type
Definition: FiniteDifferenceSolver.H:181
amrex::Gpu::DeviceVector< amrex::Real > m_stencil_coefs_r
Definition: FiniteDifferenceSolver.H:189
amrex::Real m_dr
Definition: FiniteDifferenceSolver.H:184
void EvolveB(std::array< std::unique_ptr< amrex::MultiFab >, 3 > &Bfield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Efield, std::unique_ptr< amrex::MultiFab > const &Gfield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &face_areas, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &area_mod, std::array< std::unique_ptr< amrex::MultiFab >, 3 > &ECTRhofield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > &Venl, std::array< std::unique_ptr< amrex::iMultiFab >, 3 > &flag_info_cell, std::array< std::unique_ptr< amrex::LayoutData< FaceInfoBox > >, 3 > &borrowing, int lev, amrex::Real dt)
Update the B field, over one timestep.
Definition: EvolveB.cpp:50
amrex::Vector< amrex::Real > m_h_stencil_coefs_r
Definition: FiniteDifferenceSolver.H:187
void CalculateCurrentAmpereCylindrical(std::array< std::unique_ptr< amrex::MultiFab >, 3 > &Jfield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Bfield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &edge_lengths, int lev)
Definition: HybridPICSolveE.cpp:60
FiniteDifferenceSolver(int fdtd_algo, std::array< amrex::Real, 3 > cell_size, short grid_type)
Initialize the finite-difference Maxwell solver (for a given refinement level)
Definition: FiniteDifferenceSolver.cpp:30
amrex::Real m_rmin
Definition: FiniteDifferenceSolver.H:184
void EvolveBCylindrical(std::array< std::unique_ptr< amrex::MultiFab >, 3 > &Bfield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Efield, int lev, amrex::Real dt)
Definition: EvolveB.cpp:372
This class contains the parameters needed to evaluate hybrid field solutions (kinetic ions with fluid...
Definition: HybridPICModel.H:30
Definition: PML.H:123
cell_size
Definition: compute_domain.py:37
float dt
Definition: stencil.py:442