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 6140 for trunk/NEMOGCM/NEMO/OPA_SRC/TRA/zpshde.F90 – NEMO

Ignore:
Timestamp:
2015-12-21T12:35:23+01:00 (8 years ago)
Author:
timgraham
Message:

Merge of branches/2015/dev_merge_2015 back into trunk. Merge excludes NEMOGCM/TOOLS/OBSTOOLS/ for now due to issues with the change of file type. Will sort these manually with further commits.

Branch merged as follows:
In the working copy of branch ran:
svn merge ​svn+ssh://forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/trunk@HEAD
Small conflicts due to bug fixes applied to trunk since the dev_merge_2015 was copied. Bug fixes were applied to the branch as well so these were easy to resolve.
Branch committed at this stage

In working copy run:
svn switch ​svn+ssh://forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/trunk
to switch working copy

Run:
svn merge --reintegrate ​svn+ssh://forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/branches/2015/dev_merge_2015
to merge the branch into the trunk and then commit - no conflicts at this stage.

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/NEMOGCM/NEMO/OPA_SRC/TRA/zpshde.F90

    r5836 r6140  
    3232 
    3333   !! * Substitutions 
    34 #  include "domzgr_substitute.h90" 
    3534#  include "vectopt_loop_substitute.h90" 
    3635   !!---------------------------------------------------------------------- 
    … …  
    111110               iku = mbku(ji,jj)   ;   ikum1 = MAX( iku - 1 , 1 )    ! last and before last ocean level at u- & v-points 
    112111               ikv = mbkv(ji,jj)   ;   ikvm1 = MAX( ikv - 1 , 1 )    ! if level first is a p-step, ik.m1=1 
    113                ze3wu = fse3w(ji+1,jj  ,iku) - fse3w(ji,jj,iku) 
    114                ze3wv = fse3w(ji  ,jj+1,ikv) - fse3w(ji,jj,ikv) 
     112!!gm BUG ? when applied to before fields, e3w_b should be used.... 
     113               ze3wu = e3w_n(ji+1,jj  ,iku) - e3w_n(ji,jj,iku) 
     114               ze3wv = e3w_n(ji  ,jj+1,ikv) - e3w_n(ji,jj,ikv) 
    115115               ! 
    116116               ! i- direction 
    117117               IF( ze3wu >= 0._wp ) THEN      ! case 1 
    118                   zmaxu =  ze3wu / fse3w(ji+1,jj,iku) 
     118                  zmaxu =  ze3wu / e3w_n(ji+1,jj,iku) 
    119119                  ! interpolated values of tracers 
    120120                  zti (ji,jj,jn) = pta(ji+1,jj,iku,jn) + zmaxu * ( pta(ji+1,jj,ikum1,jn) - pta(ji+1,jj,iku,jn) ) 
    … …  
    122122                  pgtu(ji,jj,jn) = umask(ji,jj,1) * ( zti(ji,jj,jn) - pta(ji,jj,iku,jn) ) 
    123123               ELSE                           ! case 2 
    124                   zmaxu = -ze3wu / fse3w(ji,jj,iku) 
     124                  zmaxu = -ze3wu / e3w_n(ji,jj,iku) 
    125125                  ! interpolated values of tracers 
    126126                  zti (ji,jj,jn) = pta(ji,jj,iku,jn) + zmaxu * ( pta(ji,jj,ikum1,jn) - pta(ji,jj,iku,jn) ) 
    … …  
    131131               ! j- direction 
    132132               IF( ze3wv >= 0._wp ) THEN      ! case 1 
    133                   zmaxv =  ze3wv / fse3w(ji,jj+1,ikv) 
     133                  zmaxv =  ze3wv / e3w_n(ji,jj+1,ikv) 
    134134                  ! interpolated values of tracers 
    135135                  ztj (ji,jj,jn) = pta(ji,jj+1,ikv,jn) + zmaxv * ( pta(ji,jj+1,ikvm1,jn) - pta(ji,jj+1,ikv,jn) ) 
    … …  
    137137                  pgtv(ji,jj,jn) = vmask(ji,jj,1) * ( ztj(ji,jj,jn) - pta(ji,jj,ikv,jn) ) 
    138138               ELSE                           ! case 2 
    139                   zmaxv =  -ze3wv / fse3w(ji,jj,ikv) 
     139                  zmaxv =  -ze3wv / e3w_n(ji,jj,ikv) 
    140140                  ! interpolated values of tracers 
    141141                  ztj (ji,jj,jn) = pta(ji,jj,ikv,jn) + zmaxv * ( pta(ji,jj,ikvm1,jn) - pta(ji,jj,ikv,jn) ) 
    … …  
    156156               iku = mbku(ji,jj) 
    157157               ikv = mbkv(ji,jj) 
    158                ze3wu  = fse3w(ji+1,jj  ,iku) - fse3w(ji,jj,iku) 
    159                ze3wv  = fse3w(ji  ,jj+1,ikv) - fse3w(ji,jj,ikv) 
    160                IF( ze3wu >= 0._wp ) THEN   ;   zhi(ji,jj) = fsdept(ji  ,jj,iku)     ! i-direction: case 1 
    161                ELSE                        ;   zhi(ji,jj) = fsdept(ji+1,jj,iku)     ! -     -      case 2 
    162                ENDIF 
    163                IF( ze3wv >= 0._wp ) THEN   ;   zhj(ji,jj) = fsdept(ji,jj  ,ikv)     ! j-direction: case 1 
    164                ELSE                        ;   zhj(ji,jj) = fsdept(ji,jj+1,ikv)     ! -     -      case 2 
     158               ze3wu  = e3w_n(ji+1,jj  ,iku) - e3w_n(ji,jj,iku) 
     159               ze3wv  = e3w_n(ji  ,jj+1,ikv) - e3w_n(ji,jj,ikv) 
     160               IF( ze3wu >= 0._wp ) THEN   ;   zhi(ji,jj) = gdept_n(ji  ,jj,iku)     ! i-direction: case 1 
     161               ELSE                        ;   zhi(ji,jj) = gdept_n(ji+1,jj,iku)     ! -     -      case 2 
     162               ENDIF 
     163               IF( ze3wv >= 0._wp ) THEN   ;   zhj(ji,jj) = gdept_n(ji,jj  ,ikv)     ! j-direction: case 1 
     164               ELSE                        ;   zhj(ji,jj) = gdept_n(ji,jj+1,ikv)     ! -     -      case 2 
    165165               ENDIF 
    166166            END DO 
    … …  
    174174               iku = mbku(ji,jj) 
    175175               ikv = mbkv(ji,jj) 
    176                ze3wu  = fse3w(ji+1,jj  ,iku) - fse3w(ji,jj,iku) 
    177                ze3wv  = fse3w(ji  ,jj+1,ikv) - fse3w(ji,jj,ikv) 
     176               ze3wu  = e3w_n(ji+1,jj  ,iku) - e3w_n(ji,jj,iku) 
     177               ze3wv  = e3w_n(ji  ,jj+1,ikv) - e3w_n(ji,jj,ikv) 
    178178               IF( ze3wu >= 0._wp ) THEN   ;   pgru(ji,jj) = umask(ji,jj,1) * ( zri(ji  ,jj    ) - prd(ji,jj,iku) )   ! i: 1 
    179179               ELSE                        ;   pgru(ji,jj) = umask(ji,jj,1) * ( prd(ji+1,jj,iku) - zri(ji,jj    ) )   ! i: 2 
    … …  
    191191      ! 
    192192   END SUBROUTINE zps_hde 
    193  
    194  
    195    SUBROUTINE zps_hde_isf( kt, kjpt, pta, pgtu , pgtv , pgtui, pgtvi,                                   & 
    196       &                              prd, pgru , pgrv , pmru , pmrv , pgzu , pgzv , pge3ru , pge3rv ,   & 
    197       &                                   pgrui, pgrvi, pmrui, pmrvi, pgzui, pgzvi, pge3rui, pge3rvi ) 
    198       !!---------------------------------------------------------------------- 
    199       !!                     ***  ROUTINE zps_hde  *** 
     193   ! 
     194   SUBROUTINE zps_hde_isf( kt, kjpt, pta, pgtu, pgtv, pgtui, pgtvi,  & 
     195      &                          prd, pgru, pgrv, pgrui, pgrvi ) 
     196      !!---------------------------------------------------------------------- 
     197      !!                     ***  ROUTINE zps_hde_isf  *** 
    200198      !!                     
    201199      !! ** Purpose :   Compute the horizontal derivative of T, S and rho 
    202200      !!      at u- and v-points with a linear interpolation for z-coordinate 
    203       !!      with partial steps. 
     201      !!      with partial steps for top (ice shelf) and bottom. 
    204202      !! 
    205203      !! ** Method  :   In z-coord with partial steps, scale factors on last  
    206204      !!      levels are different for each grid point, so that T, S and rd  
    207205      !!      points are not at the same depth as in z-coord. To have horizontal 
    208       !!      gradients again, we interpolate T and S at the good depth :  
     206      !!      gradients again, we interpolate T and S at the good depth : 
     207      !!      For the bottom case: 
    209208      !!      Linear interpolation of T, S    
    210209      !!         Computation of di(tb) and dj(tb) by vertical interpolation: 
    … …  
    235234      !!          di(rho) = rd~ - rd(i,j,k)   or   rd(i+1,j,k) - rd~ 
    236235      !! 
     236      !!      For the top case (ice shelf): As for the bottom case but upside down 
     237      !! 
    237238      !! ** Action  : compute for top and bottom interfaces 
    238239      !!              - pgtu, pgtv, pgtui, pgtvi: horizontal gradient of tracer at u- & v-points 
    239240      !!              - pgru, pgrv, pgrui, pgtvi: horizontal gradient of rho (if present) at u- & v-points 
    240       !!              - pmru, pmrv, pmrui, pmrvi: horizontal sum of rho at u- & v- point (used in dynhpg with vvl) 
    241       !!              - pgzu, pgzv, pgzui, pgzvi: horizontal gradient of z at u- and v- point (used in dynhpg with vvl) 
    242       !!              - pge3ru, pge3rv, pge3rui, pge3rvi: horizontal gradient of rho weighted by local e3w at u- & v-points  
    243       !!---------------------------------------------------------------------- 
    244       INTEGER                              , INTENT(in   )           ::  kt                ! ocean time-step index 
    245       INTEGER                              , INTENT(in   )           ::  kjpt              ! number of tracers 
    246       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   )           ::  pta               ! 4D tracers fields 
    247       !                                                              !!  u-point ! v-point ! 
    248       REAL(wp), DIMENSION(jpi,jpj,    kjpt), INTENT(  out)           ::  pgtu    , pgtv    ! bottom GRADh( ptra )   
    249       REAL(wp), DIMENSION(jpi,jpj,    kjpt), INTENT(  out)           ::  pgtui   , pgtvi   ! top    GRADh( ptra ) 
    250       REAL(wp), DIMENSION(jpi,jpj,jpk     ), INTENT(in   ), OPTIONAL ::  prd               ! 3D density anomaly fields 
    251       !                                                              !!  u-point ! v-point ! 
    252       REAL(wp), DIMENSION(jpi,jpj         ), INTENT(  out), OPTIONAL ::  pgru    , pgrv    ! bottom GRADh( prd  ) 
    253       REAL(wp), DIMENSION(jpi,jpj         ), INTENT(  out), OPTIONAL ::  pmru    , pmrv    ! bottom SUM  ( prd  ) 
    254       REAL(wp), DIMENSION(jpi,jpj         ), INTENT(  out), OPTIONAL ::  pgzu    , pgzv    ! bottom GRADh( z    )  
    255       REAL(wp), DIMENSION(jpi,jpj         ), INTENT(  out), OPTIONAL ::  pge3ru  , pge3rv  ! bottom GRADh( prd  ) weighted by e3w 
    256       REAL(wp), DIMENSION(jpi,jpj         ), INTENT(  out), OPTIONAL ::  pgrui   , pgrvi   ! top    GRADh( prd  )  
    257       REAL(wp), DIMENSION(jpi,jpj         ), INTENT(  out), OPTIONAL ::  pmrui   , pmrvi   ! top    SUM  ( prd  )  
    258       REAL(wp), DIMENSION(jpi,jpj         ), INTENT(  out), OPTIONAL ::  pgzui   , pgzvi   ! top    GRADh( z    )  
    259       REAL(wp), DIMENSION(jpi,jpj         ), INTENT(  out), OPTIONAL ::  pge3rui , pge3rvi ! top    GRADh( prd  ) weighted by e3w 
     241      !!---------------------------------------------------------------------- 
     242      INTEGER                              , INTENT(in   )           ::  kt           ! ocean time-step index 
     243      INTEGER                              , INTENT(in   )           ::  kjpt         ! number of tracers 
     244      REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   )           ::  pta          ! 4D tracers fields 
     245      REAL(wp), DIMENSION(jpi,jpj,    kjpt), INTENT(  out)           ::  pgtu, pgtv   ! hor. grad. of ptra at u- & v-pts  
     246      REAL(wp), DIMENSION(jpi,jpj,    kjpt), INTENT(  out)           ::  pgtui, pgtvi ! hor. grad. of stra at u- & v-pts (ISF) 
     247      REAL(wp), DIMENSION(jpi,jpj,jpk     ), INTENT(in   ), OPTIONAL ::  prd          ! 3D density anomaly fields 
     248      REAL(wp), DIMENSION(jpi,jpj         ), INTENT(  out), OPTIONAL ::  pgru, pgrv   ! hor. grad of prd at u- & v-pts (bottom) 
     249      REAL(wp), DIMENSION(jpi,jpj         ), INTENT(  out), OPTIONAL ::  pgrui, pgrvi ! hor. grad of prd at u- & v-pts (top) 
    260250      ! 
    261251      INTEGER  ::   ji, jj, jn      ! Dummy loop indices 
    262252      INTEGER  ::   iku, ikv, ikum1, ikvm1,ikup1, ikvp1   ! partial step level (ocean bottom level) at u- and v-points 
    263       REAL(wp) ::  ze3wu, ze3wv, zmaxu, zmaxv, zdzwu, zdzwv, zdzwuip1, zdzwvjp1  ! temporary scalars 
     253      REAL(wp) ::  ze3wu, ze3wv, zmaxu, zmaxv             ! temporary scalars 
    264254      REAL(wp), DIMENSION(jpi,jpj)      ::  zri, zrj, zhi, zhj   ! NB: 3rd dim=1 to use eos 
    265255      REAL(wp), DIMENSION(jpi,jpj,kjpt) ::  zti, ztj             !  
    … …  
    277267         DO jj = 1, jpjm1 
    278268            DO ji = 1, jpim1 
    279                iku = mbku(ji,jj)   ;   ikum1 = MAX( iku - 1 , 1 )    ! last and before last ocean level at u- & v-points 
    280                ikv = mbkv(ji,jj)   ;   ikvm1 = MAX( ikv - 1 , 1 )    ! if level first is a p-step, ik.m1=1 
     269 
     270               iku = mbku(ji,jj); ikum1 = MAX( iku - 1 , 1 )    ! last and before last ocean level at u- & v-points 
     271               ikv = mbkv(ji,jj); ikvm1 = MAX( ikv - 1 , 1 )    ! if level first is a p-step, ik.m1=1 
     272               ze3wu = gdept_n(ji+1,jj,iku) - gdept_n(ji,jj,iku) 
     273               ze3wv = gdept_n(ji,jj+1,ikv) - gdept_n(ji,jj,ikv) 
     274               ! 
     275               ! i- direction 
     276               IF( ze3wu >= 0._wp ) THEN      ! case 1 
     277                  zmaxu =  ze3wu / e3w_n(ji+1,jj,iku) 
     278                  ! interpolated values of tracers 
     279                  zti (ji,jj,jn) = pta(ji+1,jj,iku,jn) + zmaxu * ( pta(ji+1,jj,ikum1,jn) - pta(ji+1,jj,iku,jn) ) 
     280                  ! gradient of  tracers 
     281                  pgtu(ji,jj,jn) = ssumask(ji,jj) * ( zti(ji,jj,jn) - pta(ji,jj,iku,jn) ) 
     282               ELSE                           ! case 2 
     283                  zmaxu = -ze3wu / e3w_n(ji,jj,iku) 
     284                  ! interpolated values of tracers 
     285                  zti (ji,jj,jn) = pta(ji,jj,iku,jn) + zmaxu * ( pta(ji,jj,ikum1,jn) - pta(ji,jj,iku,jn) ) 
     286                  ! gradient of tracers 
     287                  pgtu(ji,jj,jn) = ssumask(ji,jj) * ( pta(ji+1,jj,iku,jn) - zti(ji,jj,jn) ) 
     288               ENDIF 
     289               ! 
     290               ! j- direction 
     291               IF( ze3wv >= 0._wp ) THEN      ! case 1 
     292                  zmaxv =  ze3wv / e3w_n(ji,jj+1,ikv) 
     293                  ! interpolated values of tracers 
     294                  ztj (ji,jj,jn) = pta(ji,jj+1,ikv,jn) + zmaxv * ( pta(ji,jj+1,ikvm1,jn) - pta(ji,jj+1,ikv,jn) ) 
     295                  ! gradient of tracers 
     296                  pgtv(ji,jj,jn) = ssvmask(ji,jj) * ( ztj(ji,jj,jn) - pta(ji,jj,ikv,jn) ) 
     297               ELSE                           ! case 2 
     298                  zmaxv =  -ze3wv / e3w_n(ji,jj,ikv) 
     299                  ! interpolated values of tracers 
     300                  ztj (ji,jj,jn) = pta(ji,jj,ikv,jn) + zmaxv * ( pta(ji,jj,ikvm1,jn) - pta(ji,jj,ikv,jn) ) 
     301                  ! gradient of tracers 
     302                  pgtv(ji,jj,jn) = ssvmask(ji,jj) * ( pta(ji,jj+1,ikv,jn) - ztj(ji,jj,jn) ) 
     303               ENDIF 
     304 
     305            END DO 
     306         END DO 
     307         CALL lbc_lnk( pgtu(:,:,jn), 'U', -1. )   ;   CALL lbc_lnk( pgtv(:,:,jn), 'V', -1. )   ! Lateral boundary cond. 
     308         ! 
     309      END DO 
     310 
     311      ! horizontal derivative of density anomalies (rd) 
     312      IF( PRESENT( prd ) ) THEN         ! depth of the partial step level 
     313         pgru(:,:)=0.0_wp   ; pgrv(:,:)=0.0_wp ;  
     314         ! 
     315         DO jj = 1, jpjm1 
     316            DO ji = 1, jpim1 
     317 
     318               iku = mbku(ji,jj) 
     319               ikv = mbkv(ji,jj) 
     320               ze3wu = gdept_n(ji+1,jj,iku) - gdept_n(ji,jj,iku) 
     321               ze3wv = gdept_n(ji,jj+1,ikv) - gdept_n(ji,jj,ikv) 
     322               ! 
     323               IF( ze3wu >= 0._wp ) THEN   ;   zhi(ji,jj) = gdept_n(ji  ,jj,iku)    ! i-direction: case 1 
     324               ELSE                        ;   zhi(ji,jj) = gdept_n(ji+1,jj,iku)    ! -     -      case 2 
     325               ENDIF 
     326               IF( ze3wv >= 0._wp ) THEN   ;   zhj(ji,jj) = gdept_n(ji,jj  ,ikv)    ! j-direction: case 1 
     327               ELSE                        ;   zhj(ji,jj) = gdept_n(ji,jj+1,ikv)    ! -     -      case 2 
     328               ENDIF 
     329 
     330            END DO 
     331         END DO 
     332 
     333         ! Compute interpolated rd from zti, ztj for the 2 cases at the depth of the partial 
     334         ! step and store it in  zri, zrj for each  case 
     335         CALL eos( zti, zhi, zri ) 
     336         CALL eos( ztj, zhj, zrj ) 
     337 
     338         DO jj = 1, jpjm1                 ! Gradient of density at the last level  
     339            DO ji = 1, jpim1 
     340               iku = mbku(ji,jj) 
     341               ikv = mbkv(ji,jj) 
     342               ze3wu = gdept_n(ji+1,jj,iku) - gdept_n(ji,jj,iku) 
     343               ze3wv = gdept_n(ji,jj+1,ikv) - gdept_n(ji,jj,ikv) 
     344 
     345               IF( ze3wu >= 0._wp ) THEN   ;   pgru(ji,jj) = ssumask(ji,jj) * ( zri(ji  ,jj    ) - prd(ji,jj,iku) )   ! i: 1 
     346               ELSE                        ;   pgru(ji,jj) = ssumask(ji,jj) * ( prd(ji+1,jj,iku) - zri(ji,jj    ) )   ! i: 2 
     347               ENDIF 
     348               IF( ze3wv >= 0._wp ) THEN   ;   pgrv(ji,jj) = ssvmask(ji,jj) * ( zrj(ji,jj      ) - prd(ji,jj,ikv) )   ! j: 1 
     349               ELSE                        ;   pgrv(ji,jj) = ssvmask(ji,jj) * ( prd(ji,jj+1,ikv) - zrj(ji,jj    ) )   ! j: 2 
     350               ENDIF 
     351 
     352            END DO 
     353         END DO 
     354 
     355         CALL lbc_lnk( pgru , 'U', -1. )   ;   CALL lbc_lnk( pgrv , 'V', -1. )   ! Lateral boundary conditions 
     356         ! 
     357      END IF 
     358      ! 
     359      !     !==  (ISH)  compute grui and gruvi  ==! 
     360      ! 
     361      DO jn = 1, kjpt      !==   Interpolation of tracers at the last ocean level   ==!            ! 
     362         DO jj = 1, jpjm1 
     363            DO ji = 1, jpim1 
     364               iku = miku(ji,jj); ikup1 = miku(ji,jj) + 1 
     365               ikv = mikv(ji,jj); ikvp1 = mikv(ji,jj) + 1 
     366               ! 
    281367               ! (ISF) case partial step top and bottom in adjacent cell in vertical 
    282368               ! cannot used e3w because if 2 cell water column, we have ps at top and bottom 
    283369               ! in this case e3w(i,j) - e3w(i,j+1) is not the distance between Tj~ and Tj 
    284370               ! the only common depth between cells (i,j) and (i,j+1) is gdepw_0 
    285                ze3wu  = (gdept_0(ji+1,jj,iku) - gdepw_0(ji+1,jj,iku)) - (gdept_0(ji,jj,iku) - gdepw_0(ji,jj,iku)) 
    286                ze3wv  = (gdept_0(ji,jj+1,ikv) - gdepw_0(ji,jj+1,ikv)) - (gdept_0(ji,jj,ikv) - gdepw_0(ji,jj,ikv)) 
    287                ! 
     371               ze3wu  =  gdept_n(ji,jj,iku) - gdept_n(ji+1,jj,iku) 
     372               ze3wv  =  gdept_n(ji,jj,ikv) - gdept_n(ji,jj+1,ikv)  
     373 
    288374               ! i- direction 
    289375               IF( ze3wu >= 0._wp ) THEN      ! case 1 
    290                   zmaxu =  ze3wu / fse3w(ji+1,jj,iku) 
    291                   ! interpolated values of tracers 
    292                   zti (ji,jj,jn) = pta(ji+1,jj,iku,jn) + zmaxu * ( pta(ji+1,jj,ikum1,jn) - pta(ji+1,jj,iku,jn) ) 
     376                  zmaxu = ze3wu / e3w_n(ji+1,jj,ikup1) 
     377                  ! interpolated values of tracers 
     378                  zti(ji,jj,jn) = pta(ji+1,jj,iku,jn) + zmaxu * ( pta(ji+1,jj,ikup1,jn) - pta(ji+1,jj,iku,jn) ) 
     379                  ! gradient of tracers 
     380                  pgtui(ji,jj,jn) = ssumask(ji,jj) * ( zti(ji,jj,jn) - pta(ji,jj,iku,jn) ) 
     381               ELSE                           ! case 2 
     382                  zmaxu = - ze3wu / e3w_n(ji,jj,ikup1) 
     383                  ! interpolated values of tracers 
     384                  zti(ji,jj,jn) = pta(ji,jj,iku,jn) + zmaxu * ( pta(ji,jj,ikup1,jn) - pta(ji,jj,iku,jn) ) 
    293385                  ! gradient of  tracers 
    294                   pgtu(ji,jj,jn) = umask(ji,jj,iku) * ( zti(ji,jj,jn) - pta(ji,jj,iku,jn) ) 
    295                ELSE                           ! case 2 
    296                   zmaxu = -ze3wu / fse3w(ji,jj,iku) 
    297                   ! interpolated values of tracers 
    298                   zti (ji,jj,jn) = pta(ji,jj,iku,jn) + zmaxu * ( pta(ji,jj,ikum1,jn) - pta(ji,jj,iku,jn) ) 
    299                   ! gradient of tracers 
    300                   pgtu(ji,jj,jn) = umask(ji,jj,iku) * ( pta(ji+1,jj,iku,jn) - zti(ji,jj,jn) ) 
     386                  pgtui(ji,jj,jn) = ssumask(ji,jj) * ( pta(ji+1,jj,iku,jn) - zti(ji,jj,jn) ) 
    301387               ENDIF 
    302388               ! 
    303389               ! j- direction 
    304390               IF( ze3wv >= 0._wp ) THEN      ! case 1 
    305                   zmaxv =  ze3wv / fse3w(ji,jj+1,ikv) 
    306                   ! interpolated values of tracers 
    307                   ztj (ji,jj,jn) = pta(ji,jj+1,ikv,jn) + zmaxv * ( pta(ji,jj+1,ikvm1,jn) - pta(ji,jj+1,ikv,jn) ) 
    308                   ! gradient of tracers 
    309                   pgtv(ji,jj,jn) = vmask(ji,jj,ikv) * ( ztj(ji,jj,jn) - pta(ji,jj,ikv,jn) ) 
    310                ELSE                           ! case 2 
    311                   zmaxv =  -ze3wv / fse3w(ji,jj,ikv) 
    312                   ! interpolated values of tracers 
    313                   ztj (ji,jj,jn) = pta(ji,jj,ikv,jn) + zmaxv * ( pta(ji,jj,ikvm1,jn) - pta(ji,jj,ikv,jn) ) 
    314                   ! gradient of tracers 
    315                   pgtv(ji,jj,jn) = vmask(ji,jj,ikv) * ( pta(ji,jj+1,ikv,jn) - ztj(ji,jj,jn) ) 
    316                ENDIF 
    317             END DO 
    318          END DO 
    319          CALL lbc_lnk( pgtu(:,:,jn), 'U', -1. )   ;   CALL lbc_lnk( pgtv(:,:,jn), 'V', -1. )   ! Lateral boundary cond. 
     391                  zmaxv =  ze3wv / e3w_n(ji,jj+1,ikvp1) 
     392                  ! interpolated values of tracers 
     393                  ztj(ji,jj,jn) = pta(ji,jj+1,ikv,jn) + zmaxv * ( pta(ji,jj+1,ikvp1,jn) - pta(ji,jj+1,ikv,jn) ) 
     394                  ! gradient of tracers 
     395                  pgtvi(ji,jj,jn) = ssvmask(ji,jj) * ( ztj(ji,jj,jn) - pta(ji,jj,ikv,jn) ) 
     396               ELSE                           ! case 2 
     397                  zmaxv =  - ze3wv / e3w_n(ji,jj,ikvp1) 
     398                  ! interpolated values of tracers 
     399                  ztj(ji,jj,jn) = pta(ji,jj,ikv,jn) + zmaxv * ( pta(ji,jj,ikvp1,jn) - pta(ji,jj,ikv,jn) ) 
     400                  ! gradient of tracers 
     401                  pgtvi(ji,jj,jn) = ssvmask(ji,jj) * ( pta(ji,jj+1,ikv,jn) - ztj(ji,jj,jn) ) 
     402               ENDIF 
     403 
     404            END DO 
     405         END DO 
     406         CALL lbc_lnk( pgtui(:,:,jn), 'U', -1. ); CALL lbc_lnk( pgtvi(:,:,jn), 'V', -1. )   ! Lateral boundary cond. 
    320407         ! 
    321408      END DO 
    … …  
    323410      IF( PRESENT( prd ) ) THEN    !==  horizontal derivative of density anomalies (rd)  ==!    (optional part) 
    324411         ! 
    325          pgru  (:,:)=0._wp   ;   pgrv  (:,:) = 0._wp 
    326          pgzu  (:,:)=0._wp   ;   pgzv  (:,:) = 0._wp  
    327          pmru  (:,:)=0._wp   ;   pmru  (:,:) = 0._wp  
    328          pge3ru(:,:)=0._wp   ;   pge3rv(:,:) = 0._wp  
    329          ! 
    330          DO jj = 1, jpjm1                 ! depth of the partial step level 
    331             DO ji = 1, jpim1 
    332                iku = mbku(ji,jj) 
    333                ikv = mbkv(ji,jj) 
    334                ze3wu  = (gdept_0(ji+1,jj,iku) - gdepw_0(ji+1,jj,iku)) - (gdept_0(ji,jj,iku) - gdepw_0(ji,jj,iku)) 
    335                ze3wv  = (gdept_0(ji,jj+1,ikv) - gdepw_0(ji,jj+1,ikv)) - (gdept_0(ji,jj,ikv) - gdepw_0(ji,jj,ikv)) 
    336                ! 
    337                IF( ze3wu >= 0._wp ) THEN   ;   zhi(ji,jj) = fsdept(ji+1,jj,iku) - ze3wu     ! i-direction: case 1 
    338                ELSE                        ;   zhi(ji,jj) = fsdept(ji  ,jj,iku) + ze3wu    ! -     -      case 2 
    339                ENDIF 
    340                IF( ze3wv >= 0._wp ) THEN   ;   zhj(ji,jj) = fsdept(ji,jj+1,ikv) - ze3wv    ! j-direction: case 1 
    341                ELSE                        ;   zhj(ji,jj) = fsdept(ji,jj  ,ikv) + ze3wv    ! -     -      case 2 
    342                ENDIF 
     412         pgrui(:,:)  =0.0_wp; pgrvi(:,:)  =0.0_wp; 
     413         DO jj = 1, jpjm1 
     414            DO ji = 1, jpim1 
     415 
     416               iku = miku(ji,jj) 
     417               ikv = mikv(ji,jj) 
     418               ze3wu  =  gdept_n(ji,jj,iku) - gdept_n(ji+1,jj,iku) 
     419               ze3wv  =  gdept_n(ji,jj,ikv) - gdept_n(ji,jj+1,ikv)  
     420               ! 
     421               IF( ze3wu >= 0._wp ) THEN   ;   zhi(ji,jj) = gdept_n(ji  ,jj,iku)    ! i-direction: case 1 
     422               ELSE                        ;   zhi(ji,jj) = gdept_n(ji+1,jj,iku)    ! -     -      case 2 
     423               ENDIF 
     424 
     425               IF( ze3wv >= 0._wp ) THEN   ;   zhj(ji,jj) = gdept_n(ji,jj  ,ikv)    ! j-direction: case 1 
     426               ELSE                        ;   zhj(ji,jj) = gdept_n(ji,jj+1,ikv)    ! -     -      case 2 
     427               ENDIF 
     428 
    343429            END DO 
    344430         END DO 
    … …  
    346432         CALL eos( zti, zhi, zri )        ! interpolated density from zti, ztj  
    347433         CALL eos( ztj, zhj, zrj )        ! at the partial step depth output in  zri, zrj  
    348  
     434         ! 
    349435         DO jj = 1, jpjm1                 ! Gradient of density at the last level  
    350436            DO ji = 1, jpim1 
    351                iku = mbku(ji,jj) ; ikum1 = MAX( iku - 1 , 1 )    ! last and before last ocean level at u- & v-points 
    352                ikv = mbkv(ji,jj) ; ikvm1 = MAX( ikv - 1 , 1 )    ! last and before last ocean level at u- & v-points 
    353                ze3wu  = (gdept_0(ji+1,jj,iku) - gdepw_0(ji+1,jj,iku)) - (gdept_0(ji,jj,iku) - gdepw_0(ji,jj,iku)) 
    354                ze3wv  = (gdept_0(ji,jj+1,ikv) - gdepw_0(ji,jj+1,ikv)) - (gdept_0(ji,jj,ikv) - gdepw_0(ji,jj,ikv)) 
    355                IF( ze3wu >= 0._wp ) THEN  
    356                   pgzu(ji,jj) = (fsde3w(ji+1,jj,iku) - ze3wu) - fsde3w(ji,jj,iku) 
    357                   pgru(ji,jj) = umask(ji,jj,iku) * ( zri(ji  ,jj) - prd(ji,jj,iku) )   ! i: 1 
    358                   pmru(ji,jj) = umask(ji,jj,iku) * ( zri(ji  ,jj) + prd(ji,jj,iku) )   ! i: 1  
    359                   pge3ru(ji,jj) = umask(ji,jj,iku)                                                                  & 
    360                                 * ( (fse3w(ji+1,jj,iku) - ze3wu )* ( zri(ji  ,jj    ) + prd(ji+1,jj,ikum1) + 2._wp) & 
    361                                    - fse3w(ji  ,jj,iku)          * ( prd(ji  ,jj,iku) + prd(ji  ,jj,ikum1) + 2._wp) )  ! j: 2 
    362                ELSE   
    363                   pgzu(ji,jj) = fsde3w(ji+1,jj,iku) - (fsde3w(ji,jj,iku) + ze3wu) 
    364                   pgru(ji,jj) = umask(ji,jj,iku) * ( prd(ji+1,jj,iku) - zri(ji,jj) )   ! i: 2 
    365                   pmru(ji,jj) = umask(ji,jj,iku) * ( prd(ji+1,jj,iku) + zri(ji,jj) )   ! i: 2 
    366                   pge3ru(ji,jj) = umask(ji,jj,iku)                                                                  & 
    367                                 * (  fse3w(ji+1,jj,iku)          * ( prd(ji+1,jj,iku) + prd(ji+1,jj,ikum1) + 2._wp) & 
    368                                    -(fse3w(ji  ,jj,iku) + ze3wu) * ( zri(ji  ,jj    ) + prd(ji  ,jj,ikum1) + 2._wp) )  ! j: 2 
    369                ENDIF 
    370                IF( ze3wv >= 0._wp ) THEN 
    371                   pgzv(ji,jj) = (fsde3w(ji,jj+1,ikv) - ze3wv) - fsde3w(ji,jj,ikv)  
    372                   pgrv(ji,jj) = vmask(ji,jj,ikv) * ( zrj(ji,jj  ) - prd(ji,jj,ikv) )   ! j: 1 
    373                   pmrv(ji,jj) = vmask(ji,jj,ikv) * ( zrj(ji,jj  ) + prd(ji,jj,ikv) )   ! j: 1 
    374                   pge3rv(ji,jj) = vmask(ji,jj,ikv)                                                                  & 
    375                                 * ( (fse3w(ji,jj+1,ikv) - ze3wv )* ( zrj(ji,jj      ) + prd(ji,jj+1,ikvm1) + 2._wp) & 
    376                                    - fse3w(ji,jj  ,ikv)          * ( prd(ji,jj  ,ikv) + prd(ji,jj  ,ikvm1) + 2._wp) )  ! j: 2 
    377                ELSE  
    378                   pgzv(ji,jj) = fsde3w(ji,jj+1,ikv) - (fsde3w(ji,jj,ikv) + ze3wv) 
    379                   pgrv(ji,jj) = vmask(ji,jj,ikv) * ( prd(ji,jj+1,ikv) - zrj(ji,jj) )   ! j: 2 
    380                   pmrv(ji,jj) = vmask(ji,jj,ikv) * ( prd(ji,jj+1,ikv) + zrj(ji,jj) )   ! j: 2 
    381                   pge3rv(ji,jj) = vmask(ji,jj,ikv)                                                                  & 
    382                                 * (  fse3w(ji,jj+1,ikv)          * ( prd(ji,jj+1,ikv) + prd(ji,jj+1,ikvm1) + 2._wp) & 
    383                                    -(fse3w(ji,jj  ,ikv) + ze3wv) * ( zrj(ji,jj      ) + prd(ji,jj  ,ikvm1) + 2._wp) )  ! j: 2 
    384                ENDIF 
    385             END DO 
    386          END DO 
    387          CALL lbc_lnk( pgru   , 'U', -1. )   ;   CALL lbc_lnk( pgrv   , 'V', -1. )   ! Lateral boundary conditions 
    388          CALL lbc_lnk( pmru   , 'U',  1. )   ;   CALL lbc_lnk( pmrv   , 'V',  1. )   ! Lateral boundary conditions 
    389          CALL lbc_lnk( pgzu   , 'U', -1. )   ;   CALL lbc_lnk( pgzv   , 'V', -1. )   ! Lateral boundary conditions 
    390          CALL lbc_lnk( pge3ru , 'U', -1. )   ;   CALL lbc_lnk( pge3rv , 'V', -1. )   ! Lateral boundary conditions 
    391          ! 
    392       END IF 
    393       ! 
    394       !     !==  (ISH)  compute grui and gruvi  ==! 
    395       ! 
    396       DO jn = 1, kjpt      !==   Interpolation of tracers at the last ocean level   ==!            ! 
    397          DO jj = 1, jpjm1 
    398             DO ji = 1, jpim1 
    399                iku = miku(ji,jj)   ;  ikup1 = miku(ji,jj) + 1 
    400                ikv = mikv(ji,jj)   ;  ikvp1 = mikv(ji,jj) + 1 
    401                ! 
    402                ! (ISF) case partial step top and bottom in adjacent cell in vertical 
    403                ! cannot used e3w because if 2 cell water column, we have ps at top and bottom 
    404                ! in this case e3w(i,j) - e3w(i,j+1) is not the distance between Tj~ and Tj 
    405                ! the only common depth between cells (i,j) and (i,j+1) is gdepw_0 
    406                ze3wu  = (gdepw_0(ji+1,jj,iku+1) - gdept_0(ji+1,jj,iku)) - (gdepw_0(ji,jj,iku+1) - gdept_0(ji,jj,iku))  
    407                ze3wv  = (gdepw_0(ji,jj+1,ikv+1) - gdept_0(ji,jj+1,ikv)) - (gdepw_0(ji,jj,ikv+1) - gdept_0(ji,jj,ikv)) 
    408                ! i- direction 
    409                IF( ze3wu >= 0._wp ) THEN      ! case 1 
    410                   zmaxu = ze3wu / fse3w(ji+1,jj,iku+1) 
    411                   ! interpolated values of tracers 
    412                   zti(ji,jj,jn) = pta(ji+1,jj,iku,jn) + zmaxu * ( pta(ji+1,jj,iku+1,jn) - pta(ji+1,jj,iku,jn) ) 
    413                   ! gradient of tracers 
    414                   pgtui(ji,jj,jn) = umask(ji,jj,iku) * ( zti(ji,jj,jn) - pta(ji,jj,iku,jn) ) 
    415                ELSE                           ! case 2 
    416                   zmaxu = - ze3wu / fse3w(ji,jj,iku+1) 
    417                   ! interpolated values of tracers 
    418                   zti(ji,jj,jn) = pta(ji,jj,iku,jn) + zmaxu * ( pta(ji,jj,iku+1,jn) - pta(ji,jj,iku,jn) ) 
    419                   ! gradient of  tracers 
    420                   pgtui(ji,jj,jn) = umask(ji,jj,iku) * ( pta(ji+1,jj,iku,jn) - zti(ji,jj,jn) ) 
    421                ENDIF 
    422                ! 
    423                ! j- direction 
    424                IF( ze3wv >= 0._wp ) THEN      ! case 1 
    425                   zmaxv =  ze3wv / fse3w(ji,jj+1,ikv+1) 
    426                   ! interpolated values of tracers 
    427                   ztj(ji,jj,jn) = pta(ji,jj+1,ikv,jn) + zmaxv * ( pta(ji,jj+1,ikv+1,jn) - pta(ji,jj+1,ikv,jn) ) 
    428                   ! gradient of tracers 
    429                   pgtvi(ji,jj,jn) = vmask(ji,jj,ikv) * ( ztj(ji,jj,jn) - pta(ji,jj,ikv,jn) ) 
    430                ELSE                           ! case 2 
    431                   zmaxv =  - ze3wv / fse3w(ji,jj,ikv+1) 
    432                   ! interpolated values of tracers 
    433                   ztj(ji,jj,jn) = pta(ji,jj,ikv,jn) + zmaxv * ( pta(ji,jj,ikv+1,jn) - pta(ji,jj,ikv,jn) ) 
    434                   ! gradient of tracers 
    435                   pgtvi(ji,jj,jn) = vmask(ji,jj,ikv) * ( pta(ji,jj+1,ikv,jn) - ztj(ji,jj,jn) ) 
    436                ENDIF 
    437             END DO!! 
    438          END DO!! 
    439          CALL lbc_lnk( pgtui(:,:,jn), 'U', -1. )   ;   CALL lbc_lnk( pgtvi(:,:,jn), 'V', -1. )   ! Lateral boundary cond. 
    440          ! 
    441       END DO 
    442  
    443       IF( PRESENT( prd ) ) THEN    !==  horizontal derivative of density anomalies (rd)  ==!    (optional part) 
    444          ! 
    445          pgrui(:,:)  =0.0_wp ; pgrvi(:,:)  =0.0_wp ; 
    446          pgzui(:,:)  =0.0_wp ; pgzvi(:,:)  =0.0_wp ; 
    447          pmrui(:,:)  =0.0_wp ; pmrui(:,:)  =0.0_wp ; 
    448          pge3rui(:,:)=0.0_wp ; pge3rvi(:,:)=0.0_wp ; 
    449          ! 
    450          DO jj = 1, jpjm1        ! depth of the partial step level 
    451             DO ji = 1, jpim1 
    452                iku = miku(ji,jj) 
    453                ikv = mikv(ji,jj) 
    454                ze3wu  = (gdepw_0(ji+1,jj,iku+1) - gdept_0(ji+1,jj,iku)) - (gdepw_0(ji,jj,iku+1) - gdept_0(ji,jj,iku)) 
    455                ze3wv  = (gdepw_0(ji,jj+1,ikv+1) - gdept_0(ji,jj+1,ikv)) - (gdepw_0(ji,jj,ikv+1) - gdept_0(ji,jj,ikv)) 
    456                ! 
    457                IF( ze3wu >= 0._wp ) THEN   ;   zhi(ji,jj) = fsdept(ji+1,jj,iku) + ze3wu    ! i-direction: case 1 
    458                ELSE                        ;   zhi(ji,jj) = fsdept(ji  ,jj,iku) - ze3wu    ! -     -      case 2 
    459                ENDIF 
    460                IF( ze3wv >= 0._wp ) THEN   ;   zhj(ji,jj) = fsdept(ji,jj+1,ikv) + ze3wv    ! j-direction: case 1 
    461                ELSE                        ;   zhj(ji,jj) = fsdept(ji,jj  ,ikv) - ze3wv    ! -     -      case 2 
    462                ENDIF 
    463             END DO 
    464          END DO 
    465          ! 
    466          CALL eos( zti, zhi, zri )        ! interpolated density from zti, ztj  
    467          CALL eos( ztj, zhj, zrj )        ! at the partial step depth output in  zri, zrj  
    468          ! 
    469          DO jj = 1, jpjm1                 ! Gradient of density at the last level  
    470             DO ji = 1, jpim1 
    471                iku = miku(ji,jj) ; ikup1 = miku(ji,jj) + 1 
    472                ikv = mikv(ji,jj) ; ikvp1 = mikv(ji,jj) + 1 
    473                ze3wu  = (gdepw_0(ji+1,jj,iku+1) - gdept_0(ji+1,jj,iku)) - (gdepw_0(ji,jj,iku+1) - gdept_0(ji,jj,iku)) 
    474                ze3wv  = (gdepw_0(ji,jj+1,ikv+1) - gdept_0(ji,jj+1,ikv)) - (gdepw_0(ji,jj,ikv+1) - gdept_0(ji,jj,ikv)) 
    475                IF( ze3wu >= 0._wp ) THEN 
    476                  pgzui  (ji,jj) = (fsde3w(ji+1,jj,iku) + ze3wu) - fsde3w(ji,jj,iku) 
    477                  pgrui  (ji,jj) = umask(ji,jj,iku)   * ( zri(ji,jj) - prd(ji,jj,iku) )          ! i: 1 
    478                  pmrui  (ji,jj) = umask(ji,jj,iku)   * ( zri(ji,jj) + prd(ji,jj,iku) )          ! i: 1  
    479                  pge3rui(ji,jj) = umask(ji,jj,iku+1)                                                                  & 
    480                     &           * ( (fse3w(ji+1,jj,iku+1) - ze3wu) * (zri(ji,jj    ) + prd(ji+1,jj,iku+1) + 2._wp)   & 
    481                     &              - fse3w(ji  ,jj,iku+1)          * (prd(ji,jj,iku) + prd(ji  ,jj,iku+1) + 2._wp)   ) ! i: 1 
    482                ELSE 
    483                  pgzui  (ji,jj) = fsde3w(ji+1,jj,iku) - (fsde3w(ji,jj,iku) - ze3wu) 
    484                  pgrui  (ji,jj) = umask(ji,jj,iku)   * ( prd(ji+1,jj,iku) - zri(ji,jj) )      ! i: 2 
    485                  pmrui  (ji,jj) = umask(ji,jj,iku)   * ( prd(ji+1,jj,iku) + zri(ji,jj) )      ! i: 2 
    486                  pge3rui(ji,jj) = umask(ji,jj,iku+1)                                                                   & 
    487                     &           * (  fse3w(ji+1,jj,iku+1)          * (prd(ji+1,jj,iku) + prd(ji+1,jj,iku+1) + 2._wp)  & 
    488                     &              -(fse3w(ji  ,jj,iku+1) + ze3wu) * (zri(ji,jj      ) + prd(ji  ,jj,iku+1) + 2._wp)  )     ! i: 2 
    489                ENDIF 
    490                IF( ze3wv >= 0._wp ) THEN 
    491                  pgzvi  (ji,jj) = (fsde3w(ji,jj+1,ikv) + ze3wv) - fsde3w(ji,jj,ikv)  
    492                  pgrvi  (ji,jj) = vmask(ji,jj,ikv)   * ( zrj(ji,jj  ) - prd(ji,jj,ikv) )        ! j: 1 
    493                  pmrvi  (ji,jj) = vmask(ji,jj,ikv)   * ( zrj(ji,jj  ) + prd(ji,jj,ikv) )        ! j: 1 
    494                  pge3rvi(ji,jj) = vmask(ji,jj,ikv+1)                                                                  &  
    495                      &           * ( (fse3w(ji,jj+1,ikv+1) - ze3wv) * ( zrj(ji,jj    ) + prd(ji,jj+1,ikv+1) + 2._wp)  & 
    496                      &                - fse3w(ji,jj  ,ikv+1)          * ( prd(ji,jj,ikv) + prd(ji,jj  ,ikv+1) + 2._wp)  ) ! j: 1 
    497                                   ! + 2 due to the formulation in density and not in anomalie in hpg sco 
    498                ELSE 
    499                  pgzvi  (ji,jj) = fsde3w(ji,jj+1,ikv) - (fsde3w(ji,jj,ikv) - ze3wv) 
    500                  pgrvi  (ji,jj) = vmask(ji,jj,ikv)   * ( prd(ji,jj+1,ikv) - zrj(ji,jj) )     ! j: 2 
    501                  pmrvi  (ji,jj) = vmask(ji,jj,ikv)   * ( prd(ji,jj+1,ikv) + zrj(ji,jj) )     ! j: 2 
    502                  pge3rvi(ji,jj) = vmask(ji,jj,ikv+1)                                                                   & 
    503                     &           * (  fse3w(ji,jj+1,ikv+1)          * ( prd(ji,jj+1,ikv) + prd(ji,jj+1,ikv+1) + 2._wp) & 
    504                     &              -(fse3w(ji,jj  ,ikv+1) + ze3wv) * ( zrj(ji,jj      ) + prd(ji,jj  ,ikv+1) + 2._wp) )  ! j: 2 
    505                ENDIF 
    506             END DO 
    507          END DO 
    508          CALL lbc_lnk( pgrui   , 'U', -1. )   ;   CALL lbc_lnk( pgrvi   , 'V', -1. )   ! Lateral boundary conditions 
    509          CALL lbc_lnk( pmrui   , 'U',  1. )   ;   CALL lbc_lnk( pmrvi   , 'V',  1. )   ! Lateral boundary conditions 
    510          CALL lbc_lnk( pgzui   , 'U', -1. )   ;   CALL lbc_lnk( pgzvi   , 'V', -1. )   ! Lateral boundary conditions 
    511          CALL lbc_lnk( pge3rui , 'U', -1. )   ;   CALL lbc_lnk( pge3rvi , 'V', -1. )   ! Lateral boundary conditions 
     437               iku = miku(ji,jj)  
     438               ikv = mikv(ji,jj)  
     439               ze3wu  =  gdept_n(ji,jj,iku) - gdept_n(ji+1,jj,iku) 
     440               ze3wv  =  gdept_n(ji,jj,ikv) - gdept_n(ji,jj+1,ikv)  
     441 
     442               IF( ze3wu >= 0._wp ) THEN ; pgrui(ji,jj) = ssumask(ji,jj) * ( zri(ji  ,jj      ) - prd(ji,jj,iku) ) ! i: 1 
     443               ELSE                      ; pgrui(ji,jj) = ssumask(ji,jj) * ( prd(ji+1,jj  ,iku) - zri(ji,jj    ) ) ! i: 2 
     444               ENDIF 
     445               IF( ze3wv >= 0._wp ) THEN ; pgrvi(ji,jj) = ssvmask(ji,jj) * ( zrj(ji  ,jj      ) - prd(ji,jj,ikv) ) ! j: 1 
     446               ELSE                      ; pgrvi(ji,jj) = ssvmask(ji,jj) * ( prd(ji  ,jj+1,ikv) - zrj(ji,jj    ) ) ! j: 2 
     447               ENDIF 
     448 
     449            END DO 
     450         END DO 
     451         CALL lbc_lnk( pgrui   , 'U', -1. ); CALL lbc_lnk( pgrvi   , 'V', -1. )   ! Lateral boundary conditions 
    512452         ! 
    513453      END IF   
Note: See TracChangeset for help on using the changeset viewer.