Skip to content

Commit

Permalink
Merge #2091
Browse files Browse the repository at this point in the history
2091: Use weak divergence in scalar advection r=simonbyrne a=simonbyrne



Co-authored-by: Simon Byrne <simonbyrne@gmail.com>
  • Loading branch information
bors[bot] and simonbyrne authored Sep 18, 2023
2 parents cb44933 + 3b48dd5 commit 5eb861b
Show file tree
Hide file tree
Showing 3 changed files with 13 additions and 14 deletions.
7 changes: 3 additions & 4 deletions docs/src/equations.md
Original file line number Diff line number Diff line change
Expand Up @@ -33,7 +33,6 @@ We make use of the following operators

### Differential operators

- ``\mathcal{D}_h`` is the [discrete horizontal spectral divergence](https://clima.github.io/ClimaCore.jl/stable/operators/#ClimaCore.Operators.Divergence).
- ``\hat{\mathcal{D}}_h`` is the [discrete horizontal spectral weak divergence](https://clima.github.io/ClimaCore.jl/stable/operators/#ClimaCore.Operators.WeakDivergence).
- ``\mathcal{D}^c_v`` is the [face-to-center vertical divergence](https://clima.github.io/ClimaCore.jl/stable/operators/#ClimaCore.Operators.DivergenceF2C).
!!! todo
Expand Down Expand Up @@ -107,7 +106,7 @@ Follows the continuity equation
This is discretized using the following
```math
\frac{\partial}{\partial t} \rho
= - \mathcal{D}_h[ \rho \bar{\boldsymbol{u}}] - \mathcal{D}^c_v \left[WI^f( J, \rho) \tilde{\boldsymbol{u}} \right]
= - \hat{\mathcal{D}}_h[ \rho \bar{\boldsymbol{u}}] - \mathcal{D}^c_v \left[WI^f( J, \rho) \tilde{\boldsymbol{u}} \right]
```

with the
Expand Down Expand Up @@ -210,7 +209,7 @@ which is stabilized with the addition of a 4th-order hyperdiffusion term on tota
This is discretized using
```math
\frac{\partial}{\partial t} \rho e \approx
- \mathcal{D}_h[ (\rho e + p) \bar{\boldsymbol{u}} ]
- \hat{\mathcal{D}}_h[ (\rho e + p) \bar{\boldsymbol{u}} ]
- \mathcal{D}^c_v \left[ WI^f(J,\rho) \, \tilde{\boldsymbol{u}} \, I^f \left(\frac{\rho e + p}{\rho} \right)
+ \boldsymbol{F}_R \right] - \nu_h \hat{\mathcal{D}}_h( \rho \mathcal{G}_h(\psi) ).
```
Expand Down Expand Up @@ -247,7 +246,7 @@ This is stabilized with the addition of a 4th-order hyperdiffusion term
This is discretized using
```math
\frac{\partial}{\partial t} \rho \chi \approx
- \mathcal{D}_h[ \rho \chi \bar{\boldsymbol{u}}]
- \hat{\mathcal{D}}_h[ \rho \chi \bar{\boldsymbol{u}}]
- \mathcal{D}^c_v \left[ WI^f(J,\rho) \, U^f\left( \tilde{\boldsymbol{u}}, \frac{\rho \chi}{\rho} \right) \right]
- \nu_\chi \hat{\mathcal{D}}_h ( \rho \, \mathcal{G}_h (\psi) )
```
Expand Down
2 changes: 1 addition & 1 deletion regression_tests/ref_counter.jl
Original file line number Diff line number Diff line change
@@ -1 +1 @@
125
126
18 changes: 9 additions & 9 deletions src/prognostic_equations/advection.jl
Original file line number Diff line number Diff line change
Expand Up @@ -16,34 +16,34 @@ NVTX.@annotate function horizontal_advection_tendency!(Yₜ, Y, p, t)
(; ᶜuʲs) = p
end

@. Yₜ.c.ρ -= divₕ(Y.c.ρ * ᶜu)
@. Yₜ.c.ρ -= wdivₕ(Y.c.ρ * ᶜu)
if p.atmos.turbconv_model isa EDMFX
for j in 1:n
@. Yₜ.c.sgsʲs.:($$j).ρa -= divₕ(Y.c.sgsʲs.:($$j).ρa * ᶜuʲs.:($$j))
@. Yₜ.c.sgsʲs.:($$j).ρa -= wdivₕ(Y.c.sgsʲs.:($$j).ρa * ᶜuʲs.:($$j))
end
end

if :ρθ in propertynames(Y.c)
@. Yₜ.c.ρθ -= divₕ(Y.c.ρθ * ᶜu)
@. Yₜ.c.ρθ -= wdivₕ(Y.c.ρθ * ᶜu)
elseif :ρe_tot in propertynames(Y.c)
(; ᶜh_tot) = p
@. Yₜ.c.ρe_tot -= divₕ(Y.c.ρ * ᶜh_tot * ᶜu)
@. Yₜ.c.ρe_tot -= wdivₕ(Y.c.ρ * ᶜh_tot * ᶜu)
end
if p.atmos.turbconv_model isa EDMFX
for j in 1:n
if :ρθ in propertynames(Y.c)
@. Yₜ.c.sgsʲs.:($$j).ρaθ -=
divₕ(Y.c.sgsʲs.:($$j).ρaθ * ᶜuʲs.:($$j))
wdivₕ(Y.c.sgsʲs.:($$j).ρaθ * ᶜuʲs.:($$j))
elseif :ρe_tot in propertynames(Y.c)
(; ᶜh_totʲs) = p
@. Yₜ.c.sgsʲs.:($$j).ρae_tot -=
divₕ(Y.c.sgsʲs.:($$j).ρa * ᶜh_totʲs.:($$j) * ᶜuʲs.:($$j))
wdivₕ(Y.c.sgsʲs.:($$j).ρa * ᶜh_totʲs.:($$j) * ᶜuʲs.:($$j))
end
end
end

if use_prognostic_tke(p.atmos.turbconv_model)
@. Yₜ.c.sgs⁰.ρatke -= divₕ(Y.c.sgs⁰.ρatke * ᶜu⁰)
@. Yₜ.c.sgs⁰.ρatke -= wdivₕ(Y.c.sgs⁰.ρatke * ᶜu⁰)
end

@. Yₜ.c.uₕ -= C12(gradₕ(ᶜp - ᶜp_ref) / Y.c.ρ + gradₕ(ᶜK + ᶜΦ))
Expand All @@ -59,15 +59,15 @@ NVTX.@annotate function horizontal_tracer_advection_tendency!(Yₜ, Y, p, t)
end

for ρχ_name in filter(is_tracer_var, propertynames(Y.c))
@. Yₜ.c.:($$ρχ_name) -= divₕ(Y.c.:($$ρχ_name) * ᶜu)
@. Yₜ.c.:($$ρχ_name) -= wdivₕ(Y.c.:($$ρχ_name) * ᶜu)
end

if p.atmos.turbconv_model isa EDMFX
for j in 1:n
for ρaχ_name in
filter(is_tracer_var, propertynames(Y.c.sgsʲs.:($j)))
@. Yₜ.c.sgsʲs.:($$j).:($$ρaχ_name) -=
divₕ(Y.c.sgsʲs.:($$j).:($$ρaχ_name) * ᶜuʲs.:($$j))
wdivₕ(Y.c.sgsʲs.:($$j).:($$ρaχ_name) * ᶜuʲs.:($$j))
end
end
end
Expand Down

0 comments on commit 5eb861b

Please sign in to comment.