Fix bug in chgres_cube subroutine search_many (#808)

LarissaReames-NOAA · web-flow · commit 5101237fb568 · 2023-04-12T10:19:57.000-04:00
The routine is called from all MPI tasks, but some arrays passed to the routine were not allocated on all tasks. This problem was only seen when compiling/running under 'Debug' mode (all regression tests crashed at the call to search_many). Fixes #797.
diff --git a/sorc/chgres_cube.fd/surface.F90 b/sorc/chgres_cube.fd/surface.F90
@@ -851,11 +851,12 @@ subroutine interp(localpet)
    if (localpet == 0) then   
      where(mask_target_one_tile == 1) mask_target_one_tile = 0
      where(mask_target_one_tile == 2) mask_target_one_tile = 1
+     call search_many(num_fields,bundle_seaice_target,tile, search_nums,localpet, &
+                    mask=mask_target_one_tile)
+   else
+     call search_many(num_fields,bundle_seaice_target, tile,search_nums,localpet)
    endif
 
-
-   call search_many(num_fields,bundle_seaice_target,data_one_tile, mask_target_one_tile,tile,search_nums,localpet, &
-                    field_data_3d=data_one_tile_3d)
  enddo
 
  deallocate(search_nums)
@@ -977,10 +978,12 @@ subroutine interp(localpet)
      allocate(water_target_one_tile(i_target,j_target))
      water_target_one_tile = 0
      where(mask_target_one_tile == 0) water_target_one_tile = 1
-   endif
 
-   call search_many(num_fields,bundle_water_target,data_one_tile, water_target_one_tile,& 
-                    tile,search_nums,localpet,latitude=latitude_one_tile)
+     call search_many(num_fields,bundle_water_target, tile,search_nums,localpet, &
+                latitude=latitude_one_tile,mask=water_target_one_tile)
+   else
+     call search_many(num_fields,bundle_water_target, tile,search_nums,localpet)
+   endif
 
    if (localpet == 0) deallocate(water_target_one_tile)
 
@@ -1068,10 +1071,12 @@ subroutine interp(localpet)
      allocate(land_target_one_tile(i_target,j_target))
      land_target_one_tile = 0
      where(mask_target_one_tile == 1) land_target_one_tile = 1
-   endif
    
-   call search_many(num_fields,bundle_allland_target,data_one_tile, land_target_one_tile,& 
-                    tile,search_nums,localpet)
+     call search_many(num_fields,bundle_allland_target, &
+                    tile,search_nums,localpet, mask=land_target_one_tile)
+   else
+     call search_many(num_fields,bundle_allland_target, tile,search_nums,localpet)
+   endif
 
    if (localpet == 0) deallocate(land_target_one_tile)   
  enddo
@@ -1202,8 +1207,12 @@ subroutine interp(localpet)
     if(ESMF_logFoundError(rcToCheck=rc,msg=ESMF_LOGERR_PASSTHRU,line=__LINE__,file=__FILE__)) &
        call error_handler("IN FieldGather", rc)
 
-   call search_many(num_fields,bundle_landice_target,data_one_tile, land_target_one_tile,& 
-                    tile,search_nums,localpet,terrain_land=data_one_tile2,field_data_3d=data_one_tile_3d)
+   if (localpet==0) then
+      call search_many(num_fields,bundle_landice_target,tile,search_nums,localpet,&
+                terrain_land=data_one_tile2,mask=land_target_one_tile)
+   else
+      call search_many(num_fields,bundle_landice_target,tile,search_nums,localpet)
+   endif
  enddo
 
  deallocate (veg_type_target_one_tile)
@@ -1416,9 +1425,12 @@ subroutine interp(localpet)
    call ESMF_FieldGather(soil_type_target_grid, data_one_tile2, rootPet=0,tile=tile, rc=rc)
    if(ESMF_logFoundError(rcToCheck=rc,msg=ESMF_LOGERR_PASSTHRU,line=__LINE__,file=__FILE__))&
       call error_handler("IN FieldGather", rc)
-      
-   call search_many(num_fields,bundle_nolandice_target,data_one_tile, mask_target_one_tile,& 
-                    tile,search_nums,localpet,soilt_climo=data_one_tile2, field_data_3d=data_one_tile_3d)
+   if (localpet==0) then 
+      call search_many(num_fields,bundle_nolandice_target,tile,search_nums,localpet, &
+                soilt_climo=data_one_tile2, mask=mask_target_one_tile)
+   else
+      call search_many(num_fields,bundle_nolandice_target, tile,search_nums,localpet)
+   endif
    
    print*,"- CALL FieldGather FOR TARGET GRID TOTAL SOIL MOISTURE, TILE: ", tile
    call ESMF_FieldGather(soilm_tot_target_grid, data_one_tile_3d, rootPet=0, tile=tile, rc=rc)
@@ -3290,20 +3302,18 @@ end subroutine regrid_many
 !!
 !! @param[in] num_field  Number of fields to process.
 !! @param[inout] bundle_target  ESMF FieldBundle holding target fields to search
-!! @param[inout] field_data_2d  A real array of size i_target,j_target to temporarily hold data for searching
-!! @param[inout] mask  An integer array of size i_target,j_target that holds masked (0) and unmasked (1) 
-!!                     values indicating where to execute search (only at unmasked points).
 !! @param[in] tile  Current cubed sphere tile.
 !! @param[inout]  search_nums  Array length num_field holding search field numbers corresponding to each field provided for searching.
 !! @param[in]  localpet  ESMF local persistent execution thread.
 !! @param[in]  latitude  (optional) A real array size i_target,j_target of latitude on the target grid 
 !! @param[in]  terrain_land  (optional) A real array size i_target,j_target of terrain height (m) on the target grid 
 !! @param[in]  soilt_climo  (optional) A real array size i_target,j_target of climatological soil type on the target grid 
-!! @param[in]  field_data_3d (optional) An empty real array of size i_target,j_target,lsoil_target to temporarily hold soil data for searching 
+!! @param[inout] mask  (optional) An integer array of size i_target,j_target that holds masked (0) and unmasked (1)
+!!                     values indicating where to execute search (only at
+!unmasked points).
 !! @author Larissa Reames, OU CIMMS/NOAA/NSSL
- subroutine search_many(num_field,bundle_target,field_data_2d,mask, tile, &
-                         search_nums,localpet,latitude,terrain_land,soilt_climo,&
-                         field_data_3d)
+ subroutine search_many(num_field,bundle_target,tile,search_nums,localpet,latitude, &
+                        terrain_land,soilt_climo, mask)
 
  use model_grid, only                  : i_target,j_target, lsoil_target
  use program_setup, only               : external_model, input_type
@@ -3313,14 +3323,14 @@ subroutine search_many(num_field,bundle_target,field_data_2d,mask, tile, &
 
  integer, intent(in)                         :: num_field
  type(esmf_fieldbundle), intent(inout)       :: bundle_target
- real(esmf_kind_r8), intent(inout)           :: field_data_2d(i_target,j_target)
- real(esmf_kind_r8), intent(inout), optional :: field_data_3d(i_target,j_target,lsoil_target) 
+
  real(esmf_kind_r8), intent(inout), optional :: latitude(i_target,j_target)
  real(esmf_kind_r8), intent(inout), optional :: terrain_land(i_target,j_target)
  real(esmf_kind_r8), intent(inout), optional :: soilt_climo(i_target,j_target)
- integer(esmf_kind_i8), intent(inout)        :: mask(i_target,j_target)
+ integer(esmf_kind_i8), intent(inout), optional  :: mask(i_target,j_target)
  
-    
+ real(esmf_kind_r8), allocatable :: field_data_2d(:,:)   
+ real(esmf_kind_r8), allocatable :: field_data_3d(:,:,:)  
  integer, intent(in)             :: tile,localpet
  integer, intent(inout)          :: search_nums(num_field)
  
@@ -3331,46 +3341,45 @@ subroutine search_many(num_field,bundle_target,field_data_2d,mask, tile, &
  integer, parameter              :: TERRAIN_FIELD_NUM= 7
  integer :: j,k, rc, ndims
 
+ 
  do k = 1,num_field
    call ESMF_FieldBundleGet(bundle_target,k,temp_field, rc=rc)
     if(ESMF_logFoundError(rcToCheck=rc,msg=ESMF_LOGERR_PASSTHRU,line=__LINE__,file=__FILE__))&
       call error_handler("IN FieldGet", rc)
    call ESMF_FieldGet(temp_field, name=fname, dimcount=ndims,rc=rc)
         if(ESMF_logFoundError(rcToCheck=rc,msg=ESMF_LOGERR_PASSTHRU,line=__LINE__,file=__FILE__))&
           call error_handler("IN FieldGet", rc)
+   if (localpet==0) then
+       allocate(field_data_2d(i_target,j_target))
+   else
+       allocate(field_data_2d(0,0))
+   endif
    if (ndims .eq. 2) then
-       print*, "processing 2d field ", trim(fname)
-       print*, "FieldGather"
        call ESMF_FieldGather(temp_field,field_data_2d,rootPet=0,tile=tile, rc=rc)
        if(ESMF_logFoundError(rcToCheck=rc,msg=ESMF_LOGERR_PASSTHRU,line=__LINE__,file=__FILE__))&
         call error_handler("IN FieldGather", rc)
      if (localpet == 0) then 
        if (present(latitude) .and. search_nums(k).eq.SST_FIELD_NUM) then
          ! Sea surface temperatures; pass latitude field to search
-         print*, "search1"
          call search(field_data_2d, mask, i_target, j_target, tile,search_nums(k),latitude=latitude)
        elseif (present(terrain_land) .and. search_nums(k) .eq. TERRAIN_FIELD_NUM) then
          ! Terrain height; pass optional climo terrain array to search
-         print*, "search2"
          call search(field_data_2d, mask, i_target, j_target, tile,search_nums(k),terrain_land=terrain_land)
        elseif (search_nums(k) .eq. SOTYP_LAND_FIELD_NUM) then
          ! Soil type over land    
          if (fname .eq. "soil_type_target_grid") then
            ! Soil type over land when interpolating input data to target grid
            ! *with* the intention of retaining interpolated data in output
-           print*, "search3"
            call search(field_data_2d, mask, i_target, j_target, tile,search_nums(k),soilt_climo=soilt_climo)
          elseif (present(soilt_climo)) then
            if (maxval(field_data_2d) > 0 .and. (trim(external_model) .ne. "GFS" .or. trim(input_type) .ne. "grib2")) then
              ! Soil type over land when interpolating input data to target grid
              ! *without* the intention of retaining data in output file
-             print*, "search4"
              call search(field_data_2d, mask, i_target, j_target, tile, search_nums(k))
            else 
              ! If no soil type field exists in input data (e.g., GFS grib2) then don't search
              ! but simply set data to the climo field. This may result in
              ! somewhat inaccurate soil moistures as no scaling will occur 
-             print*, "search5"
              field_data_2d = soilt_climo
            endif !check field value   
          endif !sotype from target grid
@@ -3384,12 +3393,17 @@ subroutine search_many(num_field,bundle_target,field_data_2d,mask, tile, &
      if(ESMF_logFoundError(rcToCheck=rc,msg=ESMF_LOGERR_PASSTHRU,line=__LINE__,file=__FILE__))&
         call error_handler("IN FieldScatter", rc)
    else
+     if (localpet==0) then
+         allocate(field_data_3d(i_target,j_target,lsoil_target))
+     else
+         allocate(field_data_3d(0,0,0))
+     endif
+ 
      ! Process 3d fields soil temperature, moisture, and liquid
-     print*, "FieldGather"
      call ESMF_FieldGather(temp_field,field_data_3d,rootPet=0,tile=tile,rc=rc)
      if(ESMF_logFoundError(rcToCheck=rc,msg=ESMF_LOGERR_PASSTHRU,line=__LINE__,file=__FILE__))&
         call error_handler("IN FieldGather", rc)
-     print*, "processing 3d field ", trim(fname)
+        
      if (localpet==0) then 
        do j = 1, lsoil_target
          field_data_2d = field_data_3d(:,:,j)
@@ -3400,7 +3414,9 @@ subroutine search_many(num_field,bundle_target,field_data_2d,mask, tile, &
      call ESMF_FieldScatter(temp_field, field_data_3d, rootPet=0, tile=tile,rc=rc)
       if(ESMF_logFoundError(rcToCheck=rc,msg=ESMF_LOGERR_PASSTHRU,line=__LINE__,file=__FILE__))&
         call error_handler("IN FieldScatter", rc)
+     deallocate(field_data_3d)
    endif !ndims
+   deallocate(field_data_2d)
  end do !fields
 
  end subroutine search_many
diff --git a/tests/chgres_cube/ftst_surface_search_many.F90 b/tests/chgres_cube/ftst_surface_search_many.F90
@@ -295,8 +295,8 @@ program surface_interp
  input_type="restart"
 
  !Call the search many routine to test search and replace
- call search_many(num_fields,bundle_search1,dummy_2d,mask_target_search,1,field_nums,localpet, &
-                  soilt_climo=soilt_climo)
+ call search_many(num_fields,bundle_search1,1,field_nums,localpet, &
+                  soilt_climo=soilt_climo,mask=mask_target_search)
 
  call ESMF_FieldBundleDestroy(bundle_search1,rc=rc)
    if(ESMF_logFoundError(rcToCheck=rc,msg=ESMF_LOGERR_PASSTHRU,line=__LINE__,file=__FILE__)) &
@@ -342,8 +342,8 @@ program surface_interp
  external_model="HRRR"
 
  !Call the search many routine to test search and replace
- call search_many(num_fields,bundle_search2,dummy_2d,mask_target_search,1,field_nums,localpet, &
-         soilt_climo=soilt_climo)
+  call search_many(num_fields,bundle_search2,1,field_nums,localpet, &
+         soilt_climo=soilt_climo,mask=mask_target_search)
 
  call ESMF_FieldBundleDestroy(bundle_search2,rc=rc)
    if(ESMF_logFoundError(rcToCheck=rc,msg=ESMF_LOGERR_PASSTHRU,line=__LINE__,file=__FILE__)) &
@@ -385,8 +385,8 @@ program surface_interp
  allocate(field_nums(num_fields))
  field_nums = (/11,7,224/)
  !Call the search many routine to test some branches of default behavior
- call search_many(num_fields,bundle_default1,dummy_2d,mask_default,1,field_nums,localpet, &
-      latitude=latitude_default,terrain_land=terrain_land,soilt_climo=soilt_climo)
+ call search_many(num_fields,bundle_default1,1,field_nums,localpet, &
+      latitude=latitude_default,terrain_land=terrain_land,soilt_climo=soilt_climo,mask=mask_default)
 
  print*,"Check results for bundle_default1."
 
@@ -441,8 +441,8 @@ program surface_interp
  input_type="grib2"
  external_model="GFS"
  !Call the search many routine to test behavior for GFS grib2 soil type
- call search_many(num_fields,bundle_default2,dummy_2d,mask_default,1,field_nums,localpet,&
-      soilt_climo=soilt_climo)
+ call search_many(num_fields,bundle_default2,1,field_nums,localpet,&
+      soilt_climo=soilt_climo,mask=mask_default)
 
  call ESMF_FieldBundleDestroy(bundle_default2,rc=rc)
    if(ESMF_logFoundError(rcToCheck=rc,msg=ESMF_LOGERR_PASSTHRU,line=__LINE__,file=__FILE__))&
@@ -478,8 +478,7 @@ program surface_interp
  field_nums(:) = (/21/)
 
  !Call the search many routine to test behavior for GFS grib2 soil type
- call search_many(num_fields,bundle_3d_search,dummy_2d,mask_target_search,1,field_nums,localpet,&
-      field_data_3d=dummy_3d)
+ call search_many(num_fields,bundle_3d_search,1,field_nums,localpet,mask=mask_target_search)
 
  call ESMF_FieldBundleDestroy(bundle_3d_search,rc=rc)
    if(ESMF_logFoundError(rcToCheck=rc,msg=ESMF_LOGERR_PASSTHRU,line=__LINE__,file=__FILE__))&
