Add generalized axis permutation tests for diag_manager

Author: abrooks1085

test_fms/diag_manager/Makefile.am

@@ -34,7 +34,7 @@ check_PROGRAMS = test_diag_manager test_diag_manager_time \
  check_time_pow check_time_rms check_subregional test_cell_measures test_var_masks \
  check_var_masks test_multiple_send_data test_diag_out_yaml test_output_every_freq \
  test_dm_weights test_prepend_date test_ens_runs test_diag_multi_file test_diag_attribute_add \
- check_new_file_freq test_zbounds_limits test_multiple_zbounds test_generalized_indicies check_generalized_indices
+ check_new_file_freq test_zbounds_limits test_multiple_zbounds test_generalized_indices check_generalized_indices
 
 # This is the source code for the test.
 test_output_every_freq_SOURCES = test_output_every_freq.F90
@@ -71,7 +71,7 @@ test_diag_attribute_add_SOURCES = test_diag_attribute_add.F90
 check_new_file_freq_SOURCES = check_new_file_freq.F90
 test_zbounds_limits_SOURCES = test_zbounds_limits.F90
 test_multiple_zbounds_SOURCES = test_multiple_zbounds.F90
-test_generalized_indicies_SOURCES = testing_utils.F90 test_generalized_indicies.F90
+test_generalized_indices_SOURCES = testing_utils.F90 test_generalized_indices.F90
 check_generalized_indices_SOURCES = testing_utils.F90 check_generalized_indices.F90
 
 TEST_EXTENSIONS = .sh

test_fms/diag_manager/check_generalized_indices.F90

@@ -16,26 +16,27 @@
 !* governing permissions and limitations under the License.
 !***********************************************************************
 
-!> @brief Checker for test_generalized_indicies output.
-!!        Verifies swapped-axis variables match identity variables under transpose:
-!!          var2_id(x,y)   == var2_swap(y,x)
-!!          var3_id(x,y,z) == var3_swap(y,x,z)
+!> @brief Checker for test_generalized_indices output.
+!!        Verifies permuted-axis variables match identity variables under axis permutations
 program check_generalized_indices
   use fms_mod,           only: fms_init, fms_end, string
+  use testing_utils,     only: check_perm
   use fms2_io_mod,       only: FmsNetcdfFile_t, read_data, open_file, close_file, get_global_attribute
   use mpp_mod,           only: mpp_error, FATAL, mpp_pe
-  use platform_mod,      only: r4_kind
+  use platform_mod,      only: r8_kind
 
   implicit none
 
   type(FmsNetcdfFile_t)           :: fileobj
   integer                         :: nx, ny, nz
   integer                         :: i
 
-  real(kind=r4_kind), allocatable :: var2_id(:,:)     ! (x,y)
-  real(kind=r4_kind), allocatable :: var2_swap(:,:)   ! (y,x)
-  real(kind=r4_kind), allocatable :: var3_id(:,:,:)   ! (x,y,z)
-  real(kind=r4_kind), allocatable :: var3_swap(:,:,:) ! (y,x,z)
+  real(kind=r8_kind), allocatable :: var2_id(:,:)     ! (x,y)
+  real(kind=r8_kind), allocatable :: var2_yx(:,:)     ! (y,x)
+  real(kind=r8_kind), allocatable :: var3_id(:,:,:)   ! (x,y,z)
+  real(kind=r8_kind), allocatable :: var3_zx(:,:,:)   ! (z,y,x)
+  real(kind=r8_kind), allocatable :: var3_yzx(:,:,:)  ! (y,z,x)
+  real(kind=r8_kind), allocatable :: var3_zxy(:,:,:)  ! (z,x,y)
 
   call fms_init()
 
@@ -48,25 +49,35 @@ program check_generalized_indices
 
   call check_global_attribute(fileobj, "test_generalized_indices")
 
-  allocate(var2_id(nx,ny), var2_swap(ny,nx))
-  allocate(var3_id(nx,ny,nz), var3_swap(ny,nx,nz))
+  allocate(var2_id(nx,ny),    var2_yx(ny,nx))
+  allocate(var3_id(nx,ny,nz), var3_zx(nz,ny,nx), var3_yzx(ny,nz,nx), var3_zxy(nz,nx,ny))
 
   ! Output every 6 hours over 48 hours => 8 records
   do i = 1, 8
-    var2_id   = -999._r4_kind
-    var2_swap = -999._r4_kind
-    var3_id   = -999._r4_kind
-    var3_swap = -999._r4_kind
-
-    print *, "Checking var2_swap vs var2_id - time_level:", string(i)
-    call read_data(fileobj, "var2_id",   var2_id,   unlim_dim_level=i)
-    call read_data(fileobj, "var2_swap", var2_swap, unlim_dim_level=i)
-    call check_var2_relation(var2_id, var2_swap)
-
-    print *, "Checking var3_swap vs var3_id - time_level:", string(i)
-    call read_data(fileobj, "var3_id",   var3_id,   unlim_dim_level=i)
-    call read_data(fileobj, "var3_swap", var3_swap, unlim_dim_level=i)
-    call check_var3_relation(var3_id, var3_swap)
+    var2_id  = -999._r8_kind
+    var2_yx  = -999._r8_kind
+    var3_id  = -999._r8_kind
+    var3_zx  = -999._r8_kind
+    var3_yzx = -999._r8_kind
+    var3_zxy = -999._r8_kind
+
+    print *, "Checking var2_yx vs var2_id - time_level:", i
+    call read_data(fileobj, "var2_id", var2_id, unlim_dim_level=i)
+    call read_data(fileobj, "var2_yx", var2_yx, unlim_dim_level=i)
+    call check_perm(var2_id, var2_yx, [2,1])
+
+    print *, "Checking var3_zx vs var3_id - time_level:", i
+    call read_data(fileobj, "var3_id", var3_id, unlim_dim_level=i)
+    call read_data(fileobj, "var3_zx", var3_zx, unlim_dim_level=i)
+    call check_perm(var3_id, var3_zx, [3,2,1])
+
+    print *, "Checking var3_yzx vs var3_id - time_level:", i
+    call read_data(fileobj, "var3_yzx", var3_yzx, unlim_dim_level=i)
+    call check_perm(var3_id, var3_yzx, [2,3,1])
+
+    print *, "Checking var3_zxy vs var3_id - time_level:", i
+    call read_data(fileobj, "var3_zxy", var3_zxy, unlim_dim_level=i)
+    call check_perm(var3_id, var3_zxy, [3,1,2])
   enddo
 
   call close_file(fileobj)
@@ -85,49 +96,4 @@ subroutine check_global_attribute(fileobj, expected_title)
       call mpp_error(FATAL, "Global attribute 'title' not expected value.")
     endif
   end subroutine check_global_attribute
-
-  subroutine check_var2_relation(v_id, v_sw)
-    real(kind=r4_kind), intent(in) :: v_id(:,:)  ! (x,y)
-    real(kind=r4_kind), intent(in) :: v_sw(:,:)  ! (y,x)
-
-    integer :: x, y
-
-    if (size(v_id,1) /= size(v_sw,2) .or. size(v_id,2) /= size(v_sw,1)) then
-      call mpp_error(FATAL, "check_var2_relation: dimension mismatch between var2_id and var2_swap")
-    endif
-
-    do x = 1, size(v_id,1)
-      do y = 1, size(v_id,2)
-        if (abs(v_id(x,y) - v_sw(y,x)) > 0) then
-          print *, mpp_pe(), "var2 mismatch at (x,y)=", x, y, " id=", v_id(x,y), " swap(y,x)=", v_sw(y,x)
-          call mpp_error(FATAL, "check_var2_relation: var2_swap != transpose(var2_id)")
-        endif
-      enddo
-    enddo
-  end subroutine check_var2_relation
-
-  subroutine check_var3_relation(v_id, v_sw)
-    real(kind=r4_kind), intent(in) :: v_id(:,:,:)  ! (x,y,z)
-    real(kind=r4_kind), intent(in) :: v_sw(:,:,:)  ! (y,x,z)
-
-    integer :: x, y, z
-
-    if (size(v_id,1) /= size(v_sw,2) .or. size(v_id,2) /= size(v_sw,1) .or. size(v_id,3) /= size(v_sw,3)) then
-      call mpp_error(FATAL, "check_var3_relation: dimension mismatch between var3_id and var3_swap")
-    endif
-
-    do x = 1, size(v_id,1)
-      do y = 1, size(v_id,2)
-        do z = 1, size(v_id,3)
-          if (abs(v_id(x,y,z) - v_sw(y,x,z)) > 0) then
-            print *, mpp_pe(), "var3 mismatch at (x,y,z)=", x, y, z, &
-                               " id=", v_id(x,y,z), " swap(y,x,z)=", v_sw(y,x,z)
-            call mpp_error(FATAL, "check_var3_relation: var3_swap != var3_id with x/y swapped")
-          endif
-        enddo
-      enddo
-    enddo
-  end subroutine check_var3_relation
-
 end program check_generalized_indices
-

…ag_manager/test_generalized_indicies.F90 …iag_manager/test_generalized_indices.F90test_fms/diag_manager/test_generalized_indicies.F90 renamed to test_fms/diag_manager/test_generalized_indices.F90 test_fms/diag_manager/test_generalized_indicies.F90 renamed to test_fms/diag_manager/test_generalized_indices.F90

@@ -16,11 +16,13 @@
 !* governing permissions and limitations under the License.
 !***********************************************************************
 
-!> @brief  Test generalized axis permutations (x/y for now) for send_data
+!> @brief  Test generalized axis permutations for send_data.
+!!         Applies predefined permutations to canonical (x,y,z,w) storage
+!!         and verifies consistency between data layout and axis metadata.
 !!         Assumes default configuration parameters: test_normal + no_mask.
 program test_generalized_indices
   use fms_mod,           only: fms_init, fms_end
-  use testing_utils,     only: allocate_buffer
+  use testing_utils,     only: allocate_buffer, permute
   use platform_mod,      only: r8_kind
   use mpp_mod,           only: FATAL, mpp_error, mpp_npes, mpp_pe, mpp_root_pe
   use time_manager_mod,  only: time_type, set_calendar_type, set_date, JULIAN, set_time, OPERATOR(+)
@@ -47,9 +49,21 @@ program test_generalized_indices
   real(r8_kind), allocatable :: cdata(:,:,:,:)   ! canonical storage: (x,y,z,w)
 
   ! Permutation test
-  integer :: p_id(3), p_swap(3)
-  integer :: id_var2_id, id_var2_swap
-  integer :: id_var3_id, id_var3_swap
+  integer, parameter :: LAYOUT_XY  = 1
+  integer, parameter :: LAYOUT_YX  = 2
+  integer, parameter :: LAYOUT_ZX  = 3
+  integer, parameter :: LAYOUT_YZX = 4
+  integer, parameter :: LAYOUT_ZXY = 5
+  integer, parameter :: PERM_TABLE(3,5) = reshape([ &
+                           1,2,3, & ! XY (identity)
+                           2,1,3, & ! YX
+                           3,2,1, & ! ZX
+                           2,3,1, & ! YZX
+                           3,1,2  & ! ZXY
+                        ], [3,5])
+
+  integer :: id_var2_id, id_var2_yx
+  integer :: id_var3_id, id_var3_zx, id_var3_yzx, id_var3_zxy
 
   call fms_init
   call set_calendar_type(JULIAN)
@@ -88,29 +102,26 @@ program test_generalized_indices
 
   missing_value = -666._r8_kind
 
-  ! Define permutations: identity (x,y,z) and swap (y,x,z)
-  p_id   = [1,2,3]
-  p_swap = [2,1,3]
-
   ! Register permuted diagnostic fields ONCE
-  id_var2_id   = register_diag_field('ocn_mod', 'var2_id',   (/axis(p_id(1)),   axis(p_id(2))/),   Time, 'Var2d id', &
-                                     'mullions', missing_value=missing_value)
-  id_var2_swap = register_diag_field('ocn_mod', 'var2_swap', (/axis(p_swap(1)), axis(p_swap(2))/), Time, 'Var2d swap', &
-                                     'mullions', missing_value=missing_value)
+  id_var2_id  = register_perm_diag_field('var2_id',  'Var2d id',  axis(1:2), LAYOUT_XY)
+  id_var3_id  = register_perm_diag_field('var3_id',  'Var3d id',  axis(1:3), LAYOUT_XY)
+
+  id_var2_yx  = register_perm_diag_field('var2_yx',  'Var2d yx',  axis(1:2), LAYOUT_YX)
+  id_var3_zx  = register_perm_diag_field('var3_zx',  'Var3d zx',  axis(1:3), LAYOUT_ZX)
 
-  id_var3_id   = register_diag_field('ocn_mod', 'var3_id',   (/axis(p_id(1)),   axis(p_id(2)),   axis(p_id(3))/), &
-                                      Time, 'Var3d id', 'mullions', missing_value=missing_value)
-  id_var3_swap = register_diag_field('ocn_mod', 'var3_swap', (/axis(p_swap(1)), axis(p_swap(2)), axis(p_swap(3))/), &
-                                      Time, 'Var3d swap', 'mullions', missing_value=missing_value)
+  id_var3_yzx = register_perm_diag_field('var3_yzx', 'Var3d yzx', axis(1:3), LAYOUT_YZX)
+  id_var3_zxy = register_perm_diag_field('var3_zxy', 'Var3d zxy', axis(1:3), LAYOUT_ZXY)
 
   if (mpp_pe() == mpp_root_pe()) then
     print *, "Testing generalized indices in default mode (test_normal + no_mask)"
     print *, "  canonical storage is (x,y,z,w)"
     print *, "  sending:"
-    print *, "    var2_id   with axes (x,y)"
-    print *, "    var2_swap with axes (y,x)"
-    print *, "    var3_id   with axes (x,y,z)"
-    print *, "    var3_swap with axes (y,x,z)"
+    print *, "    var2_id with axes (x,y)"
+    print *, "    var2_yx with axes (y,x)"
+    print *, "    var3_id with axes (x,y,z)"
+    print *, "    var3_zx with axes (z,y,x)"
+    print *, "    var3_yzx with axes (y,z,x)"
+    print *, "    var3_zxy with axes (z,x,y)"
   end if
 
   call diag_manager_set_time_end(set_date(2,1,3,0,0,0))
@@ -120,12 +131,16 @@ program test_generalized_indices
     call set_buffer(cdata, i)
 
     ! Identity: axes (x,y) / (x,y,z) with canonical storage
-    call send_var2_perm(id_var2_id,   cdata, p_id,   Time)
-    call send_var3_perm(id_var3_id,   cdata, p_id,   Time)
+    call send_perm_data(id_var2_id, cdata, PERM_TABLE(1:2, LAYOUT_XY), Time)
+    call send_perm_data(id_var3_id, cdata, PERM_TABLE(1:3, LAYOUT_XY), Time)
+
+    ! Swap: axes (y,x) / (z,y,x) while canonical storage remains (x,y,...) -> pack to temp and send
+    call send_perm_data(id_var2_yx, cdata, PERM_TABLE(1:2, LAYOUT_YX), Time)
+    call send_perm_data(id_var3_zx, cdata, PERM_TABLE(1:3, LAYOUT_ZX), Time)
 
-    ! Swap: axes (y,x) / (y,x,z) while canonical storage remains (x,y,...) -> pack to temp and send
-    call send_var2_perm(id_var2_swap, cdata, p_swap, Time)
-    call send_var3_perm(id_var3_swap, cdata, p_swap, Time)
+    ! Cyclic: axes (y,z,x) / (z,x,y) while canonical storage remains (x, y, ...) -> pack to temp and send
+    call send_perm_data(id_var3_yzx, cdata, PERM_TABLE(1:3, LAYOUT_YZX), Time)
+    call send_perm_data(id_var3_zxy, cdata, PERM_TABLE(1:3, LAYOUT_ZXY), Time)
 
     call diag_send_complete(Time_step)
     call diag_send_complete(Time_step)
@@ -136,63 +151,75 @@ program test_generalized_indices
 
 contains
 
-  subroutine send_var2_perm(id_field, buf, p, Time_in)
-    integer, intent(in)         :: id_field
-    real(r8_kind), intent(in)   :: buf(:,:,:,:)   ! canonical (x,y,z,w)
-    integer, intent(in)         :: p(3)
-    type(time_type), intent(in) :: Time_in
-
-    logical :: used_local
-    real(r8_kind), allocatable :: tmp2(:,:)
-
-    ! Support only identity (1,2,*) and xy-swap (2,1,*) for 2D
-    if (p(1)==1 .and. p(2)==2) then
-      used_local = send_data(id_field, buf(:,:,1,1), Time_in)
-    else if (p(1)==2 .and. p(2)==1) then
-      allocate(tmp2(size(buf,2), size(buf,1)))
-      tmp2 = transpose(buf(:,:,1,1))
-      used_local = send_data(id_field, tmp2, Time_in)
-      deallocate(tmp2)
-    else
-      call mpp_error(FATAL, 'send_var2_perm: only p=(1,2,*) or (2,1,*) implemented')
-    end if
-  end subroutine send_var2_perm
-
-
-  subroutine send_var3_perm(id_field, buf, p, Time_in)
+  !> @brief Apply a predefined permutation to an axis array.
+  !> Maps canonical axis ordering to a permuted layout using PERM_TABLE.
+  !> Supports rank-2 and rank-3 axis subsets.
+  subroutine permute_axis(axis_in, perm_id, axis_out)
+    integer, intent(in)  :: axis_in(:)
+    integer, intent(in)  :: perm_id
+    integer, intent(out) :: axis_out(:)
+
+    integer :: order(3)
+
+    order = PERM_TABLE(:, perm_id)
+
+    if (any(order(1:size(axis_out)) > size(axis_in))) then
+       call mpp_error(FATAL, "permute_axis: invalid permutation for given rank")
+    endif
+
+    axis_out = axis_in(order(1:size(axis_out)))
+  end subroutine permute_axis
+
+  !> @brief Register a diagnostic field with permuted axes.
+  !> Applies axis permutation before calling register_diag_field.
+  function register_perm_diag_field(var_name, long_name, axis, perm_id) result(id_var)
+    character(len=*), intent(in)  :: var_name, long_name
+    integer,          intent(in)  :: axis(:)
+    integer,          intent(in)  :: perm_id
+
+    integer              :: id_var
+    integer, allocatable :: axis_perm(:)
+
+    allocate(axis_perm(size(axis)))
+    call permute_axis(axis, perm_id, axis_perm)
+    id_var = register_diag_field('ocn_mod', var_name, axis_perm, Time, long_name, &
+                                 'mullions', missing_value=missing_value)
+    deallocate(axis_perm)
+  end function register_perm_diag_field
+
+  !> @brief Send data with optional axis permutation.
+  !> Applies 2D or 3D permutation to canonical (x,y,z,w) buffers before send_data.
+  !> Skips permutation for identity mappings.
+  subroutine send_perm_data(id_field, buf, order, Time_in)
     integer, intent(in)         :: id_field
     real(r8_kind), intent(in)   :: buf(:,:,:,:)   ! canonical (x,y,z,w)
-    integer, intent(in)         :: p(3)
+    integer, intent(in)         :: order(:)
     type(time_type), intent(in) :: Time_in
 
     logical :: used_local
-    integer :: nxloc, nyloc, nzloc, k
-    real(r8_kind), allocatable :: tmp3(:,:,:)
-
-    ! For now, support only keeping z as z
-    if (p(3) /= 3) call mpp_error(FATAL, 'send_var3_perm: only permutations with p(3)=3 implemented')
-
-    if (p(1)==1 .and. p(2)==2) then
-      used_local = send_data(id_field, buf(:,:,:,1), Time_in)
-
-    else if (p(1)==2 .and. p(2)==1) then
-      nxloc = size(buf,1)
-      nyloc = size(buf,2)
-      nzloc = size(buf,3)
-
-      allocate(tmp3(nyloc, nxloc, nzloc))
-      do k = 1, nzloc
-        tmp3(:,:,k) = transpose(buf(:,:,k,1))
-      end do
-
-      used_local = send_data(id_field, tmp3, Time_in)
-      deallocate(tmp3)
+    real(r8_kind), allocatable :: tmp2(:,:), tmp3(:,:,:)
+
+    if (size(order) == 2) then
+       if (all(order == [1,2])) then
+          used_local = send_data(id_field, buf(:,:,1,1), Time_in)
+       else
+          tmp2 = permute(buf(:,:,1,1), order)
+          used_local = send_data(id_field, tmp2, Time_in)
+       endif
+
+    else if (size(order) == 3) then
+       if (all(order == [1,2,3])) then
+          used_local = send_data(id_field, buf(:,:,:,1), Time_in)
+       else
+          tmp3 = permute(buf(:,:,:,1), order)
+          used_local = send_data(id_field, tmp3, Time_in)
+       endif
 
     else
-      call mpp_error(FATAL, 'send_var3_perm: only p=(1,2,3) or (2,1,3) implemented')
-    end if
-  end subroutine send_var3_perm
+       call mpp_error(FATAL, "send_var_perm: unsupported permutation rank")
 
+    endif
+  end subroutine send_perm_data
 
   !> @brief initialized the buffer based on the starting/ending indices
   subroutine init_buffer(buffer, is, ie, js, je, nhalo)
@@ -215,7 +242,6 @@ subroutine init_buffer(buffer, is, ie, js, je, nhalo)
     end do
   end subroutine init_buffer
 
-
   !> @brief Set the buffer based on the time_index
   subroutine set_buffer(buffer, time_index)
     real(r8_kind), intent(inout) :: buffer(:,:,:,:)