From aa737aafac4e9560b2eb54bda2cc5022f3cd6ac2 Mon Sep 17 00:00:00 2001 From: JorgSchwinger Date: Thu, 3 Aug 2023 15:46:12 +0200 Subject: [PATCH] Add capability to run C-isotopes with sediment code switched on (tested technically but not validated, build script will still refuse to run this combination of options) --- hamocc/aufr_bgc.F90 | 33 ++++++++++++++++++++++++++------- hamocc/aufw_bgc.F90 | 21 +++++++++++++++++++++ hamocc/carchm.F90 | 11 +++++++++++ hamocc/sedshi.F90 | 37 +++++++++++++++++++++---------------- 4 files changed, 79 insertions(+), 23 deletions(-) diff --git a/hamocc/aufr_bgc.F90 b/hamocc/aufr_bgc.F90 index 842720ad..5c245142 100644 --- a/hamocc/aufr_bgc.F90 +++ b/hamocc/aufr_bgc.F90 @@ -482,11 +482,14 @@ SUBROUTINE AUFR_BGC(kpie,kpje,kpke,ntr,ntrbgc,itrbgc,trc, & CALL read_netcdf_var(ncid,'powasi',powtra2(1,1,1,ipowasi),2*ks,0,iotype) #ifdef cisonew IF(lread_iso) THEN - ! Burial fields for c-isotopes still missing CALL read_netcdf_var(ncid,'ssso13',sedlay2(1,1,1,issso13),2*ks,0,iotype) CALL read_netcdf_var(ncid,'ssso14',sedlay2(1,1,1,issso14),2*ks,0,iotype) CALL read_netcdf_var(ncid,'sssc13',sedlay2(1,1,1,isssc13),2*ks,0,iotype) CALL read_netcdf_var(ncid,'sssc14',sedlay2(1,1,1,isssc14),2*ks,0,iotype) + CALL read_netcdf_var(ncid,'bur_o13',burial2(1,1,1,issso13),2,0,iotype) + CALL read_netcdf_var(ncid,'bur_o14',burial2(1,1,1,issso14),2,0,iotype) + CALL read_netcdf_var(ncid,'bur_c13',burial2(1,1,1,isssc13),2,0,iotype) + CALL read_netcdf_var(ncid,'bur_c14',burial2(1,1,1,isssc14),2,0,iotype) CALL read_netcdf_var(ncid,'powc13',powtra2(1,1,1,ipowc13),2*ks,0,iotype) CALL read_netcdf_var(ncid,'powc14',powtra2(1,1,1,ipowc14),2*ks,0,iotype) ENDIF @@ -586,18 +589,34 @@ SUBROUTINE AUFR_BGC(kpie,kpje,kpke,ntr,ntrbgc,itrbgc,trc, & ENDDO #ifndef sedbypass ! Burial fields for c-isotopes still missing + !JT added burial loop below 20.06.2023 DO k=1,2*ks DO j=1,kpje DO i=1,kpie IF(omask(i,j) .GT. 0.5) THEN rco213=ocetra(i,j,kbo(i,j),isco213)/(ocetra(i,j,kbo(i,j),isco212)+safediv) rco214=ocetra(i,j,kbo(i,j),isco214)/(ocetra(i,j,kbo(i,j),isco212)+safediv) - powtra2(i,j,k,ipowc13)=powtra2(i,j,k,ipowc)*rco213*bifr13 - powtra2(i,j,k,ipowc14)=powtra2(i,j,k,ipowc)*rco214*bifr14 - sedlay2(i,j,k,issso13)=sedlay2(i,j,k,issso)*rco213*bifr13 - sedlay2(i,j,k,issso14)=sedlay2(i,j,k,issso)*rco214*bifr14 - sedlay2(i,j,k,isssc13)=sedlay2(i,j,k,isssc)*rco213 - sedlay2(i,j,k,isssc14)=sedlay2(i,j,k,isssc)*rco214 + powtra2(i,j,k,ipowc13)=powtra2(i,j,k,ipowaic)*rco213 + powtra2(i,j,k,ipowc14)=powtra2(i,j,k,ipowaic)*rco214 + sedlay2(i,j,k,issso13)=sedlay2(i,j,k,issso12)*rco213*bifr13 + sedlay2(i,j,k,issso14)=sedlay2(i,j,k,issso12)*rco214*bifr14 + sedlay2(i,j,k,isssc13)=sedlay2(i,j,k,isssc12)*rco213 + sedlay2(i,j,k,isssc14)=sedlay2(i,j,k,isssc12)*rco214 + ENDIF + ENDDO + ENDDO + ENDDO + + DO k=1,2 + DO j=1,kpje + DO i=1,kpie + IF(omask(i,j) .GT. 0.5) THEN + rco213=ocetra(i,j,kbo(i,j),isco213)/(ocetra(i,j,kbo(i,j),isco212)+safediv) + rco214=ocetra(i,j,kbo(i,j),isco214)/(ocetra(i,j,kbo(i,j),isco212)+safediv) + burial2(i,j,k,issso13)=burial2(i,j,k,issso12)*rco213*bifr13 + burial2(i,j,k,issso14)=burial2(i,j,k,issso12)*rco214*bifr14 + burial2(i,j,k,isssc13)=burial2(i,j,k,isssc12)*rco213 + burial2(i,j,k,isssc14)=burial2(i,j,k,isssc12)*rco214 ENDIF ENDDO ENDDO diff --git a/hamocc/aufw_bgc.F90 b/hamocc/aufw_bgc.F90 index 90949313..ee2132b6 100644 --- a/hamocc/aufw_bgc.F90 +++ b/hamocc/aufw_bgc.F90 @@ -733,6 +733,23 @@ SUBROUTINE AUFW_BGC(kpie,kpje,kpke,ntr,ntrbgc,itrbgc,trc, & & 7,'kg/m**2',20,'Burial layer of clay', & & rmissing,93,io_stdo_bgc) +#ifdef cisonew + CALL NETCDF_DEF_VARDB(ncid,8,'bur_o13',3,ncdimst,ncvarid, & + & 9,'kmol/m**2',27,'Burial layer of organic 13C', & + & rmissing,94,io_stdo_bgc) + + CALL NETCDF_DEF_VARDB(ncid,8,'bur_o14',3,ncdimst,ncvarid, & + & 9,'kmol/m**2',27,'Burial layer of organic 14C', & + & rmissing,95,io_stdo_bgc) + + CALL NETCDF_DEF_VARDB(ncid,8,'bur_c13',3,ncdimst,ncvarid, & + & 9,'kmol/m**2',23,'Burial layer of Ca13CO3', & + & rmissing,96,io_stdo_bgc) + + CALL NETCDF_DEF_VARDB(ncid,8,'bur_c14',3,ncdimst,ncvarid, & + & 9,'kmol/m**2',23,'Burial layer of Ca14CO3', & + & rmissing,97,io_stdo_bgc) +#endif #endif /* sedbypass */ ! @@ -897,6 +914,10 @@ SUBROUTINE AUFW_BGC(kpie,kpje,kpke,ntr,ntrbgc,itrbgc,trc, & CALL write_netcdf_var(ncid,'ssso14',sedlay2(1,1,1,issso14),2*ks,0) CALL write_netcdf_var(ncid,'sssc13',sedlay2(1,1,1,isssc13),2*ks,0) CALL write_netcdf_var(ncid,'sssc14',sedlay2(1,1,1,isssc14),2*ks,0) + CALL write_netcdf_var(ncid,'bur_o13',burial2(1,1,1,issso13),2,0) + CALL write_netcdf_var(ncid,'bur_o14',burial2(1,1,1,issso14),2,0) + CALL write_netcdf_var(ncid,'bur_c13',burial2(1,1,1,isssc13),2,0) + CALL write_netcdf_var(ncid,'bur_c14',burial2(1,1,1,isssc14),2,0) CALL write_netcdf_var(ncid,'powc13',powtra2(1,1,1,ipowc13),2*ks,0) CALL write_netcdf_var(ncid,'powc14',powtra2(1,1,1,ipowc14),2*ks,0) #endif diff --git a/hamocc/carchm.F90 b/hamocc/carchm.F90 index f0563983..15263ab4 100644 --- a/hamocc/carchm.F90 +++ b/hamocc/carchm.F90 @@ -637,6 +637,17 @@ SUBROUTINE CARCHM(kpie,kpje,kpke,kbnd, & enddo !$OMP END PARALLEL DO enddo + +!$OMP PARALLEL DO PRIVATE(i) + do j=1,kpje + do i=1,kpie + if(omask(i,j).gt.0.5) then + burial(i,j,issso14) = burial(i,j,issso14)*c14dec + burial(i,j,isssc14) = burial(i,j,isssc14)*c14dec + endif + enddo + enddo +!$OMP END PARALLEL DO #endif #endif diff --git a/hamocc/sedshi.F90 b/hamocc/sedshi.F90 index 44058447..4d032d3b 100644 --- a/hamocc/sedshi.F90 +++ b/hamocc/sedshi.F90 @@ -230,6 +230,16 @@ SUBROUTINE SEDSHI(kpie,kpje,omask) sedlay(i,j,ks,issster)=sedlay(i,j,ks,issster) & & +refill*burial(i,j,issster)/frac +#ifdef cisonew + sedlay(i,j,ks,issso13)=sedlay(i,j,ks,issso13) & + & +refill*burial(i,j,issso13)/frac + sedlay(i,j,ks,isssc13)=sedlay(i,j,ks,isssc13) & + & +refill*burial(i,j,isssc13)/frac + sedlay(i,j,ks,issso14)=sedlay(i,j,ks,issso14) & + & +refill*burial(i,j,issso14)/frac + sedlay(i,j,ks,isssc14)=sedlay(i,j,ks,isssc14) & + & +refill*burial(i,j,isssc14)/frac +#endif ! account for losses in buried sediment burial(i,j,issso12) = burial(i,j,issso12) & & - refill*burial(i,j,issso12) @@ -239,7 +249,16 @@ SUBROUTINE SEDSHI(kpie,kpje,omask) & - refill*burial(i,j,issssil) burial(i,j,issster) = burial(i,j,issster) & & - refill*burial(i,j,issster) - +#ifdef cisonew + burial(i,j,issso13) = burial(i,j,issso13) & + & - refill*burial(i,j,issso13) + burial(i,j,isssc13) = burial(i,j,isssc13) & + & - refill*burial(i,j,isssc13) + burial(i,j,issso14) = burial(i,j,issso14) & + & - refill*burial(i,j,issso14) + burial(i,j,isssc14) = burial(i,j,isssc14) & + & - refill*burial(i,j,isssc14) +#endif endif enddo !end i-loop enddo !end j-loop @@ -262,7 +281,7 @@ SUBROUTINE SEDSHI(kpie,kpje,omask) enddo !end j-loop !$OMP END PARALLEL DO - do iv=1,4 + do iv=1,nsedtra !$OMP PARALLEL DO PRIVATE(i,uebers,frac) do j=1,kpje do i=1,kpie @@ -272,20 +291,6 @@ SUBROUTINE SEDSHI(kpie,kpje,omask) frac=porsol(i,j,k)*seddw(k)/(porsol(i,j,k-1)*seddw(k-1)) sedlay(i,j,k,iv)=sedlay(i,j,k,iv)-uebers sedlay(i,j,k-1,iv)=sedlay(i,j,k-1,iv)+uebers*frac -#ifdef cisonew - if(iv.eq.issso12)then - sedlay(i,j,k,issso13) =sedlay(i,j,k,issso13)-uebers - sedlay(i,j,k-1,issso13)=sedlay(i,j,k-1,issso13)+uebers*frac - sedlay(i,j,k,issso14) =sedlay(i,j,k,issso14)-uebers - sedlay(i,j,k-1,issso14)=sedlay(i,j,k-1,issso14)+uebers*frac - endif - if(iv.eq.isssc12)then - sedlay(i,j,k,isssc13) =sedlay(i,j,k,isssc13)-uebers - sedlay(i,j,k-1,isssc13)=sedlay(i,j,k-1,isssc13)+uebers*frac - sedlay(i,j,k,isssc14) =sedlay(i,j,k,isssc14)-uebers - sedlay(i,j,k-1,isssc14)=sedlay(i,j,k-1,isssc14)+uebers*frac - endif -#endif endif enddo !end i-loop enddo !end j-loop