diff --git a/src/tracer/MOM_lateral_boundary_diffusion.F90 b/src/tracer/MOM_lateral_boundary_diffusion.F90 index 443b9108d2..d06127d0d5 100644 --- a/src/tracer/MOM_lateral_boundary_diffusion.F90 +++ b/src/tracer/MOM_lateral_boundary_diffusion.F90 @@ -45,6 +45,8 @@ module MOM_lateral_boundary_diffusion integer :: deg !< Degree of polynomial reconstruction integer :: surface_boundary_scheme !< Which boundary layer scheme to use !! 1. ePBL; 2. KPP + logical :: limiter !< Controls wether a flux limiter is applied. + !! Only valid when method = 1. type(remapping_CS) :: remap_CS !< Control structure to hold remapping configuration type(KPP_CS), pointer :: KPP_CSp => NULL() !< KPP control structure needed to get BLD type(energetic_PBL_CS), pointer :: energetic_PBL_CSp => NULL() !< ePBL control structure needed to get BLD @@ -103,6 +105,11 @@ logical function lateral_boundary_diffusion_init(Time, G, param_file, diag, diab "Determine how to apply boundary lateral diffusion of tracers: \n"//& "1. Bulk layer approach \n"//& "2. Along layer approach", default=1) + if (CS%method == 1) then + call get_param(param_file, mdl, "APPLY_LIMITER", CS%limiter, & + "If True, apply a flux limiter in the LBD. This is only available \n"//& + "when LATERAL_BOUNDARY_METHOD=1.", default=.false.) + endif call get_param(param_file, mdl, "LBD_BOUNDARY_EXTRAP", boundary_extrap, & "Use boundary extrapolation in LBD code", & default=.false.) @@ -187,7 +194,7 @@ subroutine lateral_boundary_diffusion(G, GV, US, h, Coef_x, Coef_y, dt, Reg, CS) call fluxes_bulk_method(SURFACE, GV%ke, CS%deg, h(I,j,:), h(I+1,j,:), hbl(I,j), hbl(I+1,j), & G%areaT(I,j), G%areaT(I+1,j), tracer%t(I,j,:), tracer%t(I+1,j,:), & ppoly0_coefs(I,j,:,:), ppoly0_coefs(I+1,j,:,:), ppoly0_E(I,j,:,:), & - ppoly0_E(I+1,j,:,:), remap_method, Coef_x(I,j), uFlx_bulk(I,j), uFlx(I,j,:)) + ppoly0_E(I+1,j,:,:), remap_method, Coef_x(I,j), uFlx_bulk(I,j), uFlx(I,j,:), CS%limiter) endif enddo enddo @@ -197,7 +204,7 @@ subroutine lateral_boundary_diffusion(G, GV, US, h, Coef_x, Coef_y, dt, Reg, CS) call fluxes_bulk_method(SURFACE, GV%ke, CS%deg, h(i,J,:), h(i,J+1,:), hbl(i,J), hbl(i,J+1), & G%areaT(i,J), G%areaT(i,J+1), tracer%t(i,J,:), tracer%t(i,J+1,:), & ppoly0_coefs(i,J,:,:), ppoly0_coefs(i,J+1,:,:), ppoly0_E(i,J,:,:), & - ppoly0_E(i,J+1,:,:), remap_method, Coef_y(i,J), vFlx_bulk(i,J), vFlx(i,J,:)) + ppoly0_E(i,J+1,:,:), remap_method, Coef_y(i,J), vFlx_bulk(i,J), vFlx(i,J,:), CS%limiter) endif enddo enddo @@ -380,7 +387,7 @@ subroutine boundary_k_range(boundary, nk, h, hbl, k_top, zeta_top, k_bot, zeta_b if (hbl == 0.) return if (hbl >= SUM(h(:))) then k_bot = nk - zeta_bot = 0. + zeta_bot = 1. return endif do k=1,nk @@ -396,12 +403,12 @@ subroutine boundary_k_range(boundary, nk, h, hbl, k_top, zeta_top, k_bot, zeta_b k_top = nk zeta_top = 1. k_bot = nk - zeta_bot = 1. + zeta_bot = 0. htot = 0. if (hbl == 0.) return if (hbl >= SUM(h(:))) then k_top = 1 - zeta_top = 0. + zeta_top = 1. return endif do k=nk,1,-1 @@ -558,8 +565,8 @@ subroutine fluxes_bulk_method(boundary, nk, deg, h_L, h_R, hbl_L, hbl_R, area_L, real, intent(in ) :: khtr_u !< Horizontal diffusivities times delta t at U-point [m^2] real, intent( out) :: F_bulk !< The bulk mixed layer lateral flux [m^3 conc] real, dimension(nk), intent( out) :: F_layer !< Layerwise diffusive flux at U-point [m^3 conc] - real, optional, dimension(nk), intent( out) :: F_limit !< The amount of flux not applied due to limiter - !! F_layer(k) - F_max [m^3 conc] + logical, optional, intent(in ) :: F_limit !< If True, apply a limiter + ! Local variables real, dimension(nk) :: h_means !< Calculate the layer-wise harmonic means [m] real :: khtr_avg !< Thickness-weighted diffusivity at the u-point [m^2 s^-1] @@ -580,8 +587,9 @@ subroutine fluxes_bulk_method(boundary, nk, deg, h_L, h_R, hbl_L, hbl_R, area_L, real :: h_work_L, h_work_R !< dummy variables real :: F_max !< The maximum amount of flux that can leave a !! cell [m^3 conc] - logical :: limited !< True if the flux limiter was applied - real :: hfrac, F_bulk_remain + logical :: limiter !< True if flux limiter should be applied + real :: hfrac !< Layer fraction wrt sum of all layers [nondim] + real :: dphi !< tracer gradient [conc m^-3] if (hbl_L == 0. .or. hbl_R == 0.) then F_bulk = 0. @@ -589,6 +597,11 @@ subroutine fluxes_bulk_method(boundary, nk, deg, h_L, h_R, hbl_L, hbl_R, area_L, return endif + limiter = .false. + if (PRESENT(F_limit)) then + limiter = F_limit + endif + ! Calculate vertical indices containing the boundary layer call boundary_k_range(boundary, nk, h_L, hbl_L, k_top_L, zeta_top_L, k_bot_L, zeta_bot_L) call boundary_k_range(boundary, nk, h_R, hbl_R, k_top_R, zeta_top_R, k_bot_R, zeta_bot_R) @@ -603,7 +616,6 @@ subroutine fluxes_bulk_method(boundary, nk, deg, h_L, h_R, hbl_L, hbl_R, area_L, ! GMM, khtr_avg should be computed once khtr is 3D heff = harmonic_mean(hbl_L, hbl_R) F_bulk = -(khtr_u * heff) * (phi_R_avg - phi_L_avg) - F_bulk_remain = F_bulk ! Calculate the layerwise sum of the vertical effective thickness. This is different than the heff calculated ! above, but is used as a way to decompose the fluxes onto the individual layers h_means(:) = 0. @@ -664,47 +676,36 @@ subroutine fluxes_bulk_method(boundary, nk, deg, h_L, h_R, hbl_L, hbl_R, area_L, if (h_means(k) > 0.) then hfrac = h_means(k)*inv_heff F_layer(k) = F_bulk * hfrac - ! limit the flux to 0.2 of the tracer *gradient* - ! Why 0.2? - ! t=0 t=inf - ! 0 .2 - ! 0 1 0 .2.2.2 - ! 0 .2 - ! - F_max = -0.2 * ((area_R*(phi_R(k)*h_R(k)))-(area_L*(phi_L(k)*h_R(k)))) - - ! check if bulk flux (or F_layer) and F_max have same direction - if ( SIGN(1.,F_bulk) == SIGN(1., F_max)) then - ! Distribute bulk flux onto layers - if ( ((boundary == SURFACE) .and. (k == k_min)) .or. ((boundary == BOTTOM) .and. (k == nk)) ) then - F_layer(k) = F_bulk_remain ! GMM, are not using F_bulk_remain for now. Should we keep it? - endif - F_bulk_remain = F_bulk_remain - F_layer(k) - ! Apply flux limiter calculated above - if (F_max >= 0.) then - limited = F_layer(k) > F_max - F_layer(k) = MIN(F_layer(k),F_max) - else - limited = F_layer(k) < F_max - F_layer(k) = MAX(F_layer(k),F_max) - endif - - ! GMM, again we are not using F_limit. Should we delete it? - if (PRESENT(F_limit)) then - if (limited) then - F_limit(k) = F_layer(k) - F_max + if (limiter) then + ! limit the flux to 0.2 of the tracer *gradient* + ! Why 0.2? + ! t=0 t=inf + ! 0 .2 + ! 0 1 0 .2.2.2 + ! 0 .2 + ! + F_max = -0.2 * ((area_R*(phi_R(k)*h_R(k)))-(area_L*(phi_L(k)*h_R(k)))) + + ! check if bulk flux (or F_layer) and F_max have same direction + if ( SIGN(1.,F_bulk) == SIGN(1., F_max)) then + ! Apply flux limiter calculated above + if (F_max >= 0.) then + F_layer(k) = MIN(F_layer(k),F_max) else - F_limit(k) = 0. + F_layer(k) = MAX(F_layer(k),F_max) endif + else + ! do not apply a flux on this layer + F_layer(k) = 0. endif else - ! do not apply a flux on this layer - F_bulk_remain = F_bulk_remain - F_layer(k) - F_layer(k) = 0. - endif - else - F_layer(k) = 0. + dphi = -(phi_R(k) - phi_L(k)) + if (.not. SIGN(1.,F_bulk) == SIGN(1., dphi)) then + ! upgradient, do not apply a flux on this layer + F_layer(k) = 0. + endif + endif ! limited endif enddo endif @@ -748,47 +749,56 @@ logical function near_boundary_unit_tests( verbose ) test_name = 'Surface boundary spans the entire top cell' h_L = (/5.,5./) call boundary_k_range(SURFACE, nk, h_L, 5., k_top, zeta_top, k_bot, zeta_bot) - near_boundary_unit_tests = test_boundary_k_range(k_top, zeta_top, k_bot, zeta_bot, 1, 0., 1, 1., test_name, verbose) + near_boundary_unit_tests = near_boundary_unit_tests .or. & + test_boundary_k_range(k_top, zeta_top, k_bot, zeta_bot, 1, 0., 1, 1., test_name, verbose) test_name = 'Surface boundary spans the entire column' h_L = (/5.,5./) call boundary_k_range(SURFACE, nk, h_L, 10., k_top, zeta_top, k_bot, zeta_bot) - near_boundary_unit_tests = test_boundary_k_range(k_top, zeta_top, k_bot, zeta_bot, 1, 0., 2, 0., test_name, verbose) + near_boundary_unit_tests = near_boundary_unit_tests .or. & + test_boundary_k_range(k_top, zeta_top, k_bot, zeta_bot, 1, 0., 2, 1., test_name, verbose) test_name = 'Bottom boundary spans the entire bottom cell' h_L = (/5.,5./) call boundary_k_range(BOTTOM, nk, h_L, 5., k_top, zeta_top, k_bot, zeta_bot) - near_boundary_unit_tests = test_boundary_k_range(k_top, zeta_top, k_bot, zeta_bot, 2, 0., 2, 1., test_name, verbose) + near_boundary_unit_tests = near_boundary_unit_tests .or. & + test_boundary_k_range(k_top, zeta_top, k_bot, zeta_bot, 2, 1., 2, 0., test_name, verbose) test_name = 'Bottom boundary spans the entire column' h_L = (/5.,5./) call boundary_k_range(BOTTOM, nk, h_L, 10., k_top, zeta_top, k_bot, zeta_bot) - near_boundary_unit_tests = test_boundary_k_range(k_top, zeta_top, k_bot, zeta_bot, 1, 0., 2, 1., test_name, verbose) + near_boundary_unit_tests = near_boundary_unit_tests .or. & + test_boundary_k_range(k_top, zeta_top, k_bot, zeta_bot, 1, 1., 2, 0., test_name, verbose) test_name = 'Surface boundary intersects second layer' h_L = (/10.,10./) call boundary_k_range(SURFACE, nk, h_L, 17.5, k_top, zeta_top, k_bot, zeta_bot) - near_boundary_unit_tests = test_boundary_k_range(k_top, zeta_top, k_bot, zeta_bot, 1, 0., 2, 0.75, test_name, verbose) + near_boundary_unit_tests = near_boundary_unit_tests .or. & + test_boundary_k_range(k_top, zeta_top, k_bot, zeta_bot, 1, 0., 2, 0.75, test_name, verbose) test_name = 'Surface boundary intersects first layer' h_L = (/10.,10./) call boundary_k_range(SURFACE, nk, h_L, 2.5, k_top, zeta_top, k_bot, zeta_bot) - near_boundary_unit_tests = test_boundary_k_range(k_top, zeta_top, k_bot, zeta_bot, 1, 0., 1, 0.25, test_name, verbose) + near_boundary_unit_tests = near_boundary_unit_tests .or. & + test_boundary_k_range(k_top, zeta_top, k_bot, zeta_bot, 1, 0., 1, 0.25, test_name, verbose) test_name = 'Surface boundary is deeper than column thickness' h_L = (/10.,10./) call boundary_k_range(SURFACE, nk, h_L, 21.0, k_top, zeta_top, k_bot, zeta_bot) - near_boundary_unit_tests = test_boundary_k_range(k_top, zeta_top, k_bot, zeta_bot, 1, 0., 2, 0., test_name, verbose) + near_boundary_unit_tests = near_boundary_unit_tests .or. & + test_boundary_k_range(k_top, zeta_top, k_bot, zeta_bot, 1, 0., 2, 1., test_name, verbose) test_name = 'Bottom boundary intersects first layer' h_L = (/10.,10./) call boundary_k_range(BOTTOM, nk, h_L, 17.5, k_top, zeta_top, k_bot, zeta_bot) - near_boundary_unit_tests = test_boundary_k_range(k_top, zeta_top, k_bot, zeta_bot, 1, 0.75, 2, 1., test_name, verbose) + near_boundary_unit_tests = near_boundary_unit_tests .or. & + test_boundary_k_range(k_top, zeta_top, k_bot, zeta_bot, 1, 0.75, 2, 0., test_name, verbose) test_name = 'Bottom boundary intersects second layer' h_L = (/10.,10./) call boundary_k_range(BOTTOM, nk, h_L, 2.5, k_top, zeta_top, k_bot, zeta_bot) - near_boundary_unit_tests = test_boundary_k_range(k_top, zeta_top, k_bot, zeta_bot, 2, 0.25, 2, 1., test_name, verbose) + near_boundary_unit_tests = near_boundary_unit_tests .or. & + test_boundary_k_range(k_top, zeta_top, k_bot, zeta_bot, 2, 0.25, 2, 0., test_name, verbose) ! All cases in this section have hbl which are equal to the column thicknesses test_name = 'Equal hbl and same layer thicknesses (gradient from right to left)' @@ -804,9 +814,17 @@ logical function near_boundary_unit_tests( verbose ) ppoly0_E_R(1,1) = 1.; ppoly0_E_R(1,2) = 1. ppoly0_E_R(2,1) = 1.; ppoly0_E_R(2,2) = 1. khtr_u = 1. + ! Without limiter call fluxes_bulk_method(SURFACE, nk, deg, h_L, h_R, hbl_L, hbl_R, area_L, area_R, phi_L, phi_R, phi_pp_L, phi_pp_R, & ppoly0_E_L, ppoly0_E_R, method, khtr_u, F_bulk, F_layer) - near_boundary_unit_tests = test_layer_fluxes( verbose, nk, test_name, F_layer, (/-5.0,-5.0/) ) + near_boundary_unit_tests = near_boundary_unit_tests .or. & + test_layer_fluxes( verbose, nk, test_name, F_layer, (/-5.0,-5.0/) ) + + ! same as above, but with limiter + call fluxes_bulk_method(SURFACE, nk, deg, h_L, h_R, hbl_L, hbl_R, area_L, area_R, phi_L, phi_R, phi_pp_L, phi_pp_R, & + ppoly0_E_L, ppoly0_E_R, method, khtr_u, F_bulk, F_layer, .true.) + near_boundary_unit_tests = near_boundary_unit_tests .or. & + test_layer_fluxes( verbose, nk, test_name, F_layer, (/-1.0,-1.0/) ) test_name = 'Equal hbl and same layer thicknesses (gradient from left to right)' hbl_L = 10.; hbl_R = 10. @@ -823,7 +841,8 @@ logical function near_boundary_unit_tests( verbose ) khtr_u = 1. call fluxes_bulk_method(SURFACE, nk, deg, h_L, h_R, hbl_L, hbl_R, area_L, area_R, phi_L, phi_R, phi_pp_L, phi_pp_R,& ppoly0_E_L, ppoly0_E_R, method, khtr_u, F_bulk, F_layer) - near_boundary_unit_tests = test_layer_fluxes( verbose, nk, test_name, F_layer, (/5.0,5.0/) ) + near_boundary_unit_tests = near_boundary_unit_tests .or. & + test_layer_fluxes( verbose, nk, test_name, F_layer, (/5.0,5.0/) ) test_name = 'Equal hbl and same layer thicknesses (no gradient)' hbl_L = 10; hbl_R = 10 @@ -833,14 +852,15 @@ logical function near_boundary_unit_tests( verbose ) phi_pp_L(2,1) = 1.; phi_pp_L(2,2) = 0. phi_pp_R(1,1) = 1.; phi_pp_R(1,2) = 0. phi_pp_R(2,1) = 1.; phi_pp_R(2,2) = 0. - ppoly0_E_L(1,1) = 1.; ppoly0_E_L(1,2) = 0. - ppoly0_E_L(2,1) = 1.; ppoly0_E_L(2,2) = 0. + ppoly0_E_L(1,1) = 1.; ppoly0_E_L(1,2) = 1. + ppoly0_E_L(2,1) = 1.; ppoly0_E_L(2,2) = 1. ppoly0_E_R(1,1) = 1.; ppoly0_E_R(1,2) = 1. ppoly0_E_R(2,1) = 1.; ppoly0_E_R(2,2) = 1. khtr_u = 1. call fluxes_bulk_method(SURFACE, nk, deg, h_L, h_R, hbl_L, hbl_R, area_L, area_R, phi_L, phi_R, phi_pp_L, phi_pp_R,& ppoly0_E_L, ppoly0_E_R, method, khtr_u, F_bulk, F_layer) - near_boundary_unit_tests = test_layer_fluxes( verbose, nk, test_name, F_layer, (/0.0,0.0/) ) + near_boundary_unit_tests = near_boundary_unit_tests .or. & + test_layer_fluxes( verbose, nk, test_name, F_layer, (/0.0,0.0/) ) test_name = 'Equal hbl and different layer thicknesses (gradient right to left)' hbl_L = 16.; hbl_R = 16. @@ -857,7 +877,8 @@ logical function near_boundary_unit_tests( verbose ) khtr_u = 1. call fluxes_bulk_method(SURFACE, nk, deg, h_L, h_R, hbl_L, hbl_R, area_L, area_R, phi_L, phi_R, phi_pp_L, phi_pp_R,& ppoly0_E_L, ppoly0_E_R, method, khtr_u, F_bulk, F_layer) - near_boundary_unit_tests = test_layer_fluxes( verbose, nk, test_name, F_layer, (/-8.0,-8.0/) ) + near_boundary_unit_tests = near_boundary_unit_tests .or. & + test_layer_fluxes( verbose, nk, test_name, F_layer, (/-8.0,-8.0/) ) test_name = 'Equal hbl and same layer thicknesses (diagonal tracer values)' hbl_L = 10.; hbl_R = 10. @@ -874,7 +895,8 @@ logical function near_boundary_unit_tests( verbose ) khtr_u = 1. call fluxes_bulk_method(SURFACE, nk, deg, h_L, h_R, hbl_L, hbl_R, area_L, area_R, phi_L, phi_R, phi_pp_L, phi_pp_R,& ppoly0_E_L, ppoly0_E_R, method, khtr_u, F_bulk, F_layer) - near_boundary_unit_tests = test_layer_fluxes( verbose, nk, test_name, F_layer, (/0.0,0.0/) ) + near_boundary_unit_tests = near_boundary_unit_tests .or. & + test_layer_fluxes( verbose, nk, test_name, F_layer, (/0.0,0.0/) ) test_name = 'Different hbl and different layer thicknesses (gradient from right to left)' hbl_L = 12; hbl_R = 20 @@ -891,7 +913,8 @@ logical function near_boundary_unit_tests( verbose ) khtr_u = 1. call fluxes_bulk_method(SURFACE, nk, deg, h_L, h_R, hbl_L, hbl_R, area_L, area_R, phi_L, phi_R, phi_pp_L, phi_pp_R,& ppoly0_E_L, ppoly0_E_R, method, khtr_u, F_bulk, F_layer) - near_boundary_unit_tests = test_layer_fluxes( verbose, nk, test_name, F_layer, (/-7.5,-7.5/) ) + near_boundary_unit_tests = near_boundary_unit_tests .or. & + test_layer_fluxes( verbose, nk, test_name, F_layer, (/-7.5,-7.5/) ) ! Cases where hbl < column thickness (polynomial coefficients specified for pseudo-linear reconstruction) @@ -903,10 +926,15 @@ logical function near_boundary_unit_tests( verbose ) phi_pp_L(2,1) = 0.; phi_pp_L(2,2) = 0. phi_pp_R(1,1) = 1.; phi_pp_R(1,2) = 0. phi_pp_R(2,1) = 1.; phi_pp_R(2,2) = 0. + ppoly0_E_L(1,1) = 0.; ppoly0_E_L(1,2) = 0. + ppoly0_E_L(2,1) = 0.; ppoly0_E_L(2,2) = 0. + ppoly0_E_R(1,1) = 1.; ppoly0_E_R(1,2) = 1. + ppoly0_E_R(2,1) = 1.; ppoly0_E_R(2,2) = 1. khtr_u = 1. call fluxes_bulk_method(SURFACE, nk, deg, h_L, h_R, hbl_L, hbl_R, area_L, area_R, phi_L, phi_R, phi_pp_L, phi_pp_R,& ppoly0_E_L, ppoly0_E_R, method, khtr_u, F_bulk, F_layer) - near_boundary_unit_tests = test_layer_fluxes( verbose, nk, test_name, F_layer, (/-1.,-1./) ) + near_boundary_unit_tests = near_boundary_unit_tests .or. & + test_layer_fluxes( verbose, nk, test_name, F_layer, (/-1.,-1./) ) test_name = 'hbl < column thickness, hbl same, linear profile right' hbl_L = 2; hbl_R = 2 @@ -923,7 +951,8 @@ logical function near_boundary_unit_tests( verbose ) ppoly0_E_R(2,1) = 1.; ppoly0_E_R(2,2) = 3. call fluxes_bulk_method(SURFACE, nk, deg, h_L, h_R, hbl_L, hbl_R, area_L, area_R, phi_L, phi_R, phi_pp_L, phi_pp_R,& ppoly0_E_L, ppoly0_E_R, method, khtr_u, F_bulk, F_layer) - near_boundary_unit_tests = test_layer_fluxes( verbose, nk, test_name, F_layer, (/-1.,-1./) ) + near_boundary_unit_tests = near_boundary_unit_tests .or. & + test_layer_fluxes( verbose, nk, test_name, F_layer, (/-1.,-1./) ) test_name = 'hbl < column thickness, hbl same, linear profile right, khtr=2' hbl_L = 2; hbl_R = 2 @@ -940,7 +969,8 @@ logical function near_boundary_unit_tests( verbose ) ppoly0_E_R(2,1) = 1.; ppoly0_E_R(2,2) = 3. call fluxes_layer_method(SURFACE, nk, deg, h_L, h_R, hbl_L, hbl_R, area_L, area_R, phi_L, phi_R, phi_pp_L, & phi_pp_R, ppoly0_E_L, ppoly0_E_R, method, khtr_u, F_layer) - near_boundary_unit_tests = test_layer_fluxes( verbose, nk, test_name, F_layer, (/-2.,-2./) ) + near_boundary_unit_tests = near_boundary_unit_tests .or. & + test_layer_fluxes( verbose, nk, test_name, F_layer, (/-1.,-3./) ) ! unit tests for layer by layer method test_name = 'Different hbl and different column thicknesses (gradient from right to left)' @@ -958,7 +988,8 @@ logical function near_boundary_unit_tests( verbose ) khtr_u = 1. call fluxes_layer_method(SURFACE, nk, deg, h_L, h_R, hbl_L, hbl_R, area_L, area_R, phi_L, phi_R, phi_pp_L, & phi_pp_R, ppoly0_E_L, ppoly0_E_R, method, khtr_u, F_layer) - near_boundary_unit_tests = test_layer_fluxes( verbose, nk, test_name, F_layer, (/-7.5,0.0/) ) + near_boundary_unit_tests = near_boundary_unit_tests .or. & + test_layer_fluxes( verbose, nk, test_name, F_layer, (/-7.5,-7.5/) ) test_name = 'Different hbl and different column thicknesses (linear profile right)' @@ -976,7 +1007,8 @@ logical function near_boundary_unit_tests( verbose ) khtr_u = 1. call fluxes_layer_method(SURFACE, nk, deg, h_L, h_R, hbl_L, hbl_R, area_L, area_R, phi_L, phi_R, phi_pp_L, & phi_pp_R, ppoly0_E_L, ppoly0_E_R, method, khtr_u, F_layer) - near_boundary_unit_tests = test_layer_fluxes( verbose, nk, test_name, F_layer, (/-3.75,0.0/) ) + near_boundary_unit_tests = near_boundary_unit_tests .or. & + test_layer_fluxes( verbose, nk, test_name, F_layer, (/-3.75,0.0/) ) end function near_boundary_unit_tests !> Returns true if output of near-boundary unit tests does not match correct computed values diff --git a/src/tracer/MOM_tracer_registry.F90 b/src/tracer/MOM_tracer_registry.F90 index 07ca30dec8..16ee280355 100644 --- a/src/tracer/MOM_tracer_registry.F90 +++ b/src/tracer/MOM_tracer_registry.F90 @@ -426,12 +426,12 @@ subroutine register_tracer_diagnostics(Reg, h, Time, diag, G, GV, US, use_ALE) trim(flux_units), v_extensive = .true., x_cell_method = 'sum', & conversion=(US%L_to_m**2)*Tr%flux_scale*US%s_to_T) Tr%id_lbd_dfx = register_diag_field("ocean_model", trim(shortnm)//"_lbd_diffx", & - diag%axesCuL, Time, trim(flux_longname)//" diffusive zonal flux from the lateral boundary diffusion "& + diag%axesCuL, Time, trim(flux_longname)//" diffusive zonal flux from the lateral boundary diffusion "//& "scheme", trim(flux_units), v_extensive = .true., y_cell_method = 'sum', & conversion=(US%L_to_m**2)*Tr%flux_scale*US%s_to_T) Tr%id_lbd_dfy = register_diag_field("ocean_model", trim(shortnm)//"_lbd_diffy", & - diag%axesCvL, Time, trim(flux_longname)//" diffusive meridional flux from the lateral boundary diffusion"& - " scheme", trim(flux_units), v_extensive = .true., x_cell_method = 'sum', & + diag%axesCvL, Time, trim(flux_longname)//" diffusive meridional flux from the lateral boundary diffusion "//& + "scheme", trim(flux_units), v_extensive = .true., x_cell_method = 'sum', & conversion=(US%L_to_m**2)*Tr%flux_scale*US%s_to_T) else Tr%id_adx = register_diag_field("ocean_model", trim(shortnm)//"_adx", &