Skip to content

Commit

Permalink
fixes for Float32 compat
Browse files Browse the repository at this point in the history
  • Loading branch information
juliasloan25 committed Dec 7, 2023
1 parent ba4cf2b commit 8fddcec
Show file tree
Hide file tree
Showing 4 changed files with 24 additions and 18 deletions.
10 changes: 5 additions & 5 deletions experiments/AMIP/modular/components/land/bucket_init.jl
Original file line number Diff line number Diff line change
Expand Up @@ -163,19 +163,19 @@ Initializes the bucket model variables.
"""
function bucket_init(
::Type{FT},
tspan::Tuple{FT, FT},
tspan::Tuple{AbstractFloat, AbstractFloat},
config::String,
albedo_type::String,
land_temperature_anomaly::String,
comms_ctx::AbstractCommsContext,
regrid_dirpath::String;
space,
dt::FT,
saveat::FT,
dt::AbstractFloat,
saveat::AbstractFloat,
area_fraction,
stepper = CTS.RK4(),
date_ref::DateTime,
t_start::FT,
t_start::AbstractFloat,
) where {FT}
if config != "sphere"
println(
Expand Down Expand Up @@ -210,7 +210,7 @@ function bucket_init(
z_0b = FT(1e-3)
κ_soil = FT(0.7)
ρc_soil = FT(2e8)
t_crit = dt # This is the timescale on which snow exponentially damps to zero, in the case where all
t_crit = FT(dt) # This is the timescale on which snow exponentially damps to zero, in the case where all
# the snow would melt in time t_crit. It prevents us from having to specially time step in cases where
# all the snow melts in a single timestep.
params = BucketModelParameters(κ_soil, ρc_soil, albedo, σS_c, W_f, z_0m, z_0b, t_crit, earth_param_set)
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -103,13 +103,15 @@ function solve_eisenman_model!(Y, Ya, p, thermo_params, Δt)

# ice thickness and mixed layer temperature changes due to atmosphereic and ocean fluxes
ice_covered = parent(h_ice)[1] > 0

FT = eltype(T_ml)
if ice_covered # ice-covered
F_base = @. C0_base * (T_ml - T_base)
ΔT_ml = @. -(F_base - ocean_qflux) * Δt / (hρc_ml)
Δh_ice = @. (F_atm - F_base - ocean_qflux) * Δt / L_ice
@. e_base .+= F_base * Δt
ΔT_ml = @. -(F_base - ocean_qflux) * FT(Δt) / (hρc_ml)
Δh_ice = @. (F_atm - F_base - ocean_qflux) * FT(Δt) / L_ice
@. e_base .+= F_base * FT(Δt)
else # ice-free
ΔT_ml = @. -(F_atm - ocean_qflux) * Δt / (hρc_ml)
ΔT_ml = @. -(F_atm - ocean_qflux) * FT(Δt) / (hρc_ml)
Δh_ice = 0
end

Expand Down Expand Up @@ -317,11 +319,13 @@ function get_field(sim::EisenmanIceSimulation, ::Val{:energy})
e_base = cache.Ya.e_base
ocean_qflux = cache.Ya.ocean_qflux

FT = eltype(sim.integrator.u)

hρc_ml = p_o.h * p_o.ρ * p_o.c

e_ml = @. p_o.h * p_o.ρ * p_o.c * sim.integrator.u.T_ml # heat
e_ice = @. p_i.L_ice * sim.integrator.u.h_ice # phase
e_qflux = @. ocean_qflux * sim.integrator.t
e_qflux = @. ocean_qflux * FT(sim.integrator.t)

return @. e_ml + e_ice + e_qflux + e_base

Expand Down
7 changes: 4 additions & 3 deletions src/ConservationChecker.jl
Original file line number Diff line number Diff line change
Expand Up @@ -100,10 +100,10 @@ function check_conservation!(
parent(coupler_sim.fields.F_radiative_TOA) .= parent(Interfacer.get_field(sim, Val(:F_radiative_TOA)))

if isempty(ccs.toa_net_source)
radiation_sources_accum = sum(coupler_sim.fields.F_radiative_TOA .* coupler_sim.Δt_cpl) # ∫ J / m^2 dA
radiation_sources_accum = sum(coupler_sim.fields.F_radiative_TOA .* FT(coupler_sim.Δt_cpl)) # ∫ J / m^2 dA
else
radiation_sources_accum =
sum(coupler_sim.fields.F_radiative_TOA .* coupler_sim.Δt_cpl) .+ ccs.toa_net_source[end] # ∫ J / m^2 dA
sum(coupler_sim.fields.F_radiative_TOA .* FT(coupler_sim.Δt_cpl)) .+ ccs.toa_net_source[end] # ∫ J / m^2 dA
end
push!(ccs.toa_net_source, radiation_sources_accum)

Expand Down Expand Up @@ -211,7 +211,8 @@ function surface_water_gain_from_rates(cs::Interfacer.CoupledSimulation)
evaporation = cs.fields.F_turb_moisture # kg / m^2 / s / layer depth
precipitation_l = cs.fields.P_liq
precipitation_s = cs.fields.P_snow
@. -(evaporation + precipitation_l + precipitation_s) * cs.Δt_cpl # kg / m^2 / layer depth
FT = eltype(evaporation)
@. -(evaporation + precipitation_l + precipitation_s) * FT(cs.Δt_cpl) # kg / m^2 / layer depth
end

# setup the GKS socket application environment as nul for better performance and to avoid GKS connection errors while plotting
Expand Down
11 changes: 6 additions & 5 deletions test/component_model_tests/eisenman_seaice_tests.jl
Original file line number Diff line number Diff line change
Expand Up @@ -20,7 +20,7 @@ for FT in (Float32, Float64)
params_ocean = EisenmanOceanParameters{FT}()
params = (; p_i = params_ice, p_o = params_ocean)

Δt = Float64(1e6)
Δt = FT(1e6)

# create a boundary space
boundary_space = TestHelper.create_space(FT)
Expand Down Expand Up @@ -240,6 +240,7 @@ for FT in (Float32, Float64)
F_atm = @. Ya.F_rad + Ya.F_turb
∂F_atm∂T_sfc = get_∂F_rad_energy∂T_sfc(T_sfc_0, params_ice) .+ Ya.∂F_turb_energy∂T_sfc

FT = eltype(F_atm)
@test all(parent(Ya.e_base) .≈ 0) # no contribution from basal fluxes
@test all(parent(Y.T_ml) .≈ params_ice.T_freeze) # ocean temperature remains at the freezing point
h_ice_new =
Expand All @@ -254,7 +255,7 @@ for FT in (Float32, Float64)
# conductive flux and q flux
@testset "Non-zero conductive flux for FT=$FT" begin
Y, Ya = state_init(params_ice, boundary_space)
Δt = Float64(100)
Δt = FT(100)

# ice covered
Y.h_ice .= 10
Expand Down Expand Up @@ -288,7 +289,7 @@ for FT in (Float32, Float64)

@testset "Nonzero Q-flux (~horizontal ocean transport) for FT=$FT" begin
Y, Ya = state_init(params_ice, boundary_space)
Δt = Float64(100)
Δt = FT(100)

# ice covered
Y.h_ice .= 10
Expand Down Expand Up @@ -323,7 +324,7 @@ for FT in (Float32, Float64)

include("../../experiments/AMIP/modular/components/slab_utils.jl")
@testset "step! update + total energy calculation for FT=$FT" begin
Δt = Float64(1000)
Δt = FT(1000)

sim = eisenman_seaice_init(
FT,
Expand All @@ -347,7 +348,7 @@ for FT in (Float32, Float64)
@test all(parent(h_ice) .≈ 0.001)
step!(sim, 2 * Δt)
h_ice = sim.integrator.u.h_ice
@test all(parent(h_ice) . 0.002)
@test all(abs.(parent(h_ice) .- 0.002) .< 10eps(FT))

total_energy_calc = (get_field(sim, Val(:energy)) .- total_energy_0)
total_energy_expeted = 300 .* ones(boundary_space) .* 2 .* Δt
Expand Down

0 comments on commit 8fddcec

Please sign in to comment.