Update stratified Couette flow verification experiment for JOSS

.travis.yml

@@ -24,6 +24,6 @@ jobs:
       julia: 1.2
       os: linux
       script:
-        - julia --project=docs/ -e 'using Pkg; Pkg.instantiate(); Pkg.add(PackageSpec(name="Documenter", rev="master")); Pkg.add("Literate"); Pkg.add("PyPlot"); Pkg.develop(PackageSpec(path=pwd()))'
+        - julia --project=docs/ -e 'using Pkg; Pkg.instantiate(); Pkg.add("Documenter"); Pkg.add("Literate"); Pkg.add("PyPlot"); Pkg.develop(PackageSpec(path=pwd()))'
         - julia --project=docs/ docs/make.jl
       after_success: skip

docs/Manifest.toml

@@ -19,11 +19,9 @@ version = "0.8.1"
 
 [[Documenter]]
 deps = ["Base64", "Dates", "DocStringExtensions", "InteractiveUtils", "JSON", "LibGit2", "Logging", "Markdown", "REPL", "Test", "Unicode"]
-git-tree-sha1 = "8a7ee80972d66b68a14dfbb9757c6b52f40efe93"
-repo-rev = "master"
-repo-url = "https://github.com/JuliaDocs/Documenter.jl.git"
+git-tree-sha1 = "0be9bf63e854a2408c2ecd3c600d68d4d87d8a73"
 uuid = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
-version = "0.24.0-DEV"
+version = "0.24.2"
 
 [[InteractiveUtils]]
 deps = ["Markdown"]
@@ -38,6 +36,9 @@ version = "0.21.0"
 [[LibGit2]]
 uuid = "76f85450-5226-5b5a-8eaa-529ad045b433"
 
+[[Libdl]]
+uuid = "8f399da3-3557-5675-b5ff-fb832c97cbdb"
+
 [[Literate]]
 deps = ["Base64", "JSON", "REPL"]
 git-tree-sha1 = "707c58359f2de555ace074313baea957c3187f2b"
@@ -56,12 +57,12 @@ uuid = "a63ad114-7e13-5084-954f-fe012c677804"
 
 [[Parsers]]
 deps = ["Dates", "Test"]
-git-tree-sha1 = "ef0af6c8601db18c282d092ccbd2f01f3f0cd70b"
+git-tree-sha1 = "0139ba59ce9bc680e2925aec5b7db79065d60556"
 uuid = "69de0a69-1ddd-5017-9359-2bf0b02dc9f0"
-version = "0.3.7"
+version = "0.3.10"
 
 [[Pkg]]
-deps = ["Dates", "LibGit2", "Markdown", "Printf", "REPL", "Random", "SHA", "UUIDs"]
+deps = ["Dates", "LibGit2", "Libdl", "Logging", "Markdown", "Printf", "REPL", "Random", "SHA", "UUIDs"]
 uuid = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f"
 
 [[Printf]]

docs/make.jl

@@ -28,7 +28,7 @@ for example in examples
 end
 
 #####
-##### Build docs
+##### Build and deploy docs
 #####
 
 makedocs(
@@ -42,57 +42,38 @@ makedocs(
   sitename = "Oceananigans.jl",
      pages = [
          "Home"   => "index.md",
-         "Manual" => [
-             "Model setup" => "manual/model_setup.md",
-             "Examples" => [
-                 "One-dimensional diffusion"        => "generated/simple_diffusion.md",
-                 "Two-dimensional turbulence"       => "generated/two_dimensional_turbulence.md",
-                 "Ocean wind mixing and convection" => "generated/ocean_wind_mixing_and_convection.md",
-                 "Ocean convection with plankton"   => "generated/ocean_convection_with_plankton.md",
-                 "Internal wave"                    => "generated/internal_wave.md"
-             ],
-             "Physics" => "manual/physics.md",
-             "Numerical implementation" => [
-                 #"Overview"               => "manual/overview.md",
-                 "Pressure decomposition" => "manual/pressure_decomposition.md",
-                 "Time stepping"          => "manual/time_stepping.md",
-                 "Finite volume method"   => "manual/finite_volume.md",
-                 "Spatial operators"      => "manual/spatial_operators.md",
-                 "Boundary conditions"    => "manual/boundary_conditions.md",
-                 "Poisson solvers"        => "manual/poisson_solvers.md",
-                 "Large eddy simulation"  => "manual/large_eddy_simulation.md"
-             ],
-             "Verification" => [
-                 "Taylor-Green vortex"     => "verification/taylor_green_vortex.md",
-                 "Stratified Couette flow" => "verification/stratified_couette_flow.md"
-             ],
-             #"Appendix" => [
-                 #"Staggered grid"         => "manual/staggered_grid.md",
-                 #"Fractional step method" => "manual/fractional_step.md",
-             #],
+         "Examples" => [
+             "One-dimensional diffusion"        => "generated/simple_diffusion.md",
+             "Two-dimensional turbulence"       => "generated/two_dimensional_turbulence.md",
+             "Ocean wind mixing and convection" => "generated/ocean_wind_mixing_and_convection.md",
+             "Ocean convection with plankton"   => "generated/ocean_convection_with_plankton.md",
+             "Internal wave"                    => "generated/internal_wave.md"
+         ],
+         "Model setup" => "manual/model_setup.md",
+         "Physics" => "manual/physics.md",
+         "Numerical implementation" => [
+             #"Overview"               => "manual/overview.md",
+             "Pressure decomposition" => "manual/pressure_decomposition.md",
+             "Time stepping"          => "manual/time_stepping.md",
+             "Finite volume method"   => "manual/finite_volume.md",
+             "Spatial operators"      => "manual/spatial_operators.md",
+             "Boundary conditions"    => "manual/boundary_conditions.md",
+             "Poisson solvers"        => "manual/poisson_solvers.md",
+             "Large eddy simulation"  => "manual/large_eddy_simulation.md"
+         ],
+         "Verification experiments" => [
+             "Taylor-Green vortex"     => "verification/taylor_green_vortex.md",
+             "Stratified Couette flow" => "verification/stratified_couette_flow.md"
          ],
          "Gallery"    => "gallery.md",
-         "Benchmarks" => "benchmarks.md",
+         "Performance benchmarks" => "benchmarks.md",
          "Library"    => "library.md",
+         "Appendix" => [
+             "Staggered grid"         => "manual/staggered_grid.md",
+             "Fractional step method" => "manual/fractional_step.md",
+         ],
          "Index"      => "subject_index.md"
      ]
 )
 
 deploydocs(repo = "github.com/climate-machine/Oceananigans.jl.git")
-
-#####
-##### Delete leftover JLD2 files.
-##### See: https://github.com/climate-machine/Oceananigans.jl/issues/509
-#####
-
-const GENERATED_DIR = joinpath(@__DIR__, "build/generated")
-
-@info "Deleting leftover JLD2 files..."
-leftovers = filter(x -> occursin(".jld2", x), readdir(GENERATED_DIR))
-for fname in leftovers
-    fpath = joinpath(GENERATED_DIR, fname)
-    rm(fpath, force=true)
-    @info "Deleted: $fpath"
-end
-
-deploydocs(repo = "github.com/climate-machine/Oceananigans.jl.git")

docs/src/verification/stratified_couette_flow.md

@@ -12,6 +12,8 @@ In this test we follow \citet{Vreugdenhil18} who use stratified Couette flow as
 model. We will compare our results to theirs and the direct numerical simulation results of \citet{Deusebio15} and
 \citet{Zhou17}.
 
+## Simulation setup
+
 In our setup replicating that of \citet{Vreugdenhil18}, a domain of size ``(L_x, L_y, L_z)/h = (4\pi, 2\pi, 2)``
 is bounded in the vertical ``z``-direction by two infinite, parallel plates or walls at ``z = \pm h`` with horizontal
 velocity ``u = \pm U_w`` and constant temperature ``\theta = \pm \Theta_w`` where ``\theta`` denotes temperature in this
@@ -39,3 +41,14 @@ From here the friction Reynolds number and the Nusselt number can be defined
 \text{Re}_\tau = \frac{u_\tau h}{\nu}, \quad \text{Nu} = \frac{q_w h}{\kappa \Theta_w}
 ```
 which can be computed and compared.
+
+## Velocity and temperature
+![Velocity and temperature profiles](plots_stratified_couette_flow_stratified_couette_flow_velocity_temperature_profiles.png)
+![Velocity and temperature horizontal slices](plots_stratified_couette_flow_stratified_couette_flow_velocity_temperature_slices.png)
+
+## ``\mathrm{Re}_\tau`` and ``\mathrm{Nu}`` values
+![Re and Nu comparison plots](plots_stratified_couette_flow_stratified_couette_flow_Re_Nu_scatter.png)
+![Re and Nu time series plots](plots_stratified_couette_flow_stratified_couette_flow_Re_Nu_timeseries.png)
+
+## LES viscosity and diffusivity
+![LES profiles](plots_stratified_couette_flow_stratified_couette_flow_LES_profiles.png)

verification/stratified_couette_flow/stratified_couette_flow.jl

@@ -2,6 +2,8 @@ using Statistics, Printf
 
 using Oceananigans
 using Oceananigans.TurbulenceClosures
+using Oceananigans.OutputWriters
+using Oceananigans.Diagnostics
 
 """ Friction velocity. See equation (16) of Vreugdenhil & Taylor (2018). """
 function uτ(model, Uavg)
@@ -79,19 +81,19 @@ number `Re, Prandtl number `Pr`, and Richardson number `Ri`.
 statement and updating the time step.
 """
 function simulate_stratified_couette_flow(; Nxy, Nz, arch=GPU(), h=1, U_wall=1, Re=4250, Pr=0.7, Ri, Ni=10, end_time=1000)
-    ####
-    #### Computed parameters
-    ####
+    #####
+    ##### Computed parameters
+    #####
 
          ν = U_wall * h / Re    # From Re = U_wall h / ν
     Θ_wall = Ri * U_wall^2 / h  # From Ri = L Θ_wall / U_wall²
          κ = ν / Pr             # From Pr = ν / κ
 
     parameters = (U_wall=U_wall, Θ_wall=Θ_wall, Re=Re, Pr=Pr, Ri=Ri)
 
-    ####
-    #### Impose boundary conditions
-    ####
+    #####
+    ##### Impose boundary conditions
+    #####
 
     Tbcs = HorizontallyPeriodicBCs(    top = BoundaryCondition(Value,  Θ_wall),
                                     bottom = BoundaryCondition(Value, -Θ_wall))
@@ -102,9 +104,9 @@ function simulate_stratified_couette_flow(; Nxy, Nz, arch=GPU(), h=1, U_wall=1,
     vbcs = HorizontallyPeriodicBCs(    top = BoundaryCondition(Value, 0),
                                     bottom = BoundaryCondition(Value, 0))
 
-    ####
-    #### Non-dimensional model setup
-    ####
+    #####
+    ##### Non-dimensional model setup
+    #####
 
     model = Model(
        architecture = arch,
@@ -115,9 +117,9 @@ boundary_conditions = HorizontallyPeriodicSolutionBCs(u=ubcs, v=vbcs, T=Tbcs),
          parameters = parameters
     )
 
-    ####
-    #### Set initial conditions
-    ####
+    #####
+    ##### Set initial conditions
+    #####
 
     # Add a bit of surface-concentrated noise to the initial condition
     ε(σ, z) = σ * randn() * z/model.grid.Lz * (1 + z/model.grid.Lz)
@@ -130,9 +132,9 @@ boundary_conditions = HorizontallyPeriodicSolutionBCs(u=ubcs, v=vbcs, T=Tbcs),
 
     set_ic!(model, u=u₀, v=v₀, w=w₀, T=T₀)
 
-    ####
-    #### Print simulation banner
-    ####
+    #####
+    ##### Print simulation banner
+    #####
 
     @printf(
         """
@@ -152,9 +154,9 @@ boundary_conditions = HorizontallyPeriodicSolutionBCs(u=ubcs, v=vbcs, T=Tbcs),
              model.grid.Lx, model.grid.Ly, model.grid.Lz,
              Re, Ri, Pr, ν, κ, U_wall, Θ_wall)
 
-    ####
-    #### Set up field output writer
-    ####
+    #####
+    ##### Set up field output writer
+    #####
 
     base_dir = @sprintf("stratified_couette_flow_data_Nxy%d_Nz%d_Ri%.2f", Nxy, Nz, Ri)
     prefix = @sprintf("stratified_couette_flow_Nxy%d_Nz%d_Ri%.2f", Nxy, Nz, Ri)
@@ -183,9 +185,9 @@ boundary_conditions = HorizontallyPeriodicSolutionBCs(u=ubcs, v=vbcs, T=Tbcs),
 
     push!(model.output_writers, field_writer)
 
-    ####
-    #### Set up profile output writer
-    ####
+    #####
+    ##### Set up profile output writer
+    #####
 
     Uavg = HorizontalAverage(model.velocities.u;       return_type=Array)
     Vavg = HorizontalAverage(model.velocities.v;       return_type=Array)
@@ -207,9 +209,9 @@ boundary_conditions = HorizontallyPeriodicSolutionBCs(u=ubcs, v=vbcs, T=Tbcs),
 
     push!(model.output_writers, profile_writer)
 
-    ####
-    #### Set up statistic output writer
-    ####
+    #####
+    ##### Set up statistic output writer
+    #####
 
     Reτ = FrictionReynoldsNumber(Uavg)
      Nu = NusseltNumber(Tavg)
@@ -223,9 +225,9 @@ boundary_conditions = HorizontallyPeriodicSolutionBCs(u=ubcs, v=vbcs, T=Tbcs),
 
     push!(model.output_writers, statistics_writer)
 
-    ####
-    #### Time stepping
-    ####
+    #####
+    ##### Time stepping
+    #####
 
     wizard = TimeStepWizard(cfl=0.02, Δt=0.0001, max_change=1.1, max_Δt=0.02)
 

verification/stratified_couette_flow/stratified_couette_flow_analysis.jl

@@ -54,8 +54,8 @@ for Ri in Ris
     t   = [statistic_files[Ri]["timeseries/t/"      * i]    for i in statistic_iters[Ri]]
     ReT = [statistic_files[Ri]["timeseries/Re_tau/" * i][1] for i in statistic_iters[Ri]]
     ReB = [statistic_files[Ri]["timeseries/Re_tau/" * i][2] for i in statistic_iters[Ri]]
-    NuT = [statistic_files[Ri]["timeseries/Nu_tau/" * i][1] for i in statistic_iters[Ri]]
-    NuB = [statistic_files[Ri]["timeseries/Nu_tau/" * i][2] for i in statistic_iters[Ri]]
+    NuT = [statistic_files[Ri]["timeseries/Nu/"     * i][1] for i in statistic_iters[Ri]]
+    NuB = [statistic_files[Ri]["timeseries/Nu/"     * i][2] for i in statistic_iters[Ri]]
 
     ax1.plot(t, ReT, color=c[Ri], linestyle="-",  label="Ri = $Ri (top wall)")
     ax1.plot(t, ReB, color=c[Ri], linestyle="--", label="Ri = $Ri (bottom wall)")