[MITgcm-support] Higher-Order Boundary condition

[Baylor Fox-Kemper]

Hi Joe (and mitgcm-support),
Well, I followed the adage, 'Use the force, read the source!', and  
have discovered the following boundary conditions are supposed to be  
applied to the momentum equations:

First, without no_slip_side, we find the following situation:

  Vector Invariant Form:

harmonic drag terms:
uD2 = viscAhD*cosFacU(j,bi,bj)*( Dij-Dmj )*recip_DXC(i,j,bi,bj) -  
viscAhZ*_recip_hFacW(i,j,k,bi,bj)*( Zip-Zij )*recip_DYG(i,j,bi,bj)
vD2 = viscAhZ*_recip_hFacS(i,j,k,bi,bj)*cosFacV(j,bi,bj)*( Zpj-Zij ) 
*recip_DXG(i,j,bi,bj)+viscAhD*( Dij-Dim )*recip_DYC(i,j,bi,bj)

where Dij is the Del . U_h and Zij is the Del x U_h.

Zij is set to zero in the boundaries, and down-weighted for partial  
cells (it is volume-weighted by hfacZ).  That is, for vertical walls,  
slip conditions are assumed as the vorticity is set to zero on the  
boundary.

Dij may or may not be down-weighted for partial cells (depending on  
the flag hDivScheme), but since the normal velocity is always zero at  
a boundary, it is not necessary to specify anything else about the  
divergence.

The magnitude of the friction has effects from the cosine factors,  
which depend on the latitude and the flag isotropic_cos_scaling, but  
the boundary condition is unaffected (I think).

biharmonic:
uD4=viscA4D*cosFacU(j,bi,bj)*( Dij-Dmj )*recip_DXC(i,j,bi,bj)- 
viscA4Z*_recip_hFacW(i,j,k,bi,bj)*( Zip-Zij )*recip_DYG(i,j,bi,bj) )

This time, though Dij is Del . Del^2 U_h and Zij is Del x Del^2 U_h.

The same boundary conditions are asserted for Del^2 U_h as for U_h.   
This is true whether the viscosity is constant or variable (e.g.,  
Smagorinsky).

Thus, the boundary conditions for no_slip_side=no_slip_bottom=.false.:

Del x U_h=0 and Del x Del^2 U_h=0.  These are called **slip** and  
**superslip** by Pedlosky.

It would be relatively easy to implement hyperslip, but this  
currently isn't done.

  Flux Form:

xViscFluxV(i,j) =_dyU(i,j,bi,bj)*drF(k)*hFacZ(i,j)*(-viscAh_Z(i,j)* 
(vFld(i,j)-vFld(i-1,j))*cosFacV(J,bi,bj)+viscA4_Z(i,j)*(del2v(i,j)- 
del2v(i-1,j))

yViscFluxU(i,j) =_dxV(i,j,bi,bj)*drF(k)*hFacZ(i,j)*(-viscAh_Z(i,j)* 
(uFld(i,j)-uFld(i,j-1))+viscA4_Z(i,j)*(del2u(i,j)-del2u(i,j-1))

Now, the viscous momentum flux is down-weighted and set to zero at  
the boundaries (where hFacZ is zero).

So, we have the same situation now, except implemented a bit  
differently, as instead of making the viscous drag (flux divergence)  
sense the boundaries by setting the vorticity and Del^2 vorticity to  
zero on the boundary, we set the flux through the boundary to zero  
which is later differentiated.  Again, though, these are equivalent  
to **slip** and **superslip**.

Without no_slip_bottom, Bottom boundary conditions are set the same  
way as in the fluxform above (whether or not the vector-invariant  
equations are used, which is a bit of a cheat, because this means the  
boundary conditions are not really 'vector-invariant', since they are  
implemented differently in z than in x and y).  Also, there is no  
biharmonic vertical viscosity, so no higher-order boundary conditions  
are possible.

When no_slip_side is used, a side drag is applied, whether for  
fluxform or vector invariant (no-slip is not implemented as a flux).

vDragTerms(i,j) =-_recip_hFacS(i,j,k,bi,bj)*recip_drF(k)*recip_rAs 
(i,j,bi,bj)*( hFacZClosedW*_dyU( i ,j,bi,bj)*_recip_dxV( i ,j,bi,bj)
      &         *( viscAh_Z(i  ,j)*vFld(i,j)*cosFacV(j,bi,bj)-viscA4_Z 
(i+1,j)*del2v(i,j)*cosFacV(j,bi,bj) ))*drF(k)*2
uDragTerms(i,j) =-_recip_hFacW(i,j,k,bi,bj)*recip_drF(k)*recip_rAw 
(i,j,bi,bj)*( hFacZClosedS*_dxV(i, j ,bi,bj)*_recip_dyU(i, j ,bi,bj)
      &         *( viscAh_Z(i,j  )*uFld(i,j)-viscA4_Z(i,j  )*del2u 
(i,j) ) )*drF(k)*2

So, when the cell is closed to the west, a drag is applied on v,  
which is 2*(Ah V-A4 Del^2 V)/dxV, which is equivalent to setting V to  
zero on the cell face to the west to zero (the 2 accounts for the  
fact that V is in the middle of the cell in x, but the east and west  
coasts are on the cell face in x).  My interpretation of this is that  
the boundary conditions are V=0 and Del^2 V=0 on the coasts, which is  
is no-slip and d^2V/dx^2=0 (super-no-slip?).

When no_slip_bottom is used, a bottom drag is applied at the bottom,  
but since there is no biharmonic vertical viscosity, only no-slip is  
applied.


So, in summary, it looks like the boundary conditions are either

1) slip (V_x=0) and super-slip (V_xxx=0) with no_slip_sides=.false.

or

2) no-slip (V=0) and super-no-slip (V_xx=0) with no_slip_sides=.true.

The bottom boundary conditions are either slip or no-slip.

It seems like these continue to hold even when the more obscure flags  
(ISOTROPIC_COS_SCALING+COSINEMETH_III) are used, which had the  
possibility of changing to something like hyperslip depending on how  
they were done.

Cheers,
    -Baylor