From d069e51136be35e5d0d2c515337e9ed87ae6a630 Mon Sep 17 00:00:00 2001 From: Scott Bachman Date: Wed, 17 Oct 2018 16:37:58 -0600 Subject: [PATCH] Changed constraint on Leith viscosity. --- .../lateral/MOM_hor_visc.F90 | 46 +++++++++++-------- 1 file changed, 26 insertions(+), 20 deletions(-) diff --git a/src/parameterizations/lateral/MOM_hor_visc.F90 b/src/parameterizations/lateral/MOM_hor_visc.F90 index a9e64eebec..227188e731 100644 --- a/src/parameterizations/lateral/MOM_hor_visc.F90 +++ b/src/parameterizations/lateral/MOM_hor_visc.F90 @@ -221,8 +221,8 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, CS, FrictWorkIntz, & ! depth integrated energy dissipated by lateral friction (W/m2) Leith_Kh_h, & ! Leith Laplacian viscosity at h-points (m2 s-1) Leith_Ah_h, & ! Leith bi-harmonic viscosity at h-points (m4 s-1) - pl ! Planetary number (nondim) - + pl_h ! Planetary number (nondim) + real, dimension(SZIB_(G),SZJB_(G)) :: & dvdx, dudy, & ! components in the shearing strain (s-1) sh_xy, & ! horizontal shearing strain (du/dy + dv/dx) (1/sec) including metric terms @@ -230,7 +230,8 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, CS, bhstr_xy, & ! A copy of str_xy that only contains the biharmonic contribution (H m2 s-2) vort_xy, & ! Vertical vorticity (dv/dx - du/dy) (s-1) Leith_Kh_q, & ! Leith Laplacian viscosity at q-points (m2 s-1) - Leith_Ah_q ! Leith bi-harmonic viscosity at q-points (m4 s-1) + Leith_Ah_q, & ! Leith bi-harmonic viscosity at q-points (m4 s-1) + pl_q ! Planetary number (nondim) real, dimension(SZIB_(G),SZJB_(G),SZK_(G)) :: & Ah_q, & ! biharmonic viscosity at corner points (m4/s) @@ -272,7 +273,7 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, CS, real :: RoScl ! The scaling function for MEKE source term real :: FatH ! abs(f) at h-point for MEKE source term (s-1) real :: local_strain ! Local variable for interpolating computed strain rates (s-1). - real :: beta, u_scale, epsilon, grid_sp_h2, pl_u, pl_v, mod_Leith_pl + real :: beta, u_scale, epsilon, grid_sp_h2, grid_sp_q2 real :: DY_dxBu, DX_dyBu logical :: rescale_Kh, legacy_bound logical :: find_FrictWork @@ -578,28 +579,24 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, CS, beta = sqrt( G%dF_dx(i,j)**2 + G%dF_dy(i,j)**2 ) u_scale = sqrt((0.5*(u(I,j,k)+u(I-1,j,k)))**2 + (0.5*(v(i,J,k)+v(i,J-1,k)))**2) grid_sp_h2 = (2.0*CS%DX2h(i,j)*CS%DY2h(i,j)) / (CS%DX2h(i,j) + CS%DY2h(i,j)) - pl(i,j) = beta * grid_sp_h2 / (u_scale + epsilon) + pl_h(i,j) = beta * grid_sp_h2 / (u_scale + epsilon) enddo; enddo - mod_Leith_pl = 1.0 + do J=js-2,Jeq+1 ; do I=is-2,Ieq+1 + beta = sqrt( (0.25*(G%dF_dx(i,j)+G%dF_dx(i+1,j)+G%dF_dx(i,j+1)+G%dF_dx(i+1,j+1))**2) + & + (0.25*(G%dF_dy(i,j)+G%dF_dy(i+1,j)+G%dF_dy(i,j+1)+G%dF_dy(i+1,j+1))**2) ) + u_scale = sqrt((0.5*(u(I,j,k)+u(I,j+1,k)))**2 + (0.5*(v(i,J,k)+v(i,J+1,k)))**2) + grid_sp_q2 = (2.0*CS%DX2q(I,J)*CS%DY2q(I,J)) / (CS%DX2q(I,J) + CS%DY2q(I,J)) + pl_q(I,J) = beta * grid_sp_q2 / (u_scale + epsilon) + enddo ; enddo + do J=js-2,Jeq+1 ; do i=is-1,Ieq+1 - pl_v = 0.5 * (pl(i,j) + pl(i,j+1)) - if (pl_v > 1.0) then - vort_xy_dx(i,J) = 0.5 * ( G%dF_dx(i,j) + G%dF_dx(i,j+1)) - mod_Leith_pl = 0.0 - else vort_xy_dx(i,J) = vort_xy_dx(i,J) + 0.5 * ( G%dF_dx(i,j) + G%dF_dx(i,j+1)) - endif enddo ; enddo do j=js-1,Jeq+1 ; do I=is-2,Ieq+1 - pl_u = 0.5 * (pl(i,j) + pl(i+1,j)) - if (pl_u > 1.0) then - vort_xy_dy(I,j) = 0.5 * ( G%dF_dy(i,j) + G%dF_dy(i+1,j)) - mod_Leith_pl = 0.0 - else vort_xy_dy(I,j) = vort_xy_dy(I,j) + 0.5 * ( G%dF_dy(i,j) + G%dF_dy(i+1,j)) - endif enddo ; enddo + endif mod_Leith = 0.; if (CS%modified_Leith) mod_Leith = 1.0 @@ -616,11 +613,15 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, CS, (sh_xy(I-1,J)*sh_xy(I-1,J) + sh_xy(I,J-1)*sh_xy(I,J-1)))) endif if ((CS%Leith_Kh) .or. (CS%Leith_Ah)) then - vert_vort_mag = sqrt( & + if (pl_h(i,j) > 1) then + vert_vort_mag = sqrt( G%dF_dx(i,j)**2 + G%dF_dy(i,j)**2 ) + else + vert_vort_mag = sqrt( & 0.5*((vort_xy_dx(i,J-1)*vort_xy_dx(i,J-1) + vort_xy_dx(i,J)*vort_xy_dx(i,J)) + & (vort_xy_dy(I-1,j)*vort_xy_dy(I-1,j) + vort_xy_dy(I,j)*vort_xy_dy(I,j))) + & mod_Leith*mod_Leith_pl*0.5*((div_xx_dx(I,j)*div_xx_dx(I,j) + div_xx_dx(I-1,j)* & div_xx_dx(I-1,j)) + (div_xx_dy(i,J)*div_xx_dy(i,J) + div_xx_dy(i,J-1)*div_xx_dy(i,J-1)))) + endif endif if (CS%better_bound_Ah .or. CS%better_bound_Kh) then hrat_min = min(1.0, min(h_u(I,j), h_u(I-1,j), h_v(i,J), h_v(i,J-1)) / & @@ -749,7 +750,12 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, CS, (sh_xx(i,j+1)*sh_xx(i,j+1) + sh_xx(i+1,j)*sh_xx(i+1,j)))) endif if ((CS%Leith_Kh) .or. (CS%Leith_Ah)) then - vert_vort_mag = sqrt( & + if (pl_q(I,J) > 1) then + vert_vort_mag = sqrt( & + (0.25*(G%dF_dx(i,j)+G%dF_dx(i+1,j)+G%dF_dx(i,j+1)+G%dF_dx(i+1,j+1))**2) + & + (0.25*(G%dF_dy(i,j)+G%dF_dy(i+1,j)+G%dF_dy(i,j+1)+G%dF_dy(i+1,j+1))**2) ) + else + vert_vort_mag = sqrt( & 0.5*((vort_xy_dx(i,J)*vort_xy_dx(i,J) + vort_xy_dx(i+1,J)*vort_xy_dx(i+1,J)) + & (vort_xy_dy(I,j)*vort_xy_dy(I,j) + vort_xy_dy(I,j+1)*vort_xy_dy(I,j+1))) + & mod_Leith*mod_Leith_pl*0.5*((div_xx_dx(I,j)*div_xx_dx(I,j) + div_xx_dx(I,j+1)* &