Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

DOM flow in Mosart #64

Draft
wants to merge 38 commits into
base: master
Choose a base branch
from
Draft

DOM flow in Mosart #64

Show file tree
Hide file tree
Changes from 34 commits
Commits
Show all changes
38 commits
Select commit Hold shift + click to select a range
9a1474f
add DOM mass balance code module
devarajun May 7, 2021
02ebfba
call DOM module
devarajun May 8, 2021
bdb1e45
Declare a type for DOM variables
devarajun May 8, 2021
f0e77af
Declare dom variables
devarajun May 8, 2021
3b0d992
add DOC/DON variables to output
devarajun May 16, 2021
dc541fb
add DOC, DON as tracer index
devarajun May 21, 2021
bd2272d
add DOM mass balance subroutines
devarajun May 24, 2021
bfd5fc2
include if condition to avoid division by zero
devarajun May 24, 2021
1089202
correction to include Type for wh
devarajun May 24, 2021
22196b0
update domRUp variable
devarajun Jun 3, 2021
22fa84c
comma addition b/w 1 and cnt
devarajun Jun 3, 2021
2233bbf
correction to eroutUp
devarajun Jun 7, 2021
691552d
add file pycache
devarajun Jun 7, 2021
f11e770
Merge branch 'dev/dom' of https://github.com/MetOs-UiO/MOSART into mo…
Oct 11, 2022
fbf614d
add doc tracer to mosart in new array first step
Nov 18, 2022
d4d9c2b
changes to mosart for DOC
Jan 5, 2023
ae92015
removed domRout and domRin
Jan 5, 2023
68e31f0
corrections to the code, typos and mistakes
Jan 6, 2023
3c18a69
changes to units , kg instead of g and corrections
Jan 10, 2023
9b245f7
comments added
Jan 10, 2023
4d1e1c4
add output concentrations of DOM, corrections of units
Jan 12, 2023
06827d4
sign corrections and fixes
Jan 16, 2023
0123953
various fixes
Jan 20, 2023
bf84861
testing
Jan 20, 2023
1b8ca3b
add surface and subsurface DOM
Jan 23, 2023
ba0506b
PRINT
Jan 31, 2023
03c7f89
unit changes
Feb 1, 2023
6439c7b
correct mass balance equations
Feb 2, 2023
9ce0d68
various changes
Feb 15, 2023
f7653d9
still too much going out of subn
Mar 13, 2023
d1fc844
various changes but one error
Mar 18, 2023
c38038c
changes
Mar 20, 2023
475c98e
changes made it worse
Mar 29, 2023
bdf66e2
corrections and description
Apr 19, 2023
7fa46ae
commented out a bunch of boundaries
May 2, 2023
dae8a26
remove the last boundary statement of DOC mass balance
May 3, 2023
bbd2aa7
added text
Jun 12, 2023
b2debee
starting to remove equations from dommas in hill and subn
Jun 20, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
Binary file added cime_config/__pycache__/buildnmlcpython-39.pyc
View file
Open in desktop
Binary file not shown.
23 changes: 21 additions & 2 deletions src/cpl/nuopc/rof_import_export.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -11,7 +11,7 @@ module rof_import_export
use shr_sys_mod , only : shr_sys_abort
use nuopc_shr_methods , only : chkerr
use RunoffMod , only : rtmCTL, TRunoff, TUnit
use RtmVar , only : iulog, nt_rtm, rtm_tracers
use RtmVar , only : iulog, nt_rtm, rtm_tracers,nt_rtm_dom,rtm_tracers_dom
use RtmSpmd , only : masterproc, mpicom_rof
use RtmTimeManager , only : get_nstep
use nuopc_shr_methods , only : chkerr
Expand Down Expand Up @@ -110,6 +110,8 @@ subroutine advertise_fields(gcomp, flds_scalar_name, do_rtmflood, rc)

call fldlist_add(fldsToRof_num, fldsToRof, trim(flds_scalar_name))
call fldlist_add(fldsToRof_num, fldsToRof, 'Flrl_rofsur')
call fldlist_add(fldsToRof_num, fldsToRof, 'Flrl_subdoc')
call fldlist_add(fldsToRof_num, fldsToRof, 'Flrl_surfdoc')
call fldlist_add(fldsToRof_num, fldsToRof, 'Flrl_rofgwl')
call fldlist_add(fldsToRof_num, fldsToRof, 'Flrl_rofsub')
call fldlist_add(fldsToRof_num, fldsToRof, 'Flrl_rofi')
Expand Down Expand Up @@ -244,7 +246,7 @@ subroutine import_fields( gcomp, rc )
type(ESMF_State) :: importState
integer :: n,nt
integer :: begr, endr
integer :: nliq, nfrz
integer :: nliq, nfrz, ndoc
character(len=*), parameter :: subname='(rof_import_export:import_fields)'
!---------------------------------------------------------------------------

Expand All @@ -267,6 +269,15 @@ subroutine import_fields( gcomp, rc )
call shr_sys_abort()
endif

ndoc = 0
do nt = 1,nt_rtm_dom
if (trim(rtm_tracers_dom(nt)) == 'DOC') ndoc = nt
enddo
if (ndoc == 0) then
write(iulog,*) trim(subname),': ERROR in rtm_tracers_dom DOC ',ndoc,rtm_tracers_dom
call shr_sys_abort()
endif

begr = rtmCTL%begr
endr = rtmCTL%endr

Expand All @@ -277,6 +288,14 @@ subroutine import_fields( gcomp, rc )
do_area_correction=.true., rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return

call state_getimport(importState, 'Flrl_surfdoc', begr, endr, rtmCTL%area, output=rtmCTL%domsur(:,ndoc), &
do_area_correction=.true., rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return

call state_getimport(importState, 'Flrl_subdoc', begr, endr, rtmCTL%area, output=rtmCTL%domsub(:,ndoc), &
do_area_correction=.true., rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return

call state_getimport(importState, 'Flrl_rofsub', begr, endr, rtmCTL%area, output=rtmCTL%qsub(:,nliq), &
do_area_correction=.true., rc=rc)
if (ChkErr(rc,__LINE__,u_FILE_u)) return
Expand Down
71 changes: 71 additions & 0 deletions src/riverroute/DommasbMod.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,71 @@
MODULE DommasbMod
!Description: core code of Dissolved Organic Matter mass balance utilizing river routing models
!Developed by Marius Lambert 02-02-2023
!This module is currently made interact with MOSART routing model
! USES:
use shr_kind_mod , only : r8 => shr_kind_r8
use shr_const_mod , only : SHR_CONST_REARTH, SHR_CONST_PI
use shr_sys_mod , only : shr_sys_abort
use RunoffMod , only : TRunoff, Tdom, TUnit
use RtmVar , only : iulog

implicit none

public hillslopeRoutingDOM
public subnetworkRoutingDOM
public mainchannelRoutingDOM
!--------------------------------------------------------------------

! ! PUBLIC MEMBER FUNCTIONS:
contains

!----------------------------------------------------------------------
subroutine hillslopeRoutingDOM(iunit,nt,ntdom,theDeltaT)
! ! DESCRIPTION: solve the ODEs with Euler algorithm for hillslope routing
implicit none
integer, intent(in) :: iunit, nt, ntdom
real(r8), intent(in) :: theDeltaT
!domsur (kg/m2*s) ,domH (kg/m2), ehout (m/s), domHout (kg/m2*s), qsur (m/s), wh (m)
Tdom%domHout(iunit,ntdom) = -TRunoff%ehout(iunit,nt) * (Tdom%domH(iunit,ntdom) + Tdom%domsur(iunit,ntdom) * theDeltaT)/(TRunoff%wh(iunit,nt)-TRunoff%dwh(iunit,nt)*theDeltaT+TRunoff%qsur(iunit,nt)*theDeltaT)
!we dont want a too high out
Tdom%domHout(iunit,ntdom) = min(-TRunoff%ehout(iunit,nt) * 0.3_r8, Tdom%domHout(iunit,ntdom))
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

ehout is overland flow [m/s] into tributaries. and domHout is concentration [kg/m2/s] going to trigutary. how does different units work out for this comparison?

!cannot be more be less than 0, lower boundary
Tdom%domHout(iunit,ntdom) = max(0._r8,Tdom%domHout(iunit,ntdom))
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

does this happen (negative domHout) because of negative flow?

!cannot be more than available carbon, upper boundary
Tdom%domHout(iunit,ntdom) = min((Tdom%domH(iunit,ntdom)+Tdom%domsur(iunit,ntdom)*theDeltaT)/theDeltaT,Tdom%domHout(iunit,ntdom))

Tdom%domH(iunit,ntdom) = max(0._r8,Tdom%domH(iunit,ntdom) + (Tdom%domsur(iunit,ntdom) - Tdom%domHout(iunit,ntdom))* theDeltaT)
Copy link

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I know why you did this (you cannot output more that you get from surface flow and storage), but this also causes balance off (but maybe ok?)

end subroutine hillslopeRoutingDOM

subroutine subnetworkRoutingDOM(iunit,nt,ntdom,theDeltaT)
! solve the ODEs with Euler algorithm for subnetwork routing
implicit none
integer, intent(in) :: iunit, nt, ntdom
real(r8), intent(in) :: theDeltaT
Tdom%domTout(iunit,ntdom) = -TRunoff%etout(iunit,nt) * (Tdom%domT(iunit,ntdom) + (Tdom%domsub(iunit,ntdom)+Tdom%domHout(iunit,ntdom)) * theDeltaT)/(TRunoff%wt(iunit,nt)-TRunoff%dwt(iunit,nt)*theDeltaT+TRunoff%etin(iunit,nt)*theDeltaT)
Tdom%domTout(iunit,ntdom) = min(-TRunoff%etout(iunit,nt) *0.3_r8,Tdom%domTout(iunit,ntdom))
Tdom%domTout(iunit,ntdom) = max(0._r8,Tdom%domTout(iunit,ntdom))
Tdom%domTout(iunit,ntdom) = min((Tdom%domT(iunit,ntdom)+(Tdom%domsub(iunit,ntdom)+Tdom%domHout(iunit,ntdom))* theDeltaT)/theDeltaT,Tdom%domTout(iunit,ntdom))


Tdom%domT(iunit,ntdom) = max(0._r8,Tdom%domT(iunit,ntdom) + ( Tdom%domsub(iunit,ntdom) + Tdom%domHout(iunit,ntdom) - Tdom%domTout(iunit,ntdom) ) * theDeltaT)
end subroutine subnetworkRoutingDOM

subroutine mainchannelRoutingDOM(iunit,nt,ntdom,theDeltaT)
! solve the ODE with Euler algorithm for main-channel routing
implicit none
integer, intent(in) :: iunit, nt, ntdom
real(r8), intent(in) :: theDeltaT
real(r8) :: temp_gwl
temp_gwl = TRunoff%qgwl(iunit,nt) * TUnit%area(iunit) * TUnit%frac(iunit)

Tdom%domRout(iunit,ntdom) = -TRunoff%erout(iunit,nt) * (Tdom%domR(iunit,ntdom) + (Tdom%domRUp(iunit,ntdom) + Tdom%domToutLat(iunit,ntdom)) * theDeltaT)/(TRunoff%wr(iunit,nt)-TRunoff%dwr(iunit,nt)*theDeltaT+(TRunoff%erlateral(iunit,nt)+TRunoff%erin(iunit,nt)+temp_gwl)*theDeltaT)
Tdom%domRout(iunit,ntdom) = min(-TRunoff%erout(iunit,nt) *0.3_r8,Tdom%domRout(iunit,ntdom))
Tdom%domRout(iunit,ntdom) = max(0._r8,Tdom%domRout(iunit,ntdom))
Tdom%domRout(iunit,ntdom) = min((Tdom%domR(iunit,ntdom)+(Tdom%domRUp(iunit,ntdom) + Tdom%domToutLat(iunit,ntdom))* theDeltaT)/theDeltaT,Tdom%domRout(iunit,ntdom))


Tdom%domR(iunit,ntdom) = max(0._r8,Tdom%domR(iunit,ntdom) + (Tdom%domRUp(iunit,ntdom) + Tdom%domToutLat(iunit,ntdom) - Tdom%domRout(iunit,ntdom)) * theDeltaT)
end subroutine mainchannelRoutingDOM
!-------------------------------------------------------------------------
end MODULE DommasbMod
Loading