Merge pull request #115 from gustavo-marques/update_meke

Updates in MEKE and hor_visc
NOAA-GFDL · Aug 2, 2019 · 35b5c76 · 35b5c76
2 parents c149c20 + 815db63
commit 35b5c76
Show file tree

Hide file tree

Showing 2 changed files with 37 additions and 54 deletions.
diff --git a/src/parameterizations/lateral/MOM_MEKE.F90 b/src/parameterizations/lateral/MOM_MEKE.F90
@@ -43,7 +43,6 @@ module MOM_MEKE
   real :: MEKE_min_gamma!< Minimum value of gamma_b^2 allowed [nondim]
   real :: MEKE_Ct       !< Coefficient in the \f$\gamma_{bt}\f$ expression [nondim]
   logical :: visc_drag  !< If true use the vertvisc_type to calculate bottom drag.
-  logical :: Jansen15_drag  !< If true use the bottom drag formulation from Jansen et al. (2015)
   logical :: MEKE_GEOMETRIC !< If true, uses the GM coefficient formulation from the GEOMETRIC
                         !! framework (Marshall et al., 2012)
   logical :: GM_src_alt !< If true, use the GM energy conversion form S^2*N^2*kappa rather
@@ -326,37 +325,22 @@ subroutine step_forward_MEKE(MEKE, h, SN_u, SN_v, visc, dt, G, GV, US, CS, hu, h
     if (use_drag_rate) then
       ! Calculate a viscous drag rate (includes BBL contributions from mean flow and eddies)
 !$OMP do
-      if (CS%Jansen15_drag) then
-        do j=js,je ; do i=is,ie
-          drag_rate(i,j) = (cdrag2/MAX(1.0,G%bathyT(i,j))) * sqrt(CS%MEKE_Uscale**2 + drag_rate_visc(i,j)**2 + &
-             2.0*bottomFac2(i,j)*MEKE%MEKE(i,j)) * 2.0 * bottomFac2(i,j)*MEKE%MEKE(i,j)
-        enddo ; enddo
-      else
-        do j=js,je ; do i=is,ie
-          drag_rate(i,j) = (Rho0 * I_mass(i,j)) * sqrt( drag_rate_visc(i,j)**2 &
-             + cdrag2 * ( max(0.0, 2.0*bottomFac2(i,j)*MEKE%MEKE(i,j)) + CS%MEKE_Uscale**2 ) )
-        enddo ; enddo
-      endif
+      do j=js,je ; do i=is,ie
+        drag_rate(i,j) = (Rho0 * I_mass(i,j)) * sqrt( drag_rate_visc(i,j)**2 &
+           + cdrag2 * ( max(0.0, 2.0*bottomFac2(i,j)*MEKE%MEKE(i,j)) + CS%MEKE_Uscale**2 ) )
+      enddo ; enddo
     endif
 
     ! First stage of Strang splitting
 !$OMP do
-    if (CS%Jansen15_drag) then
-      do j=js,je ; do i=is,ie
-        ldamping = CS%MEKE_damping + drag_rate(i,j)
-        MEKE%MEKE(i,j) =  MEKE%MEKE(i,j) - MIN(MEKE%MEKE(i,j),sdt_damp*drag_rate(i,j))
-        MEKE_decay(i,j) = ldamping*G%mask2dT(i,j)
-      enddo ; enddo
-    else
-      do j=js,je ; do i=is,ie
-        ldamping = CS%MEKE_damping + drag_rate(i,j) * bottomFac2(i,j)
-        if (MEKE%MEKE(i,j)<0.) ldamping = 0.
-        ! notice that the above line ensures a damping only if MEKE is positive,
-        ! while leaving MEKE unchanged if it is negative
-        MEKE%MEKE(i,j) =  MEKE%MEKE(i,j) / (1.0 + sdt_damp*ldamping)
-        MEKE_decay(i,j) = ldamping*G%mask2dT(i,j)
-      enddo ; enddo
-    endif
+    do j=js,je ; do i=is,ie
+      ldamping = CS%MEKE_damping + drag_rate(i,j) * bottomFac2(i,j)
+      if (MEKE%MEKE(i,j)<0.) ldamping = 0.
+      ! notice that the above line ensures a damping only if MEKE is positive,
+      ! while leaving MEKE unchanged if it is negative
+      MEKE%MEKE(i,j) =  MEKE%MEKE(i,j) / (1.0 + sdt_damp*ldamping)
+      MEKE_decay(i,j) = ldamping*G%mask2dT(i,j)
+    enddo ; enddo
 !$OMP end parallel
     if (CS%MEKE_KH >= 0.0 .or. CS%KhMEKE_FAC > 0.0 .or. CS%MEKE_K4 >= 0.0) then
       ! Update halos for lateral or bi-harmonic diffusion
@@ -502,26 +486,18 @@ subroutine step_forward_MEKE(MEKE, h, SN_u, SN_v, visc, dt, G, GV, US, CS, hu, h
         ! Recalculate the drag rate, since MEKE has changed.
         if (use_drag_rate) then
 !$OMP do
-          if (CS%Jansen15_drag) then
-            do j=js,je ; do i=is,ie
-              ldamping = CS%MEKE_damping + drag_rate(i,j)
-              MEKE%MEKE(i,j) =  MEKE%MEKE(i,j) -sdt_damp*drag_rate(i,j)
-              MEKE_decay(i,j) = ldamping*G%mask2dT(i,j)
-            enddo ; enddo
-          else
-            do j=js,je ; do i=is,ie
-              drag_rate(i,j) = (Rho0 * I_mass(i,j)) * sqrt( drag_rate_visc(i,j)**2 &
-              + cdrag2 * ( max(0.0, 2.0*bottomFac2(i,j)*MEKE%MEKE(i,j)) + CS%MEKE_Uscale**2 ) )
-            enddo ; enddo
-            do j=js,je ; do i=is,ie
-              ldamping = CS%MEKE_damping + drag_rate(i,j) * bottomFac2(i,j)
-              if (MEKE%MEKE(i,j)<0.) ldamping = 0.
-              ! notice that the above line ensures a damping only if MEKE is positive,
-              ! while leaving MEKE unchanged if it is negative
-              MEKE%MEKE(i,j) =  MEKE%MEKE(i,j) / (1.0 + sdt_damp*ldamping)
-              MEKE_decay(i,j) = ldamping*G%mask2dT(i,j)
-            enddo ; enddo
-          endif
+          do j=js,je ; do i=is,ie
+            drag_rate(i,j) = (Rho0 * I_mass(i,j)) * sqrt( drag_rate_visc(i,j)**2 &
+            + cdrag2 * ( max(0.0, 2.0*bottomFac2(i,j)*MEKE%MEKE(i,j)) + CS%MEKE_Uscale**2 ) )
+          enddo ; enddo
+          do j=js,je ; do i=is,ie
+            ldamping = CS%MEKE_damping + drag_rate(i,j) * bottomFac2(i,j)
+            if (MEKE%MEKE(i,j)<0.) ldamping = 0.
+            ! notice that the above line ensures a damping only if MEKE is positive,
+            ! while leaving MEKE unchanged if it is negative
+            MEKE%MEKE(i,j) =  MEKE%MEKE(i,j) / (1.0 + sdt_damp*ldamping)
+            MEKE_decay(i,j) = ldamping*G%mask2dT(i,j)
+          enddo ; enddo
         endif
 !$OMP do
       endif
@@ -1029,9 +1005,6 @@ logical function MEKE_init(Time, G, param_file, diag, CS, MEKE, restart_CS)
   call get_param(param_file, mdl, "MEKE_USCALE", CS%MEKE_Uscale, &
                  "The background velocity that is combined with MEKE to "//&
                  "calculate the bottom drag.", units="m s-1", default=0.0)
-  call get_param(param_file, mdl, "MEKE_JANSEN15_DRAG", CS%Jansen15_drag, &
-                 "If true, use the bottom drag formulation from Jansen et al. (2015) "//&
-                 "to calculate the drag acting on MEKE.", default=.false.)
   call get_param(param_file, mdl, "MEKE_GM_SRC_ALT", CS%GM_src_alt, &
                  "If true, use the GM energy conversion form S^2*N^2*kappa rather "//&
                  "than the streamfunction for the MEKE GM source term.", default=.false.)

diff --git a/src/parameterizations/lateral/MOM_hor_visc.F90 b/src/parameterizations/lateral/MOM_hor_visc.F90
@@ -1288,8 +1288,12 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US,
       if (k==1) then
         do j=js,je ; do i=is,ie
           MEKE%mom_src(i,j) = 0.
-          MEKE%GME_snk(i,j) = 0.
         enddo ; enddo
+        if (associated(MEKE%GME_snk)) then
+          do j=js,je ; do i=is,ie
+            MEKE%GME_snk(i,j) = 0.
+          enddo ; enddo
+        endif
       endif
       if (MEKE%backscatter_Ro_c /= 0.) then
         do j=js,je ; do i=is,ie
@@ -1321,12 +1325,18 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US,
                      (u(I,j,k)-u(I,j-1,k))*G%IdyBu(I,J-1)                             &
                     +(v(i+1,J-1,k)-v(i,J-1,k))*G%IdxBu(I,J-1) )) ) )
         enddo ; enddo
-      else
+      endif ! MEKE%backscatter
+
+      if (CS%use_GME) then
         do j=js,je ; do i=is,ie
           ! MEKE%mom_src now is sign definite because it only uses the dissipation
           MEKE%mom_src(i,j) = MEKE%mom_src(i,j) + MAX(FrictWork_diss(i,j,k), FrictWorkMax(i,j,k))
         enddo ; enddo
-      endif ! MEKE%backscatter
+      else
+        do j=js,je ; do i=is,ie
+          MEKE%mom_src(i,j) = MEKE%mom_src(i,j) + FrictWork(i,j,k)
+        enddo ; enddo
+      endif ! CS%use_GME
 
       if (CS%use_GME .and. associated(MEKE)) then
         if (associated(MEKE%GME_snk)) then