From 63c8f970abacab0821a3d642df5cbb81f3e016f5 Mon Sep 17 00:00:00 2001
From: Simone Silvestri <silvestri.simone0@gmail.com>
Date: Sun, 19 Jan 2025 12:22:08 -0500
Subject: [PATCH] try out correct settings

---
 docs/Project.toml                 |  2 +-
 examples/one_degree_simulation.jl | 34 +++++++++++++++++++++----------
 2 files changed, 24 insertions(+), 12 deletions(-)

diff --git a/docs/Project.toml b/docs/Project.toml
index 3bf76c58..eb3d12a1 100644
--- a/docs/Project.toml
+++ b/docs/Project.toml
@@ -12,6 +12,6 @@ SeawaterPolynomials = "d496a93d-167e-4197-9f49-d3af4ff8fe40"
 CairoMakie = "0.10.12, 0.11, 0.12"
 DataDeps = "0.7"
 Documenter = "1"
-Oceananigans = "0.95.5 - 0.99"
+Oceananigans = "0.95.7 - 0.99"
 OrthogonalSphericalShellGrids = "0.2.1"
 SeawaterPolynomials = "0.3.5"
diff --git a/examples/one_degree_simulation.jl b/examples/one_degree_simulation.jl
index 0d216226..5354b506 100644
--- a/examples/one_degree_simulation.jl
+++ b/examples/one_degree_simulation.jl
@@ -54,7 +54,7 @@ closure = (eddy_closure, vertical_mixing)
 ##### Restoring
 #####
 
-restoring_rate  = 1 / 2days
+restoring_rate  = 1 / 10days
 z_below_surface = z_faces[end-1]
 
 mask = LinearlyTaperedPolarMask(southern=(-80, -70), northern=(70, 90), z=(z_below_surface, 0))
@@ -71,19 +71,31 @@ forcing = (T=FT, S=FS)
 ##### Ocean simulation
 ##### 
 
-# No numerical-only closure
-momentum_advection = WENOVectorInvariant(vorticity_order=5)
-tracer_advection   = WENO()
+momentum_advection = WENOVectorInvariant(vorticity_order=3)
+tracer_advection   = Centered()
+
+using Oceananigans.TurbulenceClosures: IsopycnalSkewSymmetricDiffusivity,
+                                       ExplicitTimeDiscretization,
+                                       DiffusiveFormulation
+
+using Oceananigans.TurbulenceClosures.TKEBasedVerticalDiffusivities: CATKEVerticalDiffusivity
+
+numerical_closure = HorizontalScalarDiffusivity(ν=5e3)
+eddy_closure = IsopycnalSkewSymmetricDiffusivity(κ_skew=1e3, κ_symmetric=1e3, skew_flux_formulation=DiffusiveFormulation())
+vertical_mixing = CATKEVerticalDiffusivity()
+
+closure = (eddy_closure, numerical_closure, vertical_mixing)
 
 # Spacings still don't work correctly?
 free_surface = SplitExplicitFreeSurface(grid; substeps=30)
 
 # Should we add a side drag since this is at a coarser resolution?
-ocean = ocean_simulation(grid; momentum_advection, tracer_advection, free_surface,
-                         closure, forcing)
-
-set!(ocean.model, T=ECCOMetadata(:temperature; dates=first(dates)),
-                  S=ECCOMetadata(:salinity;    dates=first(dates)))
+ocean = ocean_simulation(grid;
+                         momentum_advection,
+                         tracer_advection,
+                         closure,
+                         forcing,
+                         free_surface)
 
 #####
 ##### Atmospheric forcing
@@ -97,7 +109,7 @@ atmosphere = JRA55PrescribedAtmosphere(arch; backend=JRA55NetCDFBackend(20))
 #####
 
 coupled_model = OceanSeaIceModel(ocean; atmosphere, radiation) 
-simulation = Simulation(coupled_model; Δt=1minutes, stop_time=30days)
+simulation = Simulation(coupled_model; Δt=5minutes, stop_time=10days)
 
 #####
 ##### Run it!
@@ -139,7 +151,7 @@ ocean.output_writers[:surface] = JLD2OutputWriter(ocean.model, outputs;
 
 run!(simulation)
 
-simulation.Δt = 25minutes
+simulation.Δt = 30minutes
 simulation.stop_time = 720days
 
 run!(simulation)