Test

.testing/Makefile

@@ -263,10 +263,12 @@ $(BUILD)/timing/Makefile: MOM_ACFLAGS += --with-driver=timing_tests
 
 
 # Build executables
+.NOTPARALLEL:$(foreach e,$(UNIT_EXECS),$(BUILD)/unit/$(e))
 $(BUILD)/unit/test_%: $(BUILD)/unit/Makefile FORCE
 	cd $(@D) && $(TIME) $(MAKE) $(@F) -j
 $(BUILD)/unit/Makefile: $(foreach e,$(UNIT_EXECS),../config_src/drivers/unit_tests/$(e).F90)
 
+.NOTPARALLEL:$(foreach e,$(TIMING_EXECS),$(BUILD)/timing/$(e))
 $(BUILD)/timing/time_%: $(BUILD)/timing/Makefile FORCE
 	cd $(@D) && $(TIME) $(MAKE) $(@F) -j
 $(BUILD)/timing/Makefile: $(foreach e,$(TIMING_EXECS),../config_src/drivers/timing_tests/$(e).F90)

config_src/drivers/timing_tests/time_MOM_remapping.F90

@@ -0,0 +1,100 @@
+program time_MOM_remapping
+
+! This file is part of MOM6. See LICENSE.md for the license.
+
+use MOM_remapping, only : remapping_CS
+use MOM_remapping, only : initialize_remapping
+use MOM_remapping, only : remapping_core_h
+
+implicit none
+
+type(remapping_CS) :: CS
+integer, parameter :: nk=75, nij=20*20, nits=10, nsamp=100, nschemes = 2
+character(len=10) :: scheme_labels(nschemes)
+real, dimension(nschemes) :: timings ! Time for nits of nij calls for each scheme [s]
+real, dimension(nschemes) :: tmean ! Mean time for a call [s]
+real, dimension(nschemes) :: tstd ! Standard deviation of time for a call [s]
+real, dimension(nschemes) :: tmin ! Shortest time for a call [s]
+real, dimension(nschemes) :: tmax ! Longest time for a call [s]
+real, dimension(:,:), allocatable :: u0, u1 ! Source/target values [arbitrary but same units as each other]
+real, dimension(:,:), allocatable :: h0, h1 ! Source target thicknesses [0..1]
+real :: start, finish ! Times [s]
+real :: h0sum, h1sum ! Totals of h0 and h1 [nondim]
+integer :: ij, k, isamp, iter, ischeme ! Indices and counters
+integer :: seed_size ! Number of integers used by seed
+integer, allocatable :: seed(:) ! Random number seed
+
+! Set seed for random numbers
+call random_seed(size=seed_size)
+allocate( seed(seed_Size) )
+seed(:) = 102030405
+call random_seed(put=seed)
+
+scheme_labels(1) = 'PCM'
+scheme_labels(2) = 'PLM'
+
+! Set up some test data (note: using k,i indexing rather than i,k)
+allocate( u0(nk,nij), h0(nk,nij), u1(nk,nij), h1(nk,nij) )
+call random_number(u0) ! In range 0-1
+call random_number(h0) ! In range 0-1
+call random_number(h1) ! In range 0-1
+do ij = 1, nij
+  h0(:,ij) = max(0., h0(:,ij) - 0.05) ! Make 5% of values equal to zero
+  h1(:,ij) = max(0., h1(:,ij) - 0.05) ! Make 5% of values equal to zero
+  h0sum = h0(1,ij)
+  h1sum = h1(1,ij)
+  do k = 2, nk
+    h0sum = h0sum + h0(k,ij)
+    h1sum = h1sum + h1(k,ij)
+  enddo
+  h0(:,ij) = h0(:,ij) / h0sum
+  h1(:,ij) = h1(:,ij) / h1sum
+enddo
+
+! Loop over many samples of timing loop to collect statistics
+tmean(:) = 0.
+tstd(:) = 0.
+tmin(:) = 1.e9
+tmax(:) = 0.
+do isamp = 1, nsamp
+  ! Time reconstruction + remapping
+  do ischeme = 1, nschemes
+    call initialize_remapping(CS, remapping_scheme=trim(scheme_labels(ischeme)))
+    call cpu_time(start)
+    do iter = 1, nits ! Make many passes to reduce sampling error
+      do ij = 1, nij ! Calling nij times to make similar to cost in MOM_ALE()
+        call remapping_core_h(CS, nk, h0(:,ij), u0(:,ij), nk, h1(:,ij), u1(:,ij))
+      enddo
+    enddo
+    call cpu_time(finish)
+    timings(ischeme) = (finish-start)/real(nits*nij) ! Average time per call
+  enddo
+  tmean(:) = tmean(:) + timings(:)
+  tstd(:) = tstd(:) + timings(:)**2 ! tstd contains sum of squares here
+  tmin(:) = min( tmin(:), timings(:) )
+  tmax(:) = max( tmax(:), timings(:) )
+enddo
+tmean(:) = tmean(:) / real(nsamp) ! convert to mean
+tstd(:) = tstd(:) / real(nsamp) ! convert to mean of squares
+tstd(:) = tstd(:) - tmean(:)**2  ! convert to variance
+tstd(:) = sqrt( tstd(:) * real(nsamp) / real(nsamp-1) ) ! convert to standard deviation
+
+
+! Display results in YAML
+write(*,'(a)') "{"
+do ischeme = 1, nschemes
+  write(*,"(2x,5a)") '"MOM_remapping remapping_core_h(remapping_scheme=', &
+                     trim(scheme_labels(ischeme)), ')": {'
+  write(*,"(4x,a,1pe11.4,',')") '"min": ',tmin(ischeme)
+  write(*,"(4x,a,1pe11.4,',')") '"mean":',tmean(ischeme)
+  write(*,"(4x,a,1pe11.4,',')") '"std": ',tstd(ischeme)
+  write(*,"(4x,a,i7,',')") '"n_samples": ',nsamp
+  if (ischeme.ne.nschemes) then
+    write(*,"(4x,a,1pe11.4,'},')") '"max": ',tmax(ischeme)
+  else
+    write(*,"(4x,a,1pe11.4,'}')") '"max": ',tmax(ischeme)
+  endif
+enddo
+write(*,'(a)') "}"
+
+end program time_MOM_remapping

config_src/drivers/unit_tests/test_MOM_remapping.F90

@@ -0,0 +1,7 @@
+program test_MOM_remapping
+
+use MOM_remapping, only : remapping_unit_tests
+
+if (remapping_unit_tests(.true.)) stop 1
+
+end program test_MOM_remapping

src/ALE/MOM_ALE.F90

@@ -178,6 +178,7 @@ subroutine ALE_init( param_file, GV, US, max_depth, CS)
   logical           :: local_logical
   logical           :: remap_boundary_extrap
   logical           :: init_boundary_extrap
+  logical           :: om4_remap_via_sub_cells
   type(hybgen_regrid_CS), pointer :: hybgen_regridCS => NULL() ! Control structure for hybgen regridding
                                                          ! for sharing parameters.
 
@@ -234,6 +235,11 @@ subroutine ALE_init( param_file, GV, US, max_depth, CS)
   call get_param(param_file, mdl, "DEFAULT_ANSWER_DATE", default_answer_date, &
                  "This sets the default value for the various _ANSWER_DATE parameters.", &
                  default=99991231)
+  call get_param(param_file, mdl, "REMAPPING_USE_OM4_SUBCELLS", om4_remap_via_sub_cells, &
+                 "This selects the remapping algorithm used in OM4 that does not use "//&
+                 "the full reconstruction for the top- and lower-most sub-layers, but instead "//&
+                 "assumes they are always vanished (untrue) and so just uses their edge values. "//&
+                 "We recommend setting this option to false.", default=.true.)
   call get_param(param_file, mdl, "REMAPPING_ANSWER_DATE", CS%answer_date, &
                  "The vintage of the expressions and order of arithmetic to use for remapping.  "//&
                  "Values below 20190101 result in the use of older, less accurate expressions "//&
@@ -247,12 +253,14 @@ subroutine ALE_init( param_file, GV, US, max_depth, CS)
                              check_reconstruction=check_reconstruction, &
                              check_remapping=check_remapping, &
                              force_bounds_in_subcell=force_bounds_in_subcell, &
+                             om4_remap_via_sub_cells=om4_remap_via_sub_cells, &
                              answer_date=CS%answer_date)
   call initialize_remapping( CS%vel_remapCS, vel_string, &
                              boundary_extrapolation=init_boundary_extrap, &
                              check_reconstruction=check_reconstruction, &
                              check_remapping=check_remapping, &
                              force_bounds_in_subcell=force_bounds_in_subcell, &
+                             om4_remap_via_sub_cells=om4_remap_via_sub_cells, &
                              answer_date=CS%answer_date)
 
   call get_param(param_file, mdl, "PARTIAL_CELL_VELOCITY_REMAP", CS%partial_cell_vel_remap, &
@@ -326,6 +334,21 @@ subroutine ALE_set_extrap_boundaries( param_file, CS)
   call remapping_set_param(CS%remapCS, boundary_extrapolation=remap_boundary_extrap)
 end subroutine ALE_set_extrap_boundaries
 
+!> Sets the remapping algorithm to that of OM4
+!!
+!! The remapping aglorithm used in OM4 made poor assumptions about the reconstructions
+!! in the top/bottom layers, namely that they were always vanished and could be
+!! represented solely by their upper/lower edge value respectively.
+!! Passing .false. here uses the full reconstruction of those top and bottom layers
+!! and properly sample those layers.
+subroutine ALE_set_OM4_remap_algorithm( CS, om4_remap_via_sub_cells )
+  type(ALE_CS), pointer :: CS !< Module control structure
+  logical, intent(in)   :: om4_remap_via_sub_cells !< If true, use OM4 remapping algorithm
+
+  call remapping_set_param(CS%remapCS, om4_remap_via_sub_cells =om4_remap_via_sub_cells )
+
+end subroutine ALE_set_OM4_remap_algorithm
+
 !> Initialize diagnostics for the ALE module.
 subroutine ALE_register_diags(Time, G, GV, US, diag, CS)
   type(time_type),target,     intent(in)  :: Time  !< Time structure