WarpX
|
Top-level class for the electromagnetic spectral solver. More...
#include <SpectralSolver.H>
Public Member Functions | |
SpectralSolver (int lev, const amrex::BoxArray &realspace_ba, const amrex::DistributionMapping &dm, int norder_x, int norder_y, int norder_z, short grid_type, const amrex::Vector< amrex::Real > &v_galilean, const amrex::Vector< amrex::Real > &v_comoving, amrex::RealVect dx, amrex::Real dt, bool pml, bool periodic_single_box, bool update_with_rho, bool fft_do_time_averaging, int psatd_solution_type, int J_in_time, int rho_in_time, bool dive_cleaning, bool divb_cleaning) | |
Constructor of the class SpectralSolver. More... | |
void | ForwardTransform (int lev, const amrex::MultiFab &mf, int field_index, int i_comp=0) |
Transform the component i_comp of the MultiFab mf to Fourier space, and store the result internally (in the spectral field specified by field_index) More... | |
void | BackwardTransform (int lev, amrex::MultiFab &mf, int field_index, const amrex::IntVect &fill_guards, int i_comp=0) |
Transform spectral field specified by field_index back to real space, and store it in the component i_comp of mf More... | |
void | pushSpectralFields () |
Update the fields in spectral space, over one timestep. More... | |
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 SpectralBaseAlgorithm from objects of class SpectralSolver. More... | |
void | CurrentCorrection () |
Public interface to call the virtual function CurrentCorrection , defined in the base class SpectralBaseAlgorithm and possibly overridden by its derived classes (e.g. PsatdAlgorithm, PsatdAlgorithmComoving, etc.), from objects of class SpectralSolver through the private unique pointer algorithm . More... | |
void | VayDeposition () |
Public interface to call the virtual function VayDeposition , declared in the base class SpectralBaseAlgorithm and defined in its derived classes, from objects of class SpectralSolver through the private unique pointer algorithm . More... | |
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 . More... | |
void | ZeroOutDataComp (const int icomp) |
Set to zero the data on component icomp of field_data.fields . More... | |
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. More... | |
Public Attributes | |
SpectralFieldIndex | m_spectral_index |
Protected Attributes | |
amrex::IntVect | m_fill_guards |
Private Member Functions | |
void | ReadParameters () |
Private Attributes | |
SpectralFieldData | field_data |
std::unique_ptr< SpectralBaseAlgorithm > | algorithm |
Top-level class for the electromagnetic spectral solver.
Stores the field in spectral space, and has member functions to Fourier-transform the fields between real space and spectral space and to update fields in spectral space over one time step.
SpectralSolver::SpectralSolver | ( | int | lev, |
const amrex::BoxArray & | realspace_ba, | ||
const amrex::DistributionMapping & | dm, | ||
int | norder_x, | ||
int | norder_y, | ||
int | norder_z, | ||
short | grid_type, | ||
const amrex::Vector< amrex::Real > & | v_galilean, | ||
const amrex::Vector< amrex::Real > & | v_comoving, | ||
amrex::RealVect | dx, | ||
amrex::Real | dt, | ||
bool | pml, | ||
bool | periodic_single_box, | ||
bool | update_with_rho, | ||
bool | fft_do_time_averaging, | ||
int | psatd_solution_type, | ||
int | J_in_time, | ||
int | rho_in_time, | ||
bool | dive_cleaning, | ||
bool | divb_cleaning | ||
) |
Constructor of the class SpectralSolver.
Select the spectral algorithm to be used, allocate the corresponding coefficients for the discrete field update equations, and prepare the structures that store the fields in spectral space.
[in] | lev | mesh refinement level |
[in] | realspace_ba | BoxArray in real space |
[in] | dm | DistributionMapping for the given BoxArray |
[in] | norder_x | spectral order along x |
[in] | norder_y | spectral order along y |
[in] | norder_z | spectral order along z |
[in] | grid_type | type of grid (collocated or not) |
[in] | v_galilean | three-dimensional vector containing the components of the Galilean velocity for the standard or averaged Galilean PSATD solvers |
[in] | v_comoving | three-dimensional vector containing the components of the comoving velocity for the comoving PSATD solver |
[in] | dx | AMREX_SPACEDIM-dimensional vector containing the cell sizes along each direction |
[in] | dt | time step for the analytical integration of Maxwell's equations |
[in] | pml | whether the boxes in the given BoxArray are PML boxes |
[in] | periodic_single_box | whether there is only one periodic single box (no domain decomposition) |
[in] | update_with_rho | whether rho is used in the field update equations |
[in] | fft_do_time_averaging | whether the time averaging algorithm is used |
[in] | psatd_solution_type | whether the PSATD equations are derived from a first-order or second-order model |
[in] | J_in_time | integer that corresponds to the time dependency of J (constant, linear) for the PSATD algorithm |
[in] | rho_in_time | integer that corresponds to the time dependency of rho (linear, quadratic) for the PSATD algorithm |
[in] | dive_cleaning | whether to use div(E) cleaning to account for errors in Gauss law (new field F in the update equations) |
[in] | divb_cleaning | whether to use div(B) cleaning to account for errors in div(B) = 0 law (new field G in the update equations) |
void SpectralSolver::BackwardTransform | ( | int | lev, |
amrex::MultiFab & | mf, | ||
int | field_index, | ||
const amrex::IntVect & | fill_guards, | ||
int | i_comp = 0 |
||
) |
Transform spectral field specified by field_index
back to real space, and store it in the component i_comp
of mf
|
inline |
Public interface to call the member function ComputeSpectralDivE of the base class SpectralBaseAlgorithm from objects of class SpectralSolver.
Copy spectral data from component src_comp
to component dest_comp
of field_data.fields
.
[in] | src_comp | component of the source FabArray from which the data are copied |
[in] | dest_comp | component of the destination FabArray where the data are copied |
|
inline |
Public interface to call the virtual function CurrentCorrection
, defined in the base class SpectralBaseAlgorithm and possibly overridden by its derived classes (e.g. PsatdAlgorithm, PsatdAlgorithmComoving, etc.), from objects of class SpectralSolver through the private unique pointer algorithm
.
void SpectralSolver::ForwardTransform | ( | int | lev, |
const amrex::MultiFab & | mf, | ||
int | field_index, | ||
int | i_comp = 0 |
||
) |
Transform the component i_comp of the MultiFab mf to Fourier space, and store the result internally (in the spectral field specified by field_index)
[in] | lev | mesh refinement level |
[in] | mf | MultiFab that is transformed to Fourier space (component i_comp) |
[in] | field_index | index of the spectral field that stores the FFT result |
[in] | i_comp | component of the MultiFab mf that is transformed to Fourier space |
void SpectralSolver::pushSpectralFields | ( | ) |
Update the fields in spectral space, over one timestep.
|
private |
|
inline |
Scale the data on component icomp
of field_data.fields
by a given scale factor.
[in] | icomp | component of the FabArray where the data are scaled |
[in] | scale_factor | scale factor to use for scaling |
|
inline |
Public interface to call the virtual function VayDeposition
, declared in the base class SpectralBaseAlgorithm and defined in its derived classes, from objects of class SpectralSolver through the private unique pointer algorithm
.
|
inline |
Set to zero the data on component icomp
of field_data.fields
.
[in] | icomp | component of the FabArray where the data are set to zero |
|
private |
|
private |
|
protected |
SpectralFieldIndex SpectralSolver::m_spectral_index |