diff --git a/trunk/NDHMS/Data_Rec/module_namelist.F b/trunk/NDHMS/Data_Rec/module_namelist.F index 310ed99c2..37c0b4c5d 100644 --- a/trunk/NDHMS/Data_Rec/module_namelist.F +++ b/trunk/NDHMS/Data_Rec/module_namelist.F @@ -713,13 +713,13 @@ subroutine rt_nlst_check(nlst) if (.not. fileExists) call hydro_stop('hydro.namelist ERROR: route_lake_f not found.') endif ! Only allow lakes to be ran with gridded routing or NWM routing - if(len(trim(nlst%route_lake_f)) .ne. 0) then - if(nlst%channel_option .ne. 3) then - if(nlst%UDMP_OPT .ne. 1) then - call hydro_stop('hydro.namelist ERROR: Currently lakes only work with gridded channel routing or UDMP=1. Please change your namelist settings.') - endif - endif - endif +! if(len(trim(nlst%route_lake_f)) .ne. 0) then +! if(nlst%channel_option .ne. 3) then +! if(nlst%UDMP_OPT .ne. 1) then +! call hydro_stop('hydro.namelist ERROR: Currently lakes only work with gridded channel routing or UDMP=1. Please change your namelist settings.') +! endif +! endif +! endif if((nlst%channel_option .eq. 3) .and. (nlst%compound_channel)) then call hydro_stop("Compound channel option not available for diffusive wave routing. ") diff --git a/trunk/NDHMS/HYDRO_drv/module_HYDRO_drv.F b/trunk/NDHMS/HYDRO_drv/module_HYDRO_drv.F index ebb4eea73..035387e45 100644 --- a/trunk/NDHMS/HYDRO_drv/module_HYDRO_drv.F +++ b/trunk/NDHMS/HYDRO_drv/module_HYDRO_drv.F @@ -449,6 +449,7 @@ subroutine HYDRO_out(did, rstflag) call output_chrtout_grd_NWM(did,nlst_rt(did)%igrid) endif endif + if (RT_DOMAIN(did)%NLAKES.gt.0) then if(nlst_rt(did)%io_form_outputs .ne. 0) then ! Output lakes in NWM format @@ -984,16 +985,16 @@ subroutine driveChannelRouting(did) RT_DOMAIN(did)%QINFLOWBASE, RT_DOMAIN(did)%CHANXI, & RT_DOMAIN(did)%CHANYJ, nlst_rt(did)%channel_option, & RT_DOMAIN(did)%RETDEP_CHAN, RT_DOMAIN(did)%NLINKSL, RT_DOMAIN(did)%LINKID, & - RT_DOMAIN(did)%node_area & + RT_DOMAIN(did)%node_area, RT_DOMAIN(did)%LAKEIDX, & #ifdef MPP_LAND - ,RT_DOMAIN(did)%lake_index,RT_DOMAIN(did)%link_location,& + RT_DOMAIN(did)%lake_index,RT_DOMAIN(did)%link_location,& RT_DOMAIN(did)%mpp_nlinks,RT_DOMAIN(did)%nlinks_index, & - RT_DOMAIN(did)%yw_mpp_nlinks & - , RT_DOMAIN(did)%LNLINKSL,RT_DOMAIN(did)%LLINKID & - , rt_domain(did)%gtoNode,rt_domain(did)%toNodeInd,rt_domain(did)%nToInd & + RT_DOMAIN(did)%yw_mpp_nlinks, & + RT_DOMAIN(did)%LNLINKSL,RT_DOMAIN(did)%LLINKID, & + rt_domain(did)%gtoNode,rt_domain(did)%toNodeInd,rt_domain(did)%nToInd, & #endif - , rt_domain(did)%CH_LNKRT_SL & - ,nlst_rt(did)%GwBaseSwCRT, gw2d(did)%ho, gw2d(did)%qgw_chanrt, & + rt_domain(did)%CH_LNKRT_SL, & + nlst_rt(did)%GwBaseSwCRT, gw2d(did)%ho, gw2d(did)%qgw_chanrt, & nlst_rt(did)%gwChanCondSw, nlst_rt(did)%gwChanCondConstIn, & nlst_rt(did)%gwChanCondConstOut,rt_domain(did)%velocity & ) diff --git a/trunk/NDHMS/Routing/module_HYDRO_io.F b/trunk/NDHMS/Routing/module_HYDRO_io.F index 7eb06853d..5fe28bb90 100644 --- a/trunk/NDHMS/Routing/module_HYDRO_io.F +++ b/trunk/NDHMS/Routing/module_HYDRO_io.F @@ -1004,7 +1004,7 @@ SUBROUTINE READ_ROUTEDIM(IXRT,JXRT,route_chan_f,route_link_f, & !!-- no longer find the lakes from the 2-d hi res grid !DJG inv do j=jxrt,1,-1 ! follwoing is modified by Wei Yu 03/24/2017 - if(UDMP_OPT .eq. 0) then + if(UDMP_OPT .eq. 0 .and. channel_option .eq. 3) then NLAKES = 0 do j=1,jxrt do i = 1,ixrt diff --git a/trunk/NDHMS/Routing/module_RT.F b/trunk/NDHMS/Routing/module_RT.F index 5b0ad441f..86f9e2018 100644 --- a/trunk/NDHMS/Routing/module_RT.F +++ b/trunk/NDHMS/Routing/module_RT.F @@ -532,7 +532,7 @@ subroutine getChanDim(did) rt_domain(did)%GNLINKSL = 1 rt_domain(did)%NLINKSL = 1 endif - if(nlst_rt(did)%UDMP_OPT .eq. 1) & + if(nlst_rt(did)%UDMP_OPT .eq. 1 .or. nlst_rt(did)%channel_option .eq. 1 .or. nlst_rt(did)%channel_option .eq. 2) & call read_NSIMLAKES(rt_domain(did)%NLAKES,nlst_rt(did)%route_lake_f) call rt_allocate(did,rt_domain(did)%ix,rt_domain(did)%jx,& @@ -588,7 +588,7 @@ subroutine getChanDim(did) endif -if(nlst_rt(did)%UDMP_OPT .eq. 1) then +if(nlst_rt(did)%UDMP_OPT .eq. 1 .or. nlst_rt(did)%channel_option .eq. 1 .or. nlst_rt(did)%channel_option .eq. 2) then call read_NSIMLAKES(rt_domain(did)%NLAKES,nlst_rt(did)%route_lake_f) endif @@ -1460,7 +1460,7 @@ subroutine LandRT_ini(did) endif !End if channel option eq 3; else; -! Initialize Lake Elevations for Gridded and NWM routing. +! Initialize Lake Elevations for Gridded and NWM routing. AND REACH Community? do j=1,rt_domain(did)%NLAKES rt_domain(did)%RESHT(j) = rt_domain(did)%ORIFICEE(j) + & ((rt_domain(did)%LAKEMAXH(j) - rt_domain(did)%ORIFICEE(j) )* rt_domain(did)%ELEVLAKE(j)) diff --git a/trunk/NDHMS/Routing/module_channel_routing.F b/trunk/NDHMS/Routing/module_channel_routing.F index 404aab25f..e1fdcc18f 100644 --- a/trunk/NDHMS/Routing/module_channel_routing.F +++ b/trunk/NDHMS/Routing/module_channel_routing.F @@ -699,7 +699,7 @@ Subroutine drive_CHANNEL(latval,lonval,KT, IXRT,JXRT, SUBRTSWCRT, & RESHT, HRZAREA, LAKEMAXH, WEIRH, WEIRC, WEIRL, ORIFICEC, ORIFICEA, & ORIFICEE, ZELEV, CVOL, NLAKES, QLAKEI, QLAKEO, LAKENODE, & dist, QINFLOWBASE, CHANXI, CHANYJ, channel_option, RETDEP_CHAN, & - NLINKSL, LINKID, node_area & + NLINKSL, LINKID, node_area, lake_lookup & #ifdef MPP_LAND , lake_index,link_location,mpp_nlinks,nlinks_index,yw_mpp_nlinks & , LNLINKSL, LLINKID & @@ -792,8 +792,10 @@ Subroutine drive_CHANNEL(latval,lonval,KT, IXRT,JXRT, SUBRTSWCRT, & REAL, DIMENSION(NLAKES) :: QLLAKE !-- lateral inflow to lake in diffusion scheme REAL*8, DIMENSION(NLAKES) :: QLLAKE8 !-- lateral inflow to lake in diffusion scheme + integer, intent(in), dimension(:) :: lake_lookup !-- inverse lake index for k->lake mapping + !-- Local Variables - INTEGER :: i,j,k,t,m,jj,kk,KRT,node + INTEGER :: i,j,k,t,m,jj,kk,KRT,node,l_idx INTEGER :: DT_STEPS !-- number of timestep in routing REAL :: Qup,Quc !--Q upstream Previous, Q Upstream Current, downstream Previous REAL :: bo !--critical depth, bnd outflow just for testing @@ -812,19 +814,28 @@ Subroutine drive_CHANNEL(latval,lonval,KT, IXRT,JXRT, SUBRTSWCRT, & real ywtmp(ixrt,jxrt) integer LNLINKSL integer, dimension(LNLINKSL) :: LLINKID - real*8, dimension(LNLINKSL) :: LQLateral -! real*4, dimension(LNLINKSL) :: LQLateral + real(kind=8), dimension(LNLINKSL) :: LQLateral integer, dimension(:) :: toNodeInd integer, dimension(:,:) :: gtoNode integer :: nToNodeInd real, dimension(nToNodeInd,2) :: gQLINK + real, allocatable,dimension(:) :: tmpQLAKEO, tmpQLAKEI, tmpRESHT #else - real*8, dimension(NLINKS) :: LQLateral !--lateral flow + real(kind=8), dimension(NLINKS) :: LQLateral !--lateral flow #endif integer flag - integer :: n, kk2, nt, nsteps ! tmp + integer :: n, kk2, nt, nsteps ! tmp +#ifdef MPP_LAND + if(my_id == io_id) then +#endif + allocate(tmpQLAKEO(NLAKES)) + allocate(tmpQLAKEI(NLAKES)) + allocate(tmpRESHT(NLAKES)) +#ifdef MPP_LAND + endif +#endif QLAKEIP = 0 QLAKEI8 = 0 HLINKTMP = 0 @@ -932,6 +943,15 @@ Subroutine drive_CHANNEL(latval,lonval,KT, IXRT,JXRT, SUBRTSWCRT, & !---------- route other reaches, with upstream inflow tmpQlink = 0.0 +#ifdef MPP_LAND + if(my_id .eq. io_id) then +#endif + tmpQLAKEO = QLAKEO + tmpQLAKEI = QLAKEI + tmpRESHT = RESHT +#ifdef MPP_LAND + endif +#endif do k = 1,NLINKSL ! if (ORDER(k) .gt. 1 ) then !-- exclude first order stream Quc = 0.0 @@ -962,29 +982,38 @@ Subroutine drive_CHANNEL(latval,lonval,KT, IXRT,JXRT, SUBRTSWCRT, & end do ! do m #endif - if(TYPEL(k) .eq. 1) then !--link is a reservoir - - ! CALL LEVELPOOL(1,QLINK(k,1), Qup, QLINK(k,1), QLINK(k,2), & - ! QLateral(k), DT, RESHT(k), HRZAREA(k), LAKEMAXH(k), & - ! WEIRC(k), WEIRL(k),ORIFICEE(k), ORIFICEC(k), ORIFICEA(k)) - - elseif (channel_option .eq. 1) then !muskingum routing + if(TYPEL(k) == 1) then !--link is a reservoir + l_idx = lake_lookup(k) + if (l_idx >= 0) then !-- -999 if not a reservoir in the lookup table (belt-and-suspenders check) + CALL LEVELPOOL(l_idx,Quc, Qup, QLINK(k,2), QLateral(k), DT, & + RESHT(l_idx), HRZAREA(l_idx), WEIRH(l_idx), LAKEMAXH(l_idx), & + WEIRC(l_idx), WEIRL(l_idx), ORIFICEE(l_idx), ORIFICEC(l_idx), ORIFICEA(l_idx)) + QLAKEO(l_idx) = QLINK(k,2) !save outflow to lake + QLAKEI(l_idx) = Quc !save inflow to lake + end if + elseif (channel_option .eq. 1) then !muskingum routing Km = MUSK(k) X = MUSX(k) tmpQLINK(k,2) = MUSKING(k,Qup,(Quc+QLateral(k)),QLINK(k,1),DTRT_CH,Km,X) !upstream plust lateral inflow - elseif (channel_option .eq. 2) then ! muskingum cunge + elseif (channel_option .eq. 2) then ! muskingum cunge call SUBMUSKINGCUNGE(tmpQLINK(k,2), velocity(k), LINKID(k), & Qup,Quc, QLINK(k,1), QLateral(k), DTRT_CH, So(k), & CHANLEN(k), MannN(k), ChSSlp(k), Bw(k), Tw(k),Tw_CC(k), n_CC(k), HLINK(k) ) - else + else print *, "FATAL ERROR: no channel option selected" call hydro_stop("In drive_CHANNEL() - no channel option selected") - endif + endif ! endif !!! order(1) .ne. 1 end do !--k links +#ifdef MPP_LAND + call updateLake_seq(RESHT,nlakes,tmpRESHT) + call updateLake_seq(QLAKEO,nlakes,tmpQLAKEO) + call updateLake_seq(QLAKEI,nlakes,tmpQLAKEI) +#endif + !yw check ! gQLINK = 0.0 ! call ReachLS_write_io(tmpQLINK(:,2), gQLINK(:,2)) @@ -997,7 +1026,7 @@ Subroutine drive_CHANNEL(latval,lonval,KT, IXRT,JXRT, SUBRTSWCRT, & ! endif do k = 1, NLINKSL - if(TYPEL(k) .ne. 1) then + if(TYPEL(k) .ne. 2) then QLINK(k,2) = tmpQLINK(k,2) endif QLINK(k,1) = QLINK(k,2) !assing link flow of current to be previous for next time step @@ -1416,8 +1445,14 @@ Subroutine drive_CHANNEL(latval,lonval,KT, IXRT,JXRT, SUBRTSWCRT, & #endif if (KT .eq. 1) KT = KT + 1 - +#ifdef MPP_LAND + if (my_id == io_id) then + if(allocated(tmpRESHT)) deallocate(tmpRESHT) + if(allocated(tmpQLAKEO)) deallocate(tmpQLAKEO) + if(allocated(tmpQLAKEI)) deallocate(tmpQLAKEI) + endif +#endif end subroutine drive_CHANNEL ! ----------------------------------------------------------------