New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
p4zopt.F90 in branches/nemo_v3_3_beta/NEMOGCM/NEMO/TOP_SRC/PISCES – NEMO

source: branches/nemo_v3_3_beta/NEMOGCM/NEMO/TOP_SRC/PISCES/p4zopt.F90 @ 2317

Last change on this file since 2317 was 2317, checked in by cetlod, 14 years ago

Improve the computation of the divergence of the downward solar irradiance in traqsr.F90, see ticket #726

  • change fsdepw & fse3t into fsdepw_0 & fse3t_0 in tra_qsr_init
  • modify tra_qsr to recompute systematically etot3 in vvl case
  • suppress the namelist parameter rn_si2 which is computed as the RGB longest depth of extinction in trc_oce.F90
  • Property svn:keywords set to Id
File size: 10.4 KB
Line 
1MODULE p4zopt
2   !!======================================================================
3   !!                         ***  MODULE p4zopt  ***
4   !! TOP - PISCES : Compute the light availability in the water column
5   !!======================================================================
6   !! History :   1.0  !  2004     (O. Aumont) Original code
7   !!             2.0  !  2007-12  (C. Ethe, G. Madec)  F90
8   !!             3.2  !  2009-04  (C. Ethe, G. Madec)  optimisaion
9   !!----------------------------------------------------------------------
10#if defined  key_pisces
11   !!----------------------------------------------------------------------
12   !!   'key_pisces'                                       PISCES bio-model
13   !!----------------------------------------------------------------------
14   !!   p4z_opt       : light availability in the water column
15   !!----------------------------------------------------------------------
16   USE trc            ! tracer variables
17   USE oce_trc        ! tracer-ocean share variables
18   USE sms_pisces     ! Source Minus Sink of PISCES
19   USE iom
20
21   IMPLICIT NONE
22   PRIVATE
23
24   PUBLIC   p4z_opt        ! called in p4zbio.F90 module
25   PUBLIC   p4z_opt_init   ! called in trcsms_pisces.F90 module
26
27   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   etot, enano, ediat   !: PAR for phyto, nano and diat
28   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   emoy                 !: averaged PAR in the mixed layer
29
30   INTEGER  ::  nksrp   ! levels below which the light cannot penetrate ( depth larger than 391 m)
31   REAL(wp) ::  parlux = 0.43 / 3.e0
32
33   REAL(wp), DIMENSION(3,61), PUBLIC ::   xkrgb  !: tabulated attenuation coefficients for RGB absorption
34   
35   !!* Substitution
36#  include "top_substitute.h90"
37   !!----------------------------------------------------------------------
38   !! NEMO/TOP 3.3 , NEMO Consortium (2010)
39   !! $Id$
40   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
41   !!----------------------------------------------------------------------
42
43CONTAINS
44
45   SUBROUTINE p4z_opt( kt, jnt )
46      !!---------------------------------------------------------------------
47      !!                     ***  ROUTINE p4z_opt  ***
48      !!
49      !! ** Purpose :   Compute the light availability in the water column
50      !!              depending on the depth and the chlorophyll concentration
51      !!
52      !! ** Method  : - ???
53      !!---------------------------------------------------------------------
54      INTEGER, INTENT(in) ::   kt, jnt ! ocean time step
55      INTEGER  ::   ji, jj, jk
56      INTEGER  ::   irgb
57      REAL(wp) ::   zchl, zxsi0r
58      REAL(wp) ::   zc0 , zc1 , zc2, zc3
59      REAL(wp), DIMENSION(jpi,jpj)     ::   zdepmoy, zetmp
60      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zekg, zekr, zekb
61      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   ze1 , ze2 , ze3, ze0
62      !!---------------------------------------------------------------------
63
64
65      !     Initialisation of variables used to compute PAR
66      !     -----------------------------------------------
67      ze1 (:,:,jpk) = 0.e0
68      ze2 (:,:,jpk) = 0.e0
69      ze3 (:,:,jpk) = 0.e0
70
71      !                                        !* attenuation coef. function of Chlorophyll and wavelength (Red-Green-Blue)
72      DO jk = 1, jpkm1                         !  --------------------------------------------------------
73!CDIR NOVERRCHK
74         DO jj = 1, jpj
75!CDIR NOVERRCHK
76            DO ji = 1, jpi
77               zchl = ( trn(ji,jj,jk,jpnch) + trn(ji,jj,jk,jpdch) + rtrn ) * 1.e6
78               zchl = MIN(  10. , MAX( 0.03, zchl )  )
79               irgb = NINT( 41 + 20.* LOG10( zchl ) + rtrn )
80               !                                                         
81               zekb(ji,jj,jk) = xkrgb(1,irgb) * fse3t(ji,jj,jk)
82               zekg(ji,jj,jk) = xkrgb(2,irgb) * fse3t(ji,jj,jk)
83               zekr(ji,jj,jk) = xkrgb(3,irgb) * fse3t(ji,jj,jk)
84            END DO
85         END DO
86      END DO
87
88!!gm  Potential BUG  must discuss with Olivier about this implementation....
89!!gm           the questions are : - PAR at T-point or mean PAR over T-level....
90!!gm                               - shallow water: no penetration of light through the bottom....
91
92
93      !                                        !* Photosynthetically Available Radiation (PAR)
94      !                                        !  --------------------------------------
95!CDIR NOVERRCHK
96      DO jj = 1, jpj
97!CDIR NOVERRCHK
98         DO ji = 1, jpi
99            zc1 = parlux * qsr(ji,jj) * EXP( -0.5 * zekb(ji,jj,1) )
100            zc2 = parlux * qsr(ji,jj) * EXP( -0.5 * zekg(ji,jj,1) )
101            zc3 = parlux * qsr(ji,jj) * EXP( -0.5 * zekr(ji,jj,1) )
102            ze1  (ji,jj,1) = zc1
103            ze2  (ji,jj,1) = zc2
104            ze3  (ji,jj,1) = zc3
105            etot (ji,jj,1) = (       zc1 +        zc2 +       zc3 )
106            enano(ji,jj,1) = ( 2.1 * zc1 + 0.42 * zc2 + 0.4 * zc3 )
107            ediat(ji,jj,1) = ( 1.6 * zc1 + 0.69 * zc2 + 0.7 * zc3 )
108         END DO
109      END DO
110
111   
112      DO jk = 2, nksrp     
113!CDIR NOVERRCHK
114         DO jj = 1, jpj
115!CDIR NOVERRCHK
116            DO ji = 1, jpi
117               zc1 = ze1(ji,jj,jk-1) * EXP( -0.5 * ( zekb(ji,jj,jk-1) + zekb(ji,jj,jk) ) )
118               zc2 = ze2(ji,jj,jk-1) * EXP( -0.5 * ( zekg(ji,jj,jk-1) + zekg(ji,jj,jk) ) )
119               zc3 = ze3(ji,jj,jk-1) * EXP( -0.5 * ( zekr(ji,jj,jk-1) + zekr(ji,jj,jk) ) )
120               ze1  (ji,jj,jk) = zc1
121               ze2  (ji,jj,jk) = zc2
122               ze3  (ji,jj,jk) = zc3
123               etot (ji,jj,jk) = (       zc1 +        zc2 +       zc3 )
124               enano(ji,jj,jk) = ( 2.1 * zc1 + 0.42 * zc2 + 0.4 * zc3 )
125               ediat(ji,jj,jk) = ( 1.6 * zc1 + 0.69 * zc2 + 0.7 * zc3 )
126            END DO
127         END DO
128      END DO
129
130      IF( ln_qsr_bio ) THEN                    !* heat flux accros w-level (used in the dynamics)
131         !                                     !  ------------------------
132         zxsi0r = 1.e0 / rn_si0
133         !
134         ze0  (:,:,1) = rn_abs * qsr(:,:)
135         ze1  (:,:,1) = parlux * qsr(:,:)             ! surface value : separation in R-G-B + near surface
136         ze2  (:,:,1) = parlux * qsr(:,:)
137         ze3  (:,:,1) = parlux * qsr(:,:)
138         etot3(:,:,1) =          qsr(:,:) * tmask(:,:,1)
139         !
140         DO jk = 2, nksrp+1
141!CDIR NOVERRCHK
142            DO jj = 1, jpj
143!CDIR NOVERRCHK
144               DO ji = 1, jpi
145                  zc0 = ze0(ji,jj,jk-1) * EXP( -fse3t(ji,jj,jk-1) * zxsi0r )
146                  zc1 = ze1(ji,jj,jk-1) * EXP( -zekb(ji,jj,jk-1 ) )
147                  zc2 = ze2(ji,jj,jk-1) * EXP( -zekg(ji,jj,jk-1 ) )
148                  zc3 = ze3(ji,jj,jk-1) * EXP( -zekr(ji,jj,jk-1 ) )
149                  ze0(ji,jj,jk) = zc0
150                  ze1(ji,jj,jk) = zc1
151                  ze2(ji,jj,jk) = zc2
152                  ze3(ji,jj,jk) = zc3
153                  etot3(ji,jj,jk) = ( zc0 + zc1 + zc2 + zc3 ) * tmask(ji,jj,jk)
154              END DO
155              !
156            END DO
157            !
158        END DO
159        !
160      ENDIF
161
162      !                                        !* Euphotic depth and level
163      neln(:,:) = 1                            !  ------------------------
164      heup(:,:) = 300.
165
166      DO jk = 2, nksrp
167         DO jj = 1, jpj
168           DO ji = 1, jpi
169              IF( etot(ji,jj,jk) >= 0.0043 * qsr(ji,jj) )  THEN
170                 neln(ji,jj) = jk+1                    ! Euphotic level : 1rst T-level strictly below Euphotic layer
171                 !                                     ! nb: ensure the compatibility with nmld_trc definition in trd_mld_trc_zint
172                 heup(ji,jj) = fsdepw(ji,jj,jk+1)      ! Euphotic layer depth
173              ENDIF
174           END DO
175        END DO
176      END DO
177 
178      heup(:,:) = MIN( 300., heup(:,:) )
179
180      !                                        !* mean light over the mixed layer
181      zdepmoy(:,:)   = 0.e0                    !  -------------------------------
182      zetmp  (:,:)   = 0.e0
183      emoy   (:,:,:) = 0.e0
184
185      DO jk = 1, nksrp
186!CDIR NOVERRCHK
187         DO jj = 1, jpj
188!CDIR NOVERRCHK
189            DO ji = 1, jpi
190               IF( fsdepw(ji,jj,jk+1) <= hmld(ji,jj) ) THEN
191                  zetmp  (ji,jj) = zetmp  (ji,jj) + etot(ji,jj,jk) * fse3t(ji,jj,jk)
192                  zdepmoy(ji,jj) = zdepmoy(ji,jj) + fse3t(ji,jj,jk)
193               ENDIF
194            END DO
195         END DO
196      END DO
197      !
198      emoy(:,:,:) = etot(:,:,:)
199      !
200      DO jk = 1, nksrp
201!CDIR NOVERRCHK
202         DO jj = 1, jpj
203!CDIR NOVERRCHK
204            DO ji = 1, jpi
205               IF( fsdepw(ji,jj,jk+1) <= hmld(ji,jj) ) &
206       &           emoy(ji,jj,jk) = zetmp(ji,jj) / ( zdepmoy(ji,jj) + rtrn )
207            END DO
208         END DO
209      END DO
210
211#if defined key_diatrc
212# if ! defined key_iomput
213      ! save for outputs
214      trc2d(:,:,  jp_pcs0_2d + 10) = heup(:,:  ) * tmask(:,:,1) 
215      trc3d(:,:,:,jp_pcs0_3d + 3)  = etot(:,:,:) * tmask(:,:,:)
216# else
217      ! write diagnostics
218      IF( jnt == nrdttrc ) then
219         CALL iom_put( "Heup", heup(:,:  ) * tmask(:,:,1) )  ! euphotic layer deptht
220         CALL iom_put( "PAR" , etot(:,:,:) * tmask(:,:,:) )  ! Photosynthetically Available Radiation
221      ENDIF
222# endif
223#endif
224      !
225   END SUBROUTINE p4z_opt
226
227   SUBROUTINE p4z_opt_init
228      !!----------------------------------------------------------------------
229      !!                  ***  ROUTINE p4z_opt_init  ***
230      !!
231      !! ** Purpose :   Initialization of tabulated attenuation coef
232      !!
233      !!
234      !!----------------------------------------------------------------------
235
236!!      CALL trc_oce_rgb( xkrgb )                  ! tabulated attenuation coefficients
237      CALL trc_oce_rgb_read( xkrgb )               ! tabulated attenuation coefficients
238      nksrp = trc_oce_ext_lev( r_si2, 0.33e2 )     ! max level of light extinction (Blue Chl=0.01)
239      IF(lwp) WRITE(numout,*) '        level of light extinction = ', nksrp, ' ref depth = ', gdepw_0(nksrp+1), ' m'
240      !
241      etot (:,:,:) = 0.e0
242      enano(:,:,:) = 0.e0
243      ediat(:,:,:) = 0.e0
244      IF( ln_qsr_bio )   etot3(:,:,:) = 0.e0
245
246      !
247   END SUBROUTINE p4z_opt_init
248#else
249   !!----------------------------------------------------------------------
250   !!  Dummy module :                                   No PISCES bio-model
251   !!----------------------------------------------------------------------
252CONTAINS
253   SUBROUTINE p4z_opt                   ! Empty routine
254   END SUBROUTINE p4z_opt
255#endif 
256
257   !!======================================================================
258END MODULE  p4zopt
Note: See TracBrowser for help on using the repository browser.