Add TOA and surface radiation diagnostics

src/diagnostics/default_diagnostics.jl

@@ -232,7 +232,7 @@ end
 # Radiation mode #
 ##################
 function default_diagnostics(::RRTMGPI.AbstractRRTMGPMode; output_writer)
-    rad_diagnostics = ["rsd", "rsu", "rld", "rlu"]
+    rad_diagnostics = ["rsd", "rsdt", "rsds", "rsu", "rld", "rlu"]
 
     return [daily_averages(rad_diagnostics...; output_writer)...]
 end
@@ -242,8 +242,18 @@ function default_diagnostics(
     ::RRTMGPI.AllSkyRadiationWithClearSkyDiagnostics;
     output_writer,
 )
-    rad_diagnostics =
-        ["rsd", "rsu", "rld", "rlu", "rsdcs", "rsucs", "rldcs", "rlucs"]
+    rad_diagnostics = [
+        "rsd",
+        "rsdt",
+        "rsds",
+        "rsu",
+        "rld",
+        "rlu",
+        "rsdcs",
+        "rsucs",
+        "rldcs",
+        "rlucs",
+    ]
 
     return [daily_averages(rad_diagnostics...; output_writer)...]
 end

src/diagnostics/radiation_diagnostics.jl

@@ -39,6 +39,93 @@ add_diagnostic_variable!(
     compute! = compute_rsd!,
 )
 
+###
+# TOA downwelling shortwave radiation (3d)
+###
+compute_rsdt!(out, state, cache, time) =
+    compute_rsdt!(out, state, cache, time, cache.atmos.radiation_mode)
+compute_rsdt!(_, _, _, _, radiation_mode::T) where {T} =
+    error_diagnostic_variable("rsdt", radiation_mode)
+
+function compute_rsdt!(
+    out,
+    state,
+    cache,
+    time,
+    radiation_mode::T,
+) where {T <: RRTMGPI.AbstractRRTMGPMode}
+    nlevels = Spaces.nlevels(axes(state.c))
+    if isnothing(out)
+        return Fields.level(
+            RRTMGPI.array2field(
+                cache.radiation.radiation_model.face_sw_flux_dn,
+                axes(state.f),
+            ),
+            nlevels + half,
+        )
+    else
+        out .= Fields.level(
+            RRTMGPI.array2field(
+                cache.radiation.radiation_model.face_sw_flux_dn,
+                axes(state.f),
+            ),
+            nlevels + half,
+        )
+    end
+end
+
+add_diagnostic_variable!(
+    short_name = "rsdt",
+    long_name = "TOA Incident Shortwave Radiation",
+    standard_name = "toa_incoming_shortwave_flux",
+    units = "W m^-2",
+    comments = "Shortwave radiation incident at the top of the atmosphere",
+    compute! = compute_rsdt!,
+)
+
+###
+# Surface downwelling shortwave radiation (3d)
+###
+compute_rsds!(out, state, cache, time) =
+    compute_rsds!(out, state, cache, time, cache.atmos.radiation_mode)
+compute_rsds!(_, _, _, _, radiation_mode::T) where {T} =
+    error_diagnostic_variable("rsds", radiation_mode)
+
+function compute_rsds!(
+    out,
+    state,
+    cache,
+    time,
+    radiation_mode::T,
+) where {T <: RRTMGPI.AbstractRRTMGPMode}
+    if isnothing(out)
+        return Fields.level(
+            RRTMGPI.array2field(
+                cache.radiation.radiation_model.face_sw_flux_dn,
+                axes(state.f),
+            ),
+            half,
+        )
+    else
+        out .= Fields.level(
+            RRTMGPI.array2field(
+                cache.radiation.radiation_model.face_sw_flux_dn,
+                axes(state.f),
+            ),
+            half,
+        )
+    end
+end
+
+add_diagnostic_variable!(
+    short_name = "rsds",
+    long_name = "Surface Downwelling Shortwave Radiation",
+    standard_name = "surface_downwelling_shortwave_flux_in_air",
+    units = "W m^-2",
+    comments = "Surface solar irradiance for UV calculations.",
+    compute! = compute_rsds!,
+)
+
 ###
 # Upwelling shortwave radiation (3d)
 ###

src/solver/type_getters.jl

@@ -940,8 +940,14 @@ function get_simulation(config::AtmosConfig)
                     )
                 else
                     # Add to the accumulator
+
+                    # We use similar + .= instead of copy because of a bug in CUDA.jl 5.2.0
+                    # with view(reshape(view)) objects. See discussion in
+                    # https://github.com/CliMA/ClimaAtmos.jl/pull/2579
                     diagnostic_accumulators[diag] =
-                        copy(diagnostic_storage[diag])
+                        similar(diagnostic_storage[diag])
+                    diagnostic_accumulators[diag] .=
+                        diagnostic_storage[diag]
                 end
             catch e
                 error("Could not compute diagnostic $(variable.long_name): $e")