Add in documentation for user defined vertical diffusivity profile in the hydrostatic model about oceananigans.jl HOT 3 OPEN

julia> using Oceananigans

help?> ScalarDiffusivity
search: ScalarDiffusivity VerticalScalarDiffusivity HorizontalScalarDiffusivity ScalarBiharmonicDiffusivity VerticalScalarBiharmonicDiffusivity HorizontalScalarBiharmonicDiffusivity

  ScalarDiffusivity(time_discretization = ExplicitTimeDiscretization(),
                    formulation = ThreeDimensionalFormulation(), FT = Float64;
                    ν = 0,
                    κ = 0,
                    discrete_form = false,
                    loc = (nothing, nothing, nothing),
                    parameters = nothing)

  Return ScalarDiffusivity turbulence closure with viscosity ν and tracer diffusivities κ for each tracer field in tracers. If a single κ is provided, it is applied to all tracers. Otherwise κ
  must be a NamedTuple with values for every tracer individually.

  Arguments
  ≡≡≡≡≡≡≡≡≡

    •  time_discretization: either ExplicitTimeDiscretization() (default) or VerticallyImplicitTimeDiscretization().

    •  formulation:
       • HorizontalFormulation() for diffusivity applied in the horizontal direction(s)
       • VerticalFormulation() for diffusivity applied in the vertical direction,
       • ThreeDimensionalFormulation() (default) for diffusivity applied isotropically to all directions

    •  FT: the float datatype (default: Float64)

  Keyword arguments
  ≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡≡

    •  ν: Viscosity. Number, AbstractArray, Field, or Function.

    •  κ: Diffusivity. Number, AbstractArray, Field, Function, or NamedTuple of diffusivities with entries for each tracer.

    •  discrete_form: Boolean; default: false.

  When prescribing the viscosities or diffusivities as functions, depending on the value of keyword argument discrete_form, the constructor expects:

    •  discrete_form = false (default): functions of the grid's native coordinates and time, e.g., (x, y, z, t) for a RectilinearGrid or (λ, φ, z, t) for a LatitudeLongitudeGrid.

    •  discrete_form = true:
       • with loc = (nothing, nothing, nothing) (default): functions of (i, j, k, grid, ℓx, ℓy, ℓz) with ℓx, ℓy, and ℓz either Face() or Center().
       • with loc = (ℓx, ℓy, ℓz) with ℓx, ℓy, and ℓz either Face() or Center(): functions of (i, j, k, grid).

    •  parameters: NamedTuple with parameters used by the functions that compute viscosity and/or diffusivity; default: nothing.

   julia> using Oceananigans

  julia> ScalarDiffusivity(ν=1000, κ=2000)
  ScalarDiffusivity{ExplicitTimeDiscretization}(ν=1000.0, κ=2000.0)

  julia> const depth_scale = 100;

  julia> @inline ν(x, y, z) = 1000 * exp(z / depth_scale)
  ν (generic function with 1 method)

  julia> ScalarDiffusivity(ν=ν)
  ScalarDiffusivity{ExplicitTimeDiscretization}(ν=ν (generic function with 1 method), κ=0.0)

  julia> using Oceananigans.Grids: znode

  julia> @inline function κ(i, j, k, grid, ℓx, ℓy, ℓz, clock, fields)
             z = znode(i, j, k, grid, ℓx, ℓy, ℓz)
             return 2000 * exp(z / depth_scale)
         end
  κ (generic function with 1 method)

  julia> ScalarDiffusivity(κ=κ, discrete_form=true)
  ScalarDiffusivity{ExplicitTimeDiscretization}(ν=0.0, κ=Oceananigans.TurbulenceClosures.DiscreteDiffusionFunction{Nothing, Nothing, Nothing, Nothing, typeof(κ)})

  julia> @inline function another_κ(i, j, k, grid, clock, fields, p)
             z = znode(i, j, k, grid, Center(), Center(), Face())
             return 2000 * exp(z / p.depth_scale)
         end
  another_κ (generic function with 1 method)

  julia> ScalarDiffusivity(κ=another_κ, discrete_form=true, loc=(Center, Center, Face), parameters=(; depth_scale = 120.0))
  ScalarDiffusivity{ExplicitTimeDiscretization}(ν=0.0, κ=Oceananigans.TurbulenceClosures.DiscreteDiffusionFunction{Center, Center, Face, @NamedTuple{depth_scale::Float64}, typeof(another_κ)})