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/2011/dev_r2787_LOCEAN3_TRA_TRP/NEMOGCM/NEMO/TOP_SRC/PISCES – NEMO

source: branches/2011/dev_r2787_LOCEAN3_TRA_TRP/NEMOGCM/NEMO/TOP_SRC/PISCES/p4zopt.F90 @ 2819

Last change on this file since 2819 was 2819, checked in by cetlod, 13 years ago

Improvment of branch dev_r2787_LOCEAN3_TRA_TRP

  • Property svn:keywords set to Id
File size: 11.5 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)  optimisation
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   PUBLIC   p4z_opt_alloc
27
28   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: etot, enano, ediat   !: PAR for phyto, nano and diat
29   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: emoy                 !: averaged PAR in the mixed layer
30
31   INTEGER  ::   nksrp   ! levels below which the light cannot penetrate ( depth larger than 391 m)
32   REAL(wp) ::   parlux = 0.43_wp / 3._wp
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 3.3 , NEMO Consortium (2010)
40   !! $Id$
41   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
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      USE wrk_nemo, ONLY:   wrk_in_use, wrk_not_released
55      USE wrk_nemo, ONLY:   zdepmoy => wrk_2d_1 , zetmp => wrk_2d_2
56      USE wrk_nemo, ONLY:   zekg    => wrk_3d_2 , zekr  => wrk_3d_3 , zekb => wrk_3d_4
57      USE wrk_nemo, ONLY:   ze0     => wrk_3d_5 , ze1   => wrk_3d_6
58      USE wrk_nemo, ONLY:   ze2     => wrk_3d_7 , ze3   => wrk_3d_8
59      !
60      INTEGER, INTENT(in) ::   kt, jnt   ! ocean time step
61      !
62      INTEGER  ::   ji, jj, jk
63      INTEGER  ::   irgb
64      REAL(wp) ::   zchl, zxsi0r
65      REAL(wp) ::   zc0 , zc1 , zc2, zc3
66      !!---------------------------------------------------------------------
67
68      IF(  wrk_in_use(2, 1,2)   .OR.   wrk_in_use(3, 2,3,4,5,6,7,8)   ) THEN
69         CALL ctl_stop('p4z_opt: requested workspace arrays unavailable')   ;   RETURN
70      ENDIF
71
72      !     Initialisation of variables used to compute PAR
73      !     -----------------------------------------------
74      ze1 (:,:,jpk) = 0._wp
75      ze2 (:,:,jpk) = 0._wp
76      ze3 (:,:,jpk) = 0._wp
77
78      !                                        !* attenuation coef. function of Chlorophyll and wavelength (Red-Green-Blue)
79      DO jk = 1, jpkm1                         !  --------------------------------------------------------
80!CDIR NOVERRCHK
81         DO jj = 1, jpj
82!CDIR NOVERRCHK
83            DO ji = 1, jpi
84               zchl = ( trn(ji,jj,jk,jpnch) + trn(ji,jj,jk,jpdch) + rtrn ) * 1.e6
85               zchl = MIN(  10. , MAX( 0.03, zchl )  )
86               irgb = NINT( 41 + 20.* LOG10( zchl ) + rtrn )
87               !                                                         
88               zekb(ji,jj,jk) = xkrgb(1,irgb) * fse3t(ji,jj,jk)
89               zekg(ji,jj,jk) = xkrgb(2,irgb) * fse3t(ji,jj,jk)
90               zekr(ji,jj,jk) = xkrgb(3,irgb) * fse3t(ji,jj,jk)
91            END DO
92         END DO
93      END DO
94
95!!gm  Potential BUG  must discuss with Olivier about this implementation....
96!!gm           the questions are : - PAR at T-point or mean PAR over T-level....
97!!gm                               - shallow water: no penetration of light through the bottom....
98
99
100      !                                        !* Photosynthetically Available Radiation (PAR)
101      !                                        !  --------------------------------------
102!CDIR NOVERRCHK
103      DO jj = 1, jpj
104!CDIR NOVERRCHK
105         DO ji = 1, jpi
106            zc1 = parlux * qsr(ji,jj) * EXP( -0.5 * zekb(ji,jj,1) )
107            zc2 = parlux * qsr(ji,jj) * EXP( -0.5 * zekg(ji,jj,1) )
108            zc3 = parlux * qsr(ji,jj) * EXP( -0.5 * zekr(ji,jj,1) )
109            ze1  (ji,jj,1) = zc1
110            ze2  (ji,jj,1) = zc2
111            ze3  (ji,jj,1) = zc3
112            etot (ji,jj,1) = (       zc1 +        zc2 +       zc3 )
113            enano(ji,jj,1) = ( 2.1 * zc1 + 0.42 * zc2 + 0.4 * zc3 )
114            ediat(ji,jj,1) = ( 1.6 * zc1 + 0.69 * zc2 + 0.7 * zc3 )
115         END DO
116      END DO
117
118   
119      DO jk = 2, nksrp     
120!CDIR NOVERRCHK
121         DO jj = 1, jpj
122!CDIR NOVERRCHK
123            DO ji = 1, jpi
124               zc1 = ze1(ji,jj,jk-1) * EXP( -0.5 * ( zekb(ji,jj,jk-1) + zekb(ji,jj,jk) ) )
125               zc2 = ze2(ji,jj,jk-1) * EXP( -0.5 * ( zekg(ji,jj,jk-1) + zekg(ji,jj,jk) ) )
126               zc3 = ze3(ji,jj,jk-1) * EXP( -0.5 * ( zekr(ji,jj,jk-1) + zekr(ji,jj,jk) ) )
127               ze1  (ji,jj,jk) = zc1
128               ze2  (ji,jj,jk) = zc2
129               ze3  (ji,jj,jk) = zc3
130               etot (ji,jj,jk) = (       zc1 +        zc2 +       zc3 )
131               enano(ji,jj,jk) = ( 2.1 * zc1 + 0.42 * zc2 + 0.4 * zc3 )
132               ediat(ji,jj,jk) = ( 1.6 * zc1 + 0.69 * zc2 + 0.7 * zc3 )
133            END DO
134         END DO
135      END DO
136
137      IF( ln_qsr_bio ) THEN                    !* heat flux accros w-level (used in the dynamics)
138         !                                     !  ------------------------
139         zxsi0r = 1.e0 / rn_si0
140         !
141         ze0  (:,:,1) = rn_abs * qsr(:,:)
142         ze1  (:,:,1) = parlux * qsr(:,:)             ! surface value : separation in R-G-B + near surface
143         ze2  (:,:,1) = parlux * qsr(:,:)
144         ze3  (:,:,1) = parlux * qsr(:,:)
145         etot3(:,:,1) =          qsr(:,:) * tmask(:,:,1)
146         !
147         DO jk = 2, nksrp+1
148!CDIR NOVERRCHK
149            DO jj = 1, jpj
150!CDIR NOVERRCHK
151               DO ji = 1, jpi
152                  zc0 = ze0(ji,jj,jk-1) * EXP( -fse3t(ji,jj,jk-1) * zxsi0r )
153                  zc1 = ze1(ji,jj,jk-1) * EXP( -zekb(ji,jj,jk-1 ) )
154                  zc2 = ze2(ji,jj,jk-1) * EXP( -zekg(ji,jj,jk-1 ) )
155                  zc3 = ze3(ji,jj,jk-1) * EXP( -zekr(ji,jj,jk-1 ) )
156                  ze0(ji,jj,jk) = zc0
157                  ze1(ji,jj,jk) = zc1
158                  ze2(ji,jj,jk) = zc2
159                  ze3(ji,jj,jk) = zc3
160                  etot3(ji,jj,jk) = ( zc0 + zc1 + zc2 + zc3 ) * tmask(ji,jj,jk)
161              END DO
162              !
163            END DO
164            !
165        END DO
166        !
167      ENDIF
168
169      !                                        !* Euphotic depth and level
170      neln(:,:) = 1                            !  ------------------------
171      heup(:,:) = 300.
172
173      DO jk = 2, nksrp
174         DO jj = 1, jpj
175           DO ji = 1, jpi
176              IF( etot(ji,jj,jk) >= 0.0043 * qsr(ji,jj) )  THEN
177                 neln(ji,jj) = jk+1                    ! Euphotic level : 1rst T-level strictly below Euphotic layer
178                 !                                     ! nb: ensure the compatibility with nmld_trc definition in trd_mld_trc_zint
179                 heup(ji,jj) = fsdepw(ji,jj,jk+1)      ! Euphotic layer depth
180              ENDIF
181           END DO
182        END DO
183      END DO
184 
185      heup(:,:) = MIN( 300., heup(:,:) )
186
187      !                                        !* mean light over the mixed layer
188      zdepmoy(:,:)   = 0.e0                    !  -------------------------------
189      zetmp  (:,:)   = 0.e0
190      emoy   (:,:,:) = 0.e0
191
192      DO jk = 1, nksrp
193!CDIR NOVERRCHK
194         DO jj = 1, jpj
195!CDIR NOVERRCHK
196            DO ji = 1, jpi
197               IF( fsdepw(ji,jj,jk+1) <= hmld(ji,jj) ) THEN
198                  zetmp  (ji,jj) = zetmp  (ji,jj) + etot(ji,jj,jk) * fse3t(ji,jj,jk)
199                  zdepmoy(ji,jj) = zdepmoy(ji,jj) + fse3t(ji,jj,jk)
200               ENDIF
201            END DO
202         END DO
203      END DO
204      !
205      emoy(:,:,:) = etot(:,:,:)
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) )   emoy(ji,jj,jk) = zetmp(ji,jj) / ( zdepmoy(ji,jj) + rtrn )
213            END DO
214         END DO
215      END DO
216
217      IF( ln_diatrc ) THEN        ! save output diagnostics
218        !
219        IF( lk_iomput ) THEN
220           IF( jnt == nrdttrc ) THEN
221              CALL iom_put( "Heup", heup(:,:  ) * tmask(:,:,1) )  ! euphotic layer deptht
222              CALL iom_put( "PAR" , etot(:,:,:) * tmask(:,:,:) )  ! Photosynthetically Available Radiation
223           ENDIF
224        ELSE
225           trc2d(:,:,  jp_pcs0_2d + 10) = heup(:,:  ) * tmask(:,:,1) 
226           trc3d(:,:,:,jp_pcs0_3d + 3)  = etot(:,:,:) * tmask(:,:,:)
227        ENDIF
228        !
229      ENDIF
230      !
231      IF(  wrk_not_released(2, 1,2)           .OR.   &
232           wrk_not_released(3, 2,3,4,5,6,7,8)   )   CALL ctl_stop('p4z_opt: failed to release workspace arrays')
233      !
234   END SUBROUTINE p4z_opt
235
236
237   SUBROUTINE p4z_opt_init
238      !!----------------------------------------------------------------------
239      !!                  ***  ROUTINE p4z_opt_init  ***
240      !!
241      !! ** Purpose :   Initialization of tabulated attenuation coef
242      !!----------------------------------------------------------------------
243      !
244      CALL trc_oce_rgb( xkrgb )                  ! tabulated attenuation coefficients
245      nksrp = trc_oce_ext_lev( r_si2, 0.33e2 )     ! max level of light extinction (Blue Chl=0.01)
246      !
247      IF(lwp) WRITE(numout,*) '        level of light extinction = ', nksrp, ' ref depth = ', gdepw_0(nksrp+1), ' m'
248      !
249                         etot (:,:,:) = 0._wp
250                         enano(:,:,:) = 0._wp
251                         ediat(:,:,:) = 0._wp
252      IF( ln_qsr_bio )   etot3(:,:,:) = 0._wp
253      !
254   END SUBROUTINE p4z_opt_init
255
256
257   INTEGER FUNCTION p4z_opt_alloc()
258      !!----------------------------------------------------------------------
259      !!                     ***  ROUTINE p4z_opt_alloc  ***
260      !!----------------------------------------------------------------------
261      ALLOCATE( etot (jpi,jpj,jpk) , enano(jpi,jpj,jpk) ,     &
262         &      ediat(jpi,jpj,jpk) , emoy (jpi,jpj,jpk) , STAT=p4z_opt_alloc )
263         !
264      IF( p4z_opt_alloc /= 0 ) CALL ctl_warn('p4z_opt_alloc : failed to allocate arrays.')
265      !
266   END FUNCTION p4z_opt_alloc
267
268#else
269   !!----------------------------------------------------------------------
270   !!  Dummy module :                                   No PISCES bio-model
271   !!----------------------------------------------------------------------
272CONTAINS
273   SUBROUTINE p4z_opt                   ! Empty routine
274   END SUBROUTINE p4z_opt
275#endif 
276
277   !!======================================================================
278END MODULE  p4zopt
Note: See TracBrowser for help on using the repository browser.