Merge pull request #3620 from CliMA/aj/cloud_effective_radius

Prescribe effective radius, add to diagnostics output

Author: trontrytel

.buildkite/Manifest-v1.11.toml

@@ -2,7 +2,7 @@
 
 julia_version = "1.11.3"
 manifest_format = "2.0"
-project_hash = "5971c6c9602125f1301878afca0ead7416797e54"
+project_hash = "435748ba4bd3e84a898785b2d03f015c42fbae0c"
 
 [[deps.ADTypes]]
 git-tree-sha1 = "e1ce448a0d7f88168ffe2eeac4549c32d45a42d1"
@@ -394,9 +394,9 @@ version = "0.2.12"
 
 [[deps.ClimaParams]]
 deps = ["TOML"]
-git-tree-sha1 = "1e17ed5997da08f1ca186df1882ef10010f70481"
+git-tree-sha1 = "e8f4589e684b5053660ee3981bb2236221d5ad27"
 uuid = "5c42b081-d73a-476f-9059-fd94b934656c"
-version = "0.10.18"
+version = "0.10.19"
 
 [[deps.ClimaReproducibilityTests]]
 deps = ["OrderedCollections", "PrettyTables"]
@@ -436,9 +436,9 @@ version = "0.1.13"
 
 [[deps.CloudMicrophysics]]
 deps = ["ClimaParams", "DocStringExtensions", "ForwardDiff", "HCubature", "LazyArtifacts", "QuadGK", "RootSolvers", "SpecialFunctions", "Thermodynamics"]
-git-tree-sha1 = "ce0a087e03585d89fa925ffdb18246309ad40e55"
+git-tree-sha1 = "98c77aaefaedaa185a7bb676af16ab19aede7ff3"
 uuid = "6a9e3e04-43cd-43ba-94b9-e8782df3c71b"
-version = "0.22.6"
+version = "0.22.8"
 
     [deps.CloudMicrophysics.extensions]
     EmulatorModelsExt = ["DataFrames", "MLJ"]
@@ -1986,9 +1986,9 @@ version = "1.11.0"
 
 [[deps.RRTMGP]]
 deps = ["Adapt", "Artifacts", "ClimaComms", "DocStringExtensions", "Random"]
-git-tree-sha1 = "bf25b44442e3bedd7ea0b487a863608c7e99a8d4"
+git-tree-sha1 = "d258e9495bb6a1b7ad2b28683c9382c0053dafc1"
 uuid = "a01a1ee8-cea4-48fc-987c-fc7878d79da1"
-version = "0.20.0"
+version = "0.20.1"
 weakdeps = ["CUDA", "ClimaParams"]
 
     [deps.RRTMGP.extensions]

Project.toml

@@ -43,10 +43,10 @@ AtmosphericProfilesLibrary = "0.1.7"
 ClimaComms = "0.6.6"
 ClimaCore = "0.14.24"
 ClimaDiagnostics = "0.2.12"
-ClimaParams = "0.10.17"
+ClimaParams = "0.10.19"
 ClimaTimeSteppers = "0.8.2"
 ClimaUtilities = "0.1.22"
-CloudMicrophysics = "0.22.4"
+CloudMicrophysics = "0.22.8"
 Dates = "1"
 DiffEqBase = "6.145"
 FastGaussQuadrature = "0.5, 1"
@@ -59,7 +59,7 @@ Logging = "1"
 NCDatasets = "0.14.2"
 NullBroadcasts = "0.1"
 NVTX = "0.3"
-RRTMGP = "0.19, 0.20"
+RRTMGP = "0.20.1"
 Random = "1"
 SciMLBase = "2.12"
 StaticArrays = "1.7"

config/model_configs/sphere_aquaplanet_rhoe_nonequilmoist_allsky.yml

@@ -15,3 +15,6 @@ insolation: "timevarying"
 rayleigh_sponge: true
 surface_setup: "DefaultMoninObukhov"
 toml: [toml/sphere_aquaplanet_1M.toml]
+diagnostics:
+  - short_name: [reffclw, reffcli]
+    period: 1hours

src/callbacks/callbacks.jl

@@ -1,6 +1,7 @@
 import ClimaCore.DataLayouts as DL
 import .RRTMGPInterface as RRTMGPI
 import Thermodynamics as TD
+import CloudMicrophysics as CM
 import LinearAlgebra
 import ClimaCore.Fields
 import ClimaComms
@@ -15,6 +16,7 @@ import ClimaCore.Fields: ColumnField
 
 import ClimaUtilities.TimeVaryingInputs: evaluate!
 
+
 include("callback_helpers.jl")
 
 function flux_accumulation!(integrator)
@@ -95,6 +97,7 @@ NVTX.@annotate function rrtmgp_model_callback!(integrator)
 
     FT = Spaces.undertype(axes(Y.c))
     thermo_params = CAP.thermodynamics_params(params)
+    cmc = CAP.microphysics_cloud_params(params)
     T_min = CAP.optics_lookup_temperature_min(params)
     T_max = CAP.optics_lookup_temperature_max(params)
 
@@ -183,6 +186,7 @@ NVTX.@annotate function rrtmgp_model_callback!(integrator)
             )
             # RRTMGP needs lwp and iwp in g/m^2
             kg_to_g_factor = 1000
+            m_to_um_factor = FT(1e6)
             cloud_liquid_water_content =
                 radiation_mode.cloud isa PrescribedCloudInRadiation ?
                 p.radiation.prescribed_clouds_field.clwc :
@@ -203,8 +207,18 @@ NVTX.@annotate function rrtmgp_model_callback!(integrator)
                 max(cloud_fraction, eps(FT))
             @. ᶜfrac = cloud_fraction
             # RRTMGP needs effective radius in microns
-            @. ᶜreliq = FT(12)
-            @. ᶜreice = FT(25)
+            @. ᶜreliq = ifelse(
+                cloud_liquid_water_content > FT(0),
+                CM.CloudDiagnostics.effective_radius_const(cmc.liquid) *
+                m_to_um_factor,
+                FT(0),
+            )
+            @. ᶜreice = ifelse(
+                cloud_ice_water_content > FT(0),
+                CM.CloudDiagnostics.effective_radius_const(cmc.ice) *
+                m_to_um_factor,
+                FT(0),
+            )
         end
     end
 

src/diagnostics/radiation_diagnostics.jl

@@ -857,3 +857,88 @@ add_diagnostic_variable!(
     comments = "Upwelling clear-sky longwave radiation at the top of the atmosphere",
     compute! = compute_rlutcs!,
 )
+
+
+###
+# Effective radius for liquid clouds (3d)
+###
+compute_reffclw!(out, state, cache, time) =
+    compute_reffclw!(out, state, cache, time, cache.atmos.radiation_mode)
+compute_reffclw!(_, _, _, _, radiation_mode::T) where {T} =
+    error_diagnostic_variable("reffclw", radiation_mode)
+
+function compute_reffclw!(
+    out,
+    state,
+    cache,
+    time,
+    radiation_mode::Union{
+        RRTMGPI.AllSkyRadiationWithClearSkyDiagnostics,
+        RRTMGPI.AllSkyRadiation,
+    },
+)
+    FT = eltype(state)
+    if isnothing(out)
+        return Fields.array2field(
+            cache.radiation.rrtmgp_model.center_cloud_liquid_effective_radius,
+            axes(state.c),
+        ) .* FT(1e-6) # RRTMGP stores r_eff in microns
+    else
+        out .=
+            Fields.array2field(
+                cache.radiation.rrtmgp_model.center_cloud_liquid_effective_radius,
+                axes(state.c),
+            ) .* FT(1e-6)
+    end
+end
+
+add_diagnostic_variable!(
+    short_name = "reffclw",
+    long_name = "Effective radius for liquid clouds",
+    standard_name = "effective_radius_of_cloud_liquid_particles",
+    units = "m",
+    comments = "In-cloud ratio of the third moment over the second moment of the particle size distribution. Set to zero outside of clouds.",
+    compute! = compute_reffclw!,
+)
+
+###
+# Effective radius for ice clouds (3d)
+###
+compute_reffcli!(out, state, cache, time) =
+    compute_reffcli!(out, state, cache, time, cache.atmos.radiation_mode)
+compute_reffcli!(_, _, _, _, radiation_mode::T) where {T} =
+    error_diagnostic_variable("reffcli", radiation_mode)
+
+function compute_reffcli!(
+    out,
+    state,
+    cache,
+    time,
+    radiation_mode::Union{
+        RRTMGPI.AllSkyRadiationWithClearSkyDiagnostics,
+        RRTMGPI.AllSkyRadiation,
+    },
+)
+    FT = eltype(state)
+    if isnothing(out)
+        return Fields.array2field(
+            cache.radiation.rrtmgp_model.center_cloud_ice_effective_radius,
+            axes(state.c),
+        ) .* FT(1e-6) # RRTMGP stores r_eff in microns
+    else
+        out .=
+            Fields.array2field(
+                cache.radiation.rrtmgp_model.center_cloud_ice_effective_radius,
+                axes(state.c),
+            ) .* FT(1e-6)
+    end
+end
+
+add_diagnostic_variable!(
+    short_name = "reffcli",
+    long_name = "Effective radius for ice clouds",
+    standard_name = "effective_radius_of_cloud_ice_particles",
+    units = "m",
+    comments = "In-cloud ratio of the third moment over the second moment of the particle size distribution. Set to zero outside of clouds.",
+    compute! = compute_reffcli!,
+)

toml/sphere_aquaplanet.toml

@@ -6,3 +6,6 @@ value = 40000.0
 
 [precipitation_timescale]
 value = 400
+
+[liquid_cloud_effective_radius]
+value = 12e-6

toml/sphere_aquaplanet_1M.toml

@@ -9,3 +9,6 @@ value = 40000.0
 
 [zd_rayleigh]
 value = 40000.0
+
+[liquid_cloud_effective_radius]
+value = 12e-6