Merge pull request #169 from nmizukami/master

v1.2.1

README.md

@@ -1,4 +1,4 @@
-[![DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.3953155.svg)](https://doi.org/10.5281/zenodo.3953155)
+[![DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.4395155.svg)](https://doi.org/10.5281/zenodo.4395155)
 [![Documentation Status](https://readthedocs.org/projects/mizuroute/badge/?version=master)](https://mizuroute.readthedocs.io/en/master/?badge=master)
 
 # mizuRoute

docs/source/Control_file.rst

@@ -1,10 +1,10 @@
 Control file
 ============
 
-Control file is a simple text file, mainly defining model control such as simulation time, file name and locations, routing option etc. 
+Control file is a simple text file, defining various model controls such as simulation time, file names and locations, routing options etc. 
 Variables in control file are read in the beginning of the code (see ``./build/src/read_control.f90``) and 
 saved in fortran variable specified by tag (inside <> in table) and as public variables (see ``./build/src/public_var.f90``) . 
-Some of such public varialbes have some default values, but most of them are not defined.
+Some of control varialbes have their default values, but most of them are not defined.
 Those undefined variables need to be defined in control file.   
 Other variables in supplement table have their default values but can be also included in control file to overwrite the values. 
 The order of variables in the control file does not matter. However, grouping variables into similar themes is recommended for readibility. 
@@ -19,6 +19,7 @@ Some of rules
 * Format: <tag>    variable    ! comments
 * tag is Fortran variable name and cannot be changed and have to be enclosed by <>
 * need ! after variable, otherwise getting error.
+* Do not leave any lines empty in control file
 
 
 The following variables (not pre-defined in the code) need to be defined in control file.
@@ -86,12 +87,6 @@ The following variables (not pre-defined in the code) need to be defined in cont
 +--------+------------------------+-------------------------------------------------------------------------------------------+
 |   2,3  | <dname_data_remap>     | dimension name for data                                                                   |
 +--------+------------------------+-------------------------------------------------------------------------------------------+
-| 1,2,3  | <restart_write>        | restart ouput timing options. N[n]ever, L[l]ast, S[s]pecified.                            | 
-+--------+------------------------+-------------------------------------------------------------------------------------------+
-| 1,2,3  | <restart_date>         | specified restart date in yyyy-mm-dd (hh:mm:ss) if <restart_write> = "Specified"          | 
-+--------+------------------------+-------------------------------------------------------------------------------------------+
-| 1,2,3  | <fname_state_in>       | input restart netCDF name. If not specified, simulation start with cold start             | 
-+--------+------------------------+-------------------------------------------------------------------------------------------+
 | 1,2,3  | <route_opt>            | option for routing schemes 0-> both, 1->IRF, 2->KWT, otherwise error                      |
 +--------+------------------------+-------------------------------------------------------------------------------------------+
 
@@ -120,6 +115,8 @@ Variables that have default values but can be overwritten
 +------------------------+------------------------+--------------------------------------------------------------------------+
 | <time_units>           | From runoff input      | specified time units <unit> since yyyy-mm-dd (hh:mm:ss). See note 4      |
 +------------------------+------------------------+--------------------------------------------------------------------------+
+| <netcdf_format>        | netcdf4                | netcdf format for output netcdf. other options: classic, 64bit_offset.   |
++------------------------+------------------------+--------------------------------------------------------------------------+
 
 1. River network subset mode. 
 
@@ -138,10 +135,42 @@ Often case, river network data has different variable names than defaults. In th
 See :doc:`River parameters <seg_hru_param>`.   
 
 
+Restart options 
+---------------------
+
+mizuRoute does not write restart netCDF as default. The following control variables are used to control restart dropoff timing and use restart file for continuous run from the previous simulations.
+The restart file is written at previous time step to the specified time. In other words, if ``Specified`` is used for <restart_write> and ``1981-01-01-00000`` is specified in <restart_date>, mizuRoute writes restart file
+at ``1980-12-31 00:00:00`` for daily time step. The restart file name uses the time stamp at user specified timing. ``Annual``, ``Monthly``, ``Daily`` options also follow This convention. 
+
+The restart file name convension:  <case_name>.r.yyyy-mm-dd-sssss.nc 
+
+
++---------------------+---------------------------------------------------------------------------------------------------------+
+| tag                 | Description                                                                                             |
++=====================+=========================================================================================================+
+| <restart_dir>       | directory for restart files. defualt is <output_dir>                                                    | 
++---------------------+---------------------------------------------------------------------------------------------------------+
+| <restart_write>     | restart ouput options. N[n]ever (default), L[l]ast, S[s]pecified, Annual, M[m]onthly, D[d]aily.         | 
++---------------------+---------------------------------------------------------------------------------------------------------+
+| <restart_date>      | restart time in yyyy-mm-dd (hh:mm:ss). required if <restart_write> = "Specified"                        | 
++---------------------+---------------------------------------------------------------------------------------------------------+
+| <restart_month>     | periodic restart month (default 1). Effective if <restart_write>="Annual"                               | 
++---------------------+---------------------------------------------------------------------------------------------------------+
+| <restart_day>       | periodic restart day (default 1). Effective if <restart_write>="Annual" or "Monthly"                    | 
++---------------------+---------------------------------------------------------------------------------------------------------+
+| <restart_hour>      | periodic restart hour (default 0). Effective if <restart_write>="Annual", "Monthly", or "Daily"         | 
++---------------------+---------------------------------------------------------------------------------------------------------+
+| <fname_state_in>    | input restart netCDF name. If not specified, simulation start with cold start                           | 
++---------------------+---------------------------------------------------------------------------------------------------------+
+
+
 Output variables
 ---------------------
 
 The following variables, besides time, basinID (RN_hru ID) and reachID can be output in netCDF. Users can control which variables are output by setting <variable_name> to T or F in control file. All the variables are set to T by default.
+The output file name includes a timie stamp at the first time step.  
+
+The output file name convension:  <case_name>.h.yyyy-mm-dd-sssss.nc
 
 
 +------------------------+------------------------------------------------------------------------------------------------+

route/build/Makefile

@@ -120,16 +120,24 @@ DATATYPES = \
     public_var.f90 \
     dataTypes.f90 \
     var_lookup.f90 \
+    time_utils.f90 \
+    datetime_data.f90 \
     globalData.f90 \
     popMetadat.f90 \
     allocation.f90
 # define utilities
 UTILS = \
     nr_utility.f90 \
     ascii_util.f90 \
-    time_utils.f90 \
     ncio_utils.f90 \
     gamma_func.f90
+# initialization
+INIT = \
+    network_topo.f90 \
+    process_param.f90 \
+    process_ntopo.f90 \
+    pfafstetter.f90 \
+    domain_decomposition.f90
 # read/write files
 IO = \
     remap.f90 \
@@ -143,16 +151,9 @@ IO = \
     read_restart.f90 \
     write_restart.f90 \
     write_simoutput.f90
-# initialization
-INIT = \
-    process_time.f90 \
-    network_topo.f90 \
-    process_param.f90 \
-    process_ntopo.f90 \
-    pfafstetter.f90 \
-    domain_decomposition.f90
 # CORE
 CORE = \
+    model_finalize.f90 \
     accum_runoff.f90 \
     basinUH.f90 \
     irf_route.f90 \
@@ -161,7 +162,7 @@ CORE = \
     model_setup.f90
 
 # concatanate model subroutines
-TEMP_MODSUB = $(DATATYPES) $(UTILS) $(IO) $(INIT) $(CORE)
+TEMP_MODSUB = $(DATATYPES) $(UTILS) $(INIT) $(IO) $(CORE)
 
 # insert appropriate directory name
 MODSUB = $(patsubst %, $(F_KORE_DIR)%, $(TEMP_MODSUB))

route/build/src/accum_runoff.f90

@@ -33,7 +33,8 @@ SUBROUTINE accum_runoff(iEns,          & ! input: index of runoff ensemble to be
  !
  ! ----------------------------------------------------------------------------------------
 
- USE dataTypes,         ONLY: subbasin_omp   ! mainstem+tributary data structures
+ USE dataTypes,      ONLY: subbasin_omp   ! mainstem+tributary data structures
+ USE model_finalize, ONLY : handle_err
 
  implicit none
  ! input
@@ -56,14 +57,9 @@ SUBROUTINE accum_runoff(iEns,          & ! input: index of runoff ensemble to be
  integer(i4b)                                   :: iTrib, ix       ! loop indices
  logical(lgt), allocatable                      :: doRoute(:)      ! logical to indicate which reaches are processed
  character(len=strLen)                          :: cmessage        ! error message from subroutines
- integer*8                                      :: cr                   ! rate
- integer*8                                      :: startTime,endTime    ! date/time for the start and end of the initialization
- real(dp)                                       :: elapsedTime          ! elapsed time for the process
 
  ierr=0; message='accum_runoff/'
- call system_clock(count_rate=cr)
 
- ! check
  if (size(NETOPO_in)/=size(RCHFLX_out(iens,:))) then
   ierr=20; message=trim(message)//'sizes of NETOPO and RCHFLX mismatch'; return
  endif
@@ -84,8 +80,6 @@ SUBROUTINE accum_runoff(iEns,          & ! input: index of runoff ensemble to be
 
  nDom = size(river_basin)
 
- call system_clock(startTime)
-
  do ix = 1,nDom
    ! 1. Route tributary reaches (parallel)
    ! compute the sum of all upstream runoff at each point in the river network
@@ -107,17 +101,13 @@ SUBROUTINE accum_runoff(iEns,          & ! input: index of runoff ensemble to be
        if (.not. doRoute(jSeg)) cycle
 
        call accum_qupstream(iens, jSeg, ixDesire, NETOPO_in, RCHFLX_out, ierr, cmessage)
-       !if(ierr/=0)then; message=trim(message)//trim(cmessage); return; endif
+       if(ierr/=0) call handle_err(ierr, trim(message)//trim(cmessage))
 
      end do
    end do
 !$OMP END PARALLEL DO
 
- end do ! looping through stream segments
-
- call system_clock(endTime)
- elapsedTime = real(endTime-startTime, kind(dp))/real(cr)
- !write(*,"(A,1PG15.7,A)") '  elapsed-time [routing/accum] = ', elapsedTime, ' s'
+ end do
 
  END SUBROUTINE accum_runoff
 
@@ -142,40 +132,44 @@ subroutine accum_qupstream(iEns,       &    ! input: index of runoff ensemble to
  integer(i4b), intent(out)                :: ierr           ! error code
  character(*), intent(out)                :: message        ! error message
  ! Local variables to
- real(dp), allocatable                    :: uprflux(:)     ! upstream Reach fluxes
+ real(dp)                                 :: q_upstream     ! upstream Reach fluxes
  integer(i4b)                             :: nUps           ! number of upstream segment
  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
 
- ! initialize error control
  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)
 
+ q_upstream = 0._dp
  if (nUps>0) then
 
-   allocate(uprflux(nUps), stat=ierr, errmsg=cmessage)
-   if(ierr/=0)then; message=trim(message)//trim(cmessage)//': uprflux'; return; endif
-
    do iUps = 1,nUps
      iRch_ups = NETOPO_in(segIndex)%UREACHI(iUps)      !  index of upstream of segIndex-th reach
-     uprflux(iUps) = RCHFLX_out(iens,iRch_ups)%UPSTREAM_QI
+     q_upstream = q_upstream + RCHFLX_out(iens,iRch_ups)%UPSTREAM_QI
    end do
 
-   RCHFLX_out(iEns,segIndex)%UPSTREAM_QI = RCHFLX_out(iEns,segIndex)%UPSTREAM_QI + sum(uprflux)
+   RCHFLX_out(iEns,segIndex)%UPSTREAM_QI = RCHFLX_out(iEns,segIndex)%UPSTREAM_QI + q_upstream
 
  endif
 
  ! check
- if(NETOPO_in(segIndex)%REACHIX == ixDesire)then
-  print*, 'CHECK ACCUM_RUNOFF'
-  print*, ' UREACHK, uprflux = ', (NETOPO_in(segIndex)%UREACHK(iUps), uprflux(iUps), iUps=1,nUps)
-  print*, ' RCHFLX_out(iEns,segIndex)%BASIN_QR(1) = ', RCHFLX_out(iEns,segIndex)%BASIN_QR(1)
-  print*, ' RCHFLX_out%UPSTREAM_QI = ', RCHFLX_out(iens,segIndex)%UPSTREAM_QI
+ 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
  endif
 
  end subroutine accum_qupstream

route/build/src/basinUH.f90

@@ -36,12 +36,8 @@ SUBROUTINE IRF_route_basin(iens,          & ! input: ensemble index
  integer(i4b)                             :: iSeg            ! reach loop indix
  logical(lgt),               allocatable  :: doRoute(:)      ! logical to indicate which reaches are processed
  character(len=strLen)                    :: cmessage        ! error message from subroutines
- integer*8                                :: cr                   ! rate
- integer*8                                :: startTime,endTime    ! date/time for the start and end of the initialization
- real(dp)                                 :: elapsedTime          ! elapsed time for the process
 
  ierr=0; message='IRF_route_basin/'
- call system_clock(count_rate=cr)
 
  nSeg = size(RCHFLX_out(iens,:))
 
@@ -56,7 +52,6 @@ SUBROUTINE IRF_route_basin(iens,          & ! input: ensemble index
   doRoute(:) = .true.
  endif
 
- call system_clock(startTime)
 !$OMP PARALLEL DO schedule(dynamic,1)   &
 !$OMP          private(iSeg)            & ! loop index
 !$OMP          private(ierr, cmessage)  & ! private for a given thread
@@ -75,10 +70,6 @@ SUBROUTINE IRF_route_basin(iens,          & ! input: ensemble index
  end do
 !$OMP END PARALLEL DO
 
- call system_clock(endTime)
- elapsedTime = real(endTime-startTime, kind(dp))/real(cr)
-! write(*,"(A,1PG15.7,A)") '  elapsed-time [routing/irf_hru] = ', elapsedTime, ' s'
-
  END SUBROUTINE IRF_route_basin