Add small (molecular?) background diffusivity on top of Smag-Lilly eddy diffusivity

[hdrake]

The documentation for the Smagorinsky-Lilly closure, which we use here, states that a background diffusivity (e.g. molecular) can be added on top of the eddy diffusivity produced by dividing the eddy viscosity produces by the scheme with a given Prandtl number: https://clima.github.io/OceananigansDocumentation/stable/physics/turbulence_closures/#Smagorinsky-Lilly-turbulence-closure

However, as far as I can tell from the source code, there is no way to specify a non-zero background diffusivity, and the eddy diffusivity just gets computed from the eddy viscosity: https://github.com/CliMA/Oceananigans.jl/blob/b50c092461dbe0604e7c21a8985b3220efa1af1c/src/TurbulenceClosures/turbulence_closure_implementations/Smagorinskys/smagorinsky.jl#L136-L138

@liuchihl, we should figure out if there is already a way to specify this background diffusivity and–if not–to add this as a new parameter in the Smagorinsky-Lilly closure.

[hdrake]

Okay, I actually think this capability already exists in the code, although it doesn't seem to be documented anywhere.

Based on this test script, it looks like we can set the closure to a tuple like:

closure = (SmagorinskyLilly(), ScalarDiffusivity(ν=1.05e-6, κ=1.46e-7))

instead of just closure = SmagorinskyLilly().

Presumably, this computes eddy viscosity and diffusivity using the SmagorinskyLilly LES closure and then adds the constant scalar diffusivities on top!

[liuchihl]

After a simple test, it seems to be working. Though probably writing it to the outputwriter would be the best way to check whether it works properly or not.

julia> simulation.model.closure
(SmagorinskyLilly: C=0.16, Cb=1.0, Pr=(b = 1.0,), ScalarDiffusivity{ExplicitTimeDiscretization}(ν=1.05e-6, κ=(b=1.46e-7,)))