From 0d4e11d099571b32ca9bfc70e743d8488a132a0f Mon Sep 17 00:00:00 2001 From: Sabin Taranu Date: Mon, 12 Dec 2022 11:38:49 +0100 Subject: [PATCH 1/8] Update to src/cpl/mct/rof_comp_mct.F90 --- src/cpl/mct/rof_comp_mct.F90 | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/src/cpl/mct/rof_comp_mct.F90 b/src/cpl/mct/rof_comp_mct.F90 index 56b4c90..f2d8c42 100644 --- a/src/cpl/mct/rof_comp_mct.F90 +++ b/src/cpl/mct/rof_comp_mct.F90 @@ -33,6 +33,11 @@ module rof_comp_mct use mosart_cpl_indices , only : index_x2r_Flrl_rofsur, index_x2r_Flrl_rofi use mosart_cpl_indices , only : index_x2r_Flrl_rofgwl, index_x2r_Flrl_rofsub use mosart_cpl_indices , only : index_x2r_Flrl_irrig + use mosart_cpl_indices , only : index_x2r_Flrl_dom_withd, index_x2r_Flrl_dom_rf + use mosart_cpl_indices , only : index_x2r_Flrl_liv_withd, index_x2r_Flrl_liv_rf + use mosart_cpl_indices , only : index_x2r_Flrl_elec_withd, index_x2r_Flrl_elec_rf + use mosart_cpl_indices , only : index_x2r_Flrl_mfc_withd, index_x2r_Flrl_mfc_rf + use mosart_cpl_indices , only : index_x2r_Flrl_min_withd, index_x2r_Flrl_min_rf use mosart_cpl_indices , only : index_r2x_Forr_rofl, index_r2x_Forr_rofi, index_r2x_Flrr_flood use mosart_cpl_indices , only : index_r2x_Flrr_volr, index_r2x_Flrr_volrmch From 57f310e40f964b629309a15bfd7e5575800153f0 Mon Sep 17 00:00:00 2001 From: Sabin Taranu Date: Mon, 12 Dec 2022 11:39:41 +0100 Subject: [PATCH 2/8] Update to src/cpl/mct/mosart_cpl_indices.F90 --- src/cpl/mct/mosart_cpl_indices.F90 | 21 +++++++++++++++++++++ 1 file changed, 21 insertions(+) diff --git a/src/cpl/mct/mosart_cpl_indices.F90 b/src/cpl/mct/mosart_cpl_indices.F90 index 403db10..f41bde7 100644 --- a/src/cpl/mct/mosart_cpl_indices.F90 +++ b/src/cpl/mct/mosart_cpl_indices.F90 @@ -20,6 +20,17 @@ module mosart_cpl_indices integer, public :: index_x2r_Flrl_rofdto = 0 ! lnd->rof liquid direct to ocean runoff integer, public :: index_x2r_Flrl_rofi = 0 ! lnd->rof ice runoff forcing from land integer, public :: index_x2r_Flrl_irrig = 0 ! lnd->rof fraction of volr to be removed for irrigation + integer, public :: index_x2r_Flrl_dom_withd = 0 ! lnd->rof fraction of volr to be removed for domestic usage + integer, public :: index_x2r_Flrl_dom_rf = 0 ! lnd->rof fraction of volr to be returned after domestic usage + integer, public :: index_x2r_Flrl_liv_withd = 0 ! lnd->rof fraction of volr to be removed for livestock usage + integer, public :: index_x2r_Flrl_liv_rf = 0 ! lnd->rof fraction of volr to be returned after livestock usage + integer, public :: index_x2r_Flrl_elec_withd = 0 ! lnd->rof fraction of volr to be removed for thermoelectric usage + integer, public :: index_x2r_Flrl_elec_rf = 0 ! lnd->rof fraction of volr to be returned after thermoelectric usage + integer, public :: index_x2r_Flrl_mfc_withd = 0 ! lnd->rof fraction of volr to be removed for manufacturing usage + integer, public :: index_x2r_Flrl_mfc_rf = 0 ! lnd->rof fraction of volr to be returned after manufacturing usage + integer, public :: index_x2r_Flrl_min_withd = 0 ! lnd->rof fraction of volr to be removed for mining usage + integer, public :: index_x2r_Flrl_min_rf = 0 ! lnd->rof fraction of volr to be returned after mining usage + integer, public :: nflds_x2r = 0 ! roff to driver (part of land for now) (optional if ROF is off) @@ -65,6 +76,16 @@ subroutine mosart_cpl_indices_set(flds_x2r, flds_r2x ) index_x2r_Flrl_rofdto = mct_avect_indexra(avtmp,'Flrl_rofdto') index_x2r_Flrl_rofi = mct_avect_indexra(avtmp,'Flrl_rofi') index_x2r_Flrl_irrig = mct_avect_indexra(avtmp,'Flrl_irrig') + index_x2r_Flrl_dom_withd = mct_avect_indexra(avtmp,'Flrl_dom_withd') + index_x2r_Flrl_dom_rf = mct_avect_indexra(avtmp,'Flrl_dom_rf') + index_x2r_Flrl_liv_withd = mct_avect_indexra(avtmp,'Flrl_liv_withd') + index_x2r_Flrl_liv_rf = mct_avect_indexra(avtmp,'Flrl_liv_rf') + index_x2r_Flrl_elec_withd = mct_avect_indexra(avtmp,'Flrl_elec_withd') + index_x2r_Flrl_elec_rf = mct_avect_indexra(avtmp,'Flrl_elec_rf') + index_x2r_Flrl_mfc_withd = mct_avect_indexra(avtmp,'Flrl_mfc_withd') + index_x2r_Flrl_mfc_rf = mct_avect_indexra(avtmp,'Flrl_mfc_rf') + index_x2r_Flrl_min_withd = mct_avect_indexra(avtmp,'Flrl_min_withd') + index_x2r_Flrl_min_rf = mct_avect_indexra(avtmp,'Flrl_min_rf') nflds_x2r = mct_avect_nRattr(avtmp) From ce4b68170554e1d56fd035099183563b6691d8f0 Mon Sep 17 00:00:00 2001 From: Sabin Taranu Date: Mon, 12 Dec 2022 11:40:33 +0100 Subject: [PATCH 3/8] Update to src/cpl/nuopc/rof_import_export.F90 --- src/cpl/nuopc/rof_import_export.F90 | 52 ++++++++++++++++++++++++++++- 1 file changed, 51 insertions(+), 1 deletion(-) diff --git a/src/cpl/nuopc/rof_import_export.F90 b/src/cpl/nuopc/rof_import_export.F90 index 30fe4fb..a8956e3 100644 --- a/src/cpl/nuopc/rof_import_export.F90 +++ b/src/cpl/nuopc/rof_import_export.F90 @@ -114,7 +114,17 @@ subroutine advertise_fields(gcomp, flds_scalar_name, do_rtmflood, rc) call fldlist_add(fldsToRof_num, fldsToRof, 'Flrl_rofsub') call fldlist_add(fldsToRof_num, fldsToRof, 'Flrl_rofi') call fldlist_add(fldsToRof_num, fldsToRof, 'Flrl_irrig') - + call fldlist_add(fldsToRof_num, fldsToRof, 'Flrl_dom_withd') + call fldlist_add(fldsToRof_num, fldsToRof, 'Flrl_dom_rf') + call fldlist_add(fldsToRof_num, fldsToRof, 'Flrl_liv_withd') + call fldlist_add(fldsToRof_num, fldsToRof, 'Flrl_liv_rf') + call fldlist_add(fldsToRof_num, fldsToRof, 'Flrl_elec_withd') + call fldlist_add(fldsToRof_num, fldsToRof, 'Flrl_elec_rf') + call fldlist_add(fldsToRof_num, fldsToRof, 'Flrl_mfc_withd') + call fldlist_add(fldsToRof_num, fldsToRof, 'Flrl_mfc_rf') + call fldlist_add(fldsToRof_num, fldsToRof, 'Flrl_min_withd') + call fldlist_add(fldsToRof_num, fldsToRof, 'Flrl_min_rf') + do n = 1,fldsToRof_num call NUOPC_Advertise(importState, standardName=fldsToRof(n)%stdname, & TransferOfferGeomObject='will provide', rc=rc) @@ -293,6 +303,46 @@ subroutine import_fields( gcomp, rc ) do_area_correction=.true., rc=rc) if (ChkErr(rc,__LINE__,u_FILE_u)) return + call state_getimport(importState, 'Flrl_dom_withd', begr, endr, rtmCTL%area, output=rtmCTL%qdom_withd(:), & + do_area_correction=.true., rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + + call state_getimport(importState, 'Flrl_dom_rf', begr, endr, rtmCTL%area, output=rtmCTL%qdom_rf(:), & + do_area_correction=.true., rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + + call state_getimport(importState, 'Flrl_liv_withd', begr, endr, rtmCTL%area, output=rtmCTL%qliv_withd(:), & + do_area_correction=.true., rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + + call state_getimport(importState, 'Flrl_liv_rf', begr, endr, rtmCTL%area, output=rtmCTL%qliv_rf(:), & + do_area_correction=.true., rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + + call state_getimport(importState, 'Flrl_elec_withd', begr, endr, rtmCTL%area, output=rtmCTL%qelec_withd(:), & + do_area_correction=.true., rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + + call state_getimport(importState, 'Flrl_elec_rf', begr, endr, rtmCTL%area, output=rtmCTL%qelec_rf(:), & + do_area_correction=.true., rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + + call state_getimport(importState, 'Flrl_mfc_withd', begr, endr, rtmCTL%area, output=rtmCTL%qmfc_withd(:), & + do_area_correction=.true., rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + + call state_getimport(importState, 'Flrl_mfc_rf', begr, endr, rtmCTL%area, output=rtmCTL%qmfc_rf(:), & + do_area_correction=.true., rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + + call state_getimport(importState, 'Flrl_min_withd', begr, endr, rtmCTL%area, output=rtmCTL%qmin_withd(:), & + do_area_correction=.true., rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + + call state_getimport(importState, 'Flrl_min_rf', begr, endr, rtmCTL%area, output=rtmCTL%qmin_rf(:), & + do_area_correction=.true., rc=rc) + if (ChkErr(rc,__LINE__,u_FILE_u)) return + rtmCTL%qsub(begr:endr, nfrz) = 0.0_r8 rtmCTL%qgwl(begr:endr, nfrz) = 0.0_r8 From 0562bf4aeb0adfd73fd88d7a4f45acf39d23bc1d Mon Sep 17 00:00:00 2001 From: Sabin Taranu Date: Mon, 12 Dec 2022 11:41:32 +0100 Subject: [PATCH 4/8] Update to src/cpl/mct/mosart_import_export.F90 --- src/cpl/mct/mosart_import_export.F90 | 26 ++++++++++++++++++++++++++ 1 file changed, 26 insertions(+) diff --git a/src/cpl/mct/mosart_import_export.F90 b/src/cpl/mct/mosart_import_export.F90 index 1ea0c88..22e733b 100644 --- a/src/cpl/mct/mosart_import_export.F90 +++ b/src/cpl/mct/mosart_import_export.F90 @@ -5,6 +5,11 @@ module mosart_import_export use mosart_cpl_indices , only : index_x2r_Flrl_rofsur, index_x2r_Flrl_rofi use mosart_cpl_indices , only : index_x2r_Flrl_rofgwl, index_x2r_Flrl_rofsub use mosart_cpl_indices , only : index_x2r_Flrl_irrig + use mosart_cpl_indices , only : index_x2r_Flrl_dom_withd, index_x2r_Flrl_dom_rf + use mosart_cpl_indices , only : index_x2r_Flrl_liv_withd, index_x2r_Flrl_liv_rf + use mosart_cpl_indices , only : index_x2r_Flrl_elec_withd, index_x2r_Flrl_elec_rf + use mosart_cpl_indices , only : index_x2r_Flrl_mfc_withd, index_x2r_Flrl_mfc_rf + use mosart_cpl_indices , only : index_x2r_Flrl_min_withd, index_x2r_Flrl_min_rf use mosart_cpl_indices , only : index_r2x_Forr_rofl, index_r2x_Forr_rofi use mosart_cpl_indices , only : index_r2x_Flrr_flood use mosart_cpl_indices , only : index_r2x_Flrr_volr, index_r2x_Flrr_volrmch @@ -69,6 +74,17 @@ subroutine mosart_import( x2r ) rtmCTL%qsur(n,nfrz) = x2r(index_x2r_Flrl_rofi ,n2) * (rtmCTL%area(n)*0.001_r8) rtmCTL%qirrig(n) = x2r(index_x2r_Flrl_irrig,n2) * (rtmCTL%area(n)*0.001_r8) + rtmCTL%qdom_withd(n) = x2r(index_x2r_Flrl_dom_withd,n2) * (rtmCTL%area(n)*0.001_r8) + rtmCTL%qdom_rf(n) = x2r(index_x2r_Flrl_dom_rf,n2) * (rtmCTL%area(n)*0.001_r8) + rtmCTL%qliv_withd(n) = x2r(index_x2r_Flrl_liv_withd,n2) * (rtmCTL%area(n)*0.001_r8) + rtmCTL%qliv_rf(n) = x2r(index_x2r_Flrl_liv_rf,n2) * (rtmCTL%area(n)*0.001_r8) + rtmCTL%qelec_withd(n) = x2r(index_x2r_Flrl_elec_withd,n2) * (rtmCTL%area(n)*0.001_r8) + rtmCTL%qelec_rf(n) = x2r(index_x2r_Flrl_elec_rf,n2) * (rtmCTL%area(n)*0.001_r8) + rtmCTL%qmfc_withd(n) = x2r(index_x2r_Flrl_mfc_withd,n2) * (rtmCTL%area(n)*0.001_r8) + rtmCTL%qmfc_rf(n) = x2r(index_x2r_Flrl_mfc_rf,n2) * (rtmCTL%area(n)*0.001_r8) + rtmCTL%qmin_withd(n) = x2r(index_x2r_Flrl_min_withd,n2) * (rtmCTL%area(n)*0.001_r8) + rtmCTL%qmin_rf(n) = x2r(index_x2r_Flrl_min_rf,n2) * (rtmCTL%area(n)*0.001_r8) + rtmCTL%qsub(n,nfrz) = 0.0_r8 rtmCTL%qgwl(n,nfrz) = 0.0_r8 @@ -81,6 +97,16 @@ subroutine mosart_import( x2r ) write(iulog,F01)'import: nstep, n, Flrl_rofgwl = ',get_nstep(),n,rtmCTL%qgwl(n,nliq) write(iulog,F01)'import: nstep, n, Flrl_rofi = ',get_nstep(),n,rtmCTL%qsur(n,nfrz) write(iulog,F01)'import: nstep, n, Flrl_irrig = ',get_nstep(),n,rtmCTL%qirrig(n) + write(iulog,F01)'import: nstep, n, Flrl_dom_withd = ',get_nstep(),n,rtmCTL%qdom_withd(n) + write(iulog,F01)'import: nstep, n, Flrl_dom_rf = ',get_nstep(),n,rtmCTL%qdom_rf(n) + write(iulog,F01)'import: nstep, n, Flrl_liv_withd = ',get_nstep(),n,rtmCTL%qliv_withd(n) + write(iulog,F01)'import: nstep, n, Flrl_liv_rf = ',get_nstep(),n,rtmCTL%qliv_rf(n) + write(iulog,F01)'import: nstep, n, Flrl_elec_withd = ',get_nstep(),n,rtmCTL%qelec_withd(n) + write(iulog,F01)'import: nstep, n, Flrl_elec_rf = ',get_nstep(),n,rtmCTL%qelec_rf(n) + write(iulog,F01)'import: nstep, n, Flrl_mfc_withd = ',get_nstep(),n,rtmCTL%qmfc_withd(n) + write(iulog,F01)'import: nstep, n, Flrl_mfc_rf = ',get_nstep(),n,rtmCTL%qmfc_rf(n) + write(iulog,F01)'import: nstep, n, Flrl_min_withd = ',get_nstep(),n,rtmCTL%qmin_withd(n) + write(iulog,F01)'import: nstep, n, Flrl_min_rf = ',get_nstep(),n,rtmCTL%qmin_rf(n) end do end if From c4086080529578a039601c505e2e9b1d52e58b53 Mon Sep 17 00:00:00 2001 From: Sabin Taranu Date: Mon, 12 Dec 2022 11:43:06 +0100 Subject: [PATCH 5/8] Update to src/riverroute/RtmHistFlds.F90 --- src/riverroute/RtmHistFlds.F90 | 80 ++++++++++++++++++++++++++++++++++ 1 file changed, 80 insertions(+) diff --git a/src/riverroute/RtmHistFlds.F90 b/src/riverroute/RtmHistFlds.F90 index 8550930..54bffae 100644 --- a/src/riverroute/RtmHistFlds.F90 +++ b/src/riverroute/RtmHistFlds.F90 @@ -126,10 +126,90 @@ subroutine RtmHistFldsInit() avgflag='A', long_name='Amount of water used for irrigation (total flux received from coupler)', & ptr_rof=rtmCTL%qirrig, default='inactive') + call RtmHistAddfld (fname='QDOM_WITHD_FROM_COUPLER', units='m3/s', & + avgflag='A', long_name='Amount of water withdrew for domestic usage', & + ptr_rof=rtmCTL%qdom_withd, default='inactive') + + call RtmHistAddfld (fname='QDOM_RF_FROM_COUPLER', units='m3/s', & + avgflag='A', long_name='Amount of water returned from domestic usage', & + ptr_rof=rtmCTL%qdom_rf, default='inactive') + + call RtmHistAddfld (fname='QLIV_WITHD_FROM_COUPLER', units='m3/s', & + avgflag='A', long_name='Amount of water withdrew for livestock usage', & + ptr_rof=rtmCTL%qliv_withd, default='inactive') + + call RtmHistAddfld (fname='QLIV_RF_FROM_COUPLER', units='m3/s', & + avgflag='A', long_name='Amount of water returned from livestock usage', & + ptr_rof=rtmCTL%qliv_rf, default='inactive') + + call RtmHistAddfld (fname='QELEC_WITHD_FROM_COUPLER', units='m3/s', & + avgflag='A', long_name='Amount of water withdrew for thermoelectric usage', & + ptr_rof=rtmCTL%qelec_withd, default='inactive') + + call RtmHistAddfld (fname='QELEC_RF_FROM_COUPLER', units='m3/s', & + avgflag='A', long_name='Amount of water returned from thermoelectric usage', & + ptr_rof=rtmCTL%qelec_rf, default='inactive') + + call RtmHistAddfld (fname='QMFC_WITHD_FROM_COUPLER', units='m3/s', & + avgflag='A', long_name='Amount of water withdrew for manufacturing usage', & + ptr_rof=rtmCTL%qmfc_withd, default='inactive') + + call RtmHistAddfld (fname='QMFC_RF_FROM_COUPLER', units='m3/s', & + avgflag='A', long_name='Amount of water returned from manufacturing usage', & + ptr_rof=rtmCTL%qmfc_rf, default='inactive') + + call RtmHistAddfld (fname='QMIN_WITHD_FROM_COUPLER', units='m3/s', & + avgflag='A', long_name='Amount of water withdrew for mining usage', & + ptr_rof=rtmCTL%qmin_withd, default='inactive') + + call RtmHistAddfld (fname='QMIN_RF_FROM_COUPLER', units='m3/s', & + avgflag='A', long_name='Amount of water returned from mining usage', & + ptr_rof=rtmCTL%qmin_rf, default='inactive') + call RtmHistAddfld (fname='QIRRIG_ACTUAL', units='m3/s', & avgflag='A', long_name='Actual irrigation (if limited by river storage)', & ptr_rof=rtmCTL%qirrig_actual, default='inactive') + call RtmHistAddfld (fname='QDOM_ACTUAL_WITHD', units='m3/s', & + avgflag='A', long_name='Actual domestic withdrawal (if limited by river storage)', & + ptr_rof=rtmCTL%qdom_actual_withd, default='inactive') + + call RtmHistAddfld (fname='QELEC_ACTUAL_WITHD', units='m3/s', & + avgflag='A', long_name='Actual thermoelectric withdrawal (if limited by river storage)', & + ptr_rof=rtmCTL%qelec_actual_withd, default='inactive') + + call RtmHistAddfld (fname='QLIV_ACTUAL_WITHD', units='m3/s', & + avgflag='A', long_name='Actual livestock withdrawal (if limited by river storage)', & + ptr_rof=rtmCTL%qliv_actual_withd, default='inactive') + + call RtmHistAddfld (fname='QMFC_ACTUAL_WITHD', units='m3/s', & + avgflag='A', long_name='Actual manufacturing withdrawal (if limited by river storage)', & + ptr_rof=rtmCTL%qmfc_actual_withd, default='inactive') + + call RtmHistAddfld (fname='QMIN_ACTUAL_WITHD', units='m3/s', & + avgflag='A', long_name='Actual mining withdrawal (if limited by river storage)', & + ptr_rof=rtmCTL%qmin_actual_withd, default='inactive') + + call RtmHistAddfld (fname='QDOM_ACTUAL_RF', units='m3/s', & + avgflag='A', long_name='Actual domestic return flow (if limited by river storage)', & + ptr_rof=rtmCTL%qdom_actual_rf, default='inactive') + + call RtmHistAddfld (fname='QELEC_ACTUAL_RF', units='m3/s', & + avgflag='A', long_name='Actual thermoelectric return flow (if limited by river storage)', & + ptr_rof=rtmCTL%qelec_actual_rf, default='inactive') + + call RtmHistAddfld (fname='QLIV_ACTUAL_RF', units='m3/s', & + avgflag='A', long_name='Actual livestock return flow (if limited by river storage)', & + ptr_rof=rtmCTL%qliv_actual_rf, default='inactive') + + call RtmHistAddfld (fname='QMFC_ACTUAL_RF', units='m3/s', & + avgflag='A', long_name='Actual manufacturing return flow (if limited by river storage)', & + ptr_rof=rtmCTL%qmfc_actual_rf, default='inactive') + + call RtmHistAddfld (fname='QMIN_ACTUAL_RF', units='m3/s', & + avgflag='A', long_name='Actual mining return flow (if limited by river storage)', & + ptr_rof=rtmCTL%qmin_actual_rf, default='inactive') + ! Print masterlist of history fields call RtmHistPrintflds() From 85cd86d05d162e3e541000e32977b220fdbf9dee Mon Sep 17 00:00:00 2001 From: Sabin Taranu Date: Mon, 12 Dec 2022 12:04:02 +0100 Subject: [PATCH 6/8] Update to src/riverroute/RunoffMod.F90 --- src/riverroute/RunoffMod.F90 | 157 ++++++++++++++++++++++++----------- 1 file changed, 109 insertions(+), 48 deletions(-) diff --git a/src/riverroute/RunoffMod.F90 b/src/riverroute/RunoffMod.F90 index 995be6c..0b4274f 100644 --- a/src/riverroute/RunoffMod.F90 +++ b/src/riverroute/RunoffMod.F90 @@ -80,6 +80,27 @@ module RunoffMod real(r8), pointer :: qirrig(:) ! coupler irrigation [m3/s] real(r8), pointer :: qirrig_actual(:) ! minimum of irrigation and available main channel storage + real(r8), pointer :: qdom_withd(:) ! domestic withdrawal from the coupler + real(r8), pointer :: qdom_rf(:) ! domestic return flow from the coupler + real(r8), pointer :: qliv_withd(:) ! livestock withdrawal from the coupler + real(r8), pointer :: qliv_rf(:) ! livestock return flow from the coupler + real(r8), pointer :: qelec_withd(:) ! thermoelectric withdrawal from the coupler + real(r8), pointer :: qelec_rf(:) ! thermoelectric return flow from the coupler + real(r8), pointer :: qmfc_withd(:) ! manufactured withdrawal from the coupler + real(r8), pointer :: qmfc_rf(:) ! manufactured return flow from the coupler + real(r8), pointer :: qmin_withd(:) ! mining withdrawal from the coupler + real(r8), pointer :: qmin_rf(:) ! mining return flow from the coupler + real(r8), pointer :: qdom_actual_withd(:) ! actual domestic withdrawal based on available main channel storage + real(r8), pointer :: qliv_actual_withd(:) ! actual livestock withdrawal based on available main channel storage + real(r8), pointer :: qelec_actual_withd(:) ! actual thermoelectric withdrawal based on available main channel storage + real(r8), pointer :: qmfc_actual_withd(:) ! actual manufacturing withdrawal based on available main channel storage + real(r8), pointer :: qmin_actual_withd(:) ! actual mining withdrawal based on available main channel storage + real(r8), pointer :: qdom_actual_rf(:) ! actual domestic return flow based on available main channel storage + real(r8), pointer :: qliv_actual_rf(:) ! actual livestock return flow based on available main channel storage + real(r8), pointer :: qelec_actual_rf(:) ! actual thermoelectric return flow based on available main channel storage + real(r8), pointer :: qmfc_actual_rf(:) ! actual manufacturing return flow based on available main channel storage + real(r8), pointer :: qmin_actual_rf(:) ! actual mining return flow based on available main channel storage + ! - history (currently needed) real(r8), pointer :: runofflnd_nt1(:) @@ -291,54 +312,74 @@ subroutine RunoffInit(begr, endr, numr) integer :: ier - allocate(rtmCTL%runoff(begr:endr,nt_rtm), & - rtmCTL%dvolrdt(begr:endr,nt_rtm), & - rtmCTL%runofflnd(begr:endr,nt_rtm), & - rtmCTL%dvolrdtlnd(begr:endr,nt_rtm), & - rtmCTL%runoffocn(begr:endr,nt_rtm), & - rtmCTL%dvolrdtocn(begr:endr,nt_rtm), & - rtmCTL%runofftot(begr:endr,nt_rtm), & - rtmCTL%area(begr:endr), & - rtmCTL%volr(begr:endr,nt_rtm), & - rtmCTL%lonc(begr:endr), & - rtmCTL%latc(begr:endr), & - rtmCTL%dsig(begr:endr), & - rtmCTL%outletg(begr:endr), & - rtmCTL%runofflnd_nt1(begr:endr), & - rtmCTL%runofflnd_nt2(begr:endr), & - rtmCTL%runoffocn_nt1(begr:endr), & - rtmCTL%runoffocn_nt2(begr:endr), & - rtmCTL%runofftot_nt1(begr:endr), & - rtmCTL%runofftot_nt2(begr:endr), & - rtmCTL%runoffdir_nt1(begr:endr), & - rtmCTL%runoffdir_nt2(begr:endr), & - rtmCTL%volr_nt1(begr:endr), & - rtmCTL%volr_nt2(begr:endr), & - rtmCTL%volr_mch(begr:endr), & - rtmCTL%dvolrdtlnd_nt1(begr:endr), & - rtmCTL%dvolrdtlnd_nt2(begr:endr), & - rtmCTL%dvolrdtocn_nt1(begr:endr), & - rtmCTL%dvolrdtocn_nt2(begr:endr), & - rtmCTL%qsur_nt1(begr:endr), & - rtmCTL%qsur_nt2(begr:endr), & - rtmCTL%qsub_nt1(begr:endr), & - rtmCTL%qsub_nt2(begr:endr), & - rtmCTL%qgwl_nt1(begr:endr), & - rtmCTL%qgwl_nt2(begr:endr), & - rtmCTL%mask(begr:endr), & - rtmCTL%gindex(begr:endr), & - rtmCTL%fthresh(begr:endr), & - rtmCTL%flood(begr:endr), & - rtmCTL%direct(begr:endr,nt_rtm), & - rtmCTL%wh(begr:endr,nt_rtm), & - rtmCTL%wt(begr:endr,nt_rtm), & - rtmCTL%wr(begr:endr,nt_rtm), & - rtmCTL%erout(begr:endr,nt_rtm), & - rtmCTL%qsur(begr:endr,nt_rtm), & - rtmCTL%qsub(begr:endr,nt_rtm), & - rtmCTL%qgwl(begr:endr,nt_rtm), & - rtmCTL%qirrig(begr:endr), & - rtmCTL%qirrig_actual(begr:endr), & + allocate(rtmCTL%runoff(begr:endr,nt_rtm), & + rtmCTL%dvolrdt(begr:endr,nt_rtm), & + rtmCTL%runofflnd(begr:endr,nt_rtm), & + rtmCTL%dvolrdtlnd(begr:endr,nt_rtm), & + rtmCTL%runoffocn(begr:endr,nt_rtm), & + rtmCTL%dvolrdtocn(begr:endr,nt_rtm), & + rtmCTL%runofftot(begr:endr,nt_rtm), & + rtmCTL%area(begr:endr), & + rtmCTL%volr(begr:endr,nt_rtm), & + rtmCTL%lonc(begr:endr), & + rtmCTL%latc(begr:endr), & + rtmCTL%dsig(begr:endr), & + rtmCTL%outletg(begr:endr), & + rtmCTL%runofflnd_nt1(begr:endr), & + rtmCTL%runofflnd_nt2(begr:endr), & + rtmCTL%runoffocn_nt1(begr:endr), & + rtmCTL%runoffocn_nt2(begr:endr), & + rtmCTL%runofftot_nt1(begr:endr), & + rtmCTL%runofftot_nt2(begr:endr), & + rtmCTL%runoffdir_nt1(begr:endr), & + rtmCTL%runoffdir_nt2(begr:endr), & + rtmCTL%volr_nt1(begr:endr), & + rtmCTL%volr_nt2(begr:endr), & + rtmCTL%volr_mch(begr:endr), & + rtmCTL%dvolrdtlnd_nt1(begr:endr), & + rtmCTL%dvolrdtlnd_nt2(begr:endr), & + rtmCTL%dvolrdtocn_nt1(begr:endr), & + rtmCTL%dvolrdtocn_nt2(begr:endr), & + rtmCTL%qsur_nt1(begr:endr), & + rtmCTL%qsur_nt2(begr:endr), & + rtmCTL%qsub_nt1(begr:endr), & + rtmCTL%qsub_nt2(begr:endr), & + rtmCTL%qgwl_nt1(begr:endr), & + rtmCTL%qgwl_nt2(begr:endr), & + rtmCTL%mask(begr:endr), & + rtmCTL%gindex(begr:endr), & + rtmCTL%fthresh(begr:endr), & + rtmCTL%flood(begr:endr), & + rtmCTL%direct(begr:endr,nt_rtm), & + rtmCTL%wh(begr:endr,nt_rtm), & + rtmCTL%wt(begr:endr,nt_rtm), & + rtmCTL%wr(begr:endr,nt_rtm), & + rtmCTL%erout(begr:endr,nt_rtm), & + rtmCTL%qsur(begr:endr,nt_rtm), & + rtmCTL%qsub(begr:endr,nt_rtm), & + rtmCTL%qgwl(begr:endr,nt_rtm), & + rtmCTL%qirrig(begr:endr), & + rtmCTL%qdom_withd(begr:endr), & + rtmCTL%qdom_rf(begr:endr), & + rtmCTL%qdom_actual_withd(begr:endr), & + rtmCTL%qdom_actual_rf(begr:endr), & + rtmCTL%qliv_withd(begr:endr), & + rtmCTL%qliv_rf(begr:endr), & + rtmCTL%qliv_actual_withd(begr:endr), & + rtmCTL%qliv_actual_rf(begr:endr), & + rtmCTL%qelec_withd(begr:endr), & + rtmCTL%qelec_rf(begr:endr), & + rtmCTL%qelec_actual_withd(begr:endr), & + rtmCTL%qelec_actual_rf(begr:endr), & + rtmCTL%qmfc_withd(begr:endr), & + rtmCTL%qmfc_rf(begr:endr), & + rtmCTL%qmfc_actual_withd(begr:endr), & + rtmCTL%qmfc_actual_rf(begr:endr), & + rtmCTL%qmin_withd(begr:endr), & + rtmCTL%qmin_rf(begr:endr), & + rtmCTL%qmin_actual_withd(begr:endr), & + rtmCTL%qmin_actual_rf(begr:endr), & + rtmCTL%qirrig_actual(begr:endr), & stat=ier) if (ier /= 0) then write(iulog,*)'Rtmini ERROR allocation of runoff local arrays' @@ -356,6 +397,26 @@ subroutine RunoffInit(begr, endr, numr) rtmCTL%flood(:) = 0._r8 rtmCTL%direct(:,:) = 0._r8 rtmCTL%qirrig(:) = 0._r8 + rtmCTL%qdom_withd(:) = 0._r8 + rtmCTL%qdom_rf(:) = 0._r8 + rtmCTL%qdom_actual_withd(:) = 0._r8 + rtmCTL%qdom_actual_rf(:) = 0._r8 + rtmCTL%qliv_withd(:) = 0._r8 + rtmCTL%qliv_rf(:) = 0._r8 + rtmCTL%qliv_actual_withd(:) = 0._r8 + rtmCTL%qliv_actual_rf(:) = 0._r8 + rtmCTL%qelec_withd(:) = 0._r8 + rtmCTL%qelec_rf(:) = 0._r8 + rtmCTL%qelec_actual_withd(:) = 0._r8 + rtmCTL%qelec_actual_rf(:) = 0._r8 + rtmCTL%qmfc_withd(:) = 0._r8 + rtmCTL%qmfc_rf(:) = 0._r8 + rtmCTL%qmfc_actual_withd(:) = 0._r8 + rtmCTL%qmfc_actual_rf(:) = 0._r8 + rtmCTL%qmin_withd(:) = 0._r8 + rtmCTL%qmin_rf(:) = 0._r8 + rtmCTL%qmin_actual_withd(:) = 0._r8 + rtmCTL%qmin_actual_rf(:) = 0._r8 rtmCTL%qirrig_actual(:)= 0._r8 rtmCTL%volr_mch(:) = 0._r8 From 2efbdaae28cad407e1b2795a304d8a560d3cd73a Mon Sep 17 00:00:00 2001 From: Sabin Taranu Date: Mon, 12 Dec 2022 12:24:27 +0100 Subject: [PATCH 7/8] Update to src/riverroute/RtmHistFlds.F90 --- src/riverroute/RtmMod.F90 | 316 ++++++++++++++++++++++++++++++-------- 1 file changed, 255 insertions(+), 61 deletions(-) diff --git a/src/riverroute/RtmMod.F90 b/src/riverroute/RtmMod.F90 index c597256..6f867d8 100644 --- a/src/riverroute/RtmMod.F90 +++ b/src/riverroute/RtmMod.F90 @@ -1391,14 +1391,14 @@ subroutine Rtmrun(rstwr,nlend,rdate) ! !LOCAL VARIABLES: !EOP integer :: i, j, n, nr, ns, nt, n2, nf ! indices - real(r8) :: budget_terms(30,nt_rtm) ! BUDGET terms + real(r8) :: budget_terms(40,nt_rtm) ! BUDGET terms ! BUDGET terms 1-10 are for volumes (m3) - ! BUDGET terms 11-30 are for flows (m3/s) + ! BUDGET terms 11-40 are for flows (m3/s) real(r8) :: budget_input, budget_output, budget_volume, budget_total, & budget_euler, budget_eroutlag real(r8),save :: budget_accum(nt_rtm) ! BUDGET accumulator over run integer ,save :: budget_accum_cnt ! counter for budget_accum - real(r8) :: budget_global(30,nt_rtm) ! global budget sum + real(r8) :: budget_global(40,nt_rtm) ! global budget sum logical :: budget_check ! do global budget check real(r8) :: volr_init ! temporary storage to compute dvolrdt real(r8),parameter :: budget_tolerance = 1.0e-6 ! budget tolerance, m3/day @@ -1419,6 +1419,16 @@ subroutine Rtmrun(rstwr,nlend,rdate) real(r8) :: river_volume_minimum ! gridcell area multiplied by average river_depth_minimum [m3] real(r8) :: qgwl_volume ! volume of runoff during time step [m3] real(r8) :: irrig_volume ! volume of irrigation demand during time step [m3] + real(r8) :: dom_withd_volume ! volume of domestic withdrawal during time step [m3] + real(r8) :: liv_withd_volume ! volume of livestock withdrawal during time step [m3] + real(r8) :: elec_withd_volume ! volume of thermoelectric withdrawal during time step [m3] + real(r8) :: mfc_withd_volume ! volume of manufacturing withdrawal during time step [m3] + real(r8) :: min_withd_volume ! volume of mining withdrawal during time step [m3] + real(r8) :: dom_rf_volume ! volume of domestic return flow during time step [m3] + real(r8) :: liv_rf_volume ! volume of livestock return flow during time step [m3] + real(r8) :: elec_rf_volume ! volume of thermoelectric return flow during time step [m3] + real(r8) :: mfc_rf_volume ! volume of manufacturing return flow time step [m3] + real(r8) :: min_rf_volume ! volume of mining return flow during time step [m3] character(len=*),parameter :: subname = '(Rtmrun) ' !----------------------------------------------------------------------- @@ -1454,6 +1464,16 @@ subroutine Rtmrun(rstwr,nlend,rdate) rtmCTL%direct = 0._r8 rtmCTL%flood = 0._r8 rtmCTL%qirrig_actual = 0._r8 + rtmCTL%qdom_actual_withd = 0._r8 + rtmCTL%qliv_actual_withd = 0._r8 + rtmCTL%qelec_actual_withd = 0._r8 + rtmCTL%qmfc_actual_withd = 0._r8 + rtmCTL%qmin_actual_withd = 0._r8 + rtmCTL%qdom_actual_rf = 0._r8 + rtmCTL%qliv_actual_rf = 0._r8 + rtmCTL%qelec_actual_rf = 0._r8 + rtmCTL%qmfc_actual_rf = 0._r8 + rtmCTL%qmin_actual_rf = 0._r8 rtmCTL%runofflnd = spval rtmCTL%runoffocn = spval rtmCTL%dvolrdt = 0._r8 @@ -1462,7 +1482,7 @@ subroutine Rtmrun(rstwr,nlend,rdate) ! BUDGET ! BUDGET terms 1-10 are for volumes (m3) - ! BUDGET terms 11-30 are for flows (m3/s) + ! BUDGET terms 11-40 are for flows (m3/s) ! if (budget_check) then call t_startf('mosartr_budget') do nt = 1,nt_rtm @@ -1474,11 +1494,24 @@ subroutine Rtmrun(rstwr,nlend,rdate) budget_terms(13,nt) = budget_terms(13,nt) + rtmCTL%qsur(nr,nt) budget_terms(14,nt) = budget_terms(14,nt) + rtmCTL%qsub(nr,nt) budget_terms(15,nt) = budget_terms(15,nt) + rtmCTL%qgwl(nr,nt) - budget_terms(17,nt) = budget_terms(17,nt) + rtmCTL%qsur(nr,nt) & + budget_terms(27,nt) = budget_terms(27,nt) + rtmCTL%qsur(nr,nt) & + rtmCTL%qsub(nr,nt)+ rtmCTL%qgwl(nr,nt) if (nt==1) then budget_terms(16,nt) = budget_terms(16,nt) + rtmCTL%qirrig(nr) - budget_terms(17,nt) = budget_terms(17,nt) + rtmCTL%qirrig(nr) + budget_terms(17,nt) = budget_terms(17,nt) + rtmCTL%qdom_withd(nr) + budget_terms(18,nt) = budget_terms(18,nt) + rtmCTL%qliv_withd(nr) + budget_terms(19,nt) = budget_terms(19,nt) + rtmCTL%qelec_withd(nr) + budget_terms(20,nt) = budget_terms(20,nt) + rtmCTL%qmfc_withd(nr) + budget_terms(21,nt) = budget_terms(21,nt) + rtmCTL%qmin_withd(nr) + budget_terms(22,nt) = budget_terms(22,nt) + rtmCTL%qdom_rf(nr) + budget_terms(23,nt) = budget_terms(23,nt) + rtmCTL%qliv_rf(nr) + budget_terms(24,nt) = budget_terms(24,nt) + rtmCTL%qelec_rf(nr) + budget_terms(25,nt) = budget_terms(25,nt) + rtmCTL%qmfc_rf(nr) + budget_terms(26,nt) = budget_terms(26,nt) + rtmCTL%qmin_rf(nr) + budget_terms(27,nt) = budget_terms(27,nt) + rtmCTL%qirrig(nr) & + + rtmCTL%qdom_withd(nr) + rtmCTL%qliv_withd(nr) + rtmCTL%qelec_withd(nr) & + + rtmCTL%qmfc_withd(nr) + rtmCTL%qmin_withd(nr) - rtmCTL%qdom_rf(nr) & + - rtmCTL%qliv_rf(nr) - rtmCTL%qelec_rf(nr) - rtmCTL%qmfc_rf(nr) - rtmCTL%qmin_rf(nr) endif enddo enddo @@ -1500,39 +1533,184 @@ subroutine Rtmrun(rstwr,nlend,rdate) ! Just consider land points and only remove liquid water !----------------------------------- - call t_startf('mosartr_irrig') + call t_startf('mosartr_irrig_and_sector_water_abstractions') nt = 1 - rtmCTL%qirrig_actual = 0._r8 + rtmCTL%qirrig_actual = 0._r8 + rtmCTL%qdom_actual_withd = 0._r8 + rtmCTL%qliv_actual_withd = 0._r8 + rtmCTL%qelec_actual_withd = 0._r8 + rtmCTL%qmfc_actual_withd = 0._r8 + rtmCTL%qmin_actual_withd = 0._r8 + rtmCTL%qdom_actual_rf = 0._r8 + rtmCTL%qliv_actual_rf = 0._r8 + rtmCTL%qelec_actual_rf = 0._r8 + rtmCTL%qmfc_actual_rf = 0._r8 + rtmCTL%qmin_actual_rf = 0._r8 do nr = rtmCTL%begr,rtmCTL%endr ! calculate volume of irrigation flux during timestep irrig_volume = -rtmCTL%qirrig(nr) * coupling_period - - ! compare irrig_volume to main channel storage; + ! calculate volume of domestic withdrawal flux during timestep + dom_withd_volume = -rtmCTL%qdom_withd(nr) * coupling_period + ! calculate volume of livestock withdrawal flux during timestep + liv_withd_volume = -rtmCTL%qliv_withd(nr) * coupling_period + ! calculate volume of thermoelectric withdrawal flux during timestep + elec_withd_volume = -rtmCTL%qelec_withd(nr) * coupling_period + ! calculate volume of manufacturing withdrawal flux during timestep + mfc_withd_volume = -rtmCTL%qmfc_withd(nr) * coupling_period + ! calculate volume of mining withdrawal flux during timestep + min_withd_volume = -rtmCTL%qmin_withd(nr) * coupling_period + + ! calculate volume of domestic return flow flux during timestep + dom_rf_volume = rtmCTL%qdom_rf(nr) * coupling_period + ! calculate volume of livestock return flow flux during timestep + liv_rf_volume = rtmCTL%qliv_rf(nr) * coupling_period + ! calculate volume of thermoelectric return flow flux during timestep + elec_rf_volume = rtmCTL%qelec_rf(nr) * coupling_period + ! calculate volume of manufacturing return flow flux during timestep + mfc_rf_volume = rtmCTL%qmfc_rf(nr) * coupling_period + ! calculate volume of mining return flow flux during timestep + min_rf_volume = rtmCTL%qmin_rf(nr) * coupling_period + + ! compare water usage volume to main channel storage; + ! priority in usage: domestic> livestock > thermoelectric > manufacturing > mining > irrigation ! add overage to subsurface runoff - if(irrig_volume > TRunoff%wr(nr,nt)) then + + if(dom_withd_volume > TRunoff%wr(nr,nt)) then + ! if water missing for domestic sector; satisfy as much as possible; + ! all other sectors get 0.0; + ! to back up potential excessive water consumption which already occured in the land component we send all consumed water to subsurface runoff + ! later in the code there is a procedure on how to deal with negative subsurface runoff, by sending it to direct flux (ocean) + ! in general we should not expect that this deffense mechanism will be required to often; + ! this is because sectoral abstractions are already adjusted in the land model based on volr from coupler; + rtmCTL%qsub(nr,nt) = rtmCTL%qsub(nr,nt) & + + (TRunoff%wr(nr,nt) - (dom_withd_volume - dom_rf_volume) - (liv_withd_volume - liv_rf_volume) & + - (elec_withd_volume - elec_rf_volume) - (mfc_withd_volume - mfc_rf_volume) - (min_withd_volume - min_rf_volume) - irrig_volume) / coupling_period + + TRunoff%qsub(nr,nt) = rtmCTL%qsub(nr,nt) + ! First update the return flow volume by keeping the same ratio between return flow and withdrawal fluxes + dom_rf_volume = TRunoff%wr(nr,nt) * (dom_rf_volume/dom_withd_volume) + liv_rf_volume = 0._r8 + elec_rf_volume = 0._r8 + mfc_rf_volume = 0._r8 + min_rf_volume = 0._r8 + ! Update the withdrawal fluxes based on water availability and sector priority + dom_withd_volume = TRunoff%wr(nr,nt) + liv_withd_volume = 0._r8 + elec_withd_volume = 0._r8 + mfc_withd_volume = 0._r8 + min_withd_volume = 0._r8 + irrig_volume = 0._r8 + + elseif(liv_withd_volume > TRunoff%wr(nr,nt) - dom_withd_volume) then + ! Defense mechanism to protect against negative runoff (see details above) + rtmCTL%qsub(nr,nt) = rtmCTL%qsub(nr,nt) & + + (TRunoff%wr(nr,nt) - (dom_withd_volume - dom_rf_volume) - (liv_withd_volume - liv_rf_volume) & + - (elec_withd_volume - elec_rf_volume) - (mfc_withd_volume - mfc_rf_volume) - (min_withd_volume - min_rf_volume) - irrig_volume) / coupling_period + + TRunoff%qsub(nr,nt) = rtmCTL%qsub(nr,nt) + ! If water missing for livestock sector after satisfying domestic sector, + ! the livestock sector takes what is left and following sectors get 0 + ! First update the return flow volume by keeping the same ratio between return flow and withdrawal fluxes + liv_rf_volume = (TRunoff%wr(nr,nt) - dom_withd_volume) * (liv_rf_volume/liv_withd_volume) + elec_rf_volume = 0._r8 + mfc_rf_volume = 0._r8 + min_rf_volume = 0._r8 + ! Update the withdrawal fluxes based on water availability and sector priority + liv_withd_volume = TRunoff%wr(nr,nt) - dom_withd_volume + elec_withd_volume = 0._r8 + mfc_withd_volume = 0._r8 + min_withd_volume = 0._r8 + irrig_volume = 0._r8 + + elseif(elec_withd_volume > TRunoff%wr(nr,nt) - dom_withd_volume - liv_withd_volume) then + ! Defense mechanism to protect against negative runoff (see details above) + rtmCTL%qsub(nr,nt) = rtmCTL%qsub(nr,nt) & + + (TRunoff%wr(nr,nt) - (dom_withd_volume - dom_rf_volume) - (liv_withd_volume - liv_rf_volume) & + - (elec_withd_volume - elec_rf_volume) - (mfc_withd_volume - mfc_rf_volume) - (min_withd_volume - min_rf_volume) - irrig_volume) / coupling_period + + TRunoff%qsub(nr,nt) = rtmCTL%qsub(nr,nt) + ! If water missing for thermoelectric sector after satisfying domestic and livestock sectors, + ! the thermoelectric sector takes what is left and following sectors get 0 + ! First update the return flow volume by keeping the same ratio between return flow and withdrawal fluxes + elec_rf_volume = (TRunoff%wr(nr,nt) - dom_withd_volume - liv_withd_volume) * (elec_rf_volume/elec_withd_volume) + mfc_rf_volume = 0._r8 + min_rf_volume = 0._r8 + ! Update the withdrawal fluxes based on water availability and sector priority + elec_withd_volume = TRunoff%wr(nr,nt) - dom_withd_volume - liv_withd_volume + mfc_withd_volume = 0._r8 + min_withd_volume = 0._r8 + irrig_volume = 0._r8 + elseif(mfc_withd_volume > TRunoff%wr(nr,nt) - dom_withd_volume - liv_withd_volume - elec_withd_volume) then + ! Defense mechanism to protect against negative runoff (see details above) + rtmCTL%qsub(nr,nt) = rtmCTL%qsub(nr,nt) & + + (TRunoff%wr(nr,nt) - (dom_withd_volume - dom_rf_volume) - (liv_withd_volume - liv_rf_volume) & + - (elec_withd_volume - elec_rf_volume) - (mfc_withd_volume - mfc_rf_volume) - (min_withd_volume - min_rf_volume) - irrig_volume) / coupling_period + + TRunoff%qsub(nr,nt) = rtmCTL%qsub(nr,nt) + ! If water missing for manufacturing sector after satisfying domestic, livestock and thermoelectric sectors, + ! the manufacturing sector takes what is left and following sectors get 0 + ! First update the return flow volume by keeping the same ratio between return flow and withdrawal fluxes + mfc_rf_volume = (TRunoff%wr(nr,nt) - dom_withd_volume - liv_withd_volume - elec_withd_volume) * (mfc_rf_volume/mfc_withd_volume) + min_rf_volume = 0._r8 + ! Update the withdrawal fluxes based on water availability and sector priority + mfc_withd_volume = TRunoff%wr(nr,nt) - dom_withd_volume - liv_withd_volume - elec_withd_volume + min_withd_volume = 0._r8 + irrig_volume = 0._r8 + + elseif(min_withd_volume > TRunoff%wr(nr,nt) - dom_withd_volume - liv_withd_volume - elec_withd_volume - mfc_withd_volume) then + ! Defense mechanism to protect against negative runoff (see details above) rtmCTL%qsub(nr,nt) = rtmCTL%qsub(nr,nt) & - + (TRunoff%wr(nr,nt) - irrig_volume) / coupling_period + + (TRunoff%wr(nr,nt) - (dom_withd_volume - dom_rf_volume) - (liv_withd_volume - liv_rf_volume) & + - (elec_withd_volume - elec_rf_volume) - (mfc_withd_volume - mfc_rf_volume) - (min_withd_volume - min_rf_volume) - irrig_volume) / coupling_period + TRunoff%qsub(nr,nt) = rtmCTL%qsub(nr,nt) - irrig_volume = TRunoff%wr(nr,nt) + ! If water missing for mining sector after satisfying domestic, livestock, thermoelectric and manufacturing sectors, + ! the mining sector takes what is left and following sectors get 0 + ! First update the return flow volume by keeping the same ratio between return flow and withdrawal fluxes + min_rf_volume = (TRunoff%wr(nr,nt) - dom_withd_volume - liv_withd_volume - elec_withd_volume - mfc_withd_volume) * (min_rf_volume/min_withd_volume) + ! Update the withdrawal fluxes based on water availability and sector priority + min_withd_volume = TRunoff%wr(nr,nt) - dom_withd_volume - liv_withd_volume - elec_withd_volume - mfc_withd_volume + irrig_volume = 0._r8 + elseif(irrig_volume > TRunoff%wr(nr,nt) - dom_withd_volume - liv_withd_volume - elec_withd_volume - mfc_withd_volume - min_withd_volume) then + ! Defense mechanism to protect against negative runoff (see details above) + rtmCTL%qsub(nr,nt) = rtmCTL%qsub(nr,nt) & + + (TRunoff%wr(nr,nt) - (dom_withd_volume - dom_rf_volume) - (liv_withd_volume - liv_rf_volume) & + - (elec_withd_volume - elec_rf_volume) - (mfc_withd_volume - mfc_rf_volume) - (min_withd_volume - min_rf_volume) - irrig_volume) / coupling_period + + TRunoff%qsub(nr,nt) = rtmCTL%qsub(nr,nt) + ! Update the withdrawal fluxes based on water availability and sector priority + irrig_volume = TRunoff%wr(nr,nt) - dom_withd_volume - liv_withd_volume - elec_withd_volume - mfc_withd_volume - min_withd_volume endif + !scs: how to deal with sink points / river outlets? ! if (rtmCTL%mask(nr) == 1) then - ! actual irrigation rate [m3/s] + ! update actual irrigation/domestic/thermoelectric/livestick/mining/manufacturing rates [m3/s] ! i.e. the rate actually removed from the main channel - ! if irrig_volume is greater than TRunoff%wr - rtmCTL%qirrig_actual(nr) = - irrig_volume / coupling_period - - ! remove irrigation from wr (main channel) - TRunoff%wr(nr,nt) = TRunoff%wr(nr,nt) - irrig_volume - - + rtmCTL%qdom_actual_withd(nr) = - dom_withd_volume / coupling_period + rtmCTL%qelec_actual_withd(nr) = - elec_withd_volume / coupling_period + rtmCTL%qliv_actual_withd(nr) = - liv_withd_volume / coupling_period + rtmCTL%qmfc_actual_withd(nr) = - mfc_withd_volume / coupling_period + rtmCTL%qmin_actual_withd(nr) = - min_withd_volume / coupling_period + rtmCTL%qirrig_actual_withd(nr) = - irrig_volume / coupling_period + + ! update the actual return flow rates for all sectors (except irrigation) + rtmCTL%qdom_actual_rf(nr) = dom_rf_volume / coupling_period + rtmCTL%qelec_actual_rf(nr) = elec_rf_volume / coupling_period + rtmCTL%qliv_actual_rf(nr) = liv_rf_volume / coupling_period + rtmCTL%qmfc_actual_rf(nr) = mfc_rf_volume / coupling_period + rtmCTL%qmin_actual_rf(nr) = min_rf_volume / coupling_period + + ! remove withdrawals for all sectors and add the return flows to wr (main channel) + ! the return flows are scaled to the actual withdrawals (irrigation have no return flow) + TRunoff%wr(nr,nt) = TRunoff%wr(nr,nt) - (dom_withd_volume - dom_rf_volume) - (elec_withd_volume - elec_rf_volume) & + - (liv_withd_volume - liv_rf_volume) - (mfc_withd_volume - mfc_rf_volume) - (min_withd_volume - min_rf_volume) - irrig_volume !scs endif enddo - call t_stopf('mosartr_irrig') + call t_stopf('mosartr_irrig_and_sector_water_abstractions') !----------------------------------- @@ -1800,9 +1978,9 @@ subroutine Rtmrun(rstwr,nlend,rdate) call t_startf('mosartr_budget') do nt = 1,nt_rtm do nr = rtmCTL%begr,rtmCTL%endr - budget_terms(20,nt) = budget_terms(20,nt) + TRunoff%qsur(nr,nt) & + budget_terms(30,nt) = budget_terms(30,nt) + TRunoff%qsur(nr,nt) & + TRunoff%qsub(nr,nt) + TRunoff%qgwl(nr,nt) - budget_terms(29,nt) = budget_terms(29,nt) + TRunoff%qgwl(nr,nt) + budget_terms(39,nt) = budget_terms(39,nt) + TRunoff%qgwl(nr,nt) enddo enddo call t_stopf('mosartr_budget') @@ -1930,24 +2108,24 @@ subroutine Rtmrun(rstwr,nlend,rdate) budget_terms( 4,nt) = budget_terms( 4,nt) + TRunoff%wt(nr,nt) budget_terms( 6,nt) = budget_terms( 6,nt) + TRunoff%wr(nr,nt) budget_terms( 8,nt) = budget_terms( 8,nt) + TRunoff%wh(nr,nt)*rtmCTL%area(nr) - budget_terms(21,nt) = budget_terms(21,nt) + rtmCTL%direct(nr,nt) + budget_terms(31,nt) = budget_terms(31,nt) + rtmCTL%direct(nr,nt) if (rtmCTL%mask(nr) >= 2) then - budget_terms(18,nt) = budget_terms(18,nt) + rtmCTL%runoff(nr,nt) - budget_terms(26,nt) = budget_terms(26,nt) - erout_prev(nr,nt) - budget_terms(27,nt) = budget_terms(27,nt) + flow(nr,nt) + budget_terms(28,nt) = budget_terms(28,nt) + rtmCTL%runoff(nr,nt) + budget_terms(36,nt) = budget_terms(36,nt) - erout_prev(nr,nt) + budget_terms(37,nt) = budget_terms(37,nt) + flow(nr,nt) else - budget_terms(23,nt) = budget_terms(23,nt) - erout_prev(nr,nt) - budget_terms(24,nt) = budget_terms(24,nt) + flow(nr,nt) + budget_terms(33,nt) = budget_terms(33,nt) - erout_prev(nr,nt) + budget_terms(34,nt) = budget_terms(34,nt) + flow(nr,nt) endif - budget_terms(25,nt) = budget_terms(25,nt) - eroutup_avg(nr,nt) - budget_terms(28,nt) = budget_terms(28,nt) - erlat_avg(nr,nt) - budget_terms(22,nt) = budget_terms(22,nt) + rtmCTL%runoff(nr,nt) + rtmCTL%direct(nr,nt) + eroutup_avg(nr,nt) + budget_terms(35,nt) = budget_terms(35,nt) - eroutup_avg(nr,nt) + budget_terms(38,nt) = budget_terms(38,nt) - erlat_avg(nr,nt) + budget_terms(32,nt) = budget_terms(32,nt) + rtmCTL%runoff(nr,nt) + rtmCTL%direct(nr,nt) + eroutup_avg(nr,nt) enddo enddo nt = 1 do nr = rtmCTL%begr,rtmCTL%endr - budget_terms(19,nt) = budget_terms(19,nt) + rtmCTL%flood(nr) - budget_terms(22,nt) = budget_terms(22,nt) + rtmCTL%flood(nr) + budget_terms(29,nt) = budget_terms(29,nt) + rtmCTL%flood(nr) + budget_terms(32,nt) = budget_terms(32,nt) + rtmCTL%flood(nr) enddo ! accumulate the budget total over the run to make sure it's decreasing on avg @@ -1955,12 +2133,18 @@ subroutine Rtmrun(rstwr,nlend,rdate) do nt = 1,nt_rtm budget_volume = (budget_terms( 2,nt) - budget_terms( 1,nt)) / delt_coupling budget_input = (budget_terms(13,nt) + budget_terms(14,nt) + & - budget_terms(15,nt) + budget_terms(16,nt)) - budget_output = (budget_terms(18,nt) + budget_terms(19,nt) + & - budget_terms(21,nt)) + budget_terms(15,nt) + budget_terms(16,nt) + & + ! Modifications for dom/elec/liv/mfc/min abstractions (both withdrawal and return flow terms) + budget_terms(17,nt) + budget_terms(18,nt) + & + budget_terms(19,nt) + budget_terms(20,nt) + & + budget_terms(21,nt) + budget_terms(22,nt) + & + budget_terms(23,nt) + budget_terms(24,nt) + & + budget_terms(25,nt) + budget_terms(26,nt)) + budget_output = (budget_terms(28,nt) + budget_terms(29,nt) + & + budget_terms(31,nt)) budget_total = budget_volume - budget_input + budget_output budget_accum(nt) = budget_accum(nt) + budget_total - budget_terms(30,nt) = budget_accum(nt)/budget_accum_cnt + budget_terms(40,nt) = budget_accum(nt)/budget_accum_cnt enddo call t_stopf('mosartr_budget') @@ -1969,7 +2153,7 @@ subroutine Rtmrun(rstwr,nlend,rdate) !--- check budget ! convert fluxes from m3/s to m3 by mult by coupling_period - budget_terms(11:30,:) = budget_terms(11:30,:) * delt_coupling + budget_terms(11:40,:) = budget_terms(11:40,:) * delt_coupling ! convert terms from m3 to million m3 budget_terms(:,:) = budget_terms(:,:) * 1.0e-6_r8 @@ -1984,11 +2168,11 @@ subroutine Rtmrun(rstwr,nlend,rdate) budget_volume = (budget_global( 2,nt) - budget_global( 1,nt)) budget_input = (budget_global(13,nt) + budget_global(14,nt) + & budget_global(15,nt)) - budget_output = (budget_global(18,nt) + budget_global(19,nt) + & - budget_global(21,nt)) + budget_output = (budget_global(28,nt) + budget_global(29,nt) + & + budget_global(31,nt)) budget_total = budget_volume - budget_input + budget_output - budget_euler = budget_volume - budget_global(20,nt) + budget_global(18,nt) - budget_eroutlag = budget_global(23,nt) - budget_global(24,nt) + budget_euler = budget_volume - budget_global(30,nt) + budget_global(28,nt) + budget_eroutlag = budget_global(33,nt) - budget_global(34,nt) write(iulog,'(2a,i4)') trim(subname),' tracer = ',nt write(iulog,'(2a,i4,f22.6)') trim(subname),' volume init = ',nt,budget_global(1,nt) write(iulog,'(2a,i4,f22.6)') trim(subname),' volume final = ',nt,budget_global(2,nt) @@ -2003,15 +2187,25 @@ subroutine Rtmrun(rstwr,nlend,rdate) write(iulog,'(2a,i4,f22.6)') trim(subname),' input subsurf = ',nt,budget_global(14,nt) write(iulog,'(2a,i4,f22.6)') trim(subname),' input gwl = ',nt,budget_global(15,nt) write(iulog,'(2a,i4,f22.6)') trim(subname),' input irrig = ',nt,budget_global(16,nt) - write(iulog,'(2a,i4,f22.6)') trim(subname),' input total = ',nt,budget_global(17,nt) - !write(iulog,'(2a,i4,f22.6)') trim(subname),' input check = ',nt,budget_input - budget_global(17,nt) - !write(iulog,'(2a,i4,f22.6)') trim(subname),' input euler = ',nt,budget_global(20,nt) + write(iulog,'(2a,i4,f22.6)') trim(subname),' input dom withd = ',nt,budget_global(17,nt) + write(iulog,'(2a,i4,f22.6)') trim(subname),' input liv withd = ',nt,budget_global(18,nt) + write(iulog,'(2a,i4,f22.6)') trim(subname),' input elec withd= ',nt,budget_global(19,nt) + write(iulog,'(2a,i4,f22.6)') trim(subname),' input mfc withd = ',nt,budget_global(20,nt) + write(iulog,'(2a,i4,f22.6)') trim(subname),' input min withd = ',nt,budget_global(21,nt) + write(iulog,'(2a,i4,f22.6)') trim(subname),' input dom rf = ',nt,budget_global(22,nt) + write(iulog,'(2a,i4,f22.6)') trim(subname),' input liv rf = ',nt,budget_global(23,nt) + write(iulog,'(2a,i4,f22.6)') trim(subname),' input elec rf = ',nt,budget_global(24,nt) + write(iulog,'(2a,i4,f22.6)') trim(subname),' input mfc rf = ',nt,budget_global(25,nt) + write(iulog,'(2a,i4,f22.6)') trim(subname),' input min rf = ',nt,budget_global(26,nt) + write(iulog,'(2a,i4,f22.6)') trim(subname),' input total = ',nt,budget_global(27,nt) + !write(iulog,'(2a,i4,f22.6)') trim(subname),' input check = ',nt,budget_input - budget_global(27,nt) + !write(iulog,'(2a,i4,f22.6)') trim(subname),' input euler = ',nt,budget_global(30,nt) !write(iulog,'(2a)') trim(subname),'----------------' - write(iulog,'(2a,i4,f22.6)') trim(subname),' output flow = ',nt,budget_global(18,nt) - write(iulog,'(2a,i4,f22.6)') trim(subname),' output direct = ',nt,budget_global(21,nt) - write(iulog,'(2a,i4,f22.6)') trim(subname),' output flood = ',nt,budget_global(19,nt) - write(iulog,'(2a,i4,f22.6)') trim(subname),' output total = ',nt,budget_global(22,nt) - !write(iulog,'(2a,i4,f22.6)') trim(subname),' output check = ',nt,budget_output - budget_global(22,nt) + write(iulog,'(2a,i4,f22.6)') trim(subname),' output flow = ',nt,budget_global(28,nt) + write(iulog,'(2a,i4,f22.6)') trim(subname),' output direct = ',nt,budget_global(31,nt) + write(iulog,'(2a,i4,f22.6)') trim(subname),' output flood = ',nt,budget_global(29,nt) + write(iulog,'(2a,i4,f22.6)') trim(subname),' output total = ',nt,budget_global(32,nt) + !write(iulog,'(2a,i4,f22.6)') trim(subname),' output check = ',nt,budget_output - budget_global(32,nt) !write(iulog,'(2a)') trim(subname),'----------------' write(iulog,'(2a,i4,f22.6)') trim(subname),' sum input = ',nt,budget_input write(iulog,'(2a,i4,f22.6)') trim(subname),' sum dvolume = ',nt,budget_volume @@ -2020,16 +2214,16 @@ subroutine Rtmrun(rstwr,nlend,rdate) write(iulog,'(2a,i4,f22.6)') trim(subname),' net (dv-i+o) = ',nt,budget_total !write(iulog,'(2a,i4,f22.6)') trim(subname),' net euler = ',nt,budget_euler write(iulog,'(2a,i4,f22.6)') trim(subname),' eul erout lag = ',nt,budget_eroutlag - !write(iulog,'(2a,i4,f22.6)') trim(subname),' accum (dv-i+o)= ',nt,budget_global(30,nt) + !write(iulog,'(2a,i4,f22.6)') trim(subname),' accum (dv-i+o)= ',nt,budget_global(40,nt) !write(iulog,'(2a)') trim(subname),'----------------' - !write(iulog,'(2a,i4,f22.6)') trim(subname),' erout_prev no= ',nt,budget_global(23,nt) - !write(iulog,'(2a,i4,f22.6)') trim(subname),' erout no= ',nt,budget_global(24,nt) - !write(iulog,'(2a,i4,f22.6)') trim(subname),' eroutup_avg = ',nt,budget_global(25,nt) - !write(iulog,'(2a,i4,f22.6)') trim(subname),' erout_prev out= ',nt,budget_global(26,nt) - !write(iulog,'(2a,i4,f22.6)') trim(subname),' erout out= ',nt,budget_global(27,nt) - !write(iulog,'(2a,i4,f22.6)') trim(subname),' erlateral = ',nt,budget_global(28,nt) - !write(iulog,'(2a,i4,f22.6)') trim(subname),' euler gwl = ',nt,budget_global(29,nt) - !write(iulog,'(2a,i4,f22.6)') trim(subname),' net main chan = ',nt,budget_global(6,nt)-budget_global(5,nt)+budget_global(24,nt)-budget_global(23,nt)+budget_global(27,nt)+budget_global(28,nt)+budget_global(29,nt) + !write(iulog,'(2a,i4,f22.6)') trim(subname),' erout_prev no= ',nt,budget_global(33,nt) + !write(iulog,'(2a,i4,f22.6)') trim(subname),' erout no= ',nt,budget_global(34,nt) + !write(iulog,'(2a,i4,f22.6)') trim(subname),' eroutup_avg = ',nt,budget_global(35,nt) + !write(iulog,'(2a,i4,f22.6)') trim(subname),' erout_prev out= ',nt,budget_global(36,nt) + !write(iulog,'(2a,i4,f22.6)') trim(subname),' erout out= ',nt,budget_global(37,nt) + !write(iulog,'(2a,i4,f22.6)') trim(subname),' erlateral = ',nt,budget_global(38,nt) + !write(iulog,'(2a,i4,f22.6)') trim(subname),' euler gwl = ',nt,budget_global(39,nt) + !write(iulog,'(2a,i4,f22.6)') trim(subname),' net main chan = ',nt,budget_global(6,nt)-budget_global(5,nt)+budget_global(34,nt)-budget_global(33,nt)+budget_global(37,nt)+budget_global(38,nt)+budget_global(39,nt) !write(iulog,'(2a)') trim(subname),'----------------' if ((budget_total-budget_eroutlag) > 1.0e-6) then From 397a5b3c71a8230f939a48a0de4f694d32e9407d Mon Sep 17 00:00:00 2001 From: Sabin Taranu Date: Tue, 13 Dec 2022 13:24:15 +0100 Subject: [PATCH 8/8] Correction to src/riverroute/RtmMod.F90 (wrong array name in one location) --- src/riverroute/RtmMod.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/riverroute/RtmMod.F90 b/src/riverroute/RtmMod.F90 index 6f867d8..4334c85 100644 --- a/src/riverroute/RtmMod.F90 +++ b/src/riverroute/RtmMod.F90 @@ -1694,7 +1694,7 @@ subroutine Rtmrun(rstwr,nlend,rdate) rtmCTL%qliv_actual_withd(nr) = - liv_withd_volume / coupling_period rtmCTL%qmfc_actual_withd(nr) = - mfc_withd_volume / coupling_period rtmCTL%qmin_actual_withd(nr) = - min_withd_volume / coupling_period - rtmCTL%qirrig_actual_withd(nr) = - irrig_volume / coupling_period + rtmCTL%qirrig_actual(nr) = - irrig_volume / coupling_period ! update the actual return flow rates for all sectors (except irrigation) rtmCTL%qdom_actual_rf(nr) = dom_rf_volume / coupling_period