From 82acbf5196637b0f68475be7a3e4dbe41190f457 Mon Sep 17 00:00:00 2001 From: yairchn Date: Mon, 9 Aug 2021 10:44:25 -0500 Subject: [PATCH] move nh_pressure to closures move nh_pressure to closures update tables --- integration_tests/ARM_SGP.jl | 24 ++++---- integration_tests/Bomex.jl | 26 ++++----- integration_tests/DYCOMS_RF01.jl | 27 +++++---- integration_tests/Nieuwstadt.jl | 16 +++--- integration_tests/Rico.jl | 26 ++++----- integration_tests/Soares.jl | 21 ++++--- integration_tests/TRMM_LBA.jl | 26 ++++----- src/Operators.jl | 3 + src/TurbulenceConvection.jl | 1 + src/Turbulence_PrognosticTKE.jl | 79 +++++++-------------------- src/closures/perturbation_pressure.jl | 56 +++++++++++++++++++ src/turbulence_functions.jl | 33 ----------- src/types.jl | 44 +-------------- 13 files changed, 165 insertions(+), 217 deletions(-) create mode 100644 src/closures/perturbation_pressure.jl diff --git a/integration_tests/ARM_SGP.jl b/integration_tests/ARM_SGP.jl index c330d5ad9..e379ec2c9 100644 --- a/integration_tests/ARM_SGP.jl +++ b/integration_tests/ARM_SGP.jl @@ -11,18 +11,18 @@ include(joinpath("utils", "compute_mse.jl")) using .NameList best_mse = OrderedDict() -best_mse["qt_mean"] = 3.0998954346529128e-01 -best_mse["updraft_area"] = 1.9456011695010488e+03 -best_mse["updraft_w"] = 3.1996277691153108e+02 -best_mse["updraft_qt"] = 1.2645380500158543e+01 -best_mse["updraft_thetal"] = 2.7684653720274000e+01 -best_mse["u_mean"] = 8.7998547277817863e+01 -best_mse["tke_mean"] = 6.2536497518678300e+02 -best_mse["temperature_mean"] = 1.3316949170699788e-04 -best_mse["ql_mean"] = 2.5314704151953435e+02 -best_mse["thetal_mean"] = 1.3123118572261322e-04 -best_mse["Hvar_mean"] = 1.0770693519961310e+03 -best_mse["QTvar_mean"] = 1.9003025309800736e+02 +best_mse["qt_mean"] = 3.4761388491700224e-01 +best_mse["updraft_area"] = 1.9655632954763928e+03 +best_mse["updraft_w"] = 3.2377093752385053e+02 +best_mse["updraft_qt"] = 1.2452570404535185e+01 +best_mse["updraft_thetal"] = 2.7685659900036992e+01 +best_mse["u_mean"] = 8.7998547277817906e+01 +best_mse["tke_mean"] = 5.9501141635391093e+02 +best_mse["temperature_mean"] = 1.4210607813065063e-04 +best_mse["ql_mean"] = 2.1890556031513017e+02 +best_mse["thetal_mean"] = 1.3953340654098429e-04 +best_mse["Hvar_mean"] = 1.1148198912643536e+03 +best_mse["QTvar_mean"] = 3.0703058500395360e+02 @testset "ARM_SGP" begin println("Running ARM_SGP...") diff --git a/integration_tests/Bomex.jl b/integration_tests/Bomex.jl index b83b1dd90..2350e5100 100644 --- a/integration_tests/Bomex.jl +++ b/integration_tests/Bomex.jl @@ -11,19 +11,19 @@ include(joinpath("utils", "compute_mse.jl")) using .NameList best_mse = OrderedDict() -best_mse["qt_mean"] = 1.1113745694252045e-01 -best_mse["updraft_area"] = 7.2538205557139236e+02 -best_mse["updraft_w"] = 2.5943804836195852e+01 -best_mse["updraft_qt"] = 4.0075397466347100e+00 -best_mse["updraft_thetal"] = 2.1549266568293614e+01 -best_mse["v_mean"] = 6.5506192741355221e+01 -best_mse["u_mean"] = 5.3295866692805944e+01 -best_mse["tke_mean"] = 3.9596900219869937e+01 -best_mse["temperature_mean"] = 4.6114814306316929e-05 -best_mse["ql_mean"] = 6.7039710795893939e+00 -best_mse["thetal_mean"] = 4.6514684598969401e-05 -best_mse["Hvar_mean"] = 2.5387994922244559e+02 -best_mse["QTvar_mean"] = 4.5843881812484575e+01 +best_mse["qt_mean"] = 1.1068396401839847e-01 +best_mse["updraft_area"] = 7.2515940450436165e+02 +best_mse["updraft_w"] = 2.5856202193762641e+01 +best_mse["updraft_qt"] = 4.0061911576818243e+00 +best_mse["updraft_thetal"] = 2.1549386071253458e+01 +best_mse["v_mean"] = 6.5501952940732622e+01 +best_mse["u_mean"] = 5.3295928431425537e+01 +best_mse["tke_mean"] = 3.9582886707067345e+01 +best_mse["temperature_mean"] = 4.5842727523070819e-05 +best_mse["ql_mean"] = 6.5924727808279036e+00 +best_mse["thetal_mean"] = 4.6222694206747763e-05 +best_mse["Hvar_mean"] = 2.5900864672158934e+02 +best_mse["QTvar_mean"] = 4.6280994372607324e+01 @testset "Bomex" begin println("Running Bomex...") diff --git a/integration_tests/DYCOMS_RF01.jl b/integration_tests/DYCOMS_RF01.jl index b0c225e80..082f31a2c 100644 --- a/integration_tests/DYCOMS_RF01.jl +++ b/integration_tests/DYCOMS_RF01.jl @@ -11,20 +11,19 @@ include(joinpath("utils", "compute_mse.jl")) using .NameList best_mse = OrderedDict() -best_mse["qt_mean"] = 1.5987989356173508e-02 -best_mse["ql_mean"] = 5.9262542639081852e+00 -best_mse["updraft_area"] = 2.3887443993964089e+02 -best_mse["updraft_w"] = 4.3317357122037770e+00 -best_mse["updraft_qt"] = 1.1983341056077728e+00 -best_mse["updraft_thetal"] = 1.2739627176096105e+01 -best_mse["v_mean"] = 3.9736878160027324e+01 -best_mse["u_mean"] = 3.7037227811674803e+01 -best_mse["tke_mean"] = 1.4955313626928097e+01 -best_mse["temperature_mean"] = 1.6122516172306611e-05 -best_mse["thetal_mean"] = 1.7630641025752975e-05 -best_mse["Hvar_mean"] = 1.2452024821778659e+04 -best_mse["QTvar_mean"] = 7.6626261862502452e+02 - +best_mse["qt_mean"] = 1.5987989288888820e-02 +best_mse["ql_mean"] = 5.9262500323990217e+00 +best_mse["updraft_area"] = 2.3887443639936907e+02 +best_mse["updraft_w"] = 4.3317426538007693e+00 +best_mse["updraft_qt"] = 1.1983341214387264e+00 +best_mse["updraft_thetal"] = 1.2739627177179765e+01 +best_mse["v_mean"] = 3.9736878190256498e+01 +best_mse["u_mean"] = 3.7037227839730683e+01 +best_mse["tke_mean"] = 1.4955295527017412e+01 +best_mse["temperature_mean"] = 1.6122517253407575e-05 +best_mse["thetal_mean"] = 1.7630642054888154e-05 +best_mse["Hvar_mean"] = 1.2452017842800817e+04 +best_mse["QTvar_mean"] = 7.6626218978684426e+02 @testset "DYCOMS_RF01" begin println("Running DYCOMS_RF01...") diff --git a/integration_tests/Nieuwstadt.jl b/integration_tests/Nieuwstadt.jl index a215d6249..bdb26eff0 100644 --- a/integration_tests/Nieuwstadt.jl +++ b/integration_tests/Nieuwstadt.jl @@ -11,14 +11,14 @@ include(joinpath("utils", "compute_mse.jl")) using .NameList best_mse = OrderedDict() -best_mse["updraft_area"] = 5.9123316810068854e+02 -best_mse["updraft_w"] = 2.6743480109271673e+01 -best_mse["updraft_thetal"] = 3.0475272195292995e+01 -best_mse["u_mean"] = 1.5245784141331271e+02 -best_mse["tke_mean"] = 7.3707701211117509e+01 -best_mse["temperature_mean"] = 1.2199756418560522e-05 -best_mse["thetal_mean"] = 1.2066937466844642e-05 -best_mse["Hvar_mean"] = 1.9140385882308689e+02 +best_mse["updraft_area"] = 5.9103576007452034e+02 +best_mse["updraft_w"] = 2.6741082264407293e+01 +best_mse["updraft_thetal"] = 3.0475272117346940e+01 +best_mse["u_mean"] = 1.5245873874565478e+02 +best_mse["tke_mean"] = 7.3707728452173882e+01 +best_mse["temperature_mean"] = 1.2199914281992930e-05 +best_mse["thetal_mean"] = 1.2067109520894632e-05 +best_mse["Hvar_mean"] = 1.9140690462510639e+02 @testset "Nieuwstadt" begin println("Running Nieuwstadt...") diff --git a/integration_tests/Rico.jl b/integration_tests/Rico.jl index 8ea01450f..f2e8d5a5a 100644 --- a/integration_tests/Rico.jl +++ b/integration_tests/Rico.jl @@ -11,19 +11,19 @@ include(joinpath("utils", "compute_mse.jl")) using .NameList best_mse = OrderedDict() -best_mse["qt_mean"] = 4.2382033011491649e-01 -best_mse["updraft_area"] = 1.6696925513643459e+03 -best_mse["updraft_w"] = 2.8274002139289530e+02 -best_mse["updraft_qt"] = 1.7817151504225862e+01 -best_mse["updraft_thetal"] = 6.3573789244725759e+01 -best_mse["v_mean"] = 1.0633863265842136e+02 -best_mse["u_mean"] = 1.1434504431562398e+02 -best_mse["tke_mean"] = 8.6963627090840509e+02 -best_mse["temperature_mean"] = 1.2597722889848628e-04 -best_mse["ql_mean"] = 8.8699692429062537e+01 -best_mse["thetal_mean"] = 1.1916582780639058e-04 -best_mse["Hvar_mean"] = 2.0710735541431495e+03 -best_mse["QTvar_mean"] = 8.8908792485557092e+02 +best_mse["qt_mean"] = 4.2359420114784679e-01 +best_mse["updraft_area"] = 1.6722742675049749e+03 +best_mse["updraft_w"] = 2.8300484466673169e+02 +best_mse["updraft_qt"] = 1.7904453430368505e+01 +best_mse["updraft_thetal"] = 6.3573804437983789e+01 +best_mse["v_mean"] = 1.0634096403180666e+02 +best_mse["u_mean"] = 1.1434656097951681e+02 +best_mse["tke_mean"] = 8.6824438228445820e+02 +best_mse["temperature_mean"] = 1.2576042829475203e-04 +best_mse["ql_mean"] = 8.9402683009257942e+01 +best_mse["thetal_mean"] = 1.1896802334339157e-04 +best_mse["Hvar_mean"] = 2.0228054123757493e+03 +best_mse["QTvar_mean"] = 8.7748179858993274e+02 @testset "Rico" begin println("Running Rico...") diff --git a/integration_tests/Soares.jl b/integration_tests/Soares.jl index 5e8b9473f..2d85c0e95 100644 --- a/integration_tests/Soares.jl +++ b/integration_tests/Soares.jl @@ -11,17 +11,16 @@ include(joinpath("utils", "compute_mse.jl")) using .NameList best_mse = OrderedDict() -best_mse["qt_mean"] = 1.9338238569368907e-01 -best_mse["updraft_area"] = 4.8852480428600916e+02 -best_mse["updraft_w"] = 3.1250380743373622e+01 -best_mse["updraft_qt"] = 1.1136647638388030e+01 -best_mse["updraft_thetal"] = 2.2395194056858262e+01 -best_mse["u_mean"] = 7.2109941452760484e+02 -best_mse["tke_mean"] = 4.7919320017479706e+01 -best_mse["temperature_mean"] = 1.4374232234123472e-05 -best_mse["thetal_mean"] = 1.4379965951428025e-05 -best_mse["Hvar_mean"] = 2.3305430832031414e+02 - +best_mse["qt_mean"] = 1.9338238625419904e-01 +best_mse["updraft_area"] = 4.8852467389339699e+02 +best_mse["updraft_w"] = 3.1250501928635639e+01 +best_mse["updraft_qt"] = 1.1136647539128006e+01 +best_mse["updraft_thetal"] = 2.2395194054335001e+01 +best_mse["u_mean"] = 7.2109940121868794e+02 +best_mse["tke_mean"] = 4.7919324815294324e+01 +best_mse["temperature_mean"] = 1.4374232147002715e-05 +best_mse["thetal_mean"] = 1.4379965874656967e-05 +best_mse["Hvar_mean"] = 2.3305425470797826e+02 @testset "Soares" begin println("Running Soares...") diff --git a/integration_tests/TRMM_LBA.jl b/integration_tests/TRMM_LBA.jl index e98149ede..0ceb08ef3 100644 --- a/integration_tests/TRMM_LBA.jl +++ b/integration_tests/TRMM_LBA.jl @@ -11,19 +11,19 @@ include(joinpath("utils", "compute_mse.jl")) using .NameList best_mse = OrderedDict() -best_mse["qt_mean"] = 1.7430143247852277e+00 -best_mse["updraft_area"] = 3.2544587990881948e+04 -best_mse["updraft_w"] = 9.2113862144403333e+02 -best_mse["updraft_qt"] = 3.0303277939652101e+01 -best_mse["updraft_thetal"] = 1.1108868468322090e+02 -best_mse["v_mean"] = 2.9291515854682558e+02 -best_mse["u_mean"] = 1.6874269993592438e+03 -best_mse["tke_mean"] = 1.5877278080585038e+03 -best_mse["temperature_mean"] = 6.3087898742812351e-04 -best_mse["ql_mean"] = 9.4714318856310877e+02 -best_mse["thetal_mean"] = 4.1952114677068139e-04 -best_mse["Hvar_mean"] = 8.1299695367399718e+03 -best_mse["QTvar_mean"] = 2.7048653627007479e+03 +best_mse["qt_mean"] = 1.7430285968163504e+00 +best_mse["updraft_area"] = 3.2545401197014067e+04 +best_mse["updraft_w"] = 9.2113782722003907e+02 +best_mse["updraft_qt"] = 3.0306781659879494e+01 +best_mse["updraft_thetal"] = 1.1108867837608678e+02 +best_mse["v_mean"] = 2.9291577336310036e+02 +best_mse["u_mean"] = 1.6874237352911189e+03 +best_mse["tke_mean"] = 1.5857027462501208e+03 +best_mse["temperature_mean"] = 6.3088676948654218e-04 +best_mse["ql_mean"] = 9.4722581733362244e+02 +best_mse["thetal_mean"] = 4.1952296308447006e-04 +best_mse["Hvar_mean"] = 8.1275664794363738e+03 +best_mse["QTvar_mean"] = 2.7043370821860535e+03 @testset "TRMM_LBA" begin println("Running TRMM_LBA...") diff --git a/src/Operators.jl b/src/Operators.jl index f8e2aa1d7..fc4f6be31 100644 --- a/src/Operators.jl +++ b/src/Operators.jl @@ -16,6 +16,9 @@ struct Extrapolate end function ∇f2c(f, grid::Grid, k::Int) return (f[k] - f[k - 1]) * grid.dzi end +function ∇f2c(f, grid::Grid, k::Int, i::Int) + return (f[i, k] - f[i, k - 1]) * grid.dzi +end function ∇f2c(f_dual::SVector, grid::Grid, k::Int; bottom = NoBCGivenError(), top = NoBCGivenError()) if is_surface_face(grid, k - 1) return ∇f2c(f_dual, grid, BottomBCTag(), bottom) diff --git a/src/TurbulenceConvection.jl b/src/TurbulenceConvection.jl index 834879f47..c7426749f 100644 --- a/src/TurbulenceConvection.jl +++ b/src/TurbulenceConvection.jl @@ -68,5 +68,6 @@ include("Radiation.jl") include("forcing_functions.jl") include("Surface.jl") include("surface_functions.jl") +include("closures/perturbation_pressure.jl") end diff --git a/src/Turbulence_PrognosticTKE.jl b/src/Turbulence_PrognosticTKE.jl index a77db67af..58bd0159a 100755 --- a/src/Turbulence_PrognosticTKE.jl +++ b/src/Turbulence_PrognosticTKE.jl @@ -24,8 +24,6 @@ function initialize_io(self::EDMF_PrognosticTKE, Stats::NetCDFIO_Stats) add_profile(Stats, "nh_pressure_drag") add_profile(Stats, "nh_pressure_b") add_profile(Stats, "asp_ratio") - add_profile(Stats, "b_coeff") - add_profile(Stats, "horiz_K_eddy") add_profile(Stats, "sorting_function") add_profile(Stats, "b_mix") @@ -99,8 +97,6 @@ function io(self::EDMF_PrognosticTKE, Stats::NetCDFIO_Stats, TS::TimeStepping) mean_nh_pressure_drag = center_field(grid) mean_nh_pressure_b = center_field(grid) mean_asp_ratio = center_field(grid) - mean_b_coeff = center_field(grid) - mean_detr_sc = center_field(grid) massflux = face_field(grid) mf_h = face_field(grid) @@ -135,8 +131,6 @@ function io(self::EDMF_PrognosticTKE, Stats::NetCDFIO_Stats, TS::TimeStepping) self.UpdVar.Area.values[i, k] * self.nh_pressure_drag[i, k] / self.UpdVar.Area.bulkvalues[k] mean_asp_ratio[k] += self.UpdVar.Area.values[i, k] * self.asp_ratio[i, k] / self.UpdVar.Area.bulkvalues[k] - mean_b_coeff[k] += self.UpdVar.Area.values[i, k] * self.b_coeff[i, k] / self.UpdVar.Area.bulkvalues[k] - mean_frac_turb_entr[k] += self.UpdVar.Area.values[i, k] * self.frac_turb_entr[i, k] / self.UpdVar.Area.bulkvalues[k] mean_horiz_K_eddy[k] += @@ -159,8 +153,6 @@ function io(self::EDMF_PrognosticTKE, Stats::NetCDFIO_Stats, TS::TimeStepping) write_profile(Stats, "nh_pressure_drag", mean_nh_pressure_drag[cinterior]) write_profile(Stats, "nh_pressure_b", mean_nh_pressure_b[cinterior]) write_profile(Stats, "asp_ratio", mean_asp_ratio[cinterior]) - write_profile(Stats, "b_coeff", mean_b_coeff[cinterior]) - write_profile(Stats, "massflux", massflux[cinterior]) write_profile(Stats, "massflux_h", mf_h[cinterior]) write_profile(Stats, "massflux_qt", mf_qt[cinterior]) @@ -931,9 +923,6 @@ function compute_updraft_closures(self::EDMF_PrognosticTKE, GMV::GridMeanVariabl input.c_ed_mf = self.entrainment_ed_mf_sigma input.chi_upd = self.updraft_mixing_frac input.tke_coef = self.entrainment_smin_tke_coeff - ret_b = pressure_buoy_struct() - ret_w = pressure_drag_struct() - input_p = pressure_in_struct() # TODO: merge this with grid.cinterior kc_surf = kc_surface(grid) @@ -948,10 +937,8 @@ function compute_updraft_closures(self::EDMF_PrognosticTKE, GMV::GridMeanVariabl @inbounds for k in real_center_indicies(grid) @inbounds for i in xrange(self.n_updrafts) - input_p.updraft_top = self.UpdVar.updraft_top[i] alen = max(length(argwhere(self.UpdVar.Area.values[i, cinterior])), 1) avals = self.UpdVar.Area.values[i, cinterior] - input_p.a_med = Statistics.median(avals[1:alen]) input.zi = self.UpdVar.cloud_base[i] # entrainment input.buoy_ed_flux = self.EnvVar.TKE.buoy[k] @@ -967,7 +954,6 @@ function compute_updraft_closures(self::EDMF_PrognosticTKE, GMV::GridMeanVariabl input.w_env = interpf2c(self.EnvVar.W.values, grid, k) input.ql_up = self.UpdVar.QL.values[i, k] input.ql_env = self.EnvVar.QL.values[k] - input.nh_pressure = self.nh_pressure[i, k] input.RH_upd = self.UpdVar.RH.values[i, k] input.RH_env = self.EnvVar.RH.values[k] @@ -1015,55 +1001,32 @@ function compute_updraft_closures(self::EDMF_PrognosticTKE, GMV::GridMeanVariabl end # pressure - input_p.a_kfull = interpc2f(self.UpdVar.Area.values, grid, k, i) - if input_p.a_kfull >= self.minimum_area - input_p.dzi = grid.dzi - - input_p.b_kfull = interpc2f(self.UpdVar.B.values, grid, k, i) - input_p.rho0_kfull = ref_state.rho0[k] - input_p.alpha1 = self.pressure_normalmode_buoy_coeff1 - input_p.alpha2 = self.pressure_normalmode_buoy_coeff2 - input_p.beta1 = self.pressure_normalmode_adv_coeff - input_p.beta2 = self.pressure_normalmode_drag_coeff - input_p.w_kfull = self.UpdVar.W.values[i, k] - input_p.w_kmfull = self.UpdVar.W.values[i, k - 1] - input_p.w_kenv = self.EnvVar.W.values[k] - - if self.asp_label == "z_dependent" - input_p.asp_ratio = input_p.updraft_top / 2.0 / sqrt(input_p.a_kfull) / input_p.rd - elseif self.asp_label == "median" - input_p.asp_ratio = input_p.updraft_top / 2.0 / sqrt(input_p.a_med) / input_p.rd - elseif self.asp_label == "const" - input_p.asp_ratio = 1.0 - end - - if input_p.a_kfull > 0.0 - ret_b = self.pressure_func_buoy(input_p) - ret_w = self.pressure_func_drag(input_p) - self.nh_pressure_b[i, k] = ret_b.nh_pressure_b - self.nh_pressure_adv[i, k] = ret_w.nh_pressure_adv - self.nh_pressure_drag[i, k] = ret_w.nh_pressure_drag - - self.b_coeff[i, k] = ret_b.b_coeff - self.asp_ratio[i, k] = input_p.asp_ratio - - else - self.nh_pressure_b[i, k] = 0.0 - self.nh_pressure_adv[i, k] = 0.0 - self.nh_pressure_drag[i, k] = 0.0 - - self.b_coeff[i, k] = 0.0 - self.asp_ratio[i, k] = 0.0 - end + a_kfull = interpc2f(self.UpdVar.Area.values, grid, k, i) + if a_kfull > 0.0 + b_kfull = interpc2f(self.UpdVar.B.values, grid, k, i) + ∇w_up = ∇f2c(self.UpdVar.W.values, grid, k, i) + asp_ratio = 1.0 + self.nh_pressure_b[i, k], self.nh_pressure_adv[i, k], self.nh_pressure_drag[i, k] = + perturbation_pressure( + self.UpdVar.updraft_top[i], + 500.0, + a_kfull, + b_kfull, + ref_state.rho0[k], + self.pressure_normalmode_buoy_coeff1, + self.pressure_normalmode_buoy_coeff2, + self.pressure_normalmode_adv_coeff, + self.pressure_normalmode_drag_coeff, + self.UpdVar.W.values[i, k], + ∇w_up, + self.EnvVar.W.values[k], + asp_ratio, + ) else self.nh_pressure_b[i, k] = 0.0 self.nh_pressure_adv[i, k] = 0.0 self.nh_pressure_drag[i, k] = 0.0 - - self.b_coeff[i, k] = 0.0 - self.asp_ratio[i, k] = 0.0 end - self.nh_pressure[i, k] = self.nh_pressure_b[i, k] + self.nh_pressure_adv[i, k] + self.nh_pressure_drag[i, k] end end diff --git a/src/closures/perturbation_pressure.jl b/src/closures/perturbation_pressure.jl new file mode 100644 index 000000000..b2fc80ae9 --- /dev/null +++ b/src/closures/perturbation_pressure.jl @@ -0,0 +1,56 @@ +""" + perturbation_pressure( + updraft_top, + min_updraft_top, + a_up, + b_up, + ρ0_k, + α₁, + α₂, + β₁, + β₂, + w_up, + ∇w_up, + w_en, + asp_ratio + ) +Returns the value of perturbation pressure gradient +for updraft i following He et al. (JAMES, 2020), given: + - `updraft_top`: the height of the updraft in the previous timestep + - `min_updraft_top`: the minimal height of the updraft to avoid zero devision + - `a_up`: updraft area + - `b_up`: updraft buoyancy + - `ρ0_k`: reference density + - `α₁`: pressure closure free parameter + - `α₂`: pressure closure free parameter + - `β₁`: pressure closure free parameter + - `β₂`: pressure closure free parameter + - `w_up`: updraft vertical velocity + - `∇w_up`: updraft divergence of vertical velocity + - `w_en`: environment vertical velocity + - esp_ratio`: the specific aspect ratio of the updraft +""" +function perturbation_pressure( + updraft_top, + min_updraft_top, + a_up, + b_up, + ρ0_k, + α₁, + α₂, + β₁, + β₂, + w_up, + ∇w_up, + w_en, + asp_ratio, +) + + nh_press_buoy = -α₁ / (1 + α₂ * asp_ratio^2) * ρ0_k * a_up * b_up + + nh_pressure_adv = ρ0_k * a_up * β₁ * w_up * ∇w_up + # drag as w_dif and account for downdrafts + nh_pressure_drag = -1.0 * ρ0_k * a_up * β₂ * (w_up - w_en) * fabs(w_up - w_en) / fmax(updraft_top, min_updraft_top) + + return nh_press_buoy, nh_pressure_adv, nh_pressure_drag +end; diff --git a/src/turbulence_functions.jl b/src/turbulence_functions.jl index f6e02d509..27754264c 100644 --- a/src/turbulence_functions.jl +++ b/src/turbulence_functions.jl @@ -142,39 +142,6 @@ function entr_detr_env_moisture_deficit_div(param_set::APS, entr_in::entr_struct return _ret end - -function pressure_normalmode_buoy(press_in::pressure_in_struct) - _ret = pressure_buoy_struct() - - _ret.b_coeff = press_in.alpha1 / (1 + press_in.alpha2 * press_in.asp_ratio^2) - _ret.nh_pressure_b = -1.0 * press_in.rho0_kfull * press_in.a_kfull * press_in.b_kfull * _ret.b_coeff - - return _ret -end - -function pressure_normalmode_drag(press_in::pressure_in_struct) - _ret = pressure_drag_struct() - - _ret.nh_pressure_adv = - press_in.rho0_kfull * - press_in.a_kfull * - press_in.beta1 * - press_in.w_kfull * - (press_in.w_kfull - press_in.w_kmfull) * - press_in.dzi - - # drag as w_dif and account for downdrafts - _ret.nh_pressure_drag = - -1.0 * - press_in.rho0_kfull * - press_in.a_kfull * - press_in.beta2 * - (press_in.w_kfull - press_in.w_kenv) * - fabs(press_in.w_kfull - press_in.w_kenv) / fmax(press_in.updraft_top, 500.0) - - return _ret -end - # convective velocity scale function get_wstar(bflux, zi) return cbrt(fmax(bflux * zi, 0.0)) diff --git a/src/types.jl b/src/types.jl index 6e6f32480..29ee0adf1 100644 --- a/src/types.jl +++ b/src/types.jl @@ -1,32 +1,3 @@ -Base.@kwdef mutable struct pressure_drag_struct - nh_pressure_adv::Float64 = 0 - nh_pressure_drag::Float64 = 0 -end - -Base.@kwdef mutable struct pressure_buoy_struct - b_coeff::Float64 = 0 - nh_pressure_b::Float64 = 0 -end - -Base.@kwdef mutable struct pressure_in_struct - updraft_top::Float64 = 0 - asp_label::Float64 = 0 - a_med::Float64 = 0 - a_kfull::Float64 = 0 - b_kfull::Float64 = 0 - rho0_kfull::Float64 = 0 - alpha1::Float64 = 0 - alpha2::Float64 = 0 - beta1::Float64 = 0 - beta2::Float64 = 0 - w_kfull::Float64 = 0 - w_kmfull::Float64 = 0 - w_kenv::Float64 = 0 - dzi::Float64 = 0 - drag_sign::Float64 = 0 - asp_ratio::Float64 = 0 -end - Base.@kwdef mutable struct entr_struct entr_sc::Float64 = 0 detr_sc::Float64 = 0 @@ -753,9 +724,7 @@ mutable struct EDMF_PrognosticTKE{PS, A1, A2} n_updrafts::Int use_const_plume_spacing::Bool entr_detr_fp::Function - pressure_func_buoy::Function drag_sign::Int - pressure_func_drag::Function asp_label extrapolate_buoyancy::Bool surface_area::Float64 @@ -799,7 +768,6 @@ mutable struct EDMF_PrognosticTKE{PS, A1, A2} nh_pressure_adv::A2 nh_pressure_drag::A2 asp_ratio::A2 - b_coeff::A2 m::A2 mixing_length::A1 horiz_K_eddy::A2 @@ -864,8 +832,6 @@ mutable struct EDMF_PrognosticTKE{PS, A1, A2} error("Turbulence--EDMF_PrognosticTKE: Entrainment rate namelist option is not recognized") end - pressure_func_buoy = pressure_normalmode_buoy - pressure_func_drag_str = parse_namelist( namelist, "turbulence", @@ -875,12 +841,10 @@ mutable struct EDMF_PrognosticTKE{PS, A1, A2} valid_options = ["normalmode", "normalmode_signdf"], ) - drag_sign = false - pressure_func_drag = if pressure_func_drag_str == "normalmode" - pressure_normalmode_drag + if pressure_func_drag_str == "normalmode" + drag_sign = false elseif pressure_func_drag_str == "normalmode_signdf" drag_sign = true - pressure_normalmode_drag end asp_label = parse_namelist( @@ -994,7 +958,6 @@ mutable struct EDMF_PrognosticTKE{PS, A1, A2} nh_pressure_adv = center_field(Gr, n_updrafts) nh_pressure_drag = center_field(Gr, n_updrafts) asp_ratio = center_field(Gr, n_updrafts) - b_coeff = center_field(Gr, n_updrafts) # Mass flux m = face_field(Gr, n_updrafts) @@ -1058,9 +1021,7 @@ mutable struct EDMF_PrognosticTKE{PS, A1, A2} n_updrafts, use_const_plume_spacing, entr_detr_fp, - pressure_func_buoy, drag_sign, - pressure_func_drag, asp_label, extrapolate_buoyancy, surface_area, @@ -1104,7 +1065,6 @@ mutable struct EDMF_PrognosticTKE{PS, A1, A2} nh_pressure_adv, nh_pressure_drag, asp_ratio, - b_coeff, m, mixing_length, horiz_K_eddy,