File Gravity.H
Typedefs
-
using GradPhiPhysBCFunct = amrex::PhysBCFunctNoOp
-
class Gravity
- #include <Gravity.H>
Public Functions
-
Gravity(amrex::Amr *Parent, int _finest_level, amrex::BCRec *_phys_bc, int _density)
Constructor
- Parameters:
Parent – Parent
Amr
object_finest_level – finest level
_phys_bc –
BCRec
object, physical boundary conditions_density – index of density component
-
~Gravity()
Deconstructor
-
Gravity(const Gravity&) = delete
Remove copy/move constructors/assignment operators.
-
Gravity(Gravity&&) = delete
-
void read_params()
Read gravity-related parameters from parameter file
-
void set_numpts_in_gravity(int numpts)
Set
Gravity
object’snumpts_at_level
variable.- Parameters:
numpts – number of points
-
void install_level(int level, amrex::AmrLevel *level_data, amrex::MultiFab &volume, amrex::MultiFab *area)
Setup gravity at level
level
.- Parameters:
level – integer, level number
level_data –
AmrLevel
object containing state data at that levelvolume – MultiFab, volume elements of cells
area – MultiFab, area elements of cells
-
void set_mass_offset(amrex::Real time, bool multi_level = true) const
Set the
mass_offset
- Parameters:
time – Current time
multi_level – boolean, do we iterate over all levels or mask off fine grids?
-
amrex::Vector<std::unique_ptr<amrex::MultiFab>> &get_grad_phi_prev(int level)
Return
grad_phi_prev
at given level- Parameters:
level – level index
-
amrex::Vector<std::unique_ptr<amrex::MultiFab>> &get_grad_phi_curr(int level)
Return
grad_phi_curr
at given level- Parameters:
level – level index
-
amrex::MultiFab *get_grad_phi_prev_comp(int level, int comp)
Return given component of
grad_phi_prev
at given level- Parameters:
level – level index
comp – component index
-
void plus_grad_phi_curr(int level, amrex::Vector<std::unique_ptr<amrex::MultiFab>> &addend)
Return
grad_phi_curr
at given level plus the given vector- Parameters:
level – level index
addend – Vector of MultiFabs to add to grad phi
-
void swapTimeLevels(int level)
Swap
grad_phi_prev
withgrad_phi_curr
at given level at set newgrad_phi_curr
to 1.e50.- Parameters:
level – level index
-
void update_max_rhs()
Calculate the maximum value of the RHS over all levels. This should only be called at a synchronization point where all Castro levels have valid new time data at the same simulation time. The RHS we will use is the density multiplied by 4*pi*G and also multiplied by the metric terms, just as it would be in a real solve.
-
void solve_for_phi(int level, amrex::MultiFab &phi, const amrex::Vector<amrex::MultiFab*> &grad_phi, int is_new)
Solve Poisson’s equation to find the gravitational potential
- Parameters:
level – level index
phi – MultiFab to store gravitational potential in
grad_phi – Vector of MultiFabs, \( \nabla \Phi \)
is_new – do we use state data at current time (1) or old time (0)?
-
void solve_for_delta_phi(int crse_level, int fine_level, const amrex::Vector<amrex::MultiFab*> &rhs, const amrex::Vector<amrex::MultiFab*> &delta_phi, const amrex::Vector<amrex::Vector<amrex::MultiFab*>> &grad_delta_phi)
Find delta phi
- Parameters:
crse_level – index of coarse level
fine_level – index of fine level
rhs – Vector of MultiFabs with right hand side source terms
delta_phi – Vector of MultiFabs delta phi will be saved to
grad_delta_phi – Vector of MultiFabs, gradient of delta phi
-
void gravity_sync(int crse_level, int fine_level, const amrex::Vector<amrex::MultiFab*> &drho, const amrex::Vector<amrex::MultiFab*> &dphi)
Sync gravity across levels
- Parameters:
crse_level – index of coarse level
fine_level – index of fine level
drho –
dphi –
-
void multilevel_solve_for_new_phi(int level, int finest_level)
Multilevel solve for new phi from base level to finest level
- Parameters:
level – Base level index
finest_level – Fine level index
-
void actual_multilevel_solve(int level, int finest_level, const amrex::Vector<amrex::Vector<amrex::MultiFab*>> &grad_phi, int is_new)
Actually do the multilevel solve for new phi from base level to finest level
- Parameters:
level – Base level index
finest_level – Fine level index
grad_phi – gradient of phi
is_new – Should we use the new state (1) or previous state (0)?
-
void create_comp_minus_level_grad_phi(int level, amrex::MultiFab &comp_phi, const amrex::Vector<amrex::MultiFab*> &comp_gphi, amrex::MultiFab &comp_minus_level_phi, amrex::Vector<std::unique_ptr<amrex::MultiFab>> &comp_minus_level_grad_phi)
Compute the difference between level and composite solves
- Parameters:
level – level index
comp_phi – MultiFab containing computed phi
comp_gphi – Vector of MultiFabs containing computed grad phi
comp_minus_level_phi – MultiFab, computed minus level phi
comp_minus_level_grad_phi – Vector of MultiFabs, computed minus level grad phi
-
void GetCrsePhi(int level, amrex::MultiFab &phi_crse, amrex::Real time)
Get coarse phi on level
level
-1- Parameters:
level – level index of fine data
phi_crse – MultiFab to contain coarse phi
time – Current time
-
void get_old_grav_vector(int level, amrex::MultiFab &grav_vector, amrex::Real time)
Get old gravity vector
- Parameters:
level – Level index
grav_vector – MultiFab containing gravity vector
time – Current time
-
void get_new_grav_vector(int level, amrex::MultiFab &grav_vector, amrex::Real time)
Get new gravity vector
- Parameters:
level – Level index
grav_vector – MultiFab containing gravity vector
time – Current time
-
void test_level_grad_phi_prev(int level)
- Parameters:
level –
-
void test_level_grad_phi_curr(int level)
- Parameters:
level –
-
void test_composite_phi(int level)
- Parameters:
level –
-
void average_fine_ec_onto_crse_ec(int level, int is_new)
- Parameters:
level – Level index
is_new – Use new state data (1) or old state data (0)?
-
void make_radial_gravity(int level, amrex::Real time, RealVector &radial_grav)
Make radial gravity
- Parameters:
level – Level index
time – Current time
radial_grav – Vector containing gravity in radial direction
-
void interpolate_monopole_grav(int level, RealVector &radial_grav, amrex::MultiFab &grav_vector) const
Interpolate monopole gravitational field in the radial direction
- Parameters:
level – Level index
radial_grav – Radial gravity
grav_vector – Gravity vector
-
void compute_radial_mass(const amrex::Box &bx, amrex::Array4<amrex::Real const> const u, RealVector &radial_mass, RealVector &radial_vol, RealVector &radial_pres, int n1d, int level) const
Integrate radially outward to find radial mass distribution
- Parameters:
bx – Box
u – Density (RHS of Poisson’s equation)
radial_mass – Radially integrated mass
radial_vol – Radially integrated volume
radial_pres – Radially integrated pressure
n1d – Number of radial points in the domain
level – Level index
-
void fill_multipole_BCs(int crse_level, int fine_level, const amrex::Vector<amrex::MultiFab*> &Rhs, amrex::MultiFab &phi)
Implement multipole boundary conditions
- Parameters:
crse_level –
fine_level –
Rhs –
phi –
-
void init_multipole_grav() const
Initialize multipole gravity
-
void fill_direct_sum_BCs(int crse_level, int fine_level, const amrex::Vector<amrex::MultiFab*> &Rhs, amrex::MultiFab &phi)
Compute and fill direct sum boundary conditions
- Parameters:
crse_level – Index of coarse level
fine_level – Index of fine level
Rhs – Vector of MultiFabs, right hand side
phi – MultiFab, phi
-
void make_mg_bc()
Make multigrid boundary conditions
-
void applyMetricTerms(int level, amrex::MultiFab &Rhs, const amrex::Vector<amrex::MultiFab*> &coeffs) const
Modify Rhs and coeffs with the appropriate metric terms.
- Parameters:
level – Level index
Rhs – MultiFab
coeffs – Vector of MultiFabs
-
void unweight_cc(int level, amrex::MultiFab &cc) const
- Parameters:
level – Index of level
cc – Cell-centered data
-
void unweight_edges(int level, const amrex::Vector<amrex::MultiFab*> &edges) const
- Parameters:
level – Index of level
edges – Edge-based data
-
void add_pointmass_to_gravity(int level, amrex::MultiFab &phi, amrex::MultiFab &grav_vector) const
- Parameters:
level – Index of level
phi – Gravitational potential
grav_vector – Gravity vector
-
amrex::Vector<std::unique_ptr<amrex::MultiFab>> get_rhs(int crse_level, int nlevs, int is_new)
Get the rhs
- Parameters:
crse_level – Index of coarse level
nlevs – Number of levels
is_new – Use new (1) or old (0) state data
-
void sanity_check(int level)
This is a sanity check on whether we are trying to fill multipole boundary conditiosn for grids at this level > 0 — this case is not currently supported. Here we shrink the domain at this level by 1 in any direction which is not symmetry or periodic, then ask if the grids at this level are contained in the shrunken domain. If not, then grids at this level touch the domain boundary and we must abort.
- Parameters:
level – Level index
-
amrex::Real actual_solve_with_mlmg(int crse_level, int fine_level, const amrex::Vector<amrex::MultiFab*> &phi, const amrex::Vector<const amrex::MultiFab*> &rhs, const amrex::Vector<std::array<amrex::MultiFab*, AMREX_SPACEDIM>> &grad_phi, const amrex::Vector<amrex::MultiFab*> &res, const amrex::MultiFab *const crse_bcdata, amrex::Real rel_eps, amrex::Real abs_eps) const
Do multigrid solve
- Parameters:
crse_level – Coarse level index
fine_level – Fine level index
phi – Gravitational potential
rhs – Right hand side
grad_phi – Grad phi
res –
crse_bcdata –
rel_eps – Relative tolerance
abs_eps – Absolute tolerance
-
amrex::Real solve_phi_with_mlmg(int crse_level, int fine_level, const amrex::Vector<amrex::MultiFab*> &phi, const amrex::Vector<amrex::MultiFab*> &rhs, const amrex::Vector<amrex::Vector<amrex::MultiFab*>> &grad_phi, const amrex::Vector<amrex::MultiFab*> &res, amrex::Real time)
Do multigrid solve to find phi
- Parameters:
crse_level – Coarse level index
fine_level – Fine level index
phi – Gravitational potential
rhs – Right hand side source term
grad_phi – Grad phi
res –
time – Current time
Public Members
-
amrex::Amr *parent
Pointers to amr,amrlevel.
-
amrex::Vector<amrex::Vector<std::unique_ptr<amrex::MultiFab>>> grad_phi_curr
Pointers to grad_phi at previous and current time
-
amrex::Real max_rhs
Maximum value of the RHS (used for obtaining absolute tolerances)
-
int Density
-
int finest_level_allocated
-
amrex::BCRec *phys_bc
-
std::array<amrex::MLLinOp::BCType, AMREX_SPACEDIM> mlmg_lobc
-
std::array<amrex::MLLinOp::BCType, AMREX_SPACEDIM> mlmg_hibc
-
int numpts_at_level
-
amrex::Vector<RealVector> radial_grav_old
-
amrex::Vector<RealVector> radial_grav_new
-
amrex::Vector<RealVector> radial_mass
-
amrex::Vector<RealVector> radial_vol
-
amrex::Vector<RealVector> radial_pres
Public Static Functions
-
static void output_job_info_params(std::ostream &jobInfoFile)
- Parameters:
jobInfoFile – std::ostream object
-
static std::string get_gravity_type()
Returns
gravity_type
-
static int get_max_solve_level()
Returns
max_solve_level
-
static int NoSync()
Returns
no_sync
-
static int DoCompositeCorrection()
Returns
do_composite_phi_correction
-
static int test_results_of_solves()
Returns
test_solves
-
static void test_residual(const amrex::Box &bx, amrex::Array4<amrex::Real> const &rhs, amrex::Array4<amrex::Real> const &ecx, amrex::Array4<amrex::Real> const &ecy, amrex::Array4<amrex::Real> const &ecz, amrex::GpuArray<amrex::Real, AMREX_SPACEDIM> dx, amrex::GpuArray<amrex::Real, AMREX_SPACEDIM> problo, int coord_type)
Test whether using the edge-based gradients to compute Div(Grad(Phi)) satisfies Lap(phi) = RHS
- Parameters:
bx – box
rhs – right-hand-side
ecx – gradients wrt x
ecy – gradients wrt y
ecz – gradients wrt z
coord_type – coordinate system
-
Gravity(amrex::Amr *Parent, int _finest_level, amrex::BCRec *_phys_bc, int _density)
-
namespace multipole
Multipole gravity data
Variables
-
const int lnum_max = 30
- AMREX_GPU_MANAGED amrex::Real volumeFactor
- AMREX_GPU_MANAGED amrex::Real parityFactor
- AMREX_GPU_MANAGED amrex::Real rmax
- AMREX_GPU_MANAGED amrex::Array1D< bool, 0, 2 > doSymmetricAddLo
- AMREX_GPU_MANAGED amrex::Array1D< bool, 0, 2 > doSymmetricAddHi
- AMREX_GPU_MANAGED bool doSymmetricAdd
- AMREX_GPU_MANAGED amrex::Array1D< bool, 0, 2 > doReflectionLo
- AMREX_GPU_MANAGED amrex::Array1D< bool, 0, 2 > doReflectionHi
- AMREX_GPU_MANAGED amrex::Array2D< amrex::Real, 0, lnum_max, 0, lnum_max > factArray
- AMREX_GPU_MANAGED amrex::Array1D< amrex::Real, 0, lnum_max > parity_q0
- AMREX_GPU_MANAGED amrex::Array2D< amrex::Real, 0, lnum_max, 0, lnum_max > parity_qC_qS
-
const int lnum_max = 30