Skip to content
Closed
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
46 changes: 42 additions & 4 deletions src/parameterizations/lateral/MOM_hor_visc.F90
Original file line number Diff line number Diff line change
Expand Up @@ -127,6 +127,8 @@ module MOM_hor_visc
logical :: bound_Coriolis ! If true & SMAGORINSKY_AH is used, the biharmonic
! viscosity is modified to include a term that
! scales quadratically with the velocity shears.
logical :: side_drag ! Use side drag (Deremble et al., 2011)
real :: Cd ! Drag coefficient for side drag, nondimensional.

real ALLOCABLE_, dimension(NIMEM_,NJMEM_) :: &
Kh_bg_xx, &! The background Laplacian viscosity at h points, in units
Expand Down Expand Up @@ -244,7 +246,8 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, CS, OBC)

real, dimension(SZIB_(G),SZJB_(G)) :: &
sh_xy, & ! horizontal shearing strain (du/dy + dv/dx) (1/sec) including metric terms
str_xy, & ! str_xy is the cross term in the stress tensor (H m2 s-2)
str_xy, & ! str_xy and str_yx are cross terms in the stress tensor (H m2 s-2)
str_yx, &
FrictWorkIntz ! depth integrated energy dissipated by lateral friction (W/m2)

real, dimension(SZIB_(G),SZJB_(G),SZK_(G)) :: &
Expand Down Expand Up @@ -277,6 +280,9 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, CS, OBC)
logical :: rescale_Kh
logical :: find_FrictWork
logical :: apply_OBC = .false.
real :: vel_gam ! drag scaling factor in H m s-1.
real :: u_loc, v_loc, uvmag ! velocity at boundary and magnitude in m s-1.
real :: drag_stress_u, drag_stress_v ! side drag stresses in H m2 s-2.
logical :: use_MEKE_Ku
integer :: is, ie, js, je, Isq, Ieq, Jsq, Jeq, nz
integer :: i, j, k
Expand Down Expand Up @@ -525,7 +531,30 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, CS, OBC)
endif ! biharmonic

str_xy(I,J) = str_xy(I,J) * (hq * G%mask2dBu(I,J) * CS%reduction_xy(I,J))
enddo ; enddo

! stress tensor is symmetric except when side drag is applied
str_yx(i,j) = str_xy(i,j)

if (CS%side_drag) then
! since drag is applied on mask2dBu, at most one of the surrounding
! velocity components is non-zero
u_loc = u(i,j,k) + u(i,j+1,k)
v_loc = v(i,j,k) + v(i+1,j,k)
uvmag = sqrt(u_loc**2 + v_loc**2)
vel_gam = hq * CS%Cd * uvmag

! calculate the stress due to drag from the wall. we need to take into
! account the orientation of the wall relative to the current cell
! corner on which the stress is being calculated
drag_stress_u = vel_gam * u_loc * (2.0 * G%mask2dT(i,j) - 1.0)
drag_stress_v = vel_gam * v_loc * (2.0 * G%mask2dT(i,j) - 1.0)

! modify the str_xy and str_yx terms to include a drag
! when next to a wall
str_xy(i,j) = str_xy(i,j) + drag_stress_u*(1.0 - G%mask2dBu(i,j))
str_yx(i,j) = str_yx(i,j) + drag_stress_v*(1.0 - G%mask2dBu(i,j))
endif
enddo; enddo

do j=js,je ; do i=isq,ieq
! Evaluate 1/h x.Div(h Grad u) or the biharmonic equivalent.
Expand All @@ -543,8 +572,8 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, CS, OBC)

! Evaluate 1/h y.Div(h Grad u) or the biharmonic equivalent.
do J=Jsq,Jeq ; do i=is,ie
diffv(i,J,k) = ((G%IdyCv(i,J)*(CS%DY2q(I-1,J)*str_xy(I-1,J) - &
CS%DY2q(I,J) *str_xy(I,J)) - &
diffv(i,J,k) = ((G%IdyCv(i,J)*(CS%DY2q(I-1,J)*str_yx(I-1,J) - &
CS%DY2q(I,J) *str_yx(I,J)) - &
G%IdxCv(i,J)*(CS%DX2h(i,j) *str_xx(i,j) - &
CS%DX2h(i,j+1)*str_xx(i,j+1))) * &
G%IareaCv(i,J)) / (0.5*(h(i,j+1,k) + h(i,j,k)) + h_neglect)
Expand Down Expand Up @@ -795,6 +824,15 @@ subroutine hor_visc_init(Time, G, param_file, diag, CS)
"is strongly encouraged, and no slip BCs are not used with \n"//&
"the biharmonic viscosity.", default=.false.)

call get_param(param_file, mod, "SIDE_DRAG", CS%side_drag, &
"If true, a drag is applied to horizontal boundaries \n"//&
"proportional to velocity through SIDE_DRAG_CD.", &
default=.false.)
if (CS%side_drag) &
call get_param(param_file, mod, "SIDE_DRAG_CD", CS%Cd, &
"Drag coefficient to use for SIDE_DRAG.", &
units="nondim", default=0.0)

if (CS%bound_Kh .or. CS%bound_Ah .or. CS%better_bound_Kh .or. CS%better_bound_Ah) &
call get_param(param_file, mod, "DT", dt, &
"The (baroclinic) dynamics time step.", units = "s", &
Expand Down