Index: branches/2012/dev_v3_4_STABLE_2012/NEMOGCM/NEMO/OPATAM_SRC/DYN/dynspg_flt_tam.F90
===================================================================
--- branches/2012/dev_v3_4_STABLE_2012/NEMOGCM/NEMO/OPATAM_SRC/DYN/dynspg_flt_tam.F90	(revision 4577)
+++ branches/2012/dev_v3_4_STABLE_2012/NEMOGCM/NEMO/OPATAM_SRC/DYN/dynspg_flt_tam.F90	(revision 4578)
@@ -46,4 +46,5 @@
    USE dom_oce
    USE solver
+   USE dynspg_flt
    USE sol_oce
    USE oce_tam
@@ -63,4 +64,5 @@
    USE lib_fortran
    USE timing
+   USE iom
 
 
@@ -149,5 +151,8 @@
          spgu_tl(:,:) = 0.0_wp                    ! surface pressure gradient (i-direction)
          spgv_tl(:,:) = 0.0_wp                    ! surface pressure gradient (j-direction)
-         !CALL solver_init( nit000 )               ! Elliptic solver initialisation
+         ! Reinitialize the solver arrays
+         gcxb_tl(:,:) = 0.e0
+         gcx_tl (:,:) = 0.e0
+         CALL sol_mat( nit000 )
       ENDIF
       ! Local constant initialization
@@ -168,4 +173,5 @@
             END DO
          END DO
+
          DO jk = 1, jpkm1              ! unweighted time stepping
             DO jj = 2, jpjm1
@@ -236,4 +242,5 @@
          END DO
       END DO
+
       ! apply the lateral boundary conditions
       IF( nn_solv == 2 .AND. MAX( jpr2di, jpr2dj ) > 0 ) CALL lbc_lnk_e( gcb_tl, c_solver_pt, 1.0_wp )
@@ -291,4 +298,5 @@
          END DO
       END DO
+
       !
       IF( nn_timing == 1 )  CALL timing_stop('dyn_spg_flt_tan')
@@ -370,26 +378,19 @@
          CALL sol_mat(kt)       ! initialize matrix
 
-      ELSEIF (   neuler == 0 .AND. kt == nit000 + 1 ) THEN
+      ELSEIF (                   kt == nitend     ) THEN
 
          z2dt = 2.0_wp * rdt    ! time step: leap-frog
          CALL sol_mat(kt)       ! reinitialize matrix
 
-      ELSEIF (                   kt == nitend     ) THEN
+      ELSE
 
          z2dt = 2.0_wp * rdt    ! time step: leap-frog
-         CALL sol_mat(kt)       ! reinitialize matrix
-
-      ELSEIF ( neuler /= 0 .AND. kt == nit000     ) THEN
-
-         z2dt = 2.0_wp * rdt    ! time step: leap-frog
-         CALL sol_mat(kt)       ! initialize matrix
-
-      ELSE
-
-         z2dt = 2.0_wp * rdt    ! time step: leap-frog
 
       ENDIF
 
       z2dtg  = grav * z2dt
+
+      ! set to zero free surface specific arrays (they are actually local variables)
+      spgu_ad(:,:) = 0.0_wp    ;      spgv_ad(:,:) = 0.0_wp
 
       ! Add the trends multiplied by z2dt to the after velocity
@@ -546,10 +547,10 @@
                   ua_ad(  ji,jj,jk) = ua_ad(ji,jj,jk) * umask(ji,jj,jk)
                   ub_ad(  ji,jj,jk) = ub_ad(ji,jj,jk) + ua_ad(ji,jj,jk)
-                  spgu_ad(ji,jj   ) = spgu_ad(ji,jj)  + ua_ad(ji,jj,jk) * z2dt
                   ua_ad(  ji,jj,jk) = ua_ad(ji,jj,jk) * z2dt
+                  spgu_ad(ji,jj   ) = spgu_ad(ji,jj)  + ua_ad(ji,jj,jk)
                   va_ad(  ji,jj,jk) = va_ad(ji,jj,jk) * vmask(ji,jj,jk)
                   vb_ad(  ji,jj,jk) = vb_ad(ji,jj,jk) + va_ad(ji,jj,jk)
-                  spgv_ad(ji,jj   ) = spgv_ad(ji,jj)  + va_ad(ji,jj,jk) * z2dt
                   va_ad(  ji,jj,jk) = va_ad(ji,jj,jk) * z2dt
+                  spgv_ad(ji,jj   ) = spgv_ad(ji,jj)  + va_ad(ji,jj,jk)
                END DO
             END DO
@@ -573,4 +574,7 @@
          spgu_ad(:,:) = 0.0_wp                    ! surface pressure gradient (i-direction)
          spgv_ad(:,:) = 0.0_wp                    ! surface pressure gradient (j-direction)
+         ! Reinitialize the solver arrays
+         gcxb_ad(:,:) = 0.e0
+         gcx_ad (:,:) = 0.e0
       ENDIF
       !
@@ -615,8 +619,4 @@
          & zua_tlout,   & ! Tangent output: ua_tl
          & zva_tlout,   & ! Tangent output: va_tl
-         & zub_tlout,   & ! Tangent output: ua_tl
-         & zvb_tlout,   & ! Tangent output: va_tl
-         & zub_adin,    & ! Tangent output: ua_ad
-         & zvb_adin,    & ! Tangent output: va_ad
          & zua_adin,    & ! Adjoint input: ua_ad
          & zva_adin,    & ! Adjoint input: va_ad
@@ -628,13 +628,10 @@
 
       REAL(wp), DIMENSION(:,:), ALLOCATABLE :: &
-         & zgcx_tlin, zgcxb_tlin, zgcb_tlin, zgcx_tlout, zgcxb_tlout, zgcb_tlout,  &
-         & zgcx_adin, zgcxb_adin, zgcb_adin, zgcx_adout, zgcxb_adout, zgcb_adout,  &
-         & zspgu_tlout, zspgv_tlout, zspgu_adin, zspgv_adin
+         & zgcx_tlin, zgcxb_tlin, zgcx_tlout, zgcxb_tlout,  &
+         & zgcx_adin, zgcxb_adin, zgcx_adout, zgcxb_adout
 
       REAL(KIND=wp), DIMENSION(:,:), ALLOCATABLE :: &
          & zsshn_tlin, & ! Tangent input: sshn_tl
          & zsshn_adout,& ! Adjoint output: sshn_ad
-         & zemp_tlin,  & ! Tangent input: emp_tl
-         & zemp_adout, & ! Adjoint output: emp_ad
          & znssh         ! 2D random field for SSH
       REAL(wp) :: &
@@ -654,5 +651,5 @@
       INTEGER ::             &
          & jpert
-      INTEGER, PARAMETER :: jpertmax = 7
+      INTEGER, PARAMETER :: jpertmax = 6
 
       ! Allocate memory
@@ -676,18 +673,9 @@
          & zsshn_tlin(jpi,jpj), &
          & zsshn_adout(jpi,jpj),&
-         & zemp_tlin(jpi,jpj),  &
-         & zemp_adout(jpi,jpj), &
          & znssh(jpi,jpj)       &
          & )
 
       ALLOCATE( zgcx_tlin (jpi,jpj), zgcx_tlout (jpi,jpj), zgcx_adin (jpi,jpj), zgcx_adout (jpi,jpj),  &
-                zgcxb_tlin(jpi,jpj), zgcxb_tlout(jpi,jpj), zgcxb_adin(jpi,jpj), zgcxb_adout(jpi,jpj),  &
-                zgcb_tlin (jpi,jpj), zgcb_tlout (jpi,jpj), zgcb_adin (jpi,jpj), zgcb_adout (jpi,jpj)   &
-               & )
-
-      ALLOCATE ( zub_tlout(jpi,jpj,jpk), zvb_tlout(jpi,jpj,jpk),  &
-                 zub_adin (jpi,jpj,jpk), zvb_adin (jpi,jpj,jpk)  )
-
-      ALLOCATE( zspgu_tlout (jpi,jpj), zspgv_tlout (jpi,jpj), zspgu_adin (jpi,jpj), zspgv_adin (jpi,jpj))
+                zgcxb_tlin(jpi,jpj), zgcxb_tlout(jpi,jpj), zgcxb_adin(jpi,jpj), zgcxb_adout(jpi,jpj)   )
 
       !=========================================================================
@@ -698,6 +686,6 @@
       ! Test for time steps nit000 and nit000 + 1 (the matrix changes)
 
-      DO jstp = nit000, nit000 + 1
-         DO jpert = 1, jpertmax
+      DO jstp = nit000, nitend, nitend-nit000
+         DO jpert = jpertmax, jpertmax
             istp = jstp
 
@@ -720,15 +708,8 @@
 
             zsshn_tlin (:,:) = 0.0_wp
-            zemp_tlin  (:,:) = 0.0_wp
             zsshn_adout(:,:) = 0.0_wp
-            zemp_adout (:,:) = 0.0_wp
-            zspgu_adin (:,:) = 0.0_wp
-            zspgv_adin (:,:) = 0.0_wp
-            zspgu_tlout(:,:) = 0.0_wp
-            zspgv_tlout(:,:) = 0.0_wp
 
             zgcx_tlout (:,:) = 0.0_wp ; zgcx_adin (:,:) = 0.0_wp ; zgcx_adout (:,:) = 0.0_wp
             zgcxb_tlout(:,:) = 0.0_wp ; zgcxb_adin(:,:) = 0.0_wp ; zgcxb_adout(:,:) = 0.0_wp
-            zgcb_tlout (:,:) = 0.0_wp ; zgcb_adin (:,:) = 0.0_wp ; zgcb_adout (:,:) = 0.0_wp
 
             ub_tl(:,:,:) = 0.0_wp
@@ -737,6 +718,4 @@
             va_tl(:,:,:) = 0.0_wp
             sshn_tl(:,:) = 0.0_wp
-            emp_tl(:,:)  = 0.0_wp
-            gcb_tl(:,:)  = 0.0_wp
             gcx_tl(:,:)  = 0.0_wp
             gcxb_tl(:,:) = 0.0_wp
@@ -748,5 +727,4 @@
             va_ad(:,:,:) = 0.0_wp
             sshn_ad(:,:) = 0.0_wp
-            emp_ad(:,:)  = 0.0_wp
             gcb_ad(:,:)  = 0.0_wp
             gcx_ad(:,:)  = 0.0_wp
@@ -807,14 +785,5 @@
             ENDIF
             IF ( (jpert == 5) .OR. (jpert == jpertmax) ) THEN
-               CALL grid_random(  znssh, 'T', 0.0_wp, stdemp )
-
-               DO jj = nldj, nlej
-                  DO ji = nldi, nlei
-                     zemp_tlin(ji,jj) = znssh(ji,jj)
-                  END DO
-               END DO
-
-            ENDIF
-            IF ( (jpert == 6) .OR. (jpert == jpertmax) ) THEN
+
                CALL grid_random(  znssh, 'T', 0.0_wp, stdssh )
                DO jj = nldj, nlej
@@ -834,16 +803,14 @@
             ub_tl(:,:,:) = zub_tlin(:,:,:)
             vb_tl(:,:,:) = zvb_tlin(:,:,:)
-            emp_tl (:,:) = zemp_tlin (:,:)
             sshn_tl(:,:) = zsshn_tlin(:,:)
 
-            gcx_tl (:,:) = 0.e0              ;   gcxb_tl(:,:) = 0.e0
             gcb_tl (:,:) = 0.e0
             gcx_tl (:,:) = zgcx_tlin (:,:)   ;   gcxb_tl(:,:) = zgcxb_tlin(:,:)
 
+            CALL sol_mat( istp ) ! for nitend, it is not called in _tan so it is still set to the nit000 case   
             CALL dyn_spg_flt_tan( istp, indic )
 
             zua_tlout(:,:,:) = ua_tl(:,:,:)   ;   zva_tlout(:,:,:) = va_tl(:,:,:)
-            zspgu_tlout(:,:) = spgu_tl(:,:)   ;   zspgv_tlout(:,:) = spgv_tl(:,:)
-            zgcb_tlout (:,:) = gcb_tl (:,:)
+            zgcxb_tlout(:,:) = gcxb_tl(:,:)   ;   zgcx_tlout (:,:) = gcx_tl (:,:)
 
             !--------------------------------------------------------------------
@@ -865,10 +832,8 @@
             DO jj = nldj, nlej
                DO ji = nldi, nlei
-                  zgcb_adin (ji,jj) = zgcb_tlout (ji,jj) &
-                     &              * e1t(ji,jj) * e2t(ji,jj) * fse3u(ji,jj,1) * tmask(ji,jj,1)
-                  zspgu_adin (ji,jj) = zspgu_tlout (ji,jj) &
-                     &              * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,1) * umask(ji,jj,1)
-                  zspgv_adin(ji,jj) = zspgv_tlout(ji,jj) &
-                     &              * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,1) * vmask(ji,jj,1)
+                  zgcx_adin (ji,jj) = zgcx_tlout (ji,jj) &
+                     &              * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,1) * tmask(ji,jj,1)
+                  zgcxb_adin(ji,jj) = zgcxb_tlout(ji,jj) &
+                     &              * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,1) * tmask(ji,jj,1)
                END DO
             END DO
@@ -878,9 +843,8 @@
             !--------------------------------------------------------------------
 
-            zsp1 =   DOT_PRODUCT( zua_tlout  , zua_adin   ) &
-               &   + DOT_PRODUCT( zgcb_tlout , zgcb_adin  ) &
-               &   + DOT_PRODUCT( zspgu_tlout , zspgu_adin  ) &
-               &   + DOT_PRODUCT( zspgv_tlout , zspgv_adin  ) &
-               &   + DOT_PRODUCT( zva_tlout  , zva_adin   )
+            zsp1 =   DOT_PRODUCT( zua_tlout   , zua_adin   ) &
+               &   + DOT_PRODUCT( zgcx_tlout  , zgcx_adin  ) &
+               &   + DOT_PRODUCT( zgcxb_tlout , zgcxb_adin ) &
+               &   + DOT_PRODUCT( zva_tlout   , zva_adin   )
 
 
@@ -892,9 +856,6 @@
             va_ad(:,:,:) = zva_adin(:,:,:)
 
-            gcx_ad (:,:)   = 0.0_wp            ;   gcxb_ad(:,:) = 0.0_wp
-            gcb_ad (:,:)   = zgcb_adin (:,:)
-            spgu_ad(:,:)   = zspgu_adin(:,:)
-            spgv_ad(:,:)   = zspgv_adin(:,:)
-            ub_ad  (:,:,:) = zub_adin  (:,:,:) ;  vb_ad  (:,:,:) = zvb_adin  (:,:,:)
+            gcx_ad (:,:)   = zgcx_adin (:,:)   ;   gcxb_ad(:,:) = zgcxb_adin (:,:)
+            ub_ad  (:,:,:) = 0.0_wp ;  vb_ad  (:,:,:) = 0.0_wp
 
             CALL dyn_spg_flt_adj( istp, indic )
@@ -919,6 +880,5 @@
                &   + DOT_PRODUCT( zgcx_tlin , zgcx_adout  ) &
                &   + DOT_PRODUCT( zgcxb_tlin, zgcxb_adout ) &
-               &   + DOT_PRODUCT( zsshn_tlin, zsshn_adout ) &
-               &   + DOT_PRODUCT( zemp_tlin , zemp_adout  )
+               &   + DOT_PRODUCT( zsshn_tlin, zsshn_adout )
 
             ! Compare the scalar products
@@ -936,11 +896,9 @@
                   cl_name = 'spg_flt  Vb T1'
                CASE(5)
-                  cl_name = 'spg_flt emp T1'
-               CASE(6)
                   cl_name = 'spg_flt ssh T1'
                CASE(jpertmax)
                   cl_name = 'dyn_spg_flt T1'
                END SELECT
-            ELSEIF ( istp == nit000 + 1 ) THEN
+            ELSEIF ( istp == nitend ) THEN
                SELECT CASE (jpert)
                CASE(1)
@@ -953,6 +911,4 @@
                   cl_name = 'spg_flt  Vb T2'
                CASE(5)
-                  cl_name = 'spg_flt emp T2'
-               CASE(6)
                   cl_name = 'spg_flt ssh T2'
                CASE(jpertmax)
@@ -965,5 +921,5 @@
       END DO
 
-!nn_nmod = kmod  ! restore initial frequency of test for the SOR solver
+      nitsor(:) = jp_it0adj  ! restore nitsor to avoid non reproducible results with or without the tests
 
       ! Deallocate memory
@@ -986,14 +942,9 @@
       DEALLOCATE( &
          & zsshn_tlin, &
-         & zemp_tlin,  &
          & zsshn_adout,&
-         & zemp_adout, &
          & znssh       &
          & )
       DEALLOCATE( zgcx_tlin , zgcx_tlout , zgcx_adin , zgcx_adout,  &
-         & zgcxb_tlin, zgcxb_tlout, zgcxb_adin, zgcxb_adout,  &
-         & zgcb_tlin , zgcb_tlout , zgcb_adin , zgcb_adout    &
-         & )
-      DEALLOCATE ( zub_tlout, zvb_tlout, zub_adin , zvb_adin )
+         & zgcxb_tlin, zgcxb_tlout, zgcxb_adin, zgcxb_adout         ) 
    END SUBROUTINE dyn_spg_flt_adj_tst
 
