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 @ 1830

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

Computation of additional diagnostics for PISCES model ( under CPP key key_diaar5 )

  • needed for AR5 outputs (vertical inventories, passive tracers at surface,... )
  • new output file with suffix dbio_T
  • Property svn:executable set to *
  • Property svn:keywords set to Id
File size: 10.2 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      !!---------------------------------------------------------------------
64
65
66      !                                        !* tabulated attenuation coef.
67      IF( kt * jnt == nittrc000 ) THEN
68         !                                ! level of light extinction
69         nksrp = trc_oce_ext_lev( rn_si2, 0.33e2 )
70         IF(lwp) THEN
71           WRITE(numout,*)
72           WRITE(numout,*) ' level max of computation of qsr = ', nksrp, ' ref depth = ', gdepw_0(nksrp+1), ' m'
73         ENDIF
74!!         CALL trc_oce_rgb( xkrgb )     ! tabulated attenuation coefficients
75         CALL trc_oce_rgb_read( xkrgb )     ! tabulated attenuation coefficients
76         etot (:,:,:) = 0.e0
77         enano(:,:,:) = 0.e0
78         ediat(:,:,:) = 0.e0
79         IF( ln_qsr_bio ) etot3(:,:,:) = 0.e0
80      ENDIF
81
82
83!     Initialisation of variables used to compute PAR
84!     -----------------------------------------------
85      ze1 (:,:,jpk) = 0.e0
86      ze2 (:,:,jpk) = 0.e0
87      ze3 (:,:,jpk) = 0.e0
88
89      !                                        !* attenuation coef. function of Chlorophyll and wavelength (Red-Green-Blue)
90      DO jk = 1, jpkm1                         !  --------------------------------------------------------
91!CDIR NOVERRCHK
92         DO jj = 1, jpj
93!CDIR NOVERRCHK
94            DO ji = 1, jpi
95               zchl = ( trn(ji,jj,jk,jpnch) + trn(ji,jj,jk,jpdch) + rtrn ) * 1.e6
96               zchl = MIN(  10. , MAX( 0.03, zchl )  )
97               irgb = NINT( 41 + 20.* LOG10( zchl ) + rtrn )
98               !                                                         
99               zekb(ji,jj,jk) = xkrgb(1,irgb) * fse3t(ji,jj,jk)
100               zekg(ji,jj,jk) = xkrgb(2,irgb) * fse3t(ji,jj,jk)
101               zekr(ji,jj,jk) = xkrgb(3,irgb) * fse3t(ji,jj,jk)
102            END DO
103         END DO
104      END DO
105
106!!gm  Potential BUG  must discuss with Olivier about this implementation....
107!!gm           the questions are : - PAR at T-point or mean PAR over T-level....
108!!gm                               - shallow water: no penetration of light through the bottom....
109
110
111      !                                        !* Photosynthetically Available Radiation (PAR)
112      !                                        !  --------------------------------------
113!CDIR NOVERRCHK
114      DO jj = 1, jpj
115!CDIR NOVERRCHK
116         DO ji = 1, jpi
117            zc1 = parlux * qsr(ji,jj) * EXP( -0.5 * zekb(ji,jj,1) )
118            zc2 = parlux * qsr(ji,jj) * EXP( -0.5 * zekg(ji,jj,1) )
119            zc3 = parlux * qsr(ji,jj) * EXP( -0.5 * zekr(ji,jj,1) )
120            ze1  (ji,jj,1) = zc1
121            ze2  (ji,jj,1) = zc2
122            ze3  (ji,jj,1) = zc3
123            etot (ji,jj,1) = (       zc1 +        zc2 +       zc3 )
124            enano(ji,jj,1) = ( 2.1 * zc1 + 0.42 * zc2 + 0.4 * zc3 )
125            ediat(ji,jj,1) = ( 1.6 * zc1 + 0.69 * zc2 + 0.7 * zc3 )
126         END DO
127      END DO
128
129   
130      DO jk = 2, nksrp     
131!CDIR NOVERRCHK
132         DO jj = 1, jpj
133!CDIR NOVERRCHK
134            DO ji = 1, jpi
135               zc1 = ze1(ji,jj,jk-1) * EXP( -0.5 * ( zekb(ji,jj,jk-1) + zekb(ji,jj,jk) ) )
136               zc2 = ze2(ji,jj,jk-1) * EXP( -0.5 * ( zekg(ji,jj,jk-1) + zekg(ji,jj,jk) ) )
137               zc3 = ze3(ji,jj,jk-1) * EXP( -0.5 * ( zekr(ji,jj,jk-1) + zekr(ji,jj,jk) ) )
138               ze1  (ji,jj,jk) = zc1
139               ze2  (ji,jj,jk) = zc2
140               ze3  (ji,jj,jk) = zc3
141               etot (ji,jj,jk) = (       zc1 +        zc2 +       zc3 )
142               enano(ji,jj,jk) = ( 2.1 * zc1 + 0.42 * zc2 + 0.4 * zc3 )
143               ediat(ji,jj,jk) = ( 1.6 * zc1 + 0.69 * zc2 + 0.7 * zc3 )
144            END DO
145         END DO
146      END DO
147
148      IF( ln_qsr_bio ) THEN                    !* heat flux accros w-level (used in the dynamics)
149         !                                     !  ------------------------
150         zxsi0r = 1.e0 / rn_si0
151         !
152         ze0  (:,:,1) = rn_abs * qsr(:,:)
153         ze1  (:,:,1) = parlux * qsr(:,:)             ! surface value : separation in R-G-B + near surface
154         ze2  (:,:,1) = parlux * qsr(:,:)
155         ze3  (:,:,1) = parlux * qsr(:,:)
156         etot3(:,:,1) =          qsr(:,:) * tmask(:,:,1)
157         !
158         DO jk = 2, nksrp+1
159!CDIR NOVERRCHK
160            DO jj = 1, jpj
161!CDIR NOVERRCHK
162               DO ji = 1, jpi
163                  zc0 = ze0(ji,jj,jk-1) * EXP( -fse3t(ji,jj,jk-1) * zxsi0r )
164                  zc1 = ze1(ji,jj,jk-1) * EXP( -zekb(ji,jj,jk-1 ) )
165                  zc2 = ze2(ji,jj,jk-1) * EXP( -zekg(ji,jj,jk-1 ) )
166                  zc3 = ze3(ji,jj,jk-1) * EXP( -zekr(ji,jj,jk-1 ) )
167                  ze0(ji,jj,jk) = zc0
168                  ze1(ji,jj,jk) = zc1
169                  ze2(ji,jj,jk) = zc2
170                  ze3(ji,jj,jk) = zc3
171                  etot3(ji,jj,jk) = ( zc0 + zc1 + zc2 + zc3 ) * tmask(ji,jj,jk)
172              END DO
173              !
174            END DO
175            !
176        END DO
177        !
178      ENDIF
179
180      !                                        !* Euphotic depth and level
181      neln(:,:) = 1                            !  ------------------------
182      heup(:,:) = 300.
183
184      DO jk = 2, nksrp
185         DO jj = 1, jpj
186           DO ji = 1, jpi
187              IF( etot(ji,jj,jk) >= 0.0043 * qsr(ji,jj) )  THEN
188                 neln(ji,jj) = jk+1                    ! Euphotic level : 1rst T-level strictly below Euphotic layer
189                 !                                     ! nb: ensure the compatibility with nmld_trc definition in trd_mld_trc_zint
190                 heup(ji,jj) = fsdepw(ji,jj,jk+1)      ! Euphotic layer depth
191              ENDIF
192           END DO
193        END DO
194      END DO
195 
196      heup(:,:) = MIN( 300., heup(:,:) )
197
198      !                                        !* mean light over the mixed layer
199      zdepmoy(:,:)   = 0.e0                    !  -------------------------------
200      zetmp  (:,:)   = 0.e0
201      emoy   (:,:,:) = 0.e0
202
203      DO jk = 1, nksrp
204!CDIR NOVERRCHK
205         DO jj = 1, jpj
206!CDIR NOVERRCHK
207            DO ji = 1, jpi
208               IF( fsdepw(ji,jj,jk+1) <= hmld(ji,jj) ) THEN
209                  zetmp  (ji,jj) = zetmp  (ji,jj) + etot(ji,jj,jk) * fse3t(ji,jj,jk)
210                  zdepmoy(ji,jj) = zdepmoy(ji,jj) + fse3t(ji,jj,jk)
211               ENDIF
212            END DO
213         END DO
214      END DO
215      !
216      emoy(:,:,:) = etot(:,:,:)
217      !
218      DO jk = 1, nksrp
219!CDIR NOVERRCHK
220         DO jj = 1, jpj
221!CDIR NOVERRCHK
222            DO ji = 1, jpi
223               IF( fsdepw(ji,jj,jk+1) <= hmld(ji,jj) ) &
224       &           emoy(ji,jj,jk) = zetmp(ji,jj) / ( zdepmoy(ji,jj) + rtrn )
225            END DO
226         END DO
227      END DO
228
229#if defined key_trc_diaadd
230# if ! defined key_iomput
231      ! save for outputs
232      trc2d(:,:,  jp_pcs0_2d + 10) = heup(:,:  ) * tmask(:,:,1) 
233      trc3d(:,:,:,jp_pcs0_3d + 3)  = etot(:,:,:) * tmask(:,:,:)
234# else
235      ! write diagnostics
236      IF( jnt == nrdttrc ) then
237         CALL iom_put( "Heup", heup(:,:  ) * tmask(:,:,1) )  ! euphotic layer deptht
238         CALL iom_put( "PAR" , etot(:,:,:) * tmask(:,:,:) )  ! Photosynthetically Available Radiation
239      ENDIF
240# endif
241#endif
242      !
243   END SUBROUTINE p4z_opt
244
245#else
246   !!----------------------------------------------------------------------
247   !!  Dummy module :                                   No PISCES bio-model
248   !!----------------------------------------------------------------------
249CONTAINS
250   SUBROUTINE p4z_opt                   ! Empty routine
251   END SUBROUTINE p4z_opt
252#endif 
253
254   !!======================================================================
255END MODULE  p4zopt
Note: See TracBrowser for help on using the repository browser.