! Program Name: ! Author(s)/Contact(s): ! Abstract: ! History Log: ! ! Usage: ! Parameters: ! Input Files: ! ! Output Files: ! ! ! Condition codes: ! ! If appropriate, descriptive troubleshooting instructions or ! likely causes for failures could be mentioned here with the ! appropriate error code ! ! User controllable options: ! This subrouting includs the data structure and tools used for NHDPlus network mapping. module module_UDMAP use module_namelist, only: nlst_rt #ifdef MPP_LAND use module_mpp_land, only: my_id, local_startx_rt, local_starty_rt, & local_endx_rt,local_endy_rt, left_id, right_id, down_id, up_id, mpp_collect_1d_int_mem, & IO_id , numprocs use module_mpp_land, only: mpp_land_bcast_int, mpp_land_bcast_real8_1d, mpp_land_bcast_int1 use module_mpp_land, only: sum_int1d, global_rt_nx, global_rt_ny, write_IO_rt_int, MPP_LAND_COM_INTEGER use MODULE_mpp_ReachLS, only : updatelinkv, ReachLS_write_io, com_write1dInt, & com_decomp1dInt, pack_decomp_int, pack_decomp_real8 #endif implicit none #ifndef MPP_LAND integer, parameter :: numprocs=1 #endif #include type userDefineMapping integer, allocatable, dimension(:) :: grid_i, grid_j real, allocatable, dimension(:) :: weight, nodeArea, cellArea integer :: ngrids integer :: myid ! for bucket model definition real, allocatable, dimension(:) :: cellWeight integer, allocatable, dimension(:) :: cell_i, cell_j integer :: ncell end type userDefineMapping TYPE ( userDefineMapping ), allocatable, DIMENSION (:) :: LUDRSL integer, allocatable, dimension(:) :: bufid real*8 , allocatable, dimension(:) :: bufw integer :: LNUMRSL ! number of local links integer :: ter_rt_flag real*8, allocatable, dimension(:) :: basns_area integer :: gnpid, lnsize integer, allocatable, dimension(:) :: bufi,bufj contains subroutine UDMP_ini(nlinksl,ixrt,jxrt,rtmask, OVRTSWCRT, SUBRTSWCRT,cell_area) !This is the driver for user defined mapping file funciton application. integer :: ixrt, jxrt, OVRTSWCRT, SUBRTSWCRT, nlinksl integer, intent(in), dimension(ixrt,jxrt):: rtmask integer :: npid !local variable. real,dimension(:,:) :: cell_area ter_rt_flag = 0 if(OVRTSWCRT .eq. 1 .or. SUBRTSWCRT .eq. 1) then ter_rt_flag = 1 endif call readUDMP(ixrt,jxrt,npid,nlinksl) call UDMP2LOCAL(npid,ixrt,jxrt,rtmask, ter_rt_flag) call getUDMP_area(cell_area) end subroutine UDMP_ini subroutine readUDMP(ixrt,jxrt,npid, nlinksl) implicit none integer :: i,j,Ndata, did, Npid, nlinksl, k, m, kk integer,allocatable,dimension(:) :: g1bufid, gbufid, linkid ,bufidflag, & bufid_tmp, nprocs_map, lnsizes, istart integer :: ix_bufid, ii, ixrt,jxrt integer, allocatable, dimension(:) :: gbufi,gbufj,bufsize real*8 , allocatable, dimension(:) :: gbufw did = 1 call get_dimension(trim(nlst_rt(did)%UDMAP_FILE), ndata, npid) #ifdef MPP_LAND gnpid = npid allocate (lnsizes(numprocs)) if(my_id .eq. io_id) then allocate (istart(numprocs)) allocate (nprocs_map(ndata)) allocate(gbufi(ndata)) allocate(gbufj(ndata)) call get1d_int(trim(nlst_rt(did)%UDMAP_FILE),"i_index",gbufi) call get1d_int(trim(nlst_rt(did)%UDMAP_FILE),"j_index",gbufj) endif call get_nprocs_map(ixrt,jxrt,gbufi,gbufj,nprocs_map,ndata) if(my_id .eq. io_id) then lnsizes = 0 do i =1 , ndata if(nprocs_map(i) .gt. 0) then lnsizes(nprocs_map(i)) = lnsizes(nprocs_map(i)) + 1 endif enddo endif call mpp_land_bcast_int(numprocs,lnsizes) if(my_id .eq. io_id ) then kk = 0 do i = 1, numprocs kk = kk + lnsizes(i) end do end if if(my_id .eq. IO_id) then ii = 1 do i = 1, numprocs istart(i) = ii if(lnsizes(i) .gt. 0) then ii = lnsizes(i) + ii else istart(i) = -999 endif end do endif if(lnsizes(my_id+1) .gt. 0) allocate(bufi(lnsizes(my_id+1) )) call pack_decomp_int(gbufi, ndata, nprocs_map, lnsizes, istart,bufi) if(my_id .eq. io_id) then if(allocated(gbufi)) deallocate(gbufi) endif if(lnsizes(my_id+1) .gt. 0) allocate(bufj(lnsizes(my_id+1) )) call pack_decomp_int(gbufj, ndata, nprocs_map, lnsizes, istart,bufj) if(my_id .eq. io_id) then if(allocated(gbufj)) deallocate(gbufj) endif ! check bufid ! check polyid and linkid allocate(linkid(nlinksl)) if(my_id .eq. io_id) then call get1d_int(trim(nlst_rt(did)%route_link_f),"link",linkid) allocate(gbufid(npid)) call get1d_int(trim(nlst_rt(did)%UDMAP_FILE),"polyid",gbufid) endif #ifdef MPP_LAND if(nlinksl .gt. 0) then call mpp_land_bcast_int(nlinksl,linkid) endif call com_decomp1dInt(gbufid,npid,bufid_tmp,ix_bufid) #endif if(ix_bufid .gt. 0) then allocate(bufidflag(ix_bufid)) bufidflag = -999 endif do i = 1, ix_bufid do j = 1, nlinksl if(bufid_tmp(i) .eq. linkid(j)) then bufidflag(i) = bufid_tmp(i) goto 102 endif end do 102 continue end do #ifdef MPP_LAND call com_write1dInt(bufidflag,ix_bufid,gbufid,npid) #endif if(ix_bufid .gt. 0) then if(allocated(bufidflag)) deallocate(bufidflag) if(allocated(bufid_tmp)) deallocate(bufid_tmp) endif if(allocated(linkid)) deallocate(linkid) if(my_id .eq. io_id) then allocate(bufsize(npid)) allocate(g1bufid(ndata)) call get1d_int(trim(nlst_rt(did)%UDMAP_FILE),"overlaps",bufsize) g1bufid = -999 i = 1 do k = 1, npid do j = 1, bufsize(k) g1bufid(i) = gbufid(k) i = i + 1 end do enddo if(allocated(bufsize)) deallocate(bufsize) endif if(my_id .eq. io_id) then if(allocated(gbufid)) deallocate(gbufid) endif if(lnsizes(my_id+1) .gt. 0) allocate(bufid(lnsizes(my_id+1) )) call pack_decomp_int(g1bufid, ndata, nprocs_map, lnsizes, istart,bufid) if(my_id .eq. io_id) then if(allocated(g1bufid)) deallocate(g1bufid) endif if(my_id .eq. io_id) then allocate(gbufw(ndata)) call get1d_real8(trim(nlst_rt(did)%UDMAP_FILE),"regridweight",gbufw) endif if(lnsizes(my_id+1) .gt. 0) allocate(bufw(lnsizes(my_id+1) )) call pack_decomp_real8(gbufw, ndata, nprocs_map, lnsizes, istart,bufw) if(my_id .eq. io_id) then if(allocated(gbufw)) deallocate(gbufw) endif if(my_id .eq. io_id) then if(allocated(nprocs_map)) deallocate (nprocs_map) if(allocated(istart)) deallocate (istart) endif lnsize = lnsizes(my_id + 1) if(allocated(lnsizes)) deallocate(lnsizes) #else call hydro_stop("FATAL ERROR in UDMP : sequential not defined.") #endif end subroutine readUDMP subroutine UDMP2LOCAL(npid,ix,jx,rtmask, ter_rt_flag) implicit none integer :: i,j,k, ngrids, ix,jx, starti,startj, endi,endj, ii,jj, npid, kk integer, intent(in), dimension(ix,jx) :: rtmask integer, dimension(lnsize) :: lndflag,gridflag , tmpgridflag integer :: ter_rt_flag, m, c ! find ngrids is 0 so that we need to mapping from subsurface runoff. #ifdef MPP_LAND if(left_id .ge. 0) then starti = local_startx_rt + 1 else starti = local_startx_rt endif if(down_id .ge. 0) then startj = local_starty_rt + 1 else startj = local_starty_rt endif if(right_id .ge. 0) then endi = local_startx_rt + ix -2 else endi = local_startx_rt + ix -1 endif if(up_id .ge. 0) then endj = local_starty_rt + jx -2 else endj = local_starty_rt + jx -1 endif #else starti = 1 startj = 1 endi = ix endj = jx #endif gridflag = 0 lndflag = 0 #ifdef MPP_LAND k = 0 do i = 1, lnsize if(bufid(i) .gt. 0) then if(bufi(i) .ge. starti .and. bufj(i) .ge. startj .and. & bufi(i) .le. endi .and. bufj(i) .le. endj) then if(k .eq. 0) then k = 1 else if(bufid(i) .ne. bufid(i-1)) k = k + 1 endif lndflag(k) = lndflag(k) + 1 if(ter_rt_flag .eq. 1) then if(rtmask(bufi(i)-local_startx_rt+1,bufj(i)-local_starty_rt+1) .ge. 0) then gridflag(k) = gridflag(k) + 1 endif endif endif endif end do ! decide how many mapping land grids on current domain ! tmpgridflag = gridflag #ifdef MPP_LAND ! call mpp_collect_1d_int_mem(npid,tmpgridflag) #endif ! decide how many user defined links on current domain kk = k LNUMRSL = 0 do k = 1, lnsize if(lndflag(k) .gt. 0) LNUMRSL = LNUMRSL + 1 enddo if(LNUMRSL .gt. 0) then allocate(LUDRSL(LNUMRSL)) allocate( basns_area(LNUMRSL) ) else write(6,*) "Warning: no routing links found." call cleanBuf() return endif kk = 0 do k = 1, lnsize if( bufid(k) .ge. 0 ) then if (bufi(k) .ge. starti .and. bufj(k) .ge. startj .and. & bufi(k) .le. endi .and. bufj(k) .le. endj ) then if(kk .eq. 0) then kk = 1 else if(bufid(k) .ne. bufid(k-1)) kk = kk + 1 endif LUDRSL(kk)%myid = bufid(k) LUDRSL(kk)%ngrids = -999 if(gridflag(kk) .gt. 0) then LUDRSL(kk)%ngrids = gridflag(kk) if(.not. allocated(LUDRSL(kk)%weight) ) then allocate( LUDRSL(kk)%weight(LUDRSL(kk)%ngrids )) allocate( LUDRSL(kk)%grid_i(LUDRSL(kk)%ngrids )) allocate( LUDRSL(kk)%grid_j(LUDRSL(kk)%ngrids )) allocate( LUDRSL(kk)%nodeArea(LUDRSL(kk)%ngrids )) endif endif ! define bucket variables LUDRSL(kk)%ncell = lndflag(kk) if(.not. allocated(LUDRSL(kk)%cellweight) ) then allocate( LUDRSL(kk)%cellweight(LUDRSL(kk)%ncell)) allocate( LUDRSL(kk)%cell_i(LUDRSL(kk)%ncell)) allocate( LUDRSL(kk)%cell_j(LUDRSL(kk)%ncell)) allocate( LUDRSL(kk)%cellArea(LUDRSL(kk)%ncell)) endif endif endif enddo ! maping grid_i, grid_j and weight kk = 0 m = 1 c = 1 do i = 1, lnsize if( (bufid(i) .ge. 0) ) then if(bufi(i) .ge. starti .and. bufj(i) .ge. startj .and. & bufi(i) .le. endi .and. bufj(i) .le. endj) then if(kk .eq. 0) then kk = 1 else if(bufid(i) .ne. bufid(i-1)) then kk = kk + 1 m = 1 c = 1 endif endif if(LUDRSL(kk)%ngrids .gt. 0) then if(rtmask(bufi(i)-local_startx_rt+1,bufj(i)-local_starty_rt+1) .ge. 0) then LUDRSL(kk)%grid_i(m) = bufi(i) - local_startx_rt+1 LUDRSL(kk)%grid_j(m) = bufj(i) - local_starty_rt+1 LUDRSL(kk)%weight(m) = bufw(i) m = m + 1 endif endif !! begin define bucket variables LUDRSL(kk)%cell_i(c) = bufi(i) - local_startx_rt+1 LUDRSL(kk)%cell_j(c) = bufj(i) - local_starty_rt+1 LUDRSL(kk)%cellWeight(c) = bufw(i) c = c + 1 !! end define bucket variables endif endif end do call cleanBuf() #else call hydro_stop("FATAL ERROR in UDMP: Sequential not work.") #endif end subroutine UDMP2LOCAL subroutine cleanBuf() if(allocated(bufi)) deallocate(bufi) if(allocated(bufj)) deallocate(bufj) if(allocated(bufw)) deallocate(bufw) if(allocated(bufid)) deallocate(bufid) end subroutine cleanBuf subroutine get_dimension(fileName, ndata,npid) implicit none character(len=*) fileName integer ncid , iret, ndata,npid, dimid #ifdef MPP_LAND if(my_id .eq. IO_id) then #endif iret = nf_open(fileName, NF_NOWRITE, ncid) if (iret /= 0) then write(*,'("FATAL ERROR: Problem opening mapping file: ''", A, "''")') & trim(fileName) call hydro_stop("In get_dimension() - Problem opening mapping file.") endif iret = nf_inq_dimid(ncid, "polyid", dimid) if (iret /= 0) then print*, "nf_inq_dimid: polyid" call hydro_stop("In get_dimension() - nf_inq_dimid: polyid") endif iret = nf_inq_dimlen(ncid, dimid, npid) iret = nf_inq_dimid(ncid, "data", dimid) if (iret /= 0) then print*, "nf_inq_dimid: data" call hydro_stop("In get_file_dimension() - nf_inq_dimid: data") endif iret = nf_inq_dimlen(ncid, dimid, ndata) iret = nf_close(ncid) #ifdef MPP_LAND endif call mpp_land_bcast_int1(ndata) call mpp_land_bcast_int1(npid) #endif return end subroutine get_dimension subroutine get1d_real8(fileName,var_name,out_buff) implicit none integer :: ivar, iret,varid,ncid real*8 out_buff(:) character(len=*), intent(in) :: var_name character(len=*), intent(in) :: fileName iret = nf_open(trim(fileName), NF_NOWRITE, ncid) if (iret .ne. 0) then print*,"failed to open the netcdf file: ",trim(fileName) call hydro_stop("In get1d_real8() - failed to open the netcdf file.") return endif ivar = nf_inq_varid(ncid,trim(var_name), varid) if(ivar .ne. 0) then write(6,*) "Read Variable Error file: ",trim(fileName) write(6,*) "Read Error: could not find ",trim(var_name) call hydro_stop("In get1d_real8() - failed to read netcdf varialbe name. ") end if iret = nf_get_var_double(ncid, varid, out_buff) iret = nf_close(ncid) end subroutine get1d_real8 subroutine get1d_int(fileName,var_name,out_buff) implicit none integer :: ivar, iret,varid,ncid integer out_buff(:) character(len=*), intent(in) :: var_name character(len=*), intent(in) :: fileName iret = nf_open(trim(fileName), NF_NOWRITE, ncid) if (iret .ne. 0) then print*,"FATAL ERROR: Failed to open the netcdf file: ",trim(fileName) call hydro_stop("In get1d_int() - Failed to open the netcdf file") return endif ivar = nf_inq_varid(ncid,trim(var_name), varid) if(ivar .ne. 0) then write(6,*) "Read Variable Error file: ",trim(fileName) write(6,*) "Read Error: could not find ",trim(var_name) call hydro_stop("In get1d_int() - failed to read netcdf variable name.") end if iret = nf_get_var_int(ncid, varid, out_buff) iret = nf_close(ncid) end subroutine get1d_int subroutine getUDMP_area(cell_area) implicit none integer i,j,k, m real, dimension(:,:) :: cell_area do k = 1, LNUMRSL if(LUDRSL(k)%ngrids .gt. 0) then do m = 1, LUDRSL(k)%ngrids LUDRSL(k)%nodeArea(m) = cell_area(LUDRSL(k)%grid_i(m),LUDRSL(k)%grid_j(m)) enddo endif do m = 1, LUDRSL(k)%ncell LUDRSL(k)%cellArea(m) = cell_area(LUDRSL(k)%cell_i(m),LUDRSL(k)%cell_j(m)) enddo basns_area(k) = 0 do m = 1, LUDRSL(k)%ncell basns_area(k) = basns_area(k) + & cell_area(LUDRSL(k)%cell_i(m),LUDRSL(k)%cell_j(m)) * LUDRSL(k)%cellWeight(m) enddo end do end subroutine getUDMP_area subroutine get_basn_area_nhd(inOut) implicit none real, dimension(:) :: inOut real, dimension(gnpid) :: buf #ifdef MPP_LAND call updateLinkV(basns_area, inOut) #else inOut = basns_area #endif end subroutine get_basn_area_nhd subroutine get_nprocs_map(ix,jx,bufi,bufj,nprocs_map,ndata) implicit none integer,dimension(:) :: bufi, bufj,nprocs_map ! integer, allocatable, dimension(:) :: lbufi,lbufj, lmap integer :: ndata, lsize, ix,jx integer, dimension(ix,jx) :: mask integer, allocatable,dimension(:,:) :: gmask integer :: i,j,k, starti,startj, endi,endj, ii,jj, npid, kk #ifdef MPP_LAND mask = my_id + 1 if(my_id .eq. IO_id) allocate(gmask(global_rt_nx, global_rt_ny)) call MPP_LAND_COM_INTEGER(mask,IX,JX,99) call write_IO_rt_int(mask, gmask) if(my_id .eq. io_id ) then nprocs_map = -999 do i = 1, ndata if( (bufi(i) .gt. 0 .and. bufi(i) .le. global_rt_nx) .and. & (bufj(i) .gt. 0 .and. bufj(i) .le. global_rt_ny) ) then nprocs_map(i) = gmask(bufi(i), bufj(i)) if( gmask(bufi(i), bufj(i)) .lt. 0) then write(6,*) "mapping error in gmask : ", bufi(i) ,bufj(i) endif else write(6,*) "no mapping for i,j : ", bufi(i) ,bufj(i) endif end do if(allocated(gmask)) deallocate(gmask) endif #else call hydro_stop("FATAL ERROR in UDMP: Sequential not work.") #endif end subroutine get_nprocs_map end module module_UDMAP