Index: /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/DIA/diaar5.F90
===================================================================
--- /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/DIA/diaar5.F90	(revision 2026)
+++ /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/DIA/diaar5.F90	(revision 2027)
@@ -23,5 +23,6 @@
    PRIVATE
 
-   PUBLIC   dia_ar5   ! routine called in step.F90 module
+   PUBLIC   dia_ar5        ! routine called in step.F90 module
+   PUBLIC   dia_ar5_init   ! routine called in opa.F90 module
 
    LOGICAL, PUBLIC, PARAMETER :: lk_diaar5 = .TRUE.   ! coupled flag
@@ -67,6 +68,4 @@
       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zrhd, zrhop
       !!--------------------------------------------------------------------
-
-      IF( kt == nit000  )   CALL dia_ar5_init   ! Initialization (first time-step only)
 
       CALL iom_put( 'cellthc', fse3t(:,:,:) )
Index: /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/DYN/dynadv.F90
===================================================================
--- /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/DYN/dynadv.F90	(revision 2026)
+++ /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/DYN/dynadv.F90	(revision 2027)
@@ -22,5 +22,6 @@
    PRIVATE
 
-   PUBLIC dyn_adv     ! routine called by step module
+   PUBLIC dyn_adv       ! routine called by step module
+   PUBLIC dyn_adv_init  ! routine called by opa module
  
    LOGICAL, PUBLIC ::   ln_dynadv_vec  = .TRUE.    ! vector form flag
@@ -56,6 +57,4 @@
       !!----------------------------------------------------------------------
       !
-      IF( kt == nit000 )   CALL dyn_adv_ctl          ! initialisation & control of options
-
       SELECT CASE ( nadv )                     ! compute advection trend and add it to general trend
       CASE ( 0 )     
@@ -77,7 +76,7 @@
 
    
-   SUBROUTINE dyn_adv_ctl
+   SUBROUTINE dyn_adv_init
       !!---------------------------------------------------------------------
-      !!                  ***  ROUTINE dyn_adv_ctl  ***
+      !!                  ***  ROUTINE dyn_adv_init  ***
       !!                
       !! ** Purpose :   Control the consistency between namelist options for 
@@ -94,5 +93,5 @@
       IF(lwp) THEN                    ! Namelist print
          WRITE(numout,*)
-         WRITE(numout,*) 'dyn_adv_ctl : choice/control of the momentum advection scheme'
+         WRITE(numout,*) 'dyn_adv_init : choice/control of the momentum advection scheme'
          WRITE(numout,*) '~~~~~~~~~~~'
          WRITE(numout,*) '       Namelist namdyn_adv : chose a advection formulation & scheme for momentum'
@@ -124,5 +123,5 @@
       ENDIF
       !
-   END SUBROUTINE dyn_adv_ctl
+   END SUBROUTINE dyn_adv_init
 
   !!======================================================================
Index: /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/DYN/dynhpg.F90
===================================================================
--- /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/DYN/dynhpg.F90	(revision 2026)
+++ /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/DYN/dynhpg.F90	(revision 2027)
@@ -18,5 +18,5 @@
    !!   dyn_hpg      : update the momentum trend with the now horizontal
    !!                  gradient of the hydrostatic pressure
-   !!       hpg_ctl  : initialisation and control of options
+   !!       hpg_init : initialisation and control of options
    !!       hpg_zco  : z-coordinate scheme
    !!       hpg_zps  : z-coordinate plus partial steps (interpolation)
@@ -39,5 +39,6 @@
    PRIVATE
 
-   PUBLIC   dyn_hpg    ! routine called by step module
+   PUBLIC   dyn_hpg        ! routine called by step module
+   PUBLIC       hpg_init   ! routine called by opa module
 
    !                                              !!* Namelist namdyn_hpg : hydrostatic pressure gradient 
@@ -81,6 +82,4 @@
       !!----------------------------------------------------------------------
    
-      IF( kt == nit000 )   CALL hpg_ctl      ! initialisation & control of options
-
       IF( l_trddyn ) THEN                    ! Temporary saving of ua and va trends (l_trddyn)
          ztrdu(:,:,:) = ua(:,:,:)  
@@ -110,7 +109,7 @@
 
 
-   SUBROUTINE hpg_ctl
-      !!----------------------------------------------------------------------
-      !!                 ***  ROUTINE hpg_ctl  ***
+   SUBROUTINE hpg_init
+      !!----------------------------------------------------------------------
+      !!                 ***  ROUTINE hpg_init  ***
       !!
       !! ** Purpose :   initializations for the hydrostatic pressure gradient
@@ -174,5 +173,5 @@
 
       !
-   END SUBROUTINE hpg_ctl
+   END SUBROUTINE hpg_init
 
 
Index: /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/DYN/dynldf.F90
===================================================================
--- /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/DYN/dynldf.F90	(revision 2026)
+++ /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/DYN/dynldf.F90	(revision 2027)
@@ -30,5 +30,6 @@
    PRIVATE
 
-   PUBLIC   dyn_ldf   ! called by step module 
+   PUBLIC   dyn_ldf       ! called by step module 
+   PUBLIC   dyn_ldf_init  ! called by opa module 
 
    INTEGER ::   nldf = 0   ! type of lateral diffusion used defined from ln_dynldf_... namlist logicals)
@@ -56,6 +57,4 @@
       !!----------------------------------------------------------------------
 
-      IF( kt == nit000 )   CALL dyn_ldf_ctl      ! initialisation & control of options
-
       IF( l_trddyn )   THEN                      ! temporary save of ta and sa trends
          ztrdu(:,:,:) = ua(:,:,:) 
@@ -69,10 +68,4 @@
       CASE ( 2 )    ;   CALL dyn_ldf_bilap  ( kt )      ! iso-level bilaplacian
       CASE ( 3 )    ;   CALL dyn_ldf_bilapg ( kt )      ! s-coord. horizontal bilaplacian
-      CASE ( 4 )                                        ! iso-level laplacian + bilaplacian
-         CALL dyn_ldf_lap    ( kt )
-         CALL dyn_ldf_bilap  ( kt )
-      CASE ( 5 )                                        ! rotated laplacian + bilaplacian (s-coord)
-         CALL dyn_ldf_iso    ( kt )
-         CALL dyn_ldf_bilapg ( kt )
       !
       CASE ( -1 )                                       ! esopa: test all possibility with control print
@@ -89,11 +82,4 @@
                         CALL prt_ctl( tab3d_1=ua, clinfo1=' ldf3 - Ua: ', mask1=umask,   &
             &                         tab3d_2=va, clinfo2=       ' Va: ', mask2=vmask, clinfo3='dyn' )
-      !
-      CASE ( -2 )                                       ! neither laplacian nor bilaplacian schemes used
-         IF( kt == nit000 ) THEN
-            IF(lwp) WRITE(numout,*)
-            IF(lwp) WRITE(numout,*) 'dyn_ldf : no lateral diffusion on momentum setup'
-            IF(lwp) WRITE(numout,*) '~~~~~~~ '
-         ENDIF
       END SELECT
 
@@ -110,7 +96,7 @@
 
 
-   SUBROUTINE dyn_ldf_ctl
+   SUBROUTINE dyn_ldf_init
       !!----------------------------------------------------------------------
-      !!                  ***  ROUTINE dyn_ldf_ctl  ***
+      !!                  ***  ROUTINE dyn_ldf_init  ***
       !! 
       !! ** Purpose :   initializations of the horizontal ocean dynamics physics
@@ -123,5 +109,5 @@
       IF(lwp) THEN                        ! Namelist print
          WRITE(numout,*)
-         WRITE(numout,*) 'dyn_ldf_ctl : Choice of the lateral diffusive operator on dynamics'
+         WRITE(numout,*) 'dyn_ldf_init : Choice of the lateral diffusive operator on dynamics'
          WRITE(numout,*) '~~~~~~~~~~~'
          WRITE(numout,*) '       Namelist nam_dynldf : set lateral mixing parameters (type, direction, coefficients)'
@@ -137,5 +123,5 @@
       IF( ln_dynldf_lap   )   ioptio = ioptio + 1
       IF( ln_dynldf_bilap )   ioptio = ioptio + 1
-      IF( ioptio <  1 ) CALL ctl_warn( '          neither laplacian nor bilaplacian operator set for dynamics' )
+      IF( ioptio /= 1 ) CALL ctl_stop( '          use ONE of the 2 lap/bilap operator type on dynamics' )
       ioptio = 0
       IF( ln_dynldf_level )   ioptio = ioptio + 1
@@ -157,5 +143,5 @@
             IF ( ln_dynldf_iso   )   nldf = 1      ! isoneutral (   rotation)
          ENDIF
-         IF ( ln_sco ) THEN             ! s-coordinate
+         IF ( ln_sco ) THEN             ! z-coordinate
             IF ( ln_dynldf_level )   nldf = 0      ! iso-level  (no rotation)
             IF ( ln_dynldf_hor   )   nldf = 1      ! horizontal (   rotation)
@@ -175,5 +161,5 @@
             IF ( ln_dynldf_iso   )   ierr = 2      ! isoneutral (   rotation)
          ENDIF
-         IF ( ln_sco ) THEN             ! s-coordinate
+         IF ( ln_sco ) THEN             ! z-coordinate
             IF ( ln_dynldf_level )   nldf = 2      ! iso-level  (no rotation)
             IF ( ln_dynldf_hor   )   nldf = 3      ! horizontal (   rotation)
@@ -182,22 +168,4 @@
       ENDIF
       
-      IF( ln_dynldf_lap .AND. ln_dynldf_bilap ) THEN  ! mixed laplacian and bilaplacian operators
-         IF ( ln_zco ) THEN                ! z-coordinate
-            IF ( ln_dynldf_level )   nldf = 4      ! iso-level  (no rotation)
-            IF ( ln_dynldf_hor   )   nldf = 4      ! horizontal (no rotation)
-            IF ( ln_dynldf_iso   )   ierr = 2      ! isoneutral (   rotation)
-         ENDIF
-         IF ( ln_zps ) THEN             ! z-coordinate
-            IF ( ln_dynldf_level )   ierr = 1      ! iso-level not allowed 
-            IF ( ln_dynldf_hor   )   nldf = 4      ! horizontal (no rotation)
-            IF ( ln_dynldf_iso   )   ierr = 2      ! isoneutral (   rotation)
-         ENDIF
-         IF ( ln_sco ) THEN             ! s-coordinate
-            IF ( ln_dynldf_level )   nldf = 4      ! iso-level  (no rotation)
-            IF ( ln_dynldf_hor   )   nldf = 5      ! horizontal (   rotation)
-            IF ( ln_dynldf_iso   )   ierr = 2      ! isoneutral (   rotation)
-         ENDIF
-      ENDIF
-
       IF( lk_esopa )                 nldf = -1     ! esopa test
 
@@ -210,15 +178,12 @@
       IF(lwp) THEN
          WRITE(numout,*)
-         IF( nldf == -2 )   WRITE(numout,*) '              neither laplacian nor bilaplacian schemes used'
          IF( nldf == -1 )   WRITE(numout,*) '              ESOPA test All scheme used'
          IF( nldf ==  0 )   WRITE(numout,*) '              laplacian operator'
-         IF( nldf ==  1 )   WRITE(numout,*) '              rotated laplacian operator'
+         IF( nldf ==  1 )   WRITE(numout,*) '              Rotated laplacian operator'
          IF( nldf ==  2 )   WRITE(numout,*) '              bilaplacian operator'
-         IF( nldf ==  3 )   WRITE(numout,*) '              rotated bilaplacian'
-         IF( nldf ==  4 )   WRITE(numout,*) '              laplacian and bilaplacian operators'
-         IF( nldf ==  5 )   WRITE(numout,*) '              rotated laplacian and bilaplacian operators'
+         IF( nldf ==  3 )   WRITE(numout,*) '              Rotated bilaplacian'
       ENDIF
       !
-   END SUBROUTINE dyn_ldf_ctl
+   END SUBROUTINE dyn_ldf_init
 
    !!======================================================================
Index: /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/DYN/dynspg.F90
===================================================================
--- /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/DYN/dynspg.F90	(revision 2026)
+++ /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/DYN/dynspg.F90	(revision 2027)
@@ -27,5 +27,6 @@
    PRIVATE
 
-   PUBLIC   dyn_spg   ! routine called by step module
+   PUBLIC   dyn_spg        ! routine called by step module
+   PUBLIC   dyn_spg_init   ! routine called by opa module
 
    INTEGER ::   nspg = 0   ! type of surface pressure gradient scheme defined from lk_dynspg_... 
@@ -74,6 +75,4 @@
 
 
-      IF( kt == nit000 )   CALL dyn_spg_ctl      ! initialisation & control of options
-
       IF( l_trddyn )   THEN                      ! temporary save of ta and sa trends
          ztrdu(:,:,:) = ua(:,:,:)
@@ -119,7 +118,7 @@
 
 
-   SUBROUTINE dyn_spg_ctl
+   SUBROUTINE dyn_spg_init
       !!---------------------------------------------------------------------
-      !!                  ***  ROUTINE dyn_spg_ctl  ***
+      !!                  ***  ROUTINE dyn_spg_init  ***
       !!                
       !! ** Purpose :   Control the consistency between cpp options for 
@@ -131,5 +130,5 @@
       IF(lwp) THEN             ! Control print
          WRITE(numout,*)
-         WRITE(numout,*) 'dyn_spg_ctl : choice of the surface pressure gradient scheme'
+         WRITE(numout,*) 'dyn_spg_init : choice of the surface pressure gradient scheme'
          WRITE(numout,*) '~~~~~~~~~~~'
          WRITE(numout,*) '     Explicit free surface                  lk_dynspg_exp = ', lk_dynspg_exp
@@ -178,5 +177,5 @@
 #endif
       !
-   END SUBROUTINE dyn_spg_ctl
+   END SUBROUTINE dyn_spg_init
 
   !!======================================================================
Index: /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/DYN/dynvor.F90
===================================================================
--- /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/DYN/dynvor.F90	(revision 2026)
+++ /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/DYN/dynvor.F90	(revision 2027)
@@ -36,5 +36,6 @@
    PRIVATE
 
-   PUBLIC   dyn_vor   ! routine called by step.F90
+   PUBLIC   dyn_vor        ! routine called by step.F90
+   PUBLIC       vor_init   ! routine called by opa.F90
 
    !                                             !!* Namelist namdyn_vor: vorticity term
@@ -74,6 +75,4 @@
       INTEGER, INTENT( in ) ::   kt   ! ocean time-step index
       !!----------------------------------------------------------------------
-
-      IF( kt == nit000 )   CALL vor_ctl          ! initialisation & control of options
 
       !                                          ! vorticity term 
@@ -637,7 +636,7 @@
 
 
-   SUBROUTINE vor_ctl
+   SUBROUTINE vor_init
       !!---------------------------------------------------------------------
-      !!                  ***  ROUTINE vor_ctl  ***
+      !!                  ***  ROUTINE vor_init  ***
       !!
       !! ** Purpose :   Control the consistency between cpp options for
@@ -653,5 +652,5 @@
       IF(lwp) THEN                    ! Namelist print
          WRITE(numout,*)
-         WRITE(numout,*) 'dyn:vor_ctl : vorticity term : read namelist and control the consistency'
+         WRITE(numout,*) 'dyn:vor_init : vorticity term : read namelist and control the consistency'
          WRITE(numout,*) '~~~~~~~~~~~'
          WRITE(numout,*) '        Namelist namdyn_vor : oice of the vorticity term scheme'
@@ -700,5 +699,5 @@
       ENDIF
       !
-   END SUBROUTINE vor_ctl
+   END SUBROUTINE vor_init
 
    !!==============================================================================
Index: /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/DYN/dynzdf.F90
===================================================================
--- /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/DYN/dynzdf.F90	(revision 2026)
+++ /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/DYN/dynzdf.F90	(revision 2027)
@@ -27,5 +27,6 @@
    PRIVATE
 
-   PUBLIC   dyn_zdf    !  routine called by step.F90
+   PUBLIC   dyn_zdf       !  routine called by step.F90
+   PUBLIC   dyn_zdf_init  !  routine called by opa.F90
 
    INTEGER  ::   nzdf = 0              ! type vertical diffusion algorithm used 
@@ -57,6 +58,4 @@
       REAL(wp), DIMENSION(jpi,jpj,jpk) ::   ztrdu, ztrdv   ! 3D workspace
       !!---------------------------------------------------------------------
-
-      IF( kt == nit000 )   CALL zdf_ctl          ! initialisation & control of options
 
       !                                          ! set time step
@@ -96,7 +95,7 @@
 
 
-   SUBROUTINE zdf_ctl
+   SUBROUTINE dyn_zdf_init
       !!----------------------------------------------------------------------
-      !!                 ***  ROUTINE zdf_ctl  ***
+      !!                 ***  ROUTINE dyn_zdf_init  ***
       !!
       !! ** Purpose :   initializations of the vertical diffusion scheme
@@ -124,5 +123,5 @@
       IF(lwp) THEN                                  ! Print the choice
          WRITE(numout,*)
-         WRITE(numout,*) 'dyn:zdf_ctl : vertical dynamics physics scheme'
+         WRITE(numout,*) 'dyn_zdf_init : vertical dynamics physics scheme'
          WRITE(numout,*) '~~~~~~~~~~~'
          IF( nzdf == -1 )   WRITE(numout,*) '              ESOPA test All scheme used'
@@ -131,5 +130,5 @@
       ENDIF
       !
-   END SUBROUTINE zdf_ctl
+   END SUBROUTINE dyn_zdf_init
 
    !!==============================================================================
Index: /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/FLO/floats.F90
===================================================================
--- /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/FLO/floats.F90	(revision 2026)
+++ /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/FLO/floats.F90	(revision 2027)
@@ -25,4 +25,5 @@
 
    PUBLIC   flo_stp    ! routine called by step.F90
+   PUBLIC   flo_init   ! routine called by opa.F90
 
    !!----------------------------------------------------------------------
@@ -52,6 +53,4 @@
          IF(lwp) WRITE(numout,*) 'flo_stp : call floats routine '
          IF(lwp) WRITE(numout,*) '~~~~~~~'
-
-         CALL flo_init           ! read the namelist of floats             
 
          CALL flo_dom            ! compute/read initial position of floats
Index: /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/LDF/ldfslp.F90
===================================================================
--- /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/LDF/ldfslp.F90	(revision 2026)
+++ /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/LDF/ldfslp.F90	(revision 2027)
@@ -31,5 +31,6 @@
    PRIVATE
 
-   PUBLIC   ldf_slp    ! routine called by step.F90
+   PUBLIC   ldf_slp         ! routine called by step.F90
+   PUBLIC   ldf_slp_init    ! routine called by opa.F90
 
    LOGICAL , PUBLIC, PARAMETER              ::   lk_ldfslp = .TRUE.   !: slopes flag
@@ -97,6 +98,4 @@
       !!----------------------------------------------------------------------
       
-      IF( kt == nit000 )   CALL ldf_slp_init      ! initialization (first time-step only)
-
       zeps  =  1.e-20                             ! Local constant initialization
       zmg   = -1.0 / grav
Index: /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/ZDF/zdfbfr.F90
===================================================================
--- /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/ZDF/zdfbfr.F90	(revision 2026)
+++ /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/ZDF/zdfbfr.F90	(revision 2027)
@@ -26,5 +26,6 @@
    PRIVATE
 
-   PUBLIC   zdf_bfr    ! called by step.F90
+   PUBLIC   zdf_bfr         ! called by step.F90
+   PUBLIC   zdf_bfr_init    ! called by opa.F90
    
    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   bfrua , bfrva   !: Bottom friction coefficients set in zdfbfr
@@ -73,7 +74,4 @@
       REAL(wp) ::   zvu, zuv, zecu, zecv   ! temporary scalars
       !!----------------------------------------------------------------------
-
-
-      IF( kt == nit000 )   CALL zdf_bfr_init   ! initialisation
 
       IF( nn_bfr == 2 ) THEN                 ! quadratic botton friction
Index: /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/ZDF/zdfddm.F90
===================================================================
--- /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/ZDF/zdfddm.F90	(revision 2026)
+++ /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/ZDF/zdfddm.F90	(revision 2027)
@@ -24,5 +24,6 @@
    PRIVATE
 
-   PUBLIC   zdf_ddm    ! called by step.F90
+   PUBLIC   zdf_ddm       ! called by step.F90
+   PUBLIC   zdf_ddm_init  ! called by opa.F90
 
    LOGICAL , PUBLIC, PARAMETER ::   lk_zdfddm = .TRUE.  !: double diffusive mixing flag
@@ -85,6 +86,4 @@
       REAL(wp), DIMENSION(jpi,jpj) ::   zmsks, zmskf, zmskd1, zmskd2, zmskd3   ! 2D workspace 
       !!----------------------------------------------------------------------
-
-      IF ( kt == nit000 )   CALL zdf_ddm_init          ! Initialization (first time-step only)
 
       !                                                ! ===============
Index: /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/ZDF/zdfkpp.F90
===================================================================
--- /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/ZDF/zdfkpp.F90	(revision 2026)
+++ /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/ZDF/zdfkpp.F90	(revision 2027)
@@ -28,11 +28,16 @@
    USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
    USE prtctl          ! Print control
-   USE trdmod          ! momentum/tracers trends 
+   USE trdmod_oce      ! ocean trends definition
+   USE trdtra          ! tracers trends
 
    IMPLICIT NONE
    PRIVATE
 
-   PUBLIC   zdf_kpp   ! routine called by step.F90
-   PUBLIC   tra_kpp   ! routine called by step.F90
+   PUBLIC   zdf_kpp       ! routine called by step.F90
+   PUBLIC   zdf_kpp_init  ! routine called by opa.F90
+   PUBLIC   tra_kpp       ! routine called by step.F90
+#if defined key_top
+   PUBLIC   trc_kpp       ! routine called by trcstp.F90
+#endif
 
    LOGICAL , PUBLIC, PARAMETER ::   lk_zdfkpp = .TRUE.    !: KPP vertical mixing flag
@@ -1237,14 +1242,14 @@
       !!
       !! ** Purpose :   compute and add to the tracer trend the non-local
-      !!      tracer flux
+      !!                tracer flux
       !!
       !! ** Method  :   ???
       !!
       !! history :
-      !!     9.0  ! 05-11 (G. Madec)  Original code
+      !!     1.0  ! 2005-11 (G. Madec)  Original code
+      !!     3.3  ! 2010-06 (C. Ethe)  Merge TRA-TRC
       !!----------------------------------------------------------------------
       !! * Modules used
-      USE oce, ONLY :    ztrdt => ua,       & ! use ua as 3D workspace
-                         ztrds => va          ! use va as 3D workspace
+      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::   ztrdt, ztrds   ! 3D workspace
       !!----------------------------------------------------------------------
       INTEGER, INTENT(in) :: kt
@@ -1252,26 +1257,24 @@
 
       IF( kt == nit000 ) THEN
-         IF(lwp) WRITE(numout,*)
+         IF(lwp) WRITE(numout,*) 
          IF(lwp) WRITE(numout,*) 'tra_kpp : KPP non-local tracer fluxes'
          IF(lwp) WRITE(numout,*) '~~~~~~~   '
       ENDIF
 
-
-      ! Save ta and sa trends
-      IF( l_trdtra )   THEN
-         ztrdt(:,:,:) = ta(:,:,:)
-         ztrds(:,:,:) = sa(:,:,:)
+      IF( l_trdtra )   THEN                    !* Save ta and sa trends
+         ALLOCATE( ztrdt(jpi,jpj,jpk) )   ;    ztrdt(:,:,:) = tsa(:,:,:,jp_tem)
+         ALLOCATE( ztrds(jpi,jpj,jpk) )   ;    ztrds(:,:,:) = tsa(:,:,:,jp_sal)
       ENDIF
 
       ! add non-local temperature and salinity flux ( in convective case only)
       DO jk = 1, jpkm1
-         DO jj = 2, jpjm1
+         DO jj = 2, jpjm1 
             DO ji = fs_2, fs_jpim1
-               ta(ji,jj,jk) =  ta(ji,jj,jk)                           &
-                  &         - ( ghats(ji,jj,jk  ) * avt(ji,jj,jk  )   &
-                  &           - ghats(ji,jj,jk+1) * avt(ji,jj,jk+1) ) * wt0(ji,jj) / fse3t(ji,jj,jk)
-               sa(ji,jj,jk) = sa(ji,jj,jk)                              &
-                  &         - ( ghats(ji,jj,jk  ) * fsavs(ji,jj,jk  )   &
-                  &           - ghats(ji,jj,jk+1) * fsavs(ji,jj,jk+1) ) * ws0(ji,jj) / fse3t(ji,jj,jk)
+               tsa(ji,jj,jk,jp_tem) =  tsa(ji,jj,jk,jp_tem)                      &
+                  &                 - (  ghats(ji,jj,jk  ) * avt  (ji,jj,jk  )   & 
+                  &                    - ghats(ji,jj,jk+1) * avt  (ji,jj,jk+1) ) * wt0(ji,jj) / fse3t(ji,jj,jk)
+               tsa(ji,jj,jk,jp_sal) =  tsa(ji,jj,jk,jp_sal)                      &
+                  &                 - (  ghats(ji,jj,jk  ) * fsavs(ji,jj,jk  )   & 
+                  &                    - ghats(ji,jj,jk+1) * fsavs(ji,jj,jk+1) ) * ws0(ji,jj) / fse3t(ji,jj,jk)
             END DO
          END DO
@@ -1280,17 +1283,84 @@
       ! save the non-local tracer flux trends for diagnostic
       IF( l_trdtra )   THEN
-         ztrdt(:,:,:) = ta(:,:,:) - ztrdt(:,:,:)
-         ztrds(:,:,:) = sa(:,:,:) - ztrds(:,:,:)
+         ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:)
+         ztrds(:,:,:) = tsa(:,:,:,jp_sal) - ztrds(:,:,:)
 !!bug gm jpttdzdf ==> jpttkpp
-         CALL trd_mod(ztrdt, ztrds, jptra_trd_zdf, 'TRA', kt)
+         CALL trd_tra( kt, 'TRA', jp_tem, jptra_trd_zdf, ztrdt )
+         CALL trd_tra( kt, 'TRA', jp_sal, jptra_trd_zdf, ztrds )
+         DEALLOCATE( ztrdt )      ;     DEALLOCATE( ztrds )
       ENDIF
 
-      IF(ln_ctl) THEN  
-         CALL prt_ctl( tab3d_1=ta, clinfo1=' kpp  - Ta: ', mask1=tmask,   &
-         &             tab3d_2=sa, clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' )
+      IF(ln_ctl) THEN
+         CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' kpp  - Ta: ', mask1=tmask,   &
+         &             tab3d_2=tsa(:,:,:,jp_sal), clinfo2=       ' Sa: ', mask2=tmask, clinfo3='tra' )
       ENDIF
 
    END SUBROUTINE tra_kpp
 
+#if defined key_top
+   !!----------------------------------------------------------------------
+   !!   'key_top'                                                TOP models
+   !!----------------------------------------------------------------------
+   SUBROUTINE trc_kpp( kt )
+      !!----------------------------------------------------------------------
+      !!                  ***  ROUTINE trc_kpp  ***
+      !!
+      !! ** Purpose :   compute and add to the tracer trend the non-local
+      !!                tracer flux
+      !!
+      !! ** Method  :   ???
+      !!
+      !! history :
+      !!            9.0  ! 2005-11 (G. Madec)  Original code
+      !!       NEMO 3.3  ! 2010-06 (C. Ethe )  Adapted to passive tracers
+      !!----------------------------------------------------------------------
+      USE trc
+      USE prtctl_trc          ! Print control
+      !! * Arguments
+      INTEGER ,                         INTENT( in    )  :: kt     ! ocean time-step index
+      !! * Local declarations
+      INTEGER  ::   ji, jj, jk, jn      ! Dummy loop indices
+      REAL(wp) ::   ztra, zflx
+      REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::   ztrtrd
+      !!----------------------------------------------------------------------
+
+      IF( kt == nit000 ) THEN
+         IF(lwp) WRITE(numout,*) 
+         IF(lwp) WRITE(numout,*) 'trc_kpp : KPP non-local tracer fluxes'
+         IF(lwp) WRITE(numout,*) '~~~~~~~   '
+      ENDIF
+
+      IF( l_trdtrc )  ALLOCATE( ztrtrd(jpi,jpj,jpk) )
+      !
+      DO jn = 1, jptra
+         !
+         IF( l_trdtrc )  ztrtrd(:,:,:)  = tra(:,:,:,jn)
+         ! add non-local on passive tracer flux ( in convective case only)
+         DO jk = 1, jpkm1
+            DO jj = 2, jpjm1 
+               DO ji = fs_2, fs_jpim1
+                  ! Surface tracer flux for non-local term 
+                  zflx = - ( emps(ji,jj) * tra(ji,jj,1,jn) * rcs ) * tmask(ji,jj,1)
+                  ! compute the trend
+                  ztra = - ( ghats(ji,jj,jk  ) * fsavs(ji,jj,jk  )   &
+                  &        - ghats(ji,jj,jk+1) * fsavs(ji,jj,jk+1) ) * zflx / fse3t(ji,jj,jk)
+                  ! add the trend to the general trend
+                  tra(ji,jj,jk,jn) = tra(ji,jj,jk,jn)  + ztra
+               END DO
+            END DO
+         END DO
+         ! save the non-local tracer flux trends for diagnostic
+         IF( l_trdtrc )  ztrtrd(:,:,:)  = tra(:,:,:,jn) - ztrtrd(:,:,:)
+         CALL trd_tra( kt, 'TRC', jn, jptra_trd_zdf, ztrtrd(:,:,:,jn) )
+         !
+      END DO
+      IF( l_trdtrc )  DEALLOCATE( ztrtrd )
+      IF( ln_ctl )   THEN
+         WRITE(charout, FMT="(' kpp')")  ;  CALL prt_ctl_trc_info(charout)
+         CALL prt_ctl_trc( tab4d=tra, mask=tmask, clinfo=clname, clinfo2='trd' )
+      ENDIF
+      !
+   END SUBROUTINE trc_kpp
+#endif
 
    SUBROUTINE zdf_kpp_init
Index: /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/ZDF/zdfric.F90
===================================================================
--- /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/ZDF/zdfric.F90	(revision 2026)
+++ /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/ZDF/zdfric.F90	(revision 2027)
@@ -28,5 +28,6 @@
    PRIVATE
 
-   PUBLIC   zdf_ric    ! called by step.F90
+   PUBLIC   zdf_ric         ! called by step.F90
+   PUBLIC   zdf_ric_init    ! called by opa.F90
 
    LOGICAL, PUBLIC, PARAMETER ::   lk_zdfric = .TRUE.   !: Richardson vertical mixing flag
@@ -82,6 +83,4 @@
       !!----------------------------------------------------------------------
 
-      IF( kt == nit000  ) CALL zdf_ric_init            ! Initialization (first time-step only)
-
       !                                                ! ===============
       DO jk = 2, jpkm1                                 ! Horizontal slab
Index: /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/ZDF/zdftke.F90
===================================================================
--- /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/ZDF/zdftke.F90	(revision 2026)
+++ /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/ZDF/zdftke.F90	(revision 2027)
@@ -54,4 +54,5 @@
 
    PUBLIC   zdf_tke    ! routine called in step module
+   PUBLIC   tke_init   ! routine called in opa module
    PUBLIC   tke_rst    ! routine called in step module
 
@@ -149,6 +150,4 @@
       !!----------------------------------------------------------------------
       !
-      IF( kt == nit000 )   CALL tke_init     ! initialisation 
-                           !
                            CALL tke_tke      ! now tke (en)
                            !
Index: /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/ZDF/zdftke_old.F90
===================================================================
--- /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/ZDF/zdftke_old.F90	(revision 2026)
+++ /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/ZDF/zdftke_old.F90	(revision 2027)
@@ -49,4 +49,5 @@
 
    PUBLIC   zdf_tke_old   ! routine called in step module
+   PUBLIC   zdf_tke_init  ! routine called in opa module
 
    LOGICAL , PUBLIC, PARAMETER              ::   lk_zdftke_old = .TRUE.  !: TKE vertical mixing flag
@@ -182,6 +183,4 @@
       !!--------------------------------------------------------------------
 
-      IF( kt == nit000  )   CALL zdf_tke_init      ! Initialization (first time-step only)
-
       !                                            ! Local constant initialization
       zbbrau =  .5 * rn_ebb / rau0
Index: /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/ZDF/zdftmx.F90
===================================================================
--- /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/ZDF/zdftmx.F90	(revision 2026)
+++ /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/ZDF/zdftmx.F90	(revision 2027)
@@ -27,5 +27,6 @@
    PRIVATE
 
-   PUBLIC   zdf_tmx    ! called in step module 
+   PUBLIC   zdf_tmx         ! called in step module 
+   PUBLIC   zdf_tmx_init    ! called in opa module 
 
    LOGICAL, PUBLIC, PARAMETER ::   lk_zdftmx = .TRUE.    !: tidal mixing flag
@@ -96,7 +97,4 @@
       REAL(wp), DIMENSION(jpi,jpj) ::   zkz   ! temporary 2D workspace
       !!----------------------------------------------------------------------
-
-      !                        
-      IF( kt == nit000  )   CALL zdf_tmx_init      ! Initialization (first time-step only)
 
       !                          ! ----------------------- !
Index: /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/cla.F90
===================================================================
--- /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/cla.F90	(revision 2026)
+++ /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/cla.F90	(revision 2027)
@@ -29,4 +29,5 @@
    !! * Routine accessibility
    PUBLIC tra_cla        ! routine called by step.F90
+   PUBLIC tra_cla_init   ! routine called by opa.F90
 
    !! * Modules variables   
@@ -76,10 +77,4 @@
       !!----------------------------------------------------------------------
  
-      ! cross land advection for straits
-
-      ! Initialization
-      IF( kt == nit000 )   CALL tra_cla_init
-
-
       ! Bab el Mandeb strait horizontal advection
 
Index: /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/opa.F90
===================================================================
--- /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/opa.F90	(revision 2026)
+++ /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/opa.F90	(revision 2027)
@@ -25,4 +25,5 @@
    !!             -   ! 2007-07  (J. Chanut, A. Sellar) Unstructured open boundaries (BDY)
    !!            3.2  ! 2009-08  (S. Masson)  open/write in the listing file in mpp
+   !!            3.3  ! 2010-10  (C. Ethe, G. Madec) reorganisation of initialisation phase
    !!----------------------------------------------------------------------
 
@@ -33,6 +34,6 @@
    !!   opa_closefile  : close remaining files
    !!----------------------------------------------------------------------
-   USE oce             ! dynamics and tracers variables
-   USE dom_oce         ! ocean space domain variables
+
+   USE step_oce        ! Time stepping module definition
    USE sbc_oce         ! surface boundary condition: ocean
    USE trdmod_oce      ! ocean variables trends
@@ -40,12 +41,7 @@
    USE mppini          ! shared/distributed memory setting (mpp_init routine)
    USE domain          ! domain initialization             (dom_init routine)
-   USE obc_par         ! open boundary cond. parameters
    USE obcini          ! open boundary cond. initialization (obc_ini routine)
-   USE bdy_par         ! unstructured open boundary cond. parameters
    USE bdyini          ! unstructured open boundary cond. initialization (bdy_init routine)
    USE istate          ! initial state setting          (istate_init routine)
-   USE sbcmod          ! surface boundary condition 
-   USE eosbn2          ! equation of state                 (eos_init routine)
-   USE dynhpg          ! hydrostatic pressure gradient
    USE ldfdyn          ! lateral viscosity setting      (ldfdyn_init routine)
    USE ldftra          ! lateral diffusivity setting    (ldftra_init routine)
@@ -53,5 +49,4 @@
    USE phycst          ! physical constant                  (par_cst routine)
    USE trdmod          ! momentum/tracers trends       (trd_mod_init routine)
-   USE diaptr          ! poleward transports           (dia_ptr_init routine)
    USE step            ! OPA time-stepping                  (stp     routine)
 #if defined key_oasis3
@@ -60,15 +55,10 @@
    USE cpl_oasis4      ! OASIS4 coupling (not working)
 #endif
-   USE dynspg_oce      ! Control choice of surface pressure gradient schemes
-   USE prtctl          ! Print control                 (prt_ctl_init routine)
    USE c1d             ! 1D configuration
-   USE dyncor_c1d      ! Coriolis factor at T-point
    USE step_c1d        ! Time stepping loop for the 1D configuration
 #if defined key_top
    USE trcini          ! passive tracer initialisation
 #endif
-   
-   USE iom
-   USE in_out_manager  ! I/O manager
+
    USE lib_mpp         ! distributed memory computing
 #if defined key_iomput
@@ -76,5 +66,4 @@
 #endif
 
-   IMPLICIT NONE
    PRIVATE
 
@@ -237,38 +226,70 @@
       ELSE                            ;   CALL mpp_init2     ! eliminate land processors
       ENDIF
-      
-      CALL phy_cst                          ! Physical constants
-      CALL eos_init                         ! Equation of state
-      CALL dom_cfg                          ! Domain configuration
-      CALL dom_init                         ! Domain
-!!gm c1d case can be moved in dom_init routine
-      IF( lk_c1d ) THEN                          ! 1D configuration 
-         CALL cor_c1d                            ! Coriolis defined at T-point
-         umask(:,:,:) = tmask(:,:,:)             ! U, V and T-points are the same
-         vmask(:,:,:) = tmask(:,:,:)             ! 
-      ENDIF
-!!gm c1d end
-
-      IF( ln_ctl )   CALL prt_ctl_init      ! Print control
-
-      IF( lk_obc )   CALL obc_init          ! Open boundaries 
-      IF( lk_bdy )   CALL bdy_init          ! Unstructured open boundaries
-
-      CALL istate_init                      ! ocean initial state (Dynamics and tracers)
+ 
+
+
+     
+
+                            CALL     phy_cst    ! Physical constants
+                            CALL     eos_init   ! Equation of state
+                            CALL     dom_cfg    ! Domain configuration
+                            CALL     dom_init   ! Domain
+
+      IF( ln_ctl        )   CALL prt_ctl_init   ! Print control
+
+      IF( lk_obc        )   CALL     obc_init   ! Open boundaries 
+      IF( lk_bdy        )   CALL     bdy_init   ! Unstructured open boundaries
+
+      IF( ln_zps        )   CALL zps_hde_init   ! Partial steps:  horizontal derivative
+                            CALL  istate_init   ! ocean initial state (Dynamics and tracers)
 
       !                                     ! Ocean physics
-      CALL sbc_init                             ! Read namsbc namelist : surface module (needed for iom_init)
-      CALL ldf_tra_init                         ! Lateral ocean tracer physics
-      CALL ldf_dyn_init                         ! Lateral ocean momentum physics
-      CALL zdf_init                             ! Vertical ocean physics
-
+                            CALL     sbc_init   ! Read namsbc namelist : surface module (needed for iom_init)
+
+      !                                     ! Vertical physics
+                            CALL     zdf_init   ! namelist read
+                            CALL zdf_bfr_init   ! bottom friction
+      IF( lk_zdfric     )   CALL zdf_ric_init   ! Richardson number dependent Kz
+      IF( lk_zdftke_old )   CALL zdf_tke_init   ! TKE closure scheme for Kz (old scheme)
+      IF( lk_zdftke     )   CALL     tke_init   ! TKE closure scheme for Kz
+      IF( lk_zdfkpp     )   CALL zdf_kpp_init   ! KPP closure scheme for Kz
+      IF( lk_zdftmx     )   CALL zdf_tmx_init   ! tidal vertical mixing
+      IF( lk_zdfddm .AND. .NOT. lk_zdfkpp )   & 
+         &                  CALL zdf_ddm_init   ! double diffusive mixing
+      !                                     ! Lateral physics
+                            CALL ldf_tra_init   ! Lateral ocean tracer physics
+                            CALL ldf_dyn_init   ! Lateral ocean momentum physics
+      IF( lk_ldfslp )       CALL ldf_slp_init   ! slope of lateral mixing
+#if defined key_traldf_c2d
+      IF( lk_traldf_eiv )   CALL ldf_eiv( kstp )      ! eddy induced velocity coefficient
+#  endif
+      !                                     ! Active tracers
+                            CALL tra_qsr_init   ! penetrative solar radiation qsr
+      IF( lk_trabbc     )   CALL tra_bbc_init   ! bottom heat flux
+      IF( lk_trabbl     )   CALL tra_bbl_init   ! advective (and/or diffusive) bottom boundary layer scheme
+      IF( lk_tradmp     )   CALL tra_dmp_init   ! internal damping trends
+                            CALL tra_adv_init   ! horizontal & vertical advection
+      IF( n_cla == 1    )   CALL tra_cla_init   ! Cross Land Advection (Update Hor. advection)
+                            CALL tra_ldf_init   ! lateral mixing
+                            CALL tra_zdf_init   ! vertical mixing and after tracer fields
+
+      !                                     ! Dynamics
+                            CALL dyn_adv_init   ! advection (vector or flux form)
+                            CALL     vor_init   ! vorticity term including Coriolis
+                            CALL dyn_ldf_init   ! lateral mixing
+                            CALL     hpg_init   ! horizontal gradient of Hydrostatic pressure
+                            CALL dyn_zdf_init   ! vertical diffusion
+                            CALL dyn_spg_init   ! surface pressure gradient
 #if defined key_top
-      CALL trc_ini                          ! Passive tracers
+      !                                     ! Passive tracers
+                            CALL     trc_ini
 #endif
 
       !                                     ! diagnostics
-      CALL iom_init                             ! iom_put initialization
-      CALL dia_ptr_init                         ! Poleward TRansports initialization
-      CALL trd_mod_init                         ! Mixed-layer/Vorticity/Integral constraints trends
+                            CALL     iom_init   ! iom_put initialization
+      IF( lk_floats    )    CALL     flo_init   ! drifting Floats
+      IF( lk_diaar5    )    CALL dia_ar5_init   ! ar5 diag
+                            CALL dia_ptr_init   ! Poleward TRansports initialization
+                            CALL trd_mod_init   ! Mixed-layer/Vorticity/Integral constraints trends
       !
    END SUBROUTINE opa_init
Index: /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/step.F90
===================================================================
--- /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/step.F90	(revision 2026)
+++ /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/step.F90	(revision 2027)
@@ -26,98 +26,5 @@
    !!   stp            : OPA system time-stepping
    !!----------------------------------------------------------------------
-   USE oce             ! ocean dynamics and tracers variables
-   USE dom_oce         ! ocean space and time domain variables 
-   USE zdf_oce         ! ocean vertical physics variables
-   USE ldftra_oce      ! ocean tracer   - trends
-   USE ldfdyn_oce      ! ocean dynamics - trends
-   USE in_out_manager  ! I/O manager
-   USE iom             !
-   USE lbclnk
-
-   USE daymod          ! calendar                         (day     routine)
-
-   USE dtatem          ! ocean temperature data           (dta_tem routine)
-   USE dtasal          ! ocean salinity    data           (dta_sal routine)
-   USE sbcmod          ! surface boundary condition       (sbc     routine)
-   USE sbcrnf          ! surface boundary condition: runoff variables
-   USE sbccpl          ! surface boundary condition: coupled formulation (call send at end of step)
-   USE cpl_oasis3, ONLY : lk_cpl
-
-#if defined key_top
-   USE trcstp          ! passive tracer time-stepping      (trc_stp routine)
-#endif
-
-   USE traqsr          ! solar radiation penetration      (tra_qsr routine)
-   USE trasbc          ! surface boundary condition       (tra_sbc routine)
-   USE trabbc          ! bottom boundary condition        (tra_bbc routine)
-   USE trabbl          ! bottom boundary layer            (tra_bbl routine)
-   USE tradmp          ! internal damping                 (tra_dmp routine)
-   USE traadv          ! advection scheme control     (tra_adv_ctl routine)
-   USE traldf          ! lateral mixing                   (tra_ldf routine)
-   USE cla             ! cross land advection             (tra_cla routine)
-   !   zdfkpp          ! KPP non-local tracer fluxes      (tra_kpp routine)
-   USE trazdf          ! vertical mixing                  (tra_zdf routine)
-   USE tranxt          ! time-stepping                    (tra_nxt routine)
-   USE tranpc          ! non-penetrative convection       (tra_npc routine)
-
-   USE eosbn2          ! equation of state                (eos_bn2 routine)
-
-   USE dynadv          ! advection                        (dyn_adv routine)
-   USE dynbfr          ! Bottom friction terms            (dyn_bfr routine)
-   USE dynvor          ! vorticity term                   (dyn_vor routine)
-   USE dynhpg          ! hydrostatic pressure grad.       (dyn_hpg routine)
-   USE dynldf          ! lateral momentum diffusion       (dyn_ldf routine)
-   USE dynzdf          ! vertical diffusion               (dyn_zdf routine)
-   USE dynspg_oce      ! surface pressure gradient        (dyn_spg routine)
-   USE dynspg          ! surface pressure gradient        (dyn_spg routine)
-   USE dynnxt          ! time-stepping                    (dyn_nxt routine)
-
-   USE obc_par         ! open boundary condition variables
-   USE obcdta          ! open boundary condition data     (obc_dta routine)
-   USE obcrst          ! open boundary cond. restart      (obc_rst routine)
-   USE obcrad          ! open boundary cond. radiation    (obc_rad routine)
-
-   USE bdy_par         ! unstructured open boundary data variables
-   USE bdydta          ! unstructured open boundary data  (bdy_dta routine)
-
-   USE sshwzv          ! vertical velocity and ssh        (ssh_wzv routine)
-
-   USE ldfslp          ! iso-neutral slopes               (ldf_slp routine)
-   USE ldfeiv          ! eddy induced velocity coef.      (ldf_eiv routine)
-
-   USE zdftmx          ! tide-induced vertical mixing     (zdf_tmx routine)
-   USE zdfbfr          ! bottom friction                  (zdf_bfr routine)
-   USE zdftke_old      ! old TKE vertical mixing      (zdf_tke_old routine)
-   USE zdftke          ! TKE vertical mixing              (zdf_tke routine)
-   USE zdfkpp          ! KPP vertical mixing              (zdf_kpp routine)
-   USE zdfddm          ! double diffusion mixing          (zdf_ddm routine)
-   USE zdfevd          ! enhanced vertical diffusion      (zdf_evd routine)
-   USE zdfric          ! Richardson vertical mixing       (zdf_ric routine)
-   USE zdfmxl          ! Mixed-layer depth                (zdf_mxl routine)
-
-   USE zpshde          ! partial step: hor. derivative     (zps_hde routine)
-
-   USE diawri          ! Standard run outputs             (dia_wri routine)
-   USE trdicp          ! Ocean momentum/tracers trends    (trd_wri routine)
-   USE trdmld          ! mixed-layer trends               (trd_mld routine)
-   USE trdmld_rst      ! restart for mixed-layer trends
-   USE trdmod_oce      ! ocean momentum/tracers trends
-   USE trdvor          ! vorticity budget                 (trd_vor routine)
-   USE diagap          ! hor. mean model-data gap         (dia_gap routine)
-   USE diahdy          ! dynamic height                   (dia_hdy routine)
-   USE diaptr          ! poleward transports              (dia_ptr routine)
-   USE diaar5          ! AR5 diagnosics                   (dia_ar5 routine)
-   USE diahth          ! thermocline depth                (dia_hth routine)
-   USE diafwb          ! freshwater budget                (dia_fwb routine)
-   USE flo_oce         ! floats variables
-   USE floats          ! floats computation               (flo_stp routine)
-
-   USE stpctl          ! time stepping control            (stp_ctl routine)
-   USE restart         ! ocean restart                    (rst_wri routine)
-   USE prtctl          ! Print control                    (prt_ctl routine)
-
-#if defined key_agrif
-   USE agrif_opa_sponge ! Momemtum and tracers sponges
-#endif
+   USE step_oce         ! time stepping definition modules 
 
    IMPLICIT NONE
@@ -265,9 +172,9 @@
                              sa(:,:,:) = 0.e0
 
+                             CALL tra_swap
                              CALL tra_sbc    ( kstp )       ! surface boundary condition
       IF( ln_traqsr      )   CALL tra_qsr    ( kstp )       ! penetrative solar radiation qsr
       IF( lk_trabbc      )   CALL tra_bbc    ( kstp )       ! bottom heat flux
-      IF( lk_trabbl_dif  )   CALL tra_bbl_dif( kstp )       ! diffusive bottom boundary layer scheme
-      IF( lk_trabbl_adv  )   CALL tra_bbl_adv( kstp )       ! advective (and/or diffusive) bottom boundary layer scheme
+      IF( lk_trabbl      )   CALL tra_bbl    ( kstp )       ! advective (and/or diffusive) bottom boundary layer scheme
       IF( lk_tradmp      )   CALL tra_dmp    ( kstp )       ! internal damping trends
                              CALL tra_adv    ( kstp )       ! horizontal & vertical advection
@@ -276,5 +183,7 @@
                              CALL tra_ldf    ( kstp )       ! lateral mixing
 #if defined key_agrif
+                             CALL tra_unswap
       IF(.NOT. Agrif_Root()) CALL Agrif_Sponge_tra          ! tracers sponge
+                             CALL tra_swap
 #endif
                              CALL tra_zdf    ( kstp )       ! vertical mixing and after tracer fields
@@ -283,4 +192,5 @@
          IF( ln_zdfnpc   )   CALL tra_npc    ( kstp )            ! update after fields by non-penetrative convection
                              CALL tra_nxt    ( kstp )            ! tracer fields at next time step
+                             CALL tra_unswap
                              CALL eos( ta, sa, rhd, rhop )       ! Time-filtered in situ density for hpg computation
          IF( ln_zps      )   CALL zps_hde( kstp, ta, sa, rhd,   &   ! Partial steps: time filtered hor. derivative
@@ -289,10 +199,13 @@
          
       ELSE                                                  ! centered hpg  (eos then time stepping)
+                             CALL tra_unswap
                              CALL eos( tn, sn, rhd, rhop )       ! now in situ density for hpg computation
          IF( ln_zps      )   CALL zps_hde( kstp, tn, sn, rhd,   &   ! Partial steps: now horizontal derivative
             &                                   gtu, gsu, gru,  &   ! of t, s, rd at the bottom ocean level
             &                                   gtv, gsv, grv )
+                             CALL tra_swap
          IF( ln_zdfnpc   )   CALL tra_npc    ( kstp )       ! update after fields by non-penetrative convection
                              CALL tra_nxt    ( kstp )       ! tracer fields at next time step
+                             CALL tra_unswap
       ENDIF 
 
Index: /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/step_oce.F90
===================================================================
--- /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/step_oce.F90	(revision 2027)
+++ /branches/DEV_r2006_merge_TRA_TRC/NEMO/OPA_SRC/step_oce.F90	(revision 2027)
@@ -0,0 +1,105 @@
+MODULE step_oce
+   !!======================================================================
+   !!                       ***  MODULE step_oce  ***
+   !! Time-stepping    : manager of the ocean, tracer and ice time stepping
+   !!----------------------------------------------------------------------
+   USE oce             ! ocean dynamics and tracers variables
+   USE dom_oce         ! ocean space and time domain variables 
+   USE zdf_oce         ! ocean vertical physics variables
+   USE ldftra_oce      ! ocean tracer   - trends
+   USE ldfdyn_oce      ! ocean dynamics - trends
+   USE in_out_manager  ! I/O manager
+   USE iom             !
+   USE lbclnk
+
+   USE daymod          ! calendar                         (day     routine)
+
+   USE dtatem          ! ocean temperature data           (dta_tem routine)
+   USE dtasal          ! ocean salinity    data           (dta_sal routine)
+   USE sbcmod          ! surface boundary condition       (sbc     routine)
+   USE sbcrnf          ! surface boundary condition: runoff variables
+   USE sbccpl          ! surface boundary condition: coupled formulation (call send at end of step)
+   USE cpl_oasis3, ONLY : lk_cpl
+
+#if defined key_top
+   USE trcstp          ! passive tracer time-stepping      (trc_stp routine)
+#endif
+
+   USE traqsr          ! solar radiation penetration      (tra_qsr routine)
+   USE trasbc          ! surface boundary condition       (tra_sbc routine)
+   USE trabbc          ! bottom boundary condition        (tra_bbc routine)
+   USE trabbl          ! bottom boundary layer            (tra_bbl routine)
+   USE tradmp          ! internal damping                 (tra_dmp routine)
+   USE traadv          ! advection scheme control     (tra_adv_ctl routine)
+   USE traldf          ! lateral mixing                   (tra_ldf routine)
+   USE cla             ! cross land advection             (tra_cla routine)
+   !   zdfkpp          ! KPP non-local tracer fluxes      (tra_kpp routine)
+   USE trazdf          ! vertical mixing                  (tra_zdf routine)
+   USE tranxt          ! time-stepping                    (tra_nxt routine)
+   USE tranpc          ! non-penetrative convection       (tra_npc routine)
+
+   USE eosbn2          ! equation of state                (eos_bn2 routine)
+
+   USE dynadv          ! advection                        (dyn_adv routine)
+   USE dynbfr          ! Bottom friction terms            (dyn_bfr routine)
+   USE dynvor          ! vorticity term                   (dyn_vor routine)
+   USE dynhpg          ! hydrostatic pressure grad.       (dyn_hpg routine)
+   USE dynldf          ! lateral momentum diffusion       (dyn_ldf routine)
+   USE dynzdf          ! vertical diffusion               (dyn_zdf routine)
+   USE dynspg_oce      ! surface pressure gradient        (dyn_spg routine)
+   USE dynspg          ! surface pressure gradient        (dyn_spg routine)
+   USE dynnxt          ! time-stepping                    (dyn_nxt routine)
+
+   USE obc_par         ! open boundary condition variables
+   USE obcdta          ! open boundary condition data     (obc_dta routine)
+   USE obcrst          ! open boundary cond. restart      (obc_rst routine)
+   USE obcrad          ! open boundary cond. radiation    (obc_rad routine)
+
+   USE bdy_par         ! unstructured open boundary data variables
+   USE bdydta          ! unstructured open boundary data  (bdy_dta routine)
+
+   USE sshwzv          ! vertical velocity and ssh        (ssh_wzv routine)
+
+   USE ldfslp          ! iso-neutral slopes               (ldf_slp routine)
+   USE ldfeiv          ! eddy induced velocity coef.      (ldf_eiv routine)
+
+   USE zdftmx          ! tide-induced vertical mixing     (zdf_tmx routine)
+   USE zdfbfr          ! bottom friction                  (zdf_bfr routine)
+   USE zdftke_old      ! old TKE vertical mixing      (zdf_tke_old routine)
+   USE zdftke          ! TKE vertical mixing              (zdf_tke routine)
+   USE zdfkpp          ! KPP vertical mixing              (zdf_kpp routine)
+   USE zdfddm          ! double diffusion mixing          (zdf_ddm routine)
+   USE zdfevd          ! enhanced vertical diffusion      (zdf_evd routine)
+   USE zdfric          ! Richardson vertical mixing       (zdf_ric routine)
+   USE zdfmxl          ! Mixed-layer depth                (zdf_mxl routine)
+
+   USE zpshde          ! partial step: hor. derivative     (zps_hde routine)
+
+   USE diawri          ! Standard run outputs             (dia_wri routine)
+   USE trdicp          ! Ocean momentum/tracers trends    (trd_wri routine)
+   USE trdmld          ! mixed-layer trends               (trd_mld routine)
+   USE trdmld_rst      ! restart for mixed-layer trends
+   USE trdmod_oce      ! ocean momentum/tracers trends
+   USE trdvor          ! vorticity budget                 (trd_vor routine)
+   USE diagap          ! hor. mean model-data gap         (dia_gap routine)
+   USE diahdy          ! dynamic height                   (dia_hdy routine)
+   USE diaptr          ! poleward transports              (dia_ptr routine)
+   USE diaar5          ! AR5 diagnosics                   (dia_ar5 routine)
+   USE diahth          ! thermocline depth                (dia_hth routine)
+   USE diafwb          ! freshwater budget                (dia_fwb routine)
+   USE flo_oce         ! floats variables
+   USE floats          ! floats computation               (flo_stp routine)
+
+   USE stpctl          ! time stepping control            (stp_ctl routine)
+   USE restart         ! ocean restart                    (rst_wri routine)
+   USE prtctl          ! Print control                    (prt_ctl routine)
+
+   USE traswp          ! Swap arrays                      (tra_swp routine)
+                       !                                (tra_unswp routine) 
+
+#if defined key_agrif
+   USE agrif_opa_sponge ! Momemtum and tracers sponges
+#endif
+
+   !!======================================================================
+END MODULE step_oce
