merge main

.buildkite/pipeline.yml

@@ -144,6 +144,10 @@ steps:
           julia --color=yes --project=examples examples/hybrid/driver.jl
           --config_file $CONFIG_PATH/single_column_hydrostatic_balance_ft64.yml
           --job_id single_column_hydrostatic_balance_ft64
+
+          julia --color=yes --project=examples reproducibility_tests/test_mse.jl
+          --job_id single_column_hydrostatic_balance_ft64
+          --out_dir single_column_hydrostatic_balance_ft64/output_active
         artifact_paths: "single_column_hydrostatic_balance_ft64/output_active/*"
 
   - group: "Box Examples"
@@ -420,7 +424,7 @@ steps:
           --config_file $CONFIG_PATH/sphere_ssp_baroclinic_wave_rhoe_equilmoist_earth.yml
           --job_id sphere_ssp_baroclinic_wave_rhoe_equilmoist_earth
         artifact_paths: "sphere_ssp_baroclinic_wave_rhoe_equilmoist_earth/output_active/*"
-      
+
       - label: ":computer: Baroclinic wave (ρe) equilmoist (Hughes2023 double mountain config)"
         command: >
           julia --color=yes --project=examples examples/hybrid/driver.jl
@@ -974,7 +978,7 @@ steps:
         agents:
           slurm_gpus: 1
           slurm_mem: 16G
-      
+
       - label: "GPU: test DYAMOND interpolated initial conditions"
         command: >
           julia --color=yes --project=examples examples/hybrid/driver.jl

config/default_configs/default_config.yml

@@ -150,7 +150,7 @@ smagorinsky_lilly:
   value: false
 bubble:
   help: "Enable bubble correction for more accurate surface areas"
-  value: true
+  value: false
 start_date:
   help: "Start date of the simulation"
   value: "20100101"
@@ -259,6 +259,9 @@ implicit_diffusion:
 approximate_linear_solve_iters:
   help: "Number of iterations for the approximate linear solve (used when `implicit_diffusion` is true)"
   value: 1
+override_precip_timescale:
+  help: "If true, sets τ_precip to dt. Otherwise, τ_precip is set to the value in the toml dictionary"
+  value: true
 output_default_diagnostics:
   help: "Output the default diagnostics associated to the selected atmospheric model"
   value: true

config/longrun_configs/amip_target_diagedmf.yml

@@ -9,6 +9,7 @@ dt_save_state_to_disk: "20days"
 moist: "equil"
 cloud_model: "quadrature_sgs"
 precip_model: "0M"
+override_precip_timescale: false
 rad: "allskywithclear"
 dt_rad: "1hours"
 dt_cloud_fraction: "1hours"

config/longrun_configs/longrun_aquaplanet_allsky_diagedmf_0M.yml

@@ -25,4 +25,5 @@ edmfx_nh_pressure: true
 edmfx_sgs_mass_flux: true
 edmfx_sgs_diffusive_flux: true
 precip_model: "0M"
+override_precip_timescale: false
 toml: [toml/longrun_aquaplanet_diagedmf.toml]

config/longrun_configs/longrun_aquaplanet_allsky_progedmf_diffonly_0M.yml

@@ -27,4 +27,5 @@ edmfx_filter: true
 edmfx_sgs_mass_flux: false
 edmfx_sgs_diffusive_flux: true
 precip_model: "0M"
+override_precip_timescale: false
 toml: [toml/longrun_aquaplanet_progedmf.toml]

config/model_configs/aquaplanet_rhoe_equil_clearsky_tvinsol_0M_slabocean.yml

@@ -12,4 +12,3 @@ insolation: "timevarying"
 dt_rad: "1hours"
 prognostic_surface: "PrognosticSurfaceTemperature"
 check_conservation: true
-bubble: false

config/model_configs/aquaplanet_rhoe_equil_clearsky_tvinsol_0M_slabocean_ft64.yml

@@ -12,4 +12,3 @@ insolation: "timevarying"
 dt_rad: "1hours"
 prognostic_surface: "PrognosticSurfaceTemperature"
 check_conservation: true
-bubble: false

config/model_configs/diagnostic_edmfx_trmm_box_0M.yml

@@ -15,6 +15,7 @@ moist: equil
 cloud_model: "quadrature_sgs"
 call_cloud_diagnostics_per_stage: true
 precip_model: "0M"
+override_precip_timescale: false
 config: box
 x_max: 1e8
 y_max: 1e8

config/model_configs/diagnostic_edmfx_trmm_stretched_box.yml

@@ -16,6 +16,7 @@ moist: equil
 cloud_model: "quadrature_sgs"
 call_cloud_diagnostics_per_stage: true
 precip_model: "1M"
+override_precip_timescale: false
 config: box
 x_max: 1e8
 y_max: 1e8

config/model_configs/rcemipii_box_diagnostic_edmfx.yml

@@ -16,6 +16,7 @@ edmfx_sgs_diffusive_flux: true
 rayleigh_sponge: true
 moist: equil
 precip_model: 0M
+override_precip_timescale: false
 dt: 30secs
 t_end: 3600secs
 dt_save_state_to_disk: 12hours

config/model_configs/rcemipii_sphere_diagnostic_edmfx.yml

@@ -15,6 +15,7 @@ edmfx_sgs_diffusive_flux: true
 rayleigh_sponge: true
 moist: equil
 precip_model: 0M
+override_precip_timescale: false
 dt: 100secs
 t_end: 12hours
 dt_save_state_to_disk: 12hours

config/model_configs/single_column_hydrostatic_balance_ft64.yml

@@ -7,4 +7,5 @@ config: "column"
 hyperdiff: false
 dt: "3hours"
 FLOAT_TYPE: "Float64"
+reproducibility_test: true
 dt_save_to_sol: "3hours"

config/model_configs/sphere_aquaplanet_rhoe_equilmoist_allsky_gw_res.yml

@@ -2,7 +2,7 @@ z_max: 60000.0
 z_elem: 31
 dz_bottom: 50.0
 rayleigh_sponge: true
-dt: "400secs"
+dt: "360secs"
 t_end: "1days"
 dt_save_state_to_disk: "24hours"
 vert_diff: "FriersonDiffusion"

reproducibility_tests/ref_counter.jl

@@ -1,4 +1,4 @@
-195
+198
 
 # **README**
 #
@@ -20,6 +20,15 @@
 
 
 #=
+198
+- Added terms to the implicit solver that result in changes in the
+aquaplanet (ρe_tot) equil allsky monin_obukhov varying insol gravity wave (gfdl_restart) high top 1-moment
+
+197
+- Added single column hydrostatic balance reproducibility test
+
+196
+- Set bubble correction to false as default
 
 195
 - Use `vanleer_limiter` as default.
@@ -50,7 +59,7 @@
 - Change the model top for a few ci cases
 
 186
-- Topography dataset has been modified to the 60 arc-second ETOPO2022 dataset. 
+- Topography dataset has been modified to the 60 arc-second ETOPO2022 dataset.
   This is behaviour changing for the gravity-wave (raw-topo) parameterization
   when computing `hmax` and `T tensor`.
 

reproducibility_tests/reproducibility_test_job_ids.jl

@@ -5,4 +5,5 @@ reproducibility_test_job_ids() = [
     "sphere_aquaplanet_rhoe_equilmoist_allsky_gw_raw_zonallyasymmetric",
     "deep_sphere_baroclinic_wave_rhoe_equilmoist",
     "diagnostic_edmfx_aquaplanet",
+    "single_column_hydrostatic_balance_ft64",
 ]

src/cache/precipitation_precomputed_quantities.jl

@@ -40,3 +40,23 @@ function set_precipitation_precomputed_quantities!(Y, p, t)
     )
     return nothing
 end
+
+
+"""
+    set_sedimentation_precomputed_quantities!(Y, p, t)
+Updates the sedimentation terminal velocity stored in `p`
+for the non-equilibrium microphysics scheme
+"""
+function set_sedimentation_precomputed_quantities!(Y, p, t)
+    @assert (p.atmos.moisture_model isa NonEquilMoistModel)
+
+    (; ᶜwₗ, ᶜwᵢ) = p.precomputed
+
+    FT = eltype(Y)
+
+    # compute the precipitation terminal velocity [m/s]
+    # TODO - the actual parameterization will be added in the next PR
+    @. ᶜwₗ = FT(0)
+    @. ᶜwᵢ = FT(0)
+    return nothing
+end

src/cache/precomputed_quantities.jl

@@ -48,6 +48,8 @@ function precomputed_quantities(Y, atmos)
         ᶜspecific = specific_gs.(Y.c),
         ᶜu = similar(Y.c, C123{FT}),
         ᶠu³ = similar(Y.f, CT3{FT}),
+        ᶜwₜqₜ = similar(Y.c, Geometry.WVector{FT}),
+        ᶜwₕhₜ = similar(Y.c, Geometry.WVector{FT}),
         ᶜK = similar(Y.c, FT),
         ᶜts = similar(Y.c, TST),
         ᶜp = similar(Y.c, FT),
@@ -59,6 +61,9 @@ function precomputed_quantities(Y, atmos)
     )
     cloud_diagnostics_tuple =
         similar(Y.c, @NamedTuple{cf::FT, q_liq::FT, q_ice::FT})
+    sedimentation_quantities =
+        atmos.moisture_model isa NonEquilMoistModel ?
+        (; ᶜwₗ = similar(Y.c, FT), ᶜwᵢ = similar(Y.c, FT)) : (;)
     moisture_sgs_quantities =
         atmos.moisture_model isa NonEquilMoistModel ?
         (;
@@ -183,6 +188,7 @@ function precomputed_quantities(Y, atmos)
         advective_sgs_quantities...,
         diagnostic_sgs_quantities...,
         vert_diff_quantities...,
+        sedimentation_quantities...,
         precipitation_quantities...,
         cloud_diagnostics_tuple,
         smagorinsky_lilly_quantities...,
@@ -526,9 +532,30 @@ NVTX.@annotate function set_precomputed_quantities!(Y, p, t)
     #     (; ᶜmixing_length) = p.precomputed
     #     compute_gm_mixing_length!(ᶜmixing_length, Y, p)
     # end
-
+    (; ᶜwₜqₜ, ᶜwₕhₜ) = p.precomputed
+    @. ᶜwₜqₜ = Geometry.WVector(0)
+    @. ᶜwₕhₜ = Geometry.WVector(0)
+    #
+    # TODO - uncomment in the next PR. Right now for the purpose of testing
+    # we want to merge with 0 sedimentation and precipitation
+    #
+    if moisture_model isa NonEquilMoistModel
+        set_sedimentation_precomputed_quantities!(Y, p, t)
+        #    (; ᶜwₗ, ᶜwᵢ) = p.precomputed
+        #    @. ᶜwₜqₜ += Geometry.WVector(ᶜwₗ * Y.c.ρq_liq + ᶜwᵢ * Y.c.ρq_ice) / Y.c.ρ
+        #    @. ᶜwₕhₜ += Geometry.WVector(
+        #        ᶜwₗ * Y.c.ρq_liq * (TD.internal_energy_liquid(thermo_params, ᶜts) + ᶜΦ + norm_sqr(Geometry.UVWVector(0, 0, -ᶜwₗ) + Geometry.UVWVector(ᶜu))/2) +
+        #        ᶜwᵢ * Y.c.ρq_ice * (TD.internal_energy_ice(thermo_params, ᶜts)    + ᶜΦ + norm_sqr(Geometry.UVWVector(0, 0, -ᶜwᵢ) + Geometry.UVWVector(ᶜu))/2)
+        #    ) / Y.c.ρ
+    end
     if precip_model isa Microphysics1Moment
         set_precipitation_precomputed_quantities!(Y, p, t)
+        #    (; ᶜwᵣ, ᶜwₛ) = p.precomputed
+        #    @. ᶜwₜqₜ += Geometry.WVector(ᶜwᵣ * Y.c.ρq_rai + ᶜwₛ * Y.c.ρq_sno) / Y.c.ρ
+        #    @. ᶜwₕhₜ += Geometry.WVector(
+        #        ᶜwᵣ * Y.c.ρq_rai * (TD.internal_energy_liquid(thermo_params, ᶜts) + ᶜΦ + norm_sqr(Geometry.UVWVector(0, 0, -ᶜwᵣ) + Geometry.UVWVector(ᶜu))/2) +
+        #        ᶜwₛ * Y.c.ρq_sno * (TD.internal_energy_ice(thermo_params, ᶜts)    + ᶜΦ + norm_sqr(Geometry.UVWVector(0, 0, -ᶜwₛ) + Geometry.UVWVector(ᶜu))/2)
+        #    ) / Y.c.ρ
     end
 
     if turbconv_model isa PrognosticEDMFX

src/parameters/create_parameters.jl

@@ -7,15 +7,6 @@ import Thermodynamics.Parameters.ThermodynamicsParameters
 import CloudMicrophysics as CM
 import StaticArrays as SA
 
-function ClimaAtmosParameters(config::AtmosConfig)
-    (; toml_dict, parsed_args) = config
-    FT = CP.float_type(toml_dict)
-    return ClimaAtmosParameters(
-        toml_dict,
-        FT(CA.time_to_seconds(parsed_args["dt"])),
-    )
-end
-
 """
     ClimaAtmosParameters(FT, dt)
     ClimaAtmosParameters(toml_dict, dt)
@@ -25,6 +16,16 @@ Construct the parameter set for any ClimaAtmos configuration.
 
 If dt is passed in, it will be used to override the `precipitation_timescale` parameter.
 """
+function ClimaAtmosParameters(config::AtmosConfig)
+    (; toml_dict, parsed_args) = config
+    FT = CP.float_type(toml_dict)
+    override_dt =
+        parsed_args["override_precip_timescale"] ?
+        FT(CA.time_to_seconds(parsed_args["dt"])) : nothing
+
+    return ClimaAtmosParameters(toml_dict, override_dt)
+end
+
 ClimaAtmosParameters(::Type{FT}, dt = nothing) where {FT} =
     ClimaAtmosParameters(CP.create_toml_dict(FT), dt)