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.

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

Precip bubble #1634

Merged
merged 2 commits into from
May 29, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
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
81 changes: 81 additions & 0 deletions Exec/RegTests/Bubble/inputs_BF02_moist_bubble_Kessler
Original file line number Diff line number Diff line change
@@ -0,0 +1,81 @@
# ------------------ INPUTS TO MAIN PROGRAM -------------------
max_step = 2000
stop_time = 3600.0

amrex.fpe_trap_invalid = 1

fabarray.mfiter_tile_size = 1024 1024 1024

# PROBLEM SIZE & GEOMETRY
geometry.prob_extent = 20000.0 400.0 10000.0
amr.n_cell = 200 4 100
geometry.is_periodic = 0 1 0
xlo.type = "SlipWall"
xhi.type = "SlipWall"
zlo.type = "SlipWall"
zhi.type = "SlipWall"

# TIME STEP CONTROL
erf.fixed_dt = 0.5
erf.fixed_mri_dt_ratio = 4
#erf.no_substepping = 1
#erf.fixed_dt = 0.1

# DIAGNOSTICS & VERBOSITY
erf.sum_interval = 1 # timesteps between computing mass
erf.v = 1 # verbosity in ERF.cpp
amr.v = 1 # verbosity in Amr.cpp

# REFINEMENT / REGRIDDING
amr.max_level = 0 # maximum level number allowed

# CHECKPOINT FILES
erf.check_file = chk # root name of checkpoint file
erf.check_int = 100 # number of timesteps between checkpoints

# PLOTFILES
erf.plot_file_1 = plt # prefix of plotfile name
erf.plot_int_1 = 100 # number of timesteps between plotfiles
erf.plot_vars_1 = density rhotheta rhoQ1 rhoQ2 rhoadv_0 x_velocity y_velocity z_velocity pressure theta scalar temp pres_hse dens_hse pert_pres pert_dens eq_pot_temp qt qv qc qrain

# SOLVER CHOICES
erf.use_gravity = true
erf.use_coriolis = false

erf.dycore_horiz_adv_type = "Upwind_3rd"
erf.dycore_vert_adv_type = "Upwind_3rd"
erf.dryscal_horiz_adv_type = "Upwind_3rd"
erf.dryscal_vert_adv_type = "Upwind_3rd"
erf.moistscal_horiz_adv_type = "Upwind_3rd"
erf.moistscal_vert_adv_type = "Upwind_3rd"

# PHYSICS OPTIONS
erf.les_type = "None"
erf.pbl_type = "None"
erf.moisture_model = "Kessler"
erf.buoyancy_type = 1
erf.use_moist_background = true

erf.molec_diff_type = "ConstantAlpha"
erf.rho0_trans = 1.0 # [kg/m^3], used to convert input diffusivities
erf.dynamicViscosity = 0.0 # [kg/(m-s)] ==> nu = 75.0 m^2/s
erf.alpha_T = 0.0 # [m^2/s]
erf.alpha_C = 0.0

# INITIAL CONDITIONS
#erf.init_type = "input_sounding"
#erf.input_sounding_file = "BF02_moist_sounding"
#erf.init_sounding_ideal = true

# PROBLEM PARAMETERS (optional)
# warm bubble input
prob.x_c = 10000.0
prob.z_c = 2000.0
prob.x_r = 2000.0
prob.z_r = 2000.0
prob.T_0 = 300.0

prob.do_moist_bubble = true
prob.theta_pert = 2.0
prob.qt_init = 0.02
prob.eq_pot_temp = 320.0
81 changes: 81 additions & 0 deletions Exec/RegTests/Bubble/inputs_BF02_moist_bubble_SAM_NoIce
Original file line number Diff line number Diff line change
@@ -0,0 +1,81 @@
# ------------------ INPUTS TO MAIN PROGRAM -------------------
max_step = 2000
stop_time = 3600.0

amrex.fpe_trap_invalid = 1

fabarray.mfiter_tile_size = 1024 1024 1024

# PROBLEM SIZE & GEOMETRY
geometry.prob_extent = 20000.0 400.0 10000.0
amr.n_cell = 200 4 100
geometry.is_periodic = 0 1 0
xlo.type = "SlipWall"
xhi.type = "SlipWall"
zlo.type = "SlipWall"
zhi.type = "SlipWall"

# TIME STEP CONTROL
erf.fixed_dt = 0.5
erf.fixed_mri_dt_ratio = 4
#erf.no_substepping = 1
#erf.fixed_dt = 0.1

# DIAGNOSTICS & VERBOSITY
erf.sum_interval = 1 # timesteps between computing mass
erf.v = 1 # verbosity in ERF.cpp
amr.v = 1 # verbosity in Amr.cpp

# REFINEMENT / REGRIDDING
amr.max_level = 0 # maximum level number allowed

# CHECKPOINT FILES
erf.check_file = chk # root name of checkpoint file
erf.check_int = 100 # number of timesteps between checkpoints

# PLOTFILES
erf.plot_file_1 = plt # prefix of plotfile name
erf.plot_int_1 = 100 # number of timesteps between plotfiles
erf.plot_vars_1 = density rhotheta rhoQ1 rhoQ2 rhoadv_0 x_velocity y_velocity z_velocity pressure theta scalar temp pres_hse dens_hse pert_pres pert_dens eq_pot_temp qt qv qc qi qp qrain qsnow qgraup

# SOLVER CHOICES
erf.use_gravity = true
erf.use_coriolis = false

erf.dycore_horiz_adv_type = "Upwind_3rd"
erf.dycore_vert_adv_type = "Upwind_3rd"
erf.dryscal_horiz_adv_type = "Upwind_3rd"
erf.dryscal_vert_adv_type = "Upwind_3rd"
erf.moistscal_horiz_adv_type = "Upwind_3rd"
erf.moistscal_vert_adv_type = "Upwind_3rd"

# PHYSICS OPTIONS
erf.les_type = "None"
erf.pbl_type = "None"
erf.moisture_model = "SAM_NoIce"
erf.buoyancy_type = 1
erf.use_moist_background = true

erf.molec_diff_type = "ConstantAlpha"
erf.rho0_trans = 1.0 # [kg/m^3], used to convert input diffusivities
erf.dynamicViscosity = 0.0 # [kg/(m-s)] ==> nu = 75.0 m^2/s
erf.alpha_T = 0.0 # [m^2/s]
erf.alpha_C = 0.0

# INITIAL CONDITIONS
#erf.init_type = "input_sounding"
#erf.input_sounding_file = "BF02_moist_sounding"
#erf.init_sounding_ideal = true

# PROBLEM PARAMETERS (optional)
# warm bubble input
prob.x_c = 10000.0
prob.z_c = 2000.0
prob.x_r = 2000.0
prob.z_r = 2000.0
prob.T_0 = 300.0

prob.do_moist_bubble = true
prob.theta_pert = 2.0
prob.qt_init = 0.02
prob.eq_pot_temp = 320.0
3 changes: 2 additions & 1 deletion Exec/RegTests/Bubble/prob.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -364,7 +364,8 @@ Problem::init_custom_pert(

if(sc.moisture_type == MoistureType::SAM) {
moisture_type = 1;
} else if(sc.moisture_type == MoistureType::SAM_NoPrecip_NoIce) {
} else if(sc.moisture_type == MoistureType::SAM_NoIce ||
sc.moisture_type == MoistureType::SAM_NoPrecip_NoIce) {
moisture_type = 2;
}

Expand Down
5 changes: 4 additions & 1 deletion Source/DataStructs/DataStruct.H
Original file line number Diff line number Diff line change
Expand Up @@ -33,7 +33,7 @@ enum struct MoistureModelType{
};

enum struct MoistureType {
Kessler, SAM, SAM_NoPrecip_NoIce, Kessler_NoRain, None
Kessler, SAM, SAM_NoIce, SAM_NoPrecip_NoIce, Kessler_NoRain, None
};

enum struct WindFarmType {
Expand Down Expand Up @@ -95,6 +95,8 @@ struct SolverChoice {
pp.query("moisture_model", moisture_model_string);
if (moisture_model_string == "SAM") {
moisture_type = MoistureType::SAM;
} else if (moisture_model_string == "SAM_NoIce") {
moisture_type = MoistureType::SAM_NoIce;
} else if (moisture_model_string == "SAM_NoPrecip_NoIce") {
moisture_type = MoistureType::SAM_NoPrecip_NoIce;
} else if (moisture_model_string == "Kessler") {
Expand All @@ -110,6 +112,7 @@ struct SolverChoice {
if (moisture_type != MoistureType::None) {
if (moisture_type == MoistureType::Kessler_NoRain ||
moisture_type == MoistureType::SAM ||
moisture_type == MoistureType::SAM_NoIce ||
moisture_type == MoistureType::SAM_NoPrecip_NoIce) {
buoyancy_type = 1; // asserted in make buoyancy
} else {
Expand Down
2 changes: 1 addition & 1 deletion Source/ERF.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -1451,7 +1451,7 @@ ERF::ReadParameters ()
lsm.SetModel<NullSurf>();
Print() << "Null land surface model!\n";
} else {
Abort("Dont know this moisture_type!") ;
Abort("Dont know this LandSurfaceType!") ;
}

if (verbose > 0) {
Expand Down
12 changes: 7 additions & 5 deletions Source/IO/Plotfile.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -77,11 +77,13 @@ ERF::setPlotVariables (const std::string& pp_plot_var_names, Vector<std::string>
for (int i = 0; i < derived_names.size(); ++i) {
if ( containerHasElement(plot_var_names, derived_names[i]) ) {
if (solverChoice.use_terrain || (derived_names[i] != "z_phys" && derived_names[i] != "detJ") ) {
if ( (solverChoice.moisture_type == MoistureType::SAM) || (derived_names[i] != "qi" &&
derived_names[i] != "qsnow" &&
derived_names[i] != "qgraup" &&
derived_names[i] != "snow_accum" &&
derived_names[i] != "graup_accum") )
if ( (solverChoice.moisture_type == MoistureType::SAM ||
solverChoice.moisture_type == MoistureType::SAM_NoIce) ||
(derived_names[i] != "qi" &&
derived_names[i] != "qsnow" &&
derived_names[i] != "qgraup" &&
derived_names[i] != "snow_accum" &&
derived_names[i] != "graup_accum") )
{
tmp_plot_names.push_back(derived_names[i]);
} // moisture_type
Expand Down
5 changes: 3 additions & 2 deletions Source/Microphysics/EulerianMicrophysics.H
Original file line number Diff line number Diff line change
Expand Up @@ -27,11 +27,12 @@ public:
{
AMREX_ASSERT( Microphysics::modelType(a_model_type) == MoistureModelType::Eulerian );
m_moist_model.resize(nlev);
if (a_model_type == MoistureType::SAM or
if (a_model_type == MoistureType::SAM ||
a_model_type == MoistureType::SAM_NoIce ||
a_model_type == MoistureType::SAM_NoPrecip_NoIce) {
SetModel<SAM>();
amrex::Print() << "SAM moisture model!\n";
} else if (a_model_type == MoistureType::Kessler or
} else if (a_model_type == MoistureType::Kessler ||
a_model_type == MoistureType::Kessler_NoRain) {
SetModel<Kessler>();
amrex::Print() << "Kessler moisture model!\n";
Expand Down
4 changes: 2 additions & 2 deletions Source/Microphysics/Kessler/Kessler.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -144,12 +144,12 @@ void Kessler::AdvanceKessler (const SolverChoice &solverChoice)

autor = 0.0;
if (qcc > qcw0) {
autor = 0.001;
autor = alphaelq;
}

accrr = 0.0;
accrr = 2.2 * std::pow(qp_array(i,j,k) , 0.875);
dq_clwater_to_rain = dtn *(accrr*qcc + autor*(qcc - 0.001));
dq_clwater_to_rain = dtn *(accrr*qcc + autor*(qcc - qcw0));

// If the amount of change is more than the amount of qc present, then dq = qc
dq_clwater_to_rain = std::min(dq_clwater_to_rain, qc_array(i,j,k));
Expand Down
1 change: 1 addition & 0 deletions Source/Microphysics/Microphysics.H
Original file line number Diff line number Diff line change
Expand Up @@ -58,6 +58,7 @@ public:
static MoistureModelType modelType (const MoistureType a_moisture_type)
{
if ( (a_moisture_type == MoistureType::SAM)
|| (a_moisture_type == MoistureType::SAM_NoIce)
|| (a_moisture_type == MoistureType::SAM_NoPrecip_NoIce)
|| (a_moisture_type == MoistureType::Kessler)
|| (a_moisture_type == MoistureType::Kessler_NoRain)
Expand Down
45 changes: 8 additions & 37 deletions Source/Microphysics/SAM/Cloud_SAM.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -8,7 +8,9 @@ using namespace amrex;
/**
* Split cloud components according to saturation pressures; source theta from latent heat.
*/
void SAM::Cloud (const SolverChoice& solverChoice) {
void
SAM::Cloud (const SolverChoice& sc)
{

constexpr Real an = 1.0/(tbgmax-tbgmin);
constexpr Real bn = tbgmin*an;
Expand All @@ -18,13 +20,9 @@ void SAM::Cloud (const SolverChoice& solverChoice) {
Real fac_fus = m_fac_fus;
Real rdOcp = m_rdOcp;

SolverChoice sc = solverChoice;

int SAM_moisture_type = 1;

if(solverChoice.moisture_type == MoistureType::SAM) {
SAM_moisture_type = 1;
} else if(solverChoice.moisture_type == MoistureType::SAM_NoPrecip_NoIce) {
if (sc.moisture_type == MoistureType::SAM_NoIce ||
sc.moisture_type == MoistureType::SAM_NoPrecip_NoIce) {
SAM_moisture_type = 2;
}

Expand Down Expand Up @@ -128,18 +126,9 @@ void SAM::Cloud (const SolverChoice& solverChoice) {
qv_array , qcl_array , qci_array,
qn_array , qt_array);

// Update theta from temperature (it is essential to do this BEFORE the pressure is updated)
// This would be inconsistent with updating the pressure as part of the iteration above.
// Empirically based on the moist bubble rise case, getting the correct theta here
// depends on using the old (unchanged by the phase changes) pressure.
// Update theta
theta_array(i,j,k) = getThgivenPandT(tabs_array(i,j,k), 100.0*pres_array(i,j,k), rdOcp);

// Update pressure to be consistent with updated theta_array
pres_array(i,j,k) = getPgivenRTh(rho_array(i,j,k)*theta_array(i,j,k),qv_array(i,j,k));

// Pressure unit conversion
pres_array(i,j,k) *= 0.01;

//
// We cannot blindly relax to qsat, but we can convert qc/qi -> qv.
// The concept here is that if we put all the moisture into qv and modify
Expand All @@ -163,18 +152,9 @@ void SAM::Cloud (const SolverChoice& solverChoice) {
// Update temperature (endothermic since we evap/sublime)
tabs_array(i,j,k) -= fac_cond * delta_qc + fac_sub * delta_qi;

// Update theta from temperature (it is essential to do this BEFORE the pressure is updated)
// This would be inconsistent with updating the pressure as part of the iteration above.
// Empirically based on the moist bubble rise case, getting the correct theta here
// depends on using the old (unchanged by the phase changes) pressure.
// Update theta
theta_array(i,j,k) = getThgivenPandT(tabs_array(i,j,k), 100.0*pres_array(i,j,k), rdOcp);

// Update pressure to be consistent with updated theta_array
pres_array(i,j,k) = getPgivenRTh(rho_array(i,j,k)*theta_array(i,j,k),qv_array(i,j,k));

// Pressure unit conversion
pres_array(i,j,k) *= 0.01;

// Verify assumption that qv > qsat does not occur
erf_qsatw(tabs_array(i,j,k), pres_array(i,j,k), qsatw);
erf_qsati(tabs_array(i,j,k), pres_array(i,j,k), qsati);
Expand All @@ -189,18 +169,9 @@ void SAM::Cloud (const SolverChoice& solverChoice) {
qv_array , qcl_array , qci_array,
qn_array , qt_array);

// Update theta from temperature (it is essential to do this BEFORE the pressure is updated)
// This would be inconsistent with updating the pressure as part of the iteration above.
// Empirically based on the moist bubble rise case, getting the correct theta here
// depends on using the old (unchanged by the phase changes) pressure.
// Update theta
theta_array(i,j,k) = getThgivenPandT(tabs_array(i,j,k), 100.0*pres_array(i,j,k), rdOcp);

// Update pressure to be consistent with updated theta_array
pres_array(i,j,k) = getPgivenRTh(rho_array(i,j,k)*theta_array(i,j,k),qv_array(i,j,k));

// Pressure unit conversion
pres_array(i,j,k) *= 0.01;

}
}
});
Expand Down
6 changes: 5 additions & 1 deletion Source/Microphysics/SAM/IceFall.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -7,7 +7,11 @@ using namespace amrex;
/**
* Sedimentation of cloud ice (A32)
*/
void SAM::IceFall () {
void SAM::IceFall (const SolverChoice& sc) {

if(sc.moisture_type == MoistureType::SAM_NoIce ||
sc.moisture_type == MoistureType::SAM_NoPrecip_NoIce)
return;

Real dz = m_geom.CellSize(2);
Real dtn = dt;
Expand Down
6 changes: 4 additions & 2 deletions Source/Microphysics/SAM/Init_SAM.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -18,7 +18,8 @@ using namespace amrex;
* @param[in] geom Geometry associated with these MultiFabs and grids
* @param[in] dt_advance Timestep for the advance
*/
void SAM::Init (const MultiFab& cons_in,
void
SAM::Init (const MultiFab& cons_in,
const BoxArray& grids,
const Geometry& geom,
const Real& dt_advance,
Expand Down Expand Up @@ -83,7 +84,8 @@ void SAM::Init (const MultiFab& cons_in,
*
* @param[in] cons_in Conserved variables input
*/
void SAM::Copy_State_to_Micro (const MultiFab& cons_in)
void
SAM::Copy_State_to_Micro (const MultiFab& cons_in)
{
// Get the temperature, density, theta, qt and qp from input
for ( MFIter mfi(cons_in); mfi.isValid(); ++mfi) {
Expand Down
Loading
Loading