From 2e71b5e787debc4090f60dc99fd156e0f1672e64 Mon Sep 17 00:00:00 2001 From: rem1776 Date: Wed, 16 Jul 2025 11:53:10 -0400 Subject: [PATCH 1/5] fix issues with horiz_interp r4 code and unit tests --- horiz_interp/horiz_interp_bilinear.F90 | 1 + horiz_interp/include/horiz_interp_bilinear.inc | 18 +++++++++++++++--- horiz_interp/include/horiz_interp_conserve.inc | 15 --------------- test_fms/horiz_interp/Makefile.am | 4 ++-- test_fms/horiz_interp/test_horiz_interp.F90 | 18 ++++++++++-------- 5 files changed, 28 insertions(+), 28 deletions(-) diff --git a/horiz_interp/horiz_interp_bilinear.F90 b/horiz_interp/horiz_interp_bilinear.F90 index d8db732b22..174680f9d0 100644 --- a/horiz_interp/horiz_interp_bilinear.F90 +++ b/horiz_interp/horiz_interp_bilinear.F90 @@ -71,6 +71,7 @@ module horiz_interp_bilinear_mod !> @{ real(r8_kind), parameter :: epsln=1.e-10_r8_kind + real(r4_kind), parameter :: epsln_r4=1.e-4_r4_kind integer, parameter :: DUMMY = -999 !! Private helper routines, interfaces for mixed real precision support diff --git a/horiz_interp/include/horiz_interp_bilinear.inc b/horiz_interp/include/horiz_interp_bilinear.inc index f998b823f7..85a782866e 100644 --- a/horiz_interp/include/horiz_interp_bilinear.inc +++ b/horiz_interp/include/horiz_interp_bilinear.inc @@ -270,7 +270,11 @@ epsln2 = real(epsln,FMS_HI_KIND_)* 1.0e5_kindl call FIND_NEIGHBOR_NEW_(Interp, lon_in, lat_in, lon_out, lat_out, src_is_modulo, no_crash) else - epsln2 = real(epsln,FMS_HI_KIND_) + if(kindl == r8_kind) then + epsln2 = epsln + else + epsln2 = epsln_r4 + endif call FIND_NEIGHBOR_(Interp, lon_in, lat_in, lon_out, lat_out, src_is_modulo) endif @@ -389,8 +393,16 @@ if (x > 1.0_kindl) x = 1.0_kindl if (y > 1.0_kindl) y = 1.0_kindl endif - if( x>1.0_kindl .or. x<0.0_kindl .or. y>1.0_kindl .or. y < 0.0_kindl) & - call mpp_error(FATAL, "horiz_interp_bilinear_mod: weight should be between 0 and 1") + !! if using 4 byte reals, allow a bit more tolerance (1e-4) on weight values + if( kindl == r8_kind ) then + if( x>1.0_kindl .or. x<0.0_kindl .or. y>1.0_kindl .or. y < 0.0_kindl) & + call mpp_error(FATAL, "horiz_interp_bilinear_mod: weight should be between 0 and 1, x=" & + //string(x)//" y="//string(y)) + else + if( x>1.0_kindl + epsln2 .or. x<0.0_kindl - epsln2 .or. y>1.0_kindl + epsln2 .or. y < 0.0_kindl - epsln2) & + call mpp_error(FATAL, "horiz_interp_bilinear_mod: weight should be between 0 and 1, x=" & + //string(x)//" y="//string(y)) + endif Interp % HI_KIND_TYPE_ % wti(m,n,1)=1.0_kindl-x Interp % HI_KIND_TYPE_ % wti(m,n,2)=x Interp % HI_KIND_TYPE_ % wtj(m,n,1)=1.0_kindl-y diff --git a/horiz_interp/include/horiz_interp_conserve.inc b/horiz_interp/include/horiz_interp_conserve.inc index 1d2212dabc..9b38725424 100644 --- a/horiz_interp/include/horiz_interp_conserve.inc +++ b/horiz_interp/include/horiz_interp_conserve.inc @@ -246,17 +246,12 @@ subroutine HORIZ_INTERP_CONSERVE_NEW_1DX1D_ ( Interp, lon_in, lat_in, lon_out, l integer :: nincrease, ndecrease logical :: flip_lat - integer :: wordsz integer(kind=1) :: one_byte(8) integer, parameter :: kindl = FMS_HI_KIND_ !< compiled kind size if(.not. module_is_initialized) call mpp_error(FATAL, & 'HORIZ_INTERP_CONSERVE_NEW_1DX2D_: horiz_interp_conserve_init is not called') - wordsz=size(transfer(lon_in(1), one_byte)) - if(wordsz .NE. 4 .AND. wordsz .NE. 8) call mpp_error(FATAL, & - 'HORIZ_INTERP_CONSERVE_NEW_1DX2D_: wordsz should be 4 or 8') - if( (size(lon_out,1) .NE. size(lat_out,1)) .OR. (size(lon_out,2) .NE. size(lat_out,2)) ) & call mpp_error(FATAL, 'horiz_interp_conserve_mod: size mismatch between lon_out and lat_out') nlon_in = size(lon_in(:)) - 1; nlat_in = size(lat_in(:)) - 1 @@ -413,17 +408,12 @@ subroutine HORIZ_INTERP_CONSERVE_NEW_1DX1D_ ( Interp, lon_in, lat_in, lon_out, l real(r8_kind), allocatable, dimension(:) :: lon_out_r8, lat_out_r8 real(r8_kind), allocatable, dimension(:,:) :: lon_in_r8, lat_in_r8 real(r8_kind), allocatable, dimension(:,:) :: lon_dst, lat_dst - integer :: wordsz integer(kind=1) :: one_byte(8) integer, parameter :: kindl = FMS_HI_KIND_ !< compiled kind size if(.not. module_is_initialized) call mpp_error(FATAL, & 'HORIZ_INTERP_CONSERVE_NEW_2DX1D_: horiz_interp_conserve_init is not called') - wordsz=size(transfer(lon_in(1,1), one_byte)) - if(wordsz .NE. 8) call mpp_error(FATAL, & - 'HORIZ_INTERP_CONSERVE_NEW_2DX1D_: currently only support 64-bit real(FMS_HI_KIND_), contact developer') - if( (size(lon_in,1) .NE. size(lat_in,1)) .OR. (size(lon_in,2) .NE. size(lat_in,2)) ) & call mpp_error(FATAL, 'horiz_interp_conserve_mod: size mismatch between lon_in and lat_in') nlon_in = size(lon_in,1) - 1; nlat_in = size(lon_in,2) - 1 @@ -524,17 +514,12 @@ subroutine HORIZ_INTERP_CONSERVE_NEW_1DX1D_ ( Interp, lon_in, lat_in, lon_out, l real(r8_kind), allocatable, dimension(:,:) :: dst_area real(r8_kind), allocatable, dimension(:,:) :: lon_in_r8, lat_in_r8 real(r8_kind), allocatable, dimension(:,:) :: lon_out_r8, lat_out_r8 - integer :: wordsz integer(kind=1) :: one_byte(8) integer, parameter :: kindl = FMS_HI_KIND_ !< compiled kind size if(.not. module_is_initialized) call mpp_error(FATAL, & 'HORIZ_INTERP_CONSERVE_NEW_2DX2D_: horiz_interp_conserve_init is not called') - wordsz=size(transfer(lon_in(1,1), one_byte)) - if(wordsz .NE. 4 .AND. wordsz .NE. 8) call mpp_error(FATAL, & - 'HORIZ_INTERP_CONSERVE_NEW_2DX2D_: wordsz should be 4 or 8') - if( (size(lon_in,1) .NE. size(lat_in,1)) .OR. (size(lon_in,2) .NE. size(lat_in,2)) ) & call mpp_error(FATAL, 'horiz_interp_conserve_mod: size mismatch between lon_in and lat_in') if( (size(lon_out,1) .NE. size(lat_out,1)) .OR. (size(lon_out,2) .NE. size(lat_out,2)) ) & diff --git a/test_fms/horiz_interp/Makefile.am b/test_fms/horiz_interp/Makefile.am index 27de50a932..8b24b804fd 100644 --- a/test_fms/horiz_interp/Makefile.am +++ b/test_fms/horiz_interp/Makefile.am @@ -36,8 +36,8 @@ test_horiz_interp_r4_SOURCES = test_horiz_interp.F90 test_horiz_interp_r8_SOURCES = test_horiz_interp.F90 test_create_xgrid_order2_r8_SOURCES = test_create_xgrid_order2.F90 -test_horiz_interp_r4_CPPFLAGS=-DHI_TEST_KIND=4 -I$(MODDIR) -test_horiz_interp_r8_CPPFLAGS=-DHI_TEST_KIND=8 -I$(MODDIR) +test_horiz_interp_r4_CPPFLAGS=-DHI_TEST_KIND_=4 -I$(MODDIR) +test_horiz_interp_r8_CPPFLAGS=-DHI_TEST_KIND_=8 -I$(MODDIR) test_create_xgrid_order2_r8_CPPFLAGS=-DHI_TEST_KIND_=8 -I$(MODDIR) TEST_EXTENSIONS = .sh diff --git a/test_fms/horiz_interp/test_horiz_interp.F90 b/test_fms/horiz_interp/test_horiz_interp.F90 index 59ccdbb230..51f796d414 100644 --- a/test_fms/horiz_interp/test_horiz_interp.F90 +++ b/test_fms/horiz_interp/test_horiz_interp.F90 @@ -22,11 +22,6 @@ !! Assignment test checks that the override is copying the data type properly !! TODO some larger tests with different data sets -!! defaults to 8 real kind, make check will compile with both 4 and 8 -#ifndef HI_TEST_KIND_ -#define HI_TEST_KIND_ 8 -#endif - program horiz_interp_test use mpp_mod, only : mpp_init, mpp_exit, mpp_error, FATAL, stdout, mpp_npes, WARNING @@ -125,7 +120,7 @@ subroutine test_horiz_interp_spherical real(HI_TEST_KIND_) :: lon_dst_beg = -280._lkind, lon_dst_end = 80._lkind real(HI_TEST_KIND_) :: lat_dst_beg = -90._lkind, lat_dst_end = 90._lkind real(HI_TEST_KIND_) :: D2R = real(PI,HI_TEST_KIND_)/180._lkind - real(HI_TEST_KIND_), parameter :: SMALL = 1.0e-10_lkind + real(HI_TEST_KIND_), parameter :: tolerance = 1.0e-5_lkind ! set up longitude and latitude of source/destination grid. dlon_src = (lon_src_end-lon_src_beg)/real(ni_src, HI_TEST_KIND_) @@ -176,7 +171,7 @@ subroutine test_horiz_interp_spherical call mpp_clock_end(id1) do i=1, ni_dst-1 do j=1, nj_dst-1 - if(data_dst(i,j) - 1.0_lkind .gt. SMALL) then + if(data_dst(i,j) - 1.0_lkind .gt. tolerance) then print *, 'data_dst(i=', i, ', j=', j, ')=', data_dst(i,j), ' Expected value: 1.0' call mpp_error(FATAL, "test_horiz_interp_spherical: "// & "invalid output data after interpolation") @@ -649,7 +644,14 @@ subroutine test_horiz_interp_bicubic real(HI_TEST_KIND_) :: lon_dst_beg = -280._lkind, lon_dst_end = 80._lkind real(HI_TEST_KIND_) :: lat_dst_beg = -90._lkind, lat_dst_end = 90._lkind real(HI_TEST_KIND_) :: D2R = real(PI,HI_TEST_KIND_)/180._lkind - real(HI_TEST_KIND_), parameter :: SMALL = 1.0e-10_lkind + real(HI_TEST_KIND_) :: SMALL + + ! adjust tolerance used in checks based on kind size + if(HI_TEST_KIND_ == r8_kind) then + small = 1.0e-10_lkind + else + small = 1.0e-1_lkind + end if ! set up longitude and latitude of source/destination grid. dlon_src = (lon_src_end-lon_src_beg)/real(ni_src, lkind) From 153e2274a92c7aea0a4bec36f898a184ec5e5804 Mon Sep 17 00:00:00 2001 From: rem1776 Date: Wed, 16 Jul 2025 12:02:11 -0400 Subject: [PATCH 2/5] adjust tolerances for bicubic test --- test_fms/horiz_interp/test_horiz_interp.F90 | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/test_fms/horiz_interp/test_horiz_interp.F90 b/test_fms/horiz_interp/test_horiz_interp.F90 index 51f796d414..5f825b4cd7 100644 --- a/test_fms/horiz_interp/test_horiz_interp.F90 +++ b/test_fms/horiz_interp/test_horiz_interp.F90 @@ -650,7 +650,7 @@ subroutine test_horiz_interp_bicubic if(HI_TEST_KIND_ == r8_kind) then small = 1.0e-10_lkind else - small = 1.0e-1_lkind + small = 1.0e-3_lkind end if ! set up longitude and latitude of source/destination grid. @@ -718,13 +718,13 @@ subroutine test_horiz_interp_bicubic if( interp_t%horizInterpReals4_type%is_allocated) then if( interp_t%horizInterpReals4_type%wti(i,j,1) * interp_t%horizInterpReals4_type%wti(i,j,2) & - interp_t%horizInterpReals4_type%wti(i,j,3) .gt. SMALL .or. & - interp_t%horizInterpReals4_type%wti(i,j,3) - (57.2958_lkind * 57.2958_lkind) .gt. SMALL) then + interp_t%horizInterpReals4_type%wti(i,j,3) - (57.2958_lkind * 57.2958_lkind) .gt. 1.0e-1_lkind) then print *, i, j, interp_t%horizInterpReals4_type%wti(i,j,:) call mpp_error(FATAL, "test_horiz_interp: bicubic test failed 1Dx1D weight calculation") endif else if( interp_t%horizInterpReals8_type%wti(i,j,1) * interp_t%horizInterpReals8_type%wti(i,j,2) & - - interp_t%horizInterpReals8_type%wti(i,j,3) .gt. SMALL .and. & + - interp_t%horizInterpReals8_type%wti(i,j,3) .gt. SMALL .or. & interp_t%horizInterpReals8_type%wti(i,j,3) - (57.2958_lkind * 57.2958_lkind) .gt. SMALL) then print *, i, j, interp_t%horizInterpReals8_type%wti(i,j,:) call mpp_error(FATAL, "test_horiz_interp: bicubic test failed 1Dx1D weight calculation") @@ -764,7 +764,7 @@ subroutine test_horiz_interp_bicubic if( interp_t%horizInterpReals4_type%is_allocated) then if( interp_t%horizInterpReals4_type%wti(i,j,1) * interp_t%horizInterpReals4_type%wti(i,j,2) & - interp_t%horizInterpReals4_type%wti(i,j,3) .gt. SMALL .or. & - interp_t%horizInterpReals4_type%wti(i,j,3) - (57.2958_lkind * 57.2958_lkind) .gt. SMALL) then + interp_t%horizInterpReals4_type%wti(i,j,3) - (57.2958_lkind * 57.2958_lkind) .gt. 1.0e-1_lkind) then print *, i, j, interp_t%horizInterpReals4_type%wti(i,j,:) call mpp_error(FATAL, "test_horiz_interp: bicubic test failed 1Dx1D weight calculation") endif From fa397b1def4c5658758aee0877a4ec1adb0cdf65 Mon Sep 17 00:00:00 2001 From: rem1776 Date: Wed, 16 Jul 2025 12:08:31 -0400 Subject: [PATCH 3/5] set tolerance based on kind size for spherical --- test_fms/horiz_interp/test_horiz_interp.F90 | 8 +++++++- 1 file changed, 7 insertions(+), 1 deletion(-) diff --git a/test_fms/horiz_interp/test_horiz_interp.F90 b/test_fms/horiz_interp/test_horiz_interp.F90 index 5f825b4cd7..a99b5e9b98 100644 --- a/test_fms/horiz_interp/test_horiz_interp.F90 +++ b/test_fms/horiz_interp/test_horiz_interp.F90 @@ -120,7 +120,13 @@ subroutine test_horiz_interp_spherical real(HI_TEST_KIND_) :: lon_dst_beg = -280._lkind, lon_dst_end = 80._lkind real(HI_TEST_KIND_) :: lat_dst_beg = -90._lkind, lat_dst_end = 90._lkind real(HI_TEST_KIND_) :: D2R = real(PI,HI_TEST_KIND_)/180._lkind - real(HI_TEST_KIND_), parameter :: tolerance = 1.0e-5_lkind + real(HI_TEST_KIND_) :: tolerance + + if(HI_TEST_KIND_ == r8_kind) then + tolerance = 1.0e-10_r8_kind + else + tolerance = 1.0e-5_r4_kind + endif ! set up longitude and latitude of source/destination grid. dlon_src = (lon_src_end-lon_src_beg)/real(ni_src, HI_TEST_KIND_) From 134df597dd98a23cad7ff50463ea26ffefbe0faa Mon Sep 17 00:00:00 2001 From: rem1776 Date: Wed, 16 Jul 2025 12:36:58 -0400 Subject: [PATCH 4/5] fix line lengths --- horiz_interp/include/horiz_interp_bilinear.inc | 3 ++- test_fms/horiz_interp/test_horiz_interp.F90 | 6 ++++-- 2 files changed, 6 insertions(+), 3 deletions(-) diff --git a/horiz_interp/include/horiz_interp_bilinear.inc b/horiz_interp/include/horiz_interp_bilinear.inc index 85a782866e..8e54ae75e4 100644 --- a/horiz_interp/include/horiz_interp_bilinear.inc +++ b/horiz_interp/include/horiz_interp_bilinear.inc @@ -399,7 +399,8 @@ call mpp_error(FATAL, "horiz_interp_bilinear_mod: weight should be between 0 and 1, x=" & //string(x)//" y="//string(y)) else - if( x>1.0_kindl + epsln2 .or. x<0.0_kindl - epsln2 .or. y>1.0_kindl + epsln2 .or. y < 0.0_kindl - epsln2) & + if( x>1.0_kindl + epsln2 .or. x<0.0_kindl - epsln2 .or. & + y>1.0_kindl + epsln2 .or. y < 0.0_kindl - epsln2) & call mpp_error(FATAL, "horiz_interp_bilinear_mod: weight should be between 0 and 1, x=" & //string(x)//" y="//string(y)) endif diff --git a/test_fms/horiz_interp/test_horiz_interp.F90 b/test_fms/horiz_interp/test_horiz_interp.F90 index a99b5e9b98..e269893923 100644 --- a/test_fms/horiz_interp/test_horiz_interp.F90 +++ b/test_fms/horiz_interp/test_horiz_interp.F90 @@ -724,7 +724,8 @@ subroutine test_horiz_interp_bicubic if( interp_t%horizInterpReals4_type%is_allocated) then if( interp_t%horizInterpReals4_type%wti(i,j,1) * interp_t%horizInterpReals4_type%wti(i,j,2) & - interp_t%horizInterpReals4_type%wti(i,j,3) .gt. SMALL .or. & - interp_t%horizInterpReals4_type%wti(i,j,3) - (57.2958_lkind * 57.2958_lkind) .gt. 1.0e-1_lkind) then + interp_t%horizInterpReals4_type%wti(i,j,3) - (57.2958_lkind * 57.2958_lkind) & + .gt. 1.0e-1_lkind) then print *, i, j, interp_t%horizInterpReals4_type%wti(i,j,:) call mpp_error(FATAL, "test_horiz_interp: bicubic test failed 1Dx1D weight calculation") endif @@ -770,7 +771,8 @@ subroutine test_horiz_interp_bicubic if( interp_t%horizInterpReals4_type%is_allocated) then if( interp_t%horizInterpReals4_type%wti(i,j,1) * interp_t%horizInterpReals4_type%wti(i,j,2) & - interp_t%horizInterpReals4_type%wti(i,j,3) .gt. SMALL .or. & - interp_t%horizInterpReals4_type%wti(i,j,3) - (57.2958_lkind * 57.2958_lkind) .gt. 1.0e-1_lkind) then + interp_t%horizInterpReals4_type%wti(i,j,3) - (57.2958_lkind * 57.2958_lkind) .gt. 1.0e-1_lkind)& + then print *, i, j, interp_t%horizInterpReals4_type%wti(i,j,:) call mpp_error(FATAL, "test_horiz_interp: bicubic test failed 1Dx1D weight calculation") endif From c22390694a6e79ce22f732bfc8b809278efa271c Mon Sep 17 00:00:00 2001 From: rem1776 Date: Wed, 16 Jul 2025 16:24:57 -0400 Subject: [PATCH 5/5] fix kind name in constant --- test_fms/horiz_interp/test_horiz_interp.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/test_fms/horiz_interp/test_horiz_interp.F90 b/test_fms/horiz_interp/test_horiz_interp.F90 index e269893923..efeddd882a 100644 --- a/test_fms/horiz_interp/test_horiz_interp.F90 +++ b/test_fms/horiz_interp/test_horiz_interp.F90 @@ -123,9 +123,9 @@ subroutine test_horiz_interp_spherical real(HI_TEST_KIND_) :: tolerance if(HI_TEST_KIND_ == r8_kind) then - tolerance = 1.0e-10_r8_kind + tolerance = 1.0e-10_lkind else - tolerance = 1.0e-5_r4_kind + tolerance = 1.0e-5_lkind endif ! set up longitude and latitude of source/destination grid.