From 5cfc4d2f3f97ad4c5fe33e3620aa531fb8e0ea51 Mon Sep 17 00:00:00 2001 From: Aubrey Dugger Date: Wed, 4 May 2022 16:39:43 -0600 Subject: [PATCH] Move overland roughness assignment out of time cycling disaggregation routine (#553) * Add RT disaggregation to read2dlsm() * Added rt=.logical. to read2dlsm to call Aubrey's regrid_lowres_to_highres while LSM data are all on rank 0. Useful for OV_ROUGH2D and possibly others in the future. Co-authored-by: Ryan Cabell Co-authored-by: Soren Rasmussen --- .travis.yml.prev | 20 ------ trunk/NDHMS/Routing/Noah_distr_routing.F | 33 ++-------- trunk/NDHMS/Routing/module_HYDRO_io.F | 82 +++++++++++++++++++++--- trunk/NDHMS/Routing/module_RT.F | 4 +- 4 files changed, 79 insertions(+), 60 deletions(-) delete mode 100644 .travis.yml.prev diff --git a/.travis.yml.prev b/.travis.yml.prev deleted file mode 100644 index 0684ac00e..000000000 --- a/.travis.yml.prev +++ /dev/null @@ -1,20 +0,0 @@ -sudo: required -language: bash -services: - - docker - -before_install: - #Setup test directory on travis - - TEST_DIR_TRAVIS=$HOME/test_dir - - sudo mkdir $TEST_DIR_TRAVIS - - sudo chmod -R 777 $TEST_DIR_TRAVIS - - cp -r $TRAVIS_BUILD_DIR $TEST_DIR_TRAVIS/candidate - - git clone https://github.com/NCAR/wrf_hydro_nwm_public.git $TEST_DIR_TRAVIS/reference - - cd $TEST_DIR_TRAVIS/reference - - git checkout $TRAVIS_BRANCH - - sudo chmod -R 777 $TEST_DIR_TRAVIS - #Get docker images - - docker pull wrfhydro/dev:modeltesting - -script: - - docker run -e TRAVIS=1 -t -v $TEST_DIR_TRAVIS/candidate:/home/docker/candidate -v $TEST_DIR_TRAVIS/reference:/home/docker/reference wrfhydro/dev:modeltesting --config nwm_ana nwm_long_range gridded reach --domain_tag dev \ No newline at end of file diff --git a/trunk/NDHMS/Routing/Noah_distr_routing.F b/trunk/NDHMS/Routing/Noah_distr_routing.F index 04bec13b5..df0f04936 100644 --- a/trunk/NDHMS/Routing/Noah_distr_routing.F +++ b/trunk/NDHMS/Routing/Noah_distr_routing.F @@ -889,7 +889,7 @@ subroutine disaggregateDomain_drv(did) RT_DOMAIN(did)%NEXP, & rt_domain(did)%overland%properties%distance_to_neighbor, & RT_DOMAIN(did)%INFXSWGT, & - RT_DOMAIN(did)%OVROUGHRTFAC,RT_DOMAIN(did)%LKSATFAC, & + RT_DOMAIN(did)%LKSATFAC, & rt_domain(did)%overland%streams_and_lakes%ch_netrt, & RT_DOMAIN(did)%SH2OWGT, & RT_DOMAIN(did)%subsurface%grid_transform%smcrefrt, & @@ -897,11 +897,9 @@ subroutine disaggregateDomain_drv(did) RT_DOMAIN(did)%subsurface%grid_transform%smcmaxrt, & RT_DOMAIN(did)%subsurface%grid_transform%smcwltrt, & rt_domain(did)%subsurface%grid_transform%smcrt, & - rt_domain(did)%overland%properties%roughness, & rt_domain(did)%overland%streams_and_lakes%lake_mask, & rt_domain(did)%subsurface%properties%lksatrt, & rt_domain(did)%subsurface%properties%nexprt, & - RT_DOMAIN(did)%OV_ROUGH2d, & RT_DOMAIN(did)%subsurface%properties%sldpth, & RT_DOMAIN(did)%soiltypRT, RT_DOMAIN(did)%soiltyp, & rt_domain(did)%ELRT, RT_DOMAIN(did)%iswater, & @@ -926,9 +924,9 @@ end subroutine disaggregateDomain_drv !=================================================================================================== subroutine disaggregateDomain(IX, JX, NSOIL, IXRT, JXRT, AGGFACTRT, & SICE, SMC, SH2OX, INFXSRT, area_lsm, SMCMAX1, SMCREF1, & - SMCWLT1, VEGTYP, LKSAT, NEXP, dist,INFXSWGT,OVROUGHRTFAC, & - LKSATFAC, CH_NETRT,SH2OWGT,SMCREFRT, INFXSUBRT,SMCMAXRT, & - SMCWLTRT,SMCRT, OVROUGHRT, LAKE_MSKRT, LKSATRT, NEXPRT, OV_ROUGH2d, & + SMCWLT1, VEGTYP, LKSAT, NEXP, dist, INFXSWGT, & + LKSATFAC, CH_NETRT, SH2OWGT, SMCREFRT, INFXSUBRT, SMCMAXRT, & + SMCWLTRT, SMCRT, LAKE_MSKRT, LKSATRT, NEXPRT, & SLDPTH, soiltypRT, soiltyp, elrt, iswater, impervfrac, imperv_adj) #ifdef MPP_LAND use module_mpp_land, only: left_id,down_id,right_id, & @@ -957,7 +955,6 @@ subroutine disaggregateDomain(IX, JX, NSOIL, IXRT, JXRT, AGGFACTRT, & real, intent(in), dimension(IX,JX) :: SMCWLT1 ! coarse grid wilting point real, intent(in), dimension(IX,JX) :: LKSAT ! coarse grid lateral ksat (m/s) real, intent(in), dimension(IX,JX) :: NEXP ! coarse grid n exponent - real, intent(in), dimension(ix,jx) :: OV_ROUGH2d ! overland roughness ! LSM states: real, intent(in), dimension(IX,JX,NSOIL) :: SMC ! total soil moisture (m3/m3) real, intent(in), dimension(IX,JX,NSOIL) :: SH2OX ! liquid soil moisture (m3/m3) @@ -967,7 +964,6 @@ subroutine disaggregateDomain(IX, JX, NSOIL, IXRT, JXRT, AGGFACTRT, & real, intent(in), dimension(IXRT,JXRT,9) :: dist ! routing grid cell distances (m) and area (m2) ! TODO: can we just pass in area since we don't need other distances? - real, intent(in), dimension(IXRT,JXRT) :: OVROUGHRTFAC ! overland roughness adj factor real, intent(in), dimension(IXRT,JXRT) :: LKSATFAC ! lateral ksat adj factor real, intent(in), dimension(IXRT,JXRT) :: elrt ! elevation grid (m) integer, intent(in), dimension(IXRT,JXRT) :: CH_NETRT ! channel network routing grid @@ -987,7 +983,6 @@ subroutine disaggregateDomain(IX, JX, NSOIL, IXRT, JXRT, AGGFACTRT, & real, intent(out), dimension(IXRT,JXRT,NSOIL) :: SMCWLTRT ! wilting point on routing grid real, intent(out), dimension(IXRT,JXRT) :: LKSATRT ! lateral ksat on the routing grid (m/s) real, intent(out), dimension(IXRT,JXRT) :: NEXPRT ! n exponent on the routing grid - real, intent(out), dimension(IXRT,JXRT) :: OVROUGHRT ! overland roughness on the routing grid ! States: real, intent(out), dimension(IX,JX,NSOIL) :: SICE ! soil ice content on coarse grid (m3/m3) real, intent(out), dimension(IXRT,JXRT,NSOIL) :: SMCRT @@ -1169,25 +1164,6 @@ subroutine disaggregateDomain(IX, JX, NSOIL, IXRT, JXRT, AGGFACTRT, & ! Now do simple grid remapping tasks - ! Overland roughness: - !DJG map ov roughness as function of land use provided in VEGPARM.TBL... - ! --- added extra check for VEGTYP for 'masked-out' locations... - ! --- out of basin locations (VEGTYP=0) have OVROUGH hardwired to 0.1 - IF (VEGTYP(I,J).LE.0) then - OVROUGHRT(IXXRT,JYYRT) = 0.1 !COWS mask test - ELSE - OVROUGHRT(IXXRT,JYYRT) = OV_ROUGH2d(i,j) - ! Modify based on impervious fraction - ! See Liong et al. 1989 for linear weighting of "smoothness" (1/roughness) - ! Assuming roughness of 0.02 for impervious and native cell roughness for pervious - if (imperv_adj .ne. 0) then - OVROUGHRT(IXXRT,JYYRT) = 1. / ((1./0.02)*impervfrac(IXXRT,JYYRT) + & ! impervious fraction - (1./OVROUGHRT(IXXRT,JYYRT))*(1.-impervfrac(IXXRT,JYYRT))) ! pervious fraction - endif - ! Apply user-supplied adjustment factor - OVROUGHRT(IXXRT,JYYRT) = OVROUGHRT(IXXRT,JYYRT)*OVROUGHRTFAC(IXXRT,JYYRT) - END IF - ! Lateral ksat !DJG 6.12.08 Map lateral hydraulic conductivity and apply distributed scaling ! --- factor that will be read in from hires terrain file @@ -1283,7 +1259,6 @@ subroutine disaggregateDomain(IX, JX, NSOIL, IXRT, JXRT, AGGFACTRT, & call MPP_LAND_COM_REAL(INFXSUBRT,IXRT,JXRT,99) call MPP_LAND_COM_REAL(LKSATRT,IXRT,JXRT,99) call MPP_LAND_COM_REAL(NEXPRT,IXRT,JXRT,99) - call MPP_LAND_COM_REAL(OVROUGHRT,IXRT,JXRT,99) call MPP_LAND_COM_INTEGER(LAKE_MSKRT,IXRT,JXRT,99) do i = 1, NSOIL call MPP_LAND_COM_REAL(SMCMAXRT(:,:,i),IXRT,JXRT,99) diff --git a/trunk/NDHMS/Routing/module_HYDRO_io.F b/trunk/NDHMS/Routing/module_HYDRO_io.F index 86201a2f7..b4312a88e 100644 --- a/trunk/NDHMS/Routing/module_HYDRO_io.F +++ b/trunk/NDHMS/Routing/module_HYDRO_io.F @@ -11171,8 +11171,8 @@ subroutine hdtbl_out(did) call hdtbl_out_nc(did,ncid, count,count_flag,"LKSAT",rt_domain(did)%LKSAT,"",ixd,jxd) call hdtbl_out_nc(did,ncid, count,count_flag,"NEXP",rt_domain(did)%NEXP,"",ixd,jxd) end do - end subroutine hdtbl_out + subroutine hdtbl_in_nc(did) implicit none integer :: did @@ -11180,7 +11180,7 @@ subroutine hdtbl_in_nc(did) call read2dlsm(did,trim(nlst(did)%hydrotbl_f),"SMCMAX1",rt_domain(did)%SMCMAX1,ierr) call read2dlsm(did,trim(nlst(did)%hydrotbl_f),"SMCREF1",rt_domain(did)%SMCREF1,ierr) call read2dlsm(did,trim(nlst(did)%hydrotbl_f),"SMCWLT1",rt_domain(did)%SMCWLT1,ierr) - call read2dlsm(did,trim(nlst(did)%hydrotbl_f),"OV_ROUGH2D",rt_domain(did)%OV_ROUGH2D,ierr) + call read2dlsm(did,trim(nlst(did)%hydrotbl_f),"OV_ROUGH2D",rt_domain(did)%overland%properties%roughness,ierr, rt=.true.) call read2dlsm(did,trim(nlst(did)%hydrotbl_f),"LKSAT",rt_domain(did)%LKSAT,ierr) call read2dlsm(did,trim(nlst(did)%hydrotbl_f),"NEXP",rt_domain(did)%NEXP,ierr) ! Letting this variable be optional and setting to global default value if not found @@ -11190,7 +11190,7 @@ subroutine hdtbl_in_nc(did) endif end subroutine hdtbl_in_nc -subroutine read2dlsm(did,file,varName,varOut,ierr) +subroutine read2dlsm(did,file,varName,varOut,ierr,rt) use module_mpp_land,only: mpp_land_bcast_int1 implicit none integer :: did, ncid , iret @@ -11198,29 +11198,91 @@ subroutine read2dlsm(did,file,varName,varOut,ierr) real,dimension(:,:) :: varOut character(len=256) :: units integer, intent(out) :: ierr -#ifdef MPP_LAND + logical, optional, intent(in) :: rt + logical :: regrid + real,allocatable,dimension(:,:) :: tmpArr + +#ifdef MPP_LAND if(my_id .eq. io_id) then +#endif allocate(tmpArr(global_nx,global_ny)) iret = nf90_open(trim(file), NF90_NOWRITE, ncid) call get_2d_netcdf(trim(varName), ncid, tmpArr, units, global_nx, global_ny, & .false., ierr) iret = nf90_close(ncid) +#ifdef MPP_LAND else allocate(tmpArr(1,1)) endif - call decompose_data_real (tmpArr,varOut) +#endif + + if (present(rt)) then + regrid = rt + else + regrid = .false. + endif + + if (regrid) then + call regrid_lowres_to_highres(did, tmpArr, varOut, rt_domain(did)%ixrt, rt_domain(did)%jxrt) + else + call decompose_data_real (tmpArr,varOut) + endif + +#ifdef MPP_LAND call mpp_land_bcast_int1(ierr) +#endif + deallocate(tmpArr) +end subroutine read2dlsm + +subroutine regrid_lowres_to_highres(did, lowres_grid, highres_grid, ixrt, jxrt) + + implicit none + integer :: did + integer :: ixrt, jxrt + real, dimension(global_nx, global_ny) :: lowres_grid + real, dimension(ixrt,jxrt) :: highres_grid + ! Local variables + integer :: i, j, irt, jrt, aggfacxrt, aggfacyrt + +#ifdef MPP_LAND + real,allocatable,dimension(:,:) :: tmpArr + if(my_id .eq. io_id) then + allocate(tmpArr(global_rt_nx, global_rt_ny)) +#endif + + do j = 1,global_ny ! Start coarse grid j loop + do i = 1,global_nx ! Start coarse grid i loop + + do aggfacyrt = nlst(did)%AGGFACTRT-1,0,-1 ! Start disagg fine grid j loop + do aggfacxrt = nlst(did)%AGGFACTRT-1,0,-1 ! Start disagg fine grid i loop + + irt = i * nlst(did)%AGGFACTRT - aggfacxrt ! Define fine grid i + jrt = j * nlst(did)%AGGFACTRT - aggfacyrt ! Define fine grid j +#ifdef MPP_LAND + ! if(left_id.ge.0) irt = irt + 1 + ! if(down_id.ge.0) jrt = jrt + 1 + tmpArr(irt,jrt) = lowres_grid(i,j) #else - iret = nf90_open(trim(file), NF90_NOWRITE, ncid) - call get_2d_netcdf(trim(varName), ncid, varOut, units, rt_domain(did)%ix, rt_domain(did)%jx, & - .false., ierr) - iret = nf90_close(ncid) + highres_grid(irt,jrt) = lowres_grid(i,j) #endif -end subroutine read2dlsm + end do + end do + + end do + end do + +#ifdef MPP_LAND + else + allocate(tmpArr(1,1)) + endif + call decompose_RT_real(tmpArr, highres_grid, global_rt_nx, global_rt_ny, ixrt, jxrt) + deallocate(tmpArr) +#endif +end subroutine regrid_lowres_to_highres subroutine read_channel_only (olddateIn, hgrid, indir, dtbl) !use module_HYDRO_io, only: read_rst_crt_reach_nc diff --git a/trunk/NDHMS/Routing/module_RT.F b/trunk/NDHMS/Routing/module_RT.F index aa16207fc..c21a7a0ad 100644 --- a/trunk/NDHMS/Routing/module_RT.F +++ b/trunk/NDHMS/Routing/module_RT.F @@ -660,6 +660,7 @@ subroutine LandRT_ini(did) use config_base, only: nlst, noah_lsm use module_RT_data, only: rt_domain use module_gw_gw2d_data, only: gw2d + use module_HYDRO_io, only: regrid_lowres_to_highres #ifdef HYDRO_D use module_HYDRO_io, only: output_lake_types #endif @@ -1215,7 +1216,7 @@ subroutine LandRT_ini(did) !Apply calibration scaling factors to sfc roughness and retention depth here... rt_domain(did)%overland%properties%retention_depth = rt_domain(did)%overland%properties%retention_depth * rt_domain(did)%RETDEPRTFAC - ! Removing roughness parameter update since it has not been populated yet... currently happens in Noah_dist_routing (AD) + rt_domain(did)%overland%properties%roughness = rt_domain(did)%overland%properties%roughness * rt_domain(did)%OVROUGHRTFAC !ADCHANGE: Moved this channel cell setting from OV_RTNG so it is outside !of overland routine (frequently called) and time loop. @@ -1244,6 +1245,7 @@ subroutine LandRT_ini(did) #ifdef MPP_LAND ! communicate the value to call MPP_LAND_COM_REAL(rt_domain(did)%overland%properties%retention_depth,rt_domain(did)%IXRT,rt_domain(did)%JXRT,99) + call MPP_LAND_COM_REAL(rt_domain(did)%overland%properties%roughness,rt_domain(did)%IXRT,rt_domain(did)%JXRT,99) call MPP_LAND_COM_REAL(rt_domain(did)%overland%properties%surface_slope_x,rt_domain(did)%IXRT,rt_domain(did)%JXRT,99) call MPP_LAND_COM_REAL(rt_domain(did)%overland%properties%surface_slope_y,rt_domain(did)%IXRT,rt_domain(did)%JXRT,99) do i = 1, 8