From 65f36d76866904f53fb9f73a498365bcd2806e8f Mon Sep 17 00:00:00 2001 From: Gustavo Marques Date: Thu, 16 Apr 2020 10:40:02 -0600 Subject: [PATCH] Add diags for Lapl. and Bihar grid Reynolds #s grid_Re_Kh = (U sqtr(dx2))/Kh grid_Re_Kh = (U dx3)/Ah where dx2 is the harmonic mean of the squares of the grid [L2], and dx3 is the harmonic mean of the squares of the grid^(3/2) [L3] --- .../lateral/MOM_hor_visc.F90 | 33 +++++++++++++++++-- 1 file changed, 30 insertions(+), 3 deletions(-) diff --git a/src/parameterizations/lateral/MOM_hor_visc.F90 b/src/parameterizations/lateral/MOM_hor_visc.F90 index f3c593819a..1458706316 100644 --- a/src/parameterizations/lateral/MOM_hor_visc.F90 +++ b/src/parameterizations/lateral/MOM_hor_visc.F90 @@ -117,8 +117,9 @@ module MOM_hor_visc Kh_Max_xx, & !< The maximum permitted Laplacian viscosity [L2 T-1 ~> m2 s-1]. Ah_Max_xx, & !< The maximum permitted biharmonic viscosity [L4 T-1 ~> m4 s-1]. n1n2_h, & !< Factor n1*n2 in the anisotropic direction tensor at h-points - n1n1_m_n2n2_h !< Factor n1**2-n2**2 in the anisotropic direction tensor at h-points - + n1n1_m_n2n2_h, & !< Factor n1**2-n2**2 in the anisotropic direction tensor at h-points + grid_sp_h2, & !< Harmonic mean of the squares of the grid [L2 ~> m2] + grid_sp_h3 !< Harmonic mean of the squares of the grid^(3/2) [L3 ~> m3] real ALLOCABLE_, dimension(NIMEMB_PTR_,NJMEMB_PTR_) :: Kh_bg_xy !< The background Laplacian viscosity at q points [L2 T-1 ~> m2 s-1]. !! The actual viscosity may be the larger of this @@ -180,6 +181,7 @@ module MOM_hor_visc !>@{ !! Diagnostic id + integer :: id_grid_Re_Ah = -1, id_grid_Re_Kh = -1 integer :: id_diffu = -1, id_diffv = -1 integer :: id_Ah_h = -1, id_Ah_q = -1 integer :: id_Kh_h = -1, id_Kh_q = -1 @@ -304,9 +306,11 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, Kh_h, & ! Laplacian viscosity at thickness points [L2 T-1 ~> m2 s-1] max_diss_rate_h, & ! maximum possible energy dissipated by lateral friction [L2 T-3 ~> m2 s-3] FrictWork, & ! work done by MKE dissipation mechanisms [R L2 T-3 ~> W m-2] - FrictWork_GME, & ! work done by GME [R L2 T-3 ~> W m-2] + FrictWork_GME, & ! work done by GME [R L2 T-3 ~> W m-2] div_xx_h ! horizontal divergence [T-1 ~> s-1] real, dimension(SZI_(G),SZJ_(G),SZK_(G)) :: & + grid_Re_Kh, & !< Grid Reynolds number for Laplacian horizontal viscosity at h points [nondim] + grid_Re_Ah, & !< Grid Reynolds number for Biharmonic horizontal viscosity at h points [nondim] GME_coeff_h !< GME coeff. at h-points [L2 T-1 ~> m2 s-1] real :: Ah ! biharmonic viscosity [L4 T-1 ~> m4 s-1] real :: Kh ! Laplacian viscosity [L2 T-1 ~> m2 s-1] @@ -842,6 +846,12 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, endif if ((CS%id_Kh_h>0) .or. find_FrictWork .or. CS%debug) Kh_h(i,j,k) = Kh + + if (CS%id_grid_Re_Kh>0) then + KE = 0.125*((u(I,j,k)+u(I-1,j,k))**2 + (v(i,J,k)+v(i,J-1,k))**2) + grid_Re_Kh(i,j,k) = (sqrt(KE) * sqrt(CS%grid_sp_h2(i,j)))/Kh + endif + if (CS%id_div_xx_h>0) div_xx_h(i,j,k) = div_xx(i,j) str_xx(i,j) = -Kh * sh_xx(i,j) @@ -890,6 +900,11 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, if ((CS%id_Ah_h>0) .or. find_FrictWork .or. CS%debug) Ah_h(i,j,k) = Ah + if (CS%id_grid_Re_Ah>0) then + KE = 0.125*((u(I,j,k)+u(I-1,j,k))**2 + (v(i,J,k)+v(i,J-1,k))**2) + grid_Re_Ah(i,j,k) = (sqrt(KE) * CS%grid_sp_h3(i,j))/Ah + endif + str_xx(i,j) = str_xx(i,j) + Ah * & (CS%DY_dxT(i,j) * (G%IdyCu(I,j)*Del2u(I,j) - G%IdyCu(I-1,j)*Del2u(I-1,j)) - & CS%DX_dyT(i,j) * (G%IdxCv(i,J)*Del2v(i,J) - G%IdxCv(i,J-1)*Del2v(i,J-1))) @@ -1295,10 +1310,12 @@ subroutine horizontal_viscosity(u, v, h, diffu, diffv, MEKE, VarMix, G, GV, US, if (CS%id_FrictWork>0) call post_data(CS%id_FrictWork, FrictWork, CS%diag) if (CS%id_FrictWork_GME>0) call post_data(CS%id_FrictWork_GME, FrictWork_GME, CS%diag) if (CS%id_Ah_h>0) call post_data(CS%id_Ah_h, Ah_h, CS%diag) + if (CS%id_grid_Re_Ah>0) call post_data(CS%id_grid_Re_Ah, grid_Re_Ah, CS%diag) if (CS%id_div_xx_h>0) call post_data(CS%id_div_xx_h, div_xx_h, CS%diag) if (CS%id_vort_xy_q>0) call post_data(CS%id_vort_xy_q, vort_xy_q, CS%diag) if (CS%id_Ah_q>0) call post_data(CS%id_Ah_q, Ah_q, CS%diag) if (CS%id_Kh_h>0) call post_data(CS%id_Kh_h, Kh_h, CS%diag) + if (CS%id_grid_Re_Kh>0) call post_data(CS%id_grid_Re_Kh, grid_Re_Kh, CS%diag) if (CS%id_Kh_q>0) call post_data(CS%id_Kh_q, Kh_q, CS%diag) if (CS%id_GME_coeff_h > 0) call post_data(CS%id_GME_coeff_h, GME_coeff_h, CS%diag) if (CS%id_GME_coeff_q > 0) call post_data(CS%id_GME_coeff_q, GME_coeff_q, CS%diag) @@ -1658,6 +1675,7 @@ subroutine hor_visc_init(Time, G, US, param_file, diag, CS, MEKE) ALLOC_(CS%dx_dyBu(IsdB:IedB,JsdB:JedB)) ; CS%dx_dyBu(:,:) = 0.0 ALLOC_(CS%dy_dxBu(IsdB:IedB,JsdB:JedB)) ; CS%dy_dxBu(:,:) = 0.0 if (CS%Laplacian) then + ALLOC_(CS%grid_sp_h2(isd:ied,jsd:jed)) ; CS%grid_sp_h2(:,:) = 0.0 ALLOC_(CS%Kh_bg_xx(isd:ied,jsd:jed)) ; CS%Kh_bg_xx(:,:) = 0.0 ALLOC_(CS%Kh_bg_xy(IsdB:IedB,JsdB:JedB)) ; CS%Kh_bg_xy(:,:) = 0.0 if (CS%bound_Kh .or. CS%better_bound_Kh) then @@ -1710,6 +1728,7 @@ subroutine hor_visc_init(Time, G, US, param_file, diag, CS, MEKE) ALLOC_(CS%Idxdy2v(isd:ied,JsdB:JedB)) ; CS%Idxdy2v(:,:) = 0.0 ALLOC_(CS%Ah_bg_xx(isd:ied,jsd:jed)) ; CS%Ah_bg_xx(:,:) = 0.0 ALLOC_(CS%Ah_bg_xy(IsdB:IedB,JsdB:JedB)) ; CS%Ah_bg_xy(:,:) = 0.0 + ALLOC_(CS%grid_sp_h3(IsdB:IedB,JsdB:JedB)); CS%grid_sp_h3(:,:) = 0.0 if (CS%bound_Ah .or. CS%better_bound_Ah) then ALLOC_(CS%Ah_Max_xx(isd:ied,jsd:jed)) ; CS%Ah_Max_xx(:,:) = 0.0 ALLOC_(CS%Ah_Max_xy(IsdB:IedB,JsdB:JedB)) ; CS%Ah_Max_xy(:,:) = 0.0 @@ -1777,6 +1796,7 @@ subroutine hor_visc_init(Time, G, US, param_file, diag, CS, MEKE) do j=Jsq,Jeq+1 ; do i=Isq,Ieq+1 ! Static factors in the Smagorinsky and Leith schemes grid_sp_h2 = (2.0*CS%dx2h(i,j)*CS%dy2h(i,j)) / (CS%dx2h(i,j) + CS%dy2h(i,j)) + CS%grid_sp_h2(i,j) = grid_sp_h2 grid_sp_h3 = grid_sp_h2*sqrt(grid_sp_h2) if (CS%Smagorinsky_Kh) CS%Laplac2_const_xx(i,j) = Smag_Lap_const * grid_sp_h2 if (CS%Leith_Kh) CS%Laplac3_const_xx(i,j) = Leith_Lap_const * grid_sp_h3 @@ -1837,6 +1857,7 @@ subroutine hor_visc_init(Time, G, US, param_file, diag, CS, MEKE) do j=Jsq,Jeq+1 ; do i=Isq,Ieq+1 grid_sp_h2 = (2.0*CS%dx2h(i,j)*CS%dy2h(i,j)) / (CS%dx2h(i,j)+CS%dy2h(i,j)) grid_sp_h3 = grid_sp_h2*sqrt(grid_sp_h2) + CS%grid_sp_h3(i,j) = grid_sp_h3 if (CS%Smagorinsky_Ah) then CS%Biharm_const_xx(i,j) = Smag_bi_const * (grid_sp_h2 * grid_sp_h2) if (CS%bound_Coriolis) then @@ -1979,6 +2000,8 @@ subroutine hor_visc_init(Time, G, US, param_file, diag, CS, MEKE) cmor_standard_name='ocean_momentum_xy_biharmonic_diffusivity') CS%id_Ah_q = register_diag_field('ocean_model', 'Ahq', diag%axesBL, Time, & 'Biharmonic Horizontal Viscosity at q Points', 'm4 s-1', conversion=US%L_to_m**4*US%s_to_T) + CS%id_grid_Re_Ah = register_diag_field('ocean_model', 'grid_Re_Ah', diag%axesTL, Time, & + 'Grid Reynolds number for the Biharmonic horizontal viscosity at h points', 'nondim') endif if (CS%Laplacian) then CS%id_Kh_h = register_diag_field('ocean_model', 'Khh', diag%axesTL, Time, & @@ -1988,6 +2011,8 @@ subroutine hor_visc_init(Time, G, US, param_file, diag, CS, MEKE) cmor_standard_name='ocean_momentum_xy_laplacian_diffusivity') CS%id_Kh_q = register_diag_field('ocean_model', 'Khq', diag%axesBL, Time, & 'Laplacian Horizontal Viscosity at q Points', 'm2 s-1', conversion=US%L_to_m**2*US%s_to_T) + CS%id_grid_Re_Kh = register_diag_field('ocean_model', 'grid_Re_Kh', diag%axesTL, Time, & + 'Grid Reynolds number for the Laplacian horizontal viscosity at h points', 'nondim') if (CS%Leith_Kh) then CS%id_vort_xy_q = register_diag_field('ocean_model', 'vort_xy_q', diag%axesBL, Time, & 'Vertical vorticity at q Points', 's-1', conversion=US%s_to_T) @@ -2105,6 +2130,7 @@ subroutine hor_visc_end(CS) endif if (CS%Laplacian) then DEALLOC_(CS%Kh_bg_xx) ; DEALLOC_(CS%Kh_bg_xy) + DEALLOC_(CS%grid_sp_h2) if (CS%bound_Kh) then DEALLOC_(CS%Kh_Max_xx) ; DEALLOC_(CS%Kh_Max_xy) endif @@ -2116,6 +2142,7 @@ subroutine hor_visc_end(CS) endif endif if (CS%biharmonic) then + DEALLOC_(CS%grid_sp_h3) DEALLOC_(CS%Idx2dyCu) ; DEALLOC_(CS%Idx2dyCv) DEALLOC_(CS%Idxdy2u) ; DEALLOC_(CS%Idxdy2v) DEALLOC_(CS%Ah_bg_xx) ; DEALLOC_(CS%Ah_bg_xy)