WarpX
SpectralSolverRZ.H
Go to the documentation of this file.
1 /* Copyright 2019-2020 David Grote
2  *
3  * This file is part of WarpX.
4  *
5  * License: BSD-3-Clause-LBNL
6  */
7 #ifndef WARPX_SPECTRAL_SOLVER_RZ_H_
8 #define WARPX_SPECTRAL_SOLVER_RZ_H_
9 
10 #include "SpectralSolverRZ_fwd.H"
11 
13 #include "SpectralFieldDataRZ.H"
14 
15 /* \brief Top-level class for the electromagnetic spectral solver
16  *
17  * Stores the field in spectral space, and has member functions
18  * to Fourier-transform the fields between real space and spectral space
19  * and to update fields in spectral space over one time step.
20  */
22 {
23  public:
24  // Inline definition of the member functions of `SpectralSolverRZ`,
25  // except the constructor (see `SpectralSolverRZ.cpp`)
26  // The body of these functions is short, since the work is done in the
27  // underlying classes `SpectralFieldData` and `PsatdAlgorithm`
28 
29  // Constructor
30  SpectralSolverRZ (int lev,
31  amrex::BoxArray const & realspace_ba,
32  amrex::DistributionMapping const & dm,
33  int n_rz_azimuthal_modes,
34  int norder_z, short grid_type,
35  const amrex::Vector<amrex::Real>& v_galilean,
36  amrex::RealVect dx, amrex::Real dt,
37  bool with_pml,
38  bool update_with_rho,
39  bool fft_do_time_averaging,
40  int J_in_time,
41  int rho_in_time,
42  bool dive_cleaning,
43  bool divb_cleaning);
44 
45  /* \brief Transform the component `i_comp` of MultiFab `field_mf`
46  * to spectral space, and store the corresponding result internally
47  * (in the spectral field specified by `field_index`) */
48  void ForwardTransform (int lev, amrex::MultiFab const & field_mf, int field_index,
49  int i_comp=0);
50 
51  /* \brief Transform the two MultiFabs `field_mf1` and `field_mf2`
52  * to spectral space, and store the corresponding results internally
53  * (in the spectral field specified by `field_index1` and `field_index2`) */
54  void ForwardTransform (int lev, amrex::MultiFab const & field_mf1, int field_index1,
55  amrex::MultiFab const & field_mf2, int field_index2);
56 
57  /* \brief Transform spectral field specified by `field_index` back to
58  * real space, and store it in the component `i_comp` of `field_mf` */
59  void BackwardTransform (int lev, amrex::MultiFab& field_mf, int field_index,
60  int i_comp=0);
61 
62  /* \brief Transform spectral fields specified by `field_index1` and `field_index2`
63  * back to real space, and store it in `field_mf1` and `field_mf2`*/
64  void BackwardTransform (int lev, amrex::MultiFab& field_mf1, int field_index1,
65  amrex::MultiFab& field_mf2, int field_index2);
66 
67  /* \brief Update the fields in spectral space, over one timestep */
68  void pushSpectralFields (bool doing_pml=false);
69 
70  /* \brief Initialize K space filtering arrays */
71  void InitFilter (amrex::IntVect const & filter_npass_each_dir,
72  bool const compensation)
73  {
74  field_data.InitFilter(filter_npass_each_dir, compensation, k_space);
75  }
76 
77  /* \brief Apply K space filtering for a scalar */
78  void ApplyFilter (const int lev, int const field_index)
79  {
80  field_data.ApplyFilter(lev, field_index);
81  }
82 
83  /* \brief Apply K space filtering for a vector */
84  void ApplyFilter (const int lev, int const field_index1,
85  int const field_index2, int const field_index3)
86  {
87  field_data.ApplyFilter(lev, field_index1, field_index2, field_index3);
88  }
89 
94  void ComputeSpectralDivE (int lev, const std::array<std::unique_ptr<amrex::MultiFab>,3>& Efield,
95  amrex::MultiFab& divE);
96 
103  void CurrentCorrection ();
104 
111  void VayDeposition ();
112 
120  void CopySpectralDataComp (const int src_comp, const int dest_comp)
121  {
122  field_data.CopySpectralDataComp(src_comp, dest_comp);
123  }
124 
130  void ZeroOutDataComp (const int icomp)
131  {
133  }
134 
142  void ScaleDataComp (const int icomp, const amrex::Real scale_factor)
143  {
144  field_data.ScaleDataComp(icomp, scale_factor);
145  }
146 
148 
149  private:
150 
151  SpectralKSpaceRZ k_space; // Save the instance to initialize filtering
152  SpectralFieldDataRZ field_data; // Store field in spectral space
153  // and perform the Fourier transforms
154  std::unique_ptr<SpectralBaseAlgorithmRZ> algorithm;
155  std::unique_ptr<SpectralBaseAlgorithmRZ> PML_algorithm;
156  // Defines field update equation in spectral space,
157  // and the associated coefficients.
158  // SpectralBaseAlgorithmRZ is a base class ; this pointer is meant
159  // to point an instance of a *sub-class* defining a specific algorithm
160 
161 };
162 
163 #endif // WARPX_SPECTRAL_SOLVER_RZ_H_
Definition: SpectralFieldDataRZ.H:23
void ZeroOutDataComp(int icomp)
Set to zero the data on component icomp of fields.
Definition: SpectralFieldDataRZ.H:100
void CopySpectralDataComp(int src_comp, int dest_comp)
Copy spectral data from component src_comp to component dest_comp of fields.
Definition: SpectralFieldDataRZ.H:80
void InitFilter(amrex::IntVect const &filter_npass_each_dir, bool compensation, SpectralKSpaceRZ const &k_space)
Definition: SpectralFieldDataRZ.cpp:756
void ScaleDataComp(int icomp, amrex::Real scale_factor)
Scale the data on component icomp of fields by a given scale factor.
Definition: SpectralFieldDataRZ.H:120
void ApplyFilter(int lev, int field_index)
Definition: SpectralFieldDataRZ.cpp:773
Definition: SpectralFieldData.H:34
Definition: SpectralKSpaceRZ.H:21
Definition: SpectralSolverRZ.H:22
SpectralFieldIndex m_spectral_index
Definition: SpectralSolverRZ.H:147
void ApplyFilter(const int lev, int const field_index1, int const field_index2, int const field_index3)
Definition: SpectralSolverRZ.H:84
void CopySpectralDataComp(const int src_comp, const int dest_comp)
Copy spectral data from component src_comp to component dest_comp of field_data.fields.
Definition: SpectralSolverRZ.H:120
void ZeroOutDataComp(const int icomp)
Set to zero the data on component icomp of field_data.fields.
Definition: SpectralSolverRZ.H:130
void ApplyFilter(const int lev, int const field_index)
Definition: SpectralSolverRZ.H:78
SpectralSolverRZ(int lev, amrex::BoxArray const &realspace_ba, amrex::DistributionMapping const &dm, int n_rz_azimuthal_modes, int norder_z, short grid_type, const amrex::Vector< amrex::Real > &v_galilean, amrex::RealVect dx, amrex::Real dt, bool with_pml, bool update_with_rho, bool fft_do_time_averaging, int J_in_time, int rho_in_time, bool dive_cleaning, bool divb_cleaning)
Definition: SpectralSolverRZ.cpp:28
void InitFilter(amrex::IntVect const &filter_npass_each_dir, bool const compensation)
Definition: SpectralSolverRZ.H:71
void ScaleDataComp(const int icomp, const amrex::Real scale_factor)
Scale the data on component icomp of field_data.fields by a given scale factor.
Definition: SpectralSolverRZ.H:142
SpectralKSpaceRZ k_space
Definition: SpectralSolverRZ.H:151
SpectralFieldDataRZ field_data
Definition: SpectralSolverRZ.H:152
void VayDeposition()
Public interface to call the virtual function VayDeposition, declared in the base class SpectralBaseA...
Definition: SpectralSolverRZ.cpp:162
void BackwardTransform(int lev, amrex::MultiFab &field_mf, int field_index, int i_comp=0)
Definition: SpectralSolverRZ.cpp:105
std::unique_ptr< SpectralBaseAlgorithmRZ > algorithm
Definition: SpectralSolverRZ.H:154
void pushSpectralFields(bool doing_pml=false)
Definition: SpectralSolverRZ.cpp:126
void ForwardTransform(int lev, amrex::MultiFab const &field_mf, int field_index, int i_comp=0)
Definition: SpectralSolverRZ.cpp:82
std::unique_ptr< SpectralBaseAlgorithmRZ > PML_algorithm
Definition: SpectralSolverRZ.H:155
void CurrentCorrection()
Public interface to call the virtual function CurrentCorrection, defined in the base class SpectralBa...
Definition: SpectralSolverRZ.cpp:156
void ComputeSpectralDivE(int lev, const std::array< std::unique_ptr< amrex::MultiFab >, 3 > &Efield, amrex::MultiFab &divE)
Public interface to call the member function ComputeSpectralDivE of the base class SpectralBaseAlgori...
Definition: SpectralSolverRZ.cpp:143
float dt
Definition: stencil.py:442
tuple dx
lab frame
Definition: stencil.py:429