Skip to content

Commit

Permalink
Add four more parameters for shear turbulence in unstable stratificat…
Browse files Browse the repository at this point in the history 
…ion (#3600)

* Add special negative Ri branch to the stability function

* Update parameters

* Improve show for CATKE
  • Loading branch information
@glwagner
glwagner committed May 16, 2024
1 parent c77e1a6 commit e8537d6
Show file tree
Hide file tree
Showing 3 changed files with 110 additions and 88 deletions.
81 changes: 43 additions & 38 deletions ...rbulence_closure_implementations/CATKEVerticalDiffusivities/CATKEVerticalDiffusivities.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -122,7 +122,7 @@ include("time_step_turbulent_kinetic_energy.jl")
minimum_turbulent_kinetic_energy = 1e-9,
minimum_convective_buoyancy_flux = 1e-11,
negative_turbulent_kinetic_energy_damping_time_scale = 1minute,
turbulent_kinetic_energy_time_step = 10)
turbulent_kinetic_energy_time_step = nothing)
Return the `CATKEVerticalDiffusivity` turbulence closure for vertical mixing by
small-scale ocean turbulence based on the prognostic evolution of subgrid
Expand Down Expand Up @@ -176,7 +176,7 @@ function CATKEVerticalDiffusivity(time_discretization::TD = VerticallyImplicitTi
minimum_turbulent_kinetic_energy = 1e-9,
minimum_convective_buoyancy_flux = 1e-11,
negative_turbulent_kinetic_energy_damping_time_scale = 1minute,
turbulent_kinetic_energy_time_step = 10) where TD
turbulent_kinetic_energy_time_step = nothing) where TD

mixing_length = convert_eltype(FT, mixing_length)
turbulent_kinetic_energy_equation = convert_eltype(FT, turbulent_kinetic_energy_equation)
Expand Down Expand Up @@ -226,8 +226,8 @@ catke_first(catke1::FlavorOfCATKE, catke2::FlavorOfCATKE) = error("Can't have tw
N² = ∂z_b(i, j, k, grid, buoyancy, tracers)
S² = ∂z_u² + ∂z_v²
Ri = N² / S²
return ifelse(N² ≤ 0, zero(grid), Ri)
#return ifelse(N² == 0, zero(grid), Ri)
#return ifelse(N² ≤ 0, zero(grid), Ri)
return ifelse(N² == 0, zero(grid), Ri)
end

for S in (:MixingLength, :TurbulentKineticEnergyEquation)
Expand Down Expand Up @@ -406,41 +406,46 @@ function Base.summary(closure::CATKEVD)
return string("CATKEVerticalDiffusivity{$TD}")
end

function Base.show(io::IO, closure::FlavorOfCATKE)
print(io, summary(closure))
function Base.show(io::IO, clo::FlavorOfCATKE)
print(io, summary(clo))
print(io, '\n')
print(io, "├── maximum_tracer_diffusivity: ", prettysummary(closure.maximum_tracer_diffusivity), '\n',
"├── maximum_tke_diffusivity: ", prettysummary(closure.maximum_tke_diffusivity), '\n',
"├── maximum_viscosity: ", prettysummary(closure.maximum_viscosity), '\n',
"├── minimum_turbulent_kinetic_energy: ", prettysummary(closure.minimum_turbulent_kinetic_energy), '\n',
"├── negative_turbulent_kinetic_energy_damping_time_scale: ", prettysummary(closure.negative_turbulent_kinetic_energy_damping_time_scale), '\n',
"├── minimum_convective_buoyancy_flux: ", prettysummary(closure.minimum_convective_buoyancy_flux), '\n',
"├── turbulent_kinetic_energy_time_step: ", prettysummary(closure.turbulent_kinetic_energy_time_step), '\n',
"├── mixing_length: ", prettysummary(closure.mixing_length), '\n',
"│ ├── Cˢ: ", prettysummary(closure.mixing_length.Cˢ), '\n',
"│ ├── Cᵇ: ", prettysummary(closure.mixing_length.Cᵇ), '\n',
"│ ├── Cᶜu: ", prettysummary(closure.mixing_length.Cᶜu), '\n',
"│ ├── Cᶜc: ", prettysummary(closure.mixing_length.Cᶜc), '\n',
"│ ├── Cᶜe: ", prettysummary(closure.mixing_length.Cᶜe), '\n',
"│ ├── Cᵉc: ", prettysummary(closure.mixing_length.Cᵉc), '\n',
"│ ├── Cᵉe: ", prettysummary(closure.mixing_length.Cᵉe), '\n',
"│ ├── Cˡᵒu: ", prettysummary(closure.mixing_length.Cˡᵒu), '\n',
"│ ├── Cˡᵒc: ", prettysummary(closure.mixing_length.Cˡᵒc), '\n',
"│ ├── Cˡᵒe: ", prettysummary(closure.mixing_length.Cˡᵒe), '\n',
"│ ├── Cʰⁱu: ", prettysummary(closure.mixing_length.Cʰⁱu), '\n',
"│ ├── Cʰⁱc: ", prettysummary(closure.mixing_length.Cʰⁱc), '\n',
"│ ├── Cʰⁱe: ", prettysummary(closure.mixing_length.Cʰⁱe), '\n',
"│ ├── Cˢᵖ: ", prettysummary(closure.mixing_length.Cˢᵖ), '\n',
"│ ├── CRiᵟ: ", prettysummary(closure.mixing_length.CRiᵟ), '\n',
"│ └── CRi⁰: ", prettysummary(closure.mixing_length.CRi⁰), '\n',
"└── turbulent_kinetic_energy_equation: ", prettysummary(closure.turbulent_kinetic_energy_equation), '\n',
" ├── CˡᵒD: ", prettysummary(closure.turbulent_kinetic_energy_equation.CˡᵒD), '\n',
" ├── CʰⁱD: ", prettysummary(closure.turbulent_kinetic_energy_equation.CʰⁱD), '\n',
" ├── CᶜD: ", prettysummary(closure.turbulent_kinetic_energy_equation.CᶜD), '\n',
" ├── CᵉD: ", prettysummary(closure.turbulent_kinetic_energy_equation.CᵉD), '\n',
" ├── Cᵂu★: ", prettysummary(closure.turbulent_kinetic_energy_equation.Cᵂu★), '\n',
" ├── CᵂwΔ: ", prettysummary(closure.turbulent_kinetic_energy_equation.CᵂwΔ), '\n',
" └── Cᵂϵ: ", prettysummary(closure.turbulent_kinetic_energy_equation.Cᵂϵ))
print(io, "├── maximum_tracer_diffusivity: ", prettysummary(clo.maximum_tracer_diffusivity), '\n',
"├── maximum_tke_diffusivity: ", prettysummary(clo.maximum_tke_diffusivity), '\n',
"├── maximum_viscosity: ", prettysummary(clo.maximum_viscosity), '\n',
"├── minimum_turbulent_kinetic_energy: ", prettysummary(clo.minimum_turbulent_kinetic_energy), '\n',
"├── negative_turbulent_kinetic_energy_damping_time_scale: ", prettysummary(clo.negative_turbulent_kinetic_energy_damping_time_scale), '\n',
"├── minimum_convective_buoyancy_flux: ", prettysummary(clo.minimum_convective_buoyancy_flux), '\n',
"├── turbulent_kinetic_energy_time_step: ", prettysummary(clo.turbulent_kinetic_energy_time_step), '\n',
"├── mixing_length: ", prettysummary(clo.mixing_length), '\n',
"│ ├── Cˢ: ", prettysummary(clo.mixing_length.Cˢ), '\n',
"│ ├── Cᵇ: ", prettysummary(clo.mixing_length.Cᵇ), '\n',
"│ ├── Cʰⁱu: ", prettysummary(clo.mixing_length.Cʰⁱu), '\n',
"│ ├── Cʰⁱc: ", prettysummary(clo.mixing_length.Cʰⁱc), '\n',
"│ ├── Cʰⁱe: ", prettysummary(clo.mixing_length.Cʰⁱe), '\n',
"│ ├── Cˡᵒu: ", prettysummary(clo.mixing_length.Cˡᵒu), '\n',
"│ ├── Cˡᵒc: ", prettysummary(clo.mixing_length.Cˡᵒc), '\n',
"│ ├── Cˡᵒe: ", prettysummary(clo.mixing_length.Cˡᵒe), '\n',
"│ ├── Cᵘⁿu: ", prettysummary(clo.mixing_length.Cᵘⁿu), '\n',
"│ ├── Cᵘⁿc: ", prettysummary(clo.mixing_length.Cᵘⁿc), '\n',
"│ ├── Cᵘⁿe: ", prettysummary(clo.mixing_length.Cᵘⁿe), '\n',
"│ ├── Cᶜu: ", prettysummary(clo.mixing_length.Cᶜu), '\n',
"│ ├── Cᶜc: ", prettysummary(clo.mixing_length.Cᶜc), '\n',
"│ ├── Cᶜe: ", prettysummary(clo.mixing_length.Cᶜe), '\n',
"│ ├── Cᵉc: ", prettysummary(clo.mixing_length.Cᵉc), '\n',
"│ ├── Cᵉe: ", prettysummary(clo.mixing_length.Cᵉe), '\n',
"│ ├── Cˢᵖ: ", prettysummary(clo.mixing_length.Cˢᵖ), '\n',
"│ ├── CRiᵟ: ", prettysummary(clo.mixing_length.CRiᵟ), '\n',
"│ └── CRi⁰: ", prettysummary(clo.mixing_length.CRi⁰), '\n',
"└── turbulent_kinetic_energy_equation: ", prettysummary(clo.turbulent_kinetic_energy_equation), '\n',
" ├── CʰⁱD: ", prettysummary(clo.turbulent_kinetic_energy_equation.CʰⁱD), '\n',
" ├── CˡᵒD: ", prettysummary(clo.turbulent_kinetic_energy_equation.CˡᵒD), '\n',
" ├── CᵘⁿD: ", prettysummary(clo.turbulent_kinetic_energy_equation.CᵘⁿD), '\n',
" ├── CᶜD: ", prettysummary(clo.turbulent_kinetic_energy_equation.CᶜD), '\n',
" ├── CᵉD: ", prettysummary(clo.turbulent_kinetic_energy_equation.CᵉD), '\n',
" ├── Cᵂu★: ", prettysummary(clo.turbulent_kinetic_energy_equation.Cᵂu★), '\n',
" ├── CᵂwΔ: ", prettysummary(clo.turbulent_kinetic_energy_equation.CᵂwΔ), '\n',
" └── Cᵂϵ: ", prettysummary(clo.turbulent_kinetic_energy_equation.Cᵂϵ))
end

end # module

88 changes: 51 additions & 37 deletions ...ceClosures/turbulence_closure_implementations/CATKEVerticalDiffusivities/mixing_length.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -13,23 +13,26 @@ using ..TurbulenceClosures:
Contains mixing length parameters for CATKE vertical diffusivity.
"""
Base.@kwdef struct MixingLength{FT}
Cˢ :: FT = 1.03 # Surface distance coefficient for shear length scale
Cˢ :: FT = 1.45 # Surface distance coefficient for shear length scale
Cᵇ :: FT = Inf # Bottom distance coefficient for shear length scale
Cᶜu :: FT = 0.0 # Convective mixing length coefficient for tracers
Cᶜc :: FT = 3.72 # Convective mixing length coefficient for tracers
Cᶜe :: FT = 1.40 # Convective mixing length coefficient for TKE
Cˢᵖ :: FT = 0.197 # Sheared convective plume coefficient
CRiᵟ :: FT = 0.301 # Stability function width
CRi⁰ :: FT = 0.614 # Stability function lower Ri
Cʰⁱu :: FT = 0.174 # Shear mixing length coefficient for momentum at high Ri
Cˡᵒu :: FT = 0.270 # Shear mixing length coefficient for momentum at low Ri
Cᵘⁿu :: FT = 0.267 # Shear mixing length coefficient for momentum at negative Ri
Cᶜu :: FT = 0.814 # Convective mixing length coefficient for tracers
Cᵉu :: FT = 0.0 # Convective penetration mixing length coefficient for tracers
Cᵉc :: FT = 0.841 # Convective penetration mixing length coefficient for tracers
Cʰⁱc :: FT = 0.072 # Shear mixing length coefficient for tracers at high Ri
Cˡᵒc :: FT = 0.246 # Shear mixing length coefficient for tracers at low Ri
Cᵘⁿc :: FT = 1.711 # Shear mixing length coefficient for tracers at negative Ri
Cᶜc :: FT = 7.328 # Convective mixing length coefficient for tracers
Cᵉc :: FT = 0.054 # Convective penetration mixing length coefficient for tracers
Cʰⁱe :: FT = 1.186 # Shear mixing length coefficient for TKE at high Ri
Cˡᵒe :: FT = 2.886 # Shear mixing length coefficient for TKE at low Ri
Cᵘⁿe :: FT = 5.498 # Shear mixing length coefficient for TKE at negative Ri
Cᶜe :: FT = 0.891 # Convective mixing length coefficient for TKE
Cᵉe :: FT = 0.0 # Convective penetration mixing length coefficient for TKE
Cˢᵖ :: FT = 1.02 # Sheared convective plume coefficient
Cˡᵒu :: FT = 0.775 # Shear mixing length coefficient for momentum at low Ri
Cʰⁱu :: FT = 0.765 # Shear mixing length coefficient for momentum at high Ri
Cˡᵒc :: FT = 0.786 # Shear mixing length coefficient for tracers at low Ri
Cʰⁱc :: FT = 0.378 # Shear mixing length coefficient for tracers at high Ri
Cˡᵒe :: FT = 2.68 # Shear mixing length coefficient for TKE at low Ri
Cʰⁱe :: FT = 2.16 # Shear mixing length coefficient for TKE at high Ri
CRiᵟ :: FT = 0.421 # Stability function width
CRi⁰ :: FT = 0.261 # Stability function lower Ri
end

#####
Expand Down Expand Up @@ -208,7 +211,11 @@ end
"""Piecewise linear function between 0 (when x < c) and 1 (when x - c > w)."""
@inline step(x, c, w) = max(zero(x), min(one(x), (x - c) / w))

@inline scale(Ri, σ⁻, σ⁺ , c, w) = σ⁻ + (σ⁺ - σ⁻) * step(Ri, c, w)
@inline function scale(Ri, σ⁻, σ⁰, σ∞, c, w)
σ⁺ = σ⁰ + (σ∞ - σ⁰) * step(Ri, c, w)
σ = σ⁻ * (Ri < 0) + σ⁺ * (Ri ≥ 0)
return σ
end

#=
@inline function scale(Ri, σ₀, σ∞ , Δσ, δRi)
Expand All @@ -219,18 +226,18 @@ end
end
=#

@inline function stability_functionᶜᶜᶠ(i, j, k, grid, closure, Cˡᵒ, Cʰⁱ, velocities, tracers, buoyancy)
@inline function stability_functionᶜᶜᶠ(i, j, k, grid, closure, Cᵘⁿ, Cˡᵒ, Cʰⁱ, velocities, tracers, buoyancy)
Ri = Riᶜᶜᶠ(i, j, k, grid, velocities, tracers, buoyancy)
CRi⁰ = closure.mixing_length.CRi⁰
CRiᵟ = closure.mixing_length.CRiᵟ
return scale(Ri, Cˡᵒ, Cʰⁱ, CRi⁰, CRiᵟ)
return scale(Ri, Cᵘⁿ, Cˡᵒ, Cʰⁱ, CRi⁰, CRiᵟ)
end

@inline function stability_functionᶜᶜᶜ(i, j, k, grid, closure, Cˡᵒ, Cʰⁱ, velocities, tracers, buoyancy)
@inline function stability_functionᶜᶜᶜ(i, j, k, grid, closure, Cᵘⁿ, Cˡᵒ, Cʰⁱ, velocities, tracers, buoyancy)
Ri = Riᶜᶜᶜ(i, j, k, grid, velocities, tracers, buoyancy)
CRi⁰ = closure.mixing_length.CRi⁰
CRiᵟ = closure.mixing_length.CRiᵟ
return scale(Ri, Cˡᵒ, Cʰⁱ, CRi⁰, CRiᵟ)
return scale(Ri, Cᵘⁿ, Cˡᵒ, Cʰⁱ, CRi⁰, CRiᵟ)
end

@inline function momentum_mixing_lengthᶜᶜᶠ(i, j, k, grid, closure, velocities, tracers, buoyancy, surface_buoyancy_flux)
Expand All @@ -239,9 +246,10 @@ end
Cˢᵖ = closure.mixing_length.Cˢᵖ
ℓʰ = convective_length_scaleᶜᶜᶠ(i, j, k, grid, closure, Cᶜ, Cᵉ, Cˢᵖ, velocities, tracers, buoyancy, surface_buoyancy_flux)

Cᵘⁿ = closure.mixing_length.Cᵘⁿu
Cˡᵒ = closure.mixing_length.Cˡᵒu
Cʰⁱ = closure.mixing_length.Cʰⁱu
σ = stability_functionᶜᶜᶠ(i, j, k, grid, closure, Cˡᵒ, Cʰⁱ, velocities, tracers, buoyancy)
σ = stability_functionᶜᶜᶠ(i, j, k, grid, closure, Cᵘⁿ, Cˡᵒ, Cʰⁱ, velocities, tracers, buoyancy)

ℓ★ = σ * stable_length_scaleᶜᶜᶠ(i, j, k, grid, closure, tracers.e, velocities, tracers, buoyancy)

Expand All @@ -259,9 +267,10 @@ end
Cˢᵖ = closure.mixing_length.Cˢᵖ
ℓʰ = convective_length_scaleᶜᶜᶠ(i, j, k, grid, closure, Cᶜ, Cᵉ, Cˢᵖ, velocities, tracers, buoyancy, surface_buoyancy_flux)

Cᵘⁿ = closure.mixing_length.Cᵘⁿc
Cˡᵒ = closure.mixing_length.Cˡᵒc
Cʰⁱ = closure.mixing_length.Cʰⁱc
σ = stability_functionᶜᶜᶠ(i, j, k, grid, closure, Cˡᵒ, Cʰⁱ, velocities, tracers, buoyancy)
σ = stability_functionᶜᶜᶠ(i, j, k, grid, closure, Cᵘⁿ, Cˡᵒ, Cʰⁱ, velocities, tracers, buoyancy)
ℓ★ = σ * stable_length_scaleᶜᶜᶠ(i, j, k, grid, closure, tracers.e, velocities, tracers, buoyancy)

ℓʰ = ifelse(isnan(ℓʰ), zero(grid), ℓʰ)
Expand All @@ -278,9 +287,10 @@ end
Cˢᵖ = closure.mixing_length.Cˢᵖ
ℓʰ = convective_length_scaleᶜᶜᶠ(i, j, k, grid, closure, Cᶜ, Cᵉ, Cˢᵖ, velocities, tracers, buoyancy, surface_buoyancy_flux)

Cᵘⁿ = closure.mixing_length.Cᵘⁿe
Cˡᵒ = closure.mixing_length.Cˡᵒe
Cʰⁱ = closure.mixing_length.Cʰⁱe
σ = stability_functionᶜᶜᶠ(i, j, k, grid, closure, Cˡᵒ, Cʰⁱ, velocities, tracers, buoyancy)
σ = stability_functionᶜᶜᶠ(i, j, k, grid, closure, Cᵘⁿ, Cˡᵒ, Cʰⁱ, velocities, tracers, buoyancy)
ℓ★ = σ * stable_length_scaleᶜᶜᶠ(i, j, k, grid, closure, tracers.e, velocities, tracers, buoyancy)

ℓʰ = ifelse(isnan(ℓʰ), zero(grid), ℓʰ)
Expand All @@ -295,19 +305,23 @@ Base.summary(::MixingLength) = "CATKEVerticalDiffusivities.MixingLength"

Base.show(io::IO, ml::MixingLength) =
print(io, "CATKEVerticalDiffusivities.MixingLength parameters:", '\n',
" Cˢ: $(ml.Cˢ)", '\n',
" Cᵇ: $(ml.Cᵇ)", '\n',
" Cᶜc: $(ml.Cᶜc)", '\n',
" Cᶜe: $(ml.Cᶜe)", '\n',
" Cᵉc: $(ml.Cᵉc)", '\n',
" Cᵉe: $(ml.Cᵉe)", '\n',
" Cˡᵒu: $(ml.Cˡᵒu)", '\n',
" Cˡᵒc: $(ml.Cˡᵒc)", '\n',
" Cˡᵒe: $(ml.Cˡᵒe)", '\n',
" Cʰⁱu: $(ml.Cʰⁱu)", '\n',
" Cʰⁱc: $(ml.Cʰⁱc)", '\n',
" Cʰⁱe: $(ml.Cʰⁱe)", '\n',
" Cˢᵖ: $(ml.Cˢᵖ)", '\n',
" CRiᵟ: $(ml.CRiᵟ)", '\n',
" CRi⁰: $(ml.CRi⁰)")
" ├── Cˢ: ", ml.Cˢ, '\n',
" ├── Cᵇ: ", ml.Cᵇ, '\n',
" ├── Cʰⁱu: ", ml.Cʰⁱu, '\n',
" ├── Cʰⁱc: ", ml.Cʰⁱc, '\n',
" ├── Cʰⁱe: ", ml.Cʰⁱe, '\n',
" ├── Cˡᵒu: ", ml.Cˡᵒu, '\n',
" ├── Cˡᵒc: ", ml.Cˡᵒc, '\n',
" ├── Cˡᵒe: ", ml.Cˡᵒe, '\n',
" ├── Cᵘⁿu: ", ml.Cᵘⁿu, '\n',
" ├── Cᵘⁿc: ", ml.Cᵘⁿc, '\n',
" ├── Cᵘⁿe: ", ml.Cᵘⁿe, '\n',
" ├── Cᶜu: ", ml.Cᶜu, '\n',
" ├── Cᶜc: ", ml.Cᶜc, '\n',
" ├── Cᶜe: ", ml.Cᶜe, '\n',
" ├── Cᵉc: ", ml.Cᵉc, '\n',
" ├── Cᵉe: ", ml.Cᵉe, '\n',
" ├── Cˢᵖ: ", ml.Cˢᵖ, '\n',
" ├── CRiᵟ: ", ml.CRiᵟ, '\n',
" └── CRi⁰: ", ml.CRi⁰)

Loading

0 comments on commit e8537d6

Please sign in to comment.