Merge pull request #55 from NOAA-GFDL/dev/gfdl

Merge in latest updates from dev/gfdl

.testing/Makefile

@@ -40,7 +40,7 @@ MKMF_TEMPLATE ?= build/mkmf/templates/ncrc-gnu.mk
 # Executables
 BUILDS = symmetric asymmetric repro openmp
 CONFIGS := $(wildcard tc*)
-TESTS = grids layouts restarts nans dims openmps
+TESTS = grids layouts restarts nans dims openmps rotations
 
 # REPRO tests enable reproducibility with optimization, and often do not match
 # the DEBUG results in older GCCs and vendor compilers, so we can optionally
@@ -186,9 +186,15 @@ test: $(foreach t,$(TESTS),test.$(t))
 
 # NOTE: We remove tc3 (OBC) from grid test since it cannot run asymmetric grids
 
+# NOTE: rotation diag chksum disabled since we cannot yet compare rotationally
+#       equivalent diagnostics
+
+# TODO: restart checksum comparison is not yet implemented
+
 .PHONY: $(foreach t,$(TESTS),test.$(t))
 test.grids: $(foreach c,$(filter-out tc3,$(CONFIGS)),$(c).grid $(c).grid.diag)
 test.layouts: $(foreach c,$(CONFIGS),$(c).layout $(c).layout.diag)
+test.rotations: $(foreach c,$(CONFIGS),$(c).rotate)
 test.restarts: $(foreach c,$(CONFIGS),$(c).restart)
 test.repros: $(foreach c,$(CONFIGS),$(c).repro $(c).repro.diag)
 test.openmps: $(foreach c,$(CONFIGS),$(c).openmp $(c).openmp.diag)
@@ -210,6 +216,7 @@ endef
 
 $(eval $(call CMP_RULE,grid,symmetric asymmetric))
 $(eval $(call CMP_RULE,layout,symmetric layout))
+$(eval $(call CMP_RULE,rotate,symmetric rotate))
 $(eval $(call CMP_RULE,repro,symmetric repro))
 $(eval $(call CMP_RULE,openmp,symmetric openmp))
 $(eval $(call CMP_RULE,nan,symmetric nan))
@@ -260,7 +267,7 @@ results/%/ocean.stats.$(1): build/$(2)/MOM6
 	cp -rL $$*/* work/$$*/$(1)
 	cd work/$$*/$(1) && if [ -f Makefile ]; then make; fi
 	mkdir -p work/$$*/$(1)/RESTART
-	echo $(4) > work/$$*/$(1)/MOM_override
+	echo -e "$(4)" > work/$$*/$(1)/MOM_override
 	cd work/$$*/$(1) && $$(call MPIRUN_CMD,$(5)) -n $(6) ../../../$$< 2> debug.out > std.out \
 		|| ! sed 's/^/$$*.$(1): /' std.out debug.out \
 		&& sed 's/^/$$*.$(1): /' std.out
@@ -282,6 +289,7 @@ $(eval $(call STAT_RULE,target,target,,,,1))
 $(eval $(call STAT_RULE,repro,repro,,,,1))
 $(eval $(call STAT_RULE,openmp,openmp,,,,1))
 $(eval $(call STAT_RULE,layout,symmetric,,LAYOUT=2$(,)1,,2))
+$(eval $(call STAT_RULE,rotate,symmetric,,ROTATE_INDEX=True\nINDEX_TURNS=1,,1))
 $(eval $(call STAT_RULE,nan,symmetric,,,MALLOC_PERTURB_=256,1))
 $(eval $(call STAT_RULE,dim.t,symmetric,,T_RESCALE_POWER=11,,1))
 $(eval $(call STAT_RULE,dim.l,symmetric,,L_RESCALE_POWER=11,,1))

.testing/tc2/MOM_input

@@ -600,3 +600,4 @@ ENERGYSAVEDAYS = 0.5            !   [days] default = 3600.0
                                 ! energies of the run and other globally summed diagnostics.
 DIAG_AS_CHKSUM = True
 DEBUG = True
+USE_GM_WORK_BUG = False

src/core/MOM_barotropic.F90

@@ -20,7 +20,8 @@ module MOM_barotropic
 use MOM_open_boundary, only : ocean_OBC_type, OBC_SIMPLE, OBC_NONE, open_boundary_query
 use MOM_open_boundary, only : OBC_DIRECTION_E, OBC_DIRECTION_W
 use MOM_open_boundary, only : OBC_DIRECTION_N, OBC_DIRECTION_S, OBC_segment_type
-use MOM_restart, only : register_restart_field, query_initialized, MOM_restart_CS
+use MOM_restart, only : register_restart_field, register_restart_pair
+use MOM_restart, only : query_initialized, MOM_restart_CS
 use MOM_tidal_forcing, only : tidal_forcing_sensitivity, tidal_forcing_CS
 use MOM_time_manager, only : time_type, real_to_time, operator(+), operator(-)
 use MOM_unit_scaling, only : unit_scale_type
@@ -1536,11 +1537,15 @@ subroutine btstep(U_in, V_in, eta_in, dt, bc_accel_u, bc_accel_v, forces, pbce,
     if (.not. use_BT_cont) then
       call uvchksum("BT Dat[uv]", Datu, Datv, CS%debug_BT_HI, haloshift=1, scale=US%L_to_m*GV%H_to_m)
     endif
-    call uvchksum("BT wt_[uv]", wt_u, wt_v, G%HI, 0, .true., .true.)
-    call uvchksum("BT frhat[uv]", CS%frhatu, CS%frhatv, G%HI, 0, .true., .true.)
+    call uvchksum("BT wt_[uv]", wt_u, wt_v, G%HI, haloshift=0, &
+                  symmetric=.true., omit_corners=.true., scalar_pair=.true.)
+    call uvchksum("BT frhat[uv]", CS%frhatu, CS%frhatv, G%HI, haloshift=0, &
+                  symmetric=.true., omit_corners=.true., scalar_pair=.true.)
     call uvchksum("BT bc_accel_[uv]", bc_accel_u, bc_accel_v, G%HI, haloshift=0, scale=US%L_T2_to_m_s2)
-    call uvchksum("BT IDat[uv]", CS%IDatu, CS%IDatv, G%HI, haloshift=0, scale=US%m_to_Z)
-    call uvchksum("BT visc_rem_[uv]", visc_rem_u, visc_rem_v, G%HI, haloshift=1)
+    call uvchksum("BT IDat[uv]", CS%IDatu, CS%IDatv, G%HI, haloshift=0, &
+                  scale=US%m_to_Z, scalar_pair=.true.)
+    call uvchksum("BT visc_rem_[uv]", visc_rem_u, visc_rem_v, G%HI, &
+                  haloshift=1, scalar_pair=.true.)
   endif
 
   if (query_averaging_enabled(CS%diag)) then
@@ -3113,9 +3118,13 @@ subroutine btcalc(h, G, GV, CS, h_u, h_v, may_use_default, OBC)
   enddo ; endif
 
   if (CS%debug) then
-    call uvchksum("btcalc frhat[uv]", CS%frhatu, CS%frhatv, G%HI, 0, .true., .true.)
+    call uvchksum("btcalc frhat[uv]", CS%frhatu, CS%frhatv, G%HI, &
+                  haloshift=0, symmetric=.true., omit_corners=.true., &
+                  scalar_pair=.true.)
     if (present(h_u) .and. present(h_v)) &
-      call uvchksum("btcalc h_[uv]", h_u, h_v, G%HI, 0, .true., .true., scale=GV%H_to_m)
+      call uvchksum("btcalc h_[uv]", h_u, h_v, G%HI, haloshift=0, &
+                    symmetric=.true., omit_corners=.true., scale=GV%H_to_m, &
+                    scalar_pair=.true.)
     call hchksum(h, "btcalc h",G%HI, haloshift=1, scale=GV%H_to_m)
   endif
 
@@ -4235,6 +4244,7 @@ subroutine barotropic_init(u, v, h, eta, Time, G, GV, US, param_file, diag, CS,
     CS%debug_BT_HI%IedB=CS%iedw
     CS%debug_BT_HI%JsdB=CS%jsdw-1
     CS%debug_BT_HI%JedB=CS%jedw
+    CS%debug_BT_HI%turns = G%HI%turns
   endif
 
   ! IareaT, IdxCu, and IdyCv need to be allocated with wide halos.
@@ -4607,6 +4617,7 @@ subroutine register_barotropic_restarts(HI, GV, param_file, CS, restart_CS)
   type(vardesc) :: vd(3)
   real :: slow_rate
   integer :: isd, ied, jsd, jed, IsdB, IedB, JsdB, JedB
+
   isd = HI%isd ; ied = HI%ied ; jsd = HI%jsd ; jed = HI%jed
   IsdB = HI%IsdB ; IedB = HI%IedB ; JsdB = HI%JsdB ; JedB = HI%JedB
 
@@ -4628,17 +4639,15 @@ subroutine register_barotropic_restarts(HI, GV, param_file, CS, restart_CS)
                 hor_grid='u', z_grid='1')
   vd(3) = var_desc("vbtav","m s-1","Time mean barotropic meridional velocity",&
                 hor_grid='v', z_grid='1')
-  call register_restart_field(CS%ubtav, vd(2), .false., restart_CS)
-  call register_restart_field(CS%vbtav, vd(3), .false., restart_CS)
+  call register_restart_pair(CS%ubtav, CS%vbtav, vd(2), vd(3), .false., restart_CS)
 
   vd(2) = var_desc("ubt_IC", "m s-1", &
               longname="Next initial condition for the barotropic zonal velocity", &
               hor_grid='u', z_grid='1')
   vd(3) = var_desc("vbt_IC", "m s-1", &
               longname="Next initial condition for the barotropic meridional velocity",&
               hor_grid='v', z_grid='1')
-  call register_restart_field(CS%ubt_IC, vd(2), .false., restart_CS)
-  call register_restart_field(CS%vbt_IC, vd(3), .false., restart_CS)
+  call register_restart_pair(CS%ubt_IC, CS%vbt_IC, vd(2), vd(3), .false., restart_CS)
 
   if (GV%Boussinesq) then
     vd(2) = var_desc("uhbt_IC", "m3 s-1", &
@@ -4655,8 +4664,7 @@ subroutine register_barotropic_restarts(HI, GV, param_file, CS, restart_CS)
                 longname="Next initial condition for the barotropic meridional transport",&
                 hor_grid='v', z_grid='1')
   endif
-  call register_restart_field(CS%uhbt_IC, vd(2), .false., restart_CS)
-  call register_restart_field(CS%vhbt_IC, vd(3), .false., restart_CS)
+  call register_restart_pair(CS%uhbt_IC, CS%vhbt_IC, vd(2), vd(3), .false., restart_CS)
 
   call register_restart_field(CS%dtbt, "DTBT", .false., restart_CS, &
                               longname="Barotropic timestep", units="seconds")

src/core/MOM_dynamics_split_RK2.F90

@@ -28,7 +28,8 @@ module MOM_dynamics_split_RK2
 use MOM_file_parser,       only : get_param, log_version, param_file_type
 use MOM_get_input,         only : directories
 use MOM_io,                only : MOM_io_init, vardesc, var_desc
-use MOM_restart,           only : register_restart_field, query_initialized, save_restart
+use MOM_restart,           only : register_restart_field, register_restart_pair
+use MOM_restart,           only : query_initialized, save_restart
 use MOM_restart,           only : restart_init, is_new_run, MOM_restart_CS
 use MOM_time_manager,      only : time_type, time_type_to_real, operator(+)
 use MOM_time_manager,      only : operator(-), operator(>), operator(*), operator(/)
@@ -886,11 +887,12 @@ subroutine register_restarts_dyn_split_RK2(HI, GV, param_file, CS, restart_CS, u
   real, dimension(SZI_(HI),SZJB_(HI),SZK_(GV)), &
                          target, intent(inout) :: vh !< merid volume/mass transport [H L2 T-1 ~> m3 s-1 or kg s-1]
 
-  type(vardesc)      :: vd
+  type(vardesc)      :: vd(2)
   character(len=40)  :: mdl = "MOM_dynamics_split_RK2" ! This module's name.
   character(len=48)  :: thickness_units, flux_units
 
   integer :: isd, ied, jsd, jed, nz, IsdB, IedB, JsdB, JedB
+
   isd  = HI%isd  ; ied  = HI%ied  ; jsd  = HI%jsd  ; jed  = HI%jed ; nz = GV%ke
   IsdB = HI%IsdB ; IedB = HI%IedB ; JsdB = HI%JsdB ; JedB = HI%JedB
 
@@ -918,32 +920,26 @@ subroutine register_restarts_dyn_split_RK2(HI, GV, param_file, CS, restart_CS, u
   flux_units = get_flux_units(GV)
 
   if (GV%Boussinesq) then
-    vd = var_desc("sfc",thickness_units,"Free surface Height",'h','1')
+    vd(1) = var_desc("sfc",thickness_units,"Free surface Height",'h','1')
   else
-    vd = var_desc("p_bot",thickness_units,"Bottom Pressure",'h','1')
+    vd(1) = var_desc("p_bot",thickness_units,"Bottom Pressure",'h','1')
   endif
-  call register_restart_field(CS%eta, vd, .false., restart_CS)
-
-  vd = var_desc("u2","m s-1","Auxiliary Zonal velocity",'u','L')
-  call register_restart_field(CS%u_av, vd, .false., restart_CS)
-
-  vd = var_desc("v2","m s-1","Auxiliary Meridional velocity",'v','L')
-  call register_restart_field(CS%v_av, vd, .false., restart_CS)
-
-  vd = var_desc("h2",thickness_units,"Auxiliary Layer Thickness",'h','L')
-  call register_restart_field(CS%h_av, vd, .false., restart_CS)
+  call register_restart_field(CS%eta, vd(1), .false., restart_CS)
 
-  vd = var_desc("uh",flux_units,"Zonal thickness flux",'u','L')
-  call register_restart_field(uh, vd, .false., restart_CS)
+  vd(1) = var_desc("u2","m s-1","Auxiliary Zonal velocity",'u','L')
+  vd(2) = var_desc("v2","m s-1","Auxiliary Meridional velocity",'v','L')
+  call register_restart_pair(CS%u_av, CS%v_av, vd(1), vd(2), .false., restart_CS)
 
-  vd = var_desc("vh",flux_units,"Meridional thickness flux",'v','L')
-  call register_restart_field(vh, vd, .false., restart_CS)
+  vd(1) = var_desc("h2",thickness_units,"Auxiliary Layer Thickness",'h','L')
+  call register_restart_field(CS%h_av, vd(1), .false., restart_CS)
 
-  vd = var_desc("diffu","m s-2","Zonal horizontal viscous acceleration",'u','L')
-  call register_restart_field(CS%diffu, vd, .false., restart_CS)
+  vd(1) = var_desc("uh",flux_units,"Zonal thickness flux",'u','L')
+  vd(2) = var_desc("vh",flux_units,"Meridional thickness flux",'v','L')
+  call register_restart_pair(uh, vh, vd(1), vd(2), .false., restart_CS)
 
-  vd = var_desc("diffv","m s-2","Meridional horizontal viscous acceleration",'v','L')
-  call register_restart_field(CS%diffv, vd, .false., restart_CS)
+  vd(1) = var_desc("diffu","m s-2","Zonal horizontal viscous acceleration",'u','L')
+  vd(2) = var_desc("diffv","m s-2","Meridional horizontal viscous acceleration",'v','L')
+  call register_restart_pair(CS%diffu, CS%diffv, vd(1), vd(2), .false., restart_CS)
 
   call register_barotropic_restarts(HI, GV, param_file, CS%barotropic_CSp, &
                                     restart_CS)