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/CNRS/dev_r6526_PISCES_GAS/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z – NEMO

source: branches/CNRS/dev_r6526_PISCES_GAS/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zopt.F90 @ 6530

Last change on this file since 6530 was 6530, checked in by cetlod, 8 years ago

1st implementation of PISCES GAS
File size: 22.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)  optimisation
9   !!             3.4  !  2011-06  (O. Aumont, C. Ethe) Improve light availability of nano & diat
10   !!----------------------------------------------------------------------
11#if defined  key_pisces
12   !!----------------------------------------------------------------------
13   !!   'key_pisces'                                       PISCES bio-model
14   !!----------------------------------------------------------------------
15   !!   p4z_opt       : light availability in the water column
16   !!----------------------------------------------------------------------
17   USE trc            ! tracer variables
18   USE oce_trc        ! tracer-ocean share variables
19   USE sms_pisces     ! Source Minus Sink of PISCES
20   USE iom            ! I/O manager
21   USE fldread         !  time interpolation
22   USE prtctl_trc      !  print control for debugging
23
24
25   IMPLICIT NONE
26   PRIVATE
27
28   PUBLIC   p4z_opt        ! called in p4zbio.F90 module
29   PUBLIC   p4z_opt_init   ! called in trcsms_pisces.F90 module
30   PUBLIC   p4z_opt_alloc
31
32   !! * Shared module variables
33
34   LOGICAL  :: ln_varpar   !: boolean for variable PAR fraction
35   REAL(wp) :: parlux      !: Fraction of shortwave as PAR
36   REAL(wp) :: xparsw                 !: parlux/3
37   REAL(wp) :: xsi0r                 !:  1. /rn_si0
38   !
39   REAL(wp) :: uvlux          ! proportion of UV radiant flux to global solar radiation
40   REAL(wp) :: xa350_a300 
41   REAL(wp) :: xa350_a400
42   REAL(wp) :: xa350_a440 
43
44   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_par      ! structure of input par
45   INTEGER , PARAMETER :: nbtimes = 365  !: maximum number of times record in a file
46   INTEGER  :: ntimes_par                ! number of time steps in a file
47   REAL(wp), ALLOCATABLE, SAVE,   DIMENSION(:,:) :: par_varsw    !: PAR fraction of shortwave
48
49   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: enano, ediat   !: PAR for phyto, nano and diat
50   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: etot_ndcy      !: PAR over 24h in case of diurnal cycle
51   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: emoy           !: averaged PAR in the mixed layer
52   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ekb, ekg, ekr  !: wavelength (Red-Green-Blue)
53#ifdef key_gas
54       REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: euv        !: averaged UV in the mixed layer
55       REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: a350       !: wavelength (350 for COS)
56#endif
57
58   INTEGER  ::   nksrp   ! levels below which the light cannot penetrate ( depth larger than 391 m)
59
60   REAL(wp), DIMENSION(3,61), PUBLIC ::   xkrgb   !: tabulated attenuation coefficients for RGB absorption
61   
62   !!* Substitution
63#  include "top_substitute.h90"
64   !!----------------------------------------------------------------------
65   !! NEMO/TOP 3.3 , NEMO Consortium (2010)
66   !! $Id: p4zopt.F90 3160 2011-11-20 14:27:18Z cetlod $
67   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
68   !!----------------------------------------------------------------------
69CONTAINS
70
71   SUBROUTINE p4z_opt( kt, knt )
72      !!---------------------------------------------------------------------
73      !!                     ***  ROUTINE p4z_opt  ***
74      !!
75      !! ** Purpose :   Compute the light availability in the water column
76      !!              depending on the depth and the chlorophyll concentration
77      !!
78      !! ** Method  : - ???
79      !!---------------------------------------------------------------------
80      !
81      INTEGER, INTENT(in) ::   kt, knt   ! ocean time step
82      !
83      INTEGER  ::   ji, jj, jk
84      INTEGER  ::   irgb
85      REAL(wp) ::   zchl
86      REAL(wp) ::   zc0 , zc1 , zc2, zc3, z1_dep
87      REAL(wp), POINTER, DIMENSION(:,:  ) :: zdepmoy, zetmp1, zetmp2, zetmp3, zetmp4, zqsr100
88      REAL(wp), POINTER, DIMENSION(:,:,:) :: zpar, ze0, ze1, ze2, ze3
89      REAL(wp) ::  za300, za400, za440, zpera440
90      !!---------------------------------------------------------------------
91      !
92      IF( nn_timing == 1 )  CALL timing_start('p4z_opt')
93      !
94      ! Allocate temporary workspace
95      CALL wrk_alloc( jpi, jpj,      zqsr100, zdepmoy, zetmp1, zetmp2, zetmp3, zetmp4 )
96      CALL wrk_alloc( jpi, jpj, jpk, zpar, ze0, ze1, ze2, ze3 )
97
98      IF( knt == 1 .AND. ln_varpar ) CALL p4z_opt_sbc( kt )
99
100      !     Initialisation of variables used to compute PAR
101      !     -----------------------------------------------
102      ze1(:,:,:) = 0._wp
103      ze2(:,:,:) = 0._wp
104      ze3(:,:,:) = 0._wp
105      !                                        !* attenuation coef. function of Chlorophyll and wavelength (Red-Green-Blue)
106      DO jk = 1, jpkm1                         !  --------------------------------------------------------
107!CDIR NOVERRCHK
108         DO jj = 1, jpj
109!CDIR NOVERRCHK
110            DO ji = 1, jpi
111               zchl = ( trb(ji,jj,jk,jpnch) + trb(ji,jj,jk,jpdch) + rtrn ) * 1.e6
112               zchl = MIN(  10. , MAX( 0.05, zchl )  )
113               irgb = NINT( 41 + 20.* LOG10( zchl ) + rtrn )
114               !                                                         
115               ekb(ji,jj,jk) = xkrgb(1,irgb) * fse3t(ji,jj,jk)
116               ekg(ji,jj,jk) = xkrgb(2,irgb) * fse3t(ji,jj,jk)
117               ekr(ji,jj,jk) = xkrgb(3,irgb) * fse3t(ji,jj,jk)
118            END DO
119         END DO
120      END DO
121      !
122      IF( lk_gas ) THEN
123         ! calculation of absorbance coefficient at 350 nm
124         DO jk = 1, jpk
125            DO jj = 1, jpj
126               DO ji = 1, jpi
127                  zchl = ( trb(ji,jj,jk,jpnch) + trb(ji,jj,jk,jpdch) + rtrn ) * 1.e6
128                  !za300 = 0.429 + 0.079 * LOG10( zchl ) ! Kettle 2001
129                  !a350(ji,jj,jk) = max(0.,za300 * xa350_a300)
130                  !
131                  !Xu 2001
132                  !a350(ji,jj,jk) = 0.0077 + 0.0115 * zchl ! Xu 2001
133                  !
134                  !!Morel
135                  !za400 = 0.065*zchl**(0.63)
136                  !a350(ji,jj,jk) = max(0.,za400 * xa350_a400)
137                  !
138                  ! Preiswerk
139                  zpera440 = (-26. * LOG10( zchl )) +26.
140                  za440 = ((zchl * 0.0448)*zpera440)/(100.-zpera440)
141                  a350(ji,jj,jk) = max(0.,za440 * xa350_a440)
142                  ! Fichot
143                  !C=log10(zchl)
144                  !a350(ji,jj,jk) =
145                  !exp((0.5346*C)-(0.0263*C*C)-(0.0036*C*C*C)+(0.0012*C*C*C*C)-1.6340)
146              END DO
147            END DO
148         END DO
149      ENDIF
150
151      !                                        !* Photosynthetically Available Radiation (PAR)
152      !                                        !  --------------------------------------
153      IF( l_trcdm2dc ) THEN                     !  diurnal cycle
154         ! 1% of qsr to compute euphotic layer
155         zqsr100(:,:) = 0.01 * qsr_mean(:,:)     !  daily mean qsr
156         !
157         CALL p4z_opt_par( kt, qsr_mean, ze1, ze2, ze3 ) 
158         !
159         DO jk = 1, nksrp     
160            etot_ndcy(:,:,jk) =        ze1(:,:,jk) +        ze2(:,:,jk) +       ze3(:,:,jk)
161            enano    (:,:,jk) =  2.1 * ze1(:,:,jk) + 0.42 * ze2(:,:,jk) + 0.4 * ze3(:,:,jk)
162            ediat    (:,:,jk) =  1.6 * ze1(:,:,jk) + 0.69 * ze2(:,:,jk) + 0.7 * ze3(:,:,jk)
163         END DO
164         !
165         CALL p4z_opt_par( kt, qsr, ze1, ze2, ze3 ) 
166         !
167         DO jk = 1, nksrp     
168            etot(:,:,jk) =  ze1(:,:,jk) + ze2(:,:,jk) + ze3(:,:,jk)
169         END DO
170         !
171      ELSE
172         ! 1% of qsr to compute euphotic layer
173         zqsr100(:,:) = 0.01 * qsr(:,:)
174         !
175         CALL p4z_opt_par( kt, qsr, ze1, ze2, ze3 ) 
176         !
177         DO jk = 1, nksrp     
178            etot (:,:,jk) =        ze1(:,:,jk) +        ze2(:,:,jk) +       ze3(:,:,jk)
179            enano(:,:,jk) =  2.1 * ze1(:,:,jk) + 0.42 * ze2(:,:,jk) + 0.4 * ze3(:,:,jk)
180            ediat(:,:,jk) =  1.6 * ze1(:,:,jk) + 0.69 * ze2(:,:,jk) + 0.7 * ze3(:,:,jk)
181         END DO
182         etot_ndcy(:,:,:) =  etot(:,:,:) 
183      ENDIF
184
185      IF( lk_gas ) THEN
186         IF( l_trcdm2dc ) THEN
187            euv(:,:,1) = uvlux * qsr_mean(:,:) * EXP( -0.5 * ze1(:,:,1) )
188         ELSE
189            euv(:,:,1) = uvlux * qsr(:,:) * EXP( -0.5 * ze1(:,:,1) )
190         ENDIF
191         DO jk = 2, nksrp
192!CDIR NOVERRCHK
193            DO jj = 1, jpj
194!CDIR NOVERRCHK
195               DO ji = 1, jpi
196                  euv(ji,jj,jk) = euv(ji,jj,jk-1) * EXP( -0.5 * ( ze1(ji,jj,jk-1) + ze1(ji,jj,jk) ) )
197               ENDDO
198            ENDDO
199         ENDDO
200      ENDIF
201
202
203      IF( ln_qsr_bio ) THEN                    !* heat flux accros w-level (used in the dynamics)
204         !                                     !  ------------------------
205         CALL p4z_opt_par( kt, qsr, ze1, ze2, ze3, pe0=ze0 )
206         !
207         etot3(:,:,1) =  qsr(:,:) * tmask(:,:,1)
208         DO jk = 2, nksrp + 1
209            etot3(:,:,jk) =  ( ze0(:,:,jk) + ze1(:,:,jk) + ze2(:,:,jk) + ze3(:,:,jk) ) * tmask(:,:,jk)
210         END DO
211         !                                     !  ------------------------
212      ENDIF
213      !                                        !* Euphotic depth and level
214      neln(:,:) = 1                            !  ------------------------
215      heup(:,:) = 300.
216
217      DO jk = 2, nksrp
218         DO jj = 1, jpj
219           DO ji = 1, jpi
220              IF( etot_ndcy(ji,jj,jk) * tmask(ji,jj,jk) >= 0.43 * zqsr100(ji,jj) )  THEN
221                 neln(ji,jj) = jk+1                    ! Euphotic level : 1rst T-level strictly below Euphotic layer
222                 !                                     ! nb: ensure the compatibility with nmld_trc definition in trd_mld_trc_zint
223                 heup(ji,jj) = fsdepw(ji,jj,jk+1)      ! Euphotic layer depth
224              ENDIF
225           END DO
226        END DO
227      END DO
228      !
229      heup(:,:) = MIN( 300., heup(:,:) )
230      !                                        !* mean light over the mixed layer
231      zdepmoy(:,:)   = 0.e0                    !  -------------------------------
232      zetmp1 (:,:)   = 0.e0
233      zetmp2 (:,:)   = 0.e0
234      zetmp3 (:,:)   = 0.e0
235      zetmp4 (:,:)   = 0.e0
236
237      DO jk = 1, nksrp
238!CDIR NOVERRCHK
239         DO jj = 1, jpj
240!CDIR NOVERRCHK
241            DO ji = 1, jpi
242               IF( fsdepw(ji,jj,jk+1) <= hmld(ji,jj) ) THEN
243                  zetmp1 (ji,jj) = zetmp1 (ji,jj) + etot     (ji,jj,jk) * fse3t(ji,jj,jk) ! remineralisation
244                  zetmp2 (ji,jj) = zetmp2 (ji,jj) + etot_ndcy(ji,jj,jk) * fse3t(ji,jj,jk) ! production
245                  zetmp3 (ji,jj) = zetmp3 (ji,jj) + enano    (ji,jj,jk) * fse3t(ji,jj,jk) ! production
246                  zetmp4 (ji,jj) = zetmp4 (ji,jj) + ediat    (ji,jj,jk) * fse3t(ji,jj,jk) ! production
247                  zdepmoy(ji,jj) = zdepmoy(ji,jj) + fse3t(ji,jj,jk)
248               ENDIF
249            END DO
250         END DO
251      END DO
252      !
253      emoy(:,:,:) = etot(:,:,:)       ! remineralisation
254      zpar(:,:,:) = etot_ndcy(:,:,:)  ! diagnostic : PAR with no diurnal cycle
255      !
256      DO jk = 1, nksrp
257!CDIR NOVERRCHK
258         DO jj = 1, jpj
259!CDIR NOVERRCHK
260            DO ji = 1, jpi
261               IF( fsdepw(ji,jj,jk+1) <= hmld(ji,jj) ) THEN
262                  z1_dep = 1. / ( zdepmoy(ji,jj) + rtrn )
263                  emoy (ji,jj,jk) = zetmp1(ji,jj) * z1_dep
264                  zpar (ji,jj,jk) = zetmp2(ji,jj) * z1_dep
265                  enano(ji,jj,jk) = zetmp3(ji,jj) * z1_dep
266                  ediat(ji,jj,jk) = zetmp4(ji,jj) * z1_dep
267               ENDIF
268            END DO
269         END DO
270      END DO
271      !
272      IF( lk_iomput ) THEN
273        IF( knt == nrdttrc ) THEN
274           IF( iom_use( "Heup"  ) ) CALL iom_put( "Heup" , heup(:,:  ) * tmask(:,:,1) )  ! euphotic layer deptht
275           IF( iom_use( "PARDM" ) ) CALL iom_put( "PARDM", zpar(:,:,:) * tmask(:,:,:) )  ! Photosynthetically Available Radiation
276           IF( iom_use( "PAR"   ) ) CALL iom_put( "PAR"  , emoy(:,:,:) * tmask(:,:,:) )  ! Photosynthetically Available Radiation
277           IF (lk_gas) THEN
278              IF( iom_use( "UV"   ) ) CALL iom_put( "UV", euv(:,:,:) * tmask(:,:,:) )  ! Photosynthetically Available Radiation
279              IF( iom_use( "A350" ) ) CALL iom_put( "A350"  , a350(:,:,:) * tmask(:,:,:) )  ! Photosynthetically Available Radiation
280           ENDIF
281        ENDIF
282      ELSE
283         IF( ln_diatrc ) THEN        ! save output diagnostics
284            trc2d(:,:,  jp_pcs0_2d + 10) = heup(:,:  ) * tmask(:,:,1)
285            trc3d(:,:,:,jp_pcs0_3d + 3)  = etot(:,:,:) * tmask(:,:,:)
286         ENDIF
287      ENDIF
288      !
289      CALL wrk_dealloc( jpi, jpj,      zqsr100, zdepmoy, zetmp1, zetmp2, zetmp3, zetmp4 )
290      CALL wrk_dealloc( jpi, jpj, jpk, zpar,  ze0, ze1, ze2, ze3 )
291      !
292      IF( nn_timing == 1 )  CALL timing_stop('p4z_opt')
293      !
294   END SUBROUTINE p4z_opt
295
296   SUBROUTINE p4z_opt_par( kt, pqsr, pe1, pe2, pe3, pe0 ) 
297      !!----------------------------------------------------------------------
298      !!                  ***  routine p4z_opt_par  ***
299      !!
300      !! ** purpose :   compute PAR of each wavelength (Red-Green-Blue)
301      !!                for a given shortwave radiation
302      !!
303      !!----------------------------------------------------------------------
304      !! * arguments
305      INTEGER, INTENT(in)                                       ::  kt            !   ocean time-step
306      REAL(wp), DIMENSION(jpi,jpj)    , INTENT(in)              ::  pqsr          !   shortwave
307      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout)           ::  pe1 , pe2 , pe3   !  PAR ( R-G-B)
308      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout), OPTIONAL ::  pe0 
309      !! * local variables
310      INTEGER    ::   ji, jj, jk     ! dummy loop indices
311      REAL(wp), DIMENSION(jpi,jpj)     ::  zqsr          !   shortwave
312      !!----------------------------------------------------------------------
313
314      !  Real shortwave
315      IF( ln_varpar ) THEN  ;  zqsr(:,:) = par_varsw(:,:) * pqsr(:,:)
316      ELSE                  ;  zqsr(:,:) = xparsw         * pqsr(:,:)
317      ENDIF
318      !
319      IF( PRESENT( pe0 ) ) THEN     !  W-level
320         !
321         pe0(:,:,1) = pqsr(:,:) - 3. * zqsr(:,:)    !   ( 1 - 3 * alpha ) * q
322         pe1(:,:,1) = zqsr(:,:)         
323         pe2(:,:,1) = zqsr(:,:)
324         pe3(:,:,1) = zqsr(:,:)
325         !
326         DO jk = 2, nksrp + 1
327!CDIR NOVERRCHK
328            DO jj = 1, jpj
329!CDIR NOVERRCHK
330               DO ji = 1, jpi
331                  pe0(ji,jj,jk) = pe0(ji,jj,jk-1) * EXP( -fse3t(ji,jj,jk-1) * xsi0r )
332                  pe1(ji,jj,jk) = pe1(ji,jj,jk-1) * EXP( -ekb(ji,jj,jk-1 ) )
333                  pe2(ji,jj,jk) = pe2(ji,jj,jk-1) * EXP( -ekg(ji,jj,jk-1 ) )
334                  pe3(ji,jj,jk) = pe3(ji,jj,jk-1) * EXP( -ekr(ji,jj,jk-1 ) )
335               END DO
336              !
337            END DO
338            !
339         END DO
340        !
341      ELSE   ! T- level
342        !
343        pe1(:,:,1) = zqsr(:,:) * EXP( -0.5 * ekb(:,:,1) )
344        pe2(:,:,1) = zqsr(:,:) * EXP( -0.5 * ekg(:,:,1) )
345        pe3(:,:,1) = zqsr(:,:) * EXP( -0.5 * ekr(:,:,1) )
346        !
347        DO jk = 2, nksrp     
348!CDIR NOVERRCHK
349           DO jj = 1, jpj
350!CDIR NOVERRCHK
351              DO ji = 1, jpi
352                 pe1(ji,jj,jk) = pe1(ji,jj,jk-1) * EXP( -0.5 * ( ekb(ji,jj,jk-1) + ekb(ji,jj,jk) ) )
353                 pe2(ji,jj,jk) = pe2(ji,jj,jk-1) * EXP( -0.5 * ( ekg(ji,jj,jk-1) + ekg(ji,jj,jk) ) )
354                 pe3(ji,jj,jk) = pe3(ji,jj,jk-1) * EXP( -0.5 * ( ekr(ji,jj,jk-1) + ekr(ji,jj,jk) ) )
355              END DO
356           END DO
357        END DO   
358        !
359      ENDIF
360      !
361   END SUBROUTINE p4z_opt_par
362
363
364   SUBROUTINE p4z_opt_sbc( kt )
365      !!----------------------------------------------------------------------
366      !!                  ***  routine p4z_opt_sbc  ***
367      !!
368      !! ** purpose :   read and interpolate the variable PAR fraction
369      !!                of shortwave radiation
370      !!
371      !! ** method  :   read the files and interpolate the appropriate variables
372      !!
373      !! ** input   :   external netcdf files
374      !!
375      !!----------------------------------------------------------------------
376      !! * arguments
377      INTEGER ,                INTENT(in) ::   kt     ! ocean time step
378
379      !! * local declarations
380      INTEGER  :: ji,jj
381      REAL(wp) :: zcoef
382      !!---------------------------------------------------------------------
383      !
384      IF( nn_timing == 1 )  CALL timing_start('p4z_optsbc')
385      !
386      ! Compute par_varsw at nit000 or only if there is more than 1 time record in par coefficient file
387      IF( ln_varpar ) THEN
388         IF( kt == nit000 .OR. ( kt /= nit000 .AND. ntimes_par > 1 ) ) THEN
389            CALL fld_read( kt, 1, sf_par )
390            par_varsw(:,:) = ( sf_par(1)%fnow(:,:,1) ) / 3.0
391         ENDIF
392      ENDIF
393      !
394      IF( nn_timing == 1 )  CALL timing_stop('p4z_optsbc')
395      !
396   END SUBROUTINE p4z_opt_sbc
397
398   SUBROUTINE p4z_opt_init
399      !!----------------------------------------------------------------------
400      !!                  ***  ROUTINE p4z_opt_init  ***
401      !!
402      !! ** Purpose :   Initialization of tabulated attenuation coef
403      !!                and of the percentage of PAR in Shortwave
404      !!
405      !! ** Input   :   external ascii and netcdf files
406      !!----------------------------------------------------------------------
407      !
408      INTEGER :: numpar
409      INTEGER :: ierr
410      INTEGER :: ios                 ! Local integer output status for namelist read
411      REAL(wp), DIMENSION(nbtimes) :: zsteps                 ! times records
412      !
413      CHARACTER(len=100) ::  cn_dir          ! Root directory for location of ssr files
414      TYPE(FLD_N) ::   sn_par                ! informations about the fields to be read
415      !
416      NAMELIST/nampisopt/cn_dir, sn_par, ln_varpar, &
417#ifdef key_gas
418               & uvlux, &
419#endif
420               & parlux
421      !!----------------------------------------------------------------------
422
423      IF( nn_timing == 1 )  CALL timing_start('p4z_opt_init')
424
425      REWIND( numnatp_ref )              ! Namelist nampisopt in reference namelist : Pisces attenuation coef. and PAR
426      READ  ( numnatp_ref, nampisopt, IOSTAT = ios, ERR = 901)
427901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampisopt in reference namelist', lwp )
428
429      REWIND( numnatp_cfg )              ! Namelist nampisopt in configuration namelist : Pisces attenuation coef. and PAR
430      READ  ( numnatp_cfg, nampisopt, IOSTAT = ios, ERR = 902 )
431902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'nampisopt in configuration namelist', lwp )
432      IF(lwm) WRITE ( numonp, nampisopt )
433
434      IF(lwp) THEN
435         WRITE(numout,*) ' '
436         WRITE(numout,*) ' namelist : nampisopt '
437         WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~ '
438         WRITE(numout,*) '    PAR as a variable fraction of SW     ln_varpar      = ', ln_varpar
439         WRITE(numout,*) '    Default value for the PAR fraction   parlux         = ', parlux
440         IF( lk_gas ) &
441       &  WRITE(numout,*) '    Proportion of UV radiant flux to global solar radiation uvlux   = ', uvlux
442      ENDIF
443      !
444      xparsw = parlux / 3.0
445      xsi0r  = 1.e0 / rn_si0
446      !
447      ! Variable PAR at the surface of the ocean
448      ! ----------------------------------------
449      IF( ln_varpar ) THEN
450         IF(lwp) WRITE(numout,*) '    initialize variable par fraction '
451         IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
452         !
453         ALLOCATE( par_varsw(jpi,jpj) )
454         !
455         ALLOCATE( sf_par(1), STAT=ierr )           !* allocate and fill sf_sst (forcing structure) with sn_sst
456         IF( ierr > 0 )   CALL ctl_stop( 'STOP', 'p4z_opt_init: unable to allocate sf_par structure' )
457         !
458         CALL fld_fill( sf_par, (/ sn_par /), cn_dir, 'p4z_opt_init', 'Variable PAR fraction ', 'nampisopt' )
459                                   ALLOCATE( sf_par(1)%fnow(jpi,jpj,1)   )
460         IF( sn_par%ln_tint )      ALLOCATE( sf_par(1)%fdta(jpi,jpj,1,2) )
461
462         CALL iom_open (  TRIM( sn_par%clname ) , numpar )
463         CALL iom_gettime( numpar, zsteps, kntime=ntimes_par)  ! get number of record in file
464      ENDIF
465      !
466      CALL trc_oce_rgb( xkrgb )                  ! tabulated attenuation coefficients
467      nksrp = trc_oce_ext_lev( r_si2, 0.33e2 )     ! max level of light extinction (Blue Chl=0.01)
468      !
469      IF(lwp) WRITE(numout,*) '        level of light extinction = ', nksrp, ' ref depth = ', gdepw_1d(nksrp+1), ' m'
470      !
471                         ekr      (:,:,:) = 0._wp
472                         ekb      (:,:,:) = 0._wp
473                         ekg      (:,:,:) = 0._wp
474                         etot     (:,:,:) = 0._wp
475                         etot_ndcy(:,:,:) = 0._wp
476                         enano    (:,:,:) = 0._wp
477                         ediat    (:,:,:) = 0._wp
478      IF( ln_qsr_bio )   etot3    (:,:,:) = 0._wp
479
480      IF( lk_gas )  THEN
481        xa350_a300 = 0.36787944117_wp   !  exp(-0.02*(350.-300.))  = exp(-1)
482        xa350_a400 = 2.459603_wp        !  exp(-0.018*(350.-400.))
483        xa350_a440 = 6.049647_wp        !  exp(-0.02*(350.-440.))
484        !
485        euv  (:,:,:) = 0._wp
486        a350 (:,:,:) = 0._wp
487      ENDIF
488      !
489      IF( nn_timing == 1 )  CALL timing_stop('p4z_opt_init')
490      !
491   END SUBROUTINE p4z_opt_init
492
493
494   INTEGER FUNCTION p4z_opt_alloc()
495      !!----------------------------------------------------------------------
496      !!                     ***  ROUTINE p4z_opt_alloc  ***
497      !!----------------------------------------------------------------------
498      INTEGER :: ierr(2)
499
500      ierr(:) = 0
501
502      ALLOCATE( ekb(jpi,jpj,jpk)      , ekr(jpi,jpj,jpk), ekg(jpi,jpj,jpk),   &
503        &       enano(jpi,jpj,jpk)    , ediat(jpi,jpj,jpk), &
504        &       etot_ndcy(jpi,jpj,jpk), emoy (jpi,jpj,jpk), ierr(1) ) 
505
506      IF( lk_gas ) ALLOCATE( euv(jpi,jpj,jpk), a350(jpi,jpj,jpk), ierr(2) )
507     
508      p4z_opt_alloc = MAXVAL( ierr )
509
510      IF( p4z_opt_alloc /= 0 ) CALL ctl_warn('p4z_opt_alloc : failed to allocate arrays.')
511      !
512   END FUNCTION p4z_opt_alloc
513
514#else
515   !!----------------------------------------------------------------------
516   !!  Dummy module :                                   No PISCES bio-model
517   !!----------------------------------------------------------------------
518CONTAINS
519   SUBROUTINE p4z_opt                   ! Empty routine
520   END SUBROUTINE p4z_opt
521#endif 
522
523   !!======================================================================
524END MODULE  p4zopt
Note: See TracBrowser for help on using the repository browser.