Merge pull request #293 from ESCOMP/develop

v1.2.2

.gitignore

@@ -27,7 +27,15 @@
 *.out
 *.app
 
+# directories
+bin
+cime
+libraries
+lib
+manage_externals
+
 # other files
 *.local
 *bak*
 *tmp*
+

route/build/Makefile

@@ -15,7 +15,7 @@ FC  =
 # Define the compiler exe, e.g., gnu=>gfortran, intel=>ifort, pgi=>pgf90
 FC_EXE =
 
-# Define the compiled executable
+# Define the executable
 EXE =
 
 # Define optional setting
@@ -35,12 +35,12 @@ F_MASTER =
 # Define the NetCDF libraries and path to include files
 ifeq "$(FC)" "gnu"
  NCDF_PATH =
-endif
-ifeq "$(FC)" "intel"
+else ifeq "$(FC)" "intel"
  NCDF_PATH =
-endif
-ifeq "$(FC)" "pgi"
+else ifeq "$(FC)" "pgi"
  NCDF_PATH =
+else
+ $(error FC is not set correctly)
 endif
 
 LIBNETCDF = -Wl,-rpath,$(NCDF_PATH)/lib \
@@ -114,9 +114,15 @@ EXE_PATH = $(F_MASTER)bin
 #========================================================================
 # Define subroutines
 #
+# define utilities
+UTILS = \
+    nrtype.f90 \
+    nr_utility.f90 \
+    ascii_util.f90 \
+    ncio_utils.f90 \
+    gamma_func.f90
 # data types
 DATATYPES = \
-    nrtype.f90 \
     public_var.f90 \
     dataTypes.f90 \
     var_lookup.f90 \
@@ -125,12 +131,6 @@ DATATYPES = \
     globalData.f90 \
     popMetadat.f90 \
     allocation.f90
-# define utilities
-UTILS = \
-    nr_utility.f90 \
-    ascii_util.f90 \
-    ncio_utils.f90 \
-    gamma_func.f90
 # initialization
 INIT = \
     network_topo.f90 \
@@ -158,11 +158,14 @@ CORE = \
     basinUH.f90 \
     irf_route.f90 \
     kwt_route.f90 \
+    dfw_route.f90 \
+    kw_route.f90 \
+    mc_route.f90 \
     main_route.f90 \
     model_setup.f90
 
 # concatanate model subroutines
-TEMP_MODSUB = $(DATATYPES) $(UTILS) $(INIT) $(IO) $(CORE)
+TEMP_MODSUB = $(UTILS) $(DATATYPES) $(INIT) $(IO) $(CORE)
 
 # insert appropriate directory name
 MODSUB = $(patsubst %, $(F_KORE_DIR)%, $(TEMP_MODSUB))
@@ -183,15 +186,18 @@ all: compile install clean
 compile:
 	$(FC_EXE) $(FLAGS) $(MODSUB) $(DRIVER) \
 	$(LIBNETCDF) $(INCNETCDF) -o $(EXE)
+	@echo "Succesfully compiled"
 
 # Remove object files
 clean:
 						rm -f *.o
 						rm -f *.lst
 						rm -f *.mod
 						rm -f *__genmod.f90
+						@echo "Successfully cleaned object files"
 #
 # Copy the executable to the bin directory
 install:
 	@mkdir -p $(EXE_PATH)
 	@mv $(EXE) $(EXE_PATH)
+	@echo "Successfully installed"

route/build/src/accum_runoff.f90

@@ -1,16 +1,19 @@
 MODULE accum_runoff_module
 
 USE nrtype
-USE public_var
-
 ! data type
-USE dataTypes,          only : STRFLX         ! fluxes in each reach
-USE dataTypes,          only : RCHTOPO        ! Network topology
+USE dataTypes,      ONLY: STRFLX         ! fluxes in each reach
+USE dataTypes,      ONLY: RCHTOPO        ! Network topology
+USE dataTypes,      ONLY: subbasin_omp   ! mainstem+tributary data structures
+! global data
+USE public_var,     ONLY: iulog          ! i/o logical unit number
+USE globalData,     ONLY: idxSUM         ! index of accumulation method
+! subroutines: general
+USE model_finalize, ONLY : handle_err
 
 implicit none
 
 private
-
 public::accum_runoff
 
 CONTAINS
@@ -33,9 +36,6 @@ SUBROUTINE accum_runoff(iEns,          & ! input: index of runoff ensemble to be
  !
  ! ----------------------------------------------------------------------------------------
 
- USE dataTypes,      ONLY: subbasin_omp   ! mainstem+tributary data structures
- USE model_finalize, ONLY : handle_err
-
  implicit none
  ! input
  integer(i4b),       intent(in)                 :: iens            ! runoff ensemble index
@@ -114,17 +114,17 @@ END SUBROUTINE accum_runoff
  ! *********************************************************************
  ! subroutine: perform accumulate immediate upstream flow
  ! *********************************************************************
- subroutine accum_qupstream(iEns,       &    ! input: index of runoff ensemble to be processed
+ SUBROUTINE accum_qupstream(iEns,       &    ! input: index of runoff ensemble to be processed
                             segIndex,   &    ! input: index of reach to be processed
                             ixDesire,   &    ! input: reachID to be checked by on-screen pringing
                             NETOPO_in,  &    ! input: reach topology data structure
                             RCHFLX_out, &    ! inout: reach flux data structure
                             ierr, message)   ! output: error control
  implicit none
  ! Input
- INTEGER(I4B), intent(IN)                 :: iEns           ! runoff ensemble to be routed
- INTEGER(I4B), intent(IN)                 :: segIndex       ! segment where routing is performed
- INTEGER(I4B), intent(IN)                 :: ixDesire       ! index of the reach for verbose output
+ integer(i4b), intent(in)                 :: iEns           ! runoff ensemble to be routed
+ integer(i4b), intent(in)                 :: segIndex       ! segment where routing is performed
+ integer(i4b), intent(in)                 :: ixDesire       ! index of the reach for verbose output
  type(RCHTOPO),intent(in),    allocatable :: NETOPO_in(:)   ! River Network topology
  ! inout
  TYPE(STRFLX), intent(inout), allocatable :: RCHFLX_out(:,:)   ! Reach fluxes (ensembles, space [reaches]) for decomposed domains
@@ -137,41 +137,40 @@ subroutine accum_qupstream(iEns,       &    ! input: index of runoff ensemble to
  integer(i4b)                             :: iUps           ! upstream reach index
  integer(i4b)                             :: iRch_ups       ! index of upstream reach in NETOPO
  character(len=strLen)                    :: fmt1,fmt2      ! format string
- character(len=strLen)                    :: cmessage       ! error message from subroutine
 
  ierr=0; message='accum_qupstream/'
 
  ! identify number of upstream segments of the reach being processed
  nUps = size(NETOPO_in(segIndex)%UREACHI)
 
- RCHFLX_out(iEns,segIndex)%UPSTREAM_QI = RCHFLX_out(iEns,segIndex)%BASIN_QR(1)
+ RCHFLX_out(iEns,segIndex)%ROUTE(idxSUM)%REACH_Q = RCHFLX_out(iEns,segIndex)%BASIN_QR(1)
 
  q_upstream = 0._dp
  if (nUps>0) then
 
    do iUps = 1,nUps
      iRch_ups = NETOPO_in(segIndex)%UREACHI(iUps)      !  index of upstream of segIndex-th reach
-     q_upstream = q_upstream + RCHFLX_out(iens,iRch_ups)%UPSTREAM_QI
+     q_upstream = q_upstream + RCHFLX_out(iens,iRch_ups)%ROUTE(idxSUM)%REACH_Q
    end do
 
-   RCHFLX_out(iEns,segIndex)%UPSTREAM_QI = RCHFLX_out(iEns,segIndex)%UPSTREAM_QI + q_upstream
+   RCHFLX_out(iEns,segIndex)%ROUTE(idxSUM)%REACH_Q = RCHFLX_out(iEns,segIndex)%ROUTE(idxSUM)%REACH_Q + q_upstream
 
  endif
 
  ! check
  if(segIndex == ixDesire)then
    write(fmt1,'(A,I5,A)') '(A,1X',nUps,'(1X,I10))'
    write(fmt2,'(A,I5,A)') '(A,1X',nUps,'(1X,F20.7))'
-   write(*,'(2a)') new_line('a'),'** Check upstream discharge accumulation **'
-   write(*,'(a,x,I10,x,I10)') ' Reach index & ID =', segIndex, NETOPO_in(segIndex)%REACHID
-   write(*,'(a)')             ' * upstream reach index (NETOPO_in%UREACH) and discharge (uprflux) [m3/s] :'
-   write(*,fmt1)              ' UREACHK =', (NETOPO_in(segIndex)%UREACHK(iUps), iUps=1,nUps)
-   write(*,fmt2)              ' prflux  =', (RCHFLX_out(iens,NETOPO_in(segIndex)%UREACHI(iUps))%UPSTREAM_QI, iUps=1,nUps)
-   write(*,'(a)')             ' * local area discharge (RCHFLX_out%BASIN_QR(1)) and final discharge (RCHFLX_out%UPSTREAM_QI) [m3/s] :'
-   write(*,'(a,x,F15.7)')     ' RCHFLX_out%BASIN_QR(1) =', RCHFLX_out(iEns,segIndex)%BASIN_QR(1)
-   write(*,'(a,x,F15.7)')     ' RCHFLX_out%UPSTREAM_QI =', RCHFLX_out(iens,segIndex)%UPSTREAM_QI
+   write(iulog,'(2a)') new_line('a'),'** Check upstream discharge accumulation **'
+   write(iulog,'(a,x,I10,x,I10)') ' Reach index & ID =', segIndex, NETOPO_in(segIndex)%REACHID
+   write(iulog,'(a)')             ' * upstream reach index (NETOPO_in%UREACH) and discharge (uprflux) [m3/s] :'
+   write(iulog,fmt1)              ' UREACHK =', (NETOPO_in(segIndex)%UREACHK(iUps), iUps=1,nUps)
+   write(iulog,fmt2)              ' prflux  =', (RCHFLX_out(iens,NETOPO_in(segIndex)%UREACHI(iUps))%ROUTE(idxSUM)%REACH_Q, iUps=1,nUps)
+   write(iulog,'(a)')             ' * local area discharge (RCHFLX_out%BASIN_QR(1)) and final discharge (RCHFLX_out%ROUTE(idxSUM)%REACH_Q) [m3/s] :'
+   write(iulog,'(a,x,F15.7)')     ' RCHFLX_out%BASIN_QR(1) =', RCHFLX_out(iEns,segIndex)%BASIN_QR(1)
+   write(iulog,'(a,x,F15.7)')     ' RCHFLX_out%ROUTE(idxSUM)%REACH_Q =', RCHFLX_out(iens,segIndex)%ROUTE(idxSUM)%REACH_Q
  endif
 
- end subroutine accum_qupstream
+ END SUBROUTINE accum_qupstream
 
 END MODULE accum_runoff_module

route/build/src/ascii_util.f90

@@ -1,11 +1,18 @@
-module ascii_util_module
+MODULE ascii_util_module
+
 USE nrtype
+
 implicit none
+
 private
+
 public::file_open
 public::split_line
 public::get_vlines
-contains
+public::lower
+public::upper
+
+CONTAINS
 
  ! **********************************************************************************************
  ! new subroutine: get unused file unit (modified from DMSL)
@@ -204,4 +211,44 @@ subroutine get_vlines(unt,vlines,err,message)
  end subroutine get_vlines
 
 
-end module ascii_util_module
+  FUNCTION upper(strIn) RESULT(strOut)
+    ! convert string to upper-case
+    ! only ASCII character code works
+    implicit none
+
+    character(*), intent(in) :: strIn
+    character(len(strIn))    :: strOut
+    integer, parameter       :: DUC = ichar('A') - ichar('a')
+    character                :: ch
+    integer                  :: i
+
+    do i = 1, len(strIn)
+      ch = strIn(i:i)
+      if (ch>='a' .and. ch<='z') ch = char(ichar(ch)+DUC)
+      strOut(i:i) = ch
+    end do
+
+  END FUNCTION upper
+
+
+  pure FUNCTION lower(strIn) RESULT(strOut)
+    ! convert string to lower-case
+    ! only ASCII character code works
+     implicit none
+
+    character(*), intent(in)  :: strIn
+    character(len(strIn))     :: strOut
+    integer, parameter        :: DUC = ichar('A') - ichar('a')
+    character                 :: ch
+    integer                   :: i
+
+    do i = 1,len(strIn)
+      ch = strIn(i:i)
+      if (ch>='A' .and. ch<='Z') ch = char(ichar(ch)-DUC)
+      strOut(i:i) = ch
+    end do
+
+  END FUNCTION lower
+
+
+END MODULE ascii_util_module

route/build/src/dataTypes.f90

@@ -202,45 +202,78 @@ module dataTypes
   LOGICAL(LGT)                               :: USRTAKE  ! .TRUE. if user takes from reach, .FALSE. otherwise
  end type RCHTOPO
 
- ! ---------- kinematic wave states (collection of particles) ---------------------------------
+ ! ---------- reach states --------------------------------------------------------------------
 
+ !---------- Lagrangian kinematic wave states (collection of particles) ---------------------------------
  ! Individual flow particles
  ! NOTE: type could possibly be private
- TYPE, public :: FPOINT
-  REAL(DP)                             :: QF       ! Flow
-  REAL(DP)                             :: QM       ! Modified flow
-  REAL(DP)                             :: TI       ! initial time of point in reach
-  REAL(DP)                             :: TR       ! time point expected to exit reach
-  LOGICAL(LGT)                         :: RF       ! routing flag (T if point has exited)
- END TYPE FPOINT
+ type, public :: FPOINT
+  real(dp)                             :: QF       ! Flow
+  real(dp)                             :: QM       ! Modified flow
+  real(dp)                             :: TI       ! initial time of point in reach
+  real(dp)                             :: TR       ! time point expected to exit reach
+  logical(lgt)                         :: RF       ! routing flag (T if point has exited)
+ end type FPOINT
 
  ! Collection of flow points within a given reach
- TYPE, public :: KREACH
-  TYPE(FPOINT),allocatable             :: KWAVE(:)
- END TYPE KREACH
+ type, public :: kwtRCH
+  type(FPOINT),allocatable             :: KWAVE(:)
+ end type kwtRCH
 
  ! ---------- irf states (future flow series ) ---------------------------------
-
  ! Future flow series
- TYPE, public :: IRFREACH
-  REAL(DP), allocatable                :: qfuture(:)    ! runoff volume in future time steps for IRF routing (m3/s)
- END TYPE IRFREACH
+ type, public :: irfRCH
+  real(dp), allocatable   :: qfuture(:)    ! runoff volume in future time steps for IRF routing (m3/s)
+ end type irfRCH
+
+ ! ---------- computational molecule ---------------------------------
+ type, public :: cMolecule
+   integer(i4b)           :: KW_ROUTE
+   integer(i4b)           :: MC_ROUTE
+   integer(i4b)           :: DW_ROUTE
+ end type cMolecule
+
+ type, public :: SUBRCH
+   real(dp), allocatable  :: Q(:)        ! Discharge at sub-reaches at current step (m3/s)
+   real(dp), allocatable  :: A(:)        ! Flow area at sub-reach at current step (m2)
+   real(dp), allocatable  :: H(:)        ! Flow height at sub-reach at current step (m)
+ end type SUBRCH
+
+ type, public :: kwRch
+   type(SUBRCH)    :: molecule
+ end type kwRCH
+
+ type, public :: mcRch
+   type(SUBRCH)    :: molecule
+ end type mcRCH
+
+ type, public :: dwRch
+   type(SUBRCH)    :: molecule
+ end type dwRCH
+
+ type, public :: STRSTA
+   type(irfRCH)    :: IRF_ROUTE
+   type(kwtRCH)    :: LKW_ROUTE
+   type(kwRCH)     :: KW_ROUTE
+   type(mcRCH)     :: MC_ROUTE
+   type(dwRCH)     :: DW_ROUTE
+ end type STRSTA
 
  ! ---------- reach fluxes --------------------------------------------------------------------
+ type, public :: fluxes
+   real(dp)        :: REACH_Q
+   real(dp)        :: REACH_VOL(0:1)
+ end type fluxes
 
  ! fluxes in each reach
  TYPE, public :: STRFLX
-  REAL(DP), allocatable                :: QFUTURE(:)        ! runoff volume in future time steps (m3/s)
-  REAL(DP), allocatable                :: QFUTURE_IRF(:)    ! runoff volume in future time steps for IRF routing (m3/s)
-  REAL(DP)                             :: BASIN_QI          ! instantaneous runoff volume from the local basin (m3/s)
-  REAL(DP)                             :: BASIN_QR(0:1)     ! routed runoff volume from the local basin (m3/s)
-  REAL(DP)                             :: BASIN_QR_IRF(0:1) ! routed runoff volume from all the upstream basin (m3/s)
-  REAL(DP)                             :: REACH_Q           ! time-step average streamflow (m3/s)
-  REAL(DP)                             :: REACH_Q_IRF       ! time-step average streamflow (m3/s) from IRF routing
-  REAL(DP)                             :: UPSTREAM_QI       ! sum of upstream streamflow (m3/s)
-  REAL(DP)                             :: REACH_VOL(0:1)    ! volume of water at a reach [m3]
-  REAL(DP)                             :: TAKE              ! average take
-  logical(lgt)                         :: CHECK_IRF         ! .true. if the reach is routed
+  real(dp), allocatable                :: QFUTURE(:)        ! runoff volume in future time steps (m3/s)
+  real(dp), allocatable                :: QFUTURE_IRF(:)    ! runoff volume in future time steps for IRF routing (m3/s)
+  real(dp)                             :: BASIN_QI          ! instantaneous runoff volume from the local basin (m3/s)
+  real(dp)                             :: BASIN_QR(0:1)     ! routed runoff volume from the local basin (m3/s)
+  real(dp)                             :: BASIN_QR_IRF(0:1) ! routed runoff volume from all the upstream basin (m3/s)
+  type(fluxes), allocatable            :: ROUTE(:)          ! reach fluxes and states for each routing method
+  real(dp)                             :: TAKE              ! average take
  ENDTYPE STRFLX
 
  ! ---------- lake data types -----------------------------------------------------------------