Merge pull request #3775 from CliMA/rc/nogw

Support for intermittent calls to update non-orographic gravity wave tendencies (Replacement for #3761)

Author: akshaysridhar

NEWS.md

@@ -7,6 +7,10 @@ main
 ### Add support for reanalysis-driven single column model with time-varying forcing
 PR [#3758](https://github.com/CliMA/ClimaAtmos.jl/pull/3758) adds support for driving single-column model (SCM) simulations with time-varying ERA5 reanalysis data. This extends the existing GCM-driven SCM interface to allow site-specific simulations that resolve the diurnal cycle and are suited for calibration against observations. Users can now run reanalysis-driven cases globally using only a date and lat/lon, thanks to integrated data handling via ClimaArtifacts.jl. See the updated “Single Column Model” docs page for details on setup, variable requirements, and how to prepare ERA5 input files.
 
+### Non-orographic gravity wave tendency as a callback
+
+PR[#3761](https://github.com/CliMA/ClimaAtmos.jl/pull/3761) introduces support for intermittent calls to update the computation of non-orographic gravity wave tendencies. This PR closes issue[#3434](https://github.com/CliMA/ClimaAtmos.jl/issues/3434).
+
 ### Remove `dt_save_to_sol`
 
 The option to save the solution to the integrator object (`dt_save_to_sol`) was

config/default_configs/default_config.yml

@@ -220,6 +220,9 @@ orographic_gravity_wave:
 non_orographic_gravity_wave:
   help: "Apply parameterization for convective gravity wave forcing on horizontal mean flow [`false` (default), `true`]"
   value: false
+dt_nogw:
+  help: "Time between calling non-orographic gravity wave update"
+  value: "1800secs"
 nh_poly:
   help: "Horizontal polynomial degree. Note: The number of quadrature points in 1D within each horizontal element is then Nq = <--nh_poly> + 1"
   value: 3

config/model_configs/single_column_nonorographic_gravity_wave.yml

@@ -6,3 +6,5 @@ hyperdiff: false
 dt: "400secs"
 disable_surface_flux_tendency: true
 non_orographic_gravity_wave: true
+dt_nogw: "400secs"
+

config/model_configs/sphere_nonorographic_gravity_wave.yml

@@ -1,4 +1,8 @@
-t_end: "1500secs"
-dt: "400secs"
+t_end: "48hours"
+dt: "450secs"
 disable_surface_flux_tendency: true
 non_orographic_gravity_wave: true
+dt_nogw: "450secs"
+diagnostics:
+  - short_name: [ta, thetaa, ua, va, wa, zg, utendnogw, vtendnogw]
+    period: 1hours

docs/src/diagnostics.md

@@ -209,11 +209,6 @@ were computed differently for `EquilMoistModel` and `NonEquilMoistModel`.
 In `ClimaAtmos`, we define some helper functions to produce error messages, so
 the above code can be written as
 ```julia
-compute_relative_humidity!(out, state, cache, time) =
-    compute_relative_humidity!(out, state, cache, time, cache.atmos.moisture_model)
-compute_relative_humidity!(_, _, _, _, model::T) where {T} =
-    error_diagnostic_variable("relative_humidity", model)
-
 function compute_relative_humidity!(
     out,
     state,
@@ -228,6 +223,11 @@ function compute_relative_humidity!(
         out .= TD.relative_humidity.(thermo_params, cache.ᶜts)
     end
 end
+
+compute_relative_humidity!(out, state, cache, time) =
+    compute_relative_humidity!(out, state, cache, time, cache.atmos.moisture_model)
+compute_relative_humidity!(_, _, _, _, model) =
+    error_diagnostic_variable("relative_humidity", model)
 ```
 
 ### The `ClimaAtmos` `DiagnosticVariable`s

docs/src/gravity_wave.md

@@ -165,8 +165,9 @@ Similar to the base flux calculation but for the 3D fields, we computed ``N`` an
 ^f N[k]^2 = \frac{g}{^f T[k]} * ( ^f \overline{\frac{dT}{dz}}[k] + \frac{g}{cp}),
 ```
 ```math
-^f V_{\tau}[k] = \max(\epsilon_0, - V[k] \cdot \frac{\tau}{|\tau|}).
+^f V_{\tau}[k] = \max(\epsilon_0, - V[k] \cdot \frac{\tau}{|\tau|}),
 ```
+where ``\epsilon_0`` denotes a measure of floating-point precision.
 
 Let ``L_1 = L_0 * \max(0.5, \min(2.0, 1.0-samp*V_\tau*d^2V_{\tau}/N^2))`` where ``samp=1.0`` is the correction for coarse sampling of ``d^2V/dz^2``, and
 ```math

reproducibility_tests/ref_counter.jl

@@ -1,4 +1,4 @@
-234
+235
 
 # **README**
 #
@@ -20,6 +20,10 @@
 
 
 #=
+235
+- Related to #3775, the computation and update of non-orographic gravity wave (NOGW) 
+  are now separated into callback and tendencies update, affecting the NOGW-related
+  tests.
 
 234
 - Move the virtual mass term in pressure closure to prognostic edmf momentum equation

src/callbacks/callbacks.jl

@@ -399,6 +399,14 @@ NVTX.@annotate function rrtmgp_model_callback!(integrator)
     return nothing
 end
 
+NVTX.@annotate function nogw_model_callback!(integrator)
+    Y = integrator.u
+    p = integrator.p
+
+    non_orographic_gravity_wave_compute_tendency!(Y, p)
+    return nothing
+end
+
 #Uniform insolation, magnitudes from Wing et al. (2018)
 #Note that the TOA downward shortwave fluxes won't be the same as the values in the paper if add_isothermal_boundary_layer is true
 function set_insolation_variables!(Y, p, t, ::RCEMIPIIInsolation)

src/callbacks/get_callbacks.jl

@@ -365,5 +365,21 @@ function get_callbacks(config, sim_info, atmos, params, Y, p)
             (callbacks..., call_every_dt(rrtmgp_model_callback!, dt_rad))
     end
 
+    if atmos.non_orographic_gravity_wave isa NonOrographicGravityWave
+        dt_nogw =
+            dt isa ITime ? ITime(time_to_seconds(parsed_args["dt_nogw"])) :
+            FT(time_to_seconds(parsed_args["dt_nogw"]))
+        dt_nogw, _, _, _ = promote(dt_nogw, t_start, dt, sim_info.t_end)
+        # We use Millisecond to support fractional seconds, eg. 0.1
+        dt_nogw_ms = Dates.Millisecond(1_000 * float(dt_nogw))
+        if parsed_args["dt_save_state_to_disk"] != "Inf" &&
+           !CA.isdivisible(dt_save_state_to_disk_dates, dt_nogw_ms)
+            @warn "Non-orographic gravity wave period ($(dt_nogw_ms)) is not an even divisor of the checkpoint frequency ($dt_save_state_to_disk_dates)"
+            @warn "This simulation will not be reproducible when restarted"
+        end
+
+        callbacks = (callbacks..., call_every_dt(nogw_model_callback!, dt_nogw))
+    end
+
     return callbacks
 end

src/diagnostics/Diagnostics.jl

@@ -41,6 +41,10 @@ import ..EDOnlyEDMFX
 import ..PrognosticEDMFX
 import ..DiagnosticEDMFX
 
+# gravitywave_models
+import ..NonOrographicGravityWave
+import ..OrographicGravityWave
+
 # surface_model
 import ..PrognosticSurfaceTemperature