From d8beaf1de0223af8cd06140da3c046da3b03a0e5 Mon Sep 17 00:00:00 2001 From: Uriel Ramirez Date: Thu, 1 Jun 2023 13:22:38 -0400 Subject: [PATCH 01/16] update nearest_index to support motonically increasing/decreasing arrays --- axis_utils/include/axis_utils2.inc | 97 +++++++++++++++---------- test_fms/axis_utils/test_axis_utils.F90 | 44 +++++++---- test_fms/axis_utils/test_axis_utils2.sh | 2 +- 3 files changed, 89 insertions(+), 54 deletions(-) diff --git a/axis_utils/include/axis_utils2.inc b/axis_utils/include/axis_utils2.inc index baad2092e5..3943fd9c0a 100644 --- a/axis_utils/include/axis_utils2.inc +++ b/axis_utils/include/axis_utils2.inc @@ -267,7 +267,7 @@ !! inputs: !! !! value = arbitrary data...same units as elements in "array" - !! array = array of data points (must be monotonically increasing) + !! array = array of data points (must be monotonic) !! ia = dimension of "array" !! !! output: @@ -295,48 +295,71 @@ !! z(k1) would be the nearest data point to 12.5 and z(k2) would !! be the nearest data point to 0.0 !! @return integer nearest_index - - - function NEAREST_INDEX_(value, array) + real(kind=FMS_AU_KIND_) , intent(in) :: value !< arbitrary data...same units as elements in "array" + real(kind=FMS_AU_KIND_), dimension(:), intent(in) :: array !< array of data points (must be monotonic) integer :: NEAREST_INDEX_ integer :: ia !< dimension of "array" - integer :: i, ii, unit - real(kind=FMS_AU_KIND_) :: value !< arbitrary data...same units as elements in "array" - real(kind=FMS_AU_KIND_), dimension(:) :: array !< array of data points (must be monotonically increasing) - logical :: keep_going - - ia = size(array(:)) - - do i = 2, ia - if (array(i) < array(i-1)) then - unit = stdout() - write (unit,*) '=> Error: "nearest_index" array must be monotonically increasing' & - & // 'when searching for nearest value to ', value - write (unit,*) ' array(i) < array(i-1) for i=',i - write (unit,*) ' array(i) for i=1..ia follows:' - do ii = 1, ia - write (unit,*) 'i=',ii, ' array(i)=', array(ii) - enddo - call mpp_error(FATAL,' "nearest_index" array must be monotonically increasing.') - endif - enddo + integer :: i !< For do loops + + ia = SIZE(array(:)) + + ! check if array is monotonous + DO i = 2, ia-1 + IF( (array(i-1)array(i+1)) .OR. (array(i-1)>array(i).AND.array(i)array NOT monotonously ordered + CALL mpp_error(FATAL, 'axis_utils2::nearest_index array is NOT monotonously ordered') + END IF + END DO + + if (array(1) < array(ia)) then + !< increasing array + + !< Check if the value is outside the range of the array + if (value < array(1)) then + NEAREST_INDEX_ = 1 + return + elseif (value > array(ia)) then + NEAREST_INDEX_ = ia + return + endif - if (value < array(1) .or. value > array(ia)) then - if (value < array(1)) NEAREST_INDEX_ = 1 - if (value > array(ia)) NEAREST_INDEX_ = ia + DO i = 1, ia-1 + IF ( (array(i)<=value).AND.(array(i+1)>= value) ) THEN + IF( value - array(i) <= array(i+1) - value ) THEN + NEAREST_INDEX_ = i + return + ELSE + NEAREST_INDEX_ = i+1 + return + ENDIF + EXIT + END IF + END DO else - i = 1 - keep_going = .true. - do while (i <= ia .and. keep_going) - i = i+1 - if (value <= array(i)) then - NEAREST_INDEX_ = i - if (array(i)-value > value-array(i-1)) NEAREST_INDEX_ = i-1 - keep_going = .false. - endif - enddo + !< Decreasing Array + + !< Check if the value is outside the range of the array + if (value < array(ia)) then + NEAREST_INDEX_ = ia + return + elseif (value > array(1)) then + NEAREST_INDEX_ = 1 + return + endif + + DO i = 1, ia-1 + IF ( (array(i)>=value).AND.(array(i+1)<= value) ) THEN + IF ( array(i)-value <= value-array(i+1) ) THEN + NEAREST_INDEX_ = i + return + ELSE + NEAREST_INDEX_ = i+1 + return + END IF + END IF + END DO endif end function NEAREST_INDEX_ diff --git a/test_fms/axis_utils/test_axis_utils.F90 b/test_fms/axis_utils/test_axis_utils.F90 index aac74de010..c2765bae4b 100644 --- a/test_fms/axis_utils/test_axis_utils.F90 +++ b/test_fms/axis_utils/test_axis_utils.F90 @@ -76,9 +76,13 @@ program test_axis_utils print "(A)", "Testing frac_index (FAILURE)" call test_frac_index_fail - case ('--nearest-index') - print "(A)", "Testing nearest_index" - call test_nearest_index + case ('--nearest-index-increasing') + print "(A)", "Testing nearest_index with a monotonically increasing array" + call test_nearest_index(.true.) + + case ('--nearest-index-decreasing') + print "(A)", "Testing nearest_index with a monotonically decreasing array" + call test_nearest_index(.false.) case ('--nearest-index-fail') print "(A)", "Testing nearest_index (FAILURE)" @@ -385,28 +389,36 @@ subroutine test_frac_index_fail ret_test = frac_index(1.5_k, values) end subroutine -subroutine test_nearest_index +subroutine test_nearest_index(increasing_array) + logical, intent(in) :: increasing_array !< .True. if test using an increasing array real(k) :: arr(5) + integer :: ans(12) arr = [5._k, 12._k, 20._k, 40._k, 100._k] + if (increasing_array) then + ans=(/1, 5, 1, 2, 3, 4, 5, 1, 2, 2, 3, 3/) + else + arr = arr(ubound(arr,dim=1)::-1) + ans=(/5, 1, 5, 4, 3, 2, 1, 5, 4, 4, 3, 3/) + endif ! Test values beyond array boundaries - call nearest_index_assert(4._k, arr, 1) - call nearest_index_assert(1000._k, arr, size(arr)) + call nearest_index_assert(4._k, arr, ans(1)) + call nearest_index_assert(1000._k, arr, ans(2)) ! Test values actually in the array - call nearest_index_assert(5._k, arr, 1) - call nearest_index_assert(12._k, arr, 2) - call nearest_index_assert(20._k, arr, 3) - call nearest_index_assert(40._k, arr, 4) - call nearest_index_assert(100._k, arr, 5) + call nearest_index_assert(5._k, arr, ans(3)) + call nearest_index_assert(12._k, arr, ans(4)) + call nearest_index_assert(20._k, arr, ans(5)) + call nearest_index_assert(40._k, arr, ans(6)) + call nearest_index_assert(100._k, arr, ans(7)) ! Test the intervals between array values - call nearest_index_assert(6._k, arr, 1) - call nearest_index_assert(11._k, arr, 2) - call nearest_index_assert(15._k, arr, 2) - call nearest_index_assert(18._k, arr, 3) - call nearest_index_assert(29._k, arr, 3) + call nearest_index_assert(6._k, arr, ans(8)) + call nearest_index_assert(11._k, arr, ans(9)) + call nearest_index_assert(15._k, arr, ans(10)) + call nearest_index_assert(18._k, arr, ans(11)) + call nearest_index_assert(29._k, arr, ans(12)) end subroutine subroutine nearest_index_assert(val, arr, ret_expected) diff --git a/test_fms/axis_utils/test_axis_utils2.sh b/test_fms/axis_utils/test_axis_utils2.sh index 1288822481..eacfb89473 100755 --- a/test_fms/axis_utils/test_axis_utils2.sh +++ b/test_fms/axis_utils/test_axis_utils2.sh @@ -27,7 +27,7 @@ # Prepare the directory to run the tests. touch input.nml -TESTS_SUCCESS='--get-axis-modulo --get-axis-modulo-times --get-axis-cart --lon-in-range --frac-index --nearest-index --axis-edges --tranlon --interp-1d-1d --interp-1d-2d --interp-1d-3d' +TESTS_SUCCESS='--get-axis-modulo --get-axis-modulo-times --get-axis-cart --lon-in-range --frac-index --nearest-index-increasing --nearest-index-decreasing --axis-edges --tranlon --interp-1d-1d --interp-1d-2d --interp-1d-3d' TESTS_FAIL='--frac-index-fail --nearest-index-fail' # TODO: Enable these tests after tranlon's memory corruption bug is fixed. From aa814a61618c056652b6a73b9dfb701b72317a63 Mon Sep 17 00:00:00 2001 From: Uriel Ramirez Date: Wed, 11 Oct 2023 14:59:16 -0400 Subject: [PATCH 02/16] Fix the flip of arr --- test_fms/axis_utils/test_axis_utils.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/test_fms/axis_utils/test_axis_utils.F90 b/test_fms/axis_utils/test_axis_utils.F90 index c2765bae4b..615051471c 100644 --- a/test_fms/axis_utils/test_axis_utils.F90 +++ b/test_fms/axis_utils/test_axis_utils.F90 @@ -394,11 +394,11 @@ subroutine test_nearest_index(increasing_array) real(k) :: arr(5) integer :: ans(12) - arr = [5._k, 12._k, 20._k, 40._k, 100._k] if (increasing_array) then + arr = [5._k, 12._k, 20._k, 40._k, 100._k] ans=(/1, 5, 1, 2, 3, 4, 5, 1, 2, 2, 3, 3/) else - arr = arr(ubound(arr,dim=1)::-1) + arr = [100._k, 40._k, 20._k, 12._k, 5._k] ans=(/5, 1, 5, 4, 3, 2, 1, 5, 4, 4, 3, 3/) endif From 7ea4befd77048600b338802ec5d62562e6042766 Mon Sep 17 00:00:00 2001 From: Uriel Ramirez Date: Mon, 23 Oct 2023 15:07:32 -0400 Subject: [PATCH 03/16] Fix axis_edges for decreasing grids --- axis_utils/include/axis_utils2.inc | 7 +++++-- 1 file changed, 5 insertions(+), 2 deletions(-) diff --git a/axis_utils/include/axis_utils2.inc b/axis_utils/include/axis_utils2.inc index 3943fd9c0a..c142c0126e 100644 --- a/axis_utils/include/axis_utils2.inc +++ b/axis_utils/include/axis_utils2.inc @@ -46,6 +46,7 @@ !! indicating to reproduce the mpp_io bug where !! the null characters were not removed after reading a string attribute integer, parameter :: lkind = FMS_AU_KIND_ + integer :: edge_index(2) ndims = get_variable_num_dimensions(fileobj, name) allocate(dim_sizes(ndims)) @@ -108,8 +109,10 @@ allocate(r2d(dim_sizes(1), dim_sizes(2))) call read_data(fileobj, buffer, r2d) - edge_data(1:dim_sizes(2)) = r2d(1,:) - edge_data(dim_sizes(2)+1) = r2d(2,dim_sizes(2)) + edge_index = (/1, 2/) + if (r2d(1,1) .gt. r2d(1,2)) edge_index = (/2, 1 /) + edge_data(1:dim_sizes(2)) = r2d(edge_index(1),:) + edge_data(dim_sizes(2)+1) = r2d(edge_index(2),dim_sizes(2)) deallocate(r2d) endif deallocate(dim_sizes) From 296f2696cd91e14a447b711ad03f2cd98647ac3a Mon Sep 17 00:00:00 2001 From: Uriel Ramirez Date: Tue, 24 Oct 2023 17:21:25 -0400 Subject: [PATCH 04/16] Add tests to test axis_eges for decreasing array --- test_fms/axis_utils/test_axis_utils.F90 | 42 +++++++++++++++++++------ test_fms/axis_utils/test_axis_utils2.sh | 4 +-- 2 files changed, 34 insertions(+), 12 deletions(-) diff --git a/test_fms/axis_utils/test_axis_utils.F90 b/test_fms/axis_utils/test_axis_utils.F90 index 615051471c..c3a431b136 100644 --- a/test_fms/axis_utils/test_axis_utils.F90 +++ b/test_fms/axis_utils/test_axis_utils.F90 @@ -88,9 +88,13 @@ program test_axis_utils print "(A)", "Testing nearest_index (FAILURE)" call test_nearest_index_fail - case ('--axis-edges') - print "(A)", "Testing axis_edges" - call test_axis_edges + case ('--axis-edges-increasing') + print "(A)", "Testing axis_edges-increasing" + call test_axis_edges(.true.) + + case ('--axis-edges-decreasing') + print "(A)", "Testing axis_edges-decreasing" + call test_axis_edges(.false.) case ('--tranlon') print "(A)", "Testing tranlon" @@ -445,24 +449,42 @@ subroutine test_nearest_index_fail ret_test = nearest_index(5._k, arr) end subroutine -subroutine test_axis_edges +subroutine test_axis_edges(increasing_array) + logical, intent(in) :: increasing_array !< .True. if test using an increasing array real(k) :: data_in_var(10) real(k) :: data_in_var_edges(2,10) real(k) :: data_in_answers(11) type(FmsNetcdfFile_t) :: fileobj real(k) :: answers(11) - integer :: i + integer :: count + integer :: count_factor + integer :: factor + + if (increasing_array) then + count = 0 + factor = 1 + count_factor = -1 + else + count = 11 + factor = -1 + count_factor = 0 + endif do i=1,10 - data_in_var(i) = real(i, k) - 0.5_k + count = count + factor + data_in_var(i) = real(count, k) - 0.5_k - data_in_var_edges(1,i) = real(i-1, k) - data_in_var_edges(2,i) = real(i, k) + data_in_var_edges(1,i) = real(count-1, k) + data_in_var_edges(2,i) = real(count, k) - data_in_answers(i) = real(i-1, k) + data_in_answers(i) = real(count + count_factor, k) enddo - data_in_answers(11) = 10._k + if (increasing_array) then + data_in_answers(11) = real(count, k) + else + data_in_answers(11) = real(count + factor, k) + endif call open_netcdf_w(fileobj) diff --git a/test_fms/axis_utils/test_axis_utils2.sh b/test_fms/axis_utils/test_axis_utils2.sh index eacfb89473..8b49c6bdd5 100755 --- a/test_fms/axis_utils/test_axis_utils2.sh +++ b/test_fms/axis_utils/test_axis_utils2.sh @@ -27,11 +27,11 @@ # Prepare the directory to run the tests. touch input.nml -TESTS_SUCCESS='--get-axis-modulo --get-axis-modulo-times --get-axis-cart --lon-in-range --frac-index --nearest-index-increasing --nearest-index-decreasing --axis-edges --tranlon --interp-1d-1d --interp-1d-2d --interp-1d-3d' +TESTS_SUCCESS='--get-axis-modulo --get-axis-modulo-times --get-axis-cart --lon-in-range --frac-index --nearest-index-increasing --nearest-index-decreasing --axis-edges-increasing --axis_edges-decreasing --tranlon --interp-1d-1d --interp-1d-2d --interp-1d-3d' TESTS_FAIL='--frac-index-fail --nearest-index-fail' # TODO: Enable these tests after tranlon's memory corruption bug is fixed. -SKIP_TESTS="test_axis_utils2.15 test_axis_utils2.16" +SKIP_TESTS="test_axis_utils2.19 test_axis_utils2.20" # Run the tests From b732e6e758eafb0c3ff456be2176c36c253fcec1 Mon Sep 17 00:00:00 2001 From: Uriel Ramirez Date: Wed, 6 Dec 2023 09:36:57 -0500 Subject: [PATCH 05/16] minor documentation, clean up unused variables --- axis_utils/include/axis_utils2.inc | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/axis_utils/include/axis_utils2.inc b/axis_utils/include/axis_utils2.inc index 315c0d2250..bac9251e07 100644 --- a/axis_utils/include/axis_utils2.inc +++ b/axis_utils/include/axis_utils2.inc @@ -46,7 +46,9 @@ !! indicating to reproduce the mpp_io bug where !! the null characters were not removed after reading a string attribute integer, parameter :: lkind = FMS_AU_KIND_ - integer :: edge_index(2) + integer :: edge_index(2) !< Index to use when reading the edges from the file + !! (/1, 2/) if the axis data is monotonically increasing + !! (/2, 1/) if the axis data is monotonically decreasing ndims = get_variable_num_dimensions(fileobj, name) allocate(dim_sizes(ndims)) @@ -302,8 +304,7 @@ integer :: NEAREST_INDEX_ integer :: ia !< dimension of "array" - integer :: i, ii, iunit - logical :: keep_going + integer :: i !< For looping through "array" ia = SIZE(array(:)) From 0800d27cb63fbc90dae5afa097b10f3e00eb6a7e Mon Sep 17 00:00:00 2001 From: Uriel Ramirez Date: Wed, 6 Dec 2023 10:01:49 -0500 Subject: [PATCH 06/16] Implement axis_utils 'nearest_index' in horiz_interp_bilinear, remove 'indp' --- horiz_interp/horiz_interp_bilinear.F90 | 7 +- .../include/horiz_interp_bilinear.inc | 87 ++++++++----------- .../include/horiz_interp_bilinear_r4.fh | 3 - .../include/horiz_interp_bilinear_r8.fh | 3 - 4 files changed, 35 insertions(+), 65 deletions(-) diff --git a/horiz_interp/horiz_interp_bilinear.F90 b/horiz_interp/horiz_interp_bilinear.F90 index 64abf15263..318d2c039b 100644 --- a/horiz_interp/horiz_interp_bilinear.F90 +++ b/horiz_interp/horiz_interp_bilinear.F90 @@ -34,6 +34,7 @@ module horiz_interp_bilinear_mod use constants_mod, only: PI use horiz_interp_type_mod, only: horiz_interp_type, stats use platform_mod, only: r4_kind, r8_kind + use axis_utils2_mod, only: nearest_index implicit none private @@ -63,12 +64,6 @@ module horiz_interp_bilinear_mod integer, parameter :: DUMMY = -999 !! Private helper routines, interfaces for mixed real precision support - - interface indp - module procedure indp_r4 - module procedure indp_r8 - end interface - interface intersect module procedure intersect_r4 module procedure intersect_r8 diff --git a/horiz_interp/include/horiz_interp_bilinear.inc b/horiz_interp/include/horiz_interp_bilinear.inc index 8f81e86f79..5bac5e4b01 100644 --- a/horiz_interp/include/horiz_interp_bilinear.inc +++ b/horiz_interp/include/horiz_interp_bilinear.inc @@ -34,12 +34,20 @@ real(FMS_HI_KIND_) :: wtw, wte, wts, wtn, lon, lat, tpi, hpi real(FMS_HI_KIND_) :: glt_min, glt_max, gln_min, gln_max, min_lon, max_lon integer,parameter :: kindl = FMS_HI_KIND_ + logical :: decreasing_lat !< .True. if latitude is monotically decreasing + logical :: decreasing_lon !< .True. if longitude is monotically decreasing warns = 0 if(present(verbose)) warns = verbose src_is_modulo = .true. if (present(src_modulo)) src_is_modulo = src_modulo + decreasing_lat = .false. + if (lat_in(1) > lat_in(2)) decreasing_lat = .true. + + decreasing_lon = .false. + if (lon_in(1) > lon_in(2)) decreasing_lon = .true. + hpi = 0.5_kindl * real(pi, FMS_HI_KIND_) tpi = 4.0_kindl * hpi glt_min = hpi @@ -95,11 +103,24 @@ glt_min = min(lat,glt_min); glt_max = max(lat,glt_max) gln_min = min(lon,gln_min); gln_max = max(lon,gln_max) - is = indp(lon, lon_in ) - if( lon_in(is) .gt. lon ) is = max(is-1,1) + is = nearest_index(lon, lon_in ) + if (decreasing_lon) then + ! Lon_in is decreasing + ! This is so that is is the lower bound. + ! For example, if the array is [10 20 30 40 50] and lon is 49, `is` is going to be 5 from `nearest_index` + ! but it needs to be 4 so that it can use the data at lon_in(4) and lon_in(5) + if( lon_in(is) .gt. lon ) is = max(is-1,1) + else + ! Lon_in is increasing + ! This is so that is is the lower bound. + ! For example, if the array is [50 40 30 20 10] and lon is 11, `is` is going to be 5 from `nearest_index` + ! but it needs to be 4 so that it can use the data at lon_in(4) and lon_in(5) + if( lon_in(is) .lt. lon ) is = max(is-1,1) + endif + if( lon_in(is) .eq. lon .and. is .eq. nlon_in) is = max(is - 1,1) ie = min(is+1,nlon_in) - if(lon_in(is) .ne. lon_in(ie) .and. lon_in(is) .le. lon) then + if(lon_in(is) .ne. lon_in(ie) .and. (decreasing_lon .or. lon_in(is) .le. lon)) then wtw = ( lon_in(ie) - lon) / (lon_in(ie) - lon_in(is) ) else ! east or west of the last data value. this could be because a @@ -115,13 +136,20 @@ endif wte = 1.0_kindl - wtw - js = indp(lat, lat_in ) + js = nearest_index(lat, lat_in ) + + if (decreasing_lat) then + ! Lat_in is decreasing + if( lat_in(js) .lt. lat ) js = max(js - 1, 1) + else + ! Lat_in is increasing + if( lat_in(js) .gt. lat ) js = max(js - 1, 1) + endif - if( lat_in(js) .gt. lat ) js = max(js - 1, 1) if( lat_in(js) .eq. lat .and. js .eq. nlat_in) js = max(js - 1, 1) je = min(js + 1, nlat_in) - if ( lat_in(js) .ne. lat_in(je) .and. lat_in(js) .le. lat) then + if ( lat_in(js) .ne. lat_in(je) .and. (decreasing_lat .or. lat_in(js) .le. lat)) then wts = ( lat_in(je) - lat )/(lat_in(je)-lat_in(js)) else ! north or south of the last data value. this could be because a @@ -1180,51 +1208,4 @@ return end subroutine - - - !####################################################################### - !> @returns index of nearest data point to "value" - !! if "value" is outside the domain of "array" then INDP_ = 1 - !! or "ia" depending on whether array(1) or array(ia) is - !! closest to "value" - function INDP_ (rval, array) - integer :: INDP_ !< index of nearest data point within "array" - !! corresponding to "value". - real(FMS_HI_KIND_), dimension(:), intent(in) :: array !< array of data points (must be monotonically increasing) - real(FMS_HI_KIND_), intent(in) :: rval !< arbitrary data, same units as elements in 'array' - - !======================================================================= - - integer i, ia, iunit - logical keep_going - ! - ia = size(array(:)) - do i=2,ia - if (array(i) .lt. array(i-1)) then - iunit = stdout() - write (iunit,*) & - ' => Error: array must be monotonically increasing in "INDP_"' , & - ' when searching for nearest element to value=',rval - write (iunit,*) ' array(i) < array(i-1) for i=',i - write (iunit,*) ' array(i) for i=1..ia follows:' - call mpp_error() - endif - enddo - if (rval .lt. array(1) .or. rval .gt. array(ia)) then - if (rval .lt. array(1)) INDP_ = 1 - if (rval .gt. array(ia)) INDP_ = ia - else - i=1 - keep_going = .true. - do while (i .le. ia .and. keep_going) - i = i+1 - if (rval .le. array(i)) then - INDP_ = i - if (array(i)-rval .gt. rval-array(i-1)) INDP_ = i-1 - keep_going = .false. - endif - enddo - endif - return - end function INDP_ !> @} diff --git a/horiz_interp/include/horiz_interp_bilinear_r4.fh b/horiz_interp/include/horiz_interp_bilinear_r4.fh index 07cf041ccc..8880914e43 100644 --- a/horiz_interp/include/horiz_interp_bilinear_r4.fh +++ b/horiz_interp/include/horiz_interp_bilinear_r4.fh @@ -45,8 +45,5 @@ #undef INTERSECT_ #define INTERSECT_ intersect_r4 -#undef INDP_ -#define INDP_ indp_r4 - #include "horiz_interp_bilinear.inc" !> @} diff --git a/horiz_interp/include/horiz_interp_bilinear_r8.fh b/horiz_interp/include/horiz_interp_bilinear_r8.fh index 526fecfb69..37a2e6920b 100644 --- a/horiz_interp/include/horiz_interp_bilinear_r8.fh +++ b/horiz_interp/include/horiz_interp_bilinear_r8.fh @@ -45,8 +45,5 @@ #undef INTERSECT_ #define INTERSECT_ intersect_r8 -#undef INDP_ -#define INDP_ indp_r8 - #include "horiz_interp_bilinear.inc" !> @} From 6b035b1d6eec3605846c74dd5fcb320c773bb991 Mon Sep 17 00:00:00 2001 From: Uriel Ramirez Date: Wed, 6 Dec 2023 11:04:40 -0500 Subject: [PATCH 07/16] add a test for decreasing grids in horiz_interp --- test_fms/horiz_interp/test_horiz_interp.F90 | 18 +++++++-- test_fms/horiz_interp/test_horiz_interp2.sh | 44 +++++++++++++++++++-- 2 files changed, 54 insertions(+), 8 deletions(-) diff --git a/test_fms/horiz_interp/test_horiz_interp.F90 b/test_fms/horiz_interp/test_horiz_interp.F90 index c4e333b786..fd0d077a91 100644 --- a/test_fms/horiz_interp/test_horiz_interp.F90 +++ b/test_fms/horiz_interp/test_horiz_interp.F90 @@ -52,10 +52,11 @@ program horiz_interp_test integer :: ni_dst = 144, nj_dst = 72 integer, parameter :: max_neighbors = 400 !! took this from spherical mod !! max amount found neighbors to loop through in spherical search - + logical :: decreasing_lat = .False. !< latitude and longitude are decreasing instead of increasing + !! This is only for the bilinear test case namelist /test_horiz_interp_nml/ test_conserve, test_bicubic, test_spherical, test_bilinear, test_assign, test_solo,& - ni_src, nj_src, ni_dst,nj_dst + ni_src, nj_src, ni_dst,nj_dst, decreasing_lat type(domain2d) :: domain @@ -196,12 +197,20 @@ subroutine test_horiz_interp_bilinear real(HI_TEST_KIND_), allocatable, dimension(:) :: lon1D_src, lat1D_src, lon1D_dst, lat1D_dst real(HI_TEST_KIND_), allocatable, dimension(:,:) :: lon2D_src, lat2d_src, lon2D_dst, lat2D_dst real(HI_TEST_KIND_), allocatable, dimension(:,:) :: data_src, data_dst - real(HI_TEST_KIND_), parameter :: lon_src_beg = 0._lkind, lon_src_end = 360.0_lkind - real(HI_TEST_KIND_), parameter :: lat_src_beg = -90._lkind, lat_src_end = 90._lkind + real(HI_TEST_KIND_) :: lon_src_beg = 0._lkind, lon_src_end = 360.0_lkind + real(HI_TEST_KIND_) :: lat_src_beg = -90._lkind, lat_src_end = 90._lkind real(HI_TEST_KIND_), parameter :: D2R = real(PI,lkind)/180._lkind + real(HI_TEST_KIND_), parameter :: R2D = 180._lkind/real(PI,lkind) type(horiz_interp_type) :: interp + if (decreasing_lat) then + lon_src_beg = 360.0_lkind + lon_src_end = 0._lkind + lat_src_beg = 90._lkind + lat_src_end = -90._lkind + endif + allocate( lon1D_src(ni_src+1), lat1D_src(nj_src+1) ) allocate( lon1D_dst(ni_src+1), lat1D_dst(nj_src+1) ) allocate( lon2d_src(ni_src,nj_src), lat2d_src(ni_src,nj_src) ) @@ -377,6 +386,7 @@ subroutine test_horiz_interp_bilinear call check_dealloc(interp) endif + if (decreasing_lat) return ! --- 2dx1d version bilinear interpolation data_dst = 0.0_lkind lon1d_dst = lon1d_src diff --git a/test_fms/horiz_interp/test_horiz_interp2.sh b/test_fms/horiz_interp/test_horiz_interp2.sh index 485be882fd..cb8868d4e6 100755 --- a/test_fms/horiz_interp/test_horiz_interp2.sh +++ b/test_fms/horiz_interp/test_horiz_interp2.sh @@ -111,10 +111,45 @@ cat <<_EOF > input.nml / _EOF -test_expect_success "bilinear method with real kind=4" ' +test_expect_success "bilinear method with real kind=4 increasing latitude/longitude" ' mpirun -n 2 ./test_horiz_interp_r4 ' -test_expect_success "bilinear method with real kind=8" ' +test_expect_success "bilinear method with real kind=8 increasing latitude/longitude" ' + mpirun -n 2 ./test_horiz_interp_r8 +' + +cat <<_EOF > input.nml +&test_horiz_interp_nml + test_bilinear= .true. + ni_src = 360 + nj_src = 180 + ni_dst = 144 + nj_dst = 72 + decreasing_lat = .true. +/ +_EOF + +test_expect_success "bilinear method with real kind=4 decreasing latitude/longitude" ' + mpirun -n 2 ./test_horiz_interp_r4 +' +test_expect_success "bilinear method with real kind=8 decreasing latitude/longitude" ' + mpirun -n 2 ./test_horiz_interp_r8 +' +cat <<_EOF > input.nml +&test_horiz_interp_nml + test_bilinear= .true. + test_solo = .true. + ni_src = 360 + nj_src = 180 + ni_dst = 144 + nj_dst = 72 +/ +_EOF + +test_expect_success "bilinear method solo wrapper with real kind=4 increasing latitude/longitude" ' + mpirun -n 2 ./test_horiz_interp_r4 +' +test_expect_success "bilinear method solo wrapper with real kind=8 increasing latitude/longitude" ' mpirun -n 2 ./test_horiz_interp_r8 ' @@ -126,13 +161,14 @@ cat <<_EOF > input.nml nj_src = 180 ni_dst = 144 nj_dst = 72 + decreasing_lat = .true. / _EOF -test_expect_success "bilinear method solo wrapper with real kind=4" ' +test_expect_success "bilinear method solo wrapper with real kind=4 decreasing latitude/longitude" ' mpirun -n 2 ./test_horiz_interp_r4 ' -test_expect_success "bilinear method solo wrapper with real kind=8" ' +test_expect_success "bilinear method solo wrapper with real kind=8 decreasing latitude/longitude" ' mpirun -n 2 ./test_horiz_interp_r8 ' From 87f9e229f876b6716d73ae346adba70a29cf2723 Mon Sep 17 00:00:00 2001 From: Uriel Ramirez Date: Wed, 6 Dec 2023 11:06:02 -0500 Subject: [PATCH 08/16] rename loop index --- test_fms/axis_utils/test_axis_utils.F90 | 11 ++++++----- 1 file changed, 6 insertions(+), 5 deletions(-) diff --git a/test_fms/axis_utils/test_axis_utils.F90 b/test_fms/axis_utils/test_axis_utils.F90 index c3a431b136..c86d2ee4c4 100644 --- a/test_fms/axis_utils/test_axis_utils.F90 +++ b/test_fms/axis_utils/test_axis_utils.F90 @@ -459,6 +459,7 @@ subroutine test_axis_edges(increasing_array) integer :: count integer :: count_factor integer :: factor + integer :: index !< For looping through the data if (increasing_array) then count = 0 @@ -470,14 +471,14 @@ subroutine test_axis_edges(increasing_array) count_factor = 0 endif - do i=1,10 + do index=1,10 count = count + factor - data_in_var(i) = real(count, k) - 0.5_k + data_in_var(index) = real(count, k) - 0.5_k - data_in_var_edges(1,i) = real(count-1, k) - data_in_var_edges(2,i) = real(count, k) + data_in_var_edges(1,index) = real(count-1, k) + data_in_var_edges(2,index) = real(count, k) - data_in_answers(i) = real(count + count_factor, k) + data_in_answers(index) = real(count + count_factor, k) enddo if (increasing_array) then From 2e592e0dda5c5d2f95a27023ef8a68b776d7a29b Mon Sep 17 00:00:00 2001 From: Uriel Ramirez Date: Fri, 8 Dec 2023 13:06:02 -0500 Subject: [PATCH 09/16] fix data_override for decreasing grid + tests (this is mostly broken) --- data_override/include/data_override.inc | 16 +- test_fms/data_override/Makefile.am | 4 +- .../data_override/test_data_override2_mono.sh | 45 ++ .../test_data_override_ongrid.F90 | 448 ++++++++++++++---- 4 files changed, 400 insertions(+), 113 deletions(-) create mode 100755 test_fms/data_override/test_data_override2_mono.sh diff --git a/data_override/include/data_override.inc b/data_override/include/data_override.inc index 305e58ae6c..ad7cc6e94a 100644 --- a/data_override/include/data_override.inc +++ b/data_override/include/data_override.inc @@ -950,14 +950,14 @@ subroutine DATA_OVERRIDE_3D_(gridname,fieldname_code,return_data,time,override,d ie_src = axis_sizes(1) js_src = 1 je_src = axis_sizes(2) - do j = 1, axis_sizes(2)+1 - if( override_array(curr_position)%lat_in(j) > lat_min ) exit - js_src = j - enddo - do j = 1, axis_sizes(2)+1 - je_src = j - if( override_array(curr_position)%lat_in(j) > lat_max ) exit - enddo + ! do j = 1, axis_sizes(2)+1 + ! if( override_array(curr_position)%lat_in(j) > lat_min ) exit + ! js_src = j + ! enddo + ! do j = 1, axis_sizes(2)+1 + ! je_src = j + ! if( override_array(curr_position)%lat_in(j) > lat_max ) exit + ! enddo !--- bicubic interpolation need one extra point in each direction. Also add !--- one more point for because lat_in is in the corner but the interpolation diff --git a/test_fms/data_override/Makefile.am b/test_fms/data_override/Makefile.am index 366c773800..8b6114b88f 100644 --- a/test_fms/data_override/Makefile.am +++ b/test_fms/data_override/Makefile.am @@ -71,10 +71,10 @@ TESTS_ENVIRONMENT= test_input_path="@TEST_INPUT_PATH@" \ parser_skip=${skipflag} # Run the test program. -TESTS = test_data_override2.sh +TESTS = test_data_override2.sh test_data_override2_mono.sh # Include these files with the distribution. -EXTRA_DIST = test_data_override2.sh +EXTRA_DIST = test_data_override2.sh test_data_override2_mono.sh # Clean up CLEANFILES = input.nml *.nc* *.out diag_table data_table data_table.yaml INPUT/* *.dpi *.spi *.dyn *.spl diff --git a/test_fms/data_override/test_data_override2_mono.sh b/test_fms/data_override/test_data_override2_mono.sh new file mode 100755 index 0000000000..63ee989d01 --- /dev/null +++ b/test_fms/data_override/test_data_override2_mono.sh @@ -0,0 +1,45 @@ +#!/bin/sh + +#*********************************************************************** +#* GNU Lesser General Public License +#* +#* This file is part of the GFDL Flexible Modeling System (FMS). +#* +#* FMS is free software: you can redistribute it and/or modify it under +#* the terms of the GNU Lesser General Public License as published by +#* the Free Software Foundation, either version 3 of the License, or (at +#* your option) any later version. +#* +#* FMS is distributed in the hope that it will be useful, but WITHOUT +#* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +#* FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +#* for more details. +#* +#* You should have received a copy of the GNU Lesser General Public +#* License along with FMS. If not, see . +#*********************************************************************** +# +# Copyright (c) 2019-2021 Ed Hartnett, Uriel Ramirez, Seth Underwood + +# Set common test settings. +. ../test-lib.sh + +output_dir +[ ! -d "INPUT" ] && mkdir -p "INPUT" + +cat <<_EOF > input.nml +&test_data_override_ongrid_nml + test_case = 2 +/ +_EOF + +cat <<_EOF > data_table +"OCN", "runoff_increasing", "runoff", "./INPUT/bilinear_increasing.nc", "bilinear" , 1.0 +"OCN", "runoff_decreasing", "runoff", "./INPUT/bilinear_decreasing.nc", "bilinear" , 1.0 +_EOF + +test_expect_success "test_data_override with monotonically increasing and decreasing data sets" ' + mpirun -n 6 ../test_data_override_ongrid_r4 + ' + +test_done diff --git a/test_fms/data_override/test_data_override_ongrid.F90 b/test_fms/data_override/test_data_override_ongrid.F90 index 3f031547fa..99140865af 100644 --- a/test_fms/data_override/test_data_override_ongrid.F90 +++ b/test_fms/data_override/test_data_override_ongrid.F90 @@ -28,11 +28,12 @@ program test_data_override_ongrid use mpp_mod, only: mpp_init, mpp_exit, mpp_pe, mpp_root_pe, mpp_error, FATAL, & input_nml_file, mpp_sync use data_override_mod, only: data_override_init, data_override -use fms2_io_mod, only: fms2_io_init +use fms2_io_mod use time_manager_mod, only: set_calendar_type, time_type, set_date, NOLEAP use netcdf, only: nf90_create, nf90_def_dim, nf90_def_var, nf90_enddef, nf90_put_var, & nf90_close, nf90_put_att, nf90_clobber, nf90_64bit_offset, nf90_char, & nf90_double, nf90_unlimited +use fms_mod, only: string implicit none @@ -41,23 +42,17 @@ program test_data_override_ongrid integer :: nlon !< Number of points in x axis integer :: nlat !< Number of points in y axis type(domain2d) :: Domain !< Domain with mask table -real(lkind), allocatable, dimension(:,:) :: runoff !< Data to be written integer :: is !< Starting x index integer :: ie !< Ending x index integer :: js !< Starting y index integer :: je !< Ending y index -type(time_type) :: Time !< Time -integer :: i !< Helper indices -integer :: ncid !< Netcdf file id -integer :: err !< Error Code -integer :: dim1d, dim2d, dim3d, dim4d !< Dimension ids -integer :: varid, varid2, varid3, varid4 !< Variable ids -real(lkind), allocatable, dimension(:,:,:) :: runoff_in !< Data to be written to file -real(lkind) :: expected_result !< Expected result from data_override integer :: nhalox=2, nhaloy=2 integer :: io_status +integer, parameter :: ongrid = 1 +integer, parameter :: bilinear = 2 +integer :: test_case = ongrid -namelist / test_data_override_ongrid_nml / nhalox, nhaloy +namelist / test_data_override_ongrid_nml / nhalox, nhaloy, test_case call mpp_init call fms2_io_init @@ -65,64 +60,7 @@ program test_data_override_ongrid read (input_nml_file, test_data_override_ongrid_nml, iostat=io_status) if (io_status > 0) call mpp_error(FATAL,'=>test_data_override_ongrid: Error reading input.nml') -!< Create some files needed by data_override! -if (mpp_pe() .eq. mpp_root_pe()) then - allocate(runoff_in(1440, 1080, 10)) - do i = 1, 10 - runoff_in(:,:,i) = real(i, lkind) - enddo - err = nf90_create('INPUT/grid_spec.nc', ior(nf90_clobber, nf90_64bit_offset), ncid) - err = nf90_def_dim(ncid, 'str', 255, dim1d) - err = nf90_def_var(ncid, 'ocn_mosaic_file', nf90_char, (/dim1d/), varid) - err = nf90_enddef(ncid) - err = nf90_put_var(ncid, varid, "ocean_mosaic.nc") - err = nf90_close(ncid) - - err = nf90_create('INPUT/ocean_mosaic.nc', ior(nf90_clobber, nf90_64bit_offset), ncid) - err = nf90_def_dim(ncid, 'str', 255, dim1d) - err = nf90_def_dim(ncid, 'ntiles', 1, dim2d) - err = nf90_def_var(ncid, 'gridfiles', nf90_char, (/dim1d, dim2d/), varid) - err = nf90_enddef(ncid) - err = nf90_put_var(ncid, varid, "ocean_hgrid.nc") - err = nf90_close(ncid) - - err = nf90_create('INPUT/ocean_hgrid.nc', ior(nf90_clobber, nf90_64bit_offset), ncid) - err = nf90_def_dim(ncid, 'nx', 2880, dim1d) - err = nf90_def_dim(ncid, 'ny', 2160, dim2d) - err = nf90_def_dim(ncid, 'nxp', 2881, dim3d) - err = nf90_def_dim(ncid, 'nyp', 2161, dim4d) - err = nf90_def_var(ncid, 'x', nf90_double, (/dim3d, dim4d/), varid) - err = nf90_def_var(ncid, 'y', nf90_double, (/dim3d, dim4d/), varid2) - err = nf90_def_var(ncid, 'area', nf90_double, (/dim1d, dim2d/), varid3) - err = nf90_enddef(ncid) - err = nf90_close(ncid) - - err = nf90_create('INPUT/runoff.daitren.clim.1440x1080.v20180328.nc', ior(nf90_clobber, nf90_64bit_offset), ncid) - err = nf90_def_dim(ncid, 'i', 1440, dim1d) - err = nf90_def_dim(ncid, 'j', 1080, dim2d) - err = nf90_def_dim(ncid, 'time', nf90_unlimited, dim3d) - - err = nf90_def_var(ncid, 'i', nf90_double, (/dim1d/), varid3) - err = nf90_put_att(ncid, varid3, "cartesian_axis", "x") - - err = nf90_def_var(ncid, 'j', nf90_double, (/dim2d/), varid4) - err = nf90_put_att(ncid, varid4, "cartesian_axis", "y") - - err = nf90_def_var(ncid, 'time', nf90_double, (/dim3d/), varid2) - err = nf90_put_att(ncid, varid2, "cartesian_axis", "T") - err = nf90_put_att(ncid, varid2, "calendar", "noleap") - err = nf90_put_att(ncid, varid2, "units", "days since 0001-01-01 00:00:00") - err = nf90_def_var(ncid, 'runoff', nf90_double, (/dim1d, dim2d, dim3d/), varid) - - err = nf90_enddef(ncid) - err = nf90_put_var(ncid, varid, runoff_in) - err = nf90_put_var(ncid, varid2, (/1., 2., 3., 5., 6., 7., 8., 9., 10., 11./)) - err = nf90_close(ncid) - - deallocate(runoff_in) - -endif !< Wait for the root PE to catch up call mpp_sync @@ -131,8 +69,8 @@ program test_data_override_ongrid call set_calendar_type(NOLEAP) -nlon = 1440 -nlat = 1080 +nlon = 360 +nlat = 180 !< Create a domain nlonXnlat with mask call mpp_domains_set_stack_size(17280000) @@ -140,70 +78,60 @@ program test_data_override_ongrid call mpp_define_io_domain(Domain, (/1,1/)) call mpp_get_data_domain(Domain, is, ie, js, je) -print *, nhalox, nhaloy - -!< Set up the data -allocate(runoff(is:ie,js:je)) +select case (test_case) +case (ongrid) + call generate_ongrid_input_file () +case (bilinear) + call generate_bilinear_input_file () +end select -runoff = 999. +call mpp_sync() !< Initiliaze data_override call data_override_init(Ocean_domain_in=Domain, mode=lkind) -!< Run it when time=3 -Time = set_date(1,1,4,0,0,0) -call data_override('OCN','runoff',runoff, Time) -!< Because you are getting the data when time=3, and this is an "ongrid" case, the expected result is just -!! equal to the data at time=3, which is 3. -expected_result = 3._lkind -call compare_data(Domain, runoff, expected_result) - -!< Run it when time=4 -runoff = 999. -Time = set_date(1,1,5,0,0,0) -call data_override('OCN','runoff',runoff, Time) -!< You are getting the data when time=4, the data at time=3 is 3. and at time=5 is 4., so the expected result -!! is the average of the 2 (because this is is an "ongrid" case and there is no horizontal interpolation). -expected_result = (3._lkind + 4._lkind) / 2._lkind -call compare_data(Domain, runoff, expected_result) - -deallocate(runoff) +select case (test_case) +case (ongrid) + call ongrid_test() +case (bilinear) + call bilinear_test() +end select call mpp_exit contains -subroutine compare_data(Domain, actual_result, expected_result) -type(domain2d), intent(in) :: Domain !< Domain with mask table +subroutine compare_data(Domain_in, actual_result, expected_result) +type(domain2d), intent(in) :: Domain_in !< Domain with mask table real(lkind), intent(in) :: expected_result !< Expected result from data_override real(lkind), dimension(:,:), intent(in) :: actual_result !< Result from data_override integer :: xsizec, ysizec !< Size of the compute domain integer :: xsized, ysized !< Size of the data domain integer :: nx, ny !< Size of acual_result -integer :: nhalox, nhaloy !< Size of the halos +integer :: nhx, nhy !< Size of the halos integer :: i, j !< Helper indices !< Data is only expected to be overriden for the compute domain -not at the halos. -call mpp_get_compute_domain(Domain, xsize=xsizec, ysize=ysizec) -call mpp_get_data_domain(Domain, xsize=xsized, ysize=ysized) +call mpp_get_compute_domain(Domain_in, xsize=xsizec, ysize=ysizec) +call mpp_get_data_domain(Domain_in, xsize=xsized, ysize=ysized) !< Note that actual_result has indices at (1:nx,1:ny) not (is:ie,js:je) -nhalox= (xsized-xsizec)/2 -nhaloy = (ysized-ysizec)/2 +nhx= (xsized-xsizec)/2 +nhy = (ysized-ysizec)/2 nx = size(actual_result, 1) ny = size(actual_result, 2) do i = 1, nx do j = 1, ny - if (i <= nhalox .or. i > (nx-nhalox) .or. j <= nhaloy .or. j > (ny-nhaloy)) then + if (i <= nhx .or. i > (nx-nhx) .or. j <= nhy .or. j > (ny-nhy)) then !< This is the result at the halos it should 999. - if (actual_result(i,j) .ne. 999.) then + if (actual_result(i,j) .ne. 999._lkind) then print *, "for i=", i, " and j=", j, " result=", actual_result(i,j) call mpp_error(FATAL, "test_data_override_ongrid: Data was overriden in the halos!!") endif else if (actual_result(i,j) .ne. expected_result) then - print *, "for i=", i, " and j=", j, " result=", actual_result(i,j) + print *, "for i=", i, " and j=", j, " result=", actual_result(i,j), " expected=", expected_result call mpp_error(FATAL, "test_data_override_ongrid: Result is different from expected answer!") endif endif @@ -212,4 +140,318 @@ subroutine compare_data(Domain, actual_result, expected_result) end subroutine +subroutine create_grid_spec_file() + type(FmsNetcdfFile_t) :: fileobj + + if (open_file(fileobj, 'INPUT/grid_spec.nc', 'overwrite')) then + call register_axis(fileobj, 'str', 255) + call register_field(fileobj, 'ocn_mosaic_file', 'char', (/'str'/)) + call write_data(fileobj, 'ocn_mosaic_file', "ocean_mosaic.nc") + call close_file(fileobj) + else + call mpp_error(FATAL, "Error opening the file: 'INPUT/grid_spec.nc' to write") + endif +end subroutine create_grid_spec_file + +subroutine create_ocean_mosaic_file() + type(FmsNetcdfFile_t) :: fileobj + character(len=10) :: dimnames(2) + + dimnames(1) = 'str' + dimnames(2) = 'ntiles' + if (open_file(fileobj, 'INPUT/ocean_mosaic.nc', 'overwrite')) then + call register_axis(fileobj, dimnames(1) , 255) + call register_axis(fileobj, dimnames(2), 1) + call register_field(fileobj, 'gridfiles', 'char', dimnames) + call write_data(fileobj, 'gridfiles', (/"ocean_hgrid.nc"/)) + call close_file(fileobj) + else + call mpp_error(FATAL, "Error opening the file: 'INPUT/ocean_mosaic.nc' to write") + endif +end subroutine create_ocean_mosaic_file + +subroutine create_ocean_hgrid_file() + type(FmsNetcdfFile_t) :: fileobj + + real(lkind), allocatable, dimension(:,:) :: xdata, ydata + integer :: nx, nxp, ny, nyp, i, j + nx = nlon*2 + nxp = nx+1 + ny = nlat*2 + nyp = ny+1 + + allocate(xdata(nxp, nyp)) + xdata(1,:) = 0.0_lkind + do i = 2, nxp + xdata(i,:) = xdata(i-1,:) + 0.5 + enddo + + allocate(ydata(nxp, nyp)) + ydata(:,1) = -90.0_lkind + do i = 2, nyp + ydata(:,i) = ydata(:, i-1) + 0.5 + enddo + + if (open_file(fileobj, 'INPUT/ocean_hgrid.nc', 'overwrite')) then + call register_axis(fileobj, "nx", nx) + call register_axis(fileobj, "ny", ny) + call register_axis(fileobj, "nxp", nxp) + call register_axis(fileobj, "nyp", nyp) + call register_field(fileobj, 'x', 'float', (/'nxp', 'nyp'/)) + call register_field(fileobj, 'y', 'float', (/'nxp', 'nyp'/)) + call register_field(fileobj, 'area', 'float', (/'nx', 'ny'/)) + call write_data(fileobj, "x", xdata) + call write_data(fileobj, "y", ydata) + call close_file(fileobj) + else + call mpp_error(FATAL, "Error opening the file: 'INPUT/ocean_hgrid.nc' to write") + endif +end subroutine create_ocean_hgrid_file + +subroutine create_ongrid_data_file() + type(FmsNetcdfFile_t) :: fileobj + character(len=10) :: dimnames(3) + real(lkind), allocatable, dimension(:,:,:) :: runoff_in + real(lkind), allocatable, dimension(:) :: time_data + integer :: i + + allocate(runoff_in(nlon, nlat, 10)) + allocate(time_data(10)) + do i = 1, 10 + runoff_in(:,:,i) = real(i, lkind) + enddo + time_data = (/1., 2., 3., 5., 6., 7., 8., 9., 10., 11./) + + dimnames(1) = 'i' + dimnames(2) = 'j' + dimnames(3) = 'time' + + if (open_file(fileobj, 'INPUT/runoff.daitren.clim.1440x1080.v20180328.nc', 'overwrite')) then + call register_axis(fileobj, "i", nlon) + call register_axis(fileobj, "j", nlat) + call register_axis(fileobj, "time", unlimited) + + call register_field(fileobj, "i", "float", (/"i"/)) + call register_variable_attribute(fileobj, "i", "cartesian_axis", "x", str_len=1) + + call register_field(fileobj, "j", "float", (/"j"/)) + call register_variable_attribute(fileobj, "j", "cartesian_axis", "y", str_len=1) + + call register_field(fileobj, "time", "float", (/"time"/)) + call register_variable_attribute(fileobj, "time", "cartesian_axis", "T", str_len=1) + call register_variable_attribute(fileobj, "time", "calendar", "noleap", str_len=6) + call register_variable_attribute(fileobj, "time", "units", "days since 0001-01-01 00:00:00", str_len=30) + + call register_field(fileobj, "runoff", "float", dimnames) + call write_data(fileobj, "runoff", runoff_in) + call write_data(fileobj, "time", time_data) + call close_file(fileobj) + else + call mpp_error(FATAL, "Error opening the file: 'INPUT/runoff.daitren.clim.1440x1080.v20180328.nc' to write") + endif + deallocate(runoff_in) +end subroutine create_ongrid_data_file + +subroutine generate_ongrid_input_file() + !< Create some files needed by data_override! + if (mpp_pe() .eq. mpp_root_pe()) then + call create_grid_spec_file() + call create_ocean_mosaic_file() + call create_ocean_hgrid_file() + call create_ongrid_data_file() + endif + call mpp_sync() +end subroutine + +!> @brief Tests ongrid data overrides. +!! In the first case there is no time interpolation +!! In the second case there is time interpolation +subroutine ongrid_test() + real(lkind) :: expected_result !< Expected result from data_override + type(time_type) :: Time !< Time + real(lkind), allocatable, dimension(:,:) :: runoff !< Data to be written + + allocate(runoff(is:ie,js:je)) + + runoff = 999._lkind + !< Run it when time=3 + Time = set_date(1,1,4,0,0,0) + call data_override('OCN','runoff',runoff, Time) + !< Because you are getting the data when time=3, and this is an "ongrid" case, the expected result is just + !! equal to the data at time=3, which is 3. + expected_result = 3._lkind + call compare_data(Domain, runoff, expected_result) + + !< Run it when time=4 + runoff = 999._lkind + Time = set_date(1,1,5,0,0,0) + call data_override('OCN','runoff',runoff, Time) + !< You are getting the data when time=4, the data at time=3 is 3. and at time=5 is 4., so the expected result + !! is the average of the 2 (because this is is an "ongrid" case and there is no horizontal interpolation). + expected_result = (3._lkind + 4._lkind) / 2._lkind + call compare_data(Domain, runoff, expected_result) + + deallocate(runoff) +end subroutine ongrid_test + +!> @brief Creates an input netcdf data file to use for the ongrid data_override test case +!! with either an increasing or decreasing lat, lon grid +subroutine create_bilinear_data_file(increasing_grid) + logical, intent(in) :: increasing_grid !< .true. if increasing a file with an increasing lat/lon + + type(FmsNetcdfFile_t) :: fileobj !< Fms2_io fileobj + character(len=10) :: dimnames(3) !< dimension names for the variable + real(lkind), allocatable :: runoff_in(:,:,:) !< Data to write + real(lkind), allocatable :: time_data(:) !< Time dimension data + real(lkind), allocatable :: lat_data(:) !< Lat dimension data + real(lkind), allocatable :: lon_data(:) !< Lon dimension data + character(len=:), allocatable :: filename !< Name of the file + integer :: factor !< This is used when creating the grid data + !! -1 if the grid is decreasing + !! +1 if the grid is increasing + integer :: i, j, k !< For looping through variables + + allocate(runoff_in(nlon, nlat, 10)) + allocate(time_data(10)) + allocate(lat_data(nlat)) + allocate(lon_data(nlon)) + + do i = 1, nlon + do j = 1, nlat + do k = 1, 10 + runoff_in(i, j, k) = real(i, kind=lkind) * 1000._lkind + real(j, kind=lkind) + real(k, kind=lkind)/100._lkind + enddo + enddo + enddo + if (.not. increasing_grid) then + filename = 'INPUT/bilinear_decreasing.nc' + lon_data(1) = 359.0_lkind + lat_data(1) = 89.0_lkind + factor = -1 + call flip_array(runoff_in) + else + filename = 'INPUT/bilinear_increasing.nc' + lon_data(1) = 0.0_lkind + lat_data(1) = -90.0_lkind + factor = 1 + endif + + do i = 2, nlon + lon_data(i) = real(lon_data(i-1) + 1*factor, lkind) + enddo + + do i = 2, nlat + lat_data(i) =real(lat_data(i-1) + 1*factor, lkind) + enddo + + time_data = (/1., 2., 3., 5., 6., 7., 8., 9., 10., 11./) + + dimnames(1) = 'i' + dimnames(2) = 'j' + dimnames(3) = 'time' + + if (open_file(fileobj, filename, 'overwrite')) then + call register_axis(fileobj, "i", nlon) + call register_axis(fileobj, "j", nlat) + call register_axis(fileobj, "time", unlimited) + + call register_field(fileobj, "i", "float", (/"i"/)) + call register_variable_attribute(fileobj, "i", "cartesian_axis", "x", str_len=1) + + call register_field(fileobj, "j", "float", (/"j"/)) + call register_variable_attribute(fileobj, "j", "cartesian_axis", "y", str_len=1) + + call register_field(fileobj, "time", "float", (/"time"/)) + call register_variable_attribute(fileobj, "time", "cartesian_axis", "T", str_len=1) + call register_variable_attribute(fileobj, "time", "calendar", "noleap", str_len=6) + call register_variable_attribute(fileobj, "time", "units", "days since 0001-01-01 00:00:00", str_len=30) + + call register_field(fileobj, "runoff", "float", dimnames) + call write_data(fileobj, "runoff", runoff_in) + call write_data(fileobj, "i", lon_data) + call write_data(fileobj, "j", lat_data) + call write_data(fileobj, "time", time_data) + call close_file(fileobj) + else + call mpp_error(FATAL, "Error opening the file: 'INPUT/bilinear_increasing.nc' to write") + endif + deallocate(runoff_in) +end subroutine create_bilinear_data_file + +!> @brief Generates the input for the bilinear data_override test_case +subroutine generate_bilinear_input_file() + if (mpp_pe() .eq. mpp_root_pe()) then + call create_grid_spec_file () + call create_ocean_mosaic_file() + call create_ocean_hgrid_file() + call create_bilinear_data_file(.true.) + call create_bilinear_data_file(.false.) + endif + call mpp_sync() +end subroutine generate_bilinear_input_file + +!> @brief Tests bilinear data_override with and increasing and decreasing grid case +!! and comares the output betweeen the cases to ensure it is correct +subroutine bilinear_test() + type(time_type) :: Time !< Time + real(lkind), allocatable, dimension(:,:) :: runoff_decreasing !< Data to be written + real(lkind), allocatable, dimension(:,:) :: runoff_increasing !< Data to be written + + integer :: i, j, k + logical :: success + + allocate(runoff_decreasing(is:ie,js:je)) + allocate(runoff_increasing(is:ie,js:je)) + + runoff_decreasing = 999_lkind + runoff_increasing = 999_lkind + Time = set_date(1,1,4,0,0,0) + call data_override('OCN','runoff_increasing',runoff_decreasing, Time, override=success) + if (.not. success) call mpp_error(FATAL, "Data override failed") + call data_override('OCN','runoff_increasing',runoff_increasing, Time, override=success) + if (.not. success) call mpp_error(FATAL, "Data override failed") + + do i = 1, size(runoff_decreasing, 1) + do j = 1, size(runoff_decreasing, 2) + if (runoff_decreasing(i,j) .ne. runoff_increasing(i,j)) & + call mpp_error(FATAL, "The data is not the same: "// & + string(i)//","//string(j)//":"// & + string(runoff_decreasing(i,j))//" vs "//string(runoff_increasing(i,j))) + enddo + enddo + deallocate(runoff_decreasing, runoff_increasing) +end subroutine bilinear_test + +!> @brief Flips an array along the x and y axes +subroutine flip_array(array) + real(lkind), intent(inout) :: array(:,:,:) !< Array to flip + + integer :: nx !< size of the x dimension + integer :: ny !< size of the y dimension + integer :: nt !< size of the t dimension + + integer :: i, j, k + + nx = size(array, 1) + ny = size(array, 2) + nt = size(array, 3) + + do k = 1, nt + do j = 1, ny + do i = 1, nx / 2 + ! Swap elements along the x-axis + array(i, j, k) = array(nx - i + 1, j, k) + array(nx - i + 1, j, k) = array(i, j, k) + end do + end do + + do i = 1, nx + do j = 1, ny / 2 + ! Swap elements along the y-axis + array(i, j, k) = array(i, ny - j + 1, k) + array(i, ny - j + 1, k) = array(i, j, k) + end do + end do + end do +end subroutine flip_array end program test_data_override_ongrid From 0f2129a0b5933ad70349a29905f79132c23f4b92 Mon Sep 17 00:00:00 2001 From: Uriel Ramirez Date: Mon, 11 Dec 2023 14:42:27 -0500 Subject: [PATCH 10/16] attempt to add a data_override unit test to test increasing and decreasing data sets --- .../include/horiz_interp_bilinear.inc | 13 ++- .../test_data_override_ongrid.F90 | 90 ++++++++----------- 2 files changed, 40 insertions(+), 63 deletions(-) diff --git a/horiz_interp/include/horiz_interp_bilinear.inc b/horiz_interp/include/horiz_interp_bilinear.inc index 5bac5e4b01..9e352d9c31 100644 --- a/horiz_interp/include/horiz_interp_bilinear.inc +++ b/horiz_interp/include/horiz_interp_bilinear.inc @@ -105,19 +105,18 @@ is = nearest_index(lon, lon_in ) if (decreasing_lon) then - ! Lon_in is decreasing + ! Lon_in is increasing ! This is so that is is the lower bound. - ! For example, if the array is [10 20 30 40 50] and lon is 49, `is` is going to be 5 from `nearest_index` + ! For example, if the array is [50 40 30 20 10] and lon is 11, `is` is going to be 5 from `nearest_index` ! but it needs to be 4 so that it can use the data at lon_in(4) and lon_in(5) - if( lon_in(is) .gt. lon ) is = max(is-1,1) + if( lon_in(is) .lt. lon ) is = max(is-1,1) else ! Lon_in is increasing ! This is so that is is the lower bound. - ! For example, if the array is [50 40 30 20 10] and lon is 11, `is` is going to be 5 from `nearest_index` + ! For example, if the array is [10 20 30 40 50] and lon is 49, `is` is going to be 5 from `nearest_index` ! but it needs to be 4 so that it can use the data at lon_in(4) and lon_in(5) - if( lon_in(is) .lt. lon ) is = max(is-1,1) + if( lon_in(is) .gt. lon ) is = max(is-1,1) endif - if( lon_in(is) .eq. lon .and. is .eq. nlon_in) is = max(is - 1,1) ie = min(is+1,nlon_in) if(lon_in(is) .ne. lon_in(ie) .and. (decreasing_lon .or. lon_in(is) .le. lon)) then @@ -137,7 +136,6 @@ wte = 1.0_kindl - wtw js = nearest_index(lat, lat_in ) - if (decreasing_lat) then ! Lat_in is decreasing if( lat_in(js) .lt. lat ) js = max(js - 1, 1) @@ -145,7 +143,6 @@ ! Lat_in is increasing if( lat_in(js) .gt. lat ) js = max(js - 1, 1) endif - if( lat_in(js) .eq. lat .and. js .eq. nlat_in) js = max(js - 1, 1) je = min(js + 1, nlat_in) diff --git a/test_fms/data_override/test_data_override_ongrid.F90 b/test_fms/data_override/test_data_override_ongrid.F90 index 99140865af..7704898396 100644 --- a/test_fms/data_override/test_data_override_ongrid.F90 +++ b/test_fms/data_override/test_data_override_ongrid.F90 @@ -311,36 +311,48 @@ subroutine create_bilinear_data_file(increasing_grid) !! +1 if the grid is increasing integer :: i, j, k !< For looping through variables - allocate(runoff_in(nlon, nlat, 10)) + integer :: nlon_data + integer :: nlat_data + + nlon_data = nlon + 1 + nlat_data = nlat - 1 + allocate(runoff_in(nlon_data, nlat_data, 10)) allocate(time_data(10)) - allocate(lat_data(nlat)) - allocate(lon_data(nlon)) + allocate(lat_data(nlat_data)) + allocate(lon_data(nlon_data)) - do i = 1, nlon - do j = 1, nlat - do k = 1, 10 - runoff_in(i, j, k) = real(i, kind=lkind) * 1000._lkind + real(j, kind=lkind) + real(k, kind=lkind)/100._lkind - enddo - enddo - enddo if (.not. increasing_grid) then filename = 'INPUT/bilinear_decreasing.nc' - lon_data(1) = 359.0_lkind + lon_data(1) = 360.0_lkind lat_data(1) = 89.0_lkind factor = -1 - call flip_array(runoff_in) + do i = 1, nlon_data + do j = 1, nlat_data + do k = 1, 10 + runoff_in(i, j, k) = real(362-i, kind=lkind) * 1000._lkind + real(180-j, kind=lkind) + real(k, kind=lkind)/100._lkind + enddo + enddo + enddo else filename = 'INPUT/bilinear_increasing.nc' lon_data(1) = 0.0_lkind - lat_data(1) = -90.0_lkind + lat_data(1) = -89.0_lkind factor = 1 + + do i = 1, nlon_data + do j = 1, nlat_data + do k = 1, 10 + runoff_in(i, j, k) = real(i, kind=lkind) * 1000._lkind + real(j, kind=lkind) + real(k, kind=lkind)/100._lkind + enddo + enddo + enddo endif - do i = 2, nlon + do i = 2, nlon_data lon_data(i) = real(lon_data(i-1) + 1*factor, lkind) enddo - do i = 2, nlat + do i = 2, nlat_data lat_data(i) =real(lat_data(i-1) + 1*factor, lkind) enddo @@ -351,8 +363,8 @@ subroutine create_bilinear_data_file(increasing_grid) dimnames(3) = 'time' if (open_file(fileobj, filename, 'overwrite')) then - call register_axis(fileobj, "i", nlon) - call register_axis(fileobj, "j", nlat) + call register_axis(fileobj, "i", nlon_data) + call register_axis(fileobj, "j", nlat_data) call register_axis(fileobj, "time", unlimited) call register_field(fileobj, "i", "float", (/"i"/)) @@ -406,52 +418,20 @@ subroutine bilinear_test() runoff_decreasing = 999_lkind runoff_increasing = 999_lkind Time = set_date(1,1,4,0,0,0) - call data_override('OCN','runoff_increasing',runoff_decreasing, Time, override=success) - if (.not. success) call mpp_error(FATAL, "Data override failed") call data_override('OCN','runoff_increasing',runoff_increasing, Time, override=success) if (.not. success) call mpp_error(FATAL, "Data override failed") + call data_override('OCN','runoff_decreasing',runoff_decreasing, Time, override=success) + if (.not. success) call mpp_error(FATAL, "Data override failed") - do i = 1, size(runoff_decreasing, 1) - do j = 1, size(runoff_decreasing, 2) - if (runoff_decreasing(i,j) .ne. runoff_increasing(i,j)) & + do i = is, ie + do j = js, je + if (abs(runoff_decreasing(i,j) - runoff_increasing(i,j)) .gt. 1) then call mpp_error(FATAL, "The data is not the same: "// & string(i)//","//string(j)//":"// & string(runoff_decreasing(i,j))//" vs "//string(runoff_increasing(i,j))) + endif enddo enddo deallocate(runoff_decreasing, runoff_increasing) end subroutine bilinear_test - -!> @brief Flips an array along the x and y axes -subroutine flip_array(array) - real(lkind), intent(inout) :: array(:,:,:) !< Array to flip - - integer :: nx !< size of the x dimension - integer :: ny !< size of the y dimension - integer :: nt !< size of the t dimension - - integer :: i, j, k - - nx = size(array, 1) - ny = size(array, 2) - nt = size(array, 3) - - do k = 1, nt - do j = 1, ny - do i = 1, nx / 2 - ! Swap elements along the x-axis - array(i, j, k) = array(nx - i + 1, j, k) - array(nx - i + 1, j, k) = array(i, j, k) - end do - end do - - do i = 1, nx - do j = 1, ny / 2 - ! Swap elements along the y-axis - array(i, j, k) = array(i, ny - j + 1, k) - array(i, ny - j + 1, k) = array(i, j, k) - end do - end do - end do -end subroutine flip_array end program test_data_override_ongrid From fd63642f072a8ae748f02b51fd578249aed6fc4d Mon Sep 17 00:00:00 2001 From: Uriel Ramirez Date: Mon, 11 Dec 2023 14:52:37 -0500 Subject: [PATCH 11/16] Fix white space and line count limit --- .../data_override/test_data_override_ongrid.F90 | 14 ++++++++------ 1 file changed, 8 insertions(+), 6 deletions(-) diff --git a/test_fms/data_override/test_data_override_ongrid.F90 b/test_fms/data_override/test_data_override_ongrid.F90 index 7704898396..53f3903f59 100644 --- a/test_fms/data_override/test_data_override_ongrid.F90 +++ b/test_fms/data_override/test_data_override_ongrid.F90 @@ -233,7 +233,7 @@ subroutine create_ongrid_data_file() call register_field(fileobj, "i", "float", (/"i"/)) call register_variable_attribute(fileobj, "i", "cartesian_axis", "x", str_len=1) - + call register_field(fileobj, "j", "float", (/"j"/)) call register_variable_attribute(fileobj, "j", "cartesian_axis", "y", str_len=1) @@ -263,7 +263,7 @@ subroutine generate_ongrid_input_file() call mpp_sync() end subroutine -!> @brief Tests ongrid data overrides. +!> @brief Tests ongrid data overrides. !! In the first case there is no time interpolation !! In the second case there is time interpolation subroutine ongrid_test() @@ -329,7 +329,8 @@ subroutine create_bilinear_data_file(increasing_grid) do i = 1, nlon_data do j = 1, nlat_data do k = 1, 10 - runoff_in(i, j, k) = real(362-i, kind=lkind) * 1000._lkind + real(180-j, kind=lkind) + real(k, kind=lkind)/100._lkind + runoff_in(i, j, k) = real(362-i, kind=lkind) * 1000._lkind + & + real(180-j, kind=lkind) + real(k, kind=lkind)/100._lkind enddo enddo enddo @@ -342,7 +343,8 @@ subroutine create_bilinear_data_file(increasing_grid) do i = 1, nlon_data do j = 1, nlat_data do k = 1, 10 - runoff_in(i, j, k) = real(i, kind=lkind) * 1000._lkind + real(j, kind=lkind) + real(k, kind=lkind)/100._lkind + runoff_in(i, j, k) = real(i, kind=lkind) * 1000._lkind + real(j, kind=lkind) + & + real(k, kind=lkind)/100._lkind enddo enddo enddo @@ -369,7 +371,7 @@ subroutine create_bilinear_data_file(increasing_grid) call register_field(fileobj, "i", "float", (/"i"/)) call register_variable_attribute(fileobj, "i", "cartesian_axis", "x", str_len=1) - + call register_field(fileobj, "j", "float", (/"j"/)) call register_variable_attribute(fileobj, "j", "cartesian_axis", "y", str_len=1) @@ -422,7 +424,7 @@ subroutine bilinear_test() if (.not. success) call mpp_error(FATAL, "Data override failed") call data_override('OCN','runoff_decreasing',runoff_decreasing, Time, override=success) if (.not. success) call mpp_error(FATAL, "Data override failed") - + do i = is, ie do j = js, je if (abs(runoff_decreasing(i,j) - runoff_increasing(i,j)) .gt. 1) then From a15c4f0e0f30f27eafaa5bf36a5b3156511cff7c Mon Sep 17 00:00:00 2001 From: Uriel Ramirez Date: Mon, 11 Dec 2023 16:01:28 -0500 Subject: [PATCH 12/16] Clean up the INPUT after done with it --- test_fms/data_override/test_data_override2_mono.sh | 1 + 1 file changed, 1 insertion(+) diff --git a/test_fms/data_override/test_data_override2_mono.sh b/test_fms/data_override/test_data_override2_mono.sh index 63ee989d01..7f17f8711f 100755 --- a/test_fms/data_override/test_data_override2_mono.sh +++ b/test_fms/data_override/test_data_override2_mono.sh @@ -42,4 +42,5 @@ test_expect_success "test_data_override with monotonically increasing and decrea mpirun -n 6 ../test_data_override_ongrid_r4 ' +rm -rf INPUT test_done From d7345dd53e7eb7bcc37e681aa66d97d35d90549a Mon Sep 17 00:00:00 2001 From: Uriel Ramirez Date: Mon, 11 Dec 2023 16:38:29 -0500 Subject: [PATCH 13/16] Move the ongrid test to its own script --- test_fms/data_override/Makefile.am | 4 +- test_fms/data_override/test_data_override2.sh | 54 +----------- .../test_data_override2_ongrid.sh | 82 +++++++++++++++++++ .../test_data_override_ongrid.F90 | 5 +- 4 files changed, 90 insertions(+), 55 deletions(-) create mode 100755 test_fms/data_override/test_data_override2_ongrid.sh diff --git a/test_fms/data_override/Makefile.am b/test_fms/data_override/Makefile.am index bf94708437..d71512a364 100644 --- a/test_fms/data_override/Makefile.am +++ b/test_fms/data_override/Makefile.am @@ -72,10 +72,10 @@ TESTS_ENVIRONMENT= test_input_path="@TEST_INPUT_PATH@" \ # Run the test program. -TESTS = test_data_override2.sh test_data_override_init.sh test_data_override2_mono.sh +TESTS = test_data_override2.sh test_data_override_init.sh test_data_override2_mono.sh test_data_override2_ongrid.sh # Include these files with the distribution. -EXTRA_DIST = test_data_override2.sh test_data_override_init.sh test_data_override2_mono.sh +EXTRA_DIST = test_data_override2.sh test_data_override_init.sh test_data_override2_mono.sh test_data_override2_ongrid.sh # Clean up CLEANFILES = input.nml *.nc* *.out diag_table data_table data_table.yaml INPUT/* *.dpi *.spi *.dyn *.spl *-files/* diff --git a/test_fms/data_override/test_data_override2.sh b/test_fms/data_override/test_data_override2.sh index ffa42c05c6..48940debf7 100755 --- a/test_fms/data_override/test_data_override2.sh +++ b/test_fms/data_override/test_data_override2.sh @@ -27,62 +27,14 @@ output_dir rm -rf data_table data_table.yaml input.nml input_base.nml -if [ ! -z $parser_skip ]; then - cat <<_EOF > input_base.nml -&data_override_nml -use_data_table_yaml=.False. -/ - -&test_data_override_ongrid_nml - nhalox=halo_size - nhaloy=halo_size -/ -_EOF - printf '"OCN", "runoff", "runoff", "./INPUT/runoff.daitren.clim.1440x1080.v20180328.nc", "none" , 1.0' | cat > data_table -else -cat <<_EOF > input_base.nml -&data_override_nml -use_data_table_yaml=.True. -/ - -&test_data_override_ongrid_nml - nhalox=halo_size - nhaloy=halo_size -/ -_EOF -cat <<_EOF > data_table.yaml -data_table: - - gridname : OCN - fieldname_code : runoff - fieldname_file : runoff - file_name : INPUT/runoff.daitren.clim.1440x1080.v20180328.nc - interpol_method : none - factor : 1.0 -_EOF -fi - -[ ! -d "INPUT" ] && mkdir -p "INPUT" -for KIND in r4 r8 -do -sed 's/halo_size/2/g' input_base.nml > input.nml -test_expect_success "data_override on grid with 2 halos in x and y (${KIND})" ' - mpirun -n 6 ../test_data_override_ongrid_${KIND} -' -sed 's/halo_size/0/g' input_base.nml > input.nml -test_expect_success "data_override on grid with 2 halos in x and y (${KIND})" ' - mpirun -n 6 ../test_data_override_ongrid_${KIND} -' - -test_expect_success "data_override get_grid_v1 (${KIND})" ' - mpirun -n 1 ../test_get_grid_v1_${KIND} -' -done - for KIND in r4 r8 do # Run tests with input if enabled # skips if built with yaml parser(tests older behavior) if test ! -z "$test_input_path" && test ! -z "$parser_skip" ; then + cat <<_EOF > input.nml +_EOF + cp -r $test_input_path/data_override/INPUT . cat <<_EOF > diag_table test_data_override diff --git a/test_fms/data_override/test_data_override2_ongrid.sh b/test_fms/data_override/test_data_override2_ongrid.sh new file mode 100755 index 0000000000..faf89302eb --- /dev/null +++ b/test_fms/data_override/test_data_override2_ongrid.sh @@ -0,0 +1,82 @@ +#!/bin/sh + +#*********************************************************************** +#* GNU Lesser General Public License +#* +#* This file is part of the GFDL Flexible Modeling System (FMS). +#* +#* FMS is free software: you can redistribute it and/or modify it under +#* the terms of the GNU Lesser General Public License as published by +#* the Free Software Foundation, either version 3 of the License, or (at +#* your option) any later version. +#* +#* FMS is distributed in the hope that it will be useful, but WITHOUT +#* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +#* FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +#* for more details. +#* +#* You should have received a copy of the GNU Lesser General Public +#* License along with FMS. If not, see . +#*********************************************************************** +# +# Copyright (c) 2019-2021 Ed Hartnett, Uriel Ramirez, Seth Underwood + +# Set common test settings. +. ../test-lib.sh + +output_dir +rm -rf data_table data_table.yaml input.nml input_base.nml + +if [ ! -z $parser_skip ]; then + cat <<_EOF > input_base.nml +&data_override_nml +use_data_table_yaml=.False. +/ + +&test_data_override_ongrid_nml + nhalox=halo_size + nhaloy=halo_size +/ +_EOF + printf '"OCN", "runoff", "runoff", "./INPUT/runoff.daitren.clim.1440x1080.v20180328.nc", "none" , 1.0' | cat > data_table +else +cat <<_EOF > input_base.nml +&data_override_nml +use_data_table_yaml=.True. +/ + +&test_data_override_ongrid_nml + nhalox=halo_size + nhaloy=halo_size +/ +_EOF +cat <<_EOF > data_table.yaml +data_table: + - gridname : OCN + fieldname_code : runoff + fieldname_file : runoff + file_name : INPUT/runoff.daitren.clim.1440x1080.v20180328.nc + interpol_method : none + factor : 1.0 +_EOF +fi + +[ ! -d "INPUT" ] && mkdir -p "INPUT" +for KIND in r4 r8 +do +sed 's/halo_size/2/g' input_base.nml > input.nml +test_expect_success "data_override on grid with 2 halos in x and y (${KIND})" ' + mpirun -n 6 ../test_data_override_ongrid_${KIND} +' +sed 's/halo_size/0/g' input_base.nml > input.nml +test_expect_success "data_override on grid with 0 halos in x and y (${KIND})" ' + mpirun -n 6 ../test_data_override_ongrid_${KIND} +' + +test_expect_success "data_override get_grid_v1 (${KIND})" ' + mpirun -n 1 ../test_get_grid_v1_${KIND} +' +done +rm -rf INPUT *.nc # remove any leftover files to reduce size + +test_done \ No newline at end of file diff --git a/test_fms/data_override/test_data_override_ongrid.F90 b/test_fms/data_override/test_data_override_ongrid.F90 index 53f3903f59..ac97f03a11 100644 --- a/test_fms/data_override/test_data_override_ongrid.F90 +++ b/test_fms/data_override/test_data_override_ongrid.F90 @@ -20,13 +20,13 @@ program test_data_override_ongrid !> @brief This programs tests data_override ability to override data for an -!! on grid case +!! on grid case and when using bilinear interpolation use platform_mod, only: r4_kind, r8_kind use mpp_domains_mod, only: mpp_define_domains, mpp_define_io_domain, mpp_get_data_domain, & mpp_domains_set_stack_size, mpp_get_compute_domain, domain2d use mpp_mod, only: mpp_init, mpp_exit, mpp_pe, mpp_root_pe, mpp_error, FATAL, & - input_nml_file, mpp_sync + input_nml_file, mpp_sync, NOTE use data_override_mod, only: data_override_init, data_override use fms2_io_mod use time_manager_mod, only: set_calendar_type, time_type, set_date, NOLEAP @@ -86,6 +86,7 @@ program test_data_override_ongrid end select call mpp_sync() +call mpp_error(NOTE, "Finished creating INPUT Files") !< Initiliaze data_override call data_override_init(Ocean_domain_in=Domain, mode=lkind) From 60bab6bd0bb554c56bc5797a6c75946104caa98b Mon Sep 17 00:00:00 2001 From: Uriel Ramirez Date: Mon, 11 Dec 2023 16:48:09 -0500 Subject: [PATCH 14/16] Clean up after running each test --- test_fms/data_override/test_data_override2_ongrid.sh | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/test_fms/data_override/test_data_override2_ongrid.sh b/test_fms/data_override/test_data_override2_ongrid.sh index faf89302eb..adeff4ee5c 100755 --- a/test_fms/data_override/test_data_override2_ongrid.sh +++ b/test_fms/data_override/test_data_override2_ongrid.sh @@ -64,15 +64,19 @@ fi [ ! -d "INPUT" ] && mkdir -p "INPUT" for KIND in r4 r8 do +rm -rf INPUT/. sed 's/halo_size/2/g' input_base.nml > input.nml test_expect_success "data_override on grid with 2 halos in x and y (${KIND})" ' mpirun -n 6 ../test_data_override_ongrid_${KIND} ' + +rm -rf INPUT/. sed 's/halo_size/0/g' input_base.nml > input.nml test_expect_success "data_override on grid with 0 halos in x and y (${KIND})" ' mpirun -n 6 ../test_data_override_ongrid_${KIND} ' +rm -rf INPUT/. test_expect_success "data_override get_grid_v1 (${KIND})" ' mpirun -n 1 ../test_get_grid_v1_${KIND} ' From 4f6c6c78f0187e7cca4cbf6bbfdde76ca7d59559 Mon Sep 17 00:00:00 2001 From: Uriel Ramirez Date: Mon, 11 Dec 2023 16:55:57 -0500 Subject: [PATCH 15/16] replace a dot with a star --- test_fms/data_override/test_data_override2_ongrid.sh | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/test_fms/data_override/test_data_override2_ongrid.sh b/test_fms/data_override/test_data_override2_ongrid.sh index adeff4ee5c..2e1d7a1b03 100755 --- a/test_fms/data_override/test_data_override2_ongrid.sh +++ b/test_fms/data_override/test_data_override2_ongrid.sh @@ -64,19 +64,19 @@ fi [ ! -d "INPUT" ] && mkdir -p "INPUT" for KIND in r4 r8 do -rm -rf INPUT/. +rm -rf INPUT/* sed 's/halo_size/2/g' input_base.nml > input.nml test_expect_success "data_override on grid with 2 halos in x and y (${KIND})" ' mpirun -n 6 ../test_data_override_ongrid_${KIND} ' -rm -rf INPUT/. +rm -rf INPUT/* sed 's/halo_size/0/g' input_base.nml > input.nml test_expect_success "data_override on grid with 0 halos in x and y (${KIND})" ' mpirun -n 6 ../test_data_override_ongrid_${KIND} ' -rm -rf INPUT/. +rm -rf INPUT/* test_expect_success "data_override get_grid_v1 (${KIND})" ' mpirun -n 1 ../test_get_grid_v1_${KIND} ' From 3ab69fabca6a359e8fba0094eac82573a5f1b51b Mon Sep 17 00:00:00 2001 From: Uriel Ramirez Date: Mon, 11 Dec 2023 17:21:22 -0500 Subject: [PATCH 16/16] Add yaml r4/r8 tests --- .../data_override/test_data_override2_mono.sh | 48 +++++++++++++++++-- 1 file changed, 45 insertions(+), 3 deletions(-) diff --git a/test_fms/data_override/test_data_override2_mono.sh b/test_fms/data_override/test_data_override2_mono.sh index 7f17f8711f..cf47a152f9 100755 --- a/test_fms/data_override/test_data_override2_mono.sh +++ b/test_fms/data_override/test_data_override2_mono.sh @@ -38,9 +38,51 @@ cat <<_EOF > data_table "OCN", "runoff_decreasing", "runoff", "./INPUT/bilinear_decreasing.nc", "bilinear" , 1.0 _EOF -test_expect_success "test_data_override with monotonically increasing and decreasing data sets" ' - mpirun -n 6 ../test_data_override_ongrid_r4 - ' +for KIND in r4 r8 +do + rm -rf INPUT/* + test_expect_success "test_data_override with monotonically increasing and decreasing data sets (${KIND})" ' + mpirun -n 6 ../test_data_override_ongrid_${KIND} + ' +done + +rm -rf data_table + +cat <<_EOF > input.nml +&test_data_override_ongrid_nml + test_case = 2 +/ +&data_override_nml + use_data_table_yaml = .True. +/ +_EOF + +cat <<_EOF > data_table.yaml +data_table: +- gridname: OCN + fieldname_code: runoff_increasing + fieldname_file: runoff + file_name: ./INPUT/bilinear_increasing.nc + interpol_method: bilinear + factor: 1.0 +- gridname: OCN + fieldname_code: runoff_decreasing + fieldname_file: runoff + file_name: ./INPUT/bilinear_decreasing.nc + interpol_method: bilinear + factor: 1.0 +_EOF + +#Repeat the test with yaml if needed +if [ -z $parser_skip ]; then + for KIND in r4 r8 + do + rm -rf INPUT/* + test_expect_success "test_data_override with monotonically increasing and decreasing data sets -yaml (${KIND})" ' + mpirun -n 6 ../test_data_override_ongrid_${KIND} + ' + done +fi rm -rf INPUT test_done