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

Lon-lat variable sediment porosity #189

Merged
merged 22 commits into from
Sep 15, 2022
Merged

Lon-lat variable sediment porosity #189

Show file tree
Hide file tree
Changes from 18 commits
Commits
Show all changes
22 commits
Select commit Hold shift + click to select a range
cb2ddf8
ADD oafx to meson
jmaerz Aug 30, 2022
21aa70d
Declaring and allocating 3-D fields for porsol, porwah, porwat
jmaerz Aug 30, 2022
61e0cec
Introducing 2D field for total solid sediment volume: solfu
jmaerz Aug 30, 2022
654b6d5
Initializing 2D and 3D fields for sediment
jmaerz Aug 30, 2022
0198811
deleted real solfu, since it became an array
jmaerz Aug 31, 2022
3f35c6b
exchange k by i,j,k for porwat,porsol,porwah
jmaerz Aug 31, 2022
55948dd
make sediment diffusion coeff spatially (lon,lat,z) dependent of poro…
jmaerz Aug 31, 2022
96311dc
merge with fixed oafx branch
jmaerz Aug 31, 2022
2fedbd9
fix import of zcoefsu
jmaerz Aug 31, 2022
afc9628
fields to 3D
jmaerz Aug 31, 2022
6cb54ab
to 3D and account for lon-lat inhomogeneous porosity
jmaerz Aug 31, 2022
9623c12
to 2D / 3D
jmaerz Aug 31, 2022
f98bb5a
1D -> 3D
jmaerz Aug 31, 2022
fe7788c
Initialize parameters for dissolution/remin already in bodensed - not…
jmaerz Aug 31, 2022
a3366dc
some clean-up
jmaerz Aug 31, 2022
86ff2c8
Introduce potential for reading and application of lon-lat variable s…
jmaerz Sep 2, 2022
cbab888
merge with master - fixed oafx
jmaerz Sep 2, 2022
9e60515
changed loop order in touched loops, where possible
jmaerz Sep 5, 2022
978a95b
Minor changes to meet some review comments
jmaerz Sep 5, 2022
e44c152
Response to reviews
jmaerz Sep 15, 2022
dce4263
merge with master
jmaerz Sep 15, 2022
1a7e6ec
removed comment in sedshi.F90
jmaerz Sep 15, 2022
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
11 changes: 11 additions & 0 deletions cime_config/buildnml
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -237,6 +237,9 @@ set OALKSCEN = "''"
set OALKFILE = "''"
set WITH_DMSPH = .false.
set PI_PH_FILE = "''"
set L_SED_POR = .false.
Copy link
Contributor

Choose a reason for hiding this comment

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

Maybe a more suggestive name would be preferable? For example L_SEDPOR_2D to indicate that this activates a 2d-spatially variable sediment porosity?

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

I agree. Thanks for the suggestion! Maybe l_3Dvarsedpor (variable sediment porosity)?!

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

Done

set SEDPORFILE = "''"


# set DIAPHY defaults
set GLB_FNAMETAG = "'hd','hm','hy'"
Expand Down Expand Up @@ -706,6 +709,7 @@ if ($OCN_GRID == tnx2v1) then
set SCFILE = "'$DIN_LOC_ROOT/ocn/blom/bndcon/sss_clim_core_tnx2v1_20130927.nc'"
set FEDEPFILE = "'$DIN_LOC_ROOT/ocn/blom/bndcon/dustdep_mhw2006_tnx2v1_20130506.nc'"
set SWACLIMFILE = "''"
set SEDPORFILE = "''"
if ($BLOM_RIVER_NUTRIENTS == TRUE) then
set RIVINFILE = "'$DIN_LOC_ROOT/ocn/blom/bndcon/river_nutrients_GNEWS2000_tnx2v1_20170915.nc'"
else
Expand All @@ -726,6 +730,7 @@ else if ($OCN_GRID == tnx1v4) then
set CCFILE = "'$DIN_LOC_ROOT/ocn/blom/bndcon/chlorophyll_concentration_tnx1v4_20170608.nc'"
set SCFILE = "'$DIN_LOC_ROOT/ocn/blom/bndcon/sss_clim_core_tnx1v4_20170604.nc'"
set FEDEPFILE = "'$DIN_LOC_ROOT/ocn/blom/bndcon/dustdep_mhw2006_tnx1v4_20171107.nc'"
set SEDPORFILE = "''"
if ($HAMOCC_VSLS == TRUE) then
set SWACLIMFILE = "'$DIN_LOC_ROOT/ocn/blom/bndcon/Annual_clim_swa_tnx1v4_20210415.nc'"
else
Expand All @@ -752,6 +757,7 @@ else if ($OCN_GRID == tnx0.25v4) then
set SCFILE = "'$DIN_LOC_ROOT/ocn/blom/bndcon/sss_clim_core_tnx0.25v4_20170623.nc'"
set FEDEPFILE = "'$DIN_LOC_ROOT/ocn/blom/bndcon/dustdep_mhw2006_tnx0.25v4_20181004.nc'"
set SWACLIMFILE = "''"
set SEDPORFILE = "''"
if ($BLOM_RIVER_NUTRIENTS == TRUE) then
set RIVINFILE = "'$DIN_LOC_ROOT/ocn/blom/bndcon/river_nutrients_GNEWS2000_tnx0.25v4_20170821.nc'"
else
Expand All @@ -773,6 +779,7 @@ else if ($OCN_GRID == tnx0.125v4) then
set SCFILE = "'$DIN_LOC_ROOT/ocn/blom/bndcon/sss_clim_core_tnx0.125v4_20200722.nc'"
set FEDEPFILE = "'$DIN_LOC_ROOT/ocn/blom/bndcon/dustdep_mhw2006_tnx0.125v4_20200722.nc'"
set SWACLIMFILE = "''"
set SEDPORFILE = "''"
if ($BLOM_RIVER_NUTRIENTS == TRUE) then
set RIVINFILE = "'$DIN_LOC_ROOT/ocn/blom/bndcon/river_nutrients_GNEWS2000_tnx0.125v4_20170821.nc'"
else
Expand Down Expand Up @@ -1445,6 +1452,8 @@ cat >>! $RUNDIR/ocn_in$inststr << EOF
! OXY, NO3, SIL, D13C, and D14C
! WITH_DMSPH : Logical switch to activate DMS calculation as function of pH
! PI_PH_FILE : File name (incl. full path) for surface PI pH input data.
! L_SED_POR : Logical switch to enable lon-lat variable sediment porosity
! SEDPORFILE : File name (incl. full path) for sediment porosity
&BGCNML
ATM_CO2 = $CCSM_CO2_PPMV
FEDEPFILE = $FEDEPFILE
Expand All @@ -1470,6 +1479,8 @@ cat >>! $RUNDIR/ocn_in$inststr << EOF
INID14C = $INID14C
WITH_DMSPH = $WITH_DMSPH
PI_PH_FILE = $PI_PH_FILE
L_SED_POR = $L_SED_POR
SEDPORFILE = $SEDPORFILE
/

Copy link
Contributor

Choose a reason for hiding this comment

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

GitHub won't let me add a comment in an untouched part of the code, but further below at line 1840 there is a file blom.input_data_list created. If L_SED_POR is activated and SEDPORFILE specified, the filename should be added to this file (used by the CESM scripting system to check at submission whether all files are available)

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

I see. Thanks for the explanation! I'll add that.

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

Done

! IO-NAMELIST FOR DIAGNOSTIC iHAMOCC OUTPUT
Expand Down
89 changes: 68 additions & 21 deletions hamocc/bodensed.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -17,7 +17,7 @@
! along with BLOM. If not, see https://www.gnu.org/licenses/.


subroutine bodensed(kpie,kpje,kpke,pddpo)
subroutine bodensed(kpie,kpje,kpke,pddpo,omask,sed_por)
!**********************************************************************
!
!**** *BODENSED* - .
Expand All @@ -44,15 +44,18 @@ subroutine bodensed(kpie,kpje,kpke,pddpo)
!**********************************************************************

use mo_sedmnt, only: calcwei,calfa,clafa,claydens,calcdens,opaldens,opalwei,oplfa,orgdens,orgfa,seddzi,porwat,porwah, &
& porsol,dzs,seddw,sedict,solfu,orgwei
use mo_control_bgc, only: dtbgc,io_stdo_bgc
& porsol,dzs,seddw,sedict,solfu,orgwei,zcoefsu,zcoeflo,disso_sil,silsat,disso_poc,sed_denit,disso_caco3
use mo_control_bgc, only: dtbgc,io_stdo_bgc,l_sed_por
use mo_apply_sedpor,only: apply_sedpor
use mo_param1_bgc, only: ks
use mod_xc, only: mnproc

implicit none

integer, intent(in) :: kpie,kpje,kpke
real, intent(in) :: pddpo(kpie,kpje,kpke)
real, intent(in) :: omask(kpie,kpje)
real, intent(in) :: sed_por(kpie,kpje,ks)

! Local variables
integer :: i,j,k
Expand All @@ -79,33 +82,60 @@ subroutine bodensed(kpie,kpje,kpke,pddpo)
write(io_stdo_bgc,*) ' '
endif

porwat(1) = 0.85
porwat(2) = 0.83
porwat(3) = 0.8
porwat(4) = 0.79
porwat(5) = 0.77
porwat(6) = 0.75
porwat(7) = 0.73
porwat(8) = 0.7
porwat(9) = 0.68
porwat(10) = 0.66
porwat(11) = 0.64
porwat(12) = 0.62
! this initialization can be done later via reading a porosity map
if (l_sed_por)then
! lon-lat variable sediment porosity from input file
call apply_sedpor(kpie,kpje,ks,omask,sed_por)
Copy link
Contributor

Choose a reason for hiding this comment

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

I find it a bit overkill to define an extra module for ~10 lines of code (essentially porwat=sedpor) that could be just placed directly here.

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

:-) - agree - I was debating with myself - I'll change it.

else
porwat(:,:,1) = 0.85
porwat(:,:,2) = 0.83
porwat(:,:,3) = 0.8
porwat(:,:,4) = 0.79
porwat(:,:,5) = 0.77
porwat(:,:,6) = 0.75
porwat(:,:,7) = 0.73
porwat(:,:,8) = 0.7
porwat(:,:,9) = 0.68
porwat(:,:,10) = 0.66
porwat(:,:,11) = 0.64
porwat(:,:,12) = 0.62
endif

if (mnproc == 1) then
write(io_stdo_bgc,*) 'Pore water in sediment: ',porwat
write(io_stdo_bgc,*) 'Pore water in sediment initialized'
endif

seddzi(1) = 500.
do k = 1, ks
porsol(k) = 1. - porwat(k)
if(k >= 2) porwah(k) = 0.5 * (porwat(k) + porwat(k-1))
if(k == 1) porwah(k) = 0.5 * (1. + porwat(1))
seddzi(k+1) = 1. / dzs(k+1)
seddw(k) = 0.5 * (dzs(k) + dzs(k+1))
do j = 1, kpje
do i = 1, kpie
porsol(i,j,k) = 1. - porwat(i,j,k)
if(k >= 2) porwah(i,j,k) = 0.5 * (porwat(i,j,k) + porwat(i,j,k-1))
if(k == 1) porwah(i,j,k) = 0.5 * (1. + porwat(i,j,1))
enddo
enddo
enddo

sedict = 1.e-9 * dtbgc ! Moecular diffusion coefficient
Copy link
Contributor

Choose a reason for hiding this comment

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

Suggested change
sedict = 1.e-9 * dtbgc ! Moecular diffusion coefficient
sedict = 1.e-9 * dtbgc ! Molecular diffusion coefficient

typo

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

Done.

! Dissolution rate constant of opal (disso) [1/(kmol Si(OH)4/m3)*1/sec]
! THIS NEEDS TO BE CHANGED TO disso=3.e-8! THIS IS ONLY KEPT FOR THE MOMENT
! FOR BACKWARDS COMPATIBILITY
!disso=3.e-8 ! (2011-01-04) EMR
!disso=1.e-6 ! test vom 03.03.04 half live sil ca. 20.000 yr
disso_sil = 1.e-6*dtbgc
Copy link
Contributor

Choose a reason for hiding this comment

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

Suggested change
! Dissolution rate constant of opal (disso) [1/(kmol Si(OH)4/m3)*1/sec]
! THIS NEEDS TO BE CHANGED TO disso=3.e-8! THIS IS ONLY KEPT FOR THE MOMENT
! FOR BACKWARDS COMPATIBILITY
!disso=3.e-8 ! (2011-01-04) EMR
!disso=1.e-6 ! test vom 03.03.04 half live sil ca. 20.000 yr
disso_sil = 1.e-6*dtbgc
! Dissolution rate constant of opal (disso) [1/(kmol Si(OH)4/m3)*1/sec]
! THIS NEEDS TO BE CHANGED TO disso=3.e-8! THIS IS ONLY KEPT FOR THE MOMENT
! FOR BACKWARDS COMPATIBILITY
!disso_sil = 3.e-8*dtbgc ! (2011-01-04) EMR
!disso_sil = 1.e-6*dtbgc ! test vom 03.03.04 half live sil ca. 20.000 yr
disso_sil = 1.e-6*dtbgc

Maybe a bit pedantic, but I think it is best if the commented lines also reflect the code change.

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

Correct, I can do this.

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

Done.

! Silicate saturation concentration is 1 mol/m3
silsat = 0.001

! Degradation rate constant of POP (disso) [1/(kmol O2/m3)*1/sec]
disso_poc = 0.01 / 86400. * dtbgc ! disso=3.e-5 was quite high
Copy link
Contributor

Choose a reason for hiding this comment

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

Suggested change
disso_poc = 0.01 / 86400. * dtbgc ! disso=3.e-5 was quite high
disso_poc = 0.01 / 86400. * dtbgc

Is this comment relevant, do we need to keep it around?

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

I can't tell - maybe @JorgSchwinger made some experiences that allow for better judgement?


sedict = 1.e-9 * dtbgc
! Denitrification rate constant of POP (disso) [1/sec]
sed_denit = 0.01/86400. * dtbgc !ik denit = 1.e-6*dtbgc
Copy link
Contributor

Choose a reason for hiding this comment

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

Suggested change
sed_denit = 0.01/86400. * dtbgc !ik denit = 1.e-6*dtbgc
sed_denit = 0.01/86400. * dtbgc

Same as above

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

@JorgSchwinger - any reason to keep the comment?

Copy link
Contributor

Choose a reason for hiding this comment

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

No, I don't think so, too.

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

Removed.


! Dissolution rate constant of CaCO3 (disso) [1/(kmol CO3--/m3)*1/sec]
disso_caco3 = 1.e-7 * dtbgc

! ******************************************************************
! densities etc. for SEDIMENT SHIFTING
Expand All @@ -131,9 +161,26 @@ subroutine bodensed(kpie,kpje,kpke,pddpo)

! determine total solid sediment volume
solfu = 0.
do i = 1, kpie
do j = 1, kpje
do k = 1, ks
solfu = solfu + seddw(k) * porsol(k)
solfu(i,j) = solfu(i,j) + seddw(k) * porsol(i,j,k)
enddo
enddo
enddo

! Initialize porosity-dependent diffusion coefficients of sediment
zcoefsu(:,:,0) = 0.0
do k = 1,ks
do j = 1, kpje
do i = 1, kpie
! sediment diffusion coefficient * 1/dz * fraction of pore water at half depths
zcoefsu(i,j,k ) = -sedict * seddzi(k) * porwah(i,j,k)
zcoeflo(i,j,k-1) = -sedict * seddzi(k) * porwah(i,j,k) ! why the same ?
enddo
enddo
enddo
zcoeflo(:,:,ks) = 0.0 ! diffusion coefficient for bottom sediment layer


end subroutine bodensed
23 changes: 7 additions & 16 deletions hamocc/dipowa.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -56,7 +56,7 @@ subroutine dipowa(kpie,kpje,kpke,omask,lspin)
!**********************************************************************

use mo_carbch, only: ocetra, sedfluxo
use mo_sedmnt, only: powtra,porwat,porwah,sedict,seddw,seddzi
use mo_sedmnt, only: powtra,porwat,porwah,sedict,seddw,seddzi,zcoefsu,zcoeflo
use mo_param1_bgc, only: ks,npowtra
use mo_vgrid, only: kbo,bolay
#ifdef cisonew
Expand All @@ -77,22 +77,13 @@ subroutine dipowa(kpie,kpje,kpke,omask,lspin)
integer :: iv_oc ! index of ocetra in powtra loop

real :: sedb1(kpie,0:ks,npowtra) ! ????
real :: zcoefsu(0:ks),zcoeflo(0:ks) ! diffusion coefficients (upper/lower)
real :: tredsy(kpie,0:kpke,3) ! redsy for 'reduced system'?
real :: aprior ! start value of oceanic tracer in bottom layer

!ik accelerated sediment
!ik needed for boundary layer ventilation in fast sediment routine
real :: bolven(kpie) ! bottom layer ventilation rate

zcoefsu(0) = 0.0
do k = 1,ks
! sediment diffusion coefficient * 1/dz * fraction of pore water at half depths
zcoefsu(k ) = -sedict * seddzi(k) * porwah(k)
zcoeflo(k-1) = -sedict * seddzi(k) * porwah(k) ! why the same ?
enddo
zcoeflo(ks) = 0.0 ! diffusion coefficient for bottom sediment layer

!$OMP PARALLEL DO &
!$OMP&PRIVATE(i,k,iv,l,bolven,tredsy,sedb1,aprior,iv_oc)
j_loop: do j=1,kpje
Expand All @@ -104,8 +95,8 @@ subroutine dipowa(kpie,kpje,kpke,omask,lspin)

k = 0
do i = 1,kpie
tredsy(i,k,1) = zcoefsu(k)
tredsy(i,k,3) = zcoeflo(k)
tredsy(i,k,1) = zcoefsu(i,j,k)
tredsy(i,k,3) = zcoeflo(i,j,k)
tredsy(i,k,2) = bolven(i)*bolay(i,j) - tredsy(i,k,1) - tredsy(i,k,3)
! dz(kbo) - diff upper - diff lower
enddo
Expand All @@ -128,17 +119,17 @@ subroutine dipowa(kpie,kpje,kpke,omask,lspin)

do k = 1,ks
do i = 1,kpie
tredsy(i,k,1) = zcoefsu(k)
tredsy(i,k,3) = zcoeflo(k)
tredsy(i,k,2) = seddw(k)*porwat(k) -tredsy(i,k,1) -tredsy(i,k,3)
tredsy(i,k,1) = zcoefsu(i,j,k)
tredsy(i,k,3) = zcoeflo(i,j,k)
tredsy(i,k,2) = seddw(k)*porwat(i,j,k) -tredsy(i,k,1) -tredsy(i,k,3)
enddo
enddo

do iv = 1,npowtra
do k = 1,ks
do i = 1,kpie
! tracer_concentration(k[1:ks]) * porewater fraction(k) * dz(k)
sedb1(i,k,iv) = powtra(i,j,k,iv) * porwat(k) * seddw(k)
sedb1(i,k,iv) = powtra(i,j,k,iv) * porwat(i,j,k) * seddw(k)
enddo
enddo
enddo
Expand Down
10 changes: 6 additions & 4 deletions hamocc/hamocc_init.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -46,7 +46,7 @@ subroutine hamocc_init(read_rest,rstfnm_hamocc)
& do_ndep,do_rivinpt,do_oalk,do_sedspinup, &
& sedspin_yr_s,sedspin_yr_e,sedspin_ncyc, &
& dtb,dtbgc,io_stdo_bgc,ldtbgc, &
& ldtrunbgc,ndtdaybgc,with_dmsph
& ldtrunbgc,ndtdaybgc,with_dmsph,l_sed_por
use mo_param1_bgc, only: ks,nsedtra,npowtra
use mo_carbch, only: alloc_mem_carbch,ocetra,atm,atm_co2
use mo_biomod, only: alloc_mem_biomod
Expand All @@ -58,6 +58,7 @@ subroutine hamocc_init(read_rest,rstfnm_hamocc)
use mo_read_ndep, only: ini_read_ndep,ndepfile
use mo_read_oafx, only: ini_read_oafx,oalkfile,oalkscen
use mo_read_pi_ph, only: ini_pi_ph,pi_ph_file
use mo_read_sedpor, only: ini_read_sedpor,sedporfile,sed_por
use mo_clim_swa, only: ini_swa_clim,swaclimfile
use mo_Gdata_read, only: inidic,inialk,inipo4,inioxy,inino3, &
& inisil,inid13c,inid14c
Expand All @@ -82,7 +83,7 @@ subroutine hamocc_init(read_rest,rstfnm_hamocc)
& sedspin_yr_e,sedspin_ncyc, &
& inidic,inialk,inipo4,inioxy,inino3,inisil, &
& inid13c,inid14c,swaclimfile, &
& with_dmsph,pi_ph_file
& with_dmsph,pi_ph_file,l_sed_por,sedporfile
!
! --- Set io units and some control parameters
!
Expand Down Expand Up @@ -171,8 +172,9 @@ subroutine hamocc_init(read_rest,rstfnm_hamocc)
call set_vgrid(idm,jdm,kdm,bgc_dp)
!
! --- Initialize sediment layering
!
CALL BODENSED(idm,jdm,kdm,bgc_dp)
! First init the porosity, then apply it in bodensed
CALL ini_read_sedpor(idm,jdm,ks,omask)
CALL BODENSED(idm,jdm,kdm,bgc_dp,omask,sed_por)
!
! --- Initialize parameters, sediment and ocean tracer.
!
Expand Down
8 changes: 4 additions & 4 deletions hamocc/inventory_bgc.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -143,7 +143,7 @@ SUBROUTINE INVENTORY_BGC(kpie,kpje,kpke,dlxp,dlyp,ddpo,omask,iogrp)
DO j=1,kpje
DO i=1,kpie
ztmp1(i,j) = ztmp1(i,j) + omask(i,j)*seddw(k) &
& *dlxp(i,j)*dlyp(i,j)*porwat(k)
& *dlxp(i,j)*dlyp(i,j)*porwat(i,j,k)
ENDDO
ENDDO
ENDDO
Expand All @@ -155,7 +155,7 @@ SUBROUTINE INVENTORY_BGC(kpie,kpje,kpke,dlxp,dlyp,ddpo,omask,iogrp)
DO k=1,ks
DO j=1,kpje
DO i=1,kpie
vol = seddw(k)*dlxp(i,j)*dlyp(i,j)*porwat(k)
vol = seddw(k)*dlxp(i,j)*dlyp(i,j)*porwat(i,j,k)
ztmp1(i,j)= ztmp1(i,j) + omask(i,j)*powtra(i,j,k,l)*vol
ENDDO
ENDDO
Expand All @@ -174,7 +174,7 @@ SUBROUTINE INVENTORY_BGC(kpie,kpje,kpke,dlxp,dlyp,ddpo,omask,iogrp)
DO k=1,ks
DO j=1,kpje
DO i=1,kpie
vol = porsol(k)*seddw(k)*dlxp(i,j)*dlyp(i,j)
vol = porsol(i,j,k)*seddw(k)*dlxp(i,j)*dlyp(i,j)
ztmp1(i,j) = ztmp1(i,j) + omask(i,j)*sedlay(i,j,k,l)*vol
ENDDO
ENDDO
Expand All @@ -187,7 +187,7 @@ SUBROUTINE INVENTORY_BGC(kpie,kpje,kpke,dlxp,dlyp,ddpo,omask,iogrp)
DO k=1,ks
DO j=1,kpje
DO i=1,kpie
vol = porsol(k)*seddw(k)*dlxp(i,j)*dlyp(i,j)
vol = porsol(i,j,k)*seddw(k)*dlxp(i,j)*dlyp(i,j)
ztmp1(i,j) = ztmp1(i,j) + omask(i,j)*sedhpl(i,j,k)*vol
ENDDO
ENDDO
Expand Down
2 changes: 2 additions & 0 deletions hamocc/meson.build
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -21,6 +21,7 @@ sources += files(
'mo_apply_ndep.F90',
'mo_apply_rivin.F90',
'mo_apply_oafx.F90',
'mo_apply_sedpor.F90',
'mo_bgcmean.F90',
'mo_biomod.F90',
'mo_carbch.F90',
Expand All @@ -34,6 +35,7 @@ sources += files(
'mo_read_pi_ph.F90',
'mo_read_rivin.F90',
'mo_read_oafx.F90',
'mo_read_sedpor.F90',
'mo_sedmnt.F90',
'mo_vgrid.F90',
'ncout_hamocc.F90',
Expand Down
Loading