+Make output reproduce across layout if DEBUG=True

This commit makes a set of changes so that the chksums output by the OM4_05 configuration with DEBUG = True reproduce between 252 and 256 PE runs (note that the solutions themselves already reproduced). This is includes adding a new optional omit_corners argument to the MOM_state_chksum and MOM_thermo_chksum routines, and revising the haloshift, omit_corners or symmetric arguments on 6 subroutine calls, to check more appropriate loop ranges. All answers are bitwise identical, but there is a new optional argument to three publicly visible debugging routines.
NOAA-GFDL · Feb 3, 2023 · cfedaa3 · cfedaa3
1 parent cab1e84
commit cfedaa3
Show file tree

Hide file tree

Showing 4 changed files with 39 additions and 27 deletions.
diff --git a/src/core/MOM.F90 b/src/core/MOM.F90
@@ -1529,9 +1529,9 @@ subroutine step_MOM_thermo(CS, G, GV, US, u, v, h, tv, fluxes, dtdia, &
       call preAle_tracer_diagnostics(CS%tracer_Reg, G, GV)
 
       if (CS%debug) then
-        call MOM_state_chksum("Pre-ALE ", u, v, h, CS%uh, CS%vh, G, GV, US)
-        call hchksum(tv%T,"Pre-ALE T", G%HI, haloshift=1, scale=US%C_to_degC)
-        call hchksum(tv%S,"Pre-ALE S", G%HI, haloshift=1, scale=US%S_to_ppt)
+        call MOM_state_chksum("Pre-ALE ", u, v, h, CS%uh, CS%vh, G, GV, US, omit_corners=.true.)
+        call hchksum(tv%T,"Pre-ALE T", G%HI, haloshift=1, omit_corners=.true., scale=US%C_to_degC)
+        call hchksum(tv%S,"Pre-ALE S", G%HI, haloshift=1, omit_corners=.true., scale=US%S_to_ppt)
         call check_redundant("Pre-ALE ", u, v, G, unscale=US%L_T_to_m_s)
       endif
       call cpu_clock_begin(id_clock_ALE)

diff --git a/src/core/MOM_checksum_packages.F90 b/src/core/MOM_checksum_packages.F90
@@ -39,7 +39,7 @@ module MOM_checksum_packages
 ! =============================================================================
 
 !> Write out chksums for the model's basic state variables, including transports.
-subroutine MOM_state_chksum_5arg(mesg, u, v, h, uh, vh, G, GV, US, haloshift, symmetric, vel_scale)
+subroutine MOM_state_chksum_5arg(mesg, u, v, h, uh, vh, G, GV, US, haloshift, symmetric, omit_corners, vel_scale)
   character(len=*),                          &
                            intent(in) :: mesg !< A message that appears on the chksum lines.
   type(ocean_grid_type),   intent(in) :: G    !< The ocean's grid structure.
@@ -60,6 +60,7 @@ subroutine MOM_state_chksum_5arg(mesg, u, v, h, uh, vh, G, GV, US, haloshift, sy
   integer,       optional, intent(in) :: haloshift !< The width of halos to check (default 0).
   logical,       optional, intent(in) :: symmetric !< If true, do checksums on the fully symmetric
                                                    !! computational domain.
+  logical,       optional, intent(in) :: omit_corners !< If true, avoid checking diagonal shifts
   real,          optional, intent(in) :: vel_scale !< The scaling factor to convert velocities to [m s-1]
 
   real :: scale_vel ! The scaling factor to convert velocities to [m s-1]
@@ -73,16 +74,17 @@ subroutine MOM_state_chksum_5arg(mesg, u, v, h, uh, vh, G, GV, US, haloshift, sy
   sym = .false. ; if (present(symmetric)) sym=symmetric
   scale_vel = US%L_T_to_m_s ; if (present(vel_scale)) scale_vel = vel_scale
 
-  call uvchksum(mesg//" [uv]", u, v, G%HI, haloshift=hs, symmetric=sym, scale=scale_vel)
-  call hchksum(h, mesg//" h", G%HI, haloshift=hs, scale=GV%H_to_m)
-  call uvchksum(mesg//" [uv]h", uh, vh, G%HI, haloshift=hs, &
-                symmetric=sym, scale=GV%H_to_m*US%L_to_m**2*US%s_to_T)
+  call uvchksum(mesg//" [uv]", u, v, G%HI, haloshift=hs, symmetric=sym, &
+                omit_corners=omit_corners, scale=scale_vel)
+  call hchksum(h, mesg//" h", G%HI, haloshift=hs, omit_corners=omit_corners, scale=GV%H_to_m)
+  call uvchksum(mesg//" [uv]h", uh, vh, G%HI, haloshift=hs, symmetric=sym, &
+                omit_corners=omit_corners, scale=GV%H_to_m*US%L_to_m**2*US%s_to_T)
 end subroutine MOM_state_chksum_5arg
 
 ! =============================================================================
 
 !> Write out chksums for the model's basic state variables.
-subroutine MOM_state_chksum_3arg(mesg, u, v, h, G, GV, US, haloshift, symmetric)
+subroutine MOM_state_chksum_3arg(mesg, u, v, h, G, GV, US, haloshift, symmetric, omit_corners)
   character(len=*),                intent(in) :: mesg !< A message that appears on the chksum lines.
   type(ocean_grid_type),           intent(in) :: G  !< The ocean's grid structure.
   type(verticalGrid_type),         intent(in) :: GV !< The ocean's vertical grid structure.
@@ -97,6 +99,7 @@ subroutine MOM_state_chksum_3arg(mesg, u, v, h, G, GV, US, haloshift, symmetric)
   integer,               optional, intent(in) :: haloshift !< The width of halos to check (default 0).
   logical,               optional, intent(in) :: symmetric !< If true, do checksums on the fully
                                                     !! symmetric computational domain.
+  logical,               optional, intent(in) :: omit_corners !< If true, avoid checking diagonal shifts
 
   integer :: hs
   logical :: sym
@@ -106,34 +109,43 @@ subroutine MOM_state_chksum_3arg(mesg, u, v, h, G, GV, US, haloshift, symmetric)
   ! and js...je as their extent.
   hs = 1 ; if (present(haloshift)) hs = haloshift
   sym = .false. ; if (present(symmetric)) sym = symmetric
-  call uvchksum(mesg//" u", u, v, G%HI, haloshift=hs, symmetric=sym, scale=US%L_T_to_m_s)
-  call hchksum(h, mesg//" h",G%HI, haloshift=hs, scale=GV%H_to_m)
+  call uvchksum(mesg//" u", u, v, G%HI, haloshift=hs, symmetric=sym, &
+                omit_corners=omit_corners, scale=US%L_T_to_m_s)
+  call hchksum(h, mesg//" h",G%HI, haloshift=hs, omit_corners=omit_corners, scale=GV%H_to_m)
 end subroutine MOM_state_chksum_3arg
 
 ! =============================================================================
 
 !> Write out chksums for the model's thermodynamic state variables.
-subroutine MOM_thermo_chksum(mesg, tv, G, US, haloshift)
+subroutine MOM_thermo_chksum(mesg, tv, G, US, haloshift, omit_corners)
   character(len=*),         intent(in) :: mesg !< A message that appears on the chksum lines.
   type(thermo_var_ptrs),    intent(in) :: tv   !< A structure pointing to various
                                                !! thermodynamic variables.
   type(ocean_grid_type),    intent(in) :: G    !< The ocean's grid structure.
   type(unit_scale_type),    intent(in) :: US   !< A dimensional unit scaling type
   integer,        optional, intent(in) :: haloshift !< The width of halos to check (default 0).
+  logical,        optional, intent(in) :: omit_corners !< If true, avoid checking diagonal shifts
 
   integer :: hs
   hs=1 ; if (present(haloshift)) hs=haloshift
 
-  if (associated(tv%T)) call hchksum(tv%T, mesg//" T", G%HI, haloshift=hs, scale=US%C_to_degC)
-  if (associated(tv%S)) call hchksum(tv%S, mesg//" S", G%HI, haloshift=hs, scale=US%S_to_ppt)
-  if (associated(tv%frazil)) call hchksum(tv%frazil, mesg//" frazil", G%HI, haloshift=hs, &
-                                          scale=US%Q_to_J_kg*US%R_to_kg_m3*US%Z_to_m)
-  if (associated(tv%salt_deficit)) call hchksum(tv%salt_deficit, mesg//" salt deficit", G%HI, haloshift=hs, &
-                                                scale=US%S_to_ppt*US%RZ_to_kg_m2)
-  if (associated(tv%varT)) call hchksum(tv%varT, mesg//" varT", G%HI, haloshift=hs, scale=US%C_to_degC**2)
-  if (associated(tv%varS)) call hchksum(tv%varS, mesg//" varS", G%HI, haloshift=hs, scale=US%S_to_ppt**2)
-  if (associated(tv%covarTS)) call hchksum(tv%covarTS, mesg//" covarTS", G%HI, haloshift=hs, &
-                                           scale=US%S_to_ppt*US%C_to_degC)
+  if (associated(tv%T)) &
+    call hchksum(tv%T, mesg//" T", G%HI, haloshift=hs, omit_corners=omit_corners, scale=US%C_to_degC)
+  if (associated(tv%S)) &
+    call hchksum(tv%S, mesg//" S", G%HI, haloshift=hs, omit_corners=omit_corners, scale=US%S_to_ppt)
+  if (associated(tv%frazil)) &
+    call hchksum(tv%frazil, mesg//" frazil", G%HI, haloshift=hs, omit_corners=omit_corners, &
+                 scale=US%Q_to_J_kg*US%R_to_kg_m3*US%Z_to_m)
+  if (associated(tv%salt_deficit)) &
+    call hchksum(tv%salt_deficit, mesg//" salt deficit", G%HI, haloshift=hs, omit_corners=omit_corners, &
+                 scale=US%S_to_ppt*US%RZ_to_kg_m2)
+  if (associated(tv%varT)) &
+    call hchksum(tv%varT, mesg//" varT", G%HI, haloshift=hs, omit_corners=omit_corners, scale=US%C_to_degC**2)
+  if (associated(tv%varS)) &
+    call hchksum(tv%varS, mesg//" varS", G%HI, haloshift=hs, omit_corners=omit_corners, scale=US%S_to_ppt**2)
+  if (associated(tv%covarTS)) &
+    call hchksum(tv%covarTS, mesg//" covarTS", G%HI, haloshift=hs, omit_corners=omit_corners, &
+                 scale=US%S_to_ppt*US%C_to_degC)
 
 end subroutine MOM_thermo_chksum
 

diff --git a/src/core/MOM_dynamics_split_RK2.F90 b/src/core/MOM_dynamics_split_RK2.F90
@@ -647,7 +647,7 @@ subroutine step_MOM_dyn_split_RK2(u, v, h, tv, visc, Time_local, dt, forces, p_s
 !   call MOM_state_chksum("Predictor 1", up, vp, h, uh, vh, G, GV, US, haloshift=1)
     call MOM_accel_chksum("Predictor accel", CS%CAu_pred, CS%CAv_pred, CS%PFu, CS%PFv, &
              CS%diffu, CS%diffv, G, GV, US, CS%pbce, CS%u_accel_bt, CS%v_accel_bt, symmetric=sym)
-    call MOM_state_chksum("Predictor 1 init", u, v, h, uh, vh, G, GV, US, haloshift=2, &
+    call MOM_state_chksum("Predictor 1 init", u, v, h, uh, vh, G, GV, US, haloshift=1, &
                           symmetric=sym)
     call check_redundant("Predictor 1 up", up, vp, G, unscale=US%L_T_to_m_s)
     call check_redundant("Predictor 1 uh", uh, vh, G, unscale=GV%H_to_m*US%L_to_m**2*US%s_to_T)
@@ -867,8 +867,8 @@ subroutine step_MOM_dyn_split_RK2(u, v, h, tv, visc, Time_local, dt, forces, p_s
   call cpu_clock_end(id_clock_mom_update)
 
   if (CS%debug) then
-    call uvchksum("Corrector 1 [uv]", u, v, G%HI,haloshift=0, symmetric=sym, scale=US%L_T_to_m_s)
-    call hchksum(h, "Corrector 1 h", G%HI, haloshift=2, scale=GV%H_to_m)
+    call uvchksum("Corrector 1 [uv]", u, v, G%HI, haloshift=0, symmetric=sym, scale=US%L_T_to_m_s)
+    call hchksum(h, "Corrector 1 h", G%HI, haloshift=1, scale=GV%H_to_m)
     call uvchksum("Corrector 1 [uv]h", uh, vh, G%HI, haloshift=2, &
                   symmetric=sym, scale=GV%H_to_m*US%L_to_m**2*US%s_to_T)
   ! call MOM_state_chksum("Corrector 1", u, v, h, uh, vh, G, GV, US, haloshift=1)

diff --git a/src/parameterizations/lateral/MOM_MEKE.F90 b/src/parameterizations/lateral/MOM_MEKE.F90
@@ -260,7 +260,7 @@ subroutine step_forward_MEKE(MEKE, h, SN_u, SN_v, visc, dt, G, GV, US, CS, hu, h
         call hchksum(MEKE%MEKE, 'MEKE MEKE', G%HI, scale=US%L_T_to_m_s**2)
       call uvchksum("MEKE SN_[uv]", SN_u, SN_v, G%HI, scale=US%s_to_T, &
                     scalar_pair=.true.)
-      call uvchksum("MEKE h[uv]", hu, hv, G%HI, haloshift=1, &
+      call uvchksum("MEKE h[uv]", hu, hv, G%HI, haloshift=0, symmetric=.true., &
                     scale=GV%H_to_m*(US%L_to_m**2))
     endif
 
@@ -293,7 +293,7 @@ subroutine step_forward_MEKE(MEKE, h, SN_u, SN_v, visc, dt, G, GV, US, CS, hu, h
         enddo ; enddo
       enddo
       if (CS%MEKE_advection_bug) then
-        ! This code obviously incorrect code reproduces a bug in the original implementation of
+        ! This obviously incorrect code reproduces a bug in the original implementation of
         ! the MEKE advection.
         do j=js,je ; do I=is-1,ie
           baroHu(I,j) = hu(I,j,nz) * GV%H_to_RZ