WarpX
Public Member Functions | Private Attributes | List of all members
FiniteDifferenceSolver Class Reference

Top-level class for the electromagnetic finite-difference solver. More...

#include <FiniteDifferenceSolver.H>

Public Member Functions

 FiniteDifferenceSolver (int const fdtd_algo, std::array< amrex::Real, 3 > cell_size, bool const do_nodal)
 Initialize the finite-difference Maxwell solver (for a given refinement level) More...
 
void EvolveB (std::array< std::unique_ptr< amrex::MultiFab >, 3 > &Bfield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Efield, amrex::Real const dt)
 Update the B field, over one timestep. More...
 
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::unique_ptr< amrex::MultiFab > const &Ffield, amrex::Real const dt)
 Update the E field, over one timestep. More...
 
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 const rhocomp, amrex::Real const dt)
 Update the F field, over one timestep. More...
 
void ComputeDivE (const std::array< std::unique_ptr< amrex::MultiFab >, 3 > &Efield, amrex::MultiFab &divE)
 Update the F field, over one timestep. More...
 
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, amrex::Real const dt, std::unique_ptr< MacroscopicProperties > const &macroscopic_properties)
 Macroscopic E-update for non-vacuum medium using the user-selected finite-difference algorithm and macroscopic sigma-method defined in WarpXAlgorithmSelection.H. More...
 
void EvolveBPML (std::array< amrex::MultiFab *, 3 > Bfield, std::array< amrex::MultiFab *, 3 > const Efield, amrex::Real const dt)
 Update the B field, over one timestep. More...
 
void EvolveEPML (std::array< amrex::MultiFab *, 3 > Efield, std::array< amrex::MultiFab *, 3 > const Bfield, std::array< amrex::MultiFab *, 3 > const Jfield, amrex::MultiFab *const Ffield, MultiSigmaBox const &sigba, amrex::Real const dt, bool pml_has_particles)
 Update the E field, over one timestep. More...
 
void EvolveFPML (amrex::MultiFab *Ffield, std::array< amrex::MultiFab *, 3 > const Efield, amrex::Real const dt)
 Update the E field, over one timestep. More...
 
template<typename T_Algo >
void EvolveBCylindrical (std::array< std::unique_ptr< amrex::MultiFab >, 3 > &Bfield, std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &Efield, amrex::Real const dt)
 
template<typename T_Algo >
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::unique_ptr< amrex::MultiFab > const &Ffield, amrex::Real const dt)
 
template<typename T_Algo >
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 const rhocomp, amrex::Real const dt)
 
template<typename T_Algo >
void ComputeDivECylindrical (const std::array< std::unique_ptr< amrex::MultiFab >, 3 > &Efield, amrex::MultiFab &divE)
 

Private Attributes

int m_fdtd_algo
 
bool m_do_nodal
 
amrex::Real m_dr
 
amrex::Real m_rmin
 
amrex::Real m_nmodes
 
amrex::Gpu::DeviceVector< amrex::Real > m_stencil_coefs_r
 
amrex::Gpu::DeviceVector< amrex::Real > m_stencil_coefs_z
 

Detailed Description

Top-level class for the electromagnetic finite-difference solver.

Stores the coefficients of the finite-difference stencils, and has member functions to update fields over one time step.

Constructor & Destructor Documentation

◆ FiniteDifferenceSolver()

FiniteDifferenceSolver::FiniteDifferenceSolver ( int const  fdtd_algo,
std::array< amrex::Real, 3 >  cell_size,
bool const  do_nodal 
)

Initialize the finite-difference Maxwell solver (for a given refinement level)

This function initializes the stencil coefficients for the chosen finite-difference algorithm

Parameters
fdtd_algoIdentifies the chosen algorithm, as defined in WarpXAlgorithmSelection.H
cell_sizeCell size along each dimension, for the chosen refinement level
do_nodalWhether the solver is applied to a nodal or staggered grid

Member Function Documentation

◆ ComputeDivE()

void FiniteDifferenceSolver::ComputeDivE ( const std::array< std::unique_ptr< amrex::MultiFab >, 3 > &  Efield,
amrex::MultiFab &  divE 
)

Update the F field, over one timestep.

◆ ComputeDivECylindrical()

template<typename T_Algo >
void FiniteDifferenceSolver::ComputeDivECylindrical ( const std::array< std::unique_ptr< amrex::MultiFab >, 3 > &  Efield,
amrex::MultiFab &  divE 
)

◆ EvolveB()

void FiniteDifferenceSolver::EvolveB ( std::array< std::unique_ptr< amrex::MultiFab >, 3 > &  Bfield,
std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &  Efield,
amrex::Real const  dt 
)

Update the B field, over one timestep.

◆ EvolveBCylindrical()

template<typename T_Algo >
void FiniteDifferenceSolver::EvolveBCylindrical ( std::array< std::unique_ptr< amrex::MultiFab >, 3 > &  Bfield,
std::array< std::unique_ptr< amrex::MultiFab >, 3 > const &  Efield,
amrex::Real const  dt 
)

◆ EvolveBPML()

void FiniteDifferenceSolver::EvolveBPML ( std::array< amrex::MultiFab *, 3 >  Bfield,
std::array< amrex::MultiFab *, 3 > const  Efield,
amrex::Real const  dt 
)

Update the B field, over one timestep.

◆ EvolveE()

void FiniteDifferenceSolver::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::unique_ptr< amrex::MultiFab > const &  Ffield,
amrex::Real const  dt 
)

Update the E field, over one timestep.

◆ EvolveECylindrical()

template<typename T_Algo >
void FiniteDifferenceSolver::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::unique_ptr< amrex::MultiFab > const &  Ffield,
amrex::Real const  dt 
)

◆ EvolveEPML()

void FiniteDifferenceSolver::EvolveEPML ( std::array< amrex::MultiFab *, 3 >  Efield,
std::array< amrex::MultiFab *, 3 > const  Bfield,
std::array< amrex::MultiFab *, 3 > const  Jfield,
amrex::MultiFab *const  Ffield,
MultiSigmaBox const &  sigba,
amrex::Real const  dt,
bool  pml_has_particles 
)

Update the E field, over one timestep.

◆ EvolveF()

void FiniteDifferenceSolver::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 const  rhocomp,
amrex::Real const  dt 
)

Update the F field, over one timestep.

◆ EvolveFCylindrical()

template<typename T_Algo >
void FiniteDifferenceSolver::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 const  rhocomp,
amrex::Real const  dt 
)

◆ EvolveFPML()

void FiniteDifferenceSolver::EvolveFPML ( amrex::MultiFab *  Ffield,
std::array< amrex::MultiFab *, 3 > const  Efield,
amrex::Real const  dt 
)

Update the E field, over one timestep.

◆ MacroscopicEvolveE()

void FiniteDifferenceSolver::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,
amrex::Real const  dt,
std::unique_ptr< MacroscopicProperties > const &  macroscopic_properties 
)

Macroscopic E-update for non-vacuum medium using the user-selected finite-difference algorithm and macroscopic sigma-method defined in WarpXAlgorithmSelection.H.

Parameters
[out]Efieldvector of electric field MultiFabs updated at a given level
[in]Bfieldvector of magnetic field MultiFabs at a given level
[in]Jfieldvector of current density MultiFabs at a given level
[in]dttimestep of the simulation
[in]macroscopic_propertiescontains user-defined properties of the medium.

Member Data Documentation

◆ m_do_nodal

bool FiniteDifferenceSolver::m_do_nodal
private

◆ m_dr

amrex::Real FiniteDifferenceSolver::m_dr
private

◆ m_fdtd_algo

int FiniteDifferenceSolver::m_fdtd_algo
private

◆ m_nmodes

amrex::Real FiniteDifferenceSolver::m_nmodes
private

◆ m_rmin

amrex::Real FiniteDifferenceSolver::m_rmin
private

◆ m_stencil_coefs_r

amrex::Gpu::DeviceVector<amrex::Real> FiniteDifferenceSolver::m_stencil_coefs_r
private

◆ m_stencil_coefs_z

amrex::Gpu::DeviceVector<amrex::Real> FiniteDifferenceSolver::m_stencil_coefs_z
private

The documentation for this class was generated from the following files: