Index: branches/2014/dev_r4650_UKMO3_masked_damping/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_muscl.F90
===================================================================
--- branches/2014/dev_r4650_UKMO3_masked_damping/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_muscl.F90	(revision 4668)
+++ branches/2014/dev_r4650_UKMO3_masked_damping/NEMOGCM/NEMO/OPA_SRC/TRA/traadv_muscl.F90	(revision 5086)
@@ -16,29 +16,28 @@
    !!----------------------------------------------------------------------
    USE oce            ! ocean dynamics and active tracers
+   USE trc_oce        ! share passive tracers/Ocean variables
    USE dom_oce        ! ocean space and time domain
-   USE trdmod_oce     ! tracers trends 
-   USE trdtra         ! tracers trends 
-   USE in_out_manager ! I/O manager
+   USE trd_oce        ! trends: ocean variables
+   USE trdtra         ! tracers trends manager
    USE dynspg_oce     ! choice/control of key cpp for surface pressure gradient
-   USE trabbl         ! tracers: bottom boundary layer
-   USE lib_mpp        ! distribued memory computing
-   USE lbclnk         ! ocean lateral boundary condition (or mpp link) 
+   USE sbcrnf          ! river runoffs
    USE diaptr         ! poleward transport diagnostics
-   USE trc_oce        ! share passive tracers/Ocean variables
+   !
    USE wrk_nemo       ! Memory Allocation
    USE timing         ! Timing
    USE lib_fortran    ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)  
-   USE eosbn2          ! equation of state
-   USE sbcrnf          ! river runoffs
+   USE in_out_manager ! I/O manager
+   USE lib_mpp        ! distribued memory computing
+   USE lbclnk         ! ocean lateral boundary condition (or mpp link) 
 
    IMPLICIT NONE
    PRIVATE
 
-   PUBLIC   tra_adv_muscl  ! routine called by step.F90
-
-   LOGICAL  :: l_trd                        ! flag to compute trends
-   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: upsmsk !: mixed upstream/centered scheme near some straits
-   !                                                             !  and in closed seas (orca 2 and 4 configurations)
-   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: xind         !: mixed upstream/centered index
+   PUBLIC   tra_adv_muscl   ! routine called by traadv.F90
+   
+   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   upsmsk   !: mixed upstream/centered scheme near some straits
+   !                                                           !  and in closed seas (orca 2 and 4 configurations)
+   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   xind     !: mixed upstream/centered index
+   
    !! * Substitutions
 #  include "domzgr_substitute.h90"
@@ -51,6 +50,6 @@
 CONTAINS
 
-   SUBROUTINE tra_adv_muscl( kt, kit000, cdtype, p2dt, pun, pvn, pwn, &
-      &                                        ptb, pta, kjpt, ld_msc_ups )
+   SUBROUTINE tra_adv_muscl( kt, kit000, cdtype, p2dt, pun, pvn, pwn,   &
+      &                                                ptb, pta, kjpt, ld_msc_ups )
       !!----------------------------------------------------------------------
       !!                    ***  ROUTINE tra_adv_muscl  ***
@@ -68,6 +67,4 @@
       !!              IPSL, Sept. 2000 (http://www.lodyc.jussieu.fr/opa)
       !!----------------------------------------------------------------------
-      USE oce     , ONLY:   zwx   => ua    , zwy   => va          ! (ua,va) used as workspace
-      !
       INTEGER                              , INTENT(in   ) ::   kt              ! ocean time-step index
       INTEGER                              , INTENT(in   ) ::   kit000          ! first time step index
@@ -79,19 +76,18 @@
       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in   ) ::   ptb             ! before tracer field
       REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) ::   pta             ! tracer trend 
-
-      !
-      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices
+      !
+      INTEGER  ::   ji, jj, jk, jn            ! dummy loop indices
+      INTEGER  ::   ierr                      ! local integer
       REAL(wp) ::   zu, z0u, zzwx, zw         ! local scalars
       REAL(wp) ::   zv, z0v, zzwy, z0w        !   -      -
       REAL(wp) ::   ztra, zbtr, zdt, zalpha   !   -      -
-      REAL(wp), POINTER, DIMENSION(:,:,:) :: zslpx, zslpy
-      INTEGER  ::   ierr
+      REAL(wp), POINTER, DIMENSION(:,:,:) ::   zslpx, zslpy   ! 3D workspace
+      REAL(wp), POINTER, DIMENSION(:,:,:) ::   zwx  , zwy     ! -      - 
       !!----------------------------------------------------------------------
       !
       IF( nn_timing == 1 )  CALL timing_start('tra_adv_muscl')
       !
-      CALL wrk_alloc( jpi, jpj, jpk, zslpx, zslpy )
-      !
-
+      CALL wrk_alloc( jpi, jpj, jpk, zslpx, zslpy, zwx, zwy )
+      !
       IF( kt == kit000 )  THEN
          IF(lwp) WRITE(numout,*)
@@ -117,25 +113,20 @@
 
          !
-         ! Upstream / centered scheme indicator
+         ! Upstream / MUSCL scheme indicator
          ! ------------------------------------
+!!gm  useless
          xind(:,:,:) = 1._wp                             ! set equal to 1 where up-stream is not needed
+!!gm
          !
          IF( ld_msc_ups )  THEN
-            DO jk = 1, jpk
-               DO jj = 1, jpj
-                  DO ji = 1, jpi
-                     xind(ji,jj,jk) = 1  - MAX (           &
-                        rnfmsk(ji,jj) * rnfmsk_z(jk),      &  ! near runoff mouths (& closed sea outflows)
-                        upsmsk(ji,jj) ) * tmask(ji,jj,jk)     ! some of some straits
-                  END DO
-               END DO
+            DO jk = 1, jpkm1
+               xind(:,:,jk) = 1._wp                              &                 ! =>1 where up-stream is not needed
+                  &         - MAX ( rnfmsk(:,:) * rnfmsk_z(jk),  &                 ! =>0 near runoff mouths (& closed sea outflows)
+                  &                 upsmsk(:,:)                ) * tmask(:,:,jk)   ! =>0 near some straits
             END DO
          ENDIF 
          !
       ENDIF 
-      !
-      l_trd = .FALSE.
-      IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR. ( cdtype == 'TRC' .AND. l_trdtrc ) ) l_trd = .TRUE.
-      
+      !      
       !                                                     ! ===========
       DO jn = 1, kjpt                                       ! tracer loop
@@ -192,6 +183,6 @@
                   zalpha = 0.5 - z0u
                   zu  = z0u - 0.5 * pun(ji,jj,jk) * zdt / ( e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,jk) )
-                  zzwx = ptb(ji+1,jj,jk,jn) + xind(ji,jj,jk) * (zu * zslpx(ji+1,jj,jk))
-                  zzwy = ptb(ji  ,jj,jk,jn) + xind(ji,jj,jk) * (zu * zslpx(ji  ,jj,jk))
+                  zzwx = ptb(ji+1,jj,jk,jn) + xind(ji,jj,jk) * zu * zslpx(ji+1,jj,jk)
+                  zzwy = ptb(ji  ,jj,jk,jn) + xind(ji,jj,jk) * zu * zslpx(ji  ,jj,jk)
                   zwx(ji,jj,jk) = pun(ji,jj,jk) * ( zalpha * zzwx + (1.-zalpha) * zzwy )
                   !
@@ -199,6 +190,6 @@
                   zalpha = 0.5 - z0v
                   zv  = z0v - 0.5 * pvn(ji,jj,jk) * zdt / ( e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,jk) )
-                  zzwx = ptb(ji,jj+1,jk,jn) + xind(ji,jj,jk) * (zv * zslpy(ji,jj+1,jk))
-                  zzwy = ptb(ji,jj  ,jk,jn) + xind(ji,jj,jk) * (zv * zslpy(ji,jj  ,jk))
+                  zzwx = ptb(ji,jj+1,jk,jn) + xind(ji,jj,jk) * zv * zslpy(ji,jj+1,jk)
+                  zzwy = ptb(ji,jj  ,jk,jn) + xind(ji,jj,jk) * zv * zslpy(ji,jj  ,jk)
                   zwy(ji,jj,jk) = pvn(ji,jj,jk) * ( zalpha * zzwx + (1.-zalpha) * zzwy )
                END DO
@@ -222,7 +213,8 @@
          END DO        
          !                                 ! trend diagnostics (contribution of upstream fluxes)
-         IF( l_trd )  THEN
-            CALL trd_tra( kt, cdtype, jn, jptra_trd_xad, zwx, pun, ptb(:,:,:,jn) )
-            CALL trd_tra( kt, cdtype, jn, jptra_trd_yad, zwy, pvn, ptb(:,:,:,jn) )
+         IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR.   &
+            &( cdtype == 'TRC' .AND. l_trdtrc )      )  THEN
+            CALL trd_tra( kt, cdtype, jn, jptra_xad, zwx, pun, ptb(:,:,:,jn) )
+            CALL trd_tra( kt, cdtype, jn, jptra_yad, zwy, pvn, ptb(:,:,:,jn) )
          END IF
          !                                 ! "Poleward" heat and salt transports (contribution of upstream fluxes)
@@ -274,6 +266,6 @@
                   zalpha = 0.5 + z0w
                   zw  = z0w - 0.5 * pwn(ji,jj,jk+1) * zdt * zbtr 
-                  zzwx = ptb(ji,jj,jk+1,jn) + xind(ji,jj,jk) * (zw * zslpx(ji,jj,jk+1))
-                  zzwy = ptb(ji,jj,jk  ,jn) + xind(ji,jj,jk) * (zw * zslpx(ji,jj,jk  ))
+                  zzwx = ptb(ji,jj,jk+1,jn) + xind(ji,jj,jk) * zw * zslpx(ji,jj,jk+1)
+                  zzwy = ptb(ji,jj,jk  ,jn) + xind(ji,jj,jk) * zw * zslpx(ji,jj,jk  )
                   zwx(ji,jj,jk+1) = pwn(ji,jj,jk+1) * ( zalpha * zzwx + (1.-zalpha) * zzwy )
                END DO 
@@ -281,9 +273,8 @@
          END DO
 
-         ! Compute & add the vertical advective trend
-         DO jk = 1, jpkm1
+         DO jk = 1, jpkm1                    ! Compute & add the vertical advective 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) )
+                  zbtr = 1. / ( e1e2t(ji,jj) * fse3t(ji,jj,jk) )
                   ! vertical advective trends 
                   ztra = - zbtr * ( zwx(ji,jj,jk) - zwx(ji,jj,jk+1) )
@@ -294,9 +285,11 @@
          END DO
          !                                 ! Save the vertical advective trends for diagnostic
-         IF( l_trd ) CALL trd_tra( kt, cdtype, jn, jptra_trd_zad, zwx, pwn, ptb(:,:,:,jn) )
-         !
-      ENDDO
-      !
-      CALL wrk_dealloc( jpi, jpj, jpk, zslpx, zslpy )
+         IF( ( cdtype == 'TRA' .AND. l_trdtra ) .OR.     &
+            &( cdtype == 'TRC' .AND. l_trdtrc )      )   &
+            CALL trd_tra( kt, cdtype, jn, jptra_zad, zwx, pwn, ptb(:,:,:,jn) )
+         !
+      END DO
+      !
+      CALL wrk_dealloc( jpi, jpj, jpk, zslpx, zslpy, zwx, zwy )
       !
       IF( nn_timing == 1 )  CALL timing_stop('tra_adv_muscl')
