pkg/dic/phos_flux.F

6d54cf9ca1 Ed H*0001 #include "DIC_OPTIONS.h"
daab022f42 Step*0002
08536d17ba Step*0003 CBOP
                0004 C !ROUTINE: PHOS_FLUX
                0005
                0006 C !INTERFACE: ==========================================================
ac331e4517 Davi*0007       SUBROUTINE PHOS_FLUX( BIOac, pflux, exportflux,
daab022f42 Step*0008      I           bi,bj,imin,imax,jmin,jmax,
2e3e8c330d Jona*0009      I           myTime, myIter, myThid )
daab022f42 Step*0010
08536d17ba Step*0011 C !DESCRIPTION:
                0012 C Calculate the PO4 flux to depth from bio activity
daab022f42 Step*0013
08536d17ba Step*0014 C !USES: ===============================================================
                0015       IMPLICIT NONE
daab022f42 Step*0016 #include "SIZE.h"
                0017 #include "DYNVARS.h"
                0018 #include "EEPARAMS.h"
                0019 #include "PARAMS.h"
                0020 #include "GRID.h"
2ef8966791 Davi*0021 #include "DIC_VARS.h"
daab022f42 Step*0022
08536d17ba Step*0023 C !INPUT PARAMETERS: ===================================================
ac331e4517 Davi*0024 C  BIOac                :: biological productivity
2e3e8c330d Jona*0025 C  myTime               :: current time
                0026 C  myIter               :: current timestep
                0027 C  myThid               :: thread number
ac331e4517 Davi*0028       _RL  BIOac(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
                0029       INTEGER imin, imax, jmin, jmax, bi, bj
2e3e8c330d Jona*0030       _RL myTime
                0031       INTEGER myIter
                0032       INTEGER myThid
08536d17ba Step*0033
                0034 C !OUTPUT PARAMETERS: ===================================================
ac331e4517 Davi*0035 C  pflux                :: changes to PO4 due to flux and reminerlization
22bb60bade Mart*0036       _RL pflux     (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
65132eb687 Step*0037       _RL exportflux(1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
daab022f42 Step*0038
c2a4adfccc Mart*0039 #if (defined ALLOW_PTRACERS && defined DIC_BIOTIC)
08536d17ba Step*0040
                0041 C !LOCAL VARIABLES: ====================================================
                0042 C  i,j,k                  :: loop indices
                0043 c  ko                     :: loop-within-loop index
ac331e4517 Davi*0044 c  bexport                :: flux of phosphorus from base each "productive"
08536d17ba Step*0045 c                            layer
c2a4adfccc Mart*0046 c  depth_l                :: depth and lower interface
ac331e4517 Davi*0047 c  flux_u, flux_l         :: flux through upper and lower interfaces
2c1f2f5f32 Mart*0048 c  reminFac               :: abbreviation
ac331e4517 Davi*0049 c  zbase                  :: depth of bottom of current productive layer
2e3e8c330d Jona*0050       INTEGER i,j,k, ko, kop1
c2a4adfccc Mart*0051       _RL zbase
                0052       _RL bexport(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
2c1f2f5f32 Mart*0053       _RL reminFac
c2a4adfccc Mart*0054       _RL depth_l
                0055       _RL flux_u (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0056       _RL flux_l
08536d17ba Step*0057 CEOP
                0058
ac331e4517 Davi*0059 C- Calculate PO4 flux from base of each layer
c2a4adfccc Mart*0060       DO k=1,nlev
51c3bf0077 Step*0061        DO j=jmin,jmax
                0062         DO i=imin,imax
c2a4adfccc Mart*0063          bexport(i,j) = 0. _d 0
22bb60bade Mart*0064          IF ( _hFacC(i,j,k,bi,bj).gt.0. _d 0 ) THEN
c2a4adfccc Mart*0065 C--   If no layer below initial layer (because of bottom or
                0066 C--   topography), then remineralize in here
                0067           IF (k.EQ.Nr) THEN
2e3e8c330d Jona*0068            pflux(i,j,k)=pflux(i,j,k)+BIOac(i,j,k)*(1. _d 0-DOPfraction)
c2a4adfccc Mart*0069           ELSEIF (hFacC(i,j,k+1,bi,bj).EQ.0. _d 0) THEN
                0070            pflux(i,j,k)=pflux(i,j,k)+BIOac(i,j,k)*(1. _d 0-DOPfraction)
                0071           ELSE
ac331e4517 Davi*0072 C- flux out of layer k
c2a4adfccc Mart*0073            bexport(i,j)=BIOac(i,j,k)*(1. _d 0-DOPfraction)
                0074      &          *drF(k) * _hFacC(i,j,k,bi,bj)
                0075           ENDIF
                0076          ENDIF
                0077         ENDDO
                0078        ENDDO
                0079 C--   If available, flux phosphate downward;
                0080 C--   calculate flux to each layer from base of k
                0081        zbase=-rF(k+1)
                0082 C--   Upper flux
                0083        DO j=jmin,jmax
                0084         DO i=imin,imax
                0085          flux_u(i,j)  = bexport(i,j)
                0086         ENDDO
                0087        ENDDO
22bb60bade Mart*0088 C     Instead of running the loop to ko=Nr and masking the last
                0089 C     flux_l, let ko reach only Nr-1 and do a special loop for ko=Nr,
                0090 C     in order to save a few expensive exp-function calls
                0091        DO ko=k+1,Nr-1
c2a4adfccc Mart*0092         kop1   = MIN(Nr,ko+1)
376d422ed2 Mart*0093 #ifndef NONLIN_FRSURF
2c1f2f5f32 Mart*0094 C     For the linear free surface, hFacC can be omitted, buying another
                0095 C     performance increase of a factor of six on a vector computer.
2e3e8c330d Jona*0096 C     For now this is not implemented via run time flags, in order to
2c1f2f5f32 Mart*0097 C     avoid making this code too complicated.
                0098         depth_l  = -rF(ko) + drF(ko)
3e096c4163 Mart*0099 C       reminFac = (depth_l/zbase)**(-Kremin)
2e3e8c330d Jona*0100 C     The following form does the same, but is faster
3e096c4163 Mart*0101         reminFac = exp(-Kremin*log(depth_l/zbase))
2c1f2f5f32 Mart*0102 #endif
c2a4adfccc Mart*0103         DO j=jmin,jmax
                0104          DO i=imin,imax
                0105           IF ( bexport(i,j) .NE. 0. _d 0 ) THEN
                0106 C--   Lower flux (no flux to ocean bottom)
376d422ed2 Mart*0107 #ifdef NONLIN_FRSURF
2c1f2f5f32 Mart*0108            depth_l  = -rF(ko) + drF(ko) * _hFacC(i,j,ko,bi,bj)
3e096c4163 Mart*0109 C          reminFac = (depth_l/zbase)**(-Kremin)
2e3e8c330d Jona*0110 C     The following form does the same, but is faster
3e096c4163 Mart*0111            reminFac = exp(-Kremin*log(depth_l/zbase))
2c1f2f5f32 Mart*0112 #endif
                0113            flux_l   = bexport(i,j)*reminFac
22bb60bade Mart*0114      &          *maskC(i,j,kop1,bi,bj)
2e3e8c330d Jona*0115 C
c2a4adfccc Mart*0116            pflux(i,j,ko)=pflux(i,j,ko) + (flux_u(i,j)-flux_l)
22bb60bade Mart*0117      &          *recip_drF(ko) * _recip_hFacC(i,j,ko,bi,bj)
c2a4adfccc Mart*0118            exportflux(i,j,ko)=exportflux(i,j,ko)+flux_u(i,j)
                0119 C--   Store flux through upper layer for the next k-level
                0120            flux_u(i,j) = flux_l
                0121 C     endif bexport .ne. 0
ac331e4517 Davi*0122           ENDIF
c2a4adfccc Mart*0123 C     i,j-loops
951926fb9b Jean*0124          ENDDO
daab022f42 Step*0125         ENDDO
c2a4adfccc Mart*0126 C     ko-loop
daab022f42 Step*0127        ENDDO
22bb60bade Mart*0128 C     now do ko = Nr
                0129        ko = Nr
                0130        flux_l = 0. _d 0
                0131        DO j=jmin,jmax
                0132         DO i=imin,imax
                0133          pflux(i,j,ko)=pflux(i,j,ko) + (flux_u(i,j)-flux_l)
                0134      &        *recip_drF(ko) * _recip_hFacC(i,j,ko,bi,bj)
                0135          exportflux(i,j,ko)=exportflux(i,j,ko)+flux_u(i,j)
                0136         ENDDO
                0137        ENDDO
c2a4adfccc Mart*0138 C     k-loop
                0139       ENDDO
daab022f42 Step*0140 c
c2a4adfccc Mart*0141 #endif /* defined ALLOW_PTRACERS && defined DIC_BIOTIC */
                0142       RETURN
                0143       END