From 643b484e81e0aeb038b3038266912ad051bce9b8 Mon Sep 17 00:00:00 2001 From: "Gregory L. Wagner" Date: Wed, 28 Feb 2024 07:09:50 -0700 Subject: [PATCH] (0.90.9) Fix a sign error in CATKE's implicit dissipative buoyancy flux (#3484) * Fix sign error in implicit buoyancy flux * Bump to 0.90.9 * Put in new calibrated parameters * Update windy convection script * Pointless commit --- Project.toml | 2 +- .../CATKEVerticalDiffusivities.jl | 2 +- .../mixing_length.jl | 26 ++++---- .../turbulent_kinetic_energy_equation.jl | 10 +-- .../column_windy_convection.jl | 66 +++++++++---------- 5 files changed, 52 insertions(+), 54 deletions(-) diff --git a/Project.toml b/Project.toml index 79a70ce4b4..2004614377 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "Oceananigans" uuid = "9e8cae18-63c1-5223-a75c-80ca9d6e9a09" authors = ["Climate Modeling Alliance and contributors"] -version = "0.90.8" +version = "0.90.9" [deps] Adapt = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" diff --git a/src/TurbulenceClosures/turbulence_closure_implementations/CATKEVerticalDiffusivities/CATKEVerticalDiffusivities.jl b/src/TurbulenceClosures/turbulence_closure_implementations/CATKEVerticalDiffusivities/CATKEVerticalDiffusivities.jl index 819b4e458d..38aaa5bd12 100644 --- a/src/TurbulenceClosures/turbulence_closure_implementations/CATKEVerticalDiffusivities/CATKEVerticalDiffusivities.jl +++ b/src/TurbulenceClosures/turbulence_closure_implementations/CATKEVerticalDiffusivities/CATKEVerticalDiffusivities.jl @@ -353,7 +353,7 @@ end # # where ω = ϵ / e ∼ √e / ℓ. - diffusivities.Lᵉ[i, j, k] = - wb_e - ω + div_Jᵉ_e + diffusivities.Lᵉ[i, j, k] = wb_e - ω + div_Jᵉ_e end end diff --git a/src/TurbulenceClosures/turbulence_closure_implementations/CATKEVerticalDiffusivities/mixing_length.jl b/src/TurbulenceClosures/turbulence_closure_implementations/CATKEVerticalDiffusivities/mixing_length.jl index 445d9dd355..1267a025a7 100644 --- a/src/TurbulenceClosures/turbulence_closure_implementations/CATKEVerticalDiffusivities/mixing_length.jl +++ b/src/TurbulenceClosures/turbulence_closure_implementations/CATKEVerticalDiffusivities/mixing_length.jl @@ -13,21 +13,21 @@ using ..TurbulenceClosures: Contains mixing length parameters for CATKE vertical diffusivity. """ Base.@kwdef struct MixingLength{FT} - Cˢ :: FT = 2.4 # Surface distance coefficient for shear length scale + Cˢ :: FT = 1.03 # Surface distance coefficient for shear length scale Cᵇ :: FT = Inf # Bottom distance coefficient for shear length scale - Cᶜc :: FT = 1.5 # Convective mixing length coefficient for tracers - Cᶜe :: FT = 1.2 # Convective mixing length coefficient for TKE - Cᵉc :: FT = 0.2 # Convective penetration 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ᵉc :: FT = 0.841 # Convective penetration mixing length coefficient for tracers Cᵉe :: FT = 0.0 # Convective penetration mixing length coefficient for TKE - Cˢᵖ :: FT = 0.14 # Sheared convective plume coefficient - Cˡᵒu :: FT = 0.19 # Shear mixing length coefficient for momentum at low Ri - Cʰⁱu :: FT = 0.086 # Shear mixing length coefficient for momentum at high Ri - Cˡᵒc :: FT = 0.2 # Shear mixing length coefficient for tracers at low Ri - Cʰⁱc :: FT = 0.045 # Shear mixing length coefficient for tracers at high Ri - Cˡᵒe :: FT = 1.9 # Shear mixing length coefficient for TKE at low Ri - Cʰⁱe :: FT = 0.57 # Shear mixing length coefficient for TKE at high Ri - CRiᵟ :: FT = 0.45 # Stability function width - CRi⁰ :: FT = 0.47 # Stability function lower Ri + 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 ##### diff --git a/src/TurbulenceClosures/turbulence_closure_implementations/CATKEVerticalDiffusivities/turbulent_kinetic_energy_equation.jl b/src/TurbulenceClosures/turbulence_closure_implementations/CATKEVerticalDiffusivities/turbulent_kinetic_energy_equation.jl index 33215c98e5..58d82af1af 100644 --- a/src/TurbulenceClosures/turbulence_closure_implementations/CATKEVerticalDiffusivities/turbulent_kinetic_energy_equation.jl +++ b/src/TurbulenceClosures/turbulence_closure_implementations/CATKEVerticalDiffusivities/turbulent_kinetic_energy_equation.jl @@ -5,12 +5,12 @@ Parameters for the evolution of oceanic turbulent kinetic energy at the O(1 m) s isotropic turbulence and diapycnal mixing. """ Base.@kwdef struct TurbulentKineticEnergyEquation{FT} - CˡᵒD :: FT = 1.1 # Dissipation length scale shear coefficient for low Ri - CʰⁱD :: FT = 0.37 # Dissipation length scale shear coefficient for high Ri - CᶜD :: FT = 0.88 # Dissipation length scale convecting layer coefficient + CˡᵒD :: FT = 2.52 # Dissipation length scale shear coefficient for low Ri + CʰⁱD :: FT = 0.614 # Dissipation length scale shear coefficient for high Ri + CᶜD :: FT = 2.55 # Dissipation length scale convecting layer coefficient CᵉD :: FT = 0.0 # Dissipation length scale penetration layer coefficient - Cᵂu★ :: FT = 1.1 # Surface shear-driven TKE flux coefficient - CᵂwΔ :: FT = 4.0 # Surface convective TKE flux coefficient + Cᵂu★ :: FT = 0.138 # Surface shear-driven TKE flux coefficient + CᵂwΔ :: FT = 0.264 # Surface convective TKE flux coefficient Cᵂϵ :: FT = 1.0 # Dissipative near-bottom TKE flux coefficient end diff --git a/validation/vertical_mixing_closures/column_windy_convection.jl b/validation/vertical_mixing_closures/column_windy_convection.jl index 9394f1e5e6..cc557a7ca7 100644 --- a/validation/vertical_mixing_closures/column_windy_convection.jl +++ b/validation/vertical_mixing_closures/column_windy_convection.jl @@ -9,27 +9,28 @@ using Oceananigans.TurbulenceClosures: ConvectiveAdjustmentVerticalDiffusivity, ExplicitTimeDiscretization -##### -##### Setup simulation -##### +# Parameters +Δz = 4 # Vertical resolution +Lz = 256 # Extent of vertical domain +Nz = Int(Lz/Δz) # Vertical resolution +f₀ = 1e-4 # Coriolis parameter (s⁻¹) +N² = 1e-6 # Buoyancy gradient (s⁻²) +Jᵇ = +1e-8 # Surface buoyancy flux (m² s⁻³) +τˣ = -2e-4 # Surface kinematic momentum flux (m s⁻¹) +stop_time = 2days convective_adjustment = ConvectiveAdjustmentVerticalDiffusivity(convective_κz=0.1, convective_νz=0.01) +catke = CATKEVerticalDiffusivity() +ri_based = RiBasedVerticalDiffusivity() -grid = RectilinearGrid(size=64, z=(-256, 0), topology=(Flat, Flat, Bounded)) -coriolis = FPlane(f=1e-4) - -N² = 1e-6 -Qᵇ = +1e-8 -Qᵘ = -2e-4 # - -b_bcs = FieldBoundaryConditions(top = FluxBoundaryCondition(Qᵇ)) -u_bcs = FieldBoundaryConditions(top = FluxBoundaryCondition(Qᵘ)) +# Set up simulation -closures_to_run = [ - CATKEVerticalDiffusivity(), - #RiBasedVerticalDiffusivity(), - #convective_adjustment, - ] +grid = RectilinearGrid(size=Nz, z=(-Lz, 0), topology=(Flat, Flat, Bounded)) +coriolis = FPlane(f=f₀) +b_bcs = FieldBoundaryConditions(top = FluxBoundaryCondition(Jᵇ)) +u_bcs = FieldBoundaryConditions(top = FluxBoundaryCondition(τˣ)) +#closures_to_run = [catke, ri_based, convective_adjustment] +closures_to_run = [catke] #, ri_based, convective_adjustment] for closure in closures_to_run @@ -38,10 +39,10 @@ for closure in closures_to_run buoyancy = BuoyancyTracer(), boundary_conditions = (; b=b_bcs, u=u_bcs)) - bᵢ(x, y, z) = N² * z + bᵢ(z) = N² * z set!(model, b=bᵢ, e=1e-6) - simulation = Simulation(model, Δt=10minutes, stop_iteration=1000) + simulation = Simulation(model; Δt=2minutes, stop_time) closurename = string(nameof(typeof(closure))) @@ -50,14 +51,13 @@ for closure in closures_to_run outputs = merge(model.velocities, model.tracers, diffusivities) - simulation.output_writers[:fields] = - JLD2OutputWriter(model, outputs, - #schedule = TimeInterval(10minutes), - schedule = IterationInterval(100), - filename = "windy_convection_" * closurename, - overwrite_existing = true) + simulation.output_writers[:fields] = JLD2OutputWriter(model, outputs, + schedule = TimeInterval(20minutes), + filename = "windy_convection_" * closurename, + overwrite_existing = true) - progress(sim) = @info string("Iter: ", iteration(sim), " t: ", prettytime(sim)) + progress(sim) = @info string("Iter: ", iteration(sim), " t: ", prettytime(sim), + ", max(b): ", maximum(model.tracers.b)) simulation.callbacks[:progress] = Callback(progress, IterationInterval(100)) @info "Running a simulation of $model..." @@ -114,16 +114,16 @@ axκ = Axis(fig[1, 4], xlabel=diffusivities_label, ylabel="z (m)") xlims!(axb, -grid.Lz * N², 0) xlims!(axu, -0.1, 0.1) -xlims!(axe, -1e-4, 2e-3) -xlims!(axκ, -1e-1, 1e1) +xlims!(axe, -1e-4, 2e-4) +xlims!(axκ, -1e-1, 5e-1) colors = [:black, :blue, :red, :orange] for (i, closure) in enumerate(closures_to_run) - bn = @lift interior(b_ts[i][$n], 1, 1, :) - un = @lift interior(u_ts[i][$n], 1, 1, :) - vn = @lift interior(v_ts[i][$n], 1, 1, :) - en = @lift interior(e_ts[i][$n], 1, 1, :) + bn = @lift interior(b_ts[i][$n], 1, 1, :) + un = @lift interior(u_ts[i][$n], 1, 1, :) + vn = @lift interior(v_ts[i][$n], 1, 1, :) + en = @lift interior(e_ts[i][$n], 1, 1, :) κᶜn = @lift interior(κᶜ_ts[i][$n], 1, 1, :) κᵘn = @lift interior(κᵘ_ts[i][$n], 1, 1, :) @@ -144,9 +144,7 @@ axislegend(axκ, position=:rb) display(fig) -#= record(fig, "windy_convection.mp4", 1:Nt, framerate=24) do nn n[] = nn end -=#