Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add four more parameters for shear turbulence in unstable stratification #3600

Merged
merged 3 commits into from
May 16, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
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

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