WarpX
PML_RZ.H
Go to the documentation of this file.
1 /* Copyright 2021 David Grote
2  *
3  *
4  * This file is part of WarpX.
5  *
6  * License: BSD-3-Clause-LBNL
7  */
8 #ifndef WARPX_PML_RZ_H_
9 #define WARPX_PML_RZ_H_
10 
11 #include "PML_RZ_fwd.H"
12 
14 
15 #ifdef WARPX_USE_PSATD
17 #endif
18 
19 #include <AMReX_MultiFab.H>
20 #include <AMReX_BoxArray.H>
21 #include <AMReX_Config.H>
22 #include <AMReX_REAL.H>
23 
24 #include <AMReX_BaseFwd.H>
25 
26 #include <array>
27 #include <optional>
28 #include <string>
29 
30 class PML_RZ
31 {
32 public:
33  PML_RZ (int lev, const amrex::BoxArray& grid_ba, const amrex::DistributionMapping& grid_dm,
34  const amrex::Geometry* geom, int ncell, int do_pml_in_domain);
35 
36  void ApplyDamping(amrex::MultiFab* Et_fp, amrex::MultiFab* Ez_fp,
37  amrex::MultiFab* Bt_fp, amrex::MultiFab* Bz_fp,
38  amrex::Real dt);
39 
40  std::array<amrex::MultiFab*,2> GetE_fp ();
41  std::array<amrex::MultiFab*,2> GetB_fp ();
42 
43 #ifdef WARPX_USE_PSATD
44  void PushPSATD (int lev);
45 #endif
46 
47  void FillBoundaryE ();
48  void FillBoundaryB ();
49  void FillBoundaryE (PatchType patch_type, std::optional<bool> nodal_sync=std::nullopt);
50  void FillBoundaryB (PatchType patch_type, std::optional<bool> nodal_sync=std::nullopt);
51 
52  void CheckPoint (const std::string& dir) const;
53  void Restart (const std::string& dir);
54 
55 private:
56 
57  const int m_ncell;
58  const int m_do_pml_in_domain;
60 
61  // Only contains Er and Et, and Br and Bt
62  std::array<std::unique_ptr<amrex::MultiFab>,2> pml_E_fp;
63  std::array<std::unique_ptr<amrex::MultiFab>,2> pml_B_fp;
64 
65 #ifdef WARPX_USE_PSATD
66  void PushPMLPSATDSinglePatchRZ ( int lev,
67  SpectralSolverRZ& solver,
68  std::array<std::unique_ptr<amrex::MultiFab>,2>& pml_E,
69  std::array<std::unique_ptr<amrex::MultiFab>,2>& pml_B);
70 #endif
71 
72 };
73 
74 #endif
PatchType
Definition: WarpXAlgorithmSelection.H:59
Definition: PML_RZ.H:31
const int m_do_pml_in_domain
Definition: PML_RZ.H:58
void PushPSATD(int lev)
Definition: PML_RZ.cpp:192
std::array< std::unique_ptr< amrex::MultiFab >, 2 > pml_E_fp
Definition: PML_RZ.H:62
void FillBoundaryB()
Definition: PML_RZ.cpp:150
void Restart(const std::string &dir)
Definition: PML_RZ.cpp:179
std::array< amrex::MultiFab *, 2 > GetB_fp()
Definition: PML_RZ.cpp:127
void ApplyDamping(amrex::MultiFab *Et_fp, amrex::MultiFab *Ez_fp, amrex::MultiFab *Bt_fp, amrex::MultiFab *Bz_fp, amrex::Real dt)
Definition: PML_RZ.cpp:63
void FillBoundaryE()
Definition: PML_RZ.cpp:133
void CheckPoint(const std::string &dir) const
Definition: PML_RZ.cpp:167
const amrex::Geometry * m_geom
Definition: PML_RZ.H:59
void PushPMLPSATDSinglePatchRZ(int lev, SpectralSolverRZ &solver, std::array< std::unique_ptr< amrex::MultiFab >, 2 > &pml_E, std::array< std::unique_ptr< amrex::MultiFab >, 2 > &pml_B)
Definition: PML_RZ.cpp:201
std::array< amrex::MultiFab *, 2 > GetE_fp()
Definition: PML_RZ.cpp:121
const int m_ncell
Definition: PML_RZ.H:57
std::array< std::unique_ptr< amrex::MultiFab >, 2 > pml_B_fp
Definition: PML_RZ.H:63
PML_RZ(int lev, const amrex::BoxArray &grid_ba, const amrex::DistributionMapping &grid_dm, const amrex::Geometry *geom, int ncell, int do_pml_in_domain)
Definition: PML_RZ.cpp:39
Definition: SpectralSolverRZ.H:22
float dt
Definition: stencil.py:442