add microphysics 2M tests

Author: sajjadazimi

.buildkite/pipeline.yml

@@ -86,6 +86,13 @@ steps:
           --config_file $CONFIG_PATH/single_column_precipitation_test.yml
           --job_id single_column_precipitation_test
         artifact_paths: "single_column_precipitation_test/output_active/*"
+      
+      - label: ":umbrella: 2-moment precipitation sanity test single column"
+        command: >
+          julia --color=yes --project=.buildkite .buildkite/ci_driver.jl
+          --config_file $CONFIG_PATH/single_column_precipitation_2M_test.yml
+          --job_id single_column_precipitation_2M_test
+        artifact_paths: "single_column_precipitation_2M_test/output_active/*"
 
   - group: "Gravity wave"
     steps:

config/model_configs/single_column_precipitation_2M_test.yml

@@ -0,0 +1,22 @@
+config: "column"
+initial_condition: "PrecipitatingColumn"
+surface_setup: "DefaultMoninObukhov"
+z_elem: 200
+z_max: 10000.0
+z_stretch: false
+dt: "10secs"
+t_end: "1500secs"
+dt_save_state_to_disk: "500secs"
+dt_cloud_fraction: "60secs"
+cloud_model: "grid_scale"
+moist: "nonequil"
+precip_model: "2M"
+vert_diff: "DecayWithHeightDiffusion"
+implicit_diffusion: true
+approximate_linear_solve_iters: 2
+reproducibility_test: false
+diagnostics:
+  - short_name: [hus, clw, cli, husra, hussn, clwnm, clinm, nmra, nmsn, ta, wa]
+    period: 500secs
+  - short_name: [pr]
+    period: 10secs

post_processing/ci_plots.jl

@@ -621,23 +621,35 @@ function make_plots(
 end
 
 function make_plots(
-    ::Val{:single_column_precipitation_test},
+    sim_type::Union{
+        Val{:single_column_precipitation_test},
+        Val{:single_column_precipitation_2M_test},
+    },
     output_paths::Vector{<:AbstractString},
 )
     simdirs = SimDir.(output_paths)
 
     # TODO: Move this plotting code into the same framework as the other ones
     simdir = simdirs[1]
 
-    short_names = ["hus", "clw", "cli", "husra", "hussn", "ta"]
+    if sim_type isa Val{:single_column_precipitation_test}
+        short_names = ["hus", "clw", "cli", "husra", "hussn", "ta"]
+        figsize = (1200, 600)
+        pr_row = 3
+    else
+        short_names = ["hus", "clw", "cli", "husra", "hussn", "ta", "clwnm", "clinm", "nmra", "nmsn"]
+        figsize = (1200, 1000)
+        pr_row = 5
+    end
+
     vars = [
         slice(get(simdir; short_name), x = 0.0, y = 0.0) for
         short_name in short_names
     ]
 
     # We first prepare the axes with all the nice labels with ClimaAnalysis, then we use
     # CairoMakie to add the additional lines.
-    fig = CairoMakie.Figure(; size = (1200, 600))
+    fig = CairoMakie.Figure(; size = figsize)
 
     p_loc = [1, 1]
 
@@ -677,7 +689,7 @@ function make_plots(
     viz.line_plot1D!(
         fig,
         slice(surface_precip, x = 0.0, y = 0.0);
-        p_loc = [3, 1:3],
+        p_loc = [pr_row, 1:3],
     )
 
     file_path = joinpath(output_paths[1], "summary.pdf")

src/diagnostics/core_diagnostics.jl

@@ -889,7 +889,7 @@ add_diagnostic_variable!(
     short_name = "clwnm",
     long_name = "Number Concentration of Cloud Liquid Water",
     standard_name = "number_concentration_of_cloud_liquid_water",
-    units = "1 m^3",
+    units = "1 m^-3",
     comments = """
     This is calculated as the number of cloud liquid water droplets in the grid 
     cell divided by the cell volume.
@@ -920,7 +920,7 @@ add_diagnostic_variable!(
     short_name = "clinm",
     long_name = "Number Concentration of Cloud Ice",
     standard_name = "number_concentration_of_cloud_ice",
-    units = "1 m^3",
+    units = "1 m^-3",
     comments = """
     This is calculated as the number of cloud ice chrystals in the grid 
     cell divided by the cell volume.
@@ -951,7 +951,7 @@ add_diagnostic_variable!(
     short_name = "nmra",
     long_name = "Number Concentration of Rain",
     standard_name = "number_concentration_of_rain",
-    units = "1 m^3",
+    units = "1 m^-3",
     comments = """
     This is calculated as the number of raindrops in the grid cell divided
     by the cell volume.
@@ -982,7 +982,7 @@ add_diagnostic_variable!(
     short_name = "nmsn",
     long_name = "Number Concentration of Snow",
     standard_name = "number_concentration_of_snow",
-    units = "1 m^3",
+    units = "1 m^-3",
     comments = """
     This is calculated as the number of snow flakes in the grid cell divided
     by the cell volume.

test/parameterized_tendencies/microphysics/precipitation.jl

@@ -168,3 +168,113 @@ end
     @test minimum(p.precomputed.cloud_diagnostics_tuple.cf) >= FT(0)
     @test maximum(p.precomputed.cloud_diagnostics_tuple.cf) <= FT(1)
 end
+
+@testset "NonEquilibrium Moisture + 2-moment precipitation RHS terms" begin
+
+    ### Boilerplate default integrator objects
+    config = CA.AtmosConfig(
+        Dict(
+            "initial_condition" => "PrecipitatingColumn",
+            "moist" => "nonequil",
+            "precip_model" => "2M",
+            "config" => "column",
+            "output_default_diagnostics" => false,
+        ),
+        job_id = "precipitation_2M",
+    )
+    (; Y, p, params) = generate_test_simulation(config)
+
+    FT = eltype(Y)
+    ᶜYₜ = zero(Y)
+
+    # Set all model choices
+    (; turbconv_model, moisture_model, precip_model) = p.atmos
+
+    # Test cache to verify expected variables exist in tendency function
+    CA.set_precipitation_velocities!(Y, p, moisture_model, precip_model)
+    CA.set_precipitation_cache!(Y, p, precip_model, turbconv_model)
+    CA.set_precipitation_surface_fluxes!(Y, p, precip_model)
+    test_varnames = (
+        :ᶜSqₗᵖ,
+        :ᶜSqᵢᵖ,
+        :ᶜSqᵣᵖ,
+        :ᶜSqₛᵖ,
+        :ᶜSnₗᵖ,
+        :ᶜSnᵢᵖ,
+        :ᶜSnᵣᵖ,
+        :ᶜSnₛᵖ,
+        :surface_rain_flux,
+        :surface_snow_flux,
+        :ᶜwₗ,
+        :ᶜwᵢ,
+        :ᶜwᵣ,
+        :ᶜwₛ,
+        :ᶜwnₗ,
+        :ᶜwnᵢ,
+        :ᶜwnᵣ,
+        :ᶜwnₛ,
+        :ᶜwₜqₜ,
+        :ᶜwₕhₜ,
+    )
+    for var_name in test_varnames
+        @test var_name ∈ propertynames(p.precomputed)
+    end
+
+    # compute source terms based on the last model state
+    CA.precipitation_tendency!(
+        ᶜYₜ,
+        Y,
+        p,
+        FT(0),
+        moisture_model,
+        precip_model,
+        turbconv_model,
+    )
+
+    # check for nans
+    @assert !any(isnan, ᶜYₜ.c.ρ)
+    @assert !any(isnan, ᶜYₜ.c.ρq_tot)
+    @assert !any(isnan, ᶜYₜ.c.ρe_tot)
+    @assert !any(isnan, ᶜYₜ.c.ρq_liq)
+    @assert !any(isnan, ᶜYₜ.c.ρq_ice)
+    @assert !any(isnan, ᶜYₜ.c.ρq_rai)
+    @assert !any(isnan, ᶜYₜ.c.ρq_sno)
+    @assert !any(isnan, ᶜYₜ.c.ρn_liq)
+    @assert !any(isnan, ᶜYₜ.c.ρn_ice)
+    @assert !any(isnan, ᶜYₜ.c.ρn_rai)
+    @assert !any(isnan, ᶜYₜ.c.ρn_sno)
+    @assert !any(isnan, p.precomputed.ᶜwₗ)
+    @assert !any(isnan, p.precomputed.ᶜwᵢ)
+    @assert !any(isnan, p.precomputed.ᶜwᵣ)
+    @assert !any(isnan, p.precomputed.ᶜwₛ)
+    @assert !any(isnan, p.precomputed.ᶜwnₗ)
+    @assert !any(isnan, p.precomputed.ᶜwnᵢ)
+    @assert !any(isnan, p.precomputed.ᶜwnᵣ)
+    @assert !any(isnan, p.precomputed.ᶜwnₛ)
+
+    # test water budget
+    @test ᶜYₜ.c.ρ == ᶜYₜ.c.ρq_tot
+    @assert iszero(ᶜYₜ.c.ρ)
+
+    # test nonequilibrium cloud condensate
+    CA.cloud_condensate_tendency!(ᶜYₜ, Y, p, moisture_model, precip_model)
+    @assert !any(isnan, ᶜYₜ.c.ρq_liq)
+    @assert !any(isnan, ᶜYₜ.c.ρq_ice)
+    @assert !any(isnan, ᶜYₜ.c.ρn_liq)
+    @assert !any(isnan, ᶜYₜ.c.ρn_ice)
+
+    # test if terminal velocity is positive
+    @test minimum(p.precomputed.ᶜwₗ) >= FT(0)
+    @test minimum(p.precomputed.ᶜwᵢ) >= FT(0)
+    @test minimum(p.precomputed.ᶜwᵣ) >= FT(0)
+    @test minimum(p.precomputed.ᶜwₛ) >= FT(0)
+    @test minimum(p.precomputed.ᶜwnₗ) >= FT(0)
+    @test minimum(p.precomputed.ᶜwnᵢ) >= FT(0)
+    @test minimum(p.precomputed.ᶜwnᵣ) >= FT(0)
+    @test minimum(p.precomputed.ᶜwnₛ) >= FT(0)
+
+    # test if cloud fraction diagnostics make sense
+    @assert !any(isnan, p.precomputed.cloud_diagnostics_tuple.cf)
+    @test minimum(p.precomputed.cloud_diagnostics_tuple.cf) >= FT(0)
+    @test maximum(p.precomputed.cloud_diagnostics_tuple.cf) <= FT(1)
+end