From 65e2bc8e4bd08078b72d3e9c83fb5f83335664eb Mon Sep 17 00:00:00 2001 From: Zhengui Wang Date: Mon, 18 Nov 2024 15:31:10 -0500 Subject: [PATCH] add 2nd-order formulation for metabolism term in ICM SAV model --- sample_inputs/icm.nml | 1 + src/ICM/icm.F90 | 16 ++++++++-------- src/ICM/icm_init.F90 | 10 +++++----- src/ICM/icm_mod.F90 | 2 +- 4 files changed, 15 insertions(+), 14 deletions(-) diff --git a/sample_inputs/icm.nml b/sample_inputs/icm.nml index 1545db2f..c8d740c4 100644 --- a/sample_inputs/icm.nml +++ b/sample_inputs/icm.nml @@ -449,6 +449,7 @@ inu_ph = 0 !nudge option for pH model !----------------------------------------------------------------------- spatch0 = -999 !region flag for SAV. (1: ON all elem.; -999: spatial) sFc = 1.0 !SAV bed fraction [0,1] +iMTs = 1 !1:1st-order metabolism; 2: 2nd-order metabolism formulation sav0 = 100.0 50.0 10.0 20.0 !init. SAV leaf/stem/root/tuber conc. (g[C].m-2) !growth coefficients diff --git a/src/ICM/icm.F90 b/src/ICM/icm.F90 index b26adf7d..2daeb018 100644 --- a/src/ICM/icm.F90 +++ b/src/ICM/icm.F90 @@ -942,7 +942,7 @@ subroutine sav_calc(id,kb,zid,shtz,tdep,sLight0) !local variables integer :: i,k - real(rkind) :: drat,srat,leafC,stemC,xT,mLight,rIK,Ns,Ps,fT,fI,fN,fP,zm,sFCPm(4) + real(rkind) :: drat,srat,leafC,stemC,xT,mLight,rIK,Ns,Ps,fT,fI,fN,fP,zm,rm,sFCPm(4) real(rkind) :: sfPN,sfPNb,sfPPb,leaf0,stem0,TB,MT0(4),BM,BMb,sfT,sfI,sfN,sfP,s real(rkind) :: efT,eKd,efI,efN,efP,efEP,eGP,eMT,ePR,efPN,zEP0,eMT0 real(rkind),dimension(nvrt) :: sLight,dz,zleaf,zstem,GP,MT1,MT2,zEP @@ -1017,14 +1017,14 @@ subroutine sav_calc(id,kb,zid,shtz,tdep,sLight0) Ps=drat*PO4d(k)/sKhP(1)+bPO4(id)/sKhP(2); fP=Ps/(1+Ps) !grwoth, metabolism,transfer - GP(k)=sGPM*fT*min(fI,fN,fP)*zleaf(k) - MT1(k)=sMTB(1)*exp(sKTMT(1)*(temp(k)-sTMT(1)))*zleaf(k) - MT2(k)=sMTB(2)*exp(sKTMT(2)*(temp(k)-sTMT(2)))*zstem(k) + GP(k)=sGPM*fT*min(fI,fN,fP)*zleaf(k); rm=max(dble(nint(iMTs)),1.d0) + MT1(k)=sMTB(1)*exp(sKTMT(1)*(temp(k)-sTMT(1)))*zleaf(k)**rm + MT2(k)=sMTB(2)*exp(sKTMT(2)*(temp(k)-sTMT(2)))*zstem(k)**rm if(k==nvrt) then - MT0(1)=sMTB(1)*exp(sKTMT(1)*(temp(nvrt)-sTMT(1)))*leaf0 - MT0(2)=sMTB(2)*exp(sKTMT(2)*(temp(nvrt)-sTMT(2)))*stem0 - MT0(3)=sMTB(3)*exp(sKTMT(3)*(temp(kb+1)-sTMT(3)))*sroot - MT0(4)=sMTB(4)*exp(sKTMT(4)*(temp(kb+1)-sTMT(4)))*tuber + MT0(1)=sMTB(1)*exp(sKTMT(1)*(temp(nvrt)-sTMT(1)))*leaf0**rm + MT0(2)=sMTB(2)*exp(sKTMT(2)*(temp(nvrt)-sTMT(2)))*stem0**rm + MT0(3)=sMTB(3)*exp(sKTMT(3)*(temp(kb+1)-sTMT(3)))*sroot**rm + MT0(4)=sMTB(4)*exp(sKTMT(4)*(temp(kb+1)-sTMT(4)))*tuber**rm endif !apply minimum and maximum conc. diff --git a/src/ICM/icm_init.F90 b/src/ICM/icm_init.F90 index 6766a028..3c10629a 100644 --- a/src/ICM/icm_init.F90 +++ b/src/ICM/icm_init.F90 @@ -53,7 +53,7 @@ subroutine read_icm_param(imode) namelist /ZB/ zGPM,zKhG,zTGP,zKTGP,zAG,zRG,zMRT,zMTB,zTMT,zKTMT,zFCP,zFNP,zFPP, & & zFSP,zFCM,zFNM,zFPM,zFSM,zKhDO,zn2c,zp2c,zs2c,z2pr,p2pr namelist /PH_ICM/ ppatch0,inu_ph,pKCACO3,pKCA,pRea - namelist /SAV_ICM/ spatch0,sFc,sav0,sGPM,sTGP,sKTGP,sFAM,sFCP,sMTB,sTMT,sKTMT,sFNM,sFPM,sFCM,sFNMb, & + namelist /SAV_ICM/ spatch0,sFc,iMTs,sav0,sGPM,sTGP,sKTGP,sFAM,sFCP,sMTB,sTMT,sKTMT,sFNM,sFPM,sFCM,sFNMb, & & sFPMb,sFCMb,sKTB,sDoy,sKhN,sKhP,salpha,sKe,shtm,s2ht,sc2dw,sn2c,sp2c,savm,EP0,eGPM,eTGP, & & eKTGP,eKe,ealpha,eMTB,eTMT,eKTMT,ePRR,eFCM,eFNM,eFPM,eFCP,eFNP,eFPP,en2c,ep2c,eKhN,eKhP,eKhE namelist /MARSH/ iNmarsh,vpatch0,vmarsh0,vGPM,vFAM,vTGP,vKTGP,vFCP,vMTB,vTMT,vKTMT,vFCM,vFNM, & @@ -132,7 +132,7 @@ subroutine read_icm_param(imode) ppatch0=0; inu_ph=0; pKCACO3=0; pKCA=0; pRea=0 !init. SAV module - spatch0=0; sFc=0; sav0=0; sGPM=0; sTGP=0; sKTGP=0; sFAM=0; sFCP=0; sMTB=0; sTMT=0; sKTMT=0; + spatch0=0; sFc=0; iMTs=1; sav0=0; sGPM=0; sTGP=0; sKTGP=0; sFAM=0; sFCP=0; sMTB=0; sTMT=0; sKTMT=0; sFNM=0; sFPM=0; sFCM=0; sFNMb=0; sFPMb=0; sFCMb=0; sKhN=0; sKhP=0; salpha=0; sKe=0; shtm=0; s2ht=0; sc2dw=0; sn2c=0; sp2c=0; savm=0; EP0=0; eGPM=0; eTGP=0; eKTGP=0; eKe=0; ealpha=0; eMTB=0; eTMT=0; eKTMT=0; ePRR=0; eFCM=0; eFNM=0; eFPM=0; eFCP=0; eFNP=0; eFPP=0; en2c=0; @@ -880,19 +880,19 @@ subroutine icm_vars_init !SAV,MARSH,BA,pH,CLAM modules if(jsav/=0) then - pname((m+1):(m+29))=& + pname((m+1):(m+30))=& & (/'spatch0','sav0 ','sGPM ','sTGP ','sKTGP ',& & 'sFAM ','sFCP ','sMTB ','sTMT ','sKTMT ',& & 'sFCM ','sFNM ','sFPM ','sFCMb ','sFNMb ',& & 'sFPMb ','sKhN ','sKhP ','salpha ','sKe ',& & 'shtm ','s2ht ','sc2dw ','sn2c ','sp2c ',& - & 'sKTB ','sDoy ','sFc ','savm '/) + & 'sKTB ','sDoy ','sFc ','savm ','iMTs '/) sp(m+1)%p=>spatch0; sp(m+2)%p1=>sav0; sp(m+3)%p=>sGPM; sp(m+4)%p=>sTGP; sp(m+5)%p1=>sKTGP; m=m+5 sp(m+1)%p=>sFAM; sp(m+2)%p1=>sFCP; sp(m+3)%p1=>sMTB; sp(m+4)%p1=>sTMT; sp(m+5)%p1=>sKTMT; m=m+5 sp(m+1)%p1=>sFCM; sp(m+2)%p1=>sFNM; sp(m+3)%p1=>sFPM; sp(m+4)%p1=>sFCMb; sp(m+5)%p1=>sFNMb; m=m+5 sp(m+1)%p1=>sFPMb; sp(m+2)%p1=>sKhN; sp(m+3)%p1=>sKhP; sp(m+4)%p=>salpha; sp(m+5)%p=>sKe; m=m+5 sp(m+1)%p1=>shtm; sp(m+2)%p1=>s2ht; sp(m+3)%p=>sc2dw; sp(m+4)%p=>sn2c; sp(m+5)%p=>sp2c; m=m+5 - sp(m+1)%p=>sKTB; sp(m+2)%p1=>sDoy; sp(m+3)%p=>sFc; sp(m+4)%p2=>savm; m=m+4 + sp(m+1)%p=>sKTB; sp(m+2)%p1=>sDoy; sp(m+3)%p=>sFc; sp(m+4)%p2=>savm; sp(m+5)%p=>iMTs; m=m+5 if(jsav==2) then pname((m+1):(m+21))=& & (/'EP0 ','eGPM ','eTGP ','eKTGP ','eKe ',& diff --git a/src/ICM/icm_mod.F90 b/src/ICM/icm_mod.F90 index 06479bbb..6d2033f7 100644 --- a/src/ICM/icm_mod.F90 +++ b/src/ICM/icm_mod.F90 @@ -100,7 +100,7 @@ module icm_mod !------------------------------------------------------------------------------- !SAV parameters and variables !------------------------------------------------------------------------------- - real(rkind),save,target :: spatch0,sFc,sav0(4) !sav patch, and init conc. + real(rkind),save,target :: spatch0,sFc,iMTs,sav0(4) !sav patch, and init conc. real(rkind),save,target :: sGPM,sTGP,sKTGP(2),sFAM,sFCP(4) !growth related coefficients real(rkind),save,target :: sMTB(4),sTMT(4),sKTMT(4) !meta. coefficients (rate, temp, temp dependence) real(rkind),save,target :: sFCM(4),sFNM(4),sFPM(4) !metabolism of leaf/stem to (RPOM,RLOM,DOM,DIM)