Implementing a Flux Boundary Condition in Discrete Mode

[kenflat2]

I am implementing a temperature relaxation boundary condition at the Northern wall of my simulation in discrete mode. Do I need to multiply by the flux by the grid size in the direction normal to the boundary? In this case, I multiply by Δy because I assume Oceananigans knows about the area element of the cell but maybe doesn't explicitly understand the whole volume element in discrete mode.

function northern_temp_flux(i, k,  grid, clock, fields)
    function Tₙ(x, z)
        H = 4kilometers
        h = 1kilometers
        y = -50
        ΔT = 8
        return ΔT * (exp(z/h) - exp(-H/h)) / (1 - exp(-H/h))
    end

    x = λnode(i, grid, Center())
    z = znode(k, grid, Center())

    Δy = 1/4
    λ = 7days

    T = @inbounds fields.T[i, grid.Ny, k]

    # ATTENTION HERE
    # The Δy below is what I'm asking about
    return Δy * (T - Tₙ(x, z)) / λ
end

Any help would be appreciated!

[glwagner]

A flux has units of velocity * quantity. So what you've written at least has the correct units.

If you want to make a flux act like a forcing in one grid cell, then you multiply by the metric in the wall-normal direction.

Note there is a function yspacing(i, j, k, grid, lx, ly, lz) you can use

[kenflat2]

The temperature flux at the surface is working, just the magnitude is wrong. I multiplied the forcing by 1000 as return 1e3 * Δz * (T - Tₐ(y)) / λ and the surface relaxed quickly to the prescribed linear temperature profile.

[glwagner]

Your code looks correct so I'm not sure what the issue is. If you are interested to see, the way a flux is added to a tendency is by multiplying the flux by the cell area and dividing by volume. Therefore with a flux of units velocity * quantity, we end up with a tendency of units quantity / time-scale, which is the units for a tendency of "quantity". Besides that we use the code extensively given the fact that a flux is velocity * quantity so I think it is very unlikely there is a bug there.

Are there other physics that might cause the surface temperature to restore to T=4? (Mixing, for example). This would explain your results, because I guess it is T=4 everywhere. So if there is a competition of time-scales, that's why multiplying the flux by 1e3 shows restoring.

[kenflat2]

We do have a vertical diffusivity parameterization with kappa=3e-3 in the mixed layer, so it could be that the forcing is not strong enough to override the mixing of temperature from depth.

[kenflat2]

It is strange that this argument doesn't apply for the Northern boundary condition though since it sets the temperature to 8 at the Northern boundary within the first week, and it has the same timescale as the surface forcing.

[glwagner]

You almost certainly figured it out! Try reducing the vertical diffusivity to 0 to see. For the northern boundary, the same argument would have to be made for the horizontal diffusivity? Do you have horizontal diffusivity?