diff --git a/route/build/src/basinUH.f90 b/route/build/src/basinUH.f90 index 9079a9cd..8daf23aa 100644 --- a/route/build/src/basinUH.f90 +++ b/route/build/src/basinUH.f90 @@ -125,13 +125,15 @@ SUBROUTINE hru_irf(iSeg, & ! input: index of segment to be processed if(ierr/=0)then; message=trim(message)//trim(cmessage); return; endif if(tracer) then - ! perform river network UH routing - call irf_conv(FRAC_FUTURE_local, & ! input: unit hydrograph - RCHFLX_out(iSeg)%BASIN_solute_inst, & ! input: upstream fluxes - RCHFLX_out(iSeg)%solute_future, & ! inout: updated solute future time series - RCHFLX_out(iSeg)%BASIN_solute, & ! inout: updated fluxes at reach - ierr, message) ! output: error control - if(ierr/=0)then; message=trim(message)//trim(cmessage); return; endif + do iTracer=1,nTracer + ! perform river network UH routing + call irf_conv(FRAC_FUTURE_local, & ! input: unit hydrograph + RCHFLX_out(iSeg)%BASIN_solute_inst(iTracer), & ! input: upstream fluxes + RCHFLX_out(iSeg)%solute_future, & ! inout: updated solute future time series + RCHFLX_out(iSeg)%BASIN_solute(iTracer), & ! inout: updated fluxes at reach + ierr, message) ! output: error control + if(ierr/=0)then; message=trim(message)//trim(cmessage); return; endif + end do end if END SUBROUTINE hru_irf diff --git a/route/build/src/dataTypes.f90 b/route/build/src/dataTypes.f90 index 40cd31f8..23ffd10f 100644 --- a/route/build/src/dataTypes.f90 +++ b/route/build/src/dataTypes.f90 @@ -168,7 +168,7 @@ MODULE dataTypes real(dp) , allocatable :: basinRunoff(:) ! remapped river network catchment runoff [depth/time] (size: number of nHRU) real(dp) , allocatable :: basinEvapo(:) ! remapped river network catchment evaporation [depth/time] (size: number of nHRU) real(dp) , allocatable :: basinPrecip(:) ! remapped river network catchment precipitation [depth/time] (size: number of nHRU) - real(dp) , allocatable :: basinSolute(:) ! remapped river network catchment solute in water [mass/time] (size: number of nHRU) + real(dp) , allocatable :: basinSolute(:,:)! remapped river network catchment solute in water [mass/time] (size: nHRU, nTracer) end type runoff type, public, extends(inputData) :: wm ! water-management @@ -357,8 +357,8 @@ MODULE dataTypes real(dp) :: REACH_WM_FLUX_actual ! water management fluxes to and from each reach [m3/s] real(dp) :: WB ! reach water balance error [m3] real(dp) :: Qerror ! simulated discharge error compared to obs [m3/s] -- only for data assimilation - real(dp) :: reach_solute_mass(0:1) ! constituent mass in channel [mg] - real(dp) :: reach_solute_flux ! constituent mass flux from reach outlet [mg/s] + real(dp),allocatable :: reach_solute_mass(0:1,:) ! mass of constituent(s) in channel [mg] + real(dp),allocatable :: reach_solute_flux(:) ! mass flux of constituent(s) from reach outlet [mg/s] end type hydraulic ! fluxes and states in each reach @@ -368,8 +368,8 @@ MODULE dataTypes real(dp), allocatable :: QPASTUP_IRF(:,:) ! runoff volume in the past time steps for lake upstream [m3/s] real(dp), allocatable :: DEMANDPAST_IRF(:,:) ! demand volume for lake [m3/s] real(dp), allocatable :: solute_future(:) ! lateral solute mass in future time steps [mg/s] - real(dp) :: BASIN_solute ! instantaneous constituent mass from the local basin [mg/s] - real(dp) :: BASIN_solute_inst ! instantaneous constituent mass from the local basin [mg/s] + real(dp), allocatable :: BASIN_solute(:) ! mass of instantaneous constituent(s) from the local basin [mg/s] + real(dp), allocatable :: BASIN_solute_inst(:) ! mass of instantaneous constituent(s) from the local basin [mg/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] type(hydraulic), allocatable :: ROUTE(:) ! reach fluxes and states for each routing method diff --git a/route/build/src/globalData.f90 b/route/build/src/globalData.f90 index 8ccf1b6d..e68de20d 100644 --- a/route/build/src/globalData.f90 +++ b/route/build/src/globalData.f90 @@ -252,8 +252,11 @@ MODULE globalData real(dp), allocatable, public :: basinEvapo_main(:) ! HRU evaporation array (m/s) for mainstem real(dp), allocatable, public :: basinPrecip_trib(:) ! HRU precipitation array (m/s) for tributaries real(dp), allocatable, public :: basinPrecip_main(:) ! HRU precipitation array (m/s) for mainstem - real(dp), allocatable, public :: basinSolute_trib(:) ! HRU constituent array (mg/s) for tributaries - real(dp), allocatable, public :: basinSolute_main(:) ! HRU constituent array (mg/s) for mainstem + + ! solute (tracer) data + integer(i4b) :: nTracer ! number of tracers + real(dp), allocatable, public :: basinSolute_trib(:,:) ! HRU constituent array [nhru, nTracer] (mg/s) for tributaries + real(dp), allocatable, public :: basinSolute_main(:,:) ! HRU constituent array [nhru, nTracer] (mg/s) for mainstem ! seg water management fluxes and target volume type(wm), public :: wm_data ! SEG flux and target vol data structure for one time step for river network diff --git a/route/build/src/main_route.f90 b/route/build/src/main_route.f90 index 3996a50b..529a35f9 100644 --- a/route/build/src/main_route.f90 +++ b/route/build/src/main_route.f90 @@ -30,7 +30,7 @@ MODULE main_route_module SUBROUTINE main_route(basinRunoff_in, & ! input: basin (i.e.,HRU) runoff (m/s) basinEvapo_in, & ! input: basin (i.e.,HRU) evaporation (m/s) basinPrecip_in, & ! input: basin (i.e.,HRU) precipitation (m/s) - basinSolute_in, & ! input: basin (i.e.,HRU) constituent mass fluw (mg/s/m2) + basinSolute_in, & ! input: basin (i.e.,HRU) mass flux of constituent(s) (mg/s/m2) reachflux_in, & ! input: reach (i.e.,reach) flux (m3/s) reachvol_in, & ! input: reach (i.e.,reach) target volume for lakes (m3) ixRchProcessed, & ! input: indices of reach to be routed @@ -52,7 +52,8 @@ SUBROUTINE main_route(basinRunoff_in, & ! input: basin (i.e.,HRU) runoff (m/s) USE globalData, ONLY: TSEC ! beginning/ending of simulation time step [sec] USE globalData, ONLY: simDatetime ! current model datetime USE globalData, ONLY: rch_routes ! - USE public_var, ONLY: doesBasinRoute + USE globalData, ONLY: nTracer ! number of active tracers + USE public_var, ONLY: doesBasinRoute ! hillslope routing option USE public_var, ONLY: is_flux_wm ! logical whether or not fluxes should be passed USE public_var, ONLY: is_vol_wm ! logical whether or not target volume should be passed USE public_var, ONLY: qmodOption ! options for streamflow modification (DA) @@ -62,7 +63,7 @@ SUBROUTINE main_route(basinRunoff_in, & ! input: basin (i.e.,HRU) runoff (m/s) real(dp), allocatable, intent(in) :: basinRunoff_in(:) ! basin (i.e.,HRU) runoff (m/s) real(dp), allocatable, intent(in) :: basinEvapo_in(:) ! basin (i.e.,HRU) evaporation (m/s) real(dp), allocatable, intent(in) :: basinPrecip_in(:) ! basin (i.e.,HRU) precipitation (m/s) - real(dp), allocatable, intent(in) :: basinSolute_in(:) ! basin (i.e.,HRU) constituent (mg/s/m2) + real(dp), allocatable, intent(in) :: basinSolute_in(:,:) ! basin (i.e.,HRU) constituent(s) (mg/s/m2) real(dp), allocatable, intent(in) :: reachflux_in(:) ! reach (i.e.,reach) flux (m3/s) real(dp), allocatable, intent(in) :: reachvol_in(:) ! reach (i.e.,reach) target volume for lakes (m3) integer(i4b), allocatable, intent(in) :: ixRchProcessed(:) ! indices of reach to be routed @@ -78,11 +79,12 @@ SUBROUTINE main_route(basinRunoff_in, & ! input: basin (i.e.,HRU) runoff (m/s) character(len=strLen) :: cmessage ! error message of downwind routine real(dp) :: qobs ! observed discharge at a time step and site real(dp), allocatable :: reachRunoff_local(:) ! reach runoff (m3/s) - real(dp), allocatable :: reachSolute_local(:) ! reach constituent (mg/s) + real(dp), allocatable :: reachSolute_local(:,:) ! reach constituent(s) (mg/s) real(dp), allocatable :: reachEvapo_local(:) ! reach evaporation (m3/s) real(dp), allocatable :: reachPrecip_local(:) ! reach precipitation (m3/s) integer(i4b) :: nSeg ! number of reach to be processed integer(i4b) :: iSeg ! index of reach + integer(i4b) :: iTracer ! index of tracers integer(i4b) :: ix,jx ! loop index integer(i4b), allocatable :: reach_ix(:) integer(i4b), parameter :: no_mod=0 @@ -149,16 +151,17 @@ SUBROUTINE main_route(basinRunoff_in, & ! input: basin (i.e.,HRU) runoff (m/s) if(ierr/=0)then; message=trim(message)//trim(cmessage); return; endif if(tracer) then - allocate(reachSolute_local(nSeg), source=0._dp, stat=ierr) + allocate(reachSolute_local(nSeg,nTracer), source=0._dp, stat=ierr) if(ierr/=0)then; message=trim(message)//'problem allocating arrays for [reachSolute_local]'; return; endif - - call basin2reach_mass(basinSolute_in, & ! input: basin total constitunet (mg/s/m2) - NETOPO_in, & ! input: reach topology - RPARAM_in, & ! input: reach parameter - reachSolute_local, & ! output:total constituent going to reach (mg/s) - ierr, cmessage, & ! output: error control - ixRchProcessed) ! optional input: indices of reach to be routed - if(ierr/=0)then; message=trim(message)//trim(cmessage); return; endif + do iTracer=1,nTracer + call basin2reach_mass(basinSolute_in(:,iTracer), & ! input: basin total constitunet (mg/s/m2) + NETOPO_in, & ! input: reach topology + RPARAM_in, & ! input: reach parameter + reachSolute_local(:,iTracer), & ! output:total constituent going to reach (mg/s) + ierr, cmessage, & ! output: error control + ixRchProcessed) ! optional input: indices of reach to be routed + if(ierr/=0)then; message=trim(message)//trim(cmessage); return; endif + end do end if if (is_lake_sim) then @@ -196,7 +199,9 @@ SUBROUTINE main_route(basinRunoff_in, & ! input: basin (i.e.,HRU) runoff (m/s) RCHFLX_out(ixRchProcessed(iSeg))%BASIN_QI = reachRunoff_local(iSeg) if (tracer) then if (RCHFLX_out(ixRchProcessed(iSeg))%BASIN_QI>0) then ! this may cause mass inbalance between input and output. so may need to feed flag to land surface model - RCHFLX_out(ixRchProcessed(iSeg))%BASIN_solute_inst = reachSolute_local(iSeg) + do iTracer=1,nTracer + RCHFLX_out(ixRchProcessed(iSeg))%BASIN_solute_inst(iTracer) = reachSolute_local(iSeg,iTracer) + end do else RCHFLX_out(ixRchProcessed(iSeg))%BASIN_solute_inst = 0._dp endif @@ -218,7 +223,9 @@ SUBROUTINE main_route(basinRunoff_in, & ! input: basin (i.e.,HRU) runoff (m/s) if (tracer) then do iSeg = 1,nSeg if (RCHFLX_out(ixRchProcessed(iSeg))%BASIN_QR(1)>0) then - RCHFLX_out(ixRchProcessed(iSeg))%BASIN_solute = reachSolute_local(iSeg) ! total constituent going to reach (mg/s) + do iTracer=1,nTracer + RCHFLX_out(ixRchProcessed(iSeg))%BASIN_solute(iTracer) = reachSolute_local(iSeg,iTracer) ! total constituent going to reach (mg/s) + end do else RCHFLX_out(ixRchProcessed(iSeg))%BASIN_solute = 0._dp endif diff --git a/route/build/src/mpi_process.f90 b/route/build/src/mpi_process.f90 index 3465d500..ac728fc2 100644 --- a/route/build/src/mpi_process.f90 +++ b/route/build/src/mpi_process.f90 @@ -107,6 +107,7 @@ SUBROUTINE comm_ntopo_data(pid, & ! input: proc id USE globalData, ONLY: vol_wm_trib ! nRch target vol holder for mainstem USE globalData, ONLY: nRch_trib ! number of tributary reaches USE globalData, ONLY: nHRU_trib ! number of tributary HRUs + USE globalData, ONLY: nTracer ! number of tracer USE globalData, ONLY: hru_per_proc ! number of hrus assigned to each proc (size = num of procs+1) USE globalData, ONLY: rch_per_proc ! number of reaches assigned to each proc (size = num of procs+1) USE globalData, ONLY: ixHRU_order ! global HRU index in the order of proc assignment (size = total number of HRUs contributing to any reaches, nContribHRU) @@ -574,10 +575,10 @@ SUBROUTINE comm_ntopo_data(pid, & ! input: proc id if(ierr/=0)then; message=trim(message)//trim(cmessage); return; endif if (tracer) then - allocate(basinSolute_main(nHRU_mainstem), source=0.0_dp, stat=ierr, errmsg=cmessage) + allocate(basinSolute_main(nHRU_mainstem,nTracer), source=0.0_dp, stat=ierr, errmsg=cmessage) if(ierr/=0)then; message=trim(message)//trim(cmessage); return; endif - allocate(basinSolute_trib(nHRU_trib), source=0.0_dp, stat=ierr, errmsg=cmessage) + allocate(basinSolute_trib(nHRU_trib,nTracer), source=0.0_dp, stat=ierr, errmsg=cmessage) if(ierr/=0)then; message=trim(message)//trim(cmessage); return; endif end if @@ -629,7 +630,7 @@ SUBROUTINE comm_ntopo_data(pid, & ! input: proc id if(ierr/=0)then; message=trim(message)//trim(cmessage); return; endif if (tracer) then - allocate(basinSolute_trib(nHRU_trib), source=0.0_dp, stat=ierr, errmsg=cmessage) + allocate(basinSolute_trib(nHRU_trib,nTracer), source=0.0_dp, stat=ierr, errmsg=cmessage) if(ierr/=0)then; message=trim(message)//trim(cmessage); return; endif end if @@ -1394,6 +1395,7 @@ SUBROUTINE scatter_runoff(nNodes, & ! mpi variables: number nodes USE globalData, ONLY: nHRU ! number of all HRUs USE globalData, ONLY: nHRU_mainstem ! number of mainstem HRUs + USE globalData, ONLY: nTracer ! number of tracers USE globalData, ONLY: runoff_data ! runoff data structure USE globalData, ONLY: hru_per_proc ! number of hrus assigned to each proc (i.e., node) USE globalData, ONLY: basinRunoff_main ! HRU runoff holder for mainstem @@ -1413,10 +1415,11 @@ SUBROUTINE scatter_runoff(nNodes, & ! mpi variables: number nodes integer(i4b), intent(out) :: ierr ! error code character(len=strLen), intent(out) :: message ! error message ! local variables + integer(i4b) :: iTracer ! tracer loop index real(dp) :: basinRunoff_local(nHRU) ! temporal basin runoff (m/s) for whole domain real(dp) :: basinEvapo_local(nHRU) ! temporal basin evaporation (m/s) for whole domain real(dp) :: basinPrecip_local(nHRU) ! temporal basin precipitation (m/s) for whole domain - real(dp) :: basinSolute_local(nHRU) ! temporal basin constituent (g) for whole domain + real(dp) :: basinSolute_local(nHRU,nTracer) ! temporal basin constituent (mg) for whole domain character(len=strLen) :: cmessage ! error message from a subroutine ierr=0; message='scatter_runoff/' @@ -1433,7 +1436,7 @@ SUBROUTINE scatter_runoff(nNodes, & ! mpi variables: number nodes end if if (tracer) then - basinSolute_local (1:nHRU) = runoff_data%basinSolute(1:nHRU) + basinSolute_local (1:nHRU,1:nTracer) = runoff_data%basinSolute(1:nHRU,1:nTracer) end if if (nHRU_mainstem>0) then @@ -1459,11 +1462,10 @@ SUBROUTINE scatter_runoff(nNodes, & ! mpi variables: number nodes if (tracer) then if (.not. allocated(basinSolute_main)) then - allocate(basinSolute_main(nHRU_mainstem), stat=ierr) - if(ierr/=0)then; message=trim(message)//'problem allocating array for [basinEvapo_main]'; return; endif + allocate(basinSolute_main(nHRU_mainstem,nTracer), source=0.0_dp, stat=ierr, errmsg=cmessage) + if(ierr/=0)then; message=trim(message)//trim(cmessage); return; endif endif - basinSolute_main (1:nHRU_mainstem) = basinSolute_local (1:nHRU_mainstem) - + basinSolute_main (1:nHRU_mainstem, 1:nTracer) = basinSolute_local(1:nHRU_mainstem,1:nTracer) end if end if @@ -1491,11 +1493,13 @@ SUBROUTINE scatter_runoff(nNodes, & ! mpi variables: number nodes end if if (tracer) then - call shr_mpi_scatterV(basinSolute_local(nHRU_mainstem+1:nHRU), & - hru_per_proc(0:nNodes-1), & - basinSolute_trib, & - ierr, cmessage) - if(ierr/=0)then; message=trim(message)//trim(cmessage); return; endif + do iTracer=1,nTracer + call shr_mpi_scatterV(basinSolute_local(nHRU_mainstem+1:nHRU,iTracer), & + hru_per_proc(0:nNodes-1), & + basinSolute_trib(:,iTracer), & + ierr, cmessage) + if(ierr/=0)then; message=trim(message)//trim(cmessage); return; endif + end do end if END SUBROUTINE scatter_runoff diff --git a/route/build/src/process_remap.f90 b/route/build/src/process_remap.f90 index d06e026b..a519e48c 100644 --- a/route/build/src/process_remap.f90 +++ b/route/build/src/process_remap.f90 @@ -438,10 +438,10 @@ SUBROUTINE basin2reach_mass(basinSolute, & ! basin constituent mass flux p USE nr_utils, ONLY : arth implicit none ! Argument variables - real(dp) , intent(in) :: basinSolute(:) ! constituent mass flux per unit area (mg/s/m2) + real(dp) , intent(in) :: basinSolute(:) ! constituent mass flux per unit area (mg/s/m2) [nhru] type(RCHTOPO), allocatable, intent(in) :: NETOPO_in(:) ! River Network topology type(RCHPRP), allocatable, intent(in) :: RPARAM_in(:) ! River (non-)physical parameters - real(dp) , intent(out) :: reachSolute(:) ! constitunet mass flux into reach (g/s) + real(dp) , intent(out) :: reachSolute(:) ! constitunet mass flux into reach (mg/s) integer(i4b) , intent(out) :: ierr ! error code character(len=strLen) , intent(out) :: message ! error message integer(i4b), optional , intent(in) :: ixSubRch(:) ! subset of reach indices to be processed @@ -482,12 +482,11 @@ SUBROUTINE basin2reach_mass(basinSolute, & ! basin constituent mass flux p ! * case where HRUs drain into the segment if(nContrib > 0)then - - ! intialize the streamflow + ! intialize the solute mass reachSolute(jSeg) = 0._dp ! loop through the HRUs do iHRU=1,nContrib - ! error check - runoff depth cannot be negative (no missing value) + ! error check - solute mass cannot be negative (no missing value) if( basinSolute( hruContribIx(iHRU) ) < 0._dp )then write(iulog,*) 'Negative solute mass flux: HRU = ', hruContribId(iHRU), ' basin solute mass flux = ', basinSolute( hruContribIx(iHRU) ) write(message,'(a,i12, g12.2)') trim(message)//'Negative solute mass flux for HRU ', hruContribId(iHRU) @@ -495,9 +494,9 @@ SUBROUTINE basin2reach_mass(basinSolute, & ! basin constituent mass flux p endif ! compute the weighted average mass (mg/s) reachSolute(jSeg) = reachSolute(jSeg) + hruWeight(iHRU)*basinSolute( hruContribIx(iHRU) )*time_conv_solute*mass_conv_solute - end do ! (looping through contributing HRUs) + end do ! (looping through contributing HRUs) - ! convert basin average runoff volume (mg/s) + ! convert basin total mass (mg/s) reachSolute(jSeg) = reachSolute(jSeg)*basArea endif diff --git a/route/build/src/public_var.f90 b/route/build/src/public_var.f90 index 380115de..b83f6ae7 100644 --- a/route/build/src/public_var.f90 +++ b/route/build/src/public_var.f90 @@ -128,7 +128,7 @@ MODULE public_var character(len=strLen),public :: vname_qsim = 'runoff' ! variable name for runoff character(len=strLen),public :: vname_evapo = 'evap' ! variable name for actual evapoartion character(len=strLen),public :: vname_precip = 'precip' ! variable name for precipitation - character(len=strLen),public :: vname_solute = charMissing ! variable name for solute mass flux + character(len=strLen),public :: vname_solute(10) = charMissing ! variable name(s) for solute mass flux(es) (maximum 10 tracers) character(len=strLen),public :: vname_time = 'time' ! variable name for time character(len=strLen),public :: vname_hruid = 'hru' ! variable name for runoff hru id character(len=strLen),public :: dname_time = 'time' ! dimension name for time diff --git a/route/build/src/read_control.f90 b/route/build/src/read_control.f90 index f6e0d2b6..778a01b7 100644 --- a/route/build/src/read_control.f90 +++ b/route/build/src/read_control.f90 @@ -37,6 +37,7 @@ SUBROUTINE read_control(ctl_fname, err, message) USE globalData, ONLY: idxSUM,idxIRF,idxKWT, & idxKW,idxMC,idxDW USE globalData, ONLY: runMode ! mizuRoute run mode: standalone, cesm-coupling + USE globalData, ONLY: nTracer ! number of active tracers ! index of named variables in each structure USE var_lookup, ONLY: ixHRU USE var_lookup, ONLY: ixHRU2SEG @@ -48,6 +49,7 @@ SUBROUTINE read_control(ctl_fname, err, message) ! external subroutines USE ascii_utils, ONLY: file_open ! open file (performs a few checks as well) USE ascii_utils, ONLY: get_vlines ! get a list of character strings from non-comment lines + USE ascii_utils, ONLY: split_line ! split string into sub-strings with delimiter USE ascii_utils, ONLY: lower ! convert string to lower case USE nr_utils, ONLY: char2int ! convert integer number to a array containing individual digits @@ -61,6 +63,7 @@ SUBROUTINE read_control(ctl_fname, err, message) character(len=strLen) :: cName,cData ! name and data from cLines(iLine) character(len=strLen) :: cLength,cTime ! length and time units character(len=strLen) :: cMass ! mass units needed only when tracer is on + character(len=strLen),allocatable :: vname_tmp(:) ! multiple variable names with semi-colon delimited logical(lgt) :: isGeneric ! temporal logical scalar logical(lgt) :: onlyOneRouting ! temporal logical scalar integer(i4b) :: ipos ! index of character string @@ -158,7 +161,7 @@ SUBROUTINE read_control(ctl_fname, err, message) case(''); vname_qsim = trim(cData) ! name of runoff variable case(''); vname_evapo = trim(cData) ! name of actual evapoartion variable case(''); vname_precip = trim(cData) ! name of precipitation variable - case(''); vname_solute = trim(cData) ! name of solute mass flux variable + case(''); vname_tmp = split_line(trim(cData),delim=';,:') ! name of solute mass flux variables (up to 10 is ok), acceptable delimiter: : , : case(''); vname_time = trim(cData) ! name of time variable in the runoff file case(''); vname_hruid = trim(cData) ! name of the HRU id case(''); dname_time = trim(cData) ! name of time variable in the runoff file @@ -382,6 +385,9 @@ SUBROUTINE read_control(ctl_fname, err, message) end do ! looping through lines in the control file ! ---------- Perform minor processing and checking control variables ---------------------------------------- + ! extract tracer names into vname_solute + nTracer=size(vname_tmp) + vname_solute(1:nTracer)=vname_tmp(:) ! ---------- directory option --------------------------------------------------------------------- if (trim(restart_dir)==charMissing) then @@ -449,7 +455,7 @@ SUBROUTINE read_control(ctl_fname, err, message) cLength = units_qsim(1:ipos-1) cTime = units_qsim(ipos+1:len_trim(units_qsim)) - ! get the conversion factor for length + ! get the conversion factor for length (to m) select case(trim(cLength)) case('m'); length_conv = 1._dp case('mm'); length_conv = 1._dp/1000._dp @@ -458,7 +464,7 @@ SUBROUTINE read_control(ctl_fname, err, message) err=81; return end select - ! get the conversion factor for time + ! get the conversion factor for time (to second) select case(trim(cTime)) case('d','day'); time_conv = 1._dp/secprday case('h','hr','hour'); time_conv = 1._dp/secprhour @@ -480,17 +486,17 @@ SUBROUTINE read_control(ctl_fname, err, message) cMass = units_cc(1:ipos-1) cTime = units_cc(ipos+1:len_trim(units_cc)) - ! get the conversion factor for length + ! get the conversion factor for mass (to mg) select case(trim(cMass)) case('mg'); mass_conv_solute = 1._dp case('g'); mass_conv_solute = 1000._dp case('kg'); mass_conv_solute = 1000000._dp case default - message=trim(message)//'expect the mass units of solute to be "mg" [units='//trim(cMass)//']' + message=trim(message)//'expect the mass units of solute to be "mg","g",or "kg" [units='//trim(cMass)//']' err=81; return end select - ! get the conversion factor for time + ! get the conversion factor for time (to second) select case(trim(cTime)) case('d','day'); time_conv_solute = 1._dp/secprday case('h','hr','hour'); time_conv_solute = 1._dp/secprhour diff --git a/route/build/src/standalone/get_basin_runoff.f90 b/route/build/src/standalone/get_basin_runoff.f90 index 9c81c8ef..3d4c583b 100644 --- a/route/build/src/standalone/get_basin_runoff.f90 +++ b/route/build/src/standalone/get_basin_runoff.f90 @@ -44,6 +44,7 @@ SUBROUTINE get_hru_runoff(ierr, message) USE public_var, ONLY: verySmall ! very small values USE public_var, ONLY: realMissing ! real missing value USE globalData, ONLY: iTime ! simulation time step index + USE globalData, ONLY: nTracer ! number of active tracers USE globalData, ONLY: runoff_data ! data structure to hru runoff data USE globalData, ONLY: wm_data ! data strcuture for water management USE globalData, ONLY: remap_data ! data structure to remap data @@ -61,6 +62,7 @@ SUBROUTINE get_hru_runoff(ierr, message) integer(i4b), intent(out) :: ierr ! error code character(*), intent(out) :: message ! error message ! local variables + integer(i4b) :: iTracer ! tracer loop index type(map_time) :: tmap_sim_ro ! time-steps mapping data type(map_time) :: tmap_sim_wm ! time-steps mapping data logical(lgt) :: remove_negatives ! flag to replace the negative values to zeros @@ -107,29 +109,31 @@ SUBROUTINE get_hru_runoff(ierr, message) if(ierr/=0) then; message=trim(message)//trim(cmessage); return; endif end if - ! Optional: reading constituent + ! Optional: reading constituents if (tracer) then - call read_forcing_data(inFileInfo_ro, & ! input: meta for runoff input files - vname_solute, & ! input: solute mass flux varname - tmap_sim_ro, & ! input: ro-sim time mapping at current simulation step - runoff_data, & ! inout: constituent conentration data structure - ierr, cmessage) ! output: error control - if(ierr/=0)then; message=trim(message)//trim(cmessage); return; endif - - ! Get river network HRU constituent into runoff_data data structure - if (is_remap) then ! remap LSM simulated flux to the HRUs in the river network - call remap_runoff(runoff_data, remap_data, runoff_data%basinSolute, ierr, cmessage) - if(ierr/=0)then; message=trim(message)//trim(cmessage); return; endif - else ! runoff is already remapped to river network HRUs - remove_negatives = .true. - call sort_flux (runoff_data%hru_id, & - runoff_data%hru_ix, & - runoff_data%sim, & - remove_negatives, & - runoff_data%basinSolute, & - ierr, cmessage) + do iTracer = 1,nTracer + call read_forcing_data(inFileInfo_ro, & ! input: meta for runoff input files + vname_solute(iTracer), & ! input: solute mass flux varname + tmap_sim_ro, & ! input: ro-sim time mapping at current simulation step + runoff_data, & ! inout: constituent conentration data structure + ierr, cmessage) ! output: error control if(ierr/=0)then; message=trim(message)//trim(cmessage); return; endif - end if + + ! Get river network HRU constituent into runoff_data data structure + if (is_remap) then ! remap LSM simulated flux to the HRUs in the river network + call remap_runoff(runoff_data, remap_data, runoff_data%basinSolute(:,iTracer), ierr, cmessage) + if(ierr/=0)then; message=trim(message)//trim(cmessage); return; endif + else ! runoff is already remapped to river network HRUs + remove_negatives = .true. + call sort_flux (runoff_data%hru_id, & + runoff_data%hru_ix, & + runoff_data%sim, & + remove_negatives, & + runoff_data%basinSolute(:,iTracer), & + ierr, cmessage) + if(ierr/=0)then; message=trim(message)//trim(cmessage); return; endif + end if + end do ! end of iTracer loop end if ! Optional: lake module on -> read actual evaporation and preciptation diff --git a/route/build/src/standalone/model_setup.f90 b/route/build/src/standalone/model_setup.f90 index ec1be8ac..9a417e32 100644 --- a/route/build/src/standalone/model_setup.f90 +++ b/route/build/src/standalone/model_setup.f90 @@ -701,6 +701,7 @@ SUBROUTINE init_forc_data(ierr, message) ! output: error control USE public_var, ONLY: is_vol_wm ! logical whether or not target volume should be read USE globalData, ONLY: nHRU ! number of HRUs over the entire domain USE globalData, ONLY: nRch ! number of reaches over the entire domain + USE globalData, ONLY: nTracer ! number of tracers USE globalData, ONLY: basinID ! basin ID USE globalData, ONLY: reachID ! reach ID USE globalData, ONLY: inFileInfo_ro ! metadata of the ro/evap/p input files @@ -769,7 +770,7 @@ SUBROUTINE init_forc_data(ierr, message) ! output: error control runoff_data%basinRunoff(:) = realMissing if (tracer) then - allocate(runoff_data%basinSolute(nHRU), source=realMissing, stat=ierr, errmsg=cmessage) + allocate(runoff_data%basinSolute(nHRU,nTracer), source=realMissing, stat=ierr, errmsg=cmessage) if(ierr/=0)then; message=trim(message)//trim(cmessage); return; endif end if