From 8d0b531fbe9fea7c8b98508b7c95ddc342084c22 Mon Sep 17 00:00:00 2001 From: Zhengui Wang Date: Thu, 7 Nov 2024 13:06:49 -0500 Subject: [PATCH] 1). expand iPR for all phytoplankons, 2). add partitioning coefficient for clam induced POM flux to sediment, 3). add parameter to allow wetland effect only on surface layers --- sample_inputs/icm.nml | 6 ++++- src/ICM/icm.F90 | 57 ++++++++++++++++++++++++------------------- src/ICM/icm_init.F90 | 24 +++++++++--------- src/ICM/icm_mod.F90 | 9 ++++--- src/ICM/icm_sfm.F90 | 9 +++---- 5 files changed, 57 insertions(+), 48 deletions(-) diff --git a/sample_inputs/icm.nml b/sample_inputs/icm.nml index 8b29965a4..1545db2fc 100644 --- a/sample_inputs/icm.nml +++ b/sample_inputs/icm.nml @@ -86,7 +86,7 @@ alpha = 8.0 8.0 8.0 !init. slope of P-I curve (g[C]*m2/g[Chl]/E), (iLight=0 !----------------------------------------------------------------------- !options of phytoplankton' predation term !----------------------------------------------------------------------- -iPR = 1 !(0: linear formulation; 1: quadratic) +iPR = 1 1 1 !(0: linear formulation; 1: quadratic) PRR = 0.1 0.2 0.05 !predation rate by higher trophic level (day-1, or day-1.g-1.m3) !----------------------------------------------------------------------- @@ -559,6 +559,7 @@ vp2c = 0.01 0.01 0.01 !phosphorus to carbon ratio [nmarsh] !vKPOM (order same as WSP): (PB1,PB2,PB3) (RPOC,LPOC,DOC) (RPON,LPON,DON,NH4,NO3) (RPOP,LPOP,DOP,PO4) (COD,DO) (SRPOM, PIP) !---------------------- vAw = 0.1 !ratio between marsh area to adjacent water volume (m-1) +vdz = 2.0 !surface layer thickness that wetland effect is applied on (m) vKNO3 = 0.01 !mass-transfer coefficient for NO3 removal (m.day-1) vKTw = 0.069 !temperature dependence for wetland effect (oC-1) vRTw = 20 !reference temperature for wetland effect (oC) @@ -623,6 +624,9 @@ cDoyh = 122 122 122 122 122 274 274 274 274 274 !date range for predation cn2c = 0.167 0.167 0.167 0.167 0.167 !nitrogen to carbon ratio cp2c = 0.011 0.011 0.011 0.011 0.011 !phosphorus to carbon ratio clamm = 0.001 0.001 0.001 0.001 0.001 1.0 1.0 1.0 1.0 1.0 !minimum and maximum clam conc. [g[C].m-2) +cFCM = 0.0 0.5 0.5 !partition of FPOC from clam to different G1-3 classes +cFNM = 0.0 0.5 0.5 !partition of FPON from clam to different G1-3 classes +cFPM = 0.0 0.5 0.5 !partition of FPOP from clam to different G1-3 classes / &ERO diff --git a/src/ICM/icm.F90 b/src/ICM/icm.F90 index 153886c84..b26adf7d4 100644 --- a/src/ICM/icm.F90 +++ b/src/ICM/icm.F90 @@ -245,7 +245,7 @@ subroutine ecosystem(it) !metabolism, predation MT(i,k)=MTR(i)*GP(i,k)+MTB(i)*exp(KTMT(i)*(temp(k)-TMT(i)))*PBS(i,k) - if(iPR==0) then + if(iPR(i)==0) then PR(i,k)=PRR(i)*exp(KTMT(i)*(temp(k)-TMT(i)))*PBS(i,k) else PR(i,k)=PRR(i)*exp(KTMT(i)*(temp(k)-TMT(i)))*PBS(i,k)*PBS(i,k) @@ -701,8 +701,8 @@ subroutine marsh_calc(id,kb,zid,vhtz,vLight) real(rkind),intent(in) :: zid(nvrt),vhtz(nmarsh),vLight(nmarsh) !local variables - integer :: i,j,k,m - real(rkind) :: fR,fT,fST,fI,fN,fP,fDO,mLight,mtemp,msalt,rIK,Kd,xT,xS,vdc,drat + integer :: i,j,k,m,kv + real(rkind) :: fR,fT,fST,fI,fN,fP,fDO,mLight,mtemp,msalt,rIK,Kd,xT,xS,vdc,drat,zv,zrat real(rkind),dimension(nmarsh) :: BMw,BMb,srat,orat real(rkind) :: tdep,dz(nvrt),GP(nmarsh),MT(3,nmarsh) real(rkind),pointer,dimension(:) :: vleaf,vstem,vroot @@ -796,24 +796,35 @@ subroutine marsh_calc(id,kb,zid,vhtz,vLight) !update canopy height call get_canopy(id) elseif(vpatch(id)==1.and.jmarsh==2.and.idry_e(id)==0) then - !currently, the sink due marsh doesn't go to the sediment (todo) for this simple model + !currently, the sink due to marsh doesn't go to the sediment (todo) for this simple model db_vdNO3(id)=0; db_vdDOX(id)=0; db_vdRPOC(id)=0; db_vdLPOC(id)=0; db_vdRPON(id)=0; db_vdLPON(id)=0 db_vdRPOP(id)=0; db_vdLPOP(id)=0; db_vdSRPOC(id)=0; db_vdSRPON(id)=0; db_vdSRPOP(id)=0; db_vdPIP(id)=0 - do k=kb+1,nvrt + + !apply wetland effect only only top vdz meters + kv=kb+1; zv=0 + do k=nvrt,kb+1,-1 + zv=zv+dz(k) + if(zv>=vdz) then + kv=k; exit + endif + enddo !k + zrat=max(tdep/max(zv,1.d-5),1.d0) !concentrate the effect on surface layers + + do k=kv,nvrt fT=exp(vKTw*(temp(k)-vRTw)); fDO=DOX(k)/(vKhDO+DOX(k)) !pre-compute - vdc=vKNO3*vAw*fT*NO3(k); vdwqc(iNO3,k)=-vdc; db_vdNO3(id)=db_vdNO3(id)+vdc*dz(k) !vdNO3: g.m-2.day - vdc=vAw*fDO*fT*vOCw; vdwqc(iDOX,k)=-vdc; db_vdDOX(id)=db_vdDOX(id)+vdc*dz(k) - vdc=vKPOM(iRPOC)*vAW*RPOC(k); vdwqc(iRPOC,k)=-vdc; db_vdRPOC(id)=db_vdRPOC(id)+vdc*dz(k) - vdc=vKPOM(iLPOC)*vAW*LPOC(k); vdwqc(iLPOC,k)=-vdc; db_vdLPOC(id)=db_vdLPOC(id)+vdc*dz(k) - vdc=vKPOM(iRPON)*vAW*RPON(k); vdwqc(iRPON,k)=-vdc; db_vdRPON(id)=db_vdRPON(id)+vdc*dz(k) - vdc=vKPOM(iLPON)*vAW*LPON(k); vdwqc(iLPON,k)=-vdc; db_vdLPON(id)=db_vdLPON(id)+vdc*dz(k) - vdc=vKPOM(iRPOP)*vAW*RPOP(k); vdwqc(iRPOP,k)=-vdc; db_vdRPOP(id)=db_vdRPOP(id)+vdc*dz(k) - vdc=vKPOM(iLPOP)*vAW*LPOP(k); vdwqc(iLPOP,k)=-vdc; db_vdLPOP(id)=db_vdLPOP(id)+vdc*dz(k) + vdc=zrat*vKNO3*vAw*fT*NO3(k); vdwqc(iNO3,k)=-vdc; db_vdNO3(id)=db_vdNO3(id)+vdc*dz(k) !vdNO3: g.m-2.day + vdc=zrat*vAw*fDO*fT*vOCw; vdwqc(iDOX,k)=-vdc; db_vdDOX(id)=db_vdDOX(id)+vdc*dz(k) + vdc=zrat*vKPOM(iRPOC)*vAw*RPOC(k); vdwqc(iRPOC,k)=-vdc; db_vdRPOC(id)=db_vdRPOC(id)+vdc*dz(k) + vdc=zrat*vKPOM(iLPOC)*vAw*LPOC(k); vdwqc(iLPOC,k)=-vdc; db_vdLPOC(id)=db_vdLPOC(id)+vdc*dz(k) + vdc=zrat*vKPOM(iRPON)*vAw*RPON(k); vdwqc(iRPON,k)=-vdc; db_vdRPON(id)=db_vdRPON(id)+vdc*dz(k) + vdc=zrat*vKPOM(iLPON)*vAw*LPON(k); vdwqc(iLPON,k)=-vdc; db_vdLPON(id)=db_vdLPON(id)+vdc*dz(k) + vdc=zrat*vKPOM(iRPOP)*vAw*RPOP(k); vdwqc(iRPOP,k)=-vdc; db_vdRPOP(id)=db_vdRPOP(id)+vdc*dz(k) + vdc=zrat*vKPOM(iLPOP)*vAw*LPOP(k); vdwqc(iLPOP,k)=-vdc; db_vdLPOP(id)=db_vdLPOP(id)+vdc*dz(k) if(iSRM==1) then - vdc=vKPOM(iSRPOC)*vAw*SRPOC(k); vdwqc(iSRPOC,k)=-vdc; db_vdSRPOC(id)=db_vdSRPOC(id)+vdc*dz(k) - vdc=vKPOM(iSRPON)*vAw*SRPON(k); vdwqc(iSRPON,k)=-vdc; db_vdSRPON(id)=db_vdSRPON(id)+vdc*dz(k) - vdc=vKPOM(iSRPOP)*vAw*SRPOP(k); vdwqc(iSRPOP,k)=-vdc; db_vdSRPOP(id)=db_vdSRPOP(id)+vdc*dz(k) - vdc=vKPOM(iPIP)*vAw*PIP(k); vdwqc(iPIP,k) =-vdc; db_vdPIP(id) =db_vdPIP(id) +vdc*dz(k) + vdc=zrat*vKPOM(iSRPOC)*vAw*SRPOC(k); vdwqc(iSRPOC,k)=-vdc; db_vdSRPOC(id)=db_vdSRPOC(id)+vdc*dz(k) + vdc=zrat*vKPOM(iSRPON)*vAw*SRPON(k); vdwqc(iSRPON,k)=-vdc; db_vdSRPON(id)=db_vdSRPON(id)+vdc*dz(k) + vdc=zrat*vKPOM(iSRPOP)*vAw*SRPOP(k); vdwqc(iSRPOP,k)=-vdc; db_vdSRPOP(id)=db_vdSRPOP(id)+vdc*dz(k) + vdc=zrat*vKPOM(iPIP)*vAw*PIP(k); vdwqc(iPIP,k) =-vdc; db_vdPIP(id) =db_vdPIP(id) +vdc*dz(k) endif enddo endif !vpatch(id) @@ -1261,16 +1272,12 @@ subroutine clam_calc(id,kb,wdz) cdwqc(iDOX, kb+1)=o2c*sum(cFc*((ATFC-GP)+MT))/wdz !interaction with sediment layer - cFPOC(id,1:2)=0; cFPON(id,1:2)=0; cFPOP(id,1:2)=0 + cFPOC(id)=0; cFPON(id)=0; cFPOP(id)=0 do i=1,nclam - cFPOC(id,1)=cFPOC(id,1)+sum(cFc(i)*((1-cIF(i))+(1.0-calpha(i,1:4))*cIF(i))*Fr(i)*CLAM(i,id)*PC(1:4))+sum(cFc(i)*(RT+PR)) - cFPON(id,1)=cFPON(id,1)+sum(cFc(i)*((1-cIF(i))+(1.0-calpha(i,1:4))*cIF(i))*Fr(i)*CLAM(i,id)*PN(1:4))+sum(cFc(i)*cn2c(i)*(RT+PR)) - cFPOP(id,1)=cFPOP(id,1)+sum(cFc(i)*((1-cIF(i))+(1.0-calpha(i,1:4))*cIF(i))*Fr(i)*CLAM(i,id)*PP(1:4))+sum(cFc(i)*cp2c(i)*(RT+PR)) + cFPOC(id)=cFPOC(id)+sum(cFc(i)*((1-cIF(i))+(1.0-calpha(i,1:5))*cIF(i))*Fr(i)*CLAM(i,id)*PC(1:5))+sum(cFc(i)*(RT+PR)) + cFPON(id)=cFPON(id)+sum(cFc(i)*((1-cIF(i))+(1.0-calpha(i,1:5))*cIF(i))*Fr(i)*CLAM(i,id)*PN(1:5))+sum(cFc(i)*cn2c(i)*(RT+PR)) + cFPOP(id)=cFPOP(id)+sum(cFc(i)*((1-cIF(i))+(1.0-calpha(i,1:5))*cIF(i))*Fr(i)*CLAM(i,id)*PP(1:5))+sum(cFc(i)*cp2c(i)*(RT+PR)) enddo !i - cFPOC(id,2)=sum(cFc*((1-cIF)+(1.0-calpha(1:nclam,5))*cIF)*Fr*CLAM(1:nclam,id)*PC(5)) - cFPON(id,2)=sum(cFc*((1-cIF)+(1.0-calpha(1:nclam,5))*cIF)*Fr*CLAM(1:nclam,id)*PN(5)) - cFPOP(id,2)=sum(cFc*((1-cIF)+(1.0-calpha(1:nclam,5))*cIF)*Fr*CLAM(1:nclam,id)*PP(5)) - endif !iClam end subroutine clam_calc diff --git a/src/ICM/icm_init.F90 b/src/ICM/icm_init.F90 index ad2c05e2c..6766a0285 100644 --- a/src/ICM/icm_init.F90 +++ b/src/ICM/icm_init.F90 @@ -57,7 +57,7 @@ subroutine read_icm_param(imode) & 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, & - & vFPM,vFW,vKhN,vKhP,valpha,vKe,vSopt,vKs,vInun,vht0,vcrit,v2ht,vc2dw,vn2c,vp2c,vAw, & + & vFPM,vFW,vKhN,vKhP,valpha,vKe,vSopt,vKs,vInun,vht0,vcrit,v2ht,vc2dw,vn2c,vp2c,vAw,vdz, & & vKNO3,vKPOM,vKTW,vRTw,vKhDO,vOCw namelist /SFM/ bdz,bVb,bdiff,bsaltc,bsaltn,bsaltp,bsolid,bKTVp,bKTVd,bVp,bVd,bTR,& & btemp0,bstc0,bSTR0,bThp0,bTox0,bNH40,bNO30,bPO40,bH2S0,bCH40,bPOS0,bSA0,& @@ -70,7 +70,7 @@ subroutine read_icm_param(imode) namelist /BAG/ gpatch0,gBA0,gGPM,gTGP,gKTGP,gMTB,gPRR,gTR,gKTR,galpha,gKSED,gKBA,gKhN,gKhP, & & gp2c,gn2c,gFCP,gFNP,gFPP namelist /CLAM_ICM/ cpatch0,cFc,clam0,cfrmax,cTFR,csalt,cKDO,cDOh,cfTSSm,cRF,cIFmax,cMTB, & - & cTMT,cKTMT,cMRT,cPRR,cHSR,cDoyp,CDoyh,cn2c,cp2c,cKTFR,cKTSS,cTSS,calpha,clamm + & cTMT,cKTMT,cMRT,cPRR,cHSR,cDoyp,CDoyh,cn2c,cp2c,cKTFR,cKTSS,cTSS,calpha,clamm,cFCM,cFNM,cFPM namelist /ERO/ ierosion,erosion,etau,eporo,efrac,ediso,dfrac,dWS_POC if(imode==0) then @@ -141,7 +141,7 @@ subroutine read_icm_param(imode) !init. MARSH module iNmarsh=1; vpatch0=0; vmarsh0=0; vGPM=0; vFAM=0; vTGP=0; vKTGP=0; vFCP=0; vMTB=0; vTMT=0; vKTMT=0; vFCM=0; vFNM=0; vFPM=0; vFW=0; vKhN=0; vKhP=0; valpha=0; vKe=0; vSopt=0; vKs=0; - vInun=0; vht0=0; vcrit=0; v2ht=0; vc2dw=0; vn2c=0; vp2c=0; vAw=0; vKNO3=0; vKPOM=0; + vInun=0; vht0=0; vcrit=0; v2ht=0; vc2dw=0; vn2c=0; vp2c=0; vAw=0; vdz=0; vKNO3=0; vKPOM=0; vKTW=0; vRTw=0; vKhDO=0; vOCw=0 !init. SFM module @@ -163,7 +163,7 @@ subroutine read_icm_param(imode) !init. CLAM module cpatch0=1; clam0=0; cfrmax=0; cTFR=0; csalt=0; cKDO=0; cDOh=0; cfTSSm=0; cRF=0; cIFmax=0 cMTB=0; cTMT=0; cKTMT=0; cMRT=0; cPRR=0; cHSR=0; cDoyp=0; cDoyh=0; cn2c=0; cp2c=0 - cKTFR=0; cKTSS=0; cTSS=0; calpha=0; cFc=1.0; clamm=0 + cKTFR=0; cKTSS=0; cTSS=0; calpha=0; cFc=1.0; clamm=0; cFCM=0; cFNM=0; cFPM=0 !init. ERO module ierosion=0; erosion=0; etau=0; eporo=0; efrac=0; ediso=0; dfrac=0; dWS_POC=0 @@ -782,7 +782,7 @@ subroutine icm_vars_init !------------------------------------------------------------------------------- allocate(cpatch(nea)); cpatch=0 if(iClam==1) then - allocate(CLAM(nclam,nea),cFPOC(nea,2),cFPON(nea,2),cFPOP(nea,2), stat=istat) + allocate(CLAM(nclam,nea),cFPOC(nea),cFPON(nea),cFPOP(nea), stat=istat) if(istat/=0) call parallel_abort('failed in alloc. CLAM') endif @@ -923,11 +923,11 @@ subroutine icm_vars_init sp(m+1)%p1=>vInun; sp(m+2)%p1=>vht0; sp(m+3)%p1=>vcrit; sp(m+4)%p2=>v2ht; sp(m+5)%p1=>vc2dw; m=m+5 sp(m+1)%p1=>vn2c; sp(m+2)%p1=>vp2c; m=m+2 elseif(jmarsh==2) then !simple formulation - pname((m+1):(m+8))=& + pname((m+1):(m+9))=& & (/'vpatch0','vAw ','vKNO3 ','vKPOM ','vKTw ',& - & 'vRTw ','vKhDO ','vOCw '/) + & 'vRTw ','vKhDO ','vOCw ','vdz '/) sp(m+1)%p=>vpatch0; sp(m+2)%p=>vAw; sp(m+3)%p=>vKNO3; sp(m+4)%p1=>vKPOM; sp(m+5)%p=>vKTw; m=m+5 - sp(m+1)%p=>vRTw; sp(m+2)%p=>vKhDO; sp(m+3)%p=>vOCw; m=m+3 + sp(m+1)%p=>vRTw; sp(m+2)%p=>vKhDO; sp(m+3)%p=>vOCw; sp(m+4)%p=>vdz; m=m+4 endif if(iBA==1) then @@ -943,19 +943,19 @@ subroutine icm_vars_init endif if(iClam==1) then - pname((m+1):(m+26))=& + pname((m+1):(m+29))=& & (/'cpatch0','clam0 ','cfrmax ','cTFR ','csalt ',& & 'cKDO ','cDOh ','cfTSSm ','cRF ','cIFmax ',& & 'cMTB ','cTMT ','cKTMT ','cMRT ','cPRR ',& & 'cHSR ','cDoyp ','cDoyh ','cn2c ','cp2c ',& & 'cKTFR ','cKTSS ','cTSS ','calpha ','cFc ',& - & 'clamm '/) + & 'clamm ','cFCM ','cFNM ','cFPM '/) sp(m+1)%p=>cpatch0; sp(m+2)%p1=>clam0; sp(m+3)%p1=>cfrmax; sp(m+4)%p1=>cTFR; sp(m+5)%p1=>csalt; m=m+5 sp(m+1)%p1=>cKDO; sp(m+2)%p1=>cDOh; sp(m+3)%p1=>cfTSSm; sp(m+4)%p1=>cRF; sp(m+5)%p1=>cIFmax; m=m+5 sp(m+1)%p1=>cMTB; sp(m+2)%p1=>cTMT; sp(m+3)%p1=>cKTMT; sp(m+4)%p1=>cMRT; sp(m+5)%p1=>cPRR; m=m+5 sp(m+1)%p1=>cHSR; sp(m+2)%p2=>cDoyp; sp(m+3)%p2=>cDoyh; sp(m+4)%p1=>cn2c; sp(m+5)%p1=>cp2c; m=m+5 sp(m+1)%p2=>cKTFR; sp(m+2)%p2=>cKTSS; sp(m+3)%p2=>cTSS; sp(m+4)%p2=>calpha;sp(m+5)%p1=>cFc; m=m+5 - sp(m+1)%p2=>clamm; m=m+1 + sp(m+1)%p2=>clamm; sp(m+2)%p1=>cFCM; sp(m+3)%p1=>cFNM; sp(m+4)%p1=>cFPM; m=m+4 endif !read spatially varying parameters @@ -1312,7 +1312,7 @@ subroutine check_icm_param() !check global swtiches if(iKe/=0.and.iKe/=1.and.iKe/=2) call parallel_abort('ICM iKe: 0/1/2') if(iLight/=0.and.iLight/=1) call parallel_abort('ICM iLight: 0/1') - if(iPR/=0.and.iPR/=1) call parallel_abort('ICM iPR: 0/1') + do i=1,3; if(iPR(i)/=0.and.iPR(i)/=1) call parallel_abort('ICM iPR: 0/1'); enddo if(iSilica/=0.and.iSilica/=1) call parallel_abort('ICM iSilica: 0/1') if(iZB/=0.and.iZB/=1) call parallel_abort('ICM iZB: 0/1') if(iPh/=0.and.iPh/=1) call parallel_abort('ICM iPh: 0/1') diff --git a/src/ICM/icm_mod.F90 b/src/ICM/icm_mod.F90 index 5771496fc..06479bbbb 100644 --- a/src/ICM/icm_mod.F90 +++ b/src/ICM/icm_mod.F90 @@ -27,7 +27,7 @@ module icm_mod !------------------------------------------------------------------------------- !global switch and variables !------------------------------------------------------------------------------- - integer,save,target :: nsub,iKe,iLight,iPR,iLimit,iSFM,iBA,iRad,isflux,ibflux,idbg(10),iout_icm,nspool_icm + integer,save,target :: nsub,iKe,iLight,iPR(3),iLimit,iSFM,iBA,iRad,isflux,ibflux,idbg(10),iout_icm,nspool_icm integer,save,target :: iSilica,iZB,iPh,iSRM,isav_icm,imarsh_icm,nmarsh,idry_icm,iClam,nclam real(rkind),save,target :: KeC,KeS,KeSalt,Ke0,tss2c,PRR(3),wqc0(29),WSP(29),WSPn(29) real(rkind),save,target,dimension(3) :: alpha @@ -126,7 +126,7 @@ module icm_mod !metabolism partition, growth limit of nutrient,light,salinity,inundation real(rkind),save,target,allocatable,dimension(:) :: vFW,vKhN,vKhP,valpha,vKe,vSopt,vKs,vInun real(rkind),save,target,allocatable :: vht0(:),vcrit(:),v2ht(:,:),vc2dw(:),vn2c(:),vp2c(:) !misc - real(rkind),save,target :: vAw,vKNO3,vKTW,vRTw,vKhDO,vOCw + real(rkind),save,target :: vAw,vdz,vKNO3,vKTW,vRTw,vKhDO,vOCw real(rkind),save,target,allocatable :: vKPOM(:) integer,save,allocatable :: vpatch(:) !marsh regions @@ -174,13 +174,14 @@ module icm_mod !------------------------------------------------------------------------------- !Clam model (CLAM) parameters and variables !------------------------------------------------------------------------------- - real(rkind),save,target :: cpatch0 + real(rkind),save,target :: cpatch0,cFCM(3),cFNM(3),cFPM(3) real(rkind),save,target,allocatable,dimension(:) :: cFc,clam0,cfrmax,cTFR,csalt,cKDO,cDOh,cfTSSm,cRF, & & cIFmax,cMTB,cTMT,cKTMT,cMRT,cPRR,cHSR,cn2c,cp2c real(rkind),save,target,allocatable,dimension(:,:) :: cKTFR,cKTSS,cTSS,calpha,cDoyp,cDoyh,clamm integer,save,allocatable,dimension(:) :: cpatch - real(rkind),save,target,allocatable,dimension(:,:) :: CLAM,cFPOC,cFPON,cFPOP + real(rkind),save,target,allocatable,dimension(:) :: cFPOC,cFPON,cFPOP + real(rkind),save,target,allocatable,dimension(:,:) :: CLAM !debug real(rkind),save,pointer,dimension(:,:) :: db_cfT,db_cfS,db_cfDO,db_cfTSS,db_cFr,db_cIF,db_cTFC,db_cATFC,db_cfN, & diff --git a/src/ICM/icm_sfm.F90 b/src/ICM/icm_sfm.F90 index 668106402..dbaa3c29e 100644 --- a/src/ICM/icm_sfm.F90 +++ b/src/ICM/icm_sfm.F90 @@ -123,13 +123,10 @@ subroutine sfm_calc(id,kb,tdep,wdz,it,isub) !CLAM effect if(iClam==1.and.cpatch(id)==1) then - FPOC(1)=FPOC(1)+cFPOC(id,1) - FPON(1)=FPON(1)+cFPON(id,1) - FPOP(1)=FPOP(1)+cFPOP(id,1) do m=1,3 - FPOC(m)=FPOC(m)+cFPOC(id,2)*bFCM(m) - FPON(m)=FPON(m)+cFPON(id,2)*bFNM(m) - FPOP(m)=FPOP(m)+cFPOP(id,2)*bFPM(m) + FPOC(m)=FPOC(m)+cFPOC(id)*cFCM(m) + FPON(m)=FPON(m)+cFPON(id)*cFNM(m) + FPOP(m)=FPOP(m)+cFPOP(id)*cFPM(m) enddo endif