|
void | computePhiIGF (amrex::MultiFab const &rho, amrex::MultiFab &phi, std::array< amrex::Real, 3 > const &cell_size, amrex::BoxArray const &ba) |
| Compute the electrostatic potential using the Integrated Green Function method as in http://dx.doi.org/10.1103/PhysRevSTAB.9.044204. More...
|
|
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real | IntegratedPotential (amrex::Real x, amrex::Real y, amrex::Real z) |
| Implements equation 2 in https://doi.org/10.1103/PhysRevSTAB.10.129901 with some modification to symmetrize the function. More...
|
|
template<typename T_BoundaryHandler , typename T_PostPhiCalculationFunctor = std::nullopt_t, typename T_FArrayBoxFactory = void> |
void | computePhi (amrex::Vector< amrex::MultiFab * > const &rho, amrex::Vector< amrex::MultiFab * > &phi, std::array< amrex::Real, 3 > const beta, amrex::Real const relative_tolerance, amrex::Real absolute_tolerance, int const max_iters, int const verbosity, amrex::Vector< amrex::Geometry > const &geom, amrex::Vector< amrex::DistributionMapping > const &dmap, amrex::Vector< amrex::BoxArray > const &grids, T_BoundaryHandler const boundary_handler, [[maybe_unused]] bool is_solver_multigrid, bool const do_single_precision_comms=false, std::optional< amrex::Vector< amrex::IntVect > > rel_ref_ratio=std::nullopt, [[maybe_unused]] T_PostPhiCalculationFunctor post_phi_calculation=std::nullopt, [[maybe_unused]] std::optional< amrex::Real const > current_time=std::nullopt, [[maybe_unused]] std::optional< amrex::Vector< T_FArrayBoxFactory const * > > eb_farray_box_factory=std::nullopt) |
|
template<typename T_BoundaryHandler , typename T_PostACalculationFunctor = std::nullopt_t, typename T_FArrayBoxFactory = void> |
void | computeVectorPotential (amrex::Vector< amrex::Array< amrex::MultiFab *, 3 > > const &curr, amrex::Vector< amrex::Array< amrex::MultiFab *, 3 > > &A, amrex::Real const relative_tolerance, amrex::Real absolute_tolerance, int const max_iters, int const verbosity, amrex::Vector< amrex::Geometry > const &geom, amrex::Vector< amrex::DistributionMapping > const &dmap, amrex::Vector< amrex::BoxArray > const &grids, T_BoundaryHandler const boundary_handler, bool const do_single_precision_comms=false, std::optional< amrex::Vector< amrex::IntVect > > rel_ref_ratio=std::nullopt, [[maybe_unused]] T_PostACalculationFunctor post_A_calculation=std::nullopt, [[maybe_unused]] std::optional< amrex::Real const > current_time=std::nullopt, [[maybe_unused]] std::optional< amrex::Vector< T_FArrayBoxFactory const * > > eb_farray_box_factory=std::nullopt) |
|
template<typename T_BoundaryHandler , typename T_PostPhiCalculationFunctor = std::nullopt_t, typename T_FArrayBoxFactory = void>
void ablastr::fields::computePhi |
( |
amrex::Vector< amrex::MultiFab * > const & |
rho, |
|
|
amrex::Vector< amrex::MultiFab * > & |
phi, |
|
|
std::array< amrex::Real, 3 > const |
beta, |
|
|
amrex::Real const |
relative_tolerance, |
|
|
amrex::Real |
absolute_tolerance, |
|
|
int const |
max_iters, |
|
|
int const |
verbosity, |
|
|
amrex::Vector< amrex::Geometry > const & |
geom, |
|
|
amrex::Vector< amrex::DistributionMapping > const & |
dmap, |
|
|
amrex::Vector< amrex::BoxArray > const & |
grids, |
|
|
T_BoundaryHandler const |
boundary_handler, |
|
|
[[maybe_unused] ] bool |
is_solver_multigrid, |
|
|
bool const |
do_single_precision_comms = false , |
|
|
std::optional< amrex::Vector< amrex::IntVect > > |
rel_ref_ratio = std::nullopt , |
|
|
[[maybe_unused] ] T_PostPhiCalculationFunctor |
post_phi_calculation = std::nullopt , |
|
|
[[maybe_unused] ] std::optional< amrex::Real const > |
current_time = std::nullopt , |
|
|
[[maybe_unused] ] std::optional< amrex::Vector< T_FArrayBoxFactory const * > > |
eb_farray_box_factory = std::nullopt |
|
) |
| |
Compute the potential phi
by solving the Poisson equation
Uses rho
as a source, assuming that the source moves at a constant speed
. This uses the AMReX solver.
More specifically, this solves the equation
- Template Parameters
-
T_BoundaryHandler | handler for boundary conditions, for example |
- See also
- ElectrostaticSolver::PoissonBoundaryHandler
- Template Parameters
-
T_PostPhiCalculationFunctor | a calculation per level directly after phi was calculated |
T_FArrayBoxFactory | usually nothing or an amrex::EBFArrayBoxFactory (EB ONLY) |
- Parameters
-
[in] | rho | The charge density a given species |
[out] | phi | The potential to be computed by this function |
[in] | beta | Represents the velocity of the source of phi |
[in] | relative_tolerance | The relative convergence threshold for the MLMG solver |
[in] | absolute_tolerance | The absolute convergence threshold for the MLMG solver |
[in] | max_iters | The maximum number of iterations allowed for the MLMG solver |
[in] | verbosity | The verbosity setting for the MLMG solver |
[in] | geom | the geometry per level (e.g., from AmrMesh) |
[in] | dmap | the distribution mapping per level (e.g., from AmrMesh) |
[in] | grids | the grids per level (e.g., from AmrMesh) |
[in] | boundary_handler | a handler for boundary conditions, for example |
- See also
- ElectrostaticSolver::PoissonBoundaryHandler
- Parameters
-
[in] | is_solver_multigrid | boolean to select the Poisson solver: 1 for Multigrid, 0 for FFT |
[in] | do_single_precision_comms | perform communications in single precision |
[in] | rel_ref_ratio | mesh refinement ratio between levels (default: 1) |
[in] | post_phi_calculation | perform a calculation per level directly after phi was calculated; required for embedded boundaries (default: none) |
[in] | current_time | the current time; required for embedded boundaries (default: none) |
[in] | eb_farray_box_factory | a factory for field data, |
- See also
- amrex::EBFArrayBoxFactory; required for embedded boundaries (default: none)
template<typename T_BoundaryHandler , typename T_PostACalculationFunctor = std::nullopt_t, typename T_FArrayBoxFactory = void>
void ablastr::fields::computeVectorPotential |
( |
amrex::Vector< amrex::Array< amrex::MultiFab *, 3 > > const & |
curr, |
|
|
amrex::Vector< amrex::Array< amrex::MultiFab *, 3 > > & |
A, |
|
|
amrex::Real const |
relative_tolerance, |
|
|
amrex::Real |
absolute_tolerance, |
|
|
int const |
max_iters, |
|
|
int const |
verbosity, |
|
|
amrex::Vector< amrex::Geometry > const & |
geom, |
|
|
amrex::Vector< amrex::DistributionMapping > const & |
dmap, |
|
|
amrex::Vector< amrex::BoxArray > const & |
grids, |
|
|
T_BoundaryHandler const |
boundary_handler, |
|
|
bool const |
do_single_precision_comms = false , |
|
|
std::optional< amrex::Vector< amrex::IntVect > > |
rel_ref_ratio = std::nullopt , |
|
|
[[maybe_unused] ] T_PostACalculationFunctor |
post_A_calculation = std::nullopt , |
|
|
[[maybe_unused] ] std::optional< amrex::Real const > |
current_time = std::nullopt , |
|
|
[[maybe_unused] ] std::optional< amrex::Vector< T_FArrayBoxFactory const * > > |
eb_farray_box_factory = std::nullopt |
|
) |
| |
Compute the vector potential A
by solving the Poisson equation
Uses J
as a source, assuming that the source moves at a constant speed
. This uses the AMReX solver.
More specifically, this solves the equation
- Template Parameters
-
T_BoundaryHandler | handler for boundary conditions, for example |
- See also
- MagnetostaticSolver::MultiPoissonBoundaryHandler
- Template Parameters
-
T_PostACalculationFunctor | a calculation per level directly after A was calculated |
T_FArrayBoxFactory | usually nothing or an amrex::EBFArrayBoxFactory (EB ONLY) |
- Parameters
-
[in] | curr | The current density a given species |
[out] | A | The vector potential to be computed by this function |
[in] | relative_tolerance | The relative convergence threshold for the MLMG solver |
[in] | absolute_tolerance | The absolute convergence threshold for the MLMG solver |
[in] | max_iters | The maximum number of iterations allowed for the MLMG solver |
[in] | verbosity | The verbosity setting for the MLMG solver |
[in] | geom | the geometry per level (e.g., from AmrMesh) |
[in] | dmap | the distribution mapping per level (e.g., from AmrMesh) |
[in] | grids | the grids per level (e.g., from AmrMesh) |
[in] | boundary_handler | a handler for boundary conditions, for example |
- See also
- MagnetostaticSolver::VectorPoissonBoundaryHandler
- Parameters
-
[in] | do_single_precision_comms | perform communications in single precision |
[in] | rel_ref_ratio | mesh refinement ratio between levels (default: 1) |
[in] | post_A_calculation | perform a calculation per level directly after A was calculated; required for embedded boundaries (default: none) |
[in] | current_time | the current time; required for embedded boundaries (default: none) |
[in] | eb_farray_box_factory | a factory for field data, |
- See also
- amrex::EBFArrayBoxFactory; required for embedded boundaries (default: none)