Skip to content

Commit

Permalink
Use more ClimaCore operators, rm upwind gradient
Browse files Browse the repository at this point in the history
  • Loading branch information
@charleskawczynski
charleskawczynski committed Dec 3, 2021
1 parent 24f361f commit edd20ad
Show file tree
Hide file tree
Showing 5 changed files with 73 additions and 79 deletions.
4 changes: 4 additions & 0 deletions driver/main.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -175,7 +175,11 @@ cent_aux_vars_edmf(FT, n_up) = (;
KM = FT(0),
KH = FT(0),
mixing_length = FT(0),
m_entr_detr = FT(0),
∇m_entr_detr = FT(0),
θ_virt = FT(0),
∇θ_virt = FT(0),
∇Ri_bulk = FT(0),
massflux_tendency_h = FT(0),
massflux_tendency_qt = FT(0),
diffusive_tendency_h = FT(0),
Expand Down
84 changes: 42 additions & 42 deletions post_processing/mse_tables.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -36,14 +36,14 @@ all_best_mse["Bomex"]["Hvar_mean"] = 12489.582832388645
all_best_mse["Bomex"]["QTvar_mean"] = 4273.993352259289
#
all_best_mse["DryBubble"] = OrderedCollections.OrderedDict()
all_best_mse["DryBubble"]["updraft_area"] = 3.637544400390376e-22
all_best_mse["DryBubble"]["updraft_w"] = 1.0527021464500098e-21
all_best_mse["DryBubble"]["updraft_thetal"] = 5.075669397381562e-29
all_best_mse["DryBubble"]["updraft_area"] = 7.521107927385887e-22
all_best_mse["DryBubble"]["updraft_w"] = 1.97248709573443e-21
all_best_mse["DryBubble"]["updraft_thetal"] = 5.071999354433348e-29
all_best_mse["DryBubble"]["u_mean"] = 0.0
all_best_mse["DryBubble"]["tke_mean"] = 3.851860429936621e-21
all_best_mse["DryBubble"]["temperature_mean"] = 1.3194250993279412e-29
all_best_mse["DryBubble"]["thetal_mean"] = 1.0194700001151513e-29
all_best_mse["DryBubble"]["Hvar_mean"] = 8.396661117492861e-22
all_best_mse["DryBubble"]["tke_mean"] = 3.7569243240496424e-21
all_best_mse["DryBubble"]["temperature_mean"] = 1.8085778190787878e-29
all_best_mse["DryBubble"]["thetal_mean"] = 1.6261264438456464e-29
all_best_mse["DryBubble"]["Hvar_mean"] = 1.2744944679485061e-21
#
all_best_mse["DYCOMS_RF01"] = OrderedCollections.OrderedDict()
all_best_mse["DYCOMS_RF01"]["qt_mean"] = 0.02321251933702707
Expand All @@ -61,30 +61,30 @@ all_best_mse["DYCOMS_RF01"]["Hvar_mean"] = 1288.1655959164523
all_best_mse["DYCOMS_RF01"]["QTvar_mean"] = 518.9269070601163
#
all_best_mse["GABLS"] = OrderedCollections.OrderedDict()
all_best_mse["GABLS"]["updraft_thetal"] = 1.83662077593314e-31
all_best_mse["GABLS"]["v_mean"] = 7.551128101422123e-28
all_best_mse["GABLS"]["u_mean"] = 4.114057070189765e-31
all_best_mse["GABLS"]["tke_mean"] = 7.493040224414088e-28
all_best_mse["GABLS"]["updraft_thetal"] = 0.0
all_best_mse["GABLS"]["v_mean"] = 0.0
all_best_mse["GABLS"]["u_mean"] = 0.0
all_best_mse["GABLS"]["tke_mean"] = 0.0
all_best_mse["GABLS"]["temperature_mean"] = 0.0
all_best_mse["GABLS"]["thetal_mean"] = 1.8293883161395876e-31
all_best_mse["GABLS"]["Hvar_mean"] = 2.4092979198976328e-27
all_best_mse["GABLS"]["QTvar_mean"] = 1.9309281437739056e-28
all_best_mse["GABLS"]["qt_mean"] = 1.3642324482912016e-28
all_best_mse["GABLS"]["thetal_mean"] = 0.0
all_best_mse["GABLS"]["Hvar_mean"] = 0.0
all_best_mse["GABLS"]["QTvar_mean"] = 0.0
all_best_mse["GABLS"]["qt_mean"] = 0.0
#
all_best_mse["life_cycle_Tan2018"] = OrderedCollections.OrderedDict()
all_best_mse["life_cycle_Tan2018"]["qt_mean"] = 1.8618303434577105e-6
all_best_mse["life_cycle_Tan2018"]["ql_mean"] = 0.004385094517276016
all_best_mse["life_cycle_Tan2018"]["updraft_area"] = 0.000624106918646745
all_best_mse["life_cycle_Tan2018"]["updraft_w"] = 0.0012219001258406016
all_best_mse["life_cycle_Tan2018"]["updraft_qt"] = 1.973756649872312e-5
all_best_mse["life_cycle_Tan2018"]["updraft_thetal"] = 1.971334038041648e-8
all_best_mse["life_cycle_Tan2018"]["v_mean"] = 9.344252895542346e-6
all_best_mse["life_cycle_Tan2018"]["u_mean"] = 3.931521320501898e-8
all_best_mse["life_cycle_Tan2018"]["tke_mean"] = 3.7187564058596185e-5
all_best_mse["life_cycle_Tan2018"]["temperature_mean"] = 9.83831506541318e-10
all_best_mse["life_cycle_Tan2018"]["thetal_mean"] = 9.689212957534404e-10
all_best_mse["life_cycle_Tan2018"]["Hvar_mean"] = 910.8153523157953
all_best_mse["life_cycle_Tan2018"]["QTvar_mean"] = 61.57011640376347
all_best_mse["life_cycle_Tan2018"]["qt_mean"] = 0.0
all_best_mse["life_cycle_Tan2018"]["ql_mean"] = 0.0
all_best_mse["life_cycle_Tan2018"]["updraft_area"] = 0.0
all_best_mse["life_cycle_Tan2018"]["updraft_w"] = 0.0
all_best_mse["life_cycle_Tan2018"]["updraft_qt"] = 0.0
all_best_mse["life_cycle_Tan2018"]["updraft_thetal"] = 0.0
all_best_mse["life_cycle_Tan2018"]["v_mean"] = 0.0
all_best_mse["life_cycle_Tan2018"]["u_mean"] = 0.0
all_best_mse["life_cycle_Tan2018"]["tke_mean"] = 0.0
all_best_mse["life_cycle_Tan2018"]["temperature_mean"] = 0.0
all_best_mse["life_cycle_Tan2018"]["thetal_mean"] = 0.0
all_best_mse["life_cycle_Tan2018"]["Hvar_mean"] = 0.0
all_best_mse["life_cycle_Tan2018"]["QTvar_mean"] = 0.0
#
all_best_mse["Nieuwstadt"] = OrderedCollections.OrderedDict()
all_best_mse["Nieuwstadt"]["updraft_area"] = 98.6759339416021
Expand All @@ -97,21 +97,21 @@ all_best_mse["Nieuwstadt"]["thetal_mean"] = 9.91803796793435e-6
all_best_mse["Nieuwstadt"]["Hvar_mean"] = 1280.1426243926762
#
all_best_mse["Rico"] = OrderedCollections.OrderedDict()
all_best_mse["Rico"]["qt_mean"] = 1.4535805405900761
all_best_mse["Rico"]["updraft_area"] = 476.2072535886542
all_best_mse["Rico"]["updraft_w"] = 104.59784031707035
all_best_mse["Rico"]["updraft_qt"] = 12.155251102710238
all_best_mse["Rico"]["updraft_thetal"] = 133.7825305740903
all_best_mse["Rico"]["v_mean"] = 0.6922243424136714
all_best_mse["Rico"]["u_mean"] = 0.6135683934962348
all_best_mse["Rico"]["tke_mean"] = 89.16536703667661
all_best_mse["Rico"]["temperature_mean"] = 0.000654009537975488
all_best_mse["Rico"]["ql_mean"] = 68.95531522977028
all_best_mse["Rico"]["qt_mean"] = 1.4426449177460303
all_best_mse["Rico"]["updraft_area"] = 476.2564502191665
all_best_mse["Rico"]["updraft_w"] = 105.30176691450832
all_best_mse["Rico"]["updraft_qt"] = 12.200860942996618
all_best_mse["Rico"]["updraft_thetal"] = 133.78341019232744
all_best_mse["Rico"]["v_mean"] = 0.6924082806112366
all_best_mse["Rico"]["u_mean"] = 0.6139068278917045
all_best_mse["Rico"]["tke_mean"] = 88.82993729340717
all_best_mse["Rico"]["temperature_mean"] = 0.00065059118149235
all_best_mse["Rico"]["ql_mean"] = 69.7760396727162
all_best_mse["Rico"]["qi_mean"] = "NA"
all_best_mse["Rico"]["qr_mean"] = 765.9380058169802
all_best_mse["Rico"]["thetal_mean"] = 0.0006437287942378575
all_best_mse["Rico"]["Hvar_mean"] = 17490.272345608864
all_best_mse["Rico"]["QTvar_mean"] = 3529.0795634065917
all_best_mse["Rico"]["qr_mean"] = 765.859308043914
all_best_mse["Rico"]["thetal_mean"] = 0.0006404563982366288
all_best_mse["Rico"]["Hvar_mean"] = 20191.942729431248
all_best_mse["Rico"]["QTvar_mean"] = 4032.6331903358705
#
all_best_mse["Soares"] = OrderedCollections.OrderedDict()
all_best_mse["Soares"]["qt_mean"] = 0.12588586401602345
Expand Down
18 changes: 0 additions & 18 deletions src/Operators.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -17,24 +17,6 @@ end
struct Extrapolate <: AbstractBC end
struct FreeBoundary <: AbstractBC end # when no BC is used (one-sided derivative at surface that takes first and second interior points)

c∇_upwind(f, grid::Grid, k; bottom = NoBCGivenError(), top = NoBCGivenError()) =
c∇_upwind(ccut_upwind(f, grid, k), grid, k; bottom, top)

function c∇_upwind(f_dual::SA.SVector, grid::Grid, k; bottom = NoBCGivenError(), top = NoBCGivenError())
if is_surface_center(grid, k)
return c∇_upwind(f_dual, grid, BottomBCTag(), bottom)
elseif is_toa_center(grid, k)
return c∇_upwind(f_dual, grid, TopBCTag(), top)
else
return c∇_upwind(f_dual, grid, InteriorTag())
end
end
c∇_upwind(f::SA.SVector, grid::Grid, ::InteriorTag) = (f[2] - f[1]) * grid.Δzi
c∇_upwind(f::SA.SVector, grid::Grid, ::TopBCTag, bc::FreeBoundary) = (f[2] - f[1]) * grid.Δzi
c∇_upwind(f::SA.SVector, grid::Grid, ::TopBCTag, bc::SetGradient) = bc.value
c∇_upwind(f::SA.SVector, grid::Grid, ::BottomBCTag, bc::SetValue) = (f[1] - bc.value) * (grid.Δzi * 2)
c∇_upwind(f::SA.SVector, grid::Grid, ::BottomBCTag, bc::SetGradient) = bc.value

# Used when traversing cell faces

interpc2f(f, grid::Grid, k::CCO.PlusHalf; bottom = NoBCGivenError(), top = NoBCGivenError()) =
Expand Down
14 changes: 10 additions & 4 deletions src/turbulence_functions.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -16,15 +16,20 @@ function get_inversion(grid::Grid, state, param_set, Ri_bulk_crit)
u = center_prog_grid_mean(state).u
v = center_prog_grid_mean(state).v
Ri_bulk = center_aux_bulk(state).Ri
∇θ_virt = center_aux_turbconv(state).∇θ_virt
∇Ri_bulk = center_aux_turbconv(state).∇Ri_bulk
θ_virt_b = θ_virt[kc_surf]
z_c = grid.zc
∇c = CCO.DivergenceF2C()
wvec = CC.Geometry.WVector

# test if we need to look at the free convective limit
if (u[kc_surf]^2 + v[kc_surf]^2) <= 0.01
kmask = map(k -> (k, θ_virt[k] > θ_virt_b), real_center_indices(grid))
k_star = first(kmask[findlast(km -> km[2], kmask)])
∇θ_virt = c∇_upwind(θ_virt, grid, k_star; bottom = SetGradient(0), top = FreeBoundary())
h = (θ_virt_b - θ_virt[k_star - 1]) / ∇θ_virt + z_c[k_star - 1]
LB = CCO.LeftBiasedC2F(; bottom = CCO.SetValue(θ_virt[kc_surf]))
@. ∇θ_virt = ∇c(wvec(LB(θ_virt)))
h = (θ_virt_b - θ_virt[k_star - 1]) / ∇θ_virt[k_star] + z_c[k_star - 1]
else
Ri_bulk_fn(k) = g * (θ_virt[k] - θ_virt_b) * z_c[k] / θ_virt_b / (u[k] * u[k] + v[k] * v[k])

Expand All @@ -33,8 +38,9 @@ function get_inversion(grid::Grid, state, param_set, Ri_bulk_crit)
end
kmask = map(k -> (k, Ri_bulk_fn(k) > Ri_bulk_crit), real_center_indices(grid))
k_star = first(kmask[findlast(km -> km[2], kmask)])
∇Ri_bulk = c∇_upwind(Ri_bulk, grid, k_star; bottom = SetGradient(0), top = FreeBoundary())
h = (Ri_bulk_crit - Ri_bulk[k_star - 1]) / ∇Ri_bulk + z_c[k_star - 1]
LB = CCO.LeftBiasedC2F(; bottom = CCO.SetValue(Ri_bulk[kc_surf]))
@. ∇Ri_bulk = ∇c(wvec(LB(Ri_bulk)))
h = (Ri_bulk_crit - Ri_bulk[k_star - 1]) / ∇Ri_bulk[k_star] + z_c[k_star - 1]
end

return h
Expand Down
32 changes: 17 additions & 15 deletions src/update_aux.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -35,6 +35,9 @@ function update_aux!(edmf::EDMF_PrognosticTKE{N_up}, gm, grid, state, Case, para
prog_up_f = face_prog_updrafts(state)
aux_en_unsat = aux_en.unsat
aux_en_sat = aux_en.sat
m_entr_detr = aux_tc.m_entr_detr
∇m_entr_detr = aux_tc.∇m_entr_detr
wvec = CC.Geometry.WVector

#####
##### center variables
Expand Down Expand Up @@ -238,20 +241,20 @@ function update_aux!(edmf::EDMF_PrognosticTKE{N_up}, gm, grid, state, Case, para
##### compute_updraft_closures
#####

@inbounds for k in real_center_indices(grid)
@inbounds for i in 1:N_up
Ic = CCO.InterpolateF2C()
∇c = CCO.DivergenceF2C()
LB = CCO.LeftBiasedC2F(; bottom = CCO.SetValue(FT(0)))

@inbounds for i in 1:N_up
# compute ∇m at cell centers
a_up = aux_up[i].area
w_up = aux_up_f[i].w
w_gm = prog_gm_f.w
@. m_entr_detr = a_up * (Ic(w_up) - Ic(w_gm))
@. ∇m_entr_detr = ∇c(wvec(LB(m_entr_detr)))
@inbounds for k in real_center_indices(grid)
# entrainment
if aux_up[i].area[k] > 0.0
# compute ∇m at cell centers
a_up_c = aux_up[i].area[k]
w_up_c = interpf2c(aux_up_f[i].w, grid, k)
w_gm_c = interpf2c(prog_gm_f.w, grid, k)
m = a_up_c * (w_up_c - w_gm_c)
a_up_cut = ccut_upwind(aux_up[i].area, grid, k)
w_up_cut = daul_f2c_upwind(aux_up_f[i].w, grid, k)
w_gm_cut = daul_f2c_upwind(prog_gm_f.w, grid, k)
m_cut = a_up_cut .* (w_up_cut .- w_gm_cut)
∇m = FT(c∇_upwind(m_cut, grid, k; bottom = SetValue(0), top = FreeBoundary()))

w_min = 0.001

Expand All @@ -267,8 +270,8 @@ function update_aux!(edmf::EDMF_PrognosticTKE{N_up}, gm, grid, state, Case, para
b_up = aux_up[i].buoy[k],
b_en = aux_en.buoy[k],
tke = aux_en.tke[k],
dMdz = ∇m,
M = m,
dMdz = ∇m_entr_detr[k],
M = m_entr_detr[k],
a_up = aux_up[i].area[k],
a_en = aux_en.area[k],
R_up = edmf.pressure_plume_spacing[i],
Expand Down Expand Up @@ -298,7 +301,6 @@ function update_aux!(edmf::EDMF_PrognosticTKE{N_up}, gm, grid, state, Case, para
end
end

wvec = CC.Geometry.WVector
@inbounds for i in 1:N_up

# pressure
Expand Down

0 comments on commit edd20ad

Please sign in to comment.