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/CMIP5_IPSL/NEMO/TOP_SRC/PISCES – NEMO

source: branches/CMIP5_IPSL/NEMO/TOP_SRC/PISCES/p4zopt.F90 @ 1808

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

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