diff --git a/docs/src/tracers.md b/docs/src/tracers.md index ec737822013..bbcced5a931 100644 --- a/docs/src/tracers.md +++ b/docs/src/tracers.md @@ -1,4 +1,18 @@ -# Aerosols +# Tracers + +## Aerosols + +`ClimaAtmos` can read aerosol profiles from the CIMP5 forcing data. The data is +described in +[`ClimaArtifacts`](https://github.com/CliMA/ClimaArtifacts/tree/main/aerosol_concentrations) +and it is only relevant for the radiation transfer, and only when RRTMGP is +used. + +To enable this mode, use the `PrescribedCMIP5Aerosols` aerosols model. Upon +construction, `PrescribedCMIP5Aerosols` takes the list of aerosol names to be +used. The names have to match what is in the file. + +## Ozone `ClimaAtmos` implements two modes for ozone profiles. These are only relevant for the radiation transfer, and only when RRTMGP is used. diff --git a/src/cache/cache.jl b/src/cache/cache.jl index abf60dd48dc..7fd97d863c7 100644 --- a/src/cache/cache.jl +++ b/src/cache/cache.jl @@ -82,7 +82,7 @@ end # The model also depends on f_plane_coriolis_frequency(params) # This is a constant Coriolis frequency that is only used if space is flat -function build_cache(Y, atmos, params, surface_setup, sim_info, aerosol_names) +function build_cache(Y, atmos, params, surface_setup, sim_info) (; dt, start_date, output_dir) = sim_info FT = eltype(params) @@ -147,7 +147,7 @@ function build_cache(Y, atmos, params, surface_setup, sim_info, aerosol_names) radiation_args = atmos.radiation_mode isa RRTMGPI.AbstractRRTMGPMode ? - (start_date, params, atmos.ozone, aerosol_names, atmos.insolation) : () + (start_date, params, atmos.ozone, atmos.aerosols, atmos.insolation) : () hyperdiff = hyperdiffusion_cache(Y, atmos) precipitation = precipitation_cache(Y, atmos) @@ -157,7 +157,7 @@ function build_cache(Y, atmos, params, surface_setup, sim_info, aerosol_names) non_orographic_gravity_wave = non_orographic_gravity_wave_cache(Y, atmos) orographic_gravity_wave = orographic_gravity_wave_cache(Y, atmos) radiation = radiation_model_cache(Y, atmos, radiation_args...) - tracers = tracer_cache(Y, atmos, aerosol_names, start_date) + tracers = tracer_cache(Y, atmos, start_date) args = ( dt, diff --git a/src/cache/tracer_cache.jl b/src/cache/tracer_cache.jl index 3c4ef831df0..83dd68667fc 100644 --- a/src/cache/tracer_cache.jl +++ b/src/cache/tracer_cache.jl @@ -23,53 +23,51 @@ function ozone_cache(::PrescribedOzone, Y, start_date) return (; o3, prescribed_o3_timevaryinginput) end -function tracer_cache(Y, atmos, prescribed_aerosol_names, start_date) - if !isempty(prescribed_aerosol_names) - target_space = axes(Y.c) +aerosols_cache(_, _, _) = (;) +function aerosols_cache(aerosols::PrescribedCMIP5Aerosols, Y, start_date) + aerosol_names = aerosols.aerosol_names + target_space = axes(Y.c) - # Take the aerosol concentration file, read the keys with names matching - # the ones passed in the prescribed_aerosol_names option, and create a - # NamedTuple that uses the same keys and has as values the TimeVaryingInputs - # for those variables. - # - # The keys in the aerosol_concentrations.nc file have to match the ones passed with the - # configuration. The file also has to be defined on the globe and provide - # time series of lon-lat-z data. - prescribed_aerosol_names_as_symbols = Symbol.(prescribed_aerosol_names) - target_space = axes(Y.c) - extrapolation_bc = (Intp.Periodic(), Intp.Flat(), Intp.Flat()) - timevaryinginputs = [ - TimeVaryingInput( - AA.aerosol_concentration_file_path(; - context = ClimaComms.context(Y.c), - ), - name, - target_space; - reference_date = start_date, - regridder_type = :InterpolationsRegridder, - regridder_kwargs = (; extrapolation_bc), - method = LinearPeriodFillingInterpolation(Year(1)), - ) for name in prescribed_aerosol_names - ] + # Take the aerosol concentration file, read the keys with names matching + # the ones passed in the prescribed_aerosol_names option, and create a + # NamedTuple that uses the same keys and has as values the TimeVaryingInputs + # for those variables. + # + # The keys in the aerosol_concentrations.nc file have to match the ones passed with the + # configuration. The file also has to be defined on the globe and provide + # time series of lon-lat-z data. + prescribed_aerosol_names_as_symbols = Symbol.(prescribed_aerosol_names) + target_space = axes(Y.c) + extrapolation_bc = (Intp.Periodic(), Intp.Flat(), Intp.Flat()) + timevaryinginputs = [ + TimeVaryingInput( + AA.aerosol_concentration_file_path(; + context = ClimaComms.context(Y.c), + ), + name, + target_space; + reference_date = start_date, + regridder_type = :InterpolationsRegridder, + regridder_kwargs = (; extrapolation_bc), + method = LinearPeriodFillingInterpolation(Year(1)), + ) for name in prescribed_aerosol_names + ] + + # Field is updated in the radiation callback + prescribed_aerosols_field = similar( + Y.c, + NamedTuple{ + prescribed_aerosol_names_as_symbols, + NTuple{length(prescribed_aerosol_names_as_symbols), eltype(Y.c.ρ)}, + }, + ) + prescribed_aerosol_timevaryinginputs = + (; zip(prescribed_aerosol_names_as_symbols, timevaryinginputs)...) + return (; prescribed_aerosols_field, prescribed_aerosol_timevaryinginputs) +end - # Field is updated in the radiation callback - prescribed_aerosols_field = similar( - Y.c, - NamedTuple{ - prescribed_aerosol_names_as_symbols, - NTuple{ - length(prescribed_aerosol_names_as_symbols), - eltype(Y.c.ρ), - }, - }, - ) - prescribed_aerosol_timevaryinginputs = - (; zip(prescribed_aerosol_names_as_symbols, timevaryinginputs)...) - aerosol_cache = - (; prescribed_aerosols_field, prescribed_aerosol_timevaryinginputs) - else - aerosol_cache = (;) - end +function tracer_cache(Y, atmos, start_date) + aero_cache = aerosols_cache(atmos.aerosols, Y, start_date) o3_cache = ozone_cache(atmos.ozone, Y, start_date) - return (; aerosol_cache..., o3_cache...) + return (; aero_cache..., o3_cache...) end diff --git a/src/callbacks/callbacks.jl b/src/callbacks/callbacks.jl index 1334d397ae3..45893417e84 100644 --- a/src/callbacks/callbacks.jl +++ b/src/callbacks/callbacks.jl @@ -187,6 +187,9 @@ NVTX.@annotate function rrtmgp_model_callback!(integrator) if !(radiation_mode isa RRTMGPI.GrayRadiation) if radiation_mode.aerosol_radiation + # This is because we hard-code the names of the symbols + @assert integrator.p.atmos.aerosols isa PrescribedCMIP5Aerosols "The only aerosol model currently supported is `PrescribedCMIP5Aerosols`" + ᶜΔz = Fields.Δz_field(Y.c) ᶜaero_conc = Fields.array2field( diff --git a/src/parameterized_tendencies/radiation/radiation.jl b/src/parameterized_tendencies/radiation/radiation.jl index 8804243bc08..1db46b18d3b 100644 --- a/src/parameterized_tendencies/radiation/radiation.jl +++ b/src/parameterized_tendencies/radiation/radiation.jl @@ -92,7 +92,7 @@ function radiation_model_cache( start_date, params, ozone, - aerosol_names, + aerosols, insolation_mode; interpolation = RRTMGPI.BestFit(), bottom_extrapolation = RRTMGPI.SameAsInterpolation(), @@ -102,11 +102,13 @@ function radiation_model_cache( device = context.device if !(radiation_mode isa RRTMGPI.GrayRadiation) (; aerosol_radiation) = radiation_mode - if aerosol_radiation && !(any( - x -> x in aerosol_names, - ["DST01", "SSLT01", "SO4", "CB1", "CB2", "OC1", "OC2"], - )) - error( + if aerosol_radiation + aerosols isa PrescribedCIMP5Aerosols || + error("Only `PrescribedCIMP5Aerosols` is currently supported") + any( + x -> x in aerosols.aerosol_names, + ["DST01", "SSLT01", "SO4", "CB1", "CB2", "OC1", "OC2"], + ) && error( "Need at least one aerosol type when aerosol radiation is turned on", ) end diff --git a/src/solver/model_getters.jl b/src/solver/model_getters.jl index 53c15833b99..9ea834d82aa 100644 --- a/src/solver/model_getters.jl +++ b/src/solver/model_getters.jl @@ -498,9 +498,10 @@ function get_surface_thermo_state_type(parsed_args) return dict[parsed_args["surface_thermo_state_type"]] end -function get_tracers(parsed_args) +function get_aerosols(parsed_args) + isnothing(parsed_args["prescribed_aerosols"]) && return nothing aerosol_names = Tuple(parsed_args["prescribed_aerosols"]) - return (; aerosol_names) + return PrescribedCMIP5Aerosols(aerosol_names) end function get_tendency_model(parsed_args) diff --git a/src/solver/type_getters.jl b/src/solver/type_getters.jl index 61b651492b4..35debd3b003 100644 --- a/src/solver/type_getters.jl +++ b/src/solver/type_getters.jl @@ -20,6 +20,7 @@ function get_atmos(config::AtmosConfig, params) precip_model = get_precipitation_model(parsed_args) cloud_model = get_cloud_model(parsed_args) ozone = get_ozone(parsed_args) + aerosols = get_aerosols(parsed_args) radiation_mode = get_radiation_mode(parsed_args, FT) forcing_type = get_forcing_type(parsed_args) call_cloud_diagnostics_per_stage = @@ -30,6 +31,10 @@ function get_atmos(config::AtmosConfig, params) ozone = IdealizedOzone() end + if !isnothing(aerosols) && !(radiation_mode isa RRTMGPI.AbstractRRTMGPMode) + @warn "$aerosols is effective only with an RRTGMP model" + end + diffuse_momentum = !(forcing_type isa HeldSuarezForcing) advection_test = parsed_args["advection_test"] @@ -56,6 +61,7 @@ function get_atmos(config::AtmosConfig, params) atmos = AtmosModel(; moisture_model, ozone, + aerosols, radiation_mode, subsidence = get_subsidence_model(parsed_args, radiation_mode, FT), ls_adv = get_large_scale_advection_model(parsed_args, FT), @@ -680,14 +686,7 @@ function get_simulation(config::AtmosConfig) tracers = get_tracers(config.parsed_args) s = @timed_str begin - p = build_cache( - Y, - atmos, - params, - surface_setup, - sim_info, - tracers.aerosol_names, - ) + p = build_cache(Y, atmos, params, surface_setup, sim_info) end @info "Allocating cache (p): $s" diff --git a/src/solver/types.jl b/src/solver/types.jl index 6059303dc09..61915ee4f32 100644 --- a/src/solver/types.jl +++ b/src/solver/types.jl @@ -103,6 +103,8 @@ struct GCMDrivenInsolation <: AbstractInsolation end AbstractOzone Describe how ozone concentration should be set. + +Currently, this is only relevant for RRTGMP. """ abstract type AbstractOzone end @@ -127,6 +129,33 @@ Refer to ClimaArtifacts for more information on how to obtain the artifact. """ struct PrescribedOzone <: AbstractOzone end +""" + AbstractAerosols + +Describe how aerosol concentrations should be set. + +Currently, this is only relevant for RRTGMP. +""" +abstract type AbstractAerosols end + +""" + PrescribedCMIP5Aerosols + +Implement a time-varying ozone profile as read from disk. + +The CMIP5 forcing dataset is used. For production runs, you should acquire the +high-resolution, multi-year `aerosol_concentrations` artifact. If this is not +available, a low resolution, single-year version will be used. + +Refer to ClimaArtifacts for more information on how to obtain the artifact. + +`aerosol_names` is the list of aerosols to use. The names have to match what is +in the artifact. +""" +struct PrescribedCMIP5Aerosols{T} <: AbstractAerosols + aerosol_names::T +end + """ AbstractCloudInRadiation @@ -456,6 +485,7 @@ Base.@kwdef struct AtmosModel{ F, S, OZ, + AERO, RM, LA, EXTFORCING, @@ -489,6 +519,9 @@ Base.@kwdef struct AtmosModel{ """What to do with ozone for radiation (when using RRTGMP)""" ozone::OZ = nothing + """What to do with aerosols for radiation (when using RRTGMP)""" + aerosol::AERO = nothing + radiation_mode::RM = nothing ls_adv::LA = nothing external_forcing::EXTFORCING = nothing