Index: branches/2015/dev_r5072_UKMO2_OBS_simplification/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_qck.F90
===================================================================
--- branches/2015/dev_r5072_UKMO2_OBS_simplification/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_qck.F90	(revision 5682)
+++ branches/2015/dev_r5072_UKMO2_OBS_simplification/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_qck.F90	(revision 5974)
@@ -20,5 +20,4 @@
    USE trd_oce         ! trends: ocean variables
    USE trdtra          ! trends manager: tracers 
-   USE dynspg_oce      ! surface pressure gradient variables
    USE diaptr          ! poleward transport diagnostics
    !
@@ -102,4 +101,5 @@
          IF(lwp) WRITE(numout,*)
       ENDIF
+      !
       l_trd = .FALSE.
       IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) )   l_trd = .TRUE.
@@ -130,7 +130,7 @@
       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pta        ! tracer trend 
       !!
-      INTEGER  :: ji, jj, jk, jn   ! dummy loop indices
-      REAL(wp) :: ztra, zbtr, zdir, zdx, zdt, zmsk   ! local scalars
-      REAL(wp), POINTER, DIMENSION(:,:,:) :: zwx, zfu, zfc, zfd
+      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices
+      REAL(wp) ::   ztra, zbtr, zdir, zdx, zdt, zmsk   ! local scalars
+      REAL(wp), POINTER, DIMENSION(:,:,:) ::   zwx, zfu, zfc, zfd
       !----------------------------------------------------------------------
       !
@@ -139,16 +139,13 @@
       DO jn = 1, kjpt                                            ! tracer loop
          !                                                       ! ===========
-         zfu(:,:,:) = 0.0     ;   zfc(:,:,:) = 0.0  
-         zfd(:,:,:) = 0.0     ;   zwx(:,:,:) = 0.0     
-         !                                                  
-         DO jk = 1, jpkm1                                
-            !                                             
-            !--- Computation of the ustream and downstream value of the tracer and the mask
+         zfu(:,:,:) = 0._wp     ;   zfc(:,:,:) = 0._wp 
+         zfd(:,:,:) = 0._wp     ;   zwx(:,:,:) = 0._wp   
+         !
+!!gm why not using a SHIFT instruction...
+         DO jk = 1, jpkm1     !--- Computation of the ustream and downstream value of the tracer and the mask
             DO jj = 2, jpjm1
                DO ji = fs_2, fs_jpim1   ! vector opt.
-                  ! Upstream in the x-direction for the tracer
-                  zfc(ji,jj,jk) = ptb(ji-1,jj,jk,jn)
-                  ! Downstream in the x-direction for the tracer
-                  zfd(ji,jj,jk) = ptb(ji+1,jj,jk,jn)
+                  zfc(ji,jj,jk) = ptb(ji-1,jj,jk,jn)        ! Upstream   in the x-direction for the tracer
+                  zfd(ji,jj,jk) = ptb(ji+1,jj,jk,jn)        ! Downstream in the x-direction for the tracer
                END DO
             END DO
@@ -159,5 +156,4 @@
          ! Horizontal advective fluxes
          ! ---------------------------
-         !
          DO jk = 1, jpkm1                             
             DO jj = 2, jpjm1
@@ -220,5 +216,5 @@
             DO jj = 2, jpjm1
                DO ji = fs_2, fs_jpim1   ! vector opt.  
-                  zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
+                  zbtr = 1. / ( e1e2t(ji,jj) * fse3t(ji,jj,jk) )
                   ! horizontal advective trends
                   ztra = - zbtr * ( zwx(ji,jj,jk) - zwx(ji-1,jj,jk) )
@@ -344,5 +340,5 @@
             DO jj = 2, jpjm1
                DO ji = fs_2, fs_jpim1   ! vector opt.  
-                  zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
+                  zbtr = 1. / ( e1e2t(ji,jj) * fse3t(ji,jj,jk) )
                   ! horizontal advective trends
                   ztra = - zbtr * ( zwy(ji,jj,jk) - zwy(ji,jj-1,jk) )
@@ -380,36 +376,41 @@
       !
       INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices
-      REAL(wp) ::   zbtr , ztra      ! local scalars
       REAL(wp), POINTER, DIMENSION(:,:,:) :: zwz
       !!----------------------------------------------------------------------
       !
-      CALL wrk_alloc( jpi, jpj, jpk, zwz )
+      CALL wrk_alloc( jpi,jpj,jpk,   zwz )
+      !
+      !                          ! surface & bottom values 
+      IF( lk_vvl )   zwz(:,:, 1 ) = 0._wp             ! set to zero one for all
+                     zwz(:,:,jpk) = 0._wp             ! except at the surface in linear free surface
+      !
       !                                                          ! ===========
       DO jn = 1, kjpt                                            ! tracer loop
          !                                                       ! ===========
-         ! 1. Bottom value : flux set to zero
-         zwz(:,:,jpk) = 0.e0             ! Bottom value : flux set to zero
-         !
-         !                                 ! Surface value
-         IF( lk_vvl ) THEN   ;   zwz(:,:, 1 ) = 0.e0                      ! Variable volume : flux set to zero
-         ELSE                ;   zwz(:,:, 1 ) = pwn(:,:,1) * ptn(:,:,1,jn)   ! Constant volume : advective flux through the surface
+         !
+         DO jk = 2, jpkm1                    !* Interior point   (w-masked 2nd order centered flux)
+            DO jj = 2, jpjm1
+               DO ji = fs_2, fs_jpim1   ! vector opt.
+                  zwz(ji,jj,jk) = 0.5 * pwn(ji,jj,jk) * ( ptn(ji,jj,jk-1,jn) + ptn(ji,jj,jk,jn) ) * wmask(ji,jj,jk)
+               END DO
+            END DO
+         END DO
+         IF(.NOT.lk_vvl ) THEN               !* top value   (only in linear free surf. as zwz is multiplied by wmask)
+            IF( ln_isfcav ) THEN                  ! ice-shelf cavities (top of the ocean)
+               DO jj = 1, jpj
+                  DO ji = 1, jpi
+                     zwz(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * ptn(ji,jj,mikt(ji,jj),jn)   ! linear free surface 
+                  END DO
+               END DO   
+            ELSE                                   ! no ice-shelf cavities (only ocean surface)
+               zwz(:,:,1) = pwn(:,:,1) * ptn(:,:,1,jn)
+            ENDIF
          ENDIF
          !
-         DO jk = 2, jpkm1                  ! Interior point: second order centered tracer flux at w-point
+         DO jk = 1, jpkm1          !==  Tracer flux divergence added to the general trend  ==!
             DO jj = 2, jpjm1
                DO ji = fs_2, fs_jpim1   ! vector opt.
-                  zwz(ji,jj,jk) = 0.5 * pwn(ji,jj,jk) * ( ptn(ji,jj,jk-1,jn) + ptn(ji,jj,jk,jn) )
-               END DO
-            END DO
-         END DO
-         !
-         DO jk = 1, jpkm1          !==  Tracer flux divergence added to the general trend  ==!
-            DO jj = 2, jpjm1
-               DO ji = fs_2, fs_jpim1   ! vector opt.
-                  zbtr = 1. / ( e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk) )
-                  ! k- vertical advective trends 
-                  ztra = - zbtr * ( zwz(ji,jj,jk) - zwz(ji,jj,jk+1) ) 
-                  ! added to the general tracer trends
-                  pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ztra
+                  pta(ji,jj,jk,jn) = pta(ji,jj,jk,jn) - ( zwz(ji,jj,jk) - zwz(ji,jj,jk+1) )   &
+                     &                                / ( e1e2t(ji,jj) * fse3t(ji,jj,jk) )
                END DO
             END DO
