Skip to content

Commit

Permalink
WIP
Browse files Browse the repository at this point in the history 
WIP

add pressure closure file

cleanups in entrainment

add entr_detr, incomplete

move nh_pressure to closures

move nh_pressure to closures

update tables

running with entr in closures

running with entr in closures

working with entr_detr

apply formatter

update tables

resolve conflicts with yc/call_press_closure

remove cleanup_covariance

add entr_detr_model + merge turb entr + cleanup

fix bug

apply formatter

WIP

WIP

add pressure closure file

cleanups in entrainment

add entr_detr, incomplete

move nh_pressure to closures

move nh_pressure to closures

update tables

running with entr in closures

running with entr in closures

working with entr_detr

apply formatter

update tables

resolve conflicts with yc/call_press_closure

remove cleanup_covariance

add entr_detr_model + merge turb entr + cleanup

fix bug

apply formatter

update tables
  • Loading branch information
yairchn committed Aug 10, 2021
1 parent 94a3f72 commit cbd3d5d
Show file tree
Hide file tree
Showing 12 changed files with 254 additions and 390 deletions.
24 changes: 12 additions & 12 deletions integration_tests/ARM_SGP.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -11,18 +11,18 @@ include(joinpath("utils", "compute_mse.jl"))
using .NameList

best_mse = OrderedDict()
best_mse["qt_mean"] = 3.4297134939930601e-01
best_mse["updraft_area"] = 1.9829008493073268e+03
best_mse["updraft_w"] = 3.1072236638537851e+02
best_mse["updraft_qt"] = 1.2698381273033208e+01
best_mse["updraft_thetal"] = 2.7686390354792074e+01
best_mse["u_mean"] = 8.7998547277817863e+01
best_mse["tke_mean"] = 6.7319699431828644e+02
best_mse["temperature_mean"] = 1.3960321509272764e-04
best_mse["ql_mean"] = 2.9854035802671558e+02
best_mse["thetal_mean"] = 1.4200441599983529e-04
best_mse["Hvar_mean"] = 1.3755159256021079e+03
best_mse["QTvar_mean"] = 2.6566877259844244e+02
best_mse["qt_mean"] = 3.4420193413687117e-01
best_mse["updraft_area"] = 1.9818317322584239e+03
best_mse["updraft_w"] = 3.1101869507553300e+02
best_mse["updraft_qt"] = 1.2779005319987142e+01
best_mse["updraft_thetal"] = 2.7687174775467490e+01
best_mse["u_mean"] = 8.7998547277817835e+01
best_mse["tke_mean"] = 6.7532125088372368e+02
best_mse["temperature_mean"] = 1.3976231155597154e-04
best_mse["ql_mean"] = 2.9406894196734129e+02
best_mse["thetal_mean"] = 1.4216021128644291e-04
best_mse["Hvar_mean"] = 1.4219930941047571e+03
best_mse["QTvar_mean"] = 2.7979236616322686e+02

@testset "ARM_SGP" begin
println("Running ARM_SGP...")
Expand Down
26 changes: 13 additions & 13 deletions integration_tests/Bomex.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -11,19 +11,19 @@ include(joinpath("utils", "compute_mse.jl"))
using .NameList

best_mse = OrderedDict()
best_mse["qt_mean"] = 9.5080504754167433e-02
best_mse["updraft_area"] = 7.2534663873051477e+02
best_mse["updraft_w"] = 2.6358909128208378e+01
best_mse["updraft_qt"] = 4.0334302720334296e+00
best_mse["updraft_thetal"] = 2.1547995304586468e+01
best_mse["v_mean"] = 6.5231550206742128e+01
best_mse["u_mean"] = 5.3292889729286458e+01
best_mse["tke_mean"] = 3.8860563196392924e+01
best_mse["temperature_mean"] = 3.8436817065445734e-05
best_mse["ql_mean"] = 3.7878522086312807e+00
best_mse["thetal_mean"] = 3.9063754163395865e-05
best_mse["Hvar_mean"] = 7.2639655951863034e+01
best_mse["QTvar_mean"] = 2.2435640414468629e+01
best_mse["qt_mean"] = 9.5377282321874049e-02
best_mse["updraft_area"] = 7.2522420468798543e+02
best_mse["updraft_w"] = 2.6282677670021386e+01
best_mse["updraft_qt"] = 4.0341139540222635e+00
best_mse["updraft_thetal"] = 2.1548164483602367e+01
best_mse["v_mean"] = 6.5247595142753596e+01
best_mse["u_mean"] = 5.3292989796462209e+01
best_mse["tke_mean"] = 3.8891129116688681e+01
best_mse["temperature_mean"] = 3.8547176978568602e-05
best_mse["ql_mean"] = 3.8134281154981489e+00
best_mse["thetal_mean"] = 3.9162071364386570e-05
best_mse["Hvar_mean"] = 7.5510623226760771e+01
best_mse["QTvar_mean"] = 2.2962565404185188e+01

@testset "Bomex" begin
println("Running Bomex...")
Expand Down
16 changes: 8 additions & 8 deletions integration_tests/Nieuwstadt.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -11,14 +11,14 @@ include(joinpath("utils", "compute_mse.jl"))
using .NameList

best_mse = OrderedDict()
best_mse["updraft_area"] = 5.9412545314262150e+02
best_mse["updraft_w"] = 2.6952552348430249e+01
best_mse["updraft_thetal"] = 3.0475271481866773e+01
best_mse["u_mean"] = 1.5247657003047263e+02
best_mse["tke_mean"] = 7.3600391357561151e+01
best_mse["temperature_mean"] = 1.1971868506884125e-05
best_mse["thetal_mean"] = 1.2117924470979204e-05
best_mse["Hvar_mean"] = 1.8622261767371862e+02
best_mse["updraft_area"] = 5.9412544437504687e+02
best_mse["updraft_w"] = 2.6952552450364237e+01
best_mse["updraft_thetal"] = 3.0475271481865867e+01
best_mse["u_mean"] = 1.5247656989308092e+02
best_mse["tke_mean"] = 7.3600391623055913e+01
best_mse["temperature_mean"] = 1.1971868653164508e-05
best_mse["thetal_mean"] = 1.2117924617884567e-05
best_mse["Hvar_mean"] = 1.8622261701198534e+02

@testset "Nieuwstadt" begin
println("Running Nieuwstadt...")
Expand Down
26 changes: 13 additions & 13 deletions integration_tests/Rico.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -11,19 +11,19 @@ include(joinpath("utils", "compute_mse.jl"))
using .NameList

best_mse = OrderedDict()
best_mse["qt_mean"] = 4.2888329142054360e-01
best_mse["updraft_area"] = 1.5804026921690113e+03
best_mse["updraft_w"] = 3.0412249337738064e+02
best_mse["updraft_qt"] = 1.8041619418228048e+01
best_mse["updraft_thetal"] = 6.3565634303932690e+01
best_mse["v_mean"] = 1.0637419570945200e+02
best_mse["u_mean"] = 1.1444385357013196e+02
best_mse["tke_mean"] = 8.6459012010773438e+02
best_mse["temperature_mean"] = 1.2577596756983134e-04
best_mse["ql_mean"] = 5.0938207504953844e+01
best_mse["thetal_mean"] = 1.2604124901901089e-04
best_mse["Hvar_mean"] = 1.4598483199237567e+03
best_mse["QTvar_mean"] = 7.0850350571186209e+02
best_mse["qt_mean"] = 4.2745465446436787e-01
best_mse["updraft_area"] = 1.5810580213880007e+03
best_mse["updraft_w"] = 3.0293830804824944e+02
best_mse["updraft_qt"] = 1.8036840808532091e+01
best_mse["updraft_thetal"] = 6.3565923109974378e+01
best_mse["v_mean"] = 1.0637284782337001e+02
best_mse["u_mean"] = 1.1444354015737491e+02
best_mse["tke_mean"] = 8.6553999307149547e+02
best_mse["temperature_mean"] = 1.2540517259642333e-04
best_mse["ql_mean"] = 5.1325163153016852e+01
best_mse["thetal_mean"] = 1.2561593164473264e-04
best_mse["Hvar_mean"] = 1.4739277424672016e+03
best_mse["QTvar_mean"] = 7.1366176286414179e+02

@testset "Rico" begin
println("Running Rico...")
Expand Down
20 changes: 10 additions & 10 deletions integration_tests/Soares.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -11,16 +11,16 @@ include(joinpath("utils", "compute_mse.jl"))
using .NameList

best_mse = OrderedDict()
best_mse["qt_mean"] = 1.9368730351926020e-01
best_mse["updraft_area"] = 4.8515130524256091e+02
best_mse["updraft_w"] = 3.1260485555848359e+01
best_mse["updraft_qt"] = 1.1137263871861222e+01
best_mse["updraft_thetal"] = 2.2395183773350311e+01
best_mse["u_mean"] = 7.2111916707127898e+02
best_mse["tke_mean"] = 4.7882808390238374e+01
best_mse["temperature_mean"] = 1.3818642271008120e-05
best_mse["thetal_mean"] = 1.4383770217873919e-05
best_mse["Hvar_mean"] = 2.3235278609470686e+02
best_mse["qt_mean"] = 1.9368730493269354e-01
best_mse["updraft_area"] = 4.8515131142491634e+02
best_mse["updraft_w"] = 3.1260485803642322e+01
best_mse["updraft_qt"] = 1.1137263873078020e+01
best_mse["updraft_thetal"] = 2.2395183773363378e+01
best_mse["u_mean"] = 7.2111916666642423e+02
best_mse["tke_mean"] = 4.7882808859338290e+01
best_mse["temperature_mean"] = 1.3818642368852120e-05
best_mse["thetal_mean"] = 1.4383770315305301e-05
best_mse["Hvar_mean"] = 2.3235278783484819e+02

@testset "Soares" begin
println("Running Soares...")
Expand Down
4 changes: 2 additions & 2 deletions integration_tests/utils/generate_namelist.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -57,7 +57,7 @@ function default_namelist(case_name::String)
namelist_defaults["turbulence"]["EDMF_PrognosticTKE"]["max_area"] = 0.9
namelist_defaults["turbulence"]["EDMF_PrognosticTKE"]["entrainment_factor"] = 0.13
namelist_defaults["turbulence"]["EDMF_PrognosticTKE"]["detrainment_factor"] = 0.51
namelist_defaults["turbulence"]["EDMF_PrognosticTKE"]["entrainment_massflux_div_factor"] = 0.4
namelist_defaults["turbulence"]["EDMF_PrognosticTKE"]["entrainment_massflux_div_factor"] = 0.0
namelist_defaults["turbulence"]["EDMF_PrognosticTKE"]["turbulent_entrainment_factor"] = 0.015
namelist_defaults["turbulence"]["EDMF_PrognosticTKE"]["entrainment_ed_mf_sigma"] = 50.0
namelist_defaults["turbulence"]["EDMF_PrognosticTKE"]["entrainment_smin_tke_coeff"] = 0.3
Expand Down Expand Up @@ -345,7 +345,7 @@ function DryBubble(namelist_defaults)
namelist["meta"]["simname"] = "DryBubble"
namelist["meta"]["casename"] = "DryBubble"

namelist_defaults["turbulence"]["EDMF_PrognosticTKE"]["entrainment"] = "moisture_deficit_div"
namelist["turbulence"]["EDMF_PrognosticTKE"]["entrainment_massflux_div_factor"] = 0.4

return namelist
end
Expand Down
2 changes: 2 additions & 0 deletions src/TurbulenceConvection.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -73,5 +73,7 @@ include("forcing_functions.jl")
include("Surface.jl")
include("surface_functions.jl")
include("closures/perturbation_pressure.jl")
include("closures/entr_detr.jl")
include("closures/nondimensional_exchange_functions.jl")

end
161 changes: 37 additions & 124 deletions src/Turbulence_PrognosticTKE.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -146,8 +146,6 @@ function io(self::EDMF_PrognosticTKE, Stats::NetCDFIO_Stats, TS::TimeStepping)
write_profile(Stats, "horiz_K_eddy", mean_horiz_K_eddy[cinterior])
write_profile(Stats, "entrainment_sc", mean_entr_sc[cinterior])
write_profile(Stats, "detrainment_sc", mean_detr_sc[cinterior])
write_profile(Stats, "sorting_function", mean_sorting_function[cinterior])
write_profile(Stats, "b_mix", mean_b_mix[cinterior])
write_profile(Stats, "nh_pressure", mean_nh_pressure[cinterior])
write_profile(Stats, "nh_pressure_adv", mean_nh_pressure_adv[cinterior])
write_profile(Stats, "nh_pressure_drag", mean_nh_pressure_drag[cinterior])
Expand Down Expand Up @@ -902,102 +900,71 @@ function get_env_covar_from_GMV(
end

function compute_updraft_closures(self::EDMF_PrognosticTKE, GMV::GridMeanVariables, Case::CasesBase)
ret = entr_struct()
input = entr_in_struct()
εδ_model = entr_detr_model()
sa = eos_struct()
quadrature_order = 3
grid = get_grid(self)
param_set = parameter_set(self)

upd_cloud_diagnostics(self.UpdVar, reference_state(self))

input.wstar = self.wstar

input.dz = get_grid(self).dz
input.zbl = compute_zbl_qt_grad(GMV)
input.sort_pow = self.sorting_power
input.c_det = self.detrainment_factor
input.c_div = self.entrainment_Mdiv_factor
input.c_mu = self.entrainment_sigma
input.c_mu0 = self.entrainment_scale
input.c_ed_mf = self.entrainment_ed_mf_sigma
input.chi_upd = self.updraft_mixing_frac
input.tke_coef = self.entrainment_smin_tke_coeff

# TODO: merge this with grid.cinterior
kc_surf = kc_surface(grid)
kc_toa = kc_top_of_atmos(grid)
cinterior = kc_surf:kc_toa

ref_state = reference_state(self)
tau = get_mixing_tau(self.zi, self.wstar)

ae = 1 .- self.UpdVar.Area.bulkvalues
l = zeros(2)

@inbounds for k in real_center_indicies(grid)
@inbounds for i in xrange(self.n_updrafts)
alen = max(length(argwhere(self.UpdVar.Area.values[i, cinterior])), 1)
avals = self.UpdVar.Area.values[i, cinterior]
input.zi = self.UpdVar.cloud_base[i]
# entrainment
input.buoy_ed_flux = self.EnvVar.TKE.buoy[k]
w_upd = interpf2c(self.UpdVar.W.values, grid, k, i)
if self.UpdVar.Area.values[i, k] > 0.0
input.b_upd = self.UpdVar.B.values[i, k]
input.w_upd = interpf2c(self.UpdVar.W.values, grid, k, i)
input.z = grid.z_half[k]
input.a_upd = self.UpdVar.Area.values[i, k]
input.a_env = 1.0 - self.UpdVar.Area.bulkvalues[k]
input.tke = self.EnvVar.TKE.values[k]
input.ql_env = self.EnvVar.QL.values[k]
input.b_env = self.EnvVar.B.values[k]
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.RH_upd = self.UpdVar.RH.values[i, k]
input.RH_env = self.EnvVar.RH.values[k]

# compute dMdz at half levels
gmv_w_k = interpf2c(GMV.W.values, grid, k)
gmv_w_km = interpf2c(GMV.W.values, grid, k - 1)
upd_w_km = interpf2c(self.UpdVar.W.values, grid, k - 1, i)
input.M = input.a_upd * (input.w_upd - gmv_w_k)
M = self.UpdVar.Area.values[i, k] * (w_upd - gmv_w_k)
Mm = self.UpdVar.Area.values[i, k - 1] * (upd_w_km - gmv_w_km)
input.dMdz = (input.M - Mm) * grid.dzi
input.tke = self.EnvVar.TKE.values[k]

ret = self.entr_detr_fp(param_set, input)
self.entr_sc[i, k] = ret.entr_sc
self.detr_sc[i, k] = ret.detr_sc
self.sorting_function[i, k] = ret.sorting_function
self.b_mix[i, k] = ret.b_mix

dMdz = (M - Mm) * grid.dzi
w_min = 0.001

εδ_model.b_upd = self.UpdVar.B.values[i, k]
εδ_model.b_env = self.EnvVar.B.values[k]
εδ_model.w_upd = interpf2c(self.UpdVar.W.values, grid, k, i)
εδ_model.w_env = interpf2c(self.EnvVar.W.values, grid, k)
εδ_model.a_upd = self.UpdVar.Area.values[i, k]
εδ_model.a_env = self.EnvVar.Area.values[k]
εδ_model.ql_upd = self.UpdVar.QL.values[i, k]
εδ_model.ql_env = self.EnvVar.QL.values[k]
εδ_model.RH_upd = self.UpdVar.RH.values[i, k]
εδ_model.RH_env = self.EnvVar.RH.values[k]
εδ_model.M = M
εδ_model.dMdz = dMdz
εδ_model.tke = self.EnvVar.TKE.values[k]
εδ_model.n_up = self.n_updrafts
εδ_model.ρ = ref_state.rho0[k]
εδ_model.R_up = self.pressure_plume_spacing[i]

self.entr_sc[i, k], self.detr_sc[i, k], self.frac_turb_entr[i, k], self.horiz_K_eddy[i, k] = entr_detr(
param_set,
w_min,
self.sorting_power,
self.detrainment_factor,
self.entrainment_Mdiv_factor,
self.entrainment_sigma,
self.entrainment_scale,
self.updraft_mixing_frac,
self.entrainment_smin_tke_coeff,
self.turbulent_entrainment_factor,
εδ_model,
)
else
self.entr_sc[i, k] = 0.0
self.detr_sc[i, k] = 0.0
self.sorting_function[i, k] = 0.0
self.b_mix[i, k] = self.EnvVar.B.values[k]
end

# horizontal eddy diffusivities
if self.UpdVar.Area.values[i, k] > 0.0
self.horiz_K_eddy[i, k] =
self.UpdVar.Area.values[i, k] *
self.turbulent_entrainment_factor *
sqrt(fmax(self.EnvVar.TKE.values[k], 0.0)) *
self.pressure_plume_spacing[i]
else
self.horiz_K_eddy[i, k] = 0.0
end

# turbulent entrainment
wu_half = interp2pt(self.UpdVar.W.values[i, k], self.UpdVar.W.values[i, k - 1])
if self.UpdVar.Area.values[i, k] * wu_half > 0.0
self.frac_turb_entr[i, k] =
(2.0 / self.pressure_plume_spacing[i]^2.0) * self.horiz_K_eddy[i, k] / wu_half /
self.UpdVar.Area.values[i, k]
else
self.frac_turb_entr[i, k] = 0.0
self.horiz_K_eddy[i, k] = 0.0
end

# pressure
Expand Down Expand Up @@ -1033,16 +1000,6 @@ function compute_updraft_closures(self::EDMF_PrognosticTKE, GMV::GridMeanVariabl
return
end

function compute_zbl_qt_grad(GMV::GridMeanVariables)
grid = get_grid(GMV)
# computes inversion height as z with max gradient of qt
z∇q_tot = map(real_center_indicies(grid)) do k
(grid.z_half[k], abs(∇c2f(GMV.QT.values, grid, k)))
end
k_star = argmax(last.(z∇q_tot))
return first.(z∇q_tot)[k_star]
end

function compute_pressure_plume_spacing(self::EDMF_PrognosticTKE, GMV::GridMeanVariables, Case::CasesBase)

@inbounds for i in xrange(self.n_updrafts)
Expand Down Expand Up @@ -1785,7 +1742,6 @@ function compute_covariance(self::EDMF_PrognosticTKE, GMV::GridMeanVariables, Ca
self.UpdVar.H,
self.UpdVar.QT,
)
cleanup_covariance(self, GMV)
return
end

Expand Down Expand Up @@ -1872,49 +1828,6 @@ function initialize_covariance(self::EDMF_PrognosticTKE, GMV::GridMeanVariables,
return
end

function cleanup_covariance(self::EDMF_PrognosticTKE, GMV::GridMeanVariables)

tmp_eps = 1e-18
grid = get_grid(self)

@inbounds for k in real_center_indicies(grid)
if GMV.TKE.values[k] < tmp_eps
GMV.TKE.values[k] = 0.0
end
if GMV.Hvar.values[k] < tmp_eps
GMV.Hvar.values[k] = 0.0
end
if GMV.QTvar.values[k] < tmp_eps
GMV.QTvar.values[k] = 0.0
end
if fabs(GMV.HQTcov.values[k]) < tmp_eps
GMV.HQTcov.values[k] = 0.0
end
if fabs(GMV.W_third_m.values[k]) < tmp_eps
GMV.W_third_m.values[k] = 0.0
end
if fabs(GMV.H_third_m.values[k]) < tmp_eps
GMV.H_third_m.values[k] = 0.0
end
if fabs(GMV.QT_third_m.values[k]) < tmp_eps
GMV.QT_third_m.values[k] = 0.0
end
if self.EnvVar.Hvar.values[k] < tmp_eps
self.EnvVar.Hvar.values[k] = 0.0
end
if self.EnvVar.TKE.values[k] < tmp_eps
self.EnvVar.TKE.values[k] = 0.0
end
if self.EnvVar.QTvar.values[k] < tmp_eps
self.EnvVar.QTvar.values[k] = 0.0
end
if fabs(self.EnvVar.HQTcov.values[k]) < tmp_eps
self.EnvVar.HQTcov.values[k] = 0.0
end
end
end


function compute_covariance_shear(
self::EDMF_PrognosticTKE,
GMV::GridMeanVariables,
Expand Down
Loading

0 comments on commit cbd3d5d

Please sign in to comment.