Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Terrain with Rayleigh damping #1399

Merged
merged 6 commits into from
Jan 25, 2024
Merged

Terrain with Rayleigh damping #1399

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
6 commits
Select commit Hold shift + click to select a range
669992e
Pass terrain heights to erf_init_rayleigh
ewquon Jan 25, 2024
333d5fd
Add terrain support to default init_rayleigh_damping function
ewquon Jan 25, 2024
6e38097
Fix indexing, parallel behavior
ewquon Jan 25, 2024
cf03e92
Enable Rayleigh damping for WoA
ewquon Jan 25, 2024
193dbd9
Update WoA input to demo Rayleigh damping with terrain
ewquon Jan 25, 2024
e177c97
Merge branch 'development' of github.com:erf-model/ERF into HEAD
ewquon Jan 25, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
3 changes: 2 additions & 1 deletion Exec/DevTests/MovingTerrain/prob.H
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -49,7 +49,8 @@ public:
amrex::Vector<amrex::Real>& vbar,
amrex::Vector<amrex::Real>& wbar,
amrex::Vector<amrex::Real>& thetabar,
amrex::Geometry const& geom) override;
amrex::Geometry const& geom,
std::unique_ptr<amrex::MultiFab>& z_phys_cc) override;

protected:
std::string name() override { return "MovingTerrain"; }
Expand Down
3 changes: 2 additions & 1 deletion Exec/DevTests/MovingTerrain/prob.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -137,7 +137,8 @@ Problem::erf_init_rayleigh(
amrex::Vector<Real>& vbar,
amrex::Vector<Real>& wbar,
amrex::Vector<Real>& thetabar,
amrex::Geometry const& geom)
amrex::Geometry const& geom,
std::unique_ptr<MultiFab>& /*z_phys_cc*/)
{
// amrex::Error("Should never get here for MovingTerrain problem");
const int khi = geom.Domain().bigEnd()[2];
Expand Down
3 changes: 2 additions & 1 deletion Exec/Radiation/prob.H
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -52,7 +52,8 @@ public:
amrex::Vector<amrex::Real>& vbar,
amrex::Vector<amrex::Real>& wbar,
amrex::Vector<amrex::Real>& thetabar,
amrex::Geometry const& geom) override;
amrex::Geometry const& geom,
std::unique_ptr<amrex::MultiFab>& z_phys_cc) override;

protected:
std::string name() override { return "Supercell"; }
Expand Down
3 changes: 2 additions & 1 deletion Exec/Radiation/prob.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -219,7 +219,8 @@ Problem::erf_init_rayleigh(
amrex::Vector<amrex::Real>& vbar,
amrex::Vector<amrex::Real>& wbar,
amrex::Vector<amrex::Real>& thetabar,
amrex::Geometry const& geom)
amrex::Geometry const& geom,
std::unique_ptr<amrex::MultiFab>& /*z_phys_cc*/)
{
const int khi = geom.Domain().bigEnd()[2];

Expand Down
3 changes: 2 additions & 1 deletion Exec/RegTests/ScalarAdvDiff/prob.H
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -58,7 +58,8 @@ public:
amrex::Vector<amrex::Real>& vbar,
amrex::Vector<amrex::Real>& wbar,
amrex::Vector<amrex::Real>& thetabar,
amrex::Geometry const& geom) override;
amrex::Geometry const& geom,
std::unique_ptr<amrex::MultiFab>& z_phys_cc) override;

protected:
std::string name() override { return "Scalar Advection/Diffusion"; }
Expand Down
3 changes: 2 additions & 1 deletion Exec/RegTests/ScalarAdvDiff/prob.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -41,7 +41,8 @@ Problem::erf_init_rayleigh(
amrex::Vector<Real>& vbar,
amrex::Vector<Real>& wbar,
amrex::Vector<Real>& thetabar,
amrex::Geometry const& geom)
amrex::Geometry const& geom,
std::unique_ptr<MultiFab>& z_phys_cc)
{
const int khi = geom.Domain().bigEnd()[2];

Expand Down
3 changes: 2 additions & 1 deletion Exec/RegTests/ScalarAdvDiff/prob.cpp.convergence
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -16,7 +16,8 @@ erf_init_rayleigh(amrex::Vector<Real>& tau,
amrex::Vector<Real>& vbar,
amrex::Vector<Real>& wbar,
amrex::Vector<Real>& thetabar,
amrex::Geometry const& geom)
amrex::Geometry const& geom,
std::unique_ptr<MultiFab>& z_phys_cc)
{
const int khi = geom.Domain().bigEnd()[2];

Expand Down
9 changes: 6 additions & 3 deletions Exec/RegTests/WitchOfAgnesi/inputs
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -41,13 +41,12 @@ erf.check_int = -57600 # number of timesteps between checkpoints

# PLOTFILES
erf.plot_file_1 = plt # prefix of plotfile name
erf.plot_int_1 = 1 # number of timesteps between plotfiles
erf.plot_int_1 = 10 # number of timesteps between plotfiles
erf.plot_vars_1 = density x_velocity y_velocity z_velocity pressure theta pres_hse dens_hse pert_pres pert_dens z_phys detJ dpdx dpdy pres_hse_x pres_hse_y

# SOLVER CHOICE
erf.use_gravity = true
erf.use_coriolis = false
erf.use_rayleigh_damping = false
erf.les_type = "None"

# TERRRAIN GRID TYPE
Expand All @@ -62,5 +61,9 @@ erf.alpha_T = 0.0 # [m^2/s]

# PROBLEM PARAMETERS (optional)
prob.T_0 = 300.0
prob.U_0 = 0.0
prob.U_0 = 1.0
prob.rho_0 = 1.16

erf.use_rayleigh_damping = true
prob.dampcoef = 2000. # ==> time scale = 0.0005 s (for testing)
prob.zdamp = 1.0
5 changes: 5 additions & 0 deletions Exec/RegTests/WitchOfAgnesi/prob.H
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -9,6 +9,10 @@

struct ProbParm : ProbParmDefaults {
amrex::Real U_0 = 0.0;
amrex::Real V_0 = 0.0;
amrex::Real W_0 = 0.0;
amrex::Real dampcoef = 0.003;
amrex::Real zdamp = 1.0;
}; // namespace ProbParm

class Problem : public ProblemBase
Expand All @@ -17,6 +21,7 @@ public:
Problem();

#include "Prob/init_density_hse_dry_terrain.H"
#include "Prob/init_rayleigh_damping.H"

void init_custom_pert (
const amrex::Box& bx,
Expand Down
4 changes: 4 additions & 0 deletions Exec/RegTests/WitchOfAgnesi/prob.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -16,6 +16,10 @@ Problem::Problem()
// Parse params
amrex::ParmParse pp("prob");
pp.query("U_0", parms.U_0);
pp.query("V_0", parms.V_0);
pp.query("W_0", parms.W_0);
pp.query("dampcoef", parms.dampcoef);
pp.query("zdamp", parms.zdamp);

init_base_parms(parms.rho_0, parms.T_0);
}
Expand Down
3 changes: 2 additions & 1 deletion Exec/SquallLine_2D/prob.H
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -67,7 +67,8 @@ public:
amrex::Vector<amrex::Real>& vbar,
amrex::Vector<amrex::Real>& wbar,
amrex::Vector<amrex::Real>& thetabar,
amrex::Geometry const& geom) override;
amrex::Geometry const& geom,
std::unique_ptr<amrex::MultiFab>& z_phys_cc) override;

amrex::Real compute_theta (amrex::Real z);

Expand Down
3 changes: 2 additions & 1 deletion Exec/SquallLine_2D/prob.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -407,7 +407,8 @@ Problem::erf_init_rayleigh (amrex::Vector<amrex::Real>& tau,
amrex::Vector<amrex::Real>& vbar,
amrex::Vector<amrex::Real>& wbar,
amrex::Vector<amrex::Real>& thetabar,
amrex::Geometry const& geom)
amrex::Geometry const& geom,
std::unique_ptr<amrex::MultiFab>& /*z_phys_cc*/)
{
const int khi = geom.Domain().bigEnd()[2];

Expand Down
3 changes: 2 additions & 1 deletion Exec/SuperCell/prob.H
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -52,7 +52,8 @@ public:
amrex::Vector<amrex::Real>& vbar,
amrex::Vector<amrex::Real>& wbar,
amrex::Vector<amrex::Real>& thetabar,
amrex::Geometry const& geom) override;
amrex::Geometry const& geom,
std::unique_ptr<amrex::MultiFab>& z_phys_cc) override;

protected:
std::string name() override { return "Supercell"; }
Expand Down
3 changes: 2 additions & 1 deletion Exec/SuperCell/prob.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -219,7 +219,8 @@ Problem::erf_init_rayleigh(
amrex::Vector<amrex::Real>& vbar,
amrex::Vector<amrex::Real>& wbar,
amrex::Vector<amrex::Real>& thetabar,
amrex::Geometry const& geom)
amrex::Geometry const& geom,
std::unique_ptr<MultiFab>& /*z_phys_cc*/)
{
const int khi = geom.Domain().bigEnd()[2];

Expand Down
3 changes: 2 additions & 1 deletion Source/Initialization/ERF_init1d.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -44,7 +44,8 @@ ERF::initRayleigh ()
d_rayleigh_thetabar[lev].resize(zlen_rayleigh, 0.0_rt);

prob->erf_init_rayleigh(h_rayleigh_tau[lev], h_rayleigh_ubar[lev], h_rayleigh_vbar[lev],
h_rayleigh_wbar[lev], h_rayleigh_thetabar[lev], geom[lev]);
h_rayleigh_wbar[lev], h_rayleigh_thetabar[lev], geom[lev],
z_phys_cc[lev]);

// Copy from host version to device version
amrex::Gpu::copy(amrex::Gpu::hostToDevice, h_rayleigh_tau[lev].begin(), h_rayleigh_tau[lev].end(),
Expand Down
42 changes: 34 additions & 8 deletions Source/Prob/init_rayleigh_damping.H
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -8,21 +8,46 @@ erf_init_rayleigh (amrex::Vector<amrex::Real>& tau,
amrex::Vector<amrex::Real>& vbar,
amrex::Vector<amrex::Real>& wbar,
amrex::Vector<amrex::Real>& thetabar,
amrex::Geometry const& geom) override
amrex::Geometry const& geom,
std::unique_ptr<amrex::MultiFab>& z_phys_cc) override
{
const auto ztop = geom.ProbHi()[2];
const int khi = geom.Domain().bigEnd()[2];
amrex::Vector<amrex::Real> zcc(khi+1);

if (z_phys_cc) {
// use_terrain=1
// calculate the damping strength based on the max height at each k
amrex::MultiArray4<const amrex::Real> const& ma_z_phys = z_phys_cc->const_arrays();
for (int k = 0; k <= khi; k++) {
zcc[k] = amrex::ParReduce(amrex::TypeList<amrex::ReduceOpMax>{},
amrex::TypeList<amrex::Real>{},
*z_phys_cc,
amrex::IntVect(0),
[=] AMREX_GPU_DEVICE (int box_no, int i, int j, int) noexcept
-> amrex::GpuTuple<amrex::Real>
{
return { ma_z_phys[box_no](i,j,k) };
});
}
amrex::ParallelDescriptor::ReduceRealMax(zcc.data(), zcc.size());

} else {
const amrex::Real* prob_lo = geom.ProbLo();
const auto dx = geom.CellSize();
for (int k = 0; k <= khi; k++)
{
zcc[k] = prob_lo[2] + (k+0.5) * dx[2];
}
}

const auto ztop = zcc[khi];
amrex::Print() << "Rayleigh damping (coef="<<parms.dampcoef<<" s^-1)"
<< " between " << ztop-parms.zdamp << " and " << ztop << " m"
<< std::endl;

const amrex::Real* prob_lo = geom.ProbLo();
const auto dx = geom.CellSize();
const int khi = geom.Domain().bigEnd()[2];

for (int k = 0; k <= khi; k++)
{
const amrex::Real z = prob_lo[2] + (k + 0.5) * dx[2];
const amrex::Real zfrac = 1 - (ztop - z) / parms.zdamp;
const amrex::Real zfrac = 1 - (ztop - zcc[k]) / parms.zdamp;
if (zfrac >= 0)
{
const amrex::Real sinefac = std::sin(PIoTwo*zfrac);
Expand All @@ -31,6 +56,7 @@ erf_init_rayleigh (amrex::Vector<amrex::Real>& tau,
vbar[k] = parms.V_0;
wbar[k] = parms.W_0;
thetabar[k] = parms.T_0;
amrex::Print() << " " << zcc[k] << " " << sinefac*sinefac << std::endl;
}
else
{
Expand Down
14 changes: 8 additions & 6 deletions Source/prob_common.H
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -150,20 +150,22 @@ public:
/**
* Function to define the quantities needed to impose Rayleigh damping
*
* @param[in] tau strength of Rayleigh damping
* @param[in] ubar reference value for x-velocity used to define Rayleigh damping
* @param[in] vbar reference value for y-velocity used to define Rayleigh damping
* @param[in] wbar reference value for z-velocity used to define Rayleigh damping
* @param[in] thetabar reference value for potential temperature used to define Rayleigh damping
* @param[out] tau strength of Rayleigh damping
* @param[out] ubar reference value for x-velocity used to define Rayleigh damping
* @param[out] vbar reference value for y-velocity used to define Rayleigh damping
* @param[out] wbar reference value for z-velocity used to define Rayleigh damping
* @param[out] thetabar reference value for potential temperature used to define Rayleigh damping
* @param[in] geom container for geometric information
* @param[in] z_phys_cc height coordinate at cell centers
*/
virtual void
erf_init_rayleigh (amrex::Vector<amrex::Real>& /*tau*/,
amrex::Vector<amrex::Real>& /*ubar*/,
amrex::Vector<amrex::Real>& /*vbar*/,
amrex::Vector<amrex::Real>& /*wbar*/,
amrex::Vector<amrex::Real>& /*thetabar*/,
amrex::Geometry const& /*geom*/)
amrex::Geometry const& /*geom*/,
std::unique_ptr<amrex::MultiFab>& /*z_phys_cc*/)
{
amrex::Error("Should never call erf_init_rayleigh for "+name()+" problem");
}
Expand Down
Loading