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.
Changeset 5965 for branches/2014/dev_r4650_UKMO14.5_SST_BIAS_CORRECTION/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zopt.F90 – NEMO

Ignore:
Timestamp:
2015-12-01T16:35:30+01:00 (8 years ago)
Author:
timgraham
Message:

Upgraded branch to r5518 of trunk (v3.6 stable revision)

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2014/dev_r4650_UKMO14.5_SST_BIAS_CORRECTION/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zopt.F90

    r4624 r5965  
    3535   REAL(wp) :: parlux      !: Fraction of shortwave as PAR 
    3636   REAL(wp) :: xparsw                 !: parlux/3 
     37   REAL(wp) :: xsi0r                 !:  1. /rn_si0 
    3738 
    3839   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_par      ! structure of input par 
     
    4243 
    4344   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: enano, ediat   !: PAR for phyto, nano and diat  
     45   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: etot_ndcy      !: PAR over 24h in case of diurnal cycle 
    4446   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: emoy           !: averaged PAR in the mixed layer 
     47   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ekb, ekg, ekr  !: wavelength (Red-Green-Blue) 
    4548 
    4649   INTEGER  ::   nksrp   ! levels below which the light cannot penetrate ( depth larger than 391 m) 
     
    5760CONTAINS 
    5861 
    59    SUBROUTINE p4z_opt( kt, jnt ) 
     62   SUBROUTINE p4z_opt( kt, knt ) 
    6063      !!--------------------------------------------------------------------- 
    6164      !!                     ***  ROUTINE p4z_opt  *** 
     
    6770      !!--------------------------------------------------------------------- 
    6871      ! 
    69       INTEGER, INTENT(in) ::   kt, jnt   ! ocean time step 
     72      INTEGER, INTENT(in) ::   kt, knt   ! ocean time step 
    7073      ! 
    7174      INTEGER  ::   ji, jj, jk 
    7275      INTEGER  ::   irgb 
    73       REAL(wp) ::   zchl, zxsi0r 
     76      REAL(wp) ::   zchl 
    7477      REAL(wp) ::   zc0 , zc1 , zc2, zc3, z1_dep 
    75       REAL(wp), POINTER, DIMENSION(:,:  ) :: zdepmoy, zetmp, zetmp1, zetmp2 
    76       REAL(wp), POINTER, DIMENSION(:,:,:) :: zekg, zekr, zekb, ze0, ze1, ze2, ze3 
     78      REAL(wp), POINTER, DIMENSION(:,:  ) :: zdepmoy, zetmp1, zetmp2, zetmp3, zetmp4, zqsr100 
     79      REAL(wp), POINTER, DIMENSION(:,:,:) :: zpar, ze0, ze1, ze2, ze3 
    7780      !!--------------------------------------------------------------------- 
    7881      ! 
     
    8083      ! 
    8184      ! Allocate temporary workspace 
    82       CALL wrk_alloc( jpi, jpj,      zdepmoy, zetmp, zetmp1, zetmp2 )  
    83       CALL wrk_alloc( jpi, jpj, jpk, zekg, zekr, zekb, ze0, ze1, ze2, ze3 ) 
    84  
    85       IF( jnt == 1 .AND. ln_varpar ) CALL p4z_optsbc( kt ) 
     85      CALL wrk_alloc( jpi, jpj,      zqsr100, zdepmoy, zetmp1, zetmp2, zetmp3, zetmp4 ) 
     86      CALL wrk_alloc( jpi, jpj, jpk, zpar, ze0, ze1, ze2, ze3 ) 
     87 
     88      IF( knt == 1 .AND. ln_varpar ) CALL p4z_opt_sbc( kt ) 
    8689 
    8790      !     Initialisation of variables used to compute PAR 
    8891      !     ----------------------------------------------- 
    89       ze1(:,:,jpk) = 0._wp 
    90       ze2(:,:,jpk) = 0._wp 
    91       ze3(:,:,jpk) = 0._wp 
    92  
     92      ze1(:,:,:) = 0._wp 
     93      ze2(:,:,:) = 0._wp 
     94      ze3(:,:,:) = 0._wp 
    9395      !                                        !* attenuation coef. function of Chlorophyll and wavelength (Red-Green-Blue) 
    9496      DO jk = 1, jpkm1                         !  -------------------------------------------------------- 
     
    9799!CDIR NOVERRCHK 
    98100            DO ji = 1, jpi 
    99                zchl = ( trn(ji,jj,jk,jpnch) + trn(ji,jj,jk,jpdch) + rtrn ) * 1.e6 
     101               zchl = ( trb(ji,jj,jk,jpnch) + trb(ji,jj,jk,jpdch) + rtrn ) * 1.e6 
    100102               zchl = MIN(  10. , MAX( 0.05, zchl )  ) 
    101103               irgb = NINT( 41 + 20.* LOG10( zchl ) + rtrn ) 
    102104               !                                                          
    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) 
     105               ekb(ji,jj,jk) = xkrgb(1,irgb) * fse3t(ji,jj,jk) 
     106               ekg(ji,jj,jk) = xkrgb(2,irgb) * fse3t(ji,jj,jk) 
     107               ekr(ji,jj,jk) = xkrgb(3,irgb) * fse3t(ji,jj,jk) 
    106108            END DO 
    107109         END DO 
    108110      END DO 
    109  
    110  
    111111      !                                        !* Photosynthetically Available Radiation (PAR) 
    112112      !                                        !  -------------------------------------- 
    113  
    114       IF( ln_varpar ) THEN 
    115          ze1(:,:,1) = par_varsw(:,:) * qsr(:,:) * EXP( -0.5 * zekb(:,:,1) ) 
    116          ze2(:,:,1) = par_varsw(:,:) * qsr(:,:) * EXP( -0.5 * zekg(:,:,1) ) 
    117          ze3(:,:,1) = par_varsw(:,:) * qsr(:,:) * EXP( -0.5 * zekr(:,:,1) ) 
     113      IF( l_trcdm2dc ) THEN                     !  diurnal cycle 
     114         ! 1% of qsr to compute euphotic layer 
     115         zqsr100(:,:) = 0.01 * qsr_mean(:,:)     !  daily mean qsr 
     116         ! 
     117         CALL p4z_opt_par( kt, qsr_mean, ze1, ze2, ze3 )  
     118         ! 
     119         DO jk = 1, nksrp       
     120            etot_ndcy(:,:,jk) =        ze1(:,:,jk) +        ze2(:,:,jk) +       ze3(:,:,jk) 
     121            enano    (:,:,jk) =  2.1 * ze1(:,:,jk) + 0.42 * ze2(:,:,jk) + 0.4 * ze3(:,:,jk) 
     122            ediat    (:,:,jk) =  1.6 * ze1(:,:,jk) + 0.69 * ze2(:,:,jk) + 0.7 * ze3(:,:,jk) 
     123         END DO 
     124         ! 
     125         CALL p4z_opt_par( kt, qsr, ze1, ze2, ze3 )  
     126         ! 
     127         DO jk = 1, nksrp       
     128            etot(:,:,jk) =  ze1(:,:,jk) + ze2(:,:,jk) + ze3(:,:,jk) 
     129         END DO 
     130         ! 
    118131      ELSE 
    119          ze1(:,:,1) = xparsw         * qsr(:,:) * EXP( -0.5 * zekb(:,:,1) ) 
    120          ze2(:,:,1) = xparsw         * qsr(:,:) * EXP( -0.5 * zekg(:,:,1) ) 
    121          ze3(:,:,1) = xparsw         * qsr(:,:) * EXP( -0.5 * zekr(:,:,1) ) 
    122       ENDIF 
    123  
    124 !CDIR NOVERRCHK 
    125       DO jj = 1, jpj 
    126 !CDIR NOVERRCHK 
    127          DO ji = 1, jpi 
    128             zc1 = ze1(ji,jj,1) 
    129             zc2 = ze2(ji,jj,1)  
    130             zc3 = ze3(ji,jj,1) 
    131             etot (ji,jj,1) = (       zc1 +        zc2 +       zc3 ) 
    132             enano(ji,jj,1) = ( 2.1 * zc1 + 0.42 * zc2 + 0.4 * zc3 ) 
    133             ediat(ji,jj,1) = ( 1.6 * zc1 + 0.69 * zc2 + 0.7 * zc3 ) 
    134          END DO 
    135       END DO 
    136  
    137      
    138       DO jk = 2, nksrp       
    139 !CDIR NOVERRCHK 
    140          DO jj = 1, jpj 
    141 !CDIR NOVERRCHK 
    142             DO ji = 1, jpi 
    143                zc1 = ze1(ji,jj,jk-1) * EXP( -0.5 * ( zekb(ji,jj,jk-1) + zekb(ji,jj,jk) ) ) 
    144                zc2 = ze2(ji,jj,jk-1) * EXP( -0.5 * ( zekg(ji,jj,jk-1) + zekg(ji,jj,jk) ) ) 
    145                zc3 = ze3(ji,jj,jk-1) * EXP( -0.5 * ( zekr(ji,jj,jk-1) + zekr(ji,jj,jk) ) ) 
    146                ze1  (ji,jj,jk) = zc1 
    147                ze2  (ji,jj,jk) = zc2 
    148                ze3  (ji,jj,jk) = zc3 
    149                etot (ji,jj,jk) = (       zc1 +        zc2 +       zc3 ) 
    150                enano(ji,jj,jk) = ( 2.1 * zc1 + 0.42 * zc2 + 0.4 * zc3 ) 
    151                ediat(ji,jj,jk) = ( 1.6 * zc1 + 0.69 * zc2 + 0.7 * zc3 ) 
    152             END DO 
    153          END DO 
    154       END DO 
     132         ! 1% of qsr to compute euphotic layer 
     133         zqsr100(:,:) = 0.01 * qsr(:,:) 
     134         ! 
     135         CALL p4z_opt_par( kt, qsr, ze1, ze2, ze3 )  
     136         ! 
     137         DO jk = 1, nksrp       
     138            etot (:,:,jk) =        ze1(:,:,jk) +        ze2(:,:,jk) +       ze3(:,:,jk) 
     139            enano(:,:,jk) =  2.1 * ze1(:,:,jk) + 0.42 * ze2(:,:,jk) + 0.4 * ze3(:,:,jk) 
     140            ediat(:,:,jk) =  1.6 * ze1(:,:,jk) + 0.69 * ze2(:,:,jk) + 0.7 * ze3(:,:,jk) 
     141         END DO 
     142         etot_ndcy(:,:,:) =  etot(:,:,:)  
     143      ENDIF 
     144 
    155145 
    156146      IF( ln_qsr_bio ) THEN                    !* heat flux accros w-level (used in the dynamics) 
    157147         !                                     !  ------------------------ 
    158          zxsi0r = 1.e0 / rn_si0 
    159          ! 
    160          ze0(:,:,1) = rn_abs * qsr(:,:) 
    161          !                                                    ! surface value : separation in R-G-B + near surface 
    162          IF( ln_varpar ) THEN 
    163             ze0(:,:,1) = ( 1. - 3. * par_varsw(:,:) ) * qsr(:,:) 
    164             ze1(:,:,1) = par_varsw(:,:)               * qsr(:,:)          
    165             ze2(:,:,1) = par_varsw(:,:)               * qsr(:,:) 
    166             ze3(:,:,1) = par_varsw(:,:)               * qsr(:,:) 
    167          ELSE 
    168             ze0(:,:,1) = ( 1. - 3. * xparsw )  * qsr(:,:) 
    169             ze1(:,:,1) = xparsw                * qsr(:,:)          
    170             ze2(:,:,1) = xparsw                * qsr(:,:) 
    171             ze3(:,:,1) = xparsw                * qsr(:,:) 
    172          ENDIF 
     148         CALL p4z_opt_par( kt, qsr, ze1, ze2, ze3, pe0=ze0 ) 
     149         ! 
    173150         etot3(:,:,1) =  qsr(:,:) * tmask(:,:,1) 
    174          ! 
    175          ! 
    176151         DO jk = 2, nksrp + 1 
    177 !CDIR NOVERRCHK 
    178             DO jj = 1, jpj 
    179 !CDIR NOVERRCHK 
    180                DO ji = 1, jpi 
    181                   zc0 = ze0(ji,jj,jk-1) * EXP( -fse3t(ji,jj,jk-1) * zxsi0r ) 
    182                   zc1 = ze1(ji,jj,jk-1) * EXP( -zekb(ji,jj,jk-1 ) ) 
    183                   zc2 = ze2(ji,jj,jk-1) * EXP( -zekg(ji,jj,jk-1 ) ) 
    184                   zc3 = ze3(ji,jj,jk-1) * EXP( -zekr(ji,jj,jk-1 ) ) 
    185                   ze0(ji,jj,jk) = zc0 
    186                   ze1(ji,jj,jk) = zc1 
    187                   ze2(ji,jj,jk) = zc2 
    188                   ze3(ji,jj,jk) = zc3 
    189                   etot3(ji,jj,jk) = ( zc0 + zc1 + zc2 + zc3 ) * tmask(ji,jj,jk) 
    190               END DO 
    191               ! 
    192             END DO 
    193             ! 
    194         END DO 
    195         ! 
    196       ENDIF 
    197  
     152            etot3(:,:,jk) =  ( ze0(:,:,jk) + ze1(:,:,jk) + ze2(:,:,jk) + ze3(:,:,jk) ) * tmask(:,:,jk) 
     153         END DO 
     154         !                                     !  ------------------------ 
     155      ENDIF 
    198156      !                                        !* Euphotic depth and level 
    199157      neln(:,:) = 1                            !  ------------------------ 
     
    203161         DO jj = 1, jpj 
    204162           DO ji = 1, jpi 
    205               IF( etot(ji,jj,jk) * tmask(ji,jj,jk) >= 0.0043 * qsr(ji,jj) )  THEN 
     163              IF( etot_ndcy(ji,jj,jk) * tmask(ji,jj,jk) >= 0.43 * zqsr100(ji,jj) )  THEN 
    206164                 neln(ji,jj) = jk+1                    ! Euphotic level : 1rst T-level strictly below Euphotic layer 
    207165                 !                                     ! nb: ensure the compatibility with nmld_trc definition in trd_mld_trc_zint 
     
    211169        END DO 
    212170      END DO 
    213   
     171      ! 
    214172      heup(:,:) = MIN( 300., heup(:,:) ) 
    215  
    216173      !                                        !* mean light over the mixed layer 
    217174      zdepmoy(:,:)   = 0.e0                    !  ------------------------------- 
    218       zetmp  (:,:)   = 0.e0 
    219175      zetmp1 (:,:)   = 0.e0 
    220176      zetmp2 (:,:)   = 0.e0 
     177      zetmp3 (:,:)   = 0.e0 
     178      zetmp4 (:,:)   = 0.e0 
    221179 
    222180      DO jk = 1, nksrp 
     
    226184            DO ji = 1, jpi 
    227185               IF( fsdepw(ji,jj,jk+1) <= hmld(ji,jj) ) THEN 
    228                   zetmp  (ji,jj) = zetmp  (ji,jj) + etot (ji,jj,jk) * fse3t(ji,jj,jk) 
    229                   zetmp1 (ji,jj) = zetmp1 (ji,jj) + enano(ji,jj,jk) * fse3t(ji,jj,jk) 
    230                   zetmp2 (ji,jj) = zetmp2 (ji,jj) + ediat(ji,jj,jk) * fse3t(ji,jj,jk) 
     186                  zetmp1 (ji,jj) = zetmp1 (ji,jj) + etot     (ji,jj,jk) * fse3t(ji,jj,jk) ! remineralisation 
     187                  zetmp2 (ji,jj) = zetmp2 (ji,jj) + etot_ndcy(ji,jj,jk) * fse3t(ji,jj,jk) ! production 
     188                  zetmp3 (ji,jj) = zetmp3 (ji,jj) + enano    (ji,jj,jk) * fse3t(ji,jj,jk) ! production 
     189                  zetmp4 (ji,jj) = zetmp4 (ji,jj) + ediat    (ji,jj,jk) * fse3t(ji,jj,jk) ! production 
    231190                  zdepmoy(ji,jj) = zdepmoy(ji,jj) + fse3t(ji,jj,jk) 
    232191               ENDIF 
     
    235194      END DO 
    236195      ! 
    237       emoy(:,:,:) = etot(:,:,:) 
     196      emoy(:,:,:) = etot(:,:,:)       ! remineralisation 
     197      zpar(:,:,:) = etot_ndcy(:,:,:)  ! diagnostic : PAR with no diurnal cycle  
    238198      ! 
    239199      DO jk = 1, nksrp 
     
    244204               IF( fsdepw(ji,jj,jk+1) <= hmld(ji,jj) ) THEN 
    245205                  z1_dep = 1. / ( zdepmoy(ji,jj) + rtrn ) 
    246                   emoy (ji,jj,jk) = zetmp (ji,jj) * z1_dep 
    247                   enano(ji,jj,jk) = zetmp1(ji,jj) * z1_dep 
    248                   ediat(ji,jj,jk) = zetmp2(ji,jj) * z1_dep 
     206                  emoy (ji,jj,jk) = zetmp1(ji,jj) * z1_dep 
     207                  zpar (ji,jj,jk) = zetmp2(ji,jj) * z1_dep 
     208                  enano(ji,jj,jk) = zetmp3(ji,jj) * z1_dep 
     209                  ediat(ji,jj,jk) = zetmp4(ji,jj) * z1_dep 
    249210               ENDIF 
    250211            END DO 
    251212         END DO 
    252213      END DO 
    253  
    254       IF( ln_diatrc ) THEN        ! save output diagnostics 
     214      ! 
     215      IF( lk_iomput ) THEN 
     216        IF( knt == nrdttrc ) THEN 
     217           IF( iom_use( "Heup"  ) ) CALL iom_put( "Heup" , heup(:,:  ) * tmask(:,:,1) )  ! euphotic layer deptht 
     218           IF( iom_use( "PARDM" ) ) CALL iom_put( "PARDM", zpar(:,:,:) * tmask(:,:,:) )  ! Photosynthetically Available Radiation 
     219           IF( iom_use( "PAR"   ) ) CALL iom_put( "PAR"  , emoy(:,:,:) * tmask(:,:,:) )  ! Photosynthetically Available Radiation 
     220        ENDIF 
     221      ELSE 
     222         IF( ln_diatrc ) THEN        ! save output diagnostics 
     223            trc2d(:,:,  jp_pcs0_2d + 10) = heup(:,:  ) * tmask(:,:,1) 
     224            trc3d(:,:,:,jp_pcs0_3d + 3)  = etot(:,:,:) * tmask(:,:,:) 
     225         ENDIF 
     226      ENDIF 
     227      ! 
     228      CALL wrk_dealloc( jpi, jpj,      zqsr100, zdepmoy, zetmp1, zetmp2, zetmp3, zetmp4 ) 
     229      CALL wrk_dealloc( jpi, jpj, jpk, zpar,  ze0, ze1, ze2, ze3 ) 
     230      ! 
     231      IF( nn_timing == 1 )  CALL timing_stop('p4z_opt') 
     232      ! 
     233   END SUBROUTINE p4z_opt 
     234 
     235   SUBROUTINE p4z_opt_par( kt, pqsr, pe1, pe2, pe3, pe0 )  
     236      !!---------------------------------------------------------------------- 
     237      !!                  ***  routine p4z_opt_par  *** 
     238      !! 
     239      !! ** purpose :   compute PAR of each wavelength (Red-Green-Blue) 
     240      !!                for a given shortwave radiation 
     241      !! 
     242      !!---------------------------------------------------------------------- 
     243      !! * arguments 
     244      INTEGER, INTENT(in)                                       ::  kt            !   ocean time-step 
     245      REAL(wp), DIMENSION(jpi,jpj)    , INTENT(in)              ::  pqsr          !   shortwave 
     246      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout)           ::  pe1 , pe2 , pe3   !  PAR ( R-G-B) 
     247      REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout), OPTIONAL ::  pe0   
     248      !! * local variables 
     249      INTEGER    ::   ji, jj, jk     ! dummy loop indices 
     250      REAL(wp), DIMENSION(jpi,jpj)     ::  zqsr          !   shortwave 
     251      !!---------------------------------------------------------------------- 
     252 
     253      !  Real shortwave 
     254      IF( ln_varpar ) THEN  ;  zqsr(:,:) = par_varsw(:,:) * pqsr(:,:) 
     255      ELSE                  ;  zqsr(:,:) = xparsw         * pqsr(:,:) 
     256      ENDIF 
     257      ! 
     258      IF( PRESENT( pe0 ) ) THEN     !  W-level 
     259         ! 
     260         pe0(:,:,1) = pqsr(:,:) - 3. * zqsr(:,:)    !   ( 1 - 3 * alpha ) * q 
     261         pe1(:,:,1) = zqsr(:,:)          
     262         pe2(:,:,1) = zqsr(:,:) 
     263         pe3(:,:,1) = zqsr(:,:) 
     264         ! 
     265         DO jk = 2, nksrp + 1 
     266!CDIR NOVERRCHK 
     267            DO jj = 1, jpj 
     268!CDIR NOVERRCHK 
     269               DO ji = 1, jpi 
     270                  pe0(ji,jj,jk) = pe0(ji,jj,jk-1) * EXP( -fse3t(ji,jj,jk-1) * xsi0r ) 
     271                  pe1(ji,jj,jk) = pe1(ji,jj,jk-1) * EXP( -ekb(ji,jj,jk-1 ) ) 
     272                  pe2(ji,jj,jk) = pe2(ji,jj,jk-1) * EXP( -ekg(ji,jj,jk-1 ) ) 
     273                  pe3(ji,jj,jk) = pe3(ji,jj,jk-1) * EXP( -ekr(ji,jj,jk-1 ) ) 
     274               END DO 
     275              ! 
     276            END DO 
     277            ! 
     278         END DO 
    255279        ! 
    256         IF( lk_iomput ) THEN 
    257            IF( jnt == nrdttrc ) THEN 
    258               CALL iom_put( "Heup", heup(:,:  ) * tmask(:,:,1) )  ! euphotic layer deptht 
    259               CALL iom_put( "PAR" , emoy(:,:,:) * tmask(:,:,:) )  ! Photosynthetically Available Radiation 
    260            ENDIF 
    261         ELSE 
    262            trc2d(:,:,  jp_pcs0_2d + 10) = heup(:,:  ) * tmask(:,:,1)   
    263            trc3d(:,:,:,jp_pcs0_3d + 3)  = etot(:,:,:) * tmask(:,:,:) 
    264         ENDIF 
     280      ELSE   ! T- level 
    265281        ! 
    266       ENDIF 
    267       ! 
    268       CALL wrk_dealloc( jpi, jpj,      zdepmoy, zetmp, zetmp1, zetmp2 ) 
    269       CALL wrk_dealloc( jpi, jpj, jpk, zekg, zekr, zekb, ze0, ze1, ze2, ze3 ) 
    270       ! 
    271       IF( nn_timing == 1 )  CALL timing_stop('p4z_opt') 
    272       ! 
    273    END SUBROUTINE p4z_opt 
    274  
    275    SUBROUTINE p4z_optsbc( kt ) 
    276       !!---------------------------------------------------------------------- 
    277       !!                  ***  routine p4z_optsbc  *** 
     282        pe1(:,:,1) = zqsr(:,:) * EXP( -0.5 * ekb(:,:,1) ) 
     283        pe2(:,:,1) = zqsr(:,:) * EXP( -0.5 * ekg(:,:,1) ) 
     284        pe3(:,:,1) = zqsr(:,:) * EXP( -0.5 * ekr(:,:,1) ) 
     285        ! 
     286        DO jk = 2, nksrp       
     287!CDIR NOVERRCHK 
     288           DO jj = 1, jpj 
     289!CDIR NOVERRCHK 
     290              DO ji = 1, jpi 
     291                 pe1(ji,jj,jk) = pe1(ji,jj,jk-1) * EXP( -0.5 * ( ekb(ji,jj,jk-1) + ekb(ji,jj,jk) ) ) 
     292                 pe2(ji,jj,jk) = pe2(ji,jj,jk-1) * EXP( -0.5 * ( ekg(ji,jj,jk-1) + ekg(ji,jj,jk) ) ) 
     293                 pe3(ji,jj,jk) = pe3(ji,jj,jk-1) * EXP( -0.5 * ( ekr(ji,jj,jk-1) + ekr(ji,jj,jk) ) ) 
     294              END DO 
     295           END DO 
     296        END DO     
     297        ! 
     298      ENDIF 
     299      !  
     300   END SUBROUTINE p4z_opt_par 
     301 
     302 
     303   SUBROUTINE p4z_opt_sbc( kt ) 
     304      !!---------------------------------------------------------------------- 
     305      !!                  ***  routine p4z_opt_sbc  *** 
    278306      !! 
    279307      !! ** purpose :   read and interpolate the variable PAR fraction 
     
    286314      !!---------------------------------------------------------------------- 
    287315      !! * arguments 
    288       INTEGER, INTENT( in  ) ::   kt   ! ocean time step 
     316      INTEGER ,                INTENT(in) ::   kt     ! ocean time step 
    289317 
    290318      !! * local declarations 
     
    299327         IF( kt == nit000 .OR. ( kt /= nit000 .AND. ntimes_par > 1 ) ) THEN 
    300328            CALL fld_read( kt, 1, sf_par ) 
    301             par_varsw(:,:) = ( sf_par(1)%fnow(:,:,1) )/3.0 
     329            par_varsw(:,:) = ( sf_par(1)%fnow(:,:,1) ) / 3.0 
    302330         ENDIF 
    303331      ENDIF 
     
    305333      IF( nn_timing == 1 )  CALL timing_stop('p4z_optsbc') 
    306334      ! 
    307    END SUBROUTINE p4z_optsbc 
     335   END SUBROUTINE p4z_opt_sbc 
    308336 
    309337   SUBROUTINE p4z_opt_init 
     
    349377      ! 
    350378      xparsw = parlux / 3.0 
     379      xsi0r  = 1.e0 / rn_si0 
    351380      ! 
    352381      ! Variable PAR at the surface of the ocean 
     
    374403      IF(lwp) WRITE(numout,*) '        level of light extinction = ', nksrp, ' ref depth = ', gdepw_1d(nksrp+1), ' m' 
    375404      ! 
    376                          etot (:,:,:) = 0._wp 
    377                          enano(:,:,:) = 0._wp 
    378                          ediat(:,:,:) = 0._wp 
    379       IF( ln_qsr_bio )   etot3(:,:,:) = 0._wp 
     405                         ekr      (:,:,:) = 0._wp 
     406                         ekb      (:,:,:) = 0._wp 
     407                         ekg      (:,:,:) = 0._wp 
     408                         etot     (:,:,:) = 0._wp 
     409                         etot_ndcy(:,:,:) = 0._wp 
     410                         enano    (:,:,:) = 0._wp 
     411                         ediat    (:,:,:) = 0._wp 
     412      IF( ln_qsr_bio )   etot3    (:,:,:) = 0._wp 
    380413      !  
    381414      IF( nn_timing == 1 )  CALL timing_stop('p4z_opt_init') 
     
    388421      !!                     ***  ROUTINE p4z_opt_alloc  *** 
    389422      !!---------------------------------------------------------------------- 
    390       ALLOCATE( enano(jpi,jpj,jpk), ediat(jpi,jpj,jpk), emoy (jpi,jpj,jpk), STAT=p4z_opt_alloc )  
     423      ALLOCATE( ekb(jpi,jpj,jpk)      , ekr(jpi,jpj,jpk), ekg(jpi,jpj,jpk),   & 
     424        &       enano(jpi,jpj,jpk)    , ediat(jpi,jpj,jpk), & 
     425        &       etot_ndcy(jpi,jpj,jpk), emoy (jpi,jpj,jpk), STAT=p4z_opt_alloc )  
    391426         ! 
    392427      IF( p4z_opt_alloc /= 0 ) CALL ctl_warn('p4z_opt_alloc : failed to allocate arrays.') 
Note: See TracChangeset for help on using the changeset viewer.