From 78a2ec8ce044aef6a638797e26985b787e6f95aa Mon Sep 17 00:00:00 2001 From: Julia Sloan Date: Tue, 5 Dec 2023 18:01:50 -0800 Subject: [PATCH] fixes for Float32 compat --- .../AMIP/modular/components/land/bucket_init.jl | 10 +++++----- .../components/ocean/eisenman_seaice_init.jl | 14 +++++++++----- .../modular/components/ocean/prescr_seaice_init.jl | 2 +- src/ConservationChecker.jl | 11 ++++++----- src/Diagnostics.jl | 6 +++--- .../component_model_tests/eisenman_seaice_tests.jl | 11 ++++++----- test/diagnostics_tests.jl | 2 +- test/field_exchanger_tests.jl | 6 ++++-- 8 files changed, 35 insertions(+), 27 deletions(-) diff --git a/experiments/AMIP/modular/components/land/bucket_init.jl b/experiments/AMIP/modular/components/land/bucket_init.jl index 370573712a..83dffc75cf 100644 --- a/experiments/AMIP/modular/components/land/bucket_init.jl +++ b/experiments/AMIP/modular/components/land/bucket_init.jl @@ -163,19 +163,19 @@ Initializes the bucket model variables. """ function bucket_init( ::Type{FT}, - tspan::Tuple{FT, FT}, + tspan::Tuple{Float64, Float64}, config::String, albedo_type::String, land_temperature_anomaly::String, comms_ctx::AbstractCommsContext, regrid_dirpath::String; space, - dt::FT, - saveat::FT, + dt::Float64, + saveat::Float64, area_fraction, stepper = CTS.RK4(), date_ref::DateTime, - t_start::FT, + t_start::Float64, ) where {FT} if config != "sphere" println( @@ -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) diff --git a/experiments/AMIP/modular/components/ocean/eisenman_seaice_init.jl b/experiments/AMIP/modular/components/ocean/eisenman_seaice_init.jl index 35488bf2d0..21fd472aa3 100644 --- a/experiments/AMIP/modular/components/ocean/eisenman_seaice_init.jl +++ b/experiments/AMIP/modular/components/ocean/eisenman_seaice_init.jl @@ -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 @@ -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 diff --git a/experiments/AMIP/modular/components/ocean/prescr_seaice_init.jl b/experiments/AMIP/modular/components/ocean/prescr_seaice_init.jl index 7b9d596882..67c60c6969 100644 --- a/experiments/AMIP/modular/components/ocean/prescr_seaice_init.jl +++ b/experiments/AMIP/modular/components/ocean/prescr_seaice_init.jl @@ -79,7 +79,7 @@ function ice_rhs!(du, u, p, _) # do not count tendencies that lead to temperatures above freezing, and mask out no-ice areas area_mask = Regridder.binary_mask.(area_fraction, threshold = eps(FT)) - unphysical = @. Regridder.binary_mask.(T_freeze - (Y.T_sfc + FT(rhs) * p.dt), threshold = FT(0)) .* area_mask + unphysical = @. Regridder.binary_mask.(T_freeze - (Y.T_sfc + FT(rhs) * FT(p.dt)), threshold = FT(0)) .* area_mask parent(dY.T_sfc) .= parent(rhs .* unphysical) @. p.q_sfc = TD.q_vap_saturation_generic.(p.thermo_params, Y.T_sfc, p.ρ_sfc, TD.Ice()) diff --git a/src/ConservationChecker.jl b/src/ConservationChecker.jl index 82fb47d458..55fd8399cd 100644 --- a/src/ConservationChecker.jl +++ b/src/ConservationChecker.jl @@ -100,16 +100,16 @@ 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) # save atmos previous = getproperty(ccs, sim_name) - current = sum(FT.(Interfacer.get_field(sim, Val(:energy)))) # # ∫ J / m^3 dV + current = sum(Interfacer.get_field(sim, Val(:energy))) # # ∫ J / m^3 dV push!(previous, current) total += current + radiation_sources_accum @@ -122,7 +122,7 @@ function check_conservation!( current = FT(0) else previous = getproperty(ccs, sim_name) - current = sum(FT.(Interfacer.get_field(sim, Val(:energy)) .* area_fraction)) # # ∫ J / m^3 dV + current = sum(Interfacer.get_field(sim, Val(:energy)) .* area_fraction) # # ∫ J / m^3 dV end push!(previous, current) total += current @@ -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 diff --git a/src/Diagnostics.jl b/src/Diagnostics.jl index 269e317956..6104ba599b 100644 --- a/src/Diagnostics.jl +++ b/src/Diagnostics.jl @@ -6,7 +6,7 @@ This module contains functions for defining, gathering and outputting online mod module Diagnostics using ClimaCore: Spaces, Fields, InputOutput -using ClimaCoupler.Interfacer: CoupledSimulation +using ClimaCoupler.Interfacer: CoupledSimulation, float_type using Dates using ClimaCoupler.TimeManager: AbstractFrequency, Monthly, EveryTimestep, trigger_callback using ClimaComms @@ -78,12 +78,12 @@ end Collects diagnostics in diags names. """ function collect_diags(cs::CoupledSimulation, dg::DiagnosticsGroup) - + FT = float_type(cs) diags = (;) diag_names = propertynames(dg.field_vector) for name in diag_names - diags = (; diags..., zip((name,), (get_var(cs, Val(name)),))...) + diags = (; diags..., zip((name,), (FT.(get_var(cs, Val(name))),))...) end return Fields.FieldVector(; diags...) diff --git a/test/component_model_tests/eisenman_seaice_tests.jl b/test/component_model_tests/eisenman_seaice_tests.jl index 5bc16cdb25..4e09dbe4ce 100644 --- a/test/component_model_tests/eisenman_seaice_tests.jl +++ b/test/component_model_tests/eisenman_seaice_tests.jl @@ -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) @@ -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 = @@ -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 @@ -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 @@ -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, @@ -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 diff --git a/test/diagnostics_tests.jl b/test/diagnostics_tests.jl index 2ba04fc9c6..0435e77be1 100644 --- a/test/diagnostics_tests.jl +++ b/test/diagnostics_tests.jl @@ -19,7 +19,7 @@ import ClimaCoupler.Diagnostics: post_save, save_time_format -get_var(cs::Interfacer.CoupledSimulation, ::Val{:x}) = float_type(cs)(1) +get_var(cs::Interfacer.CoupledSimulation, ::Val{:x}) = 1 for FT in (Float32, Float64) @testset "init_diagnostics for FT=$FT" begin diff --git a/test/field_exchanger_tests.jl b/test/field_exchanger_tests.jl index 8fcad85cf5..5eb57521c6 100644 --- a/test/field_exchanger_tests.jl +++ b/test/field_exchanger_tests.jl @@ -74,7 +74,10 @@ end struct TestSurfaceSimulationLand{C} <: SurfaceModelSimulation cache::C end -get_field(::TestSurfaceSimulationLand, ::Val{:area_fraction}) = 0.5 +function get_field(sim::TestSurfaceSimulationLand, ::Val{:area_fraction}) + FT = eltype(sim.cache.turbulent_energy_flux) + return FT(0.5) +end function update_field!(sim::TestSurfaceSimulationLand, ::Val{:turbulent_energy_flux}, field) parent(sim.cache.turbulent_energy_flux) .= parent(field) end @@ -110,7 +113,6 @@ for FT in (Float32, Float64) end end - @testset "import_combined_surface_fields! for FT=$FT" begin # coupler cache setup boundary_space = TestHelper.create_space(FT)