Merge pull request #3988 from CliMA/zs/diagedmf_massflux

add mass flux to microphysic tracers in diagnostic EDMF

Author: imreddyTeja

Project.toml

@@ -64,7 +64,7 @@ NVTX = "0.3, 1"
 NullBroadcasts = "0.1"
 RRTMGP = "0.21.5"
 Random = "1"
-SciMLBase = "2.108"
+SciMLBase = "~2.108, ~2.109, ~2.110, ~2.111, ~2.112"
 SparseMatrixColorings = "0.4.20"
 StaticArrays = "1.9"
 Statistics = "1"

src/prognostic_equations/edmfx_sgs_flux.jl

@@ -269,7 +269,60 @@ function edmfx_sgs_mass_flux_tendency!(
                 )
                 @. Yₜ.c.ρq_tot += vtt
             end
-            # TODO: add environment flux?
+        end
+
+        # Microphysics tracers fluxes
+        if p.atmos.moisture_model isa NonEquilMoistModel && (
+            p.atmos.microphysics_model isa Microphysics1Moment ||
+            p.atmos.microphysics_model isa Microphysics2Moment
+        )
+            microphysics_tracers = (
+                (@name(c.ρq_liq), @name(ᶜq_liqʲs.:(1))),
+                (@name(c.ρq_ice), @name(ᶜq_iceʲs.:(1))),
+                (@name(c.ρq_rai), @name(ᶜq_raiʲs.:(1))),
+                (@name(c.ρq_sno), @name(ᶜq_snoʲs.:(1))),
+                (@name(c.ρn_liq), @name(ᶜn_liqʲs.:(1))),
+                (@name(c.ρn_rai), @name(ᶜn_raiʲs.:(1))),
+            )
+            for j in 1:n
+                # TODO using unrolled_foreach here allocates! (breaks the flame tests
+                # even though they use 0M microphysics)
+                # MatrixFields.unrolled_foreach(
+                #     microphysics_tracers,
+                # ) do (ρχ_name, χʲ_name)
+                for (ρχ_name, χʲ_name) in microphysics_tracers
+                    MatrixFields.has_field(Y, ρχ_name) || continue
+
+                    ᶜχʲ = MatrixFields.get_field(p.precomputed, χʲ_name)
+                    ᶜρχ = MatrixFields.get_field(Y, ρχ_name)
+                    ᶜχ = (@. lazy(specific(ᶜρχ, Y.c.ρ)))
+
+                    @. ᶠu³_diff = ᶠu³ʲs.:($$j) - ᶠu³
+                    # @. ᶜa_scalar =
+                    #     (ᶜχʲ - ᶜχ) *
+                    #     draft_area(ᶜρaʲs.:($$j), ᶜρʲs.:($$j))
+                    # TODO: remove this filter when mass flux is treated implicitly
+                    @. ᶜa_scalar =
+                        (ᶜχʲ - ᶜχ) * min(
+                            min(draft_area(ᶜρaʲs.:($$j), ᶜρʲs.:($$j)), a_max),
+                            FT(0.02) / max(
+                                Geometry.WVector(
+                                    ᶜinterp(ᶠu³_diff),
+                                ).components.data.:1,
+                                eps(FT),
+                            ),
+                        )
+                    vtt = vertical_transport_precip_massflux(
+                        ᶜρʲs.:($j),
+                        ᶠu³_diff,
+                        ᶜa_scalar,
+                        dt,
+                        edmfx_sgsflux_upwinding,
+                    )
+                    ᶜρχₜ = MatrixFields.get_field(Yₜ, ρχ_name)
+                    @. ᶜρχₜ += vtt
+                end
+            end
         end
         # TODO: the following adds the environment flux to the tendency
         # Make active and test later