Skip to content

Commit

Permalink
Merge branch 'development' into RemoveQKE
Browse files Browse the repository at this point in the history
  • Loading branch information
@AMLattanzi
AMLattanzi authored Oct 25, 2024
2 parents 92c1687 + cee2934 commit e37dcc7
Show file tree
Hide file tree
Showing 20 changed files with 162 additions and 95 deletions.
17 changes: 14 additions & 3 deletions Docs/sphinx_doc/Inputs.rst
View file
Original file line number Diff line number Diff line change
Expand Up @@ -23,11 +23,21 @@ Governing Equations
| | equations (instead of | | |
| | the compressible equations) | | |
+--------------------------+-----------------------------+---------------+-------------+
| **erf.use_fft** | use FFT rather than | true / false | false |
| | multigrid to solve the | | |
| | the Poisson equations | | |
+--------------------------+-----------------------------+---------------+-------------+
| **erf.mg_v** | verbosity of the multigrid | Integer >= 0 | 0 |
| | solver if used | | |
| | the Poisson equations | | |
+--------------------------+-----------------------------+---------------+-------------+

.. note::

To solve the anelastic equations, you must set ERF_USE_POISSON_SOLVE = TRUE if using
gmake or ERF_ENABLE_POISSON_SOLVE if using cmake.
To solve the anelastic equations, you must set USE_POISSON_SOLVE = TRUE if using
gmake or ERF_ENABLE_POISSON_SOLVE if using cmake. This will enable use of the
AMReX-based Poisson solver. To optionally use the FFT solver, you must additionally
set USE_FFT = TRUE if using gmake.

Problem Geometry
================
Expand Down Expand Up @@ -1221,7 +1231,8 @@ integrating the hydrostatic equation from the surface.

If **erf.init_type = custom** or **erf.init_type = input_sounding**, ``erf.nc_init_file`` and ``erf.nc_bdy_file`` do not need to be set.

Setting **erf.project_initial_velocity = 1** will have no effect if the code is not built with **ERF_USE_POISSON_SOLVE** defined.
Setting **erf.project_initial_velocity = 1** will have no effect if the code is not built with **ERF_USE_POISSON_SOLVE** defined
if using cmake or with **USE_POISSON_SOLVE = TRUE** if using gmake.

Map Scale Factors
=================
Expand Down
2 changes: 1 addition & 1 deletion Source/ERF.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -348,7 +348,7 @@ public:
const std::string &mfPrefix) const;

void erf_enforce_hse (int lev,
amrex::MultiFab& dens, amrex::MultiFab& pres, amrex::MultiFab& pi,
amrex::MultiFab& dens, amrex::MultiFab& pres, amrex::MultiFab& pi, amrex::MultiFab& th,
std::unique_ptr<amrex::MultiFab>& z_cc);

#ifdef ERF_USE_NETCDF
Expand Down
10 changes: 5 additions & 5 deletions Source/ERF.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -552,7 +552,7 @@ ERF::post_timestep (int nstep, Real time, Real dt_lev0)
// Copy z_phs_nd and detJ_cc at end of timestep
MultiFab::Copy(*z_phys_nd[lev], *z_phys_nd_new[lev], 0, 0, 1, z_phys_nd[lev]->nGrowVect());
MultiFab::Copy( *detJ_cc[lev], *detJ_cc_new[lev], 0, 0, 1, detJ_cc[lev]->nGrowVect());
MultiFab::Copy(base_state[lev],base_state_new[lev],0,0,3,1);
MultiFab::Copy(base_state[lev],base_state_new[lev],0,0,BaseState::num_comps,base_state[lev].nGrowVect());

make_zcc(geom[lev],*z_phys_nd[lev],*z_phys_cc[lev]);
}
Expand Down Expand Up @@ -925,8 +925,8 @@ ERF::InitData_post ()
base_state[lev].FillBoundary(geom[lev].periodicity());

// For moving terrain only
if (solverChoice.terrain_type != TerrainType::Static) {
MultiFab::Copy(base_state_new[lev],base_state[lev],0,0,3,1);
if (solverChoice.terrain_type == TerrainType::Moving) {
MultiFab::Copy(base_state_new[lev],base_state[lev],0,0,BaseState::num_comps,base_state[lev].nGrowVect());
base_state_new[lev].FillBoundary(geom[lev].periodicity());
}

Expand Down Expand Up @@ -1662,7 +1662,7 @@ ERF::MakeHorizontalAverages ()
fab_arr(i, j, k, 1) = cons_arr(i, j, k, RhoTheta_comp) / dens;
if (!use_moisture) {
if (is_anelastic) {
fab_arr(i,j,k,2) = hse_arr(i,j,k,1);
fab_arr(i,j,k,2) = hse_arr(i,j,k,BaseState::p0_comp);
} else {
fab_arr(i,j,k,2) = getPgivenRTh(cons_arr(i,j,k,RhoTheta_comp));
}
Expand All @@ -1682,7 +1682,7 @@ ERF::MakeHorizontalAverages ()
ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) {
Real dens = cons_arr(i, j, k, Rho_comp);
if (is_anelastic) {
fab_arr(i,j,k,2) = hse_arr(i,j,k,1);
fab_arr(i,j,k,2) = hse_arr(i,j,k,BaseState::p0_comp);
} else {
Real qv = cons_arr(i, j, k, RhoQ1_comp) / dens;
fab_arr(i, j, k, 2) = getPgivenRTh(cons_arr(i, j, k, RhoTheta_comp), qv);
Expand Down
11 changes: 11 additions & 0 deletions Source/ERF_IndexDefines.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -57,6 +57,17 @@
#define PrimQ5_comp (RhoQ5_comp-1)
#define PrimQ6_comp (RhoQ6_comp-1)

// Base state variables
namespace BaseState {
enum {
r0_comp = 0,
p0_comp,
pi0_comp,
th0_comp,
num_comps
};
}

// NOTE: we still use this indexing even if no moisture
// NOTE: We assume a single boundary condition for all passive scalars
// NOTE: and a single boundary condition for all components of the base state
Expand Down
10 changes: 5 additions & 5 deletions Source/ERF_make_new_arrays.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -28,11 +28,11 @@ ERF::init_stuff (int lev, const BoxArray& ba, const DistributionMapping& dm,
// ********************************************************************************************
// Base state holds r_0, pres_0, pi_0 (in that order)
// ********************************************************************************************
tmp_base_state.define(ba,dm,3,1);
tmp_base_state.define(ba,dm,BaseState::num_comps,3);
tmp_base_state.setVal(0.);

if (solverChoice.use_terrain && solverChoice.terrain_type != TerrainType::Static) {
base_state_new[lev].define(ba,dm,3,1);
if (solverChoice.use_terrain && solverChoice.terrain_type == TerrainType::Moving) {
base_state_new[lev].define(ba,dm,BaseState::num_comps,3);
base_state_new[lev].setVal(0.);
}

Expand All @@ -42,7 +42,7 @@ ERF::init_stuff (int lev, const BoxArray& ba, const DistributionMapping& dm,
if (solverChoice.use_terrain) {
z_phys_cc[lev] = std::make_unique<MultiFab>(ba,dm,1,1);

if (solverChoice.terrain_type != TerrainType::Static)
if (solverChoice.terrain_type == TerrainType::Moving)
{
detJ_cc_new[lev] = std::make_unique<MultiFab>(ba,dm,1,1);
detJ_cc_src[lev] = std::make_unique<MultiFab>(ba,dm,1,1);
Expand All @@ -65,7 +65,7 @@ ERF::init_stuff (int lev, const BoxArray& ba, const DistributionMapping& dm,
int ngrow = ComputeGhostCells(solverChoice.advChoice, solverChoice.use_NumDiff) + 2;
tmp_zphys_nd = std::make_unique<MultiFab>(ba_nd,dm,1,IntVect(ngrow,ngrow,ngrow));

if (solverChoice.terrain_type != TerrainType::Static) {
if (solverChoice.terrain_type == TerrainType::Moving) {
z_phys_nd_new[lev] = std::make_unique<MultiFab>(ba_nd,dm,1,IntVect(ngrow,ngrow,ngrow));
z_phys_nd_src[lev] = std::make_unique<MultiFab>(ba_nd,dm,1,IntVect(ngrow,ngrow,ngrow));
}
Expand Down
2 changes: 1 addition & 1 deletion Source/ERF_make_new_level.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -232,7 +232,7 @@ ERF::MakeNewLevelFromCoarse (int lev, Real time, const BoxArray& ba,
//
InterpFromCoarseLevel(base_state[lev], base_state[lev].nGrowVect(),
IntVect(0,0,0), // do not fill ghost cells outside the domain
base_state[lev-1], 0, 0, 3,
base_state[lev-1], 0, 0, base_state[lev].nComp(),
geom[lev-1], geom[lev],
refRatio(lev-1), &cell_cons_interp,
domain_bcs_type, BCVars::cons_bc);
Expand Down
10 changes: 4 additions & 6 deletions Source/IO/ERF_Plotfile.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -348,17 +348,15 @@ ERF::WritePlotFile (int which, Vector<std::string> plot_var_names)
mf_comp += 1;
}

MultiFab r_hse(base_state[lev], make_alias, 0, 1); // r_0 is first component
MultiFab p_hse(base_state[lev], make_alias, 1, 1); // p_0 is second component
MultiFab r_hse(base_state[lev], make_alias, BaseState::r0_comp, 1);
MultiFab p_hse(base_state[lev], make_alias, BaseState::p0_comp, 1);
if (containerHasElement(plot_var_names, "pres_hse"))
{
// p_0 is second component of base_state
MultiFab::Copy(mf[lev],p_hse,0,mf_comp,1,0);
mf_comp += 1;
}
if (containerHasElement(plot_var_names, "dens_hse"))
{
// r_0 is first component of base_state
MultiFab::Copy(mf[lev],r_hse,0,mf_comp,1,0);
mf_comp += 1;
}
Expand Down Expand Up @@ -530,7 +528,7 @@ ERF::WritePlotFile (int which, Vector<std::string> plot_var_names)
const Array4<Real const>& S_arr = vars_new[lev][Vars::cons].const_array(mfi);
if (solverChoice.anelastic[lev] == 1) {
ParallelFor(gbx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
p_arr(i,j,k) = hse_arr(i,j,k,1);
p_arr(i,j,k) = hse_arr(i,j,k,BaseState::p0_comp);
});
} else {
ParallelFor(gbx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
Expand Down Expand Up @@ -624,7 +622,7 @@ ERF::WritePlotFile (int which, Vector<std::string> plot_var_names)
const Array4<Real const>& S_arr = vars_new[lev][Vars::cons].const_array(mfi);
if (solverChoice.anelastic[lev] == 1) {
ParallelFor(gbx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
p_arr(i,j,k) = hse_arr(i,j,k,1);
p_arr(i,j,k) = hse_arr(i,j,k,BaseState::p0_comp);
});
} else {
ParallelFor(gbx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept {
Expand Down
2 changes: 1 addition & 1 deletion Source/IO/ERF_Write1DProfiles.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -263,7 +263,7 @@ ERF::derive_diag_profiles(Real /*time*/,
int nvars = vars_new[lev][Vars::cons].nComp();
MultiFab mf_cons(vars_new[lev][Vars::cons], make_alias, 0, nvars);

MultiFab p_hse (base_state[lev], make_alias, 1, 1); // p_0 is second component
MultiFab p_hse (base_state[lev], make_alias, BaseState::p0_comp, 1);

bool use_moisture = (solverChoice.moisture_type != MoistureType::None);

Expand Down
2 changes: 1 addition & 1 deletion Source/IO/ERF_Write1DProfiles_stag.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -345,7 +345,7 @@ ERF::derive_diag_profiles_stag (Real /*time*/,
int nvars = vars_new[lev][Vars::cons].nComp();
MultiFab mf_cons(vars_new[lev][Vars::cons], make_alias, 0, nvars);

MultiFab p_hse (base_state[lev], make_alias, 1, 1); // p_0 is second component
MultiFab p_hse (base_state[lev], make_alias, BaseState::p0_comp, 1);

bool use_moisture = (solverChoice.moisture_type != MoistureType::None);

Expand Down
2 changes: 1 addition & 1 deletion Source/IO/ERF_WriteScalarProfiles.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -39,7 +39,7 @@ ERF::sum_integrated_quantities (Real time)
MultiFab pert_dens(vars_new[lev][Vars::cons].boxArray(),
vars_new[lev][Vars::cons].DistributionMap(),
1,0);
MultiFab r_hse (base_state[lev], make_alias, 0, 1); // r_0 is first component
MultiFab r_hse (base_state[lev], make_alias, BaseState::r0_comp, 1);
for ( MFIter mfi(pert_dens,TilingIfNotGPU()); mfi.isValid(); ++mfi)
{
const Box& bx = mfi.tilebox();
Expand Down
1 change: 1 addition & 0 deletions Source/Initialization/ERF_Metgrid_utils.H
View file
Original file line number Diff line number Diff line change
Expand Up @@ -94,6 +94,7 @@ init_base_state_from_metgrid (const bool use_moisture,
amrex::FArrayBox& r_hse_fab,
amrex::FArrayBox& p_hse_fab,
amrex::FArrayBox& pi_hse_fab,
amrex::FArrayBox& th_hse_fab,
amrex::FArrayBox& z_phys_cc_fab,
const amrex::Vector<amrex::FArrayBox>& NC_ght_fab,
const amrex::Vector<amrex::FArrayBox>& NC_psfc_fab,
Expand Down
Loading

0 comments on commit e37dcc7

Please sign in to comment.