From 8897481ce0ae4f7ec3aaad9d62d97b750230423b Mon Sep 17 00:00:00 2001 From: Gabriele Bozzola Date: Thu, 22 Aug 2024 08:48:21 -0700 Subject: [PATCH] Use ClimaUtilities inputs in run_amip ClimaUtilities comes with a rich array of tools to read input data. Among these are the `SpaceVaryingInput` and the `TimeVaryingInput` to read data from NetCDF files (and more). These functions are more efficient, versatile, and come with various interpolation options. This commit moves AMIP to using those. In the process, it moves from using ClimaCoreTempestRemap to Interpolations. --- experiments/ClimaEarth/Project.toml | 2 +- .../components/ocean/prescr_seaice.jl | 6 - .../ClimaEarth/components/ocean/slab_ocean.jl | 7 -- experiments/ClimaEarth/run_amip.jl | 119 ++++++------------ perf/Manifest.toml | 6 +- perf/Project.toml | 2 + 6 files changed, 46 insertions(+), 96 deletions(-) diff --git a/experiments/ClimaEarth/Project.toml b/experiments/ClimaEarth/Project.toml index e7b8544f0..85e3d925d 100644 --- a/experiments/ClimaEarth/Project.toml +++ b/experiments/ClimaEarth/Project.toml @@ -68,7 +68,7 @@ Colors = "0.12" Dierckx = "0.5" ForwardDiff = "0.10" Glob = "1" -HDF5_jll = "1" +Interpolations = "0.14, 0.15" IntervalSets = "0.5, 0.6, 0.7" JSON = "0.21" MPI = "0.20" diff --git a/experiments/ClimaEarth/components/ocean/prescr_seaice.jl b/experiments/ClimaEarth/components/ocean/prescr_seaice.jl index 4b46bbc93..df08c4721 100644 --- a/experiments/ClimaEarth/components/ocean/prescr_seaice.jl +++ b/experiments/ClimaEarth/components/ocean/prescr_seaice.jl @@ -148,12 +148,6 @@ end ### ### Sea ice model-specific functions (not explicitly required by ClimaCoupler.jl) ### -""" - scale_sic(SIC, _info) -Ensures that the space of the SIC struct matches that of the mask, and converts the units from area % to area fraction. -""" -scale_sic(SIC, _info) = Utilities.swap_space!(axes(_info.land_fraction), SIC) ./ BCReader.float_type_bcf(_info)(100.0) - # setting that SIC < 0.5 is counted as ocean if binary remapping. get_ice_fraction(h_ice::FT, mono::Bool, threshold = 0.5) where {FT} = mono ? h_ice : Regridder.binary_mask(h_ice, threshold) diff --git a/experiments/ClimaEarth/components/ocean/slab_ocean.jl b/experiments/ClimaEarth/components/ocean/slab_ocean.jl index fa1277a8f..2989e2f90 100644 --- a/experiments/ClimaEarth/components/ocean/slab_ocean.jl +++ b/experiments/ClimaEarth/components/ocean/slab_ocean.jl @@ -168,13 +168,6 @@ end ### ### Slab ocean model-specific functions (not explicitly required by ClimaCoupler.jl) ### -""" - scale_sst(SST::FT, _info) -Ensures that the space of the SST struct matches that of the land_fraction, and converts the units to Kelvin (N.B.: this is dataset specific) -""" -scale_sst(SST, _info) = - (Utilities.swap_space!(axes(_info.land_fraction), SST) .+ BCReader.float_type_bcf(_info)(273.15)) - # ode function slab_ocean_rhs!(dY, Y, cache, t) p, F_turb_energy, F_radiative, area_fraction = cache diff --git a/experiments/ClimaEarth/run_amip.jl b/experiments/ClimaEarth/run_amip.jl index d70cb4c4a..ac62a6905 100644 --- a/experiments/ClimaEarth/run_amip.jl +++ b/experiments/ClimaEarth/run_amip.jl @@ -48,7 +48,6 @@ import ClimaCore as CC # ## Coupler specific imports import ClimaCoupler import ClimaCoupler: - BCReader, ConservationChecker, Checkpointer, Diagnostics, @@ -59,6 +58,10 @@ import ClimaCoupler: TimeManager, Utilities +import ClimaUtilities.SpaceVaryingInputs: SpaceVaryingInput +import ClimaUtilities.TimeVaryingInputs: TimeVaryingInput, evaluate! +import Interpolations + pkg_dir = pkgdir(ClimaCoupler) #= @@ -168,12 +171,11 @@ end The data files are downloaded from the `ClimaCoupler` artifacts directory. If the data files are not present, they are downloaded from the original sources. =# - include(joinpath(pkgdir(ClimaCoupler), "artifacts", "artifact_funcs.jl")) -sst_data = artifact_data(sst_dataset_path(), "sst", "SST", dir_paths.regrid, date0, t_start, t_end, comms_ctx) -sic_data = artifact_data(sic_dataset_path(), "sic", "SEAICE", dir_paths.regrid, date0, t_start, t_end, comms_ctx) -co2_data = artifact_data(co2_dataset_path(), "mauna_loa_co2", "co2", dir_paths.regrid, date0, t_start, t_end, comms_ctx) -land_mask_data = artifact_data(mask_dataset_path(), "seamask") +sst_data = joinpath(sst_dataset_path(), "sst.nc") +sic_data = joinpath(sic_dataset_path(), "sic.nc") +co2_data = joinpath(co2_dataset_path(), "mauna_loa_co2.nc") +land_mask_data = joinpath(mask_dataset_path(), "seamask.nc") #= ## Component Model Initialization @@ -215,18 +217,7 @@ Note that land-sea area fraction is different to the land-sea mask, which is a b (masks are used internally by the coupler to indicate passive cells that are not populated by a given component model). =# -land_area_fraction = - FT.( - Regridder.land_fraction( - FT, - dir_paths.regrid, - comms_ctx, - land_mask_data, - "LSMASK", - boundary_space, - mono = mono_surface, - ) - ) +land_area_fraction = SpaceVaryingInput(land_mask_data, "LSMASK", boundary_space) Utilities.show_memory_usage(comms_ctx) #= @@ -268,22 +259,17 @@ if mode_name == "amip" ) ## ocean stub - SST_info = BCReader.bcfile_info_init( - FT, - dir_paths.regrid, + SST_timevaryinginput = TimeVaryingInput( sst_data, "SST", boundary_space, - comms_ctx, - interpolate_daily = true, - scaling_function = scale_sst, ## convert to Kelvin - land_fraction = land_area_fraction, - date0 = date0, - mono = mono_surface, + reference_date = date0, + file_reader_kwargs = (; preprocess_func = (data) -> data + FT(273.15),), ## convert to Kelvin ) - BCReader.update_midmonth_data!(date0, SST_info) - SST_init = BCReader.interpolate_midmonth_to_daily(date0, SST_info) + SST_init = CC.Fields.zeros(boundary_space) + evaluate!(SST_init, SST_timevaryinginput, t_start) + ocean_sim = Interfacer.SurfaceStub((; T_sfc = SST_init, ρ_sfc = CC.Fields.zeros(boundary_space), @@ -298,21 +284,17 @@ if mode_name == "amip" )) ## sea ice model - SIC_info = BCReader.bcfile_info_init( - FT, - dir_paths.regrid, + SIC_timevaryinginput = TimeVaryingInput( sic_data, "SEAICE", boundary_space, - comms_ctx, - interpolate_daily = true, - scaling_function = scale_sic, ## convert to fraction - land_fraction = land_area_fraction, - date0 = date0, - mono = mono_surface, + reference_date = date0, + file_reader_kwargs = (; preprocess_func = (data) -> data / 100,), ## convert to fraction ) - BCReader.update_midmonth_data!(date0, SIC_info) - SIC_init = BCReader.interpolate_midmonth_to_daily(date0, SIC_info) + + SIC_init = CC.Fields.zeros(boundary_space) + evaluate!(SIC_init, SIC_timevaryinginput, t_start) + ice_fraction = get_ice_fraction.(SIC_init, mono_surface) ice_sim = ice_init( FT; @@ -325,24 +307,18 @@ if mode_name == "amip" ) ## CO2 concentration from temporally varying file - CO2_info = BCReader.bcfile_info_init( - FT, - dir_paths.regrid, - co2_data, - "co2", - boundary_space, - comms_ctx, - interpolate_daily = true, - land_fraction = ones(boundary_space), - date0 = date0, - mono = mono_surface, - ) + CO2_timevaryinginput = TimeVaryingInput(co2_data, "co2", boundary_space, reference_date = date0) - BCReader.update_midmonth_data!(date0, CO2_info) - CO2_init = BCReader.interpolate_midmonth_to_daily(date0, CO2_info) - Interfacer.update_field!(atmos_sim, Val(:co2), CO2_init) + CO2_init = CC.Fields.zeros(boundary_space) + evaluate!(CO2_init, CO2_timevaryinginput, t_start) + CO2_field = Interfacer.update_field!(atmos_sim, Val(:co2), CO2_init) - mode_specifics = (; name = mode_name, SST_info = SST_info, SIC_info = SIC_info, CO2_info = CO2_info) + mode_specifics = (; + name = mode_name, + SST_timevaryinginput = SST_timevaryinginput, + SIC_timevaryinginput = SIC_timevaryinginput, + CO2_timevaryinginput = CO2_timevaryinginput, + ) Utilities.show_memory_usage(comms_ctx) elseif mode_name in ("slabplanet", "slabplanet_aqua", "slabplanet_terra") @@ -394,7 +370,7 @@ elseif mode_name in ("slabplanet", "slabplanet_aqua", "slabplanet_terra") thermo_params = thermo_params, )) - mode_specifics = (; name = mode_name, SST_info = nothing, SIC_info = nothing) + mode_specifics = (; name = mode_name, SST_timevaryinginput = nothing, SIC_timevaryinginput = nothing) Utilities.show_memory_usage(comms_ctx) elseif mode_name == "slabplanet_eisenman" @@ -438,7 +414,7 @@ elseif mode_name == "slabplanet_eisenman" thermo_params = thermo_params, ) - mode_specifics = (; name = mode_name, SST_info = nothing, SIC_info = nothing) + mode_specifics = (; name = mode_name, SST_timevaryinginput = nothing, SIC_timevaryinginput = nothing) Utilities.show_memory_usage(comms_ctx) end @@ -712,28 +688,13 @@ function solve_coupler!(cs) if cs.mode.name == "amip" - ## update values of SST, SIC, and CO2 for this timestep - if cs.dates.date[1] >= BCReader.next_date_in_file(cs.mode.SST_info) - BCReader.update_midmonth_data!(cs.dates.date[1], cs.mode.SST_info) - end - SST_current = BCReader.interpolate_midmonth_to_daily(cs.dates.date[1], cs.mode.SST_info) - Interfacer.update_field!(ocean_sim, Val(:surface_temperature), SST_current) + evaluate!(Interfacer.get_field(ocean_sim, Val(:surface_temperature)), cs.mode.SST_timevaryinginput, t) + evaluate!(Interfacer.get_field(ice_sim, Val(:area_fraction)), cs.mode.SIC_timevaryinginput, t) - if cs.dates.date[1] >= BCReader.next_date_in_file(cs.mode.SIC_info) - BCReader.update_midmonth_data!(cs.dates.date[1], cs.mode.SIC_info) - end - SIC_current = - get_ice_fraction.( - BCReader.interpolate_midmonth_to_daily(cs.dates.date[1], cs.mode.SIC_info), - cs.mode.SIC_info.mono, - ) - Interfacer.update_field!(ice_sim, Val(:area_fraction), SIC_current) - - if cs.dates.date[1] >= BCReader.next_date_in_file(cs.mode.CO2_info) - BCReader.update_midmonth_data!(cs.dates.date[1], cs.mode.CO2_info) - end - CO2_current = BCReader.interpolate_midmonth_to_daily(cs.dates.date[1], cs.mode.CO2_info) - Interfacer.update_field!(atmos_sim, Val(:co2), CO2_current) + # TODO: get_field with :co2 is not implemented, so this is a little awkward + current_CO2 = CC.Fields.zeros(boundary_space) + evaluate!(current_CO2, cs.mode.CO2_timevaryinginput, t) + Interfacer.update_field!(atmos_sim, Val(:co2), current_CO2) ## calculate and accumulate diagnostics at each timestep, if we're using diagnostics in this run if !isempty(cs.diagnostics) diff --git a/perf/Manifest.toml b/perf/Manifest.toml index ba6bcb4d9..29b513c34 100644 --- a/perf/Manifest.toml +++ b/perf/Manifest.toml @@ -2,7 +2,7 @@ julia_version = "1.10.4" manifest_format = "2.0" -project_hash = "302422306467d2326750e939162edbd82fdbb61b" +project_hash = "a9717df068fda1bfc383500582ffde650f049077" [[deps.ADTypes]] git-tree-sha1 = "99a6f5d0ce1c7c6afdb759daa30226f71c54f6b0" @@ -353,9 +353,9 @@ version = "0.7.33" [[deps.ClimaUtilities]] deps = ["Artifacts", "Dates"] -git-tree-sha1 = "d572737e623699c37cdc8bb63d789fedc3dd3a95" +git-tree-sha1 = "286cf36603ec68538dda4b0333af694fe70dd628" uuid = "b3f4f4ca-9299-4f7f-bd9b-81e1242a7513" -version = "0.1.13" +version = "0.1.14" weakdeps = ["Adapt", "CUDA", "ClimaComms", "ClimaCore", "ClimaCoreTempestRemap", "Interpolations", "NCDatasets"] [deps.ClimaUtilities.extensions] diff --git a/perf/Project.toml b/perf/Project.toml index 9117033b5..591fdc2b9 100644 --- a/perf/Project.toml +++ b/perf/Project.toml @@ -13,6 +13,7 @@ ClimaCoupler = "4ade58fe-a8da-486c-bd89-46df092ec0c7" ClimaLand = "08f4d4ce-cf43-44bb-ad95-9d2d5f413532" ClimaParams = "5c42b081-d73a-476f-9059-fd94b934656c" ClimaTimeSteppers = "595c0a79-7f3d-439a-bc5a-b232dc3bde79" +ClimaUtilities = "b3f4f4ca-9299-4f7f-bd9b-81e1242a7513" CloudMicrophysics = "6a9e3e04-43cd-43ba-94b9-e8782df3c71b" Colors = "5ae59095-9a9b-59fe-a467-6f913c188581" CommonDataModel = "1fbeeb36-5f17-413c-809b-666fb144f157" @@ -22,6 +23,7 @@ ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210" GR = "28b8d3ca-fb5f-59d9-8090-bfdbd6d07a71" Glob = "c27321d9-0574-5035-807b-f59d2c89b15c" HDF5_jll = "0234f1f7-429e-5d53-9886-15a909be8d59" +Interpolations = "a98d9a8b-a2ab-59e6-89dd-64a1c18fca59" IntervalSets = "8197267c-284f-5f27-9208-e0e47529a953" JLD2 = "033835bb-8acc-5ee8-8aae-3f567f8a3819" JSON = "682c06a0-de6a-54ab-a142-c8b1cf79cde6"