Index: /branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/DIA/diaptr.F90
===================================================================
--- /branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/DIA/diaptr.F90	(revision 8969)
+++ /branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/DIA/diaptr.F90	(revision 8970)
@@ -397,5 +397,5 @@
       REWIND( numnam_cfg )              ! Namelist namptr in configuration namelist : Poleward transport
       READ  ( numnam_cfg, namptr, IOSTAT = ios, ERR = 902 )
-902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namptr in configuration namelist', lwp )
+902   IF( ios >  0 ) CALL ctl_nam ( ios , 'namptr in configuration namelist', lwp )
       IF(lwm) WRITE ( numond, namptr )
 
Index: /branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90
===================================================================
--- /branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90	(revision 8969)
+++ /branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90	(revision 8970)
@@ -239,20 +239,20 @@
          IF( nn_print >= 1 ) THEN
             WRITE(numout,*)
-            WRITE(numout,*) '          conversion local  ==> global i-index domain'
+            WRITE(numout,*) '          conversion local  ==> global i-index domain (mig)'
             WRITE(numout,25)              (mig(ji),ji = 1,jpi)
             WRITE(numout,*)
             WRITE(numout,*) '          conversion global ==> local  i-index domain'
-            WRITE(numout,*) '             starting index'
+            WRITE(numout,*) '             starting index (mi0)'
             WRITE(numout,25)              (mi0(ji),ji = 1,jpiglo)
-            WRITE(numout,*) '             ending index'
+            WRITE(numout,*) '             ending index (mi1)'
             WRITE(numout,25)              (mi1(ji),ji = 1,jpiglo)
             WRITE(numout,*)
-            WRITE(numout,*) '          conversion local  ==> global j-index domain'
+            WRITE(numout,*) '          conversion local  ==> global j-index domain (mjg)'
             WRITE(numout,25)              (mjg(jj),jj = 1,jpj)
             WRITE(numout,*)
             WRITE(numout,*) '          conversion global ==> local  j-index domain'
-            WRITE(numout,*) '             starting index'
+            WRITE(numout,*) '             starting index (mj0)'
             WRITE(numout,25)              (mj0(jj),jj = 1,jpjglo)
-            WRITE(numout,*) '             ending index'
+            WRITE(numout,*) '             ending index (mj1)'
             WRITE(numout,25)              (mj1(jj),jj = 1,jpjglo)
          ENDIF
Index: /branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/STO/stopar.F90
===================================================================
--- /branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/STO/stopar.F90	(revision 8969)
+++ /branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/STO/stopar.F90	(revision 8970)
@@ -159,5 +159,5 @@
       INTEGER  :: ji, jj, jk, jsto, jflt
       REAL(wp) :: stomax
-
+      !!----------------------------------------------------------------------
       !
       ! Update 2D stochastic arrays
@@ -235,5 +235,5 @@
          CALL lbc_lnk( sto3d(:,:,:,jsto), sto3d_typ(jsto), sto3d_sgn(jsto) )
       END DO
-
+      !
    END SUBROUTINE sto_par
 
@@ -267,6 +267,15 @@
       REWIND( numnam_cfg )              ! Namelist namdyn_adv in configuration namelist : Momentum advection scheme
       READ  ( numnam_cfg, namsto, IOSTAT = ios, ERR = 902 )
-902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsto in configuration namelist', lwp )
+902   IF( ios >  0 ) CALL ctl_nam ( ios , 'namsto in configuration namelist', lwp )
       IF(lwm) WRITE ( numond, namsto )
+
+      IF( .NOT.ln_rststo ) THEN   ! no use of stochastic parameterization
+         IF(lwp) THEN
+            WRITE(numout,*)
+            WRITE(numout,*) 'sto_par_init : NO use of stochastic parameterization'
+            WRITE(numout,*) '~~~~~~~~~~~~'
+         ENDIF
+         RETURN
+      ENDIF
 
       !IF(ln_ens_rst_in) cn_storst_in = cn_mem//cn_storst_in
@@ -673,8 +682,6 @@
       !!                  ***  ROUTINE sto_rst_read  ***
       !!
-
       !! ** Purpose :   read stochastic parameters from restart file
       !!----------------------------------------------------------------------
-
       INTEGER  :: jsto, jseed
       INTEGER(KIND=8)     ::   ziseed(4)           ! RNG seeds in integer type
@@ -683,4 +690,5 @@
       CHARACTER(LEN=9)    ::   clsto3d='sto3d_000' ! stochastic parameter variable name
       CHARACTER(LEN=10)   ::   clseed='seed0_0000' ! seed variable name
+      !!----------------------------------------------------------------------
 
       IF ( jpsto2d > 0 .OR. jpsto3d > 0 ) THEN
Index: /branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/TRA/trabbc.F90
===================================================================
--- /branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/TRA/trabbc.F90	(revision 8969)
+++ /branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/TRA/trabbc.F90	(revision 8970)
@@ -160,9 +160,9 @@
          !
          CASE ( 1 )                          !* constant flux
-            IF(lwp) WRITE(numout,*) '      *** constant heat flux  =   ', rn_geoflx_cst
+            IF(lwp) WRITE(numout,*) '      ===>>   constant heat flux  =   ', rn_geoflx_cst
             qgh_trd0(:,:) = r1_rau0_rcp * rn_geoflx_cst
             !
          CASE ( 2 )                          !* variable geothermal heat flux : read the geothermal fluxes in mW/m2
-            IF(lwp) WRITE(numout,*) '      *** variable geothermal heat flux'
+            IF(lwp) WRITE(numout,*) '      ===>>   variable geothermal heat flux'
             !
             ALLOCATE( sf_qgh(1), STAT=ierror )
@@ -172,5 +172,5 @@
             ENDIF
             ALLOCATE( sf_qgh(1)%fnow(jpi,jpj,1)   )
-            IF( sn_qgh%ln_tint )ALLOCATE( sf_qgh(1)%fdta(jpi,jpj,1,2) )
+            IF( sn_qgh%ln_tint )   ALLOCATE( sf_qgh(1)%fdta(jpi,jpj,1,2) )
             ! fill sf_chl with sn_chl and control print
             CALL fld_fill( sf_qgh, (/ sn_qgh /), cn_dir, 'tra_bbc_init',   &
@@ -186,5 +186,5 @@
          !
       ELSE
-         IF(lwp) WRITE(numout,*) '      *** no geothermal heat flux'
+         IF(lwp) WRITE(numout,*) '      ===>>   no geothermal heat flux'
       ENDIF
       !
Index: /branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfgls.F90
===================================================================
--- /branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfgls.F90	(revision 8969)
+++ /branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfgls.F90	(revision 8970)
@@ -19,5 +19,4 @@
    USE dom_oce        ! ocean space and time domain
    USE domvvl         ! ocean space and time domain : variable volume layer
-   USE zdf_oce        ! ocean vertical physics
    USE zdfdrg  , ONLY : r_z0_top , r_z0_bot   ! top/bottom roughness
    USE zdfdrg  , ONLY : rCdU_top , rCdU_bot   ! top/bottom friction
@@ -133,4 +132,6 @@
       !!              coefficients using the GLS turbulent closure scheme.
       !!----------------------------------------------------------------------
+      USE zdf_oce , ONLY : en, avtb, avmb   ! ocean vertical physics
+      !!
       INTEGER                   , INTENT(in   ) ::   kt             ! ocean time step
       REAL(wp), DIMENSION(:,:,:), INTENT(in   ) ::   p_sh2          ! shear production term
@@ -580,5 +581,5 @@
          zdep (:,:)   = ((zhsro(:,:) + gdept_n(:,:,1)) / zhsro(:,:))**(rmm*ra_sf)
          zflxs(:,:)   = (rnn + rsbc_tke1 * (rnn + rmm*ra_sf) * zdep(:,:))*(1._wp + rsbc_tke1*zdep(:,:))**(2._wp*rmm/3._wp-1_wp)
-         zdep (:,:)   = rsbc_psi1 * (zwall_psi(:,:,1)*avm(:,:,1)+zwall_psi(:,:,2)*avm(:,:,2)) * &
+         zdep (:,:)   = rsbc_psi1 * (zwall_psi(:,:,1)*p_avm(:,:,1)+zwall_psi(:,:,2)*p_avm(:,:,2)) * &
             &           ustar2_surf(:,:)**rmm * zkar(:,:)**rnn * (zhsro(:,:) + gdept_n(:,:,1))**(rnn-1.)
          zflxs(:,:)   = zdep(:,:) * zflxs(:,:)
@@ -815,9 +816,9 @@
          END DO
       END DO
-      avt(:,:,1)  = 0._wp
+      p_avt(:,:,1) = 0._wp
       !
       IF(ln_ctl) THEN
-         CALL prt_ctl( tab3d_1=en , clinfo1=' gls  - e: ', tab3d_2=avt, clinfo2=' t: ', ovlap=1, kdim=jpk)
-         CALL prt_ctl( tab3d_1=avm, clinfo1=' gls  - m: ', ovlap=1, kdim=jpk )
+         CALL prt_ctl( tab3d_1=en   , clinfo1=' gls  - e: ', tab3d_2=p_avt, clinfo2=' t: ', ovlap=1, kdim=jpk)
+         CALL prt_ctl( tab3d_1=p_avm, clinfo1=' gls  - m: ', ovlap=1, kdim=jpk )
       ENDIF
       !
@@ -1144,4 +1145,6 @@
       !!                set to rn_emin or recomputed (nn_igls/=0)
       !!----------------------------------------------------------------------
+      USE zdf_oce , ONLY : en, avt_k, avm_k   ! ocean vertical physics
+      !!
       INTEGER         , INTENT(in) ::   kt     ! ocean time-step
       CHARACTER(len=*), INTENT(in) ::   cdrw   ! "READ"/"WRITE" flag
Index: /branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftke.F90
===================================================================
--- /branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftke.F90	(revision 8969)
+++ /branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftke.F90	(revision 8970)
@@ -43,5 +43,4 @@
    USE domvvl         ! domain: variable volume layer
    USE sbc_oce        ! surface boundary condition: ocean
-   USE zdf_oce        ! vertical physics: ocean variables
    USE zdfdrg         ! vertical physics: top/bottom drag coef.
    USE zdfmxl         ! vertical physics: mixed layer
@@ -185,4 +184,6 @@
       !! ** Action  : - en : now turbulent kinetic energy)
       !! ---------------------------------------------------------------------
+      USE zdf_oce , ONLY : en   ! ocean vertical physics
+      !!
       REAL(wp), DIMENSION(:,:,:) , INTENT(in   ) ::   pdepw          ! depth of w-points
       REAL(wp), DIMENSION(:,:,:) , INTENT(in   ) ::   p_e3t, p_e3w   ! level thickness (t- & w-points)
@@ -469,4 +470,6 @@
       !! ** Action  : - avt, avm : now vertical eddy diffusivity and viscosity (w-point)
       !!----------------------------------------------------------------------
+      USE zdf_oce , ONLY : en, avtb, avmb, avtb_2d   ! ocean vertical physics
+      !!
       REAL(wp), DIMENSION(:,:,:), INTENT(in   ) ::   pdepw          ! depth (w-points)
       REAL(wp), DIMENSION(:,:,:), INTENT(in   ) ::   p_e3t, p_e3w   ! level thickness (t- & w-points)
@@ -618,6 +621,6 @@
       !
       IF(ln_ctl) THEN
-         CALL prt_ctl( tab3d_1=en , clinfo1=' tke  - e: ', tab3d_2=avt, clinfo2=' t: ', ovlap=1, kdim=jpk)
-         CALL prt_ctl( tab3d_1=avm, clinfo1=' tke  - m: ', ovlap=1, kdim=jpk )
+         CALL prt_ctl( tab3d_1=en   , clinfo1=' tke  - e: ', tab3d_2=p_avt, clinfo2=' t: ', ovlap=1, kdim=jpk)
+         CALL prt_ctl( tab3d_1=p_avm, clinfo1=' tke  - m: ', ovlap=1, kdim=jpk )
       ENDIF
       !
@@ -641,10 +644,12 @@
       !! ** Action  :   Increase by 1 the nstop flag is setting problem encounter
       !!----------------------------------------------------------------------
+      USE zdf_oce , ONLY : ln_zdfiwm   ! Internal Wave Mixing flag
+      !!
       INTEGER ::   ji, jj, jk   ! dummy loop indices
       INTEGER ::   ios
       !!
-      NAMELIST/namzdf_tke/ rn_ediff, rn_ediss , rn_ebb , rn_emin  ,   &
-         &                 rn_emin0, rn_bshear, nn_mxl , ln_mxl0  ,   &
-         &                 rn_mxl0 , nn_pdl   , ln_drg , ln_lc    , rn_lc    ,   &
+      NAMELIST/namzdf_tke/ rn_ediff, rn_ediss , rn_ebb , rn_emin  ,          &
+         &                 rn_emin0, rn_bshear, nn_mxl , ln_mxl0  ,          &
+         &                 rn_mxl0 , nn_pdl   , ln_drg , ln_lc    , rn_lc,   &
          &                 nn_etau , nn_htau  , rn_efr   
       !!----------------------------------------------------------------------
@@ -711,5 +716,5 @@
       IF( nn_htau < 0   .OR.  nn_htau > 1 )   CALL ctl_stop( 'bad flag: nn_htau is 0, 1 or 2 ' )
       IF( nn_etau == 3 .AND. .NOT. ln_cpl )   CALL ctl_stop( 'nn_etau == 3 : HF taum only known in coupled mode' )
-
+      !
       IF( ln_mxl0 ) THEN
          IF(lwp) WRITE(numout,*) '   use a surface mixing length = F(stress) :   set rn_mxl0 = rmxl_min'
@@ -728,12 +733,6 @@
          END SELECT
       ENDIF
-      !                               !* set vertical eddy coef. to the background value
-      DO jk = 1, jpk
-         avt(:,:,jk) = avtb(jk) * wmask(:,:,jk)
-         avm(:,:,jk) = avmb(jk) * wmask(:,:,jk)
-      END DO
-      dissl(:,:,:) = 1.e-12_wp
-      !                              
-      CALL tke_rst( nit000, 'READ' )  !* read or initialize all required files
+      !                                !* read or initialize all required files
+      CALL tke_rst( nit000, 'READ' )      ! (en, avt_k, avm_k, dissl) 
       !
    END SUBROUTINE zdf_tke_init
@@ -750,4 +749,6 @@
       !!                set to rn_emin or recomputed 
       !!----------------------------------------------------------------------
+      USE zdf_oce , ONLY : en, avt_k, avm_k   ! ocean vertical physics
+      !!
       INTEGER         , INTENT(in) ::   kt     ! ocean time-step
       CHARACTER(len=*), INTENT(in) ::   cdrw   ! "READ"/"WRITE" flag
@@ -772,10 +773,12 @@
             ELSE                                          ! start TKE from rest
                IF(lwp) WRITE(numout,*) '   ==>>   previous run without TKE scheme, set en to background values'
-               en(:,:,:) = rn_emin * wmask(:,:,:)
+               en   (:,:,:) = rn_emin * wmask(:,:,:)
+               dissl(:,:,:) = 1.e-12_wp
                ! avt_k, avm_k already set to the background value in zdf_phy_init
             ENDIF
          ELSE                                   !* Start from rest
             IF(lwp) WRITE(numout,*) '   ==>>   start from rest: set en to the background value'
-            en(:,:,:) = rn_emin * wmask(:,:,:)
+            en   (:,:,:) = rn_emin * wmask(:,:,:)
+            dissl(:,:,:) = 1.e-12_wp
             ! avt_k, avm_k already set to the background value in zdf_phy_init
          ENDIF
Index: /branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/nemogcm.F90
===================================================================
--- /branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/nemogcm.F90	(revision 8969)
+++ /branches/2017/dev_CNRS_2017/NEMOGCM/NEMO/OPA_SRC/nemogcm.F90	(revision 8970)
@@ -458,5 +458,5 @@
 
       !                                      ! Active tracers
-                           CALL tra_qsr_init      ! penetrative solar radiation qsr
+      IF( ln_traqsr    )   CALL tra_qsr_init      ! penetrative solar radiation qsr
                            CALL tra_bbc_init      ! bottom heat flux
       IF( ln_trabbl    )   CALL tra_bbl_init      ! advective (and/or diffusive) bottom boundary layer scheme
