Index: trunk/NEMO/OPA_SRC/step.F90
===================================================================
--- trunk/NEMO/OPA_SRC/step.F90	(revision 445)
+++ trunk/NEMO/OPA_SRC/step.F90	(revision 467)
@@ -4,5 +4,19 @@
    !! Time-stepping    : manager of the ocean, tracer and ice time stepping
    !!======================================================================
-
+   !! History :
+   !!        !  91-03  ()  Original code
+   !!        !  91-11  (G. Madec)
+   !!        !  92-06  (M. Imbard)  add a first output record
+   !!        !  96-04  (G. Madec)  introduction of dynspg
+   !!        !  96-04  (M.A. Foujols)  introduction of passive tracer
+   !!   8.0  !  97-06  (G. Madec)  new architecture of call
+   !!   8.2  !  97-06  (G. Madec, M. Imbard, G. Roullet)  free surface
+   !!   8.2  !  99-02  (G. Madec, N. Grima)  hpg implicit
+   !!   8.2  !  00-07  (J-M Molines, M. Imbard)  Open Bondary Conditions
+   !!   9.0  !  02-06  (G. Madec)  free form, suppress macro-tasking
+   !!    "   !  04-08  (C. Talandier) New trends organization
+   !!    "   !  05-01  (C. Ethe) Add the KPP closure scheme
+   !!    "   !  05-11  (V. Garnier) Surface pressure gradient organization
+   !!    "   !  05-11  (G. Madec)  Reorganisation of tra and dyn calls
    !!----------------------------------------------------------------------
    !!   stp            : OPA system time-stepping
@@ -12,6 +26,6 @@
    USE dom_oce         ! ocean space and time domain variables 
    USE zdf_oce         ! ocean vertical physics variables
-   USE ldftra_oce
-   USE ldfdyn_oce
+   USE ldftra_oce      ! ocean tracer   - trends
+   USE ldfdyn_oce      ! ocean dynamics - trends
    USE cpl_oce         ! coupled ocean-atmosphere variables
    USE in_out_manager  ! I/O manager
@@ -34,45 +48,28 @@
    USE trcstp          ! passive tracer time-stepping      (trc_stp routine)
 
-   USE dynhpg          ! hydrostatic pressure grad.       (dyn_hpg routine)
-   USE dynhpg_atsk     ! hydrostatic pressure grad.  (dyn_hpg_atsk routine)
-   USE dynspg_oce      ! surface pressure gradient        (dyn_spg routine)
-   USE dynspg          ! surface pressure gradient        (dyn_spg routine)
-   USE dynkeg          ! kinetic energy gradient          (dyn_keg routine)
-   USE dynvor          ! vorticity term              (dyn_vor_... routines)
-   USE dynzad          ! vertical advection               (dyn_adv routine)
-   USE dynldf_bilapg   ! lateral mixing            (dyn_ldf_bilapg routine)
-   USE dynldf_bilap    ! lateral mixing             (dyn_ldf_bilap routine)
-   USE dynldf_iso      ! lateral mixing               (dyn_ldf_iso routine)
-   USE dynldf_lap      ! lateral mixing               (dyn_ldf_lap routine)
-   USE dynzdf_imp      ! vertical diffusion: implicit     (dyn_zdf routine)
-   USE dynzdf_imp_atsk ! vertical diffusion: implicit     (dyn_zdf routine)
-   USE dynzdf_iso      ! vertical diffusion: isopycnal    (dyn_zdf routine)
-   USE dynzdf_exp      ! vertical diffusion: explicit (dyn_zdf_exp routine)
-   USE dynnxt          ! time-stepping                    (dyn_nxt routine)
-
+   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 traldf_bilapg   ! lateral mixing            (tra_ldf_bilapg routine)
-   USE traldf_bilap    ! lateral mixing             (tra_ldf_bilap routine)
-   USE traldf_iso      ! lateral mixing               (tra_ldf_iso routine)
-   USE traldf_iso_zps  ! lateral mixing           (tra_ldf_iso_zps routine)
-   USE traldf_lap      ! lateral mixing               (tra_ldf_lap routine)
-   USE traqsr          ! solar radiation penetration      (tra_qsr 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 tranxt          ! time-stepping                    (tra_nxt routine)
-   USE traadv_ctl      ! advection scheme control     (tra_adv_ctl routine)
-   USE traadv_cen2     ! 2nd order centered scheme   (tra_adv_cen2 routine)
-   USE traadv_tvd      ! TVD scheme                (tra_adv_tvd    routine)
-   USE traadv_muscl    ! MUSCL scheme              (tra_adv_muscl  routine)
-   USE traadv_muscl2   ! MUSCL2 scheme             (tra_adv_muscl2 routine)
-   USE cla             ! cross land advection             (tra_cla routine)
-   USE trazdf_exp      ! vertical diffusion: explicit (tra_zdf_exp routine)
-   USE trazdf_imp      ! vertical diffusion: implicit (tra_zdf_imp routine)
-   USE trazdf_iso      ! vertical diffusion           (tra_zdf_exp routine)
-   USE trazdf_iso_vopt ! vertical diffusion           (tra_zdf_exp routine)
-   USE trasbc          ! surface boundary condition       (tra_sbc routine)
 
    USE eosbn2          ! equation of state                (eos_bn2 routine)
+
+   USE dynhpg          ! hydrostatic pressure grad.       (dyn_hpg routine)
+   USE dynkeg          ! kinetic energy gradient          (dyn_keg routine)
+   USE dynvor          ! vorticity term                   (dyn_vor routine)
+   USE dynzad          ! vertical advection               (dyn_adv 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
@@ -91,4 +88,5 @@
    USE zdfbfr          ! bottom friction                  (zdf_bfr routine)
    USE zdftke          ! TKE vertical mixing              (zdf_tke routine)
+   USE zdftke_jki      ! TKE vertical mixing              (zdf_tke routine)
    USE zdfkpp          ! KPP vertical mixing              (zdf_kpp routine)
    USE zdfddm          ! double diffusion mixing          (zdf_ddm routine)
@@ -140,9 +138,10 @@
 CONTAINS
 
-   SUBROUTINE stp( &
 #if !defined key_agrif
-   kstp &
-#endif   
-   )      !!----------------------------------------------------------------------
+   SUBROUTINE stp( kstp )
+#else
+   SUBROUTINE stp( )
+#endif
+      !!----------------------------------------------------------------------
       !!                     ***  ROUTINE stp  ***
       !!                      
@@ -160,25 +159,11 @@
       !!              -8- Outputs and diagnostics
       !!
-      !! History :
-      !!        !  91-03  ()  Original code
-      !!        !  91-11  (G. Madec)
-      !!        !  92-06  (M. Imbard)  add a first output record
-      !!        !  96-04  (G. Madec)  introduction of dynspg
-      !!        !  96-04  (M.A. Foujols)  introduction of passive tracer
-      !!   8.0  !  97-06  (G. Madec)  new architecture of call
-      !!   8.2  !  97-06  (G. Madec, M. Imbard, G. Roullet)  free surface
-      !!   8.2  !  99-02  (G. Madec, N. Grima)  hpg implicit
-      !!   8.2  !  00-07  (J-M Molines, M. Imbard)  Open Bondary Conditions
-      !!   9.0  !  02-06  (G. Madec)  free form, suppress macro-tasking
-      !!    "   !  04-08  (C. Talandier) New trends organization
-      !!    "   !  05-01  (C. Ethe) Add the KPP closure scheme
-      !!    "   !  05-11  (V. Garnier) Surface pressure gradient organization
       !!----------------------------------------------------------------------
       !! * Arguments
-      INTEGER &
 #if !defined key_agrif   
-      , INTENT( in ) &
+      INTEGER, INTENT( in ) :: kstp   ! ocean time-step index
+#else
+      INTEGER               :: kstp   ! ocean time-step index
 #endif      
-      ::   kstp   ! ocean time-step index
 
       !! * local declarations
@@ -206,5 +191,5 @@
       IF( lk_dtasst  )   CALL dta_sst( kstp )         ! Sea Surface Temperature data
 
-      IF( lk_dtasss  )   CALL dta_sss( kstp )         ! Sea Surface salinity data
+      IF( lk_dtasss  )   CALL dta_sss( kstp )         ! Sea Surface Salinity data
 
       IF( lk_obc     )   CALL obc_dta( kstp )         ! update dynamic and tracer data at open boundaries
@@ -258,6 +243,11 @@
       !                                                     ! Vertical eddy viscosity and diffusivity coefficients
       IF( lk_zdfric )   CALL zdf_ric( kstp )                       ! Richardson number dependent Kz
+#if defined key_mpp_omp
+      IF( lk_zdftke )   CALL zdf_tke_jki( kstp )                   ! TKE closure scheme for Kz - j-k-i loops
+#else
       IF( lk_zdftke )   CALL zdf_tke( kstp )                       ! TKE closure scheme for Kz
+#endif
       IF( lk_zdfkpp )   CALL zdf_kpp( kstp )                       ! KPP closure scheme for Kz
+
       IF( lk_zdfcst )   avt (:,:,:) = avt0 * tmask(:,:,:)          ! Constant Kz (reset avt to the background value)
 
@@ -266,4 +256,6 @@
             CASE ( 05 )                         ! ORCA R2 configuration
                avt  (:,:,2) = avt  (:,:,2) + 1.e-3 * upsrnfh(:,:)   ! increase diffusivity of rivers mouths
+            CASE ( 025 )                         ! ORCA R025 configuration
+               avt  (:,:,2) = avt  (:,:,2) + 2.e-3 * upsrnfh(:,:)   ! increase diffusivity of rivers mouths
          END SELECT
       ENDIF
@@ -310,58 +302,45 @@
       !-----------------------------------------------------------------------
 
-                               ta(:,:,:) = 0.e0               ! set tracer trends to zero
-                               sa(:,:,:) = 0.e0
-
-                               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_tradmp        )   CALL tra_dmp( kstp )           ! internal damping trends
-
-      !                                                       ! horizontal & vertical advection
-      IF( kstp == nit000   )   CALL tra_adv_ctl                    ! chose/control the scheme used
-      IF( ln_traadv_cen2   )   CALL tra_adv_cen2  ( kstp )         ! 2nd order centered scheme
-      IF( ln_traadv_tvd    )   CALL tra_adv_tvd   ( kstp )         ! TVD scheme
-      IF( ln_traadv_muscl  )   CALL tra_adv_muscl ( kstp )         ! MUSCL scheme
-      IF( ln_traadv_muscl2 )   CALL tra_adv_muscl2( kstp )         ! MUSCL2 scheme
-
-      IF( n_cla == 1       )   CALL tra_cla( kstp )           ! Cross Land Advection (Update Hor. advection)
-
-      !                                                       ! lateral mixing 
-      IF( l_traldf_lap     )   CALL tra_ldf_lap    ( kstp )           ! iso-level laplacian
-      IF( l_traldf_bilap   )   CALL tra_ldf_bilap  ( kstp )           ! iso-level bilaplacian 
-      IF( l_traldf_bilapg  )   CALL tra_ldf_bilapg ( kstp )           ! s-coord. horizontal bilaplacian
-      IF( l_traldf_iso     )   CALL tra_ldf_iso    ( kstp )           ! iso-neutral/geopot. laplacian 
-      IF( l_traldf_iso_zps )   CALL tra_ldf_iso_zps( kstp )           ! partial step iso-neutral/geopot. laplacian
-
+                             ta(:,:,:) = 0.e0               ! set tracer trends to zero
+                             sa(:,:,:) = 0.e0
+
+                             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_tradmp      )   CALL tra_dmp    ( kstp )       ! internal damping trends
+
+                             CALL tra_adv    ( kstp )       ! horizontal & vertical advection
+
+      IF( n_cla == 1     )   CALL tra_cla    ( kstp )       ! Cross Land Advection (Update Hor. advection)
+
+      IF( lk_zdfkpp )        CALL tra_kpp    ( kstp )       ! KPP non-local tracer fluxes
+
+                             CALL tra_ldf    ( kstp )       ! lateral mixing
 #if defined key_agrif
       IF (.NOT. Agrif_Root())  CALL Agrif_Sponge_tra( kstp )          ! tracers sponge
 #endif
-      !                                                       ! vertical diffusion
-      IF( l_trazdf_exp     )   CALL tra_zdf_exp     ( kstp )          ! explicit time stepping (time splitting scheme)
-      IF( l_trazdf_imp     )   CALL tra_zdf_imp     ( kstp )          ! implicit time stepping (euler backward)
-      IF( l_trazdf_iso     )   CALL tra_zdf_iso     ( kstp )          ! isopycnal
-      IF( l_trazdf_iso_vo  )   CALL tra_zdf_iso_vopt( kstp )          ! vector opt. isopycnal
-
-                               CALL tra_nxt( kstp )           ! tracer fields at next time step
-
-      IF( ln_zdfnpc        )   CALL tra_npc( kstp )           ! update the new (t,s) fields by non
-      !                                                       ! penetrative convective adjustment
-
-      IF( ln_dynhpg_imp    ) THEN                             ! semi-implicit hpg 
-                                  CALL eos( ta, sa, rhd, rhop )   ! Time-filtered in situ density used in dynhpg module
-         IF( lk_zps    )          CALL zps_hde( kstp, ta, sa, rhd,  & ! Partial steps: time filtered hor. gradient 
-            &                                        gtu, gsu, gru, & ! of t, s, rd at the bottom ocean level
-            &                                        gtv, gsv, grv )  
-      ELSE                                                    ! centered hpg (default case)
-                                  CALL eos( tb, sb, rhd, rhop )       ! now (swap=before) in situ density for dynhpg module
-         IF( lk_zps    )          CALL zps_hde( kstp, tb, sb, rhd,  & ! Partial steps: now horizontal gradient
-            &                                        gtu, gsu, gru, & ! of t, s, rd at the bottom ocean level
-            &                                        gtv, gsv, grv )  
+                             CALL tra_zdf    ( kstp )       ! vertical mixing
+
+                             CALL tra_nxt( kstp )           ! tracer fields at next time step
+
+      IF( ln_zdfnpc      )   CALL tra_npc( kstp )           ! update the new (t,s) fields by non
+      !                                                     ! penetrative convective adjustment
+
+      IF( ln_dynhpg_imp  ) THEN                             ! semi-implicit hpg
+                               CALL eos( ta, sa, rhd, rhop )          ! Time-filtered in situ density used in dynhpg module
+         IF( ln_zps    )       CALL zps_hde( kstp, ta, sa, rhd,   &   ! Partial steps: time filtered hor. gradient
+            &                                     gtu, gsu, gru,  &   ! of t, s, rd at the bottom ocean level
+            &                                     gtv, gsv, grv )
+      ELSE                                                  ! centered hpg (default case)
+                               CALL eos( tb, sb, rhd, rhop )          ! now (swap=before) in situ density for dynhpg module
+         IF( ln_zps    )       CALL zps_hde( kstp, tb, sb, rhd,   &   ! Partial steps: now horizontal gradient
+            &                                     gtu, gsu, gru,  &   ! of t, s, rd at the bottom ocean level
+            &                                     gtv, gsv, grv )
       ENDIF
 
@@ -371,4 +350,5 @@
       ! N.B. ta, sa arrays are used as workspace in this section 
       !-----------------------------------------------------------------------
+
 
                                ua(:,:,:) = 0.e0               ! set dynamics trends to zero
@@ -377,36 +357,17 @@
                                CALL dyn_keg( kstp )           ! horizontal gradient of kinetic energy
 
-      !                                                       ! vorticity term including Coriolis
-      IF( kstp == nit000   )   CALL dyn_vor_ctl                      ! chose/control the scheme used
-      IF( ln_dynvor_ens    )   CALL dyn_vor_enstrophy( kstp )        ! enstrophy conserving scheme
-      IF( ln_dynvor_ene    )   CALL dyn_vor_energy   ( kstp )        ! energy conserving scheme
-      IF( ln_dynvor_mix    )   CALL dyn_vor_mixed    ( kstp )        ! mixed energy/enstrophy conserving scheme
-      IF( ln_dynvor_een    )   CALL dyn_vor_ene_ens  ( kstp )        ! combined energy/enstrophy conserving scheme
-      
-      !                                                       ! lateral mixing 
-      IF( l_dynldf_lap     )   CALL dyn_ldf_lap    ( kstp )          ! iso-level laplacian
-      IF( l_dynldf_bilap   )   CALL dyn_ldf_bilap  ( kstp )          ! iso-level bilaplacian 
-      IF( l_dynldf_bilapg  )   CALL dyn_ldf_bilapg ( kstp )          ! s-coord. horizontal bilaplacian
-      IF( l_dynldf_iso     )   CALL dyn_ldf_iso    ( kstp )          ! iso-neutral laplacian 
-
+                               CALL dyn_vor( kstp )           ! vorticity term including Coriolis
+
+                               CALL dyn_ldf( kstp )           ! lateral mixing
 #if defined key_agrif
       IF (.NOT. Agrif_Root())  CALL Agrif_Sponge_dyn( kstp )         ! momemtum sponge
 #endif
-      !                                                       ! horizontal gradient of Hydrostatic pressure 
-      IF ( lk_jki ) THEN
-                               CALL dyn_hpg_atsk( kstp )             ! autotask case (j-k-i loop)
-      ELSE
-                               CALL dyn_hpg     ( kstp )             ! default case  (k-j-i loop)
-      ENDIF
-
-                               CALL dyn_zad    ( kstp )       ! vertical advection       
-
-      !                                                       ! vertical diffusion
-      IF( l_dynzdf_exp     )   CALL dyn_zdf_exp    ( kstp )          ! explicit time stepping (time splitting scheme)
-      IF( l_dynzdf_imp     )   CALL dyn_zdf_imp    ( kstp )          ! implicit time stepping (euler backward)
-      IF( l_dynzdf_imp_tsk )   CALL dyn_zdf_imp_tsk( kstp )          ! autotask implicit time stepping (euler backward)
-      IF( l_dynzdf_iso     )   CALL dyn_zdf_iso    ( kstp )          ! iso-neutral case
-
-      IF( lk_dynspg_rl ) THEN 
+                               CALL dyn_hpg( kstp )           ! horizontal gradient of Hydrostatic pressure
+
+                               CALL dyn_zad( kstp )           ! vertical advection
+
+                               CALL dyn_zdf( kstp )           ! vertical diffusion
+
+      IF( lk_dynspg_rl ) THEN
          IF( lk_obc    )       CALL obc_spg( kstp )           ! surface pressure gradient at open boundaries
       ENDIF
@@ -417,5 +378,5 @@
                                CALL dyn_spg( kstp, indic )    ! surface pressure gradient
 
-                               CALL dyn_nxt( kstp )           ! velocity at next time step 
+                               CALL dyn_nxt( kstp )           ! lateral velocity at next time step
 
 
