Skip to content

Commit

Permalink
Apply CA PR #1453 surface fluxes API changes
Browse files Browse the repository at this point in the history
  • Loading branch information
valeriabarra committed May 5, 2023
1 parent 7183ae7 commit 878c8a4
Show file tree
Hide file tree
Showing 11 changed files with 32 additions and 38 deletions.
6 changes: 2 additions & 4 deletions experiments/AMIP/modular/components/flux_calculator.jl
Original file line number Diff line number Diff line change
Expand Up @@ -54,10 +54,8 @@ function calculate_surface_fluxes_atmos_grid!(integrator, info_sfc)
Fields.bycolumn(axes(Y.c.uₕ)) do colidx
get_surface_fluxes!(Y, p, colidx)
# corrections (accounting for inhomogeneous surfaces)
@. p.dif_flux_energy_bc[colidx] = # todo: get rid - shouldn't make any difference anyway
Geometry.WVector(correct_e_over_ice(p.surface_conditions[colidx], ice_fraction[colidx]))
@. p.dif_flux_ρq_tot_bc[colidx] =
Geometry.WVector(correct_q_over_ice(p.surface_conditions[colidx], ice_fraction[colidx]))
@. p.ρ_dif_flux_q_tot[colidx] = # checking right quantity in ClimaAtmos v0.11.0
-Geometry.WVector(correct_e_over_ice(p.surface_conditions[colidx], ice_fraction[colidx]))

end
end
Expand Down
3 changes: 1 addition & 2 deletions experiments/AMIP/modular/components/push_pull.jl
Original file line number Diff line number Diff line change
Expand Up @@ -6,8 +6,7 @@ updates F_A, F_R, P_liq, and F_E in place based on values used in the atmos_sim
function atmos_push!(cs)
atmos_sim = cs.model_sims.atmos_sim
csf = cs.fields
dummmy_remap!(csf.F_A, .-atmos_sim.integrator.p.dif_flux_energy_bc)
dummmy_remap!(csf.F_E, .-atmos_sim.integrator.p.dif_flux_ρq_tot_bc)
dummmy_remap!(csf.F_E, atmos_sim.integrator.p.ρ_dif_flux_q_tot)
dummmy_remap!(csf.F_R, level(atmos_sim.integrator.p.ᶠradiation_flux, half))
dummmy_remap!(csf.P_liq, atmos_sim.integrator.p.col_integrated_rain .+ atmos_sim.integrator.p.col_integrated_snow)
end
Expand Down
14 changes: 7 additions & 7 deletions experiments/AMIP/moist_mpi_aquaplanet/atmos/coupled_atmos.jl
Original file line number Diff line number Diff line change
Expand Up @@ -3,7 +3,7 @@ using ClimaCore.Geometry: ⊗
function vertical_diffusion_boundary_layer_coupled_tendency!(Yₜ, Y, p, t)
ᶜρ = Y.c.ρ
(; ᶜp, ᶠv_a, ᶠz_a, ᶠK_E) = p # assume ᶜts and ᶜp have been updated
(; dif_flux_energy, dif_flux_ρq_tot, dif_flux_uₕ) = p
(; ρ_dif_flux_h_tot, dif_flux_ρq_tot, ρ_dif_flux_uₕ) = p
ᶠgradᵥ = Operators.GradientC2F() # apply BCs to ᶜdivᵥ, which wraps ᶠgradᵥ

Fields.field_values(ᶠv_a) .= Fields.field_values(Spaces.level(Y.c.uₕ, 1)) .* one.(Fields.field_values(ᶠz_a)) # TODO: fix VIJFH copyto! to remove this
Expand All @@ -13,7 +13,7 @@ function vertical_diffusion_boundary_layer_coupled_tendency!(Yₜ, Y, p, t)
if :ρe in propertynames(Y.c)
ᶜdivᵥ = Operators.DivergenceF2C(
top = Operators.SetValue(Geometry.WVector(FT(0))),
bottom = Operators.SetValue(.-dif_flux_energy),
bottom = Operators.SetValue(.-ρ_dif_flux_h_tot),
)
@. Yₜ.c.ρe += ᶜdivᵥ(ᶠK_E * ᶠinterp(ᶜρ) * ᶠgradᵥ((Y.c.ρe + ᶜp) / ᶜρ))
end
Expand All @@ -32,7 +32,7 @@ function vertical_diffusion_boundary_layer_coupled_tendency!(Yₜ, Y, p, t)
if :uₕ in propertynames(Y.c)
ᶜdivᵥ = Operators.DivergenceF2C(
top = Operators.SetValue(Geometry.Contravariant3Vector(FT(0)) Geometry.Covariant12Vector(FT(0), FT(0))),
bottom = Operators.SetValue(.-dif_flux_uₕ),
bottom = Operators.SetValue(.-ρ_dif_flux_uₕ),
)

@. Yₜ.c.uₕ += ᶜdivᵥ(ᶠK_E * ᶠgradᵥ(Y.c.uₕ))
Expand All @@ -52,8 +52,8 @@ function vertical_diffusion_boundary_layer_coupled_cache(Y; Cd = FT(0.0044), Ch
dif_flux_ρq_tot = Ref(Geometry.WVector(FT(0)))
end

#dif_flux_uₕ = similar(z_bottom, Geometry.Contravariant3Vector{FT}) .⊗ similar(z_bottom, Geometry.Covariant12Vector{FT}) # this breaks
dif_flux_uₕ =
#ρ_dif_flux_uₕ = similar(z_bottom, Geometry.Contravariant3Vector{FT}) .⊗ similar(z_bottom, Geometry.Covariant12Vector{FT}) # this breaks
ρ_dif_flux_uₕ =
Geometry.Contravariant3Vector.(zeros(axes(z_bottom))) .⊗
Geometry.Covariant12Vector.(zeros(axes(z_bottom)), zeros(axes(z_bottom)))

Expand All @@ -75,9 +75,9 @@ function vertical_diffusion_boundary_layer_coupled_cache(Y; Cd = FT(0.0044), Ch
ᶠz_a,
ᶠK_E = similar(Y.f, FT),
flux_coefficients = similar(z_bottom, coef_type),
dif_flux_energy = similar(z_bottom, Geometry.WVector{FT}),
ρ_dif_flux_h_tot = similar(z_bottom, Geometry.WVector{FT}),
dif_flux_ρq_tot,
dif_flux_uₕ,
ρ_dif_flux_uₕ,
Cd,
Ch,
)
Expand Down
6 changes: 3 additions & 3 deletions experiments/AMIP/moist_mpi_aquaplanet/coupler_driver.jl
Original file line number Diff line number Diff line change
Expand Up @@ -48,8 +48,8 @@ walltime = @elapsed for t in (tspan[1]:Δt_cpl:tspan[end])
T_S .= FT(0)
dummmy_remap!(T_S, slab_sim.integrator.u.T_sfc)

#atmos_sim.integrator.p.dif_flux_energy .= ClimaCore.Geometry.WVector.(ClimaCore.Fields.zeros(axes(atmos_sim.integrator.p.dif_flux_energy)))
parent(atmos_sim.integrator.p.dif_flux_energy) .= FT(0)
#atmos_sim.integrator.p.ρ_dif_flux_h_tot .= ClimaCore.Geometry.WVector.(ClimaCore.Fields.zeros(axes(atmos_sim.integrator.p.ρ_dif_flux_h_tot)))
parent(atmos_sim.integrator.p.ρ_dif_flux_h_tot) .= FT(0)
calculate_surface_fluxes_atmos_grid!(atmos_sim.integrator, T_S)

# run
Expand All @@ -58,7 +58,7 @@ walltime = @elapsed for t in (tspan[1]:Δt_cpl:tspan[end])
## Slab
# pre: get accumulated flux from atmos
F_S .= ClimaCore.Fields.zeros(boundary_space)
dummmy_remap!(F_S, atmos_sim.integrator.p.dif_flux_energy)
dummmy_remap!(F_S, atmos_sim.integrator.p.ρ_dif_flux_h_tot)

# save the accumulated flux
slab_F_sfc = slab_sim.integrator.p.F_sfc
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@ calculate_surface_fluxes_atmos_grid!(integrator)
"""
function calculate_surface_fluxes_atmos_grid!(integrator, T_sfc)
p = integrator.p
(; ᶜts, dif_flux_energy, dif_flux_ρq_tot, dif_flux_uₕ, params, Cd, Ch) = p
(; ᶜts, ρ_dif_flux_h_tot, dif_flux_ρq_tot, ρ_dif_flux_uₕ, params, Cd, Ch) = p

Y = integrator.u

Expand All @@ -32,7 +32,7 @@ function calculate_surface_fluxes_atmos_grid!(integrator, T_sfc)

# Total energy flux
if :ρe in propertynames(Y.c)
@. dif_flux_energy = Geometry.WVector(tsf.shf + tsf.shf + Rn)
@. ρ_dif_flux_h_tot = Geometry.WVector(tsf.shf + tsf.shf + Rn)
end

# Moisture mass flux
Expand All @@ -46,7 +46,7 @@ function calculate_surface_fluxes_atmos_grid!(integrator, T_sfc)
normal = Geometry.WVector.(ones(u_space)) # TODO: this will need to change for topography
ρ_1 = Fields.Field(Fields.field_values(Fields.level(Y.c.ρ, 1)), u_space) # TODO: delete when "space not the same instance" error is dealt with
if :uₕ in propertynames(Y.c)
parent(dif_flux_uₕ) .= # TODO: remove parent when "space not the same instance" error is dealt with
parent(ρ_dif_flux_uₕ) .= # TODO: remove parent when "space not the same instance" error is dealt with
parent(
Geometry.Contravariant3Vector.(normal) .⊗
Geometry.Covariant12Vector.(Geometry.UVVector.(tsf.ρτxz ./ ρ_1, tsf.ρτyz ./ ρ_1)),
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -109,8 +109,8 @@ function ke_dissipation(sim)
drag_uv =
.-Geometry.UVVector.(
Geometry.Covariant12Vector.(
sim.integrator.p.dif_flux_uₕ_bc.components.data.:1,
sim.integrator.p.dif_flux_uₕ_bc.components.data.:2,
sim.integrator.p.ρ_dif_flux_uₕ.components.data.:1,
sim.integrator.p.ρ_dif_flux_uₕ.components.data.:2,
)
)
dot.(Geometry.UVVector.(level(sim.integrator.u.c.uₕ, 1)), drag_uv) .* level(sim.integrator.u.c.ρ, 1)
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -53,10 +53,8 @@ function calculate_surface_fluxes_atmos_grid!(integrator, info_sfc)
Fields.bycolumn(axes(Y.c.uₕ)) do colidx
get_surface_fluxes!(Y, p, colidx)
# corrections (accounting for inhomogeneous surfaces)
@. p.dif_flux_energy_bc[colidx] = # todo: get rid - shouldn't make any difference anyway
Geometry.WVector(correct_e_over_ice(p.surface_conditions[colidx], ice_mask[colidx]))
@. p.dif_flux_ρq_tot_bc[colidx] =
Geometry.WVector(correct_q_over_ice(p.surface_conditions[colidx], ice_mask[colidx]))
@. p.ρ_dif_flux_q_tot[colidx] =
-Geometry.WVector(correct_q_over_ice(p.surface_conditions[colidx], ice_mask[colidx]))

end
end
Expand Down
3 changes: 1 addition & 2 deletions experiments/AMIP/moist_mpi_earth/push_pull.jl
Original file line number Diff line number Diff line change
Expand Up @@ -6,8 +6,7 @@ updates F_A, F_R, P_liq, and F_E in place based on values used in the atmos_sim
function atmos_push!(cs)
atmos_sim = cs.model_sims.atmos_sim
csf = cs.fields
dummmy_remap!(csf.F_A, .-atmos_sim.integrator.p.dif_flux_energy_bc)
dummmy_remap!(csf.F_E, .-atmos_sim.integrator.p.dif_flux_ρq_tot_bc)
dummmy_remap!(csf.F_E, atmos_sim.integrator.p.ρ_dif_flux_q_tot)
dummmy_remap!(csf.F_R, level(atmos_sim.integrator.p.ᶠradiation_flux, half))
dummmy_remap!(csf.P_liq, atmos_sim.integrator.p.col_integrated_rain .+ atmos_sim.integrator.p.col_integrated_snow)
end
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,7 @@ using ClimaCore.Geometry: ⊗
function vertical_diffusion_boundary_layer_coupled_tendency!(Yₜ, Y, p, t)
ᶜρ = Y.c.ρ
(; ᶜp, ᶠv_a, ᶠz_a, ᶠK_E) = p # assume ᶜts and ᶜp have been updated
(; dif_flux_energy, dif_flux_ρq_tot, dif_flux_uₕ) = p
(; ρ_dif_flux_h_tot, dif_flux_ρq_tot, ρ_dif_flux_uₕ) = p
ᶠgradᵥ = Operators.GradientC2F() # apply BCs to ᶜdivᵥ, which wraps ᶠgradᵥ

Fields.field_values(ᶠv_a) .= Fields.field_values(Spaces.level(Y.c.uₕ, 1)) .* one.(Fields.field_values(ᶠz_a)) # TODO: fix VIJFH copyto! to remove this
Expand All @@ -14,7 +14,7 @@ function vertical_diffusion_boundary_layer_coupled_tendency!(Yₜ, Y, p, t)
if :ρe in propertynames(Y.c)
ᶜdivᵥ = Operators.DivergenceF2C(
top = Operators.SetValue(Geometry.WVector(FT(0))),
bottom = Operators.SetValue(.-dif_flux_energy),
bottom = Operators.SetValue(.-ρ_dif_flux_h_tot),
)
@. Yₜ.c.ρe += ᶜdivᵥ(ᶠK_E * ᶠinterp(ᶜρ) * ᶠgradᵥ((Y.c.ρe + ᶜp) / ᶜρ))
end
Expand All @@ -33,7 +33,7 @@ function vertical_diffusion_boundary_layer_coupled_tendency!(Yₜ, Y, p, t)
if :uₕ in propertynames(Y.c)
ᶜdivᵥ = Operators.DivergenceF2C(
top = Operators.SetValue(Geometry.Contravariant3Vector(FT(0)) Geometry.Covariant12Vector(FT(0), FT(0))),
bottom = Operators.SetValue(.-dif_flux_uₕ),
bottom = Operators.SetValue(.-ρ_dif_flux_uₕ),
)

@. Yₜ.c.uₕ += ᶜdivᵥ(ᶠK_E * ᶠgradᵥ(Y.c.uₕ))
Expand All @@ -53,7 +53,7 @@ function vertical_diffusion_boundary_layer_coupled_cache(Y; Cd = FT(0.0014), Ch
dif_flux_ρq_tot = Ref(Geometry.WVector(FT(0)))
end

dif_flux_uₕ =
ρ_dif_flux_uₕ =
Geometry.Contravariant3Vector.(zeros(axes(z_bottom))) .⊗
Geometry.Covariant12Vector.(zeros(axes(z_bottom)), zeros(axes(z_bottom)))

Expand All @@ -75,9 +75,9 @@ function vertical_diffusion_boundary_layer_coupled_cache(Y; Cd = FT(0.0014), Ch
ᶠz_a,
ᶠK_E = similar(Y.f, FT),
flux_coefficients = similar(z_bottom, coef_type),
dif_flux_energy = similar(z_bottom, Geometry.WVector{FT}),
ρ_dif_flux_h_tot = similar(z_bottom, Geometry.WVector{FT}),
dif_flux_ρq_tot,
dif_flux_uₕ,
ρ_dif_flux_uₕ,
∂F_aero∂T_sfc = zeros(axes(z_bottom)),
Cd,
Ch,
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -188,7 +188,7 @@ walltime = @elapsed for t in (tspan[1]:Δt_cpl:tspan[end])

# coupler_push!: get accumulated fluxes from atmos in the surface fields
F_A .= ClimaCore.Fields.zeros(boundary_space)
dummmy_remap!(F_A, atmos_sim.integrator.p.dif_flux_energy)
dummmy_remap!(F_A, atmos_sim.integrator.p.ρ_dif_flux_h_tot)
F_R .= ClimaCore.Fields.zeros(boundary_space)
parsed_args["rad"] == "gray" ? dummmy_remap!(F_R, level(atmos_sim.integrator.p.ᶠradiation_flux, half)) : nothing # TODO: albedo hard coded...
dF_A .= ClimaCore.Fields.zeros(boundary_space)
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -8,7 +8,7 @@ calculate_surface_fluxes_atmos_grid!(integrator)
"""
function calculate_surface_fluxes_atmos_grid!(integrator, info_sfc)
p = integrator.p
(; ᶜts, dif_flux_energy, dif_flux_ρq_tot, dif_flux_uₕ, ∂F_aero∂T_sfc, params, Cd, Ch) = p
(; ᶜts, ρ_dif_flux_h_tot, dif_flux_ρq_tot, ρ_dif_flux_uₕ, ∂F_aero∂T_sfc, params, Cd, Ch) = p

(; T_sfc, z0m, z0b, ice_mask) = info_sfc
Y = integrator.u
Expand All @@ -31,9 +31,9 @@ function calculate_surface_fluxes_atmos_grid!(integrator, info_sfc)
# Total energy flux
if :ρe in propertynames(Y.c)

flux_energy = ones(axes(dif_flux_energy))
flux_energy = ones(axes(ρ_dif_flux_h_tot))
parent(flux_energy) .= parent(tsf.shf .+ tsf.lhf .* swap_space!(zeros(axes(tsf.shf)), abs.(ice_mask .- FT(1)))) # only SHF above sea ice
@. dif_flux_energy = Geometry.WVector(flux_energy) #Geometry.WVector.(swap_space!(zeros(axes(dif_flux_energy)), tsf.shf .+ tsf.lhf) )
@. ρ_dif_flux_h_tot = Geometry.WVector(flux_energy) #Geometry.WVector.(swap_space!(zeros(axes(ρ_dif_flux_h_tot)), tsf.shf .+ tsf.lhf) )
end

# Moisture mass flux
Expand All @@ -48,7 +48,7 @@ function calculate_surface_fluxes_atmos_grid!(integrator, info_sfc)
normal = Geometry.WVector.(ones(u_space)) # TODO: this will need to change for topography
ρ_1 = Fields.Field(Fields.field_values(Fields.level(Y.c.ρ, 1)), u_space) # TODO: delete when "space not the same instance" error is dealt with
if :uₕ in propertynames(Y.c)
parent(dif_flux_uₕ) .= # TODO: remove parent when "space not the same instance" error is dealt with
parent(ρ_dif_flux_uₕ) .= # TODO: remove parent when "space not the same instance" error is dealt with
parent(
Geometry.Contravariant3Vector.(normal) .⊗
Geometry.Covariant12Vector.(Geometry.UVVector.(tsf.ρτxz ./ ρ_1, tsf.ρτyz ./ ρ_1)),
Expand Down

0 comments on commit 878c8a4

Please sign in to comment.