pkg/bling/bling_light.F

c0d1c06c15 Matt*0001 #include "BLING_OPTIONS.h"
a284455135 Matt*0002 #ifdef ALLOW_AUTODIFF
                0003 # include "AUTODIFF_OPTIONS.h"
                0004 #endif
c0d1c06c15 Matt*0005
                0006 CBOP
4dea327916 aver*0007 C !ROUTINE: BLING_LIGHT
                0008
                0009 C !INTERFACE: ==========================================================
a284455135 Matt*0010       SUBROUTINE BLING_LIGHT(
e0f9a7ba0b Matt*0011      I           mld,
                0012      U           irr_inst, irr_eff,
2a2b7d0c36 Mart*0013      I           bi, bj, iMin, iMax, jMin, jMax,
e0f9a7ba0b Matt*0014      I           myTime, myIter, myThid)
                0015
4dea327916 aver*0016 C !DESCRIPTION:
                0017 C     o calculate effective light for phytoplankton growth
                0018 C       There are multiple types of light.
                0019 C     - irr_inst is the instantaneous irradiance field.
                0020 C     - irr_mix is the same, but with the irr_inst averaged throughout
                0021 C       the mixed layer. This quantity is intended to represent the
                0022 C       light to which phytoplankton subject to turbulent transport in
                0023 C       the mixed-layer would be exposed.
                0024 C     - irr_eff is the effective irradiance for photosynthesis,
                0025 C       given either by irr_inst or irr_mix, depending on model
                0026 C       options and location.
                0027 C     o instantaneous light is calculated either from
                0028 C     - date and latitude, then exponentially attenuated down the
                0029 C       water column, or
                0030 C     - short-wave radiation read from external forcing file,
                0031 C       attenuated down the water column, or
                0032 C     - short-wave radiation distributed through the water column
                0033 C       according to SWFRAC routine
                0034
                0035 C !USES: ===============================================================
e0f9a7ba0b Matt*0036       IMPLICIT NONE
                0037
c0d1c06c15 Matt*0038 C     === Global variables ===
                0039 #include "SIZE.h"
                0040 #include "EEPARAMS.h"
                0041 #include "PARAMS.h"
                0042 #include "GRID.h"
                0043 #include "BLING_VARS.h"
e0f9a7ba0b Matt*0044 #ifdef USE_QSW
2a2b7d0c36 Mart*0045 # include "FFIELDS.h"
e0f9a7ba0b Matt*0046 #endif
a284455135 Matt*0047 #ifdef ALLOW_AUTODIFF_TAMC
c0d1c06c15 Matt*0048 # include "tamc.h"
                0049 #endif
                0050
4dea327916 aver*0051 C !INPUT PARAMETERS: ===================================================
c0d1c06c15 Matt*0052 C     bi,bj         :: tile indices
                0053 C     iMin,iMax     :: computation domain: 1rst index range
                0054 C     jMin,jMax     :: computation domain: 2nd  index range
                0055 C     myTime        :: current time
                0056 C     myIter        :: current timestep
                0057 C     myThid        :: thread Id. number
2a2b7d0c36 Mart*0058       INTEGER bi, bj, iMin, iMax, jMin, jMax
c0d1c06c15 Matt*0059       INTEGER myThid
                0060       INTEGER myIter
                0061       _RL     myTime
                0062       _RL mld       (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
4dea327916 aver*0063
                0064 C !OUTPUT PARAMETERS: ==================================================
e37161e05a Jean*0065 C     irr_inst      :: instantaneous light
                0066 C     irr_eff       :: effective light for photosynthesis
c0d1c06c15 Matt*0067       _RL irr_inst  (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
                0068       _RL irr_eff   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
                0069
4dea327916 aver*0070 C !LOCAL VARIABLES: ====================================================
e37161e05a Jean*0071       INTEGER i,j,k
                0072       LOGICAL QSW_underice
                0073 #ifdef ALLOW_CAL
                0074       INTEGER mydate(4)
                0075 #endif
                0076       _RL localTime
                0077       _RL utcTime, diffutc
                0078       _RL sat_atten
                0079       _RL sat_atten_sum(1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0080       _RL chl_sat_sum  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
c0d1c06c15 Matt*0081       _RL atten
                0082       _RL irr_surf  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
e0f9a7ba0b Matt*0083 #ifdef ML_MEAN_LIGHT
c0d1c06c15 Matt*0084       _RL irr_mix   (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
4dea327916 aver*0085       _RL SumMLIrr  (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
                0086       _RL tmp_ML    (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
c0d1c06c15 Matt*0087 #endif
                0088 #ifndef USE_QSW
2a2b7d0c36 Mart*0089       _RL sfac      (1-OLx:sNx+OLx,1-OLy:sNy+OLy)
c0d1c06c15 Matt*0090 #endif
00fa2d4ddd mmaz*0091 #ifdef PHYTO_SELF_SHADING
                0092       _RL k0_rd, chi_rd, e_rd
                0093       _RL k0_bg, chi_bg, e_bg
a284455135 Matt*0094       _RL kChl_rd   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
                0095       _RL kChl_bg   (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
00fa2d4ddd mmaz*0096       _RL atten_rd
                0097       _RL atten_bg
a284455135 Matt*0098       _RL irr_rd    (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
                0099       _RL irr_bg    (1-OLx:sNx+OLx,1-OLy:sNy+OLy,Nr)
4dea327916 aver*0100 #endif /* PHYTO_SELF_SHADING */
4e4ad91a39 Jean*0101 #ifdef ALLOW_AUTODIFF_TAMC
edb6656069 Mart*0102 C     tkey :: tape key (tile dependent)
4dea327916 aver*0103 C     kkey :: tape key (tile and level dependent)
                0104       INTEGER tkey, kkey
00fa2d4ddd mmaz*0105 #endif
c0d1c06c15 Matt*0106 CEOP
                0107
e9828181c3 Yixi*0108 C  Remove light under ice
                0109 C  If using Qsw and seaice/thsice, then ice fraction is already
                0110 C  taken into account
4dea327916 aver*0111       QSW_underice = .FALSE.
                0112 #ifdef USE_QSW
                0113       IF ( useSEAICE ) QSW_underice = .TRUE.
                0114       IF ( useThSIce ) QSW_underice = .TRUE.
                0115 #endif
a284455135 Matt*0116
e37161e05a Jean*0117       DO j=1-OLy,sNy+OLy
                0118        DO i=1-OLx,sNx+OLx
                0119         chl_sat_sum(i,j)   = 0. _d 0
                0120         sat_atten_sum(i,j) = 0. _d 0
                0121 #ifdef ML_MEAN_LIGHT
                0122         SumMLIrr(i,j) = 0. _d 0
                0123         tmp_ML(i,j)   = 0. _d 0
                0124 #endif
2a2b7d0c36 Mart*0125 #ifndef USE_QSW
                0126         sfac(i,j)     = 0. _d 0
                0127 #endif
e37161e05a Jean*0128        ENDDO
                0129       ENDDO
4dea327916 aver*0130       DO k=1,Nr
e37161e05a Jean*0131        DO j=1-OLy,sNy+OLy
                0132         DO i=1-OLx,sNx+OLx
4dea327916 aver*0133          irr_eff(i,j,k) = 0. _d 0
a284455135 Matt*0134 #ifdef PHYTO_SELF_SHADING
4dea327916 aver*0135          irr_rd(i,j,k)        = 0. _d 0
                0136          irr_bg(i,j,k)        = 0. _d 0
a284455135 Matt*0137 #endif
c0d1c06c15 Matt*0138         ENDDO
                0139        ENDDO
4dea327916 aver*0140       ENDDO
e37161e05a Jean*0141
4dea327916 aver*0142 #ifdef PHYTO_SELF_SHADING
e9828181c3 Yixi*0143 C  Specify co-efficients for bio-optical model {kChl = k0 +chi[chl]^e}
                0144 C  in red and blue-green fractions (Morel 1988; Foujols et al. 2000)
4dea327916 aver*0145       k0_rd  = 0.225 _d 0
                0146       k0_bg  = 0.0232 _d 0
                0147       chi_rd = 0.037 _d 0
                0148       chi_bg = 0.074 _d 0
                0149       e_rd   = 0.629 _d 0
                0150       e_bg   = 0.674 _d 0
                0151 #endif
c0d1c06c15 Matt*0152
e9828181c3 Yixi*0153 C ---------------------------------------------------------------------
                0154 C  instantaneous light, mixed layer averaged light
c0d1c06c15 Matt*0155
e37161e05a Jean*0156 #ifdef USE_QSW
2a2b7d0c36 Mart*0157 C  Convert insolation to photosynthetically-available radiations (PAR)
                0158       DO j=jMin,jMax
                0159        DO i=iMin,iMax
                0160         irr_surf(i,j) = MAX( epsln,
                0161      &                  -parfrac*Qsw(i,j,bi,bj)*maskInC(i,j,bi,bj))
                0162        ENDDO
                0163       ENDDO
e37161e05a Jean*0164 #else
2a2b7d0c36 Mart*0165 C  Find surface insolation (light) as function of date and latitude
                0166       CALL GCHEM_INSOLATION(
                0167      O     sfac,
                0168      I     iMin, iMax, jMin, jMax, bi, bj,
                0169      I     myTime, myIter, myThid)
                0170
                0171 C Photosynthetically active solar radiation (PAR) just below surface
                0172       DO j=jMin,jMax
                0173        DO i=iMin,iMax
                0174         irr_surf(i,j) = MAX(1. _d -5,sfac(i,j)*parfrac)
e37161e05a Jean*0175        ENDDO
                0176       ENDDO
2a2b7d0c36 Mart*0177 #endif /* USE_QSW */
                0178
                0179 C  Remove light under ice
                0180       IF ( .NOT. QSW_underice ) THEN
                0181        DO j=jMin,jMax
                0182         DO i=iMin,iMax
                0183          irr_surf(i,j) = irr_surf(i,j)*(1. _d 0 - FIce(i,j,bi,bj))
                0184         ENDDO
                0185        ENDDO
                0186       ENDIF
e37161e05a Jean*0187
4dea327916 aver*0188 #ifdef ALLOW_AUTODIFF_TAMC
                0189       tkey = bi + (bj - 1)*nSx + (ikey_dynamics - 1)*nSx*nSy
2a2b7d0c36 Mart*0190 CADJ STORE irr_surf = comlev1_bibj, key=tkey, kind=isbyte
4dea327916 aver*0191 #endif /* ALLOW_AUTODIFF_TAMC */
                0192
                0193       DO k=1,Nr
                0194
                0195 #ifdef ALLOW_AUTODIFF_TAMC
                0196        kkey = k + (tkey-1)*Nr
                0197 # ifdef ML_MEAN_LIGHT
                0198 CADJ STORE tmp_ML = comlev1_bibj_k, key=kkey, kind=isbyte
                0199 # endif /* ML_MEAN_LIGHT */
                0200 #endif
                0201
                0202 C Top layer
                0203        IF ( k.EQ.1) THEN
                0204
2a2b7d0c36 Mart*0205         DO j=jMin,jMax
                0206          DO i=iMin,iMax
4dea327916 aver*0207
                0208           IF ( maskC(i,j,k,bi,bj).EQ.oneRS ) THEN
e0f9a7ba0b Matt*0209
00fa2d4ddd mmaz*0210 #ifdef PHYTO_SELF_SHADING
e9828181c3 Yixi*0211 C  Use bio-optical model of Manizza et al. (2005) to account for
                0212 C  effect of self-shading on ligt available for phytoplankton
                0213 C  growth. As written this DOES NOT feedback onto the absorption
                0214 C  of shortwave radiation calculated in the physical model, which
                0215 C  is instead calculated in the subroutine swfrac
00fa2d4ddd mmaz*0216
e9828181c3 Yixi*0217 C  Attenuation coefficient adjusted to chlorophyll in top layer
4dea327916 aver*0218 #ifdef ALLOW_AUTODIFF
                0219            IF ( chl(i,j,1,bi,bj) .GT. 0. _d 0 ) THEN
                0220 #endif
                0221             kChl_rd(i,j,1) = k0_rd + chi_rd*(chl(i,j,1,bi,bj)**e_rd)
                0222             kChl_bg(i,j,1) = k0_bg + chi_bg*(chl(i,j,1,bi,bj)**e_bg)
                0223 #ifdef ALLOW_AUTODIFF
                0224            ELSE
                0225             kChl_rd(i,j,1) = k0_rd
                0226             kChl_bg(i,j,1) = k0_bg
                0227            ENDIF
                0228 #endif
e9828181c3 Yixi*0229 C  Light attenuation in middle of top layer
4dea327916 aver*0230            atten_rd = kChl_rd(i,j,1)*drF(1)/2. _d 0*hFacC(i,j,1,bi,bj)
                0231            atten_bg = kChl_bg(i,j,1)*drF(1)/2. _d 0*hFacC(i,j,1,bi,bj)
e9828181c3 Yixi*0232 C  Irradiance in middle of top layer
4dea327916 aver*0233            irr_rd(i,j,1) = irr_surf(i,j) * exp(-atten_rd) * 0.5 _d 0
                0234            irr_bg(i,j,1) = irr_surf(i,j) * exp(-atten_bg) * 0.5 _d 0
                0235            irr_inst(i,j,1) = irr_rd(i,j,1) + irr_bg(i,j,1)
e37161e05a Jean*0236 #else /* PHYTO_SELF_SHADING */
4dea327916 aver*0237 C SW radiation attenuated exponentially
e9828181c3 Yixi*0238 C  Light attenuation in middle of top layer
                0239            atten = k0_2d(i,j,bi,bj)*halfRL*drF(1)*hFacC(i,j,1,bi,bj)
4dea327916 aver*0240            irr_inst(i,j,1) = irr_surf(i,j)*exp(-atten)
                0241
e37161e05a Jean*0242 #endif /* PHYTO_SELF_SHADING */
4dea327916 aver*0243
                0244           ENDIF
                0245          ENDDO
                0246         ENDDO
e37161e05a Jean*0247
4dea327916 aver*0248 C k>1: below surface layer
                0249        ELSE
                0250
                0251 #ifdef ALLOW_AUTODIFF_TAMC
                0252 # ifdef PHYTO_SELF_SHADING
                0253 CADJ STORE irr_bg(:,:,k-1) = comlev1_bibj_k, key=kkey, kind=isbyte
                0254 CADJ STORE irr_rd(:,:,k-1) = comlev1_bibj_k, key=kkey, kind=isbyte
                0255 # endif
                0256 #endif
                0257
2a2b7d0c36 Mart*0258         DO j=jMin,jMax
                0259          DO i=iMin,iMax
4dea327916 aver*0260
                0261           IF ( maskC(i,j,k,bi,bj).EQ.oneRS ) THEN
00fa2d4ddd mmaz*0262
4dea327916 aver*0263 #ifdef PHYTO_SELF_SHADING
e9828181c3 Yixi*0264 C  Attenuation coefficient adjusted to chlorophyll in kth layer
4dea327916 aver*0265 #ifdef ALLOW_AUTODIFF
                0266            IF ( chl(i,j,k,bi,bj) .GT. 0. _d 0 ) THEN
                0267 #endif
                0268             kChl_rd(i,j,k) = k0_rd + chi_rd*(chl(i,j,k,bi,bj)**e_rd)
                0269             kChl_bg(i,j,k) = k0_bg + chi_bg*(chl(i,j,k,bi,bj)**e_bg)
                0270 #ifdef ALLOW_AUTODIFF
                0271            ELSE
                0272             kChl_rd(i,j,k) = k0_rd
                0273             kChl_bg(i,j,k) = k0_bg
                0274            ENDIF
                0275 #endif
e9828181c3 Yixi*0276 C  Light attenuation from one more layer
4dea327916 aver*0277            atten_rd = kChl_rd(i,j,k)*drF(k)/2. _d 0*hFacC(i,j,k,bi,bj)
                0278      &        + kChl_rd(i,j,k-1)*drF(k-1)/2. _d 0*hFacC(i,j,k-1,bi,bj)
                0279            atten_bg = kChl_bg(i,j,k)*drF(k)/2. _d 0*hFacC(i,j,k,bi,bj)
                0280      &        + kChl_bg(i,j,k-1)*drF(k-1)/2. _d 0*hFacC(i,j,k-1,bi,bj)
e9828181c3 Yixi*0281 C  Irradiance in middle of layer k
4dea327916 aver*0282            irr_rd(i,j,k) = irr_rd(i,j,k-1)*exp(-atten_rd)
                0283            irr_bg(i,j,k) = irr_bg(i,j,k-1)*exp(-atten_bg)
                0284            irr_inst(i,j,k) = irr_rd(i,j,k) + irr_bg(i,j,k)
00fa2d4ddd mmaz*0285
e37161e05a Jean*0286 #else /* PHYTO_SELF_SHADING */
e0f9a7ba0b Matt*0287 C SW radiation attenuated exponentially
e9828181c3 Yixi*0288 C  Attenuation from one more layer
                0289            atten = k0_2d(i,j,bi,bj)*halfRL*( drF(k)*hFacC(i,j,k,bi,bj)
                0290      &                                 + drF(k-1)*hFacC(i,j,k-1,bi,bj) )
e37161e05a Jean*0291            irr_inst(i,j,k) = irr_inst(i,j,k-1)*exp(-atten)
4dea327916 aver*0292
e37161e05a Jean*0293 #endif /* PHYTO_SELF_SHADING */
4dea327916 aver*0294
                0295           ENDIF
                0296          ENDDO
                0297         ENDDO
c0d1c06c15 Matt*0298
e37161e05a Jean*0299        ENDIF /* if k=1 then, else */
e0f9a7ba0b Matt*0300
82e538d851 aver*0301 C Satellite chl
2a2b7d0c36 Mart*0302        DO j=jMin,jMax
                0303         DO i=iMin,iMax
4dea327916 aver*0304          IF ( maskC(i,j,k,bi,bj).EQ.oneRS ) THEN
                0305
e37161e05a Jean*0306           IF ( irr_surf(i,j).GT.zeroRL ) THEN
82e538d851 aver*0307 c           sat_atten = irr_inst(i,j,k)/irr_surf(i,j)
4dea327916 aver*0308 #ifdef PHYTO_SELF_SHADING
e37161e05a Jean*0309            sat_atten = exp(-2. _d 0 * k0_bg * (-rC(k)))
82e538d851 aver*0310 #else
e9828181c3 Yixi*0311            sat_atten = exp(-2. _d 0 * k0_2d(i,j,bi,bj) * (-rC(k)))
82e538d851 aver*0312 #endif
4dea327916 aver*0313            chl_sat_sum(i,j) = chl_sat_sum(i,j)
82e538d851 aver*0314      &                 + chl(i,j,k,bi,bj)*sat_atten
4dea327916 aver*0315            sat_atten_sum(i,j) = sat_atten_sum(i,j) + sat_atten
                0316           ENDIF
82e538d851 aver*0317
c0d1c06c15 Matt*0318 #ifdef ML_MEAN_LIGHT
e9828181c3 Yixi*0319 C  Mean irradiance in the mixed layer
e37161e05a Jean*0320           IF ( (-rF(k+1).LE. mld(i,j)) .AND.
                0321      &         (-rF(k+1).LT.MLmix_max) ) THEN
4dea327916 aver*0322            SumMLIrr(i,j) = SumMLIrr(i,j)+drF(k)*irr_inst(i,j,k)
                0323            tmp_ML(i,j) = tmp_ML(i,j) + drF(k)
                0324            irr_mix(i,j) = SumMLIrr(i,j)/tmp_ML(i,j)
                0325           ENDIF
c0d1c06c15 Matt*0326 #endif
                0327
e37161e05a Jean*0328          ENDIF
                0329         ENDDO
                0330        ENDDO
                0331
                0332 C     end first k loop
                0333       ENDDO
                0334
                0335 C Satellite chlorophyll
2a2b7d0c36 Mart*0336
                0337 C get time (in h) within the day:
                0338       utcTime = MOD( myTime/3600. _d 0, 24. _d 0 )
                0339 #ifdef ALLOW_CAL
                0340 C mydate is utc time
                0341       IF ( useCAL ) THEN
                0342        CALL CAL_GETDATE( myIter, myTime, mydate, myThid )
                0343        i = mydate(2)/10000
                0344        j = mydate(2)/100
                0345        j = MOD(j,100)
                0346        k = MOD(mydate(2),100)
                0347        utcTime = i + j/60. _d 0 + k/3600. _d 0
                0348       ENDIF
                0349 #endif
82e538d851 aver*0350 C Update diagnostic only if ~13:30 local time, when satellite observes
2a2b7d0c36 Mart*0351       DO j=jMin,jMax
                0352        DO i=iMin,iMax
e37161e05a Jean*0353         IF ( usingSphericalPolarGrid .OR. usingCurvilinearGrid ) THEN
                0354 C       local-time difference (in h) from UTC time (note: 15 = 360/24)
                0355          diffutc = XC(i,j,bi,bj)/15. _d 0
                0356         ELSE
                0357 C       for other grid (e.g., cartesian), assumes no difference in time
                0358          diffutc = 0. _d 0
                0359         ENDIF
                0360         localTime = utcTime + diffutc + 24. _d 0
                0361         localTime = MOD( localTime, 24. _d 0 )
                0362         IF ( localTime.GT.chlsat_locTimWindow(1) .AND.
                0363      &       localTime.LT.chlsat_locTimWindow(2) ) THEN
                0364          chl_sat(i,j,bi,bj) = chl_sat_sum(i,j)
                0365      &                      / (sat_atten_sum(i,j) + epsln)
                0366         ENDIF
                0367        ENDDO
                0368       ENDDO
c0d1c06c15 Matt*0369
e37161e05a Jean*0370       DO k=1,Nr
2a2b7d0c36 Mart*0371        DO j=jMin,jMax
                0372         DO i=iMin,iMax
e37161e05a Jean*0373          IF ( maskC(i,j,k,bi,bj).EQ.oneRS ) THEN
e0f9a7ba0b Matt*0374
e37161e05a Jean*0375           irr_eff(i,j,k) = irr_inst(i,j,k)
c0d1c06c15 Matt*0376 #ifdef ML_MEAN_LIGHT
e9828181c3 Yixi*0377 C  Inside mixed layer, effective light is set to mean mixed layer light
e37161e05a Jean*0378           IF ( (-rF(k+1).LE. mld(i,j)) .AND.
                0379      &         (-rF(k+1).LT.MLmix_max) ) THEN
c0d1c06c15 Matt*0380            irr_eff(i,j,k) = irr_mix(i,j)
                0381           ENDIF
e0f9a7ba0b Matt*0382 #endif
c0d1c06c15 Matt*0383
                0384          ENDIF
                0385         ENDDO
                0386        ENDDO
                0387       ENDDO
e0f9a7ba0b Matt*0388
82e538d851 aver*0389 #ifdef ALLOW_DIAGNOSTICS
                0390       IF ( useDiagnostics ) THEN
4dea327916 aver*0391        CALL DIAGNOSTICS_FILL(chl_sat,'BLGCHLSA',0,1,1,bi,bj,myThid)
82e538d851 aver*0392       ENDIF
                0393 #endif /* ALLOW_DIAGNOSTICS */
                0394
c0d1c06c15 Matt*0395       RETURN
                0396       END