From caef59e82294704df8ae3daff183db9e92b3b457 Mon Sep 17 00:00:00 2001 From: Ryan Mulhall <35538242+rem1776@users.noreply.github.com> Date: Wed, 9 Oct 2024 11:28:16 -0400 Subject: [PATCH] fix: add allocation checks to `horiz_interp_type_eq` (#1584) --- horiz_interp/horiz_interp_type.F90 | 161 ++++++++++++++------ test_fms/horiz_interp/test_horiz_interp.F90 | 52 +++++-- 2 files changed, 157 insertions(+), 56 deletions(-) diff --git a/horiz_interp/horiz_interp_type.F90 b/horiz_interp/horiz_interp_type.F90 index e87870698..a2bc90a82 100644 --- a/horiz_interp/horiz_interp_type.F90 +++ b/horiz_interp/horiz_interp_type.F90 @@ -164,58 +164,131 @@ subroutine horiz_interp_type_eq(horiz_interp_out, horiz_interp_in) call mpp_error(FATAL,'horiz_interp_type_eq: horiz_interp_type variable on right hand side is unassigned') endif - horiz_interp_out%ilon = horiz_interp_in%ilon - horiz_interp_out%jlat = horiz_interp_in%jlat - horiz_interp_out%i_lon = horiz_interp_in%i_lon - horiz_interp_out%j_lat = horiz_interp_in%j_lat - horiz_interp_out%found_neighbors = horiz_interp_in%found_neighbors - horiz_interp_out%num_found = horiz_interp_in%num_found - horiz_interp_out%nlon_src = horiz_interp_in%nlon_src - horiz_interp_out%nlat_src = horiz_interp_in%nlat_src - horiz_interp_out%nlon_dst = horiz_interp_in%nlon_dst - horiz_interp_out%nlat_dst = horiz_interp_in%nlat_dst + if( allocated(horiz_interp_in%ilon )) & + horiz_interp_out%ilon = horiz_interp_in%ilon + + if( allocated(horiz_interp_in%jlat )) & + horiz_interp_out%jlat = horiz_interp_in%jlat + + if( allocated(horiz_interp_in%i_lon )) & + horiz_interp_out%i_lon = horiz_interp_in%i_lon + + if( allocated(horiz_interp_in%j_lat )) & + horiz_interp_out%j_lat = horiz_interp_in%j_lat + + if( allocated(horiz_interp_in%found_neighbors )) & + horiz_interp_out%found_neighbors = horiz_interp_in%found_neighbors + + if( allocated(horiz_interp_in%num_found )) & + horiz_interp_out%num_found = horiz_interp_in%num_found + + if( allocated(horiz_interp_in%i_src )) & + horiz_interp_out%i_src = horiz_interp_in%i_src + + if( allocated(horiz_interp_in%j_src )) & + horiz_interp_out%j_src = horiz_interp_in%j_src + + if( allocated(horiz_interp_in%i_dst )) & + horiz_interp_out%i_dst = horiz_interp_in%i_dst + + if( allocated(horiz_interp_in%j_dst )) & + horiz_interp_out%j_dst = horiz_interp_in%j_dst + + horiz_interp_out%nlon_src = horiz_interp_in%nlon_src + horiz_interp_out%nlat_src = horiz_interp_in%nlat_src + horiz_interp_out%nlon_dst = horiz_interp_in%nlon_dst + horiz_interp_out%nlat_dst = horiz_interp_in%nlat_dst horiz_interp_out%interp_method = horiz_interp_in%interp_method horiz_interp_out%I_am_initialized = .true. - horiz_interp_out%i_src = horiz_interp_in%i_src - horiz_interp_out%j_src = horiz_interp_in%j_src - horiz_interp_out%i_dst = horiz_interp_in%i_dst - horiz_interp_out%j_dst = horiz_interp_in%j_dst if(horiz_interp_in%horizInterpReals8_type%is_allocated) then - horiz_interp_out%horizInterpReals8_type%faci = horiz_interp_in%horizInterpReals8_type%faci - horiz_interp_out%horizInterpReals8_type%facj = horiz_interp_in%horizInterpReals8_type%facj - horiz_interp_out%horizInterpReals8_type%area_src = horiz_interp_in%horizInterpReals8_type%area_src - horiz_interp_out%horizInterpReals8_type%area_dst = horiz_interp_in%horizInterpReals8_type%area_dst - horiz_interp_out%horizInterpReals8_type%wti = horiz_interp_in%horizInterpReals8_type%wti - horiz_interp_out%horizInterpReals8_type%wtj = horiz_interp_in%horizInterpReals8_type%wtj - horiz_interp_out%horizInterpReals8_type%src_dist = horiz_interp_in%horizInterpReals8_type%src_dist - horiz_interp_out%horizInterpReals8_type%rat_x = horiz_interp_in%horizInterpReals8_type%rat_x - horiz_interp_out%horizInterpReals8_type%rat_y = horiz_interp_in%horizInterpReals8_type%rat_y - horiz_interp_out%horizInterpReals8_type%lon_in = horiz_interp_in%horizInterpReals8_type%lon_in - horiz_interp_out%horizInterpReals8_type%lat_in = horiz_interp_in%horizInterpReals8_type%lat_in - horiz_interp_out%horizInterpReals8_type%area_frac_dst = horiz_interp_in%horizInterpReals8_type%area_frac_dst - horiz_interp_out%horizInterpReals8_type%max_src_dist = horiz_interp_in%horizInterpReals8_type%max_src_dist - horiz_interp_out%horizInterpReals8_type%is_allocated = .true. + + if( allocated(horiz_interp_in%horizInterpReals8_type%faci)) & + horiz_interp_out%horizInterpReals8_type%faci = horiz_interp_in%horizInterpReals8_type%faci + + if( allocated( horiz_interp_in%horizInterpReals8_type%facj)) & + horiz_interp_out%horizInterpReals8_type%facj = horiz_interp_in%horizInterpReals8_type%facj + + if( allocated( horiz_interp_in%horizInterpReals8_type%area_src)) & + horiz_interp_out%horizInterpReals8_type%area_src = horiz_interp_in%horizInterpReals8_type%area_src + + if( allocated( horiz_interp_in%horizInterpReals8_type%area_dst)) & + horiz_interp_out%horizInterpReals8_type%area_dst = horiz_interp_in%horizInterpReals8_type%area_dst + + if( allocated( horiz_interp_in%horizInterpReals8_type%wti)) & + horiz_interp_out%horizInterpReals8_type%wti = horiz_interp_in%horizInterpReals8_type%wti + + if( allocated( horiz_interp_in%horizInterpReals8_type%wtj)) & + horiz_interp_out%horizInterpReals8_type%wtj = horiz_interp_in%horizInterpReals8_type%wtj + + if( allocated( horiz_interp_in%horizInterpReals8_type%src_dist)) & + horiz_interp_out%horizInterpReals8_type%src_dist = horiz_interp_in%horizInterpReals8_type%src_dist + + if( allocated( horiz_interp_in%horizInterpReals8_type%rat_x)) & + horiz_interp_out%horizInterpReals8_type%rat_x = horiz_interp_in%horizInterpReals8_type%rat_x + + if( allocated( horiz_interp_in%horizInterpReals8_type%rat_y)) & + horiz_interp_out%horizInterpReals8_type%rat_y = horiz_interp_in%horizInterpReals8_type%rat_y + + if( allocated( horiz_interp_in%horizInterpReals8_type%lon_in)) & + horiz_interp_out%horizInterpReals8_type%lon_in = horiz_interp_in%horizInterpReals8_type%lon_in + + if( allocated( horiz_interp_in%horizInterpReals8_type%lat_in)) & + horiz_interp_out%horizInterpReals8_type%lat_in = horiz_interp_in%horizInterpReals8_type%lat_in + + if( allocated( horiz_interp_in%horizInterpReals8_type%area_frac_dst)) & + horiz_interp_out%horizInterpReals8_type%area_frac_dst = horiz_interp_in%horizInterpReals8_type%area_frac_dst + + horiz_interp_out%horizInterpReals8_type%max_src_dist = horiz_interp_in%horizInterpReals8_type%max_src_dist + + horiz_interp_out%horizInterpReals8_type%is_allocated = .true. ! this was left out previous to mixed mode - horiz_interp_out%horizInterpReals8_type%mask_in = horiz_interp_in%horizInterpReals8_type%mask_in + if( allocated(horiz_interp_in%horizInterpReals8_type%mask_in)) & + horiz_interp_out%horizInterpReals8_type%mask_in = horiz_interp_in%horizInterpReals8_type%mask_in else if (horiz_interp_in%horizInterpReals4_type%is_allocated) then - horiz_interp_out%horizInterpReals4_type%faci = horiz_interp_in%horizInterpReals4_type%faci - horiz_interp_out%horizInterpReals4_type%facj = horiz_interp_in%horizInterpReals4_type%facj - horiz_interp_out%horizInterpReals4_type%area_src = horiz_interp_in%horizInterpReals4_type%area_src - horiz_interp_out%horizInterpReals4_type%area_dst = horiz_interp_in%horizInterpReals4_type%area_dst - horiz_interp_out%horizInterpReals4_type%wti = horiz_interp_in%horizInterpReals4_type%wti - horiz_interp_out%horizInterpReals4_type%wtj = horiz_interp_in%horizInterpReals4_type%wtj - horiz_interp_out%horizInterpReals4_type%src_dist = horiz_interp_in%horizInterpReals4_type%src_dist - horiz_interp_out%horizInterpReals4_type%rat_x = horiz_interp_in%horizInterpReals4_type%rat_x - horiz_interp_out%horizInterpReals4_type%rat_y = horiz_interp_in%horizInterpReals4_type%rat_y - horiz_interp_out%horizInterpReals4_type%lon_in = horiz_interp_in%horizInterpReals4_type%lon_in - horiz_interp_out%horizInterpReals4_type%lat_in = horiz_interp_in%horizInterpReals4_type%lat_in - horiz_interp_out%horizInterpReals4_type%area_frac_dst = horiz_interp_in%horizInterpReals4_type%area_frac_dst - horiz_interp_out%horizInterpReals4_type%max_src_dist = horiz_interp_in%horizInterpReals4_type%max_src_dist - horiz_interp_out%horizInterpReals4_type%is_allocated = .true. + if( allocated(horiz_interp_in%horizInterpReals4_type%faci)) & + horiz_interp_out%horizInterpReals4_type%faci = horiz_interp_in%horizInterpReals4_type%faci + + if( allocated( horiz_interp_in%horizInterpReals4_type%facj)) & + horiz_interp_out%horizInterpReals4_type%facj = horiz_interp_in%horizInterpReals4_type%facj + + if( allocated( horiz_interp_in%horizInterpReals4_type%area_src)) & + horiz_interp_out%horizInterpReals4_type%area_src = horiz_interp_in%horizInterpReals4_type%area_src + + if( allocated( horiz_interp_in%horizInterpReals4_type%area_dst)) & + horiz_interp_out%horizInterpReals4_type%area_dst = horiz_interp_in%horizInterpReals4_type%area_dst + + if( allocated( horiz_interp_in%horizInterpReals4_type%wti)) & + horiz_interp_out%horizInterpReals4_type%wti = horiz_interp_in%horizInterpReals4_type%wti + + if( allocated( horiz_interp_in%horizInterpReals4_type%wtj)) & + horiz_interp_out%horizInterpReals4_type%wtj = horiz_interp_in%horizInterpReals4_type%wtj + + if( allocated( horiz_interp_in%horizInterpReals4_type%src_dist)) & + horiz_interp_out%horizInterpReals4_type%src_dist = horiz_interp_in%horizInterpReals4_type%src_dist + + if( allocated( horiz_interp_in%horizInterpReals4_type%rat_x)) & + horiz_interp_out%horizInterpReals4_type%rat_x = horiz_interp_in%horizInterpReals4_type%rat_x + + if( allocated( horiz_interp_in%horizInterpReals4_type%rat_y)) & + horiz_interp_out%horizInterpReals4_type%rat_y = horiz_interp_in%horizInterpReals4_type%rat_y + + if( allocated( horiz_interp_in%horizInterpReals4_type%lon_in)) & + horiz_interp_out%horizInterpReals4_type%lon_in = horiz_interp_in%horizInterpReals4_type%lon_in + + if( allocated( horiz_interp_in%horizInterpReals4_type%lat_in)) & + horiz_interp_out%horizInterpReals4_type%lat_in = horiz_interp_in%horizInterpReals4_type%lat_in + + if( allocated( horiz_interp_in%horizInterpReals4_type%area_frac_dst)) & + horiz_interp_out%horizInterpReals4_type%area_frac_dst = horiz_interp_in%horizInterpReals4_type%area_frac_dst + + horiz_interp_out%horizInterpReals4_type%max_src_dist = horiz_interp_in%horizInterpReals4_type%max_src_dist + + horiz_interp_out%horizInterpReals4_type%is_allocated = .true. ! this was left out previous to mixed mode - horiz_interp_out%horizInterpReals4_type%mask_in = horiz_interp_in%horizInterpReals4_type%mask_in + if( allocated(horiz_interp_in%horizInterpReals4_type%mask_in)) & + horiz_interp_out%horizInterpReals4_type%mask_in = horiz_interp_in%horizInterpReals4_type%mask_in else call mpp_error(FATAL, "horiz_interp_type_eq: cannot assign unallocated real values from horiz_interp_in") diff --git a/test_fms/horiz_interp/test_horiz_interp.F90 b/test_fms/horiz_interp/test_horiz_interp.F90 index eb8afba07..31f0e3178 100644 --- a/test_fms/horiz_interp/test_horiz_interp.F90 +++ b/test_fms/horiz_interp/test_horiz_interp.F90 @@ -37,7 +37,7 @@ program horiz_interp_test use mpp_domains_mod, only : mpp_domains_init, domain2d use fms_mod, only : check_nml_error, fms_init use horiz_interp_mod, only : horiz_interp_init, horiz_interp_new, horiz_interp_del -use horiz_interp_mod, only : horiz_interp, horiz_interp_type +use horiz_interp_mod, only : horiz_interp, horiz_interp_type, assignment(=) use horiz_interp_type_mod, only: SPHERICAL use constants_mod, only : constants_init, PI use horiz_interp_bilinear_mod, only: horiz_interp_bilinear_new @@ -111,7 +111,7 @@ program horiz_interp_test subroutine test_horiz_interp_spherical !! grid data real(HI_TEST_KIND_), allocatable, dimension(:,:) :: lat_in_2D, lon_in_2D - type(horiz_interp_type) :: interp_t + type(horiz_interp_type) :: interp_t, interp_copy !! input data real(HI_TEST_KIND_), allocatable, dimension(:,:) :: data_src, data_dst !! output data @@ -125,7 +125,6 @@ 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_) :: R2D = 180._lkind/real(PI,HI_TEST_KIND_) real(HI_TEST_KIND_), parameter :: SMALL = 1.0e-10_lkind ! set up longitude and latitude of source/destination grid. @@ -170,6 +169,7 @@ subroutine test_horiz_interp_spherical call horiz_interp_new(interp_t, lon_in_2d, lat_in_2d, lon_out_2d, lon_out_2d, interp_method="spherical") call horiz_interp(interp_t, data_src, data_dst) call horiz_interp_spherical_wght(interp_t, wghts, verbose=1) + interp_copy = interp_t else call horiz_interp(data_src, lon_in_2D, lat_in_2D, lon_out_2D, lat_out_2D, data_dst, interp_method="spherical") endif @@ -185,7 +185,9 @@ subroutine test_horiz_interp_spherical enddo if(.not. test_solo) then call horiz_interp_del(interp_t) + call horiz_interp_del(interp_copy) call check_dealloc(interp_t) + call check_dealloc(interp_copy) endif deallocate(data_src, data_dst) deallocate(lat_in_2D, lon_in_2D) @@ -203,9 +205,8 @@ subroutine test_horiz_interp_bilinear 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 + type(horiz_interp_type) :: interp, interp_copy if (decreasing_lat) then lon_src_beg = 360.0_lkind @@ -256,6 +257,7 @@ subroutine test_horiz_interp_bilinear if (.not. test_solo) then call horiz_interp_new(interp, lon1D_src, lat1D_src, lon1D_dst, lat1D_dst, interp_method = "bilinear") call horiz_interp(interp, data_src, data_dst) + interp_copy = interp else call horiz_interp(data_src, lon1D_src, lat1D_src, lon1D_dst, lat1D_dst, data_dst, interp_method = "bilinear") endif @@ -313,7 +315,9 @@ subroutine test_horiz_interp_bilinear end do if(.not. test_solo) then call horiz_interp_del(interp) + call horiz_interp_del(interp_copy) call check_dealloc(interp) + call check_dealloc(interp_copy) endif ! --- 1dx2d version bilinear interpolation @@ -329,6 +333,7 @@ subroutine test_horiz_interp_bilinear if(.not. test_solo) then call horiz_interp_new(interp, lon1D_src, lat1D_src, lon2D_dst, lat2D_dst, interp_method = "bilinear") call horiz_interp(interp, data_src, data_dst) + interp_copy = interp else call horiz_interp(data_src, lon1D_src, lat1D_src, lon2D_dst, lat2D_dst, data_dst,interp_method="bilinear") endif @@ -386,7 +391,9 @@ subroutine test_horiz_interp_bilinear end do if(.not. test_solo) then call horiz_interp_del(interp) + call horiz_interp_del(interp_copy) call check_dealloc(interp) + call check_dealloc(interp_copy) endif if (decreasing_lat) return @@ -405,6 +412,7 @@ subroutine test_horiz_interp_bilinear call horiz_interp_new(interp,lon2D_src,lat2D_src,lon1D_dst(1:ni_src),lat1D_dst(1:nj_src), & interp_method = "bilinear") call horiz_interp(interp, data_src, data_dst) + interp_copy = interp else call horiz_interp(data_src, lon2D_src, lat2d_src, lon1D_dst(1:ni_src),lat1D_dst(1:nj_src), data_dst, & interp_method="bilinear") @@ -502,7 +510,9 @@ subroutine test_horiz_interp_bilinear end do if(.not. test_solo) then call horiz_interp_del(interp) + call horiz_interp_del(interp_copy) call check_dealloc(interp) + call check_dealloc(interp_copy) endif ! --- 2dx2d version bilinear interpolation @@ -514,6 +524,7 @@ subroutine test_horiz_interp_bilinear if(.not. test_solo) then call horiz_interp_new(interp, lon2D_src, lat2D_src, lon2D_dst, lat2D_dst, interp_method = "bilinear") call horiz_interp(interp, data_src, data_dst) + interp_copy = interp else call horiz_interp(data_src, lon2D_src, lat2d_src, lon2D_dst, lat2D_dst, data_dst, interp_method="bilinear") endif @@ -601,7 +612,9 @@ subroutine test_horiz_interp_bilinear endif if(.not. test_solo) then call horiz_interp_del(interp) + call horiz_interp_del(interp_copy) call check_dealloc(interp) + call check_dealloc(interp_copy) endif !check that data are equal do j=1, nj_src @@ -620,8 +633,7 @@ end subroutine test_horiz_interp_bilinear subroutine test_horiz_interp_bicubic !! grid data real(HI_TEST_KIND_), allocatable, dimension(:) :: lat_in_1D, lon_in_1D - real(HI_TEST_KIND_), allocatable, dimension(:,:) :: lat_in_2D, lon_in_2D - type(horiz_interp_type) :: interp_t + type(horiz_interp_type) :: interp_t, interp_copy !! input data real(HI_TEST_KIND_), allocatable, dimension(:,:) :: data_src, data_dst !! output data @@ -637,7 +649,6 @@ 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_) :: R2D = 180._lkind/real(PI,HI_TEST_KIND_) real(HI_TEST_KIND_), parameter :: SMALL = 1.0e-10_lkind ! set up longitude and latitude of source/destination grid. @@ -691,6 +702,7 @@ subroutine test_horiz_interp_bicubic if(.not. test_solo) then call horiz_interp_new(interp_t, lon_in_1d, lat_in_1d, lon_out_1d, lat_out_1d, interp_method="bicubic") call horiz_interp(interp_t, data_src, data_dst) + interp_copy = interp_t else call horiz_interp(data_src, lon_in_1D, lat_in_1D, lon_out_1D, lat_out_1D, data_dst, interp_method="bicubic") endif @@ -719,7 +731,9 @@ subroutine test_horiz_interp_bicubic enddo enddo call horiz_interp_del(interp_t) + call horiz_interp_del(interp_copy) call check_dealloc(interp_t) + call check_dealloc(interp_copy) endif do i=isc, iec do j=jsc, jec @@ -737,6 +751,7 @@ subroutine test_horiz_interp_bicubic if(.not. test_solo) then call horiz_interp_new(interp_t, lon_in_1d, lat_in_1d, lon_out_2d, lat_out_2d, interp_method="bicubic") call horiz_interp(interp_t, data_src, data_dst) + interp_copy = interp_t else call horiz_interp(data_src, lon_in_1D, lat_in_1D, lon_out_2D, lat_out_2D, data_dst, interp_method="bicubic") endif @@ -762,7 +777,9 @@ subroutine test_horiz_interp_bicubic enddo enddo call horiz_interp_del(interp_t) + call horiz_interp_del(interp_copy) call check_dealloc(interp_t) + call check_dealloc(interp_copy) endif do i=isc, iec do j=jsc, jec @@ -782,14 +799,13 @@ subroutine test_horiz_interp_conserve 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, data1_dst, data2_dst, data3_dst, data4_dst - real(HI_TEST_KIND_), allocatable, dimension(:,:) :: data1_solo, data2_solo, data3_solo, data4_solo real(HI_TEST_KIND_) :: lon_src_beg = 0._lkind, lon_src_end = 360._lkind real(HI_TEST_KIND_) :: lat_src_beg = -90._lkind, lat_src_end = 90._lkind 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 - type(horiz_interp_type) :: interp_conserve + type(horiz_interp_type) :: interp_conserve, interp_copy allocate(lon2D_src(ni_src+1, nj_src+1), lat2D_src(ni_src+1, nj_src+1) ) allocate(lon1D_src(ni_src+1), lat1D_src(nj_src+1), data_src(ni_src, nj_src) ) @@ -861,22 +877,29 @@ subroutine test_horiz_interp_conserve call horiz_interp_new(interp_conserve, lon1D_src, lat1D_src, lon1D_dst, lat1D_dst, & interp_method = "conservative") call horiz_interp(interp_conserve, data_src, data1_dst) + interp_copy = interp_conserve call horiz_interp_del(interp_conserve) + call horiz_interp_del(interp_copy) call check_dealloc(interp_conserve) + call check_dealloc(interp_copy) else call horiz_interp(data_src, lon1D_src, lat1D_src, lon1D_dst, lat1D_dst, data1_dst, & interp_method="conservative") endif call mpp_clock_end(id1) + ! --- 1dx2d version conservative interpolation call mpp_clock_begin(id2) if(.not. test_solo) then call horiz_interp_new(interp_conserve, lon1D_src, lat1D_src, lon2D_dst, lat2D_dst, & interp_method = "conservative") call horiz_interp(interp_conserve, data_src, data2_dst) + interp_copy = interp_conserve call horiz_interp_del(interp_conserve) + call horiz_interp_del(interp_copy) call check_dealloc(interp_conserve) + call check_dealloc(interp_copy) else call horiz_interp(data_src, lon1D_src, lat1D_src, lon2D_dst, lat2D_dst, data2_dst, & interp_method="conservative") @@ -889,8 +912,11 @@ subroutine test_horiz_interp_conserve call horiz_interp_new(interp_conserve, lon2D_src, lat2D_src, lon1D_dst, lat1D_dst, & interp_method = "conservative") call horiz_interp(interp_conserve, data_src, data3_dst) + interp_copy = interp_conserve call horiz_interp_del(interp_conserve) + call horiz_interp_del(interp_copy) call check_dealloc(interp_conserve) + call check_dealloc(interp_copy) else call horiz_interp(data_src, lon2D_src, lat2D_src, lon1D_dst, lat1D_dst, data3_dst, & interp_method="conservative") @@ -903,8 +929,11 @@ subroutine test_horiz_interp_conserve call horiz_interp_new(interp_conserve, lon2D_src, lat2D_src, lon2D_dst, lat2D_dst, & interp_method = "conservative") call horiz_interp(interp_conserve, data_src, data4_dst) + interp_copy = interp_conserve call horiz_interp_del(interp_conserve) + call horiz_interp_del(interp_copy) call check_dealloc(interp_conserve) + call check_dealloc(interp_copy) else call horiz_interp(data_src, lon2D_src, lat2D_src, lon2D_dst, lat2D_dst, data4_dst, & interp_method="conservative") @@ -963,7 +992,7 @@ subroutine test_horiz_interp_conserve !! Also tests creating the types via the method-specific *_new routines to ensure !! they can be created/deleted without allocation errors. subroutine test_assignment() - type(horiz_interp_type) :: Interp_new1, Interp_new2, Interp_cp, intp_3 + type(horiz_interp_type) :: Interp_new1, Interp_new2, Interp_cp real(HI_TEST_KIND_), allocatable, dimension(:) :: lat_in_1D, lon_in_1D !< 1D grid data points real(HI_TEST_KIND_), allocatable, dimension(:,:) :: lat_in_2D, lon_in_2D !< 2D grid data points real(HI_TEST_KIND_), allocatable, dimension(:) :: lat_out_1D, lon_out_1D !< 1D grid output points @@ -980,7 +1009,6 @@ subroutine test_assignment() real(HI_TEST_KIND_) :: lat_dst_beg = -90._lkind, lat_dst_end = 90._lkind !< destination grid !! starting/ending latitudes real(HI_TEST_KIND_) :: D2R = real(PI,HI_TEST_KIND_)/180._lkind !< radians per degree - real(HI_TEST_KIND_) :: R2D = 180._lkind/real(PI,HI_TEST_KIND_) !< degrees per radian real(HI_TEST_KIND_), allocatable :: lon_src_1d(:), lat_src_1d(:) !< src data used for bicubic test real(HI_TEST_KIND_), allocatable :: lon_dst_1d(:), lat_dst_1d(:) !< destination data used for bicubic test integer :: icount !< index for setting the output array when taking midpoints for bilinear