pkg/dic/car_flux_omega_top.F

8e447d7023 Step*0001 #include "DIC_OPTIONS.h"
                0002
                0003 CBOP
2e3e8c330d Jona*0004 C !ROUTINE: CAR_FLUX_OMEGA_TOP
8e447d7023 Step*0005
                0006 C !INTERFACE: ==========================================================
                0007       SUBROUTINE CAR_FLUX_OMEGA_TOP( bioac, cflux,
                0008      I           bi,bj,imin,imax,jmin,jmax,
2e3e8c330d Jona*0009      I           myTime, myIter, myThid )
8e447d7023 Step*0010
                0011 C !DESCRIPTION:
                0012 C  Calculate carbonate fluxes
                0013 C  HERE ONLY HAVE DISSOLUTION WHEN OMEGA < 1.0
                0014 C  Karsten Friis and Mick Follows Sep 2004
                0015
                0016 C !USES: ===============================================================
                0017       IMPLICIT NONE
                0018 #include "SIZE.h"
                0019 #include "DYNVARS.h"
                0020 #include "EEPARAMS.h"
                0021 #include "PARAMS.h"
                0022 #include "GRID.h"
2ef8966791 Davi*0023 #include "DIC_VARS.h"
8e447d7023 Step*0024
                0025 C !INPUT PARAMETERS: ===================================================
                0026 C  bioac                :: biological productivity
2e3e8c330d Jona*0027 C  myTime               :: current time
                0028 C  myIter               :: current timestep
                0029 C  myThid               :: thread number
                0030       _RL  bioac(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
                0031       INTEGER imin, imax, jmin, jmax, bi, bj
8e447d7023 Step*0032       _RL myTime
2e3e8c330d Jona*0033       INTEGER myIter
8e447d7023 Step*0034       INTEGER myThid
                0035
                0036 C !OUTPUT PARAMETERS: ===================================================
                0037 C cflux                :: carbonate flux
                0038       _RL  cflux(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
                0039
2e3e8c330d Jona*0040 #if ( defined DIC_BIOTIC && defined DIC_CALCITE_SAT )
8e447d7023 Step*0041 C !LOCAL VARIABLES: ====================================================
                0042 C  i,j,k                  :: loop indices
2e3e8c330d Jona*0043 C  ko                     :: loop-within-loop index
                0044 C caexport                :: flux of carbonate from base each "productive"
                0045 C                            layer
                0046 C depth_u, depth_l        :: depths of upper and lower interfaces
                0047 C flux_u, flux_l          :: flux through upper and lower interfaces
8e447d7023 Step*0048        _RL caexport(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
2e3e8c330d Jona*0049        INTEGER i,j,k, ko, kBottom
8e447d7023 Step*0050        _RL flux_u, flux_l
2e3e8c330d Jona*0051 C variables for calcium carbonate dissolution
8e447d7023 Step*0052        _RL DissolutionRate
                0053        _RL dumrate
                0054
2e3e8c330d Jona*0055 C diagnostics
8e447d7023 Step*0056 c     _RL   exp_tot
                0057 c     _RL   flx_tot
                0058 c     integer knum
                0059 c     _RL   omeg_bot
                0060 c     _RL   tmp
                0061
                0062 CEOP
2e3e8c330d Jona*0063 C initialize fluxes
                0064       flux_u = 0. _d 0
                0065       flux_l = 0. _d 0
                0066
                0067 C to either remineralize in bottom or top layer if flux reaches bottom layer
                0068 C      selectCalciteBottomRemin =0 : in bottom layer, =1 : in top layer
                0069 C set some nominal particulate sinking rate, try 100m/day:
                0070 C       WsinkPIC = 100. _d 0/86400. _d 0
                0071 C calculate carbonate flux from base of each nlev
51c3bf0077 Step*0072        DO j=jmin,jmax
                0073         DO i=imin,imax
8e447d7023 Step*0074 c        exp_tot=0
2e3e8c330d Jona*0075          do k=1,Nr
c5830e1a8b Jona*0076             cflux(i,j,k)=0. _d 0
8e447d7023 Step*0077          enddo
2e3e8c330d Jona*0078
c2435f7f83 Jona*0079          kBottom   = kLowC(i,j,bi,bj)
2e3e8c330d Jona*0080
8e447d7023 Step*0081          DO k=1,nLev
c5830e1a8b Jona*0082           if (hFacC(i,j,k,bi,bj).gt.0. _d 0) then
75e97f1e14 Davi*0083            caexport(i,j)= R_CP*rain_ratio(i,j,bi,bj)*bioac(i,j,k)*
c5830e1a8b Jona*0084      &           (1. _d 0-DOPfraction)*drF(k)*hFacC(i,j,k,bi,bj)
8e447d7023 Step*0085 c          exp_tot=exp_tot+caexport(i,j)
2e3e8c330d Jona*0086 C calculate flux to each layer from base of k
a2c35952f2 Jona*0087            DO ko=k+1,Nr
c5830e1a8b Jona*0088             if (hFacC(i,j,ko,bi,bj).gt.0. _d 0) then
8e447d7023 Step*0089               if (ko .eq. k+1) then
                0090                 flux_u = caexport(i,j)
                0091               else
                0092                 flux_u = flux_l
                0093               endif
                0094
                0095 C flux through lower face of cell
c5830e1a8b Jona*0096               if (omegaC(i,j,ko,bi,bj) .gt. 1. _d 0) then
8e447d7023 Step*0097                 flux_l = flux_u
                0098
2e3e8c330d Jona*0099 C if at bottom, remineralize remaining flux
c2435f7f83 Jona*0100                 if (ko.eq.kBottom) then
2e3e8c330d Jona*0101                   if (selectCalciteBottomRemin.EQ.1) then
                0102 C ... at surface
8e447d7023 Step*0103                      cflux(i,j,1)=cflux(i,j,1)+
                0104      &                  ( (flux_l)/(drF(1)*hFacC(i,j,1,bi,bj)) )
                0105                   else
                0106
2e3e8c330d Jona*0107 C ... at bottom
c5830e1a8b Jona*0108                      flux_l=0. _d 0
8e447d7023 Step*0109                   endif
                0110                 endif
                0111               else
a2c35952f2 Jona*0112 C if dissolution
                0113                IF (selectCalciteDissolution .eq. 0) THEN
                0114 C Use constant dissolution rate
                0115                 flux_l = flux_u
                0116      &                 *(1.0-calciteDissolRate(1)*drF(k)/WsinkPIC)
                0117                ELSEIF (selectCalciteDissolution .eq. 1) THEN
                0118 C Use Micks version of dissolution rate based on vertical sinking/remin balance
                0119                 DissolutionRate = calciteDissolRate(1)*(
                0120      &             (
                0121      &              1. _d 0-omegaC(i,j,ko,bi,bj)
                0122      &             )**calciteDissolExp(1)
                0123      &             )/(86400. _d 0)
                0124
                0125                 dumrate = -1. _d 0*DissolutionRate*drF(ko)*
55d8dbd247 Step*0126      &                       hFacC(i,j,ko,bi,bj)/WsinkPIC
a2c35952f2 Jona*0127                 flux_l = flux_u*exp(dumrate)
                0128                ELSEIF (selectCalciteDissolution .eq. 2) THEN
                0129 C Use Karstens version of dissolution rate based on from Keir (1980) Geochem.
                0130 C  Cosmochem. Acta, and bugfix H/T Oliver Jahn (Keir,s dissolution rate in log(%))
                0131                 DissolutionRate = EXP(calciteDissolRate(1))*(
                0132      &             (
                0133      &              1. _d 0-omegaC(i,j,ko,bi,bj)
                0134      &             )**calciteDissolExp(1)
                0135      &             )/(100. _d 0 * 86400. _d 0)
                0136
                0137                 flux_l = flux_u*(1.0-DissolutionRate*drF(k)/WsinkPIC)
                0138                ELSEIF (selectCalciteDissolution .eq. 3) THEN
                0139 C Use Naviaux et al. 2019, Marine Chemistry dissolution rates
                0140 C  The value of 0.8272 differs slightly from the paper (0.8)), but reduces the
                0141 C  large discontinuity in dissolution rates between the two states.
                0142                 IF (omegaC(i,j,ko,bi,bj) .GT. 0.8272 _d 0) THEN
                0143                  DissolutionRate = calciteDissolRate(1) * (
                0144      &             (
                0145      &              1. _d 0-omegaC(i,j,ko,bi,bj)
                0146      &             )**calciteDissolExp(1))
                0147                 ELSE
                0148                  DissolutionRate = calciteDissolRate(2) * (
                0149      &             (
                0150      &              1. _d 0-omegaC(i,j,ko,bi,bj)
                0151      &             )**calciteDissolExp(2))
                0152                 ENDIF
                0153                 flux_l = flux_u*(1.0-DissolutionRate*drF(k)/WsinkPIC)
                0154                ENDIF
                0155
2e3e8c330d Jona*0156 C TEST ............................
8e447d7023 Step*0157 c           if(i .eq. 76 .and. j .eq. 36)then
                0158 c            write(6,*)'k,flux_l/flux_u',ko,(flux_l/flux_u)
a2c35952f2 Jona*0159 c            write(6,*)'K, DissolutionRate, drF(k), drF(ko), WsinkPIC,OmegaC'
                0160 c            write(6,*)ko,DissolutionRate,drF(k),drF(ko),WsinkPIC,
8e447d7023 Step*0161 c    &            omegaC(i,j,ko,bi,bj)
                0162 c           endif
2e3e8c330d Jona*0163 C TEST ............................
                0164 C no flux to ocean bottom
c2435f7f83 Jona*0165                  if (ko.eq.kBottom)
c5830e1a8b Jona*0166      &                      flux_l=0. _d 0
8e447d7023 Step*0167               endif
                0168
2e3e8c330d Jona*0169 C flux divergence
a2c35952f2 Jona*0170               cflux(i,j,ko)=cflux(i,j,ko) +
8e447d7023 Step*0171      &          ( (flux_u-flux_l)/(drF(ko)*hFacC(i,j,ko,bi,bj)) )
2e3e8c330d Jona*0172 C TEST ............................
8e447d7023 Step*0173 c            if(i .eq. 76 .and. j .eq. 36)then
                0174 c               write(6,*)'k,flux_l/flux_u',ko,(flux_l/flux_u)
                0175 c              write(6,*)'k,flux_l,cflux ',ko,flux_l,cflux(i,j,ko)
                0176 c            endif
2e3e8c330d Jona*0177 C TEST ............................
a2c35952f2 Jona*0178             else
2e3e8c330d Jona*0179 C if no layer below initial layer, remineralize
8e447d7023 Step*0180                if (ko.eq.k+1) then
2e3e8c330d Jona*0181                 if ( selectCalciteBottomRemin.EQ.1 .AND.
                0182      &                  omegaC(i,j,k,bi,bj) .GT. 1. _d 0 ) then
                0183 C ... at surface
8e447d7023 Step*0184                    cflux(i,j,1)=cflux(i,j,1)
c5830e1a8b Jona*0185      &                  +bioac(i,j,k)*(1. _d 0-DOPfraction)*
75e97f1e14 Davi*0186      &                    R_CP*rain_ratio(i,j,bi,bj)
8e447d7023 Step*0187      &                   *drF(k)*hFacC(i,j,k,bi,bj)/
                0188      &                    (drF(1)*hFacC(i,j,1,bi,bj) )
                0189                 else
2e3e8c330d Jona*0190 C ... at bottom
8e447d7023 Step*0191                   cflux(i,j,k)=cflux(i,j,k)
c5830e1a8b Jona*0192      &                  +bioac(i,j,k)*(1. _d 0-DOPfraction)*
75e97f1e14 Davi*0193      &                    R_CP*rain_ratio(i,j,bi,bj)
8e447d7023 Step*0194                 endif
                0195                endif
a2c35952f2 Jona*0196             endif
                0197            ENDDO
8e447d7023 Step*0198
                0199           endif
                0200          ENDDO
2e3e8c330d Jona*0201 C diagnostic
8e447d7023 Step*0202 c        flx_tot=0
                0203 c        k=0
                0204 c        do k=1,nR
                0205 c          flx_tot=flx_tot+cflux(i,j,k)*drF(k)*hFacC(i,j,k,bi,bj)
c5830e1a8b Jona*0206 c          if (hFacC(i,j,k,bi,bj).gt.0 _d 0) then
8e447d7023 Step*0207 c             knum=k
                0208 c             omeg_bot=omegaC(i,j,k,bi,bj)
                0209 c          endif
                0210 c        enddo
c5830e1a8b Jona*0211 c        if (hFacC(i,j,k,bi,bj).gt. 0. _d 0) then
8e447d7023 Step*0212 c         tmp=abs(exp_tot-flx_tot)
c5830e1a8b Jona*0213 c         if (tmp>1 _d -20) then
8e447d7023 Step*0214 c          print*,'QQ car_flux', knum,
                0215 c    &                 omeg_bot, exp_tot, flx_tot, exp_tot-flx_tot
                0216 c         endif
                0217 c        endif
2e3e8c330d Jona*0218 C end diagnostic
8e447d7023 Step*0219         ENDDO
                0220        ENDDO
                0221 c
2e3e8c330d Jona*0222 #endif /* DIC_BIOTIC and DIC_CALCITE_SAT */
8e447d7023 Step*0223        RETURN
                0224        END