From 1b79c172f1f8c2b0cdbb8d57d1c82a3c455781c2 Mon Sep 17 00:00:00 2001 From: akshaysridhar Date: Wed, 6 Mar 2024 17:42:22 -0800 Subject: [PATCH 1/3] new file: coupler_component_init.jl modified: coupler_driver_calibration.jl new file: coupler_driver_init.jl new file: coupler_parse_args.jl modified: ../src/Regridder.jl --- calibration/coupler_component_init.jl | 257 +++++++++++++++ calibration/coupler_driver_calibration.jl | 385 +--------------------- calibration/coupler_driver_init.jl | 77 +++++ calibration/coupler_parse_args.jl | 47 +++ src/Regridder.jl | 2 +- 5 files changed, 385 insertions(+), 383 deletions(-) create mode 100644 calibration/coupler_component_init.jl create mode 100644 calibration/coupler_driver_init.jl create mode 100644 calibration/coupler_parse_args.jl diff --git a/calibration/coupler_component_init.jl b/calibration/coupler_component_init.jl new file mode 100644 index 0000000000..c90e8625a7 --- /dev/null +++ b/calibration/coupler_component_init.jl @@ -0,0 +1,257 @@ + +# get the paths to the necessary data files: land-sea mask, sst map, sea ice concentration +include(joinpath(pkgdir(ClimaCoupler), "artifacts", "artifact_funcs.jl")) +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") + +atmos_sim = atmos_init(FT, config_dict_atmos); +thermo_params = get_thermo_params(atmos_sim) # TODO: this should be shared by all models + +#= +We use a common `Space` for all global surfaces. This enables the MPI processes to operate on the same columns in both +the atmospheric and surface components, so exchanges are parallelized. Note this is only possible when the +atmosphere and surface are of the same horizontal resolution. +=# +## init a 2D boundary space at the surface +boundary_space = Spaces.horizontal_space(atmos_sim.domain.face_space) + +# init land-sea fraction +land_fraction = + FT.( + Regridder.land_fraction( + FT, + REGRID_DIR, + comms_ctx, + land_mask_data, + "LSMASK", + boundary_space, + mono = mono_surface, + ) + ) + +@info mode_name +if mode_name == "amip" + @info "AMIP boundary conditions - do not expect energy conservation" + + ## land + land_sim = bucket_init( + FT, + tspan, + config_dict["land_domain_type"], + config_dict["land_albedo_type"], + config_dict["land_temperature_anomaly"], + comms_ctx, + REGRID_DIR; + dt = Δt_cpl, + space = boundary_space, + saveat = saveat, + area_fraction = land_fraction, + date_ref = date0, + t_start = t_start, + ) + + ## ocean + SST_info = bcfile_info_init( + FT, + REGRID_DIR, + sst_data, + "SST", + boundary_space, + comms_ctx, + interpolate_daily = true, + scaling_function = clean_sst, ## convert to Kelvin + land_fraction = land_fraction, + date0 = date0, + mono = mono_surface, + ) + + update_midmonth_data!(date0, SST_info) + SST_init = interpolate_midmonth_to_daily(date0, SST_info) + ocean_sim = SurfaceStub((; + T_sfc = SST_init, + ρ_sfc = ClimaCore.Fields.zeros(boundary_space), + z0m = FT(1e-3), + z0b = FT(1e-3), + beta = FT(1), + α = FT(0.06), + area_fraction = (FT(1) .- land_fraction), + phase = TD.Liquid(), + thermo_params = thermo_params, + )) + + ## sea ice + SIC_info = bcfile_info_init( + FT, + REGRID_DIR, + sic_data, + "SEAICE", + boundary_space, + comms_ctx, + interpolate_daily = true, + scaling_function = clean_sic, ## convert to fraction + land_fraction = land_fraction, + date0 = date0, + mono = mono_surface, + ) + update_midmonth_data!(date0, SIC_info) + SIC_init = interpolate_midmonth_to_daily(date0, SIC_info) + ice_fraction = get_ice_fraction.(SIC_init, mono_surface) + ice_sim = ice_init( + FT; + tspan = tspan, + dt = Δt_cpl, + space = boundary_space, + saveat = saveat, + area_fraction = ice_fraction, + thermo_params = thermo_params, + ) + + ## CO2 concentration + CO2_info = bcfile_info_init( + FT, + REGRID_DIR, + co2_data, + "co2", + boundary_space, + comms_ctx, + interpolate_daily = true, + land_fraction = ones(boundary_space), + date0 = date0, + mono = mono_surface, + ) + + update_midmonth_data!(date0, CO2_info) + CO2_init = interpolate_midmonth_to_daily(date0, CO2_info) + update_field!(atmos_sim, Val(:co2_gm), CO2_init) + + mode_specifics = (; name = mode_name, SST_info = SST_info, SIC_info = SIC_info, CO2_info = CO2_info) + +elseif mode_name in ("slabplanet", "slabplanet_aqua", "slabplanet_terra") + + land_fraction = mode_name == "slabplanet_aqua" ? land_fraction .* 0 : land_fraction + land_fraction = mode_name == "slabplanet_terra" ? land_fraction .* 0 .+ 1 : land_fraction + + ## land + land_sim = bucket_init( + FT, + tspan, + config_dict["land_domain_type"], + config_dict["land_albedo_type"], + config_dict["land_temperature_anomaly"], + comms_ctx, + REGRID_DIR; + dt = Δt_cpl, + space = boundary_space, + saveat = saveat, + area_fraction = land_fraction, + date_ref = date0, + t_start = t_start, + ) + + ## ocean + ocean_sim = ocean_init( + FT; + tspan = tspan, + dt = Δt_cpl, + space = boundary_space, + saveat = saveat, + area_fraction = (FT(1) .- land_fraction), ## NB: this ocean fraction includes areas covered by sea ice (unlike the one contained in the cs) + thermo_params = thermo_params, + evolving = evolving_ocean, + ) + + ## sea ice (here set to zero area coverage) + ice_sim = SurfaceStub((; + T_sfc = ClimaCore.Fields.ones(boundary_space), + ρ_sfc = ClimaCore.Fields.zeros(boundary_space), + z0m = FT(0), + z0b = FT(0), + beta = FT(1), + α = FT(1), + area_fraction = ClimaCore.Fields.zeros(boundary_space), + phase = TD.Ice(), + thermo_params = thermo_params, + )) + + mode_specifics = (; name = mode_name, SST_info = nothing, SIC_info = nothing) +end + +#= +## Coupler Initialization +The coupler needs to contain exchange information, manage the calendar and be able to access all component models. It can also optionally +save online diagnostics. These are all initialized here and saved in a global `CouplerSimulation` struct, `cs`. +=# + +## coupler exchange fields +coupler_field_names = ( + :T_S, + :z0m_S, + :z0b_S, + :ρ_sfc, + :q_sfc, + :albedo, + :beta, + :F_turb_energy, + :F_turb_moisture, + :F_turb_ρτxz, + :F_turb_ρτyz, + :F_radiative, + :P_liq, + :P_snow, + :F_radiative_TOA, + :P_net, +) +coupler_fields = + NamedTuple{coupler_field_names}(ntuple(i -> ClimaCore.Fields.zeros(boundary_space), length(coupler_field_names))) + +## model simulations +model_sims = (atmos_sim = atmos_sim, ice_sim = ice_sim, land_sim = land_sim, ocean_sim = ocean_sim); + +## dates +dates = (; date = [date], date0 = [date0], date1 = [Dates.firstdayofmonth(date0)], new_month = [false]) + +#= +### Online Diagnostics +User can write custom diagnostics in the `user_diagnostics.jl`. +=# +monthly_3d_diags = init_diagnostics( + (:T, :u, :q_tot, :q_liq_ice), + atmos_sim.domain.center_space; + save = Monthly(), + operations = (; accumulate = TimeMean([Int(0)])), + output_dir = COUPLER_OUTPUT_DIR, + name_tag = "monthly_mean_3d_", +) + +monthly_2d_diags = init_diagnostics( + (:precipitation_rate, :toa_fluxes, :T_sfc, :tubulent_energy_fluxes), + boundary_space; + save = Monthly(), + operations = (; accumulate = TimeMean([Int(0)])), + output_dir = COUPLER_OUTPUT_DIR, + name_tag = "monthly_mean_2d_", +) + +diagnostics = (monthly_3d_diags, monthly_2d_diags) + +#= +## Initialize Conservation Checks +=# +## init conservation info collector +conservation_checks = nothing +if energy_check + @assert( + mode_name[1:10] == "slabplanet" && !CA.is_distributed(ClimaComms.context(boundary_space)), + "Only non-distributed slabplanet allowable for energy_check" + ) + conservation_checks = (; energy = EnergyConservationCheck(model_sims), water = WaterConservationCheck(model_sims)) +end + +dir_paths = (; output = COUPLER_OUTPUT_DIR, artifacts = COUPLER_ARTIFACTS_DIR) +checkpoint_cb = + HourlyCallback(dt = FT(480), func = checkpoint_sims, ref_date = [dates.date[1]], active = hourly_checkpoint) # 20 days +update_firstdayofmonth!_cb = + MonthlyCallback(dt = FT(1), func = update_firstdayofmonth!, ref_date = [dates.date1[1]], active = true) # for BCReader +callbacks = (; checkpoint = checkpoint_cb, update_firstdayofmonth! = update_firstdayofmonth!_cb) diff --git a/calibration/coupler_driver_calibration.jl b/calibration/coupler_driver_calibration.jl index 62f1930f6a..ea7149398e 100644 --- a/calibration/coupler_driver_calibration.jl +++ b/calibration/coupler_driver_calibration.jl @@ -1,385 +1,6 @@ -redirect_stderr(IOContext(stderr, :stacktrace_types_limited => Ref(false))) - -using ClimaComms -comms_ctx = ClimaComms.context() -const pid, nprocs = ClimaComms.init(comms_ctx) - - -import SciMLBase: ODEProblem, solve, step!, init, reinit! -using LinearAlgebra -import Test: @test -using Dates -using Plots -using Statistics: mean -import ClimaAtmos as CA -import YAML - -using ClimaCore.Utilities: half, PlusHalf -using ClimaCore: InputOutput, Fields -import ClimaCore.Spaces as Spaces - -## coupler specific imports -import ClimaCoupler -import ClimaCoupler.Regridder -import ClimaCoupler.Regridder: - update_surface_fractions!, combine_surfaces!, combine_surfaces_from_sol!, dummmy_remap!, binary_mask -import ClimaCoupler.ConservationChecker: - EnergyConservationCheck, WaterConservationCheck, check_conservation!, plot_global_conservation -import ClimaCoupler.Utilities: swap_space! -import ClimaCoupler.BCReader: - bcfile_info_init, float_type_bcf, update_midmonth_data!, next_date_in_file, interpolate_midmonth_to_daily -import ClimaCoupler.TimeManager: - current_date, - datetime_to_strdate, - trigger_callback, - Monthly, - EveryTimestep, - HourlyCallback, - MonthlyCallback, - update_firstdayofmonth!, - trigger_callback! -import ClimaCoupler.Diagnostics: get_var, init_diagnostics, accumulate_diagnostics!, save_diagnostics, TimeMean -import ClimaCoupler.PostProcessor: postprocess - -import ClimaCoupler.Interfacer: - CoupledSimulation, - float_type, - AtmosModelSimulation, - SurfaceModelSimulation, - SurfaceStub, - SeaIceModelSimulation, - LandModelSimulation, - OceanModelSimulation, - get_field, - update_field! -import ClimaCoupler.FluxCalculator: - PartitionedStateFluxes, - CombinedStateFluxes, - combined_turbulent_fluxes!, - MoninObukhovScheme, - partitioned_turbulent_fluxes! -import ClimaCoupler.FieldExchanger: - import_atmos_fields!, - import_combined_surface_fields!, - update_sim!, - update_model_sims!, - reinit_model_sims!, - step_model_sims! -import ClimaCoupler.Checkpointer: checkpoint_model_state, get_model_state_vector, restart_model_state! - -## helpers for component models -include("../experiments/AMIP/components/atmosphere/climaatmos_init.jl") -include("../experiments/AMIP/components/land/bucket_init.jl") -include("../experiments/AMIP/components/land/bucket_utils.jl") -include("../experiments/AMIP/components/ocean/slab_ocean_init.jl") -include("../experiments/AMIP/components/ocean/prescr_seaice_init.jl") -include("../experiments/AMIP/user_io/user_diagnostics.jl") -include("../experiments/AMIP/user_io/user_logging.jl") - -## coupler defaults -# get component model dictionaries -include("../experiments/AMIP/cli_options.jl") -parsed_args = parse_commandline(argparse_settings()) -config_dict = YAML.load_file("./experiments/amip_coupled/coupler_config.yml") -config_dict = YAML.load_file(joinpath(experiment_dir, "coupler_config.yml")); -config_dict["t_end"] = "150secs"; -config_dict["output_dir"] = output_dir; -config_dict = merge(parsed_args, config_dict) -config_dict_atmos = get_atmos_config(config_dict) - -# merge dictionaries of command line arguments, coupler dictionary and component model dictionaries -# (if there are common keys, the last dictorionary in the `merge` arguments takes precedence) -config_dict = merge(config_dict_atmos, config_dict) - - -## read in some parsed command line arguments -mode_name = config_dict["mode_name"] -run_name = config_dict["run_name"] -energy_check = config_dict["energy_check"] -FT = config_dict["FLOAT_TYPE"] == "Float64" ? Float64 : Float32 -land_sim_name = "bucket" -t_end = Float64(time_to_seconds(config_dict["t_end"])) -t_start = 0.0 -tspan = (t_start, t_end) -Δt_cpl = Float64(config_dict["dt_cpl"]) -saveat = Float64(time_to_seconds(config_dict["dt_save_to_sol"])) -date0 = date = DateTime(config_dict["start_date"], dateformat"yyyymmdd") -mono_surface = config_dict["mono_surface"] -hourly_checkpoint = config_dict["hourly_checkpoint"] -restart_dir = config_dict["restart_dir"] -restart_t = Int(config_dict["restart_t"]) -evolving_ocean = config_dict["evolving_ocean"] -config_dict["print_config_dict"] = false - -## I/O directory setup -COUPLER_OUTPUT_DIR = "/Users/akshaysridhar/Research/Codes/ClimaCoupler.jl/calibration/output/amip/" -mkpath(COUPLER_OUTPUT_DIR) - -REGRID_DIR = joinpath(COUPLER_OUTPUT_DIR, "regrid_tmp/") -mkpath(REGRID_DIR) - -COUPLER_ARTIFACTS_DIR = COUPLER_OUTPUT_DIR * "_artifacts" -isdir(COUPLER_ARTIFACTS_DIR) ? nothing : mkpath(COUPLER_ARTIFACTS_DIR) - -config_dict["print_config_dict"] ? @info(config_dict) : nothing - -# get the paths to the necessary data files: land-sea mask, sst map, sea ice concentration -include(joinpath(pkgdir(ClimaCoupler), "artifacts", "artifact_funcs.jl")) -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") - -config_dict_atmos["output_dir"] = COUPLER_OUTPUT_DIR -atmos_sim = atmos_init(FT, config_dict_atmos); -thermo_params = get_thermo_params(atmos_sim) # TODO: this should be shared by all models - -#= -We use a common `Space` for all global surfaces. This enables the MPI processes to operate on the same columns in both -the atmospheric and surface components, so exchanges are parallelized. Note this is only possible when the -atmosphere and surface are of the same horizontal resolution. -=# -## init a 2D boundary space at the surface -boundary_space = Spaces.horizontal_space(atmos_sim.domain.face_space) - -# init land-sea fraction -land_fraction = - FT.( - Regridder.land_fraction( - FT, - REGRID_DIR, - comms_ctx, - land_mask_data, - "LSMASK", - boundary_space, - mono = mono_surface, - ) - ) - -@info mode_name -if mode_name == "amip" - @info "AMIP boundary conditions - do not expect energy conservation" - - ## land - land_sim = bucket_init( - FT, - tspan, - config_dict["land_domain_type"], - config_dict["land_albedo_type"], - config_dict["land_temperature_anomaly"], - comms_ctx, - REGRID_DIR; - dt = Δt_cpl, - space = boundary_space, - saveat = saveat, - area_fraction = land_fraction, - date_ref = date0, - t_start = t_start, - ) - - ## ocean - SST_info = bcfile_info_init( - FT, - REGRID_DIR, - sst_data, - "SST", - boundary_space, - comms_ctx, - interpolate_daily = true, - scaling_function = clean_sst, ## convert to Kelvin - land_fraction = land_fraction, - date0 = date0, - mono = mono_surface, - ) - - update_midmonth_data!(date0, SST_info) - SST_init = interpolate_midmonth_to_daily(date0, SST_info) - ocean_sim = SurfaceStub((; - T_sfc = SST_init, - ρ_sfc = ClimaCore.Fields.zeros(boundary_space), - z0m = FT(1e-3), - z0b = FT(1e-3), - beta = FT(1), - α = FT(0.06), - area_fraction = (FT(1) .- land_fraction), - phase = TD.Liquid(), - thermo_params = thermo_params, - )) - - ## sea ice - SIC_info = bcfile_info_init( - FT, - REGRID_DIR, - sic_data, - "SEAICE", - boundary_space, - comms_ctx, - interpolate_daily = true, - scaling_function = clean_sic, ## convert to fraction - land_fraction = land_fraction, - date0 = date0, - mono = mono_surface, - ) - update_midmonth_data!(date0, SIC_info) - SIC_init = interpolate_midmonth_to_daily(date0, SIC_info) - ice_fraction = get_ice_fraction.(SIC_init, mono_surface) - ice_sim = ice_init( - FT; - tspan = tspan, - dt = Δt_cpl, - space = boundary_space, - saveat = saveat, - area_fraction = ice_fraction, - thermo_params = thermo_params, - ) - - ## CO2 concentration - CO2_info = bcfile_info_init( - FT, - REGRID_DIR, - co2_data, - "co2", - boundary_space, - comms_ctx, - interpolate_daily = true, - land_fraction = ones(boundary_space), - date0 = date0, - mono = mono_surface, - ) - - update_midmonth_data!(date0, CO2_info) - CO2_init = interpolate_midmonth_to_daily(date0, CO2_info) - update_field!(atmos_sim, Val(:co2_gm), CO2_init) - - mode_specifics = (; name = mode_name, SST_info = SST_info, SIC_info = SIC_info, CO2_info = CO2_info) - -elseif mode_name in ("slabplanet", "slabplanet_aqua", "slabplanet_terra") - - land_fraction = mode_name == "slabplanet_aqua" ? land_fraction .* 0 : land_fraction - land_fraction = mode_name == "slabplanet_terra" ? land_fraction .* 0 .+ 1 : land_fraction - - ## land - land_sim = bucket_init( - FT, - tspan, - config_dict["land_domain_type"], - config_dict["land_albedo_type"], - config_dict["land_temperature_anomaly"], - comms_ctx, - REGRID_DIR; - dt = Δt_cpl, - space = boundary_space, - saveat = saveat, - area_fraction = land_fraction, - date_ref = date0, - t_start = t_start, - ) - - ## ocean - ocean_sim = ocean_init( - FT; - tspan = tspan, - dt = Δt_cpl, - space = boundary_space, - saveat = saveat, - area_fraction = (FT(1) .- land_fraction), ## NB: this ocean fraction includes areas covered by sea ice (unlike the one contained in the cs) - thermo_params = thermo_params, - evolving = evolving_ocean, - ) - - ## sea ice (here set to zero area coverage) - ice_sim = SurfaceStub((; - T_sfc = ClimaCore.Fields.ones(boundary_space), - ρ_sfc = ClimaCore.Fields.zeros(boundary_space), - z0m = FT(0), - z0b = FT(0), - beta = FT(1), - α = FT(1), - area_fraction = ClimaCore.Fields.zeros(boundary_space), - phase = TD.Ice(), - thermo_params = thermo_params, - )) - - mode_specifics = (; name = mode_name, SST_info = nothing, SIC_info = nothing) -end - -#= -## Coupler Initialization -The coupler needs to contain exchange information, manage the calendar and be able to access all component models. It can also optionally -save online diagnostics. These are all initialized here and saved in a global `CouplerSimulation` struct, `cs`. -=# - -## coupler exchange fields -coupler_field_names = ( - :T_S, - :z0m_S, - :z0b_S, - :ρ_sfc, - :q_sfc, - :albedo, - :beta, - :F_turb_energy, - :F_turb_moisture, - :F_turb_ρτxz, - :F_turb_ρτyz, - :F_radiative, - :P_liq, - :P_snow, - :F_radiative_TOA, - :P_net, -) -coupler_fields = - NamedTuple{coupler_field_names}(ntuple(i -> ClimaCore.Fields.zeros(boundary_space), length(coupler_field_names))) - -## model simulations -model_sims = (atmos_sim = atmos_sim, ice_sim = ice_sim, land_sim = land_sim, ocean_sim = ocean_sim); - -## dates -dates = (; date = [date], date0 = [date0], date1 = [Dates.firstdayofmonth(date0)], new_month = [false]) - -#= -### Online Diagnostics -User can write custom diagnostics in the `user_diagnostics.jl`. -=# -monthly_3d_diags = init_diagnostics( - (:T, :u, :q_tot, :q_liq_ice), - atmos_sim.domain.center_space; - save = Monthly(), - operations = (; accumulate = TimeMean([Int(0)])), - output_dir = COUPLER_OUTPUT_DIR, - name_tag = "monthly_mean_3d_", -) - -monthly_2d_diags = init_diagnostics( - (:precipitation_rate, :toa_fluxes, :T_sfc, :tubulent_energy_fluxes), - boundary_space; - save = Monthly(), - operations = (; accumulate = TimeMean([Int(0)])), - output_dir = COUPLER_OUTPUT_DIR, - name_tag = "monthly_mean_2d_", -) - -diagnostics = (monthly_3d_diags, monthly_2d_diags) - -#= -## Initialize Conservation Checks -=# -## init conservation info collector -conservation_checks = nothing -if energy_check - @assert( - mode_name[1:10] == "slabplanet" && !CA.is_distributed(ClimaComms.context(boundary_space)), - "Only non-distributed slabplanet allowable for energy_check" - ) - conservation_checks = (; energy = EnergyConservationCheck(model_sims), water = WaterConservationCheck(model_sims)) -end - -dir_paths = (; output = COUPLER_OUTPUT_DIR, artifacts = COUPLER_ARTIFACTS_DIR) -checkpoint_cb = - HourlyCallback(dt = FT(480), func = checkpoint_sims, ref_date = [dates.date[1]], active = hourly_checkpoint) # 20 days -update_firstdayofmonth!_cb = - MonthlyCallback(dt = FT(1), func = update_firstdayofmonth!, ref_date = [dates.date1[1]], active = true) # for BCReader -callbacks = (; checkpoint = checkpoint_cb, update_firstdayofmonth! = update_firstdayofmonth!_cb) +include("coupler_driver_init.jl") +include("coupler_parse_args.jl") +include("coupler_component_init.jl") ## coupler simulation cs = CoupledSimulation{FT}( diff --git a/calibration/coupler_driver_init.jl b/calibration/coupler_driver_init.jl new file mode 100644 index 0000000000..806a83ae1c --- /dev/null +++ b/calibration/coupler_driver_init.jl @@ -0,0 +1,77 @@ +redirect_stderr(IOContext(stderr, :stacktrace_types_limited => Ref(false))) + +using ClimaComms +comms_ctx = ClimaComms.context() +const pid, nprocs = ClimaComms.init(comms_ctx) + + +import SciMLBase: ODEProblem, solve, step!, init, reinit! +using LinearAlgebra +import Test: @test +using Dates +using Plots +using Statistics: mean +import ClimaAtmos as CA +import YAML + +using ClimaCore.Utilities: half, PlusHalf +using ClimaCore: InputOutput, Fields +import ClimaCore.Spaces as Spaces + +## coupler specific imports +import ClimaCoupler +import ClimaCoupler.Regridder +import ClimaCoupler.Regridder: + update_surface_fractions!, combine_surfaces!, combine_surfaces_from_sol!, dummmy_remap!, binary_mask +import ClimaCoupler.ConservationChecker: + EnergyConservationCheck, WaterConservationCheck, check_conservation!, plot_global_conservation +import ClimaCoupler.Utilities: swap_space! +import ClimaCoupler.BCReader: + bcfile_info_init, float_type_bcf, update_midmonth_data!, next_date_in_file, interpolate_midmonth_to_daily +import ClimaCoupler.TimeManager: + current_date, + datetime_to_strdate, + trigger_callback, + Monthly, + EveryTimestep, + HourlyCallback, + MonthlyCallback, + update_firstdayofmonth!, + trigger_callback! +import ClimaCoupler.Diagnostics: get_var, init_diagnostics, accumulate_diagnostics!, save_diagnostics, TimeMean +import ClimaCoupler.PostProcessor: postprocess + +import ClimaCoupler.Interfacer: + CoupledSimulation, + float_type, + AtmosModelSimulation, + SurfaceModelSimulation, + SurfaceStub, + SeaIceModelSimulation, + LandModelSimulation, + OceanModelSimulation, + get_field, + update_field! +import ClimaCoupler.FluxCalculator: + PartitionedStateFluxes, + CombinedStateFluxes, + combined_turbulent_fluxes!, + MoninObukhovScheme, + partitioned_turbulent_fluxes! +import ClimaCoupler.FieldExchanger: + import_atmos_fields!, + import_combined_surface_fields!, + update_sim!, + update_model_sims!, + reinit_model_sims!, + step_model_sims! +import ClimaCoupler.Checkpointer: checkpoint_model_state, get_model_state_vector, restart_model_state! + +## helpers for component models +include("../experiments/AMIP/components/atmosphere/climaatmos_init.jl") +include("../experiments/AMIP/components/land/bucket_init.jl") +include("../experiments/AMIP/components/land/bucket_utils.jl") +include("../experiments/AMIP/components/ocean/slab_ocean_init.jl") +include("../experiments/AMIP/components/ocean/prescr_seaice_init.jl") +include("../experiments/AMIP/user_io/user_diagnostics.jl") +include("../experiments/AMIP/user_io/user_logging.jl") diff --git a/calibration/coupler_parse_args.jl b/calibration/coupler_parse_args.jl new file mode 100644 index 0000000000..1aa6121891 --- /dev/null +++ b/calibration/coupler_parse_args.jl @@ -0,0 +1,47 @@ +## coupler defaults +# get component model dictionaries +include("../experiments/AMIP/cli_options.jl") +parsed_args = parse_commandline(argparse_settings()) +config_dict = YAML.load_file("./experiments/amip_coupled/coupler_config.yml") +config_dict = YAML.load_file(joinpath(experiment_dir, "coupler_config.yml")); +config_dict["t_end"] = "150secs"; +config_dict["output_dir"] = output_dir; +config_dict = merge(parsed_args, config_dict) +config_dict_atmos = get_atmos_config(config_dict) + +# merge dictionaries of command line arguments, coupler dictionary and component model dictionaries +# (if there are common keys, the last dictorionary in the `merge` arguments takes precedence) +config_dict = merge(config_dict_atmos, config_dict) + + +## read in some parsed command line arguments +mode_name = config_dict["mode_name"] +run_name = config_dict["run_name"] +energy_check = config_dict["energy_check"] +FT = config_dict["FLOAT_TYPE"] == "Float64" ? Float64 : Float32 +land_sim_name = "bucket" +t_end = Float64(time_to_seconds(config_dict["t_end"])) +t_start = 0.0 +tspan = (t_start, t_end) +Δt_cpl = Float64(config_dict["dt_cpl"]) +saveat = Float64(time_to_seconds(config_dict["dt_save_to_sol"])) +date0 = date = DateTime(config_dict["start_date"], dateformat"yyyymmdd") +mono_surface = config_dict["mono_surface"] +hourly_checkpoint = config_dict["hourly_checkpoint"] +restart_dir = config_dict["restart_dir"] +restart_t = Int(config_dict["restart_t"]) +evolving_ocean = config_dict["evolving_ocean"] +config_dict["print_config_dict"] = false + +## I/O directory setup +COUPLER_OUTPUT_DIR = "/Users/akshaysridhar/Research/Codes/ClimaCoupler.jl/calibration/output/amip/" +mkpath(COUPLER_OUTPUT_DIR) + +REGRID_DIR = joinpath(COUPLER_OUTPUT_DIR, "regrid_tmp/") +mkpath(REGRID_DIR) + +COUPLER_ARTIFACTS_DIR = COUPLER_OUTPUT_DIR * "_artifacts" +isdir(COUPLER_ARTIFACTS_DIR) ? nothing : mkpath(COUPLER_ARTIFACTS_DIR) + +config_dict["print_config_dict"] ? @info(config_dict) : nothing +config_dict_atmos["output_dir"] = COUPLER_OUTPUT_DIR diff --git a/src/Regridder.jl b/src/Regridder.jl index ce166a9323..3fd81c2918 100644 --- a/src/Regridder.jl +++ b/src/Regridder.jl @@ -207,7 +207,7 @@ function hdwrite_regridfile_rll_to_cgll( # read the remapped file with sparse matrices offline_outvector, coords = NCDataset(datafile_cgll, "r") do ds_wt ( - offline_outvector = Array(ds_wt[varname])[:, :, :], # ncol, z, times + offline_outvector = nomissing(Array(ds_wt[varname])[:, :]), # ncol, times coords = get_coords(ds_wt, space), ) end From 1018f51b4d3db0d4f444d1887d497553cca2b200 Mon Sep 17 00:00:00 2001 From: akshaysridhar Date: Wed, 6 Mar 2024 17:58:28 -0800 Subject: [PATCH 2/3] modified: coupler_driver_calibration.jl modified: coupler_parse_args.jl new file: experiments/amip_coupled/truth_simulation/target_amip_n1_shortrun.yml new file: experiments/amip_coupled/truth_simulation/target_amip_n1_shortrun_parameters.toml modified: generate_observations.sbatch --- calibration/coupler_driver_calibration.jl | 19 +- calibration/coupler_parse_args.jl | 1 - .../target_amip_n1_shortrun.yml | 106 ++++ .../target_amip_n1_shortrun_parameters.toml | 487 ++++++++++++++++++ calibration/generate_observations.sbatch | 5 +- 5 files changed, 598 insertions(+), 20 deletions(-) create mode 100644 calibration/experiments/amip_coupled/truth_simulation/target_amip_n1_shortrun.yml create mode 100644 calibration/experiments/amip_coupled/truth_simulation/target_amip_n1_shortrun_parameters.toml diff --git a/calibration/coupler_driver_calibration.jl b/calibration/coupler_driver_calibration.jl index ea7149398e..cf423ead4a 100644 --- a/calibration/coupler_driver_calibration.jl +++ b/calibration/coupler_driver_calibration.jl @@ -1,6 +1,4 @@ -include("coupler_driver_init.jl") -include("coupler_parse_args.jl") -include("coupler_component_init.jl") + ## coupler simulation cs = CoupledSimulation{FT}( @@ -21,21 +19,6 @@ cs = CoupledSimulation{FT}( dir_paths, ); -#= -## Restart component model states if specified -=# -#if restart_dir !== "unspecified" -# for sim in cs.model_sims -# if get_model_state_vector(sim) !== nothing -# @skipping restart -# restart_model_state!(sim, comms_ctx, restart_t; input_dir = restart_dir) -# end -# end -#end - -#= -## Initialize Component Model Exchange -=# turbulent_fluxes = nothing if config_dict["turb_flux_partition"] == "PartitionedStateFluxes" turbulent_fluxes = PartitionedStateFluxes() diff --git a/calibration/coupler_parse_args.jl b/calibration/coupler_parse_args.jl index 1aa6121891..de943a63da 100644 --- a/calibration/coupler_parse_args.jl +++ b/calibration/coupler_parse_args.jl @@ -34,7 +34,6 @@ evolving_ocean = config_dict["evolving_ocean"] config_dict["print_config_dict"] = false ## I/O directory setup -COUPLER_OUTPUT_DIR = "/Users/akshaysridhar/Research/Codes/ClimaCoupler.jl/calibration/output/amip/" mkpath(COUPLER_OUTPUT_DIR) REGRID_DIR = joinpath(COUPLER_OUTPUT_DIR, "regrid_tmp/") diff --git a/calibration/experiments/amip_coupled/truth_simulation/target_amip_n1_shortrun.yml b/calibration/experiments/amip_coupled/truth_simulation/target_amip_n1_shortrun.yml new file mode 100644 index 0000000000..ed7aa7ffb0 --- /dev/null +++ b/calibration/experiments/amip_coupled/truth_simulation/target_amip_n1_shortrun.yml @@ -0,0 +1,106 @@ +surface_thermo_state_type: "GCMSurfaceThermoState" +topo_smoothing: false +warn_allocations_diagnostics: false +hyperdiff: "ClimaHyperdiffusion" +dt: "150secs" +output_dir: "experiments/amip_coupled/truth_simulation" +prognostic_tke: false +override_τ_precip: true +use_newton_rtol: false +netcdf_output_at_levels: false +device: "auto" +t_end: "150secs" +dz_top: 3000.0 +y_elem: 6 +z_stretch: false +bubble: true +ode_algo: "ARS343" +max_newton_iters_ode: 1 +start_date: "19790101" +check_precipitation: false +forcing: ~ +edmfx_nh_pressure: false +scalar_hyperdiffusion_coefficient: 1.5 +prognostic_surface: "false" +test_dycore_consistency: false +moist: "equil" +perf_mode: "PerfStandard" +edmf_coriolis: ~ +rad: "gray" +rayleigh_sponge: true +initial_condition: "DecayingProfile" +cloud_model: "quadrature" +krylov_rtol: 0.1 +divergence_damping_factor: 1.0 +edmfx_entr_model: ~ +eisenstat_walker_forcing_alpha: 2.0 +dt_cloud_fraction: "3hours" +smoothing_order: 3 +idealized_h2o: false +surface_setup: "PrescribedSurface" +perturb_initstate: true +jvp_step_adjustment: 1.0 +discrete_hydrostatic_balance: false +netcdf_interpolate_z_over_msl: false +log_progress: true +dz_bottom: 30.0 +h_elem: 16 +dt_save_state_to_disk: "Inf" +netcdf_interpolation_num_points: ~ +advection_test: false +z_max: 30000.0 +apply_limiter: false +topography: "NoWarp" +reference_job_id: ~ +precip_model: "0M" +perf_summary: false +vorticity_hyperdiffusion_coefficient: 1.5 +viscous_sponge: false +surface_temperature: "ZonallySymmetric" +diagnostics: + - short_name: + - "pfull" + - "wa" + - "va" + - "rv" + period: "150secs" + reduction: "average" +job_id: "target_amip_n1_shortrun" +orographic_gravity_wave: ~ +dt_rad: "1hours" +approximate_linear_solve_iters: 1 +edmfx_upwinding: "none" +tracer_upwinding: "none" +nh_poly: 3 +edmfx_sgs_diffusive_flux: false +y_max: 300000.0 +non_orographic_gravity_wave: false +use_reference_state: true +config: "sphere" +energy_upwinding: "none" +FLOAT_TYPE: "Float64" +updraft_number: 1 +split_ode: true +regression_test: false +check_conservation: false +ls_adv: ~ +output_default_diagnostics: true +implicit_diffusion: false +x_max: 300000.0 +edmfx_sgs_mass_flux: false +z_elem: 50 +newton_rtol: 1.0e-5 +fps: 5 +edmfx_sgsflux_upwinding: "none" +turbconv: ~ +x_elem: 6 +idealized_clouds: false +vert_diff: "true" +use_krylov_method: false +subsidence: ~ +use_dynamic_krylov_rtol: false +idealized_insolation: true +toml: + - "/Users/akshaysridhar/.julia/packages/ClimaCoupler/utcpx/toml/default_coarse.toml" +edmfx_detr_model: ~ +dt_save_to_sol: "1days" diff --git a/calibration/experiments/amip_coupled/truth_simulation/target_amip_n1_shortrun_parameters.toml b/calibration/experiments/amip_coupled/truth_simulation/target_amip_n1_shortrun_parameters.toml new file mode 100644 index 0000000000..1e2ece0d34 --- /dev/null +++ b/calibration/experiments/amip_coupled/truth_simulation/target_amip_n1_shortrun_parameters.toml @@ -0,0 +1,487 @@ +[f_plane_coriolis_frequency] +used_in = ["ClimaAtmos"] +value = 0 +type = "float" + +[orbit_eccentricity_at_epoch] +used_in = ["Insolation"] +value = 0.016708634 +type = "float" + +[potential_temperature_reference_pressure] +used_in = ["Thermodynamics"] +value = 100000 +type = "float" +description = "Reference pressure used in potential temperature definition" + +[prandtl_number_0_businger] +used_in = ["SurfaceFluxes"] +value = 0.74 +type = "float" +description = "Pr_0 for Businger universal functions. From Businger et al, 1971. DOI: 10.1175/1520-0469(1971)028<0181:FPRITA>2.0.CO;2." + +[mixing_length_tke_surf_scale] +used_in = ["ClimaAtmos"] +value = 3.75 +type = "float" +description = "Ratio of turbulence kinetic energy to squared friction velocity in the surface layer for the EDMF mixing length closure; denoted κ_*². See Lopez-Gomez et al. (2020) [https://doi.org/10.1029/2020MS002162], Table 1. Note: the square root, i.e. κ_*, is listed in the reference." + +[temperature_triple_point] +used_in = ["Thermodynamics"] +value = 273.16 +type = "float" + +[length_orbit_semi_major] +used_in = ["Insolation"] +value = 149597870000 +type = "float" +description = "derived: 1 * [astronomical_unit]" + +[isobaric_specific_heat_ice] +used_in = ["Thermodynamics"] +value = 2100 +type = "float" + +[mixing_length_smin_rm] +used_in = ["ClimaAtmos"] +value = 1.5 +type = "float" +description = "Upper ratio limit for smooth minimum function in mixing length closure. See Lopez-Gomez et al. (2020) Eq 40 [https://doi.org/10.1029/2020MS002161]." + +[C_E] +used_in = ["ClimaAtmos"] +value = 0.044 +type = "float" +description = "vertical diffusion coefficient" + +[latent_heat_sublimation_at_reference] +used_in = ["Thermodynamics"] +value = 2834400 +type = "float" + +[day] +used_in = ["Insolation"] +value = 86400 +type = "float" + +[isobaric_specific_heat_vapor] +used_in = ["Thermodynamics"] +value = 1859 +type = "float" + +[EDMF_surface_area] +used_in = ["ClimaAtmos"] +value = 0.1 +type = "float" +description = "Combined updraft surface area fraction; used to compute boundary conditions for prognostic updraft variables. The surface area for each updraft is `surface_area / N_updrafts`. See Cohen et al. (2020) [https://doi.org/10.1029/2020MS002162], Table 2." + +[c_smag] +used_in = ["ClimaAtmos"] +value = 0.2 +type = "float" +description = "Smagorinsky coefficient" + +[mixing_length_smin_ub] +used_in = ["ClimaAtmos"] +value = 0.1 +type = "float" +description = "Lower limit for smooth minimum function in mixing length closure. See Lopez-Gomez et al. (2020) Eq 40 [https://doi.org/10.1029/2020MS002161]." + +[molar_mass_water] +used_in = ["Thermodynamics", "RRTMGP"] +value = 0.01801528 +type = "float" + +[most_stability_exponent_businger] +used_in = ["SurfaceFluxes"] +value = 4.42 +type = "float" +description = "γ for Businger universal functions. From Businger et al, 1971. DOI: 10.1175/1520-0469(1971)028<0181:FPRITA>2.0.CO;2." + +[temperature_saturation_adjustment_min] +used_in = ["Thermodynamics"] +value = 150 +type = "float" + +[min_area_limiter_power] +used_in = ["ClimaAtmos"] +value = 10 +type = "float" +description = "Constant coefficient for the exponent in the minimum area limiter term in entrainment. Parameter not described in the literature." + +[pressure_normalmode_drag_coeff] +used_in = ["ClimaAtmos"] +value = 10.0 +type = "float" +description = "Updraft pressure drag coefficent in perturbation pressure closure. See He et al. 2022 Eq 34 [https://doi.org/10.1002/essoar.10505084.2]." + +[held_suarez_T_equator_dry] +used_in = ["ClimaAtmos"] +value = 315 +type = "float" +description = "Equator temperature. See Held and Suarez (1994) https://doi.org/10.1175/1520-0477(1994)075%3C1825:APFTIO%3E2.0.CO;2" + +[mixing_length_Prandtl_number_0] +used_in = ["ClimaAtmos"] +value = 0.74 +type = "float" +description = "Turbulent Prandtl number in neutral conditions; denoted Pr_{t,0}. See Lopez-Gomez et al. (2020) [https://doi.org/10.1029/2020MS002162], Table 1 and Eq 36." + +[longitude_perihelion_at_epoch] +used_in = ["Insolation"] +value = 4.938188299449 +type = "float" +description = "(282.937348 degrees) in radians" + +[stefan_boltzmann_constant] +used_in = ["RRTMGP"] +value = 5.67e-8 +type = "float" + +[temperature_min_at_reference] +used_in = ["Thermodynamics"] +value = 220 +type = "float" + +[entropy_water_vapor] +used_in = ["Thermodynamics"] +value = 10513.6 +type = "float" + +[equator_pole_temperature_gradient_wet] +used_in = ["ClimaAtmos"] +value = 65 +type = "float" +description = "Temperature gradient between equator and pole for moist adiabatic atmosphere. See Held and Suarez (1994) https://doi.org/10.1175/1520-0477(1994)075%3C1825:APFTIO%3E2.0.CO;2" + +[held_suarez_T_equator_wet] +used_in = ["ClimaAtmos"] +value = 294 +type = "float" +description = "Equator temperature. See Held and Suarez (1994) https://doi.org/10.1175/1520-0477(1994)075%3C1825:APFTIO%3E2.0.CO;2" + +[entropy_dry_air] +used_in = ["Thermodynamics"] +value = 6864.8 +type = "float" + +[zd_rayleigh] +used_in = ["ClimaAtmos"] +value = 15000.0 +type = "float" +description = "rayleigh sponge height" + +[coefficient_a_m_businger] +used_in = ["SurfaceFluxes"] +value = 4.7 +type = "float" +description = "a_m for Businger momentum universal functions. From Businger et al, 1971. DOI: 10.1175/1520-0469(1971)028<0181:FPRITA>2.0.CO;2." + +[temperature_saturation_adjustment_max] +used_in = ["Thermodynamics"] +value = 1000 +type = "float" + +[mean_anomalistic_at_epoch] +used_in = ["Insolation"] +value = 6.24006014121 +type = "float" +description = "(357.52911 degrees) in radians" + +[mixing_length_diss_coeff] +used_in = ["ClimaAtmos"] +value = 0.22 +type = "float" +description = "Turbulence kinetic energy dissipation coefficient for the EDMF mixing length closure; denoted c_d. See Lopez-Gomez et al. (2020) [https://doi.org/10.1029/2020MS002162], Table 1." + +[potential_temp_vertical_gradient] +used_in = ["ClimaAtmos"] +value = 10 +type = "float" +description = "Potential temperature gradient with height. See Held and Suarez (1994) https://doi.org/10.1175/1520-0477(1994)075%3C1825:APFTIO%3E2.0.CO;2" + +[EDMF_max_area] +used_in = ["ClimaAtmos"] +value = 0.9 +type = "float" +description = "Maximum area fraction per updraft. Parameter not described in the literature." + +[entr_coeff] +used_in = ["ClimaAtmos"] +value = 1 +type = "float" +description = "TODO: Remove this. Constant entrainment coefficient used for testing EDMF" + +[mixing_length_Ri_crit] +used_in = ["ClimaAtmos"] +value = 0.25 +type = "float" +description = "Critical gradient Richardson number. It is an upper limit to the gradient Richardson number . See Li (2019) [https://doi.org/10.1016/j.atmosres.2018.09.015], Section 6.2 for details." + +[gas_constant] +used_in = ["Thermodynamics", "RRTMGP"] +value = 8.3144598 +type = "float" + +[orbit_obliquity_at_epoch] +used_in = ["Insolation"] +value = 0.408979125113246 +type = "float" +description = "(23.432777778 degrees) in radians" + +[drag_layer_vertical_extent] +used_in = ["ClimaAtmos"] +value = 0.7 +type = "float" +description = "Vertical extend of drag layer. See Held and Suarez (1994) https://doi.org/10.1175/1520-0477(1994)075%3C1825:APFTIO%3E2.0.CO;2" + +[temperature_mean_at_reference] +used_in = ["Thermodynamics"] +value = 290 +type = "float" + +[detr_inv_tau] +used_in = ["ClimaAtmos"] +value = 900 +type = "float" +description = "Detrainment timescale" + +[entr_inv_tau] +used_in = ["ClimaAtmos"] +value = 900 +type = "float" +description = "Entrainment timescale" + +[temperature_homogenous_nucleation] +used_in = ["Thermodynamics"] +value = 233 +type = "float" + +[pressure_normalmode_buoy_coeff1] +used_in = ["ClimaAtmos"] +value = 0.12 +type = "float" +description = "Pressure buoyancy coefficient (encapsulating virtual mass loading effect) in perturbation pressure closure. See He et al. 2022 Eq 34 [https://doi.org/10.1002/essoar.10505084.2]." + +[mixing_length_static_stab_coeff] +used_in = ["ClimaAtmos"] +value = 0.4 +type = "float" +description = "Static stability coefficient for the EDMF mixing length closure; denoted c_b. See Lopez-Gomez et al. (2020) [https://doi.org/10.1029/2020MS002162], Table 1." + +[held_suarez_minimum_temperature] +used_in = ["ClimaAtmos"] +value = 200 +type = "float" +description = "Minimum temperature. See Held and Suarez (1994) https://doi.org/10.1175/1520-0477(1994)075%3C1825:APFTIO%3E2.0.CO;2" + +[angular_velocity_planet_rotation] +used_in = ["ClimaAtmos"] +value = 7.2921159e-5 +type = "float" + +[min_area_limiter_scale] +used_in = ["ClimaAtmos"] +value = 0.001 +type = "float" +description = "Constant coefficient that scales the minimum area limiter term in entrainment. Parameter not described in the literature." + +[entropy_reference_temperature] +used_in = ["Thermodynamics"] +value = 298.15 +type = "float" + +[most_stability_parameter_businger] +used_in = ["SurfaceFluxes"] +value = 2.5 +type = "float" +description = "ζ_a for Businger universal functions. From Businger et al, 1971. DOI: 10.1175/1520-0469(1971)028<0181:FPRITA>2.0.CO;2." + +[molar_mass_dry_air] +used_in = ["Thermodynamics", "RRTMGP"] +value = 0.02897 +type = "float" + +[zd_viscous] +used_in = ["ClimaAtmos"] +value = 15000.0 +type = "float" +description = "viscous sponge height" + +[planet_radius] +used_in = ["ClimaAtmos"] +value = 6371000 +type = "float" + +[adiabatic_exponent_dry_air] +used_in = ["Thermodynamics", "RRTMGP"] +value = 0.28571428571 +type = "float" +description = "(2/7)" + +[epoch_time] +used_in = ["Insolation"] +value = 211813488000 +type = "float" +description = "derived: 2451545.0 * [day]" + +[detr_buoy_coeff] +used_in = ["ClimaAtmos"] +value = 0.12 +type = "float" +description = "Coefficient for the b/w^2 term in the detrainment closure. See Tan et al. (2018) [https://doi.org/10.1002/2017MS001162], Eq 27." + +[precipitation_timescale] +used_in = ["CloudMicrophysics"] +value = 150.0 +type = "float" + +[C_H] +used_in = ["ClimaAtmos"] +value = 0.0044 +type = "float" +description = "bulk transfer coefficient" + +[EDMF_min_area] +used_in = ["ClimaAtmos"] +value = 1.0e-5 +type = "float" +description = "Minimum area fraction per updraft. Parameter not described in the literature." + +[isobaric_specific_heat_liquid] +used_in = ["Thermodynamics"] +value = 4181 +type = "float" + +[latent_heat_vaporization_at_reference] +used_in = ["Thermodynamics"] +value = 2500800 +type = "float" + +[detr_vertdiv_coeff] +used_in = ["ClimaAtmos"] +value = 1 +type = "float" +description = "Coefficient for the vertical divergence term in the detrainment closure. Parameter not described in the literature." + +[gravitational_acceleration] +used_in = ["Thermodynamics", "RRTMGP"] +value = 9.81 +type = "float" + +[equator_pole_temperature_gradient_dry] +used_in = ["ClimaAtmos"] +value = 60 +type = "float" +description = "Temperature gradient between equator and pole for dry adiabatic atmosphere. See Held and Suarez (1994) https://doi.org/10.1175/1520-0477(1994)075%3C1825:APFTIO%3E2.0.CO;2" + +[coefficient_a_h_businger] +used_in = ["SurfaceFluxes"] +value = 4.7 +type = "float" +description = "a_h for Businger heat universal functions. From Businger et al, 1971. DOI: 10.1175/1520-0469(1971)028<0181:FPRITA>2.0.CO;2." + +[pow_icenuc] +used_in = ["Thermodynamics"] +value = 1 +type = "float" + +[max_area_limiter_scale] +used_in = ["ClimaAtmos"] +value = 0.1 +type = "float" +description = "Constant coefficient that scales the maximum area limiter term in detrainment. Parameter not described in the literature." + +[avogadro_constant] +used_in = ["RRTMGP"] +value = 6.02214076e23 +type = "float" + +[pressure_triple_point] +used_in = ["Thermodynamics"] +value = 611.657 +type = "float" + +[thermodynamics_temperature_reference] +used_in = ["Thermodynamics"] +value = 273.16 +type = "float" + +[temperature_water_freeze] +used_in = ["Thermodynamics"] +value = 273.15 +type = "float" + +[total_solar_irradiance] +used_in = ["Insolation"] +value = 1362 +type = "float" + +[detr_coeff] +used_in = ["ClimaAtmos"] +value = 0.001 +type = "float" +description = "TODO: Remove this. Constant detrainment coefficient used for testing EDMF" + +[minimum_updraft_top] +used_in = ["ClimaAtmos"] +value = 500.0 +type = "float" +description = "Minimum updraft height limiter to avoid zero division in pressure drag and entrainment/detrainment closures. Parameter not described in the literature." + +[von_karman_constant] +used_in = ["SurfaceFluxes"] +value = 0.4 +type = "float" + +[mixing_length_eddy_viscosity_coefficient] +used_in = ["ClimaAtmos"] +value = 0.14 +type = "float" +description = "Turbulence kinetic energy diffusivity coefficient for the EDMF mixing length closure; denoted c_m. See Lopez-Gomez et al. (2020) [https://doi.org/10.1029/2020MS002162], Table 1." + +[alpha_rayleigh_w] +used_in = ["ClimaAtmos"] +value = 1.0 +type = "float" +description = "rayleigh sponge vert velocity coeff" + +[mean_sea_level_pressure] +used_in = ["Thermodynamics"] +value = 101325 +type = "float" + +[max_area_limiter_power] +used_in = ["ClimaAtmos"] +value = 10 +type = "float" +description = "Constant coefficient for the exponent in the maximum area limiter term in detrainment. Parameter not described in the literature." + +[anomalistic_year_length] +used_in = ["Insolation"] +value = 31558464 +type = "float" +description = "derived: 365.25 * [day]" + +[alpha_rayleigh_uh] +used_in = ["ClimaAtmos"] +value = 0.0 +type = "float" +description = "rayleigh sponge horizontal velocity coefficient" + +[kappa_2_sponge] +used_in = ["ClimaAtmos"] +value = 1.0e6 +type = "float" +description = "viscous sponge coefficient" + +[astronomical_unit] +used_in = ["ClimaAtmos"] +value = 149597870000 +type = "float" + +[mixing_length_Prandtl_number_scale] +used_in = ["ClimaAtmos"] +value = 4.076923076923077 +type = "float" +description = "Cospectral budget factor for turbulent Prandtl number for the EDMF mixing length closure, denoted ω_pr. In Lopez-Gomez et al. (2020) [https://doi.org/10.1029/2020MS002162], Eq. 36, it is described as a phenomenological constant, denoted by ω₂ and set to 53/13 ≈ 4.0769..." diff --git a/calibration/generate_observations.sbatch b/calibration/generate_observations.sbatch index 00af05b297..4ef003a605 100644 --- a/calibration/generate_observations.sbatch +++ b/calibration/generate_observations.sbatch @@ -19,8 +19,11 @@ using NCDatasets import JLD2 experiment_dir = joinpath("experiments", "amip_coupled") -output_dir = joinpath(experiment_dir, "truth_simulation") +COUPLER_OUTPUT_DIR = joinpath(experiment_dir, "truth_simulation") +include("coupler_driver_init.jl") +include("coupler_parse_args.jl") +include("coupler_component_init.jl") include("coupler_driver_calibration.jl"); solve_coupler!(cs); # Integrate the coupled model From 1368b9096531bf357f4b269f47d5d7ee2ecff059 Mon Sep 17 00:00:00 2001 From: akshaysridhar Date: Wed, 6 Mar 2024 18:20:44 -0800 Subject: [PATCH 3/3] renamed: ../../ekp_config.yml -> ekp_config.yml modified: ../../generate_observations.sbatch modified: ../../observation_map.jl --- .../amip_coupled}/ekp_config.yml | 0 calibration/generate_observations.sbatch | 3 +- calibration/observation_map.jl | 34 +++++++++++++++++++ 3 files changed, 36 insertions(+), 1 deletion(-) rename calibration/{ => experiments/amip_coupled}/ekp_config.yml (100%) diff --git a/calibration/ekp_config.yml b/calibration/experiments/amip_coupled/ekp_config.yml similarity index 100% rename from calibration/ekp_config.yml rename to calibration/experiments/amip_coupled/ekp_config.yml diff --git a/calibration/generate_observations.sbatch b/calibration/generate_observations.sbatch index 4ef003a605..068be41bba 100644 --- a/calibration/generate_observations.sbatch +++ b/calibration/generate_observations.sbatch @@ -27,7 +27,8 @@ include("coupler_component_init.jl") include("coupler_driver_calibration.jl"); solve_coupler!(cs); # Integrate the coupled model -ta = ncread(joinpath(output_dir, "", "ta_150s_average.nc"), "ta") +testdir = "/Users/akshaysridhar/Research/Codes/ClimaCoupler.jl/calibration/experiments/amip_coupled/truth_simulation/" +wa = NCDataset(joinpath(testdir, "", "wa_inst.nc"))["wa"] include(joinpath(experiment_dir, "observation_map.jl")) (; observation, variance) = process_member_data(ta; output_variance = true) JLD2.save_object(joinpath(experiment_dir, "obs_mean.jld2"), observation) diff --git a/calibration/observation_map.jl b/calibration/observation_map.jl index 71318c8d57..5c2a236a3b 100644 --- a/calibration/observation_map.jl +++ b/calibration/observation_map.jl @@ -1,2 +1,36 @@ ### Place holder for NCEP data from the ClimaCoupler outputs # + +function observation_map(iteration) + experiment_id = "amip_coupled" + config = YAML.load_file(joinpath("experiments", experiment_id, "ekp_config.yml")) + output_dir = config["output_dir"] + ensemble_size = config["ensemble_size"] + model_output = "wa_inst.nc" + dims = 1 + G_ensemble = Array{Float64}(undef, dims..., ensemble_size) + for m in 1:ensemble_size + member_path = + TOMLInterface.path_to_ensemble_member(output_dir, iteration, m) + ta = ncread(joinpath(member_path, model_output), "wa") + G_ensemble[:, m] = process_member_data(ta) + end + return G_ensemble +end + +function process_member_data(wa; output_variance = false) + # Cut off first 120 days to get equilibrium, take second level slice + level_slice = 2 + wa_second_height = wa[3:size(wa)[1], :, :, level_slice] + # Average over long and latitude + area_avg_wa_second_height = + longitudinal_avg(latitudinal_avg(wa_second_height)) + observation = Float64[area_avg_wa_second_height[3]] + if !(output_variance) + return observation + else + variance = Matrix{Float64}(undef, 1, 1) + variance[1] = var(area_avg_wa_second_height) + return (; observation, variance) + end +end \ No newline at end of file