Index: NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/tests/BENCH/MY_SRC/diawri.F90
===================================================================
--- NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/tests/BENCH/MY_SRC/diawri.F90	(revision 10345)
+++ NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/tests/BENCH/MY_SRC/diawri.F90	(revision 10358)
@@ -52,4 +52,5 @@
 
 #if defined key_si3
+   USE ice
    USE icewri 
 #endif
@@ -119,5 +120,5 @@
       ! Output the initial state and forcings
       IF( ninist == 1 ) THEN                       
-         CALL dia_wri_state( 'output.init', kt )
+         CALL dia_wri_state( 'output.init' )
          ninist = 0
       ENDIF
@@ -410,5 +411,5 @@
 
       IF( ninist == 1 ) THEN     !==  Output the initial state and forcings  ==!
-         CALL dia_wri_state( 'output.init', kt )
+         CALL dia_wri_state( 'output.init' )
          ninist = 0
       ENDIF
@@ -418,5 +419,5 @@
 #endif
 
-   SUBROUTINE dia_wri_state( cdfile_name, kt )
+   SUBROUTINE dia_wri_state( cdfile_name )
       !!---------------------------------------------------------------------
       !!                 ***  ROUTINE dia_wri_state  ***
@@ -432,158 +433,56 @@
       !!----------------------------------------------------------------------
       CHARACTER (len=* ), INTENT( in ) ::   cdfile_name      ! name of the file created
-      INTEGER           , INTENT( in ) ::   kt               ! ocean time-step index
-      !! 
-      CHARACTER (len=32) :: clname
-      CHARACTER (len=40) :: clop
-      INTEGER  ::   id_i , nz_i, nh_i       
-      INTEGER, DIMENSION(1) ::   idex             ! local workspace
-      REAL(wp) ::   zsto, zout, zmax, zjulian
+      !!
+      INTEGER :: inum
       !!----------------------------------------------------------------------
       ! 
-      ! 0. Initialisation
-      ! -----------------
-
-      ! Define name, frequency of output and means
-      clname = cdfile_name
-      IF( .NOT. Agrif_Root() ) clname = TRIM(Agrif_CFixed())//'_'//TRIM(clname)
-      zsto = rdt
-      clop = "inst(x)"           ! no use of the mask value (require less cpu time)
-      zout = rdt
-      zmax = ( nitend - nit000 + 1 ) * rdt
-
       IF(lwp) WRITE(numout,*)
       IF(lwp) WRITE(numout,*) 'dia_wri_state : single instantaneous ocean state'
       IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~   and forcing fields file created '
-      IF(lwp) WRITE(numout,*) '                and named :', clname, '.nc'
-
-
-      ! 1. Define NETCDF files and fields at beginning of first time step
-      ! -----------------------------------------------------------------
-
-      ! Compute julian date from starting date of the run
-      CALL ymds2ju( nyear, nmonth, nday, rdt, zjulian )         ! time axis 
-      zjulian = zjulian - adatrj   !   set calendar origin to the beginning of the experiment
-      CALL histbeg( clname, jpi, glamt, jpj, gphit,   &
-          1, jpi, 1, jpj, nit000-1, zjulian, rdt, nh_i, id_i, domain_id=nidom, snc4chunks=snc4set ) ! Horizontal grid : glamt and gphit
-      CALL histvert( id_i, "deptht", "Vertical T levels",   &    ! Vertical grid : gdept
-          "m", jpk, gdept_1d, nz_i, "down")
-
-      ! Declare all the output fields as NetCDF variables
-
-      CALL histdef( id_i, "vosaline", "Salinity"              , "PSU"    ,   &   ! salinity
-         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
-      CALL histdef( id_i, "votemper", "Temperature"           , "C"      ,   &   ! temperature
-         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
-      CALL histdef( id_i, "sossheig", "Sea Surface Height"    , "m"      ,   &  ! ssh
-         &          jpi, jpj, nh_i, 1  , 1, 1  , nz_i, 32, clop, zsto, zout )
-      CALL histdef( id_i, "vozocrtx", "Zonal Current"         , "m/s"    ,   &   ! zonal current
-         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
-      CALL histdef( id_i, "vomecrty", "Meridional Current"    , "m/s"    ,   &   ! meridonal current
-         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
-      CALL histdef( id_i, "vovecrtz", "Vertical Velocity"     , "m/s"    ,   &   ! vertical current
-         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
-         !
+      IF(lwp) WRITE(numout,*) '                and named :', cdfile_name, '...nc'
+
+#if defined key_si3
+     CALL iom_open( TRIM(cdfile_name), inum, ldwrt = .TRUE., kiolib = jprstlib, kdlev = jpl )
+#else
+     CALL iom_open( TRIM(cdfile_name), inum, ldwrt = .TRUE., kiolib = jprstlib )
+#endif
+
+      CALL iom_rstput( 0, 0, inum, 'votemper', tsn(:,:,:,jp_tem) )    ! now temperature
+      CALL iom_rstput( 0, 0, inum, 'vosaline', tsn(:,:,:,jp_sal) )    ! now salinity
+      CALL iom_rstput( 0, 0, inum, 'sossheig', sshn              )    ! sea surface height
+      CALL iom_rstput( 0, 0, inum, 'vozocrtx', un                )    ! now i-velocity
+      CALL iom_rstput( 0, 0, inum, 'vomecrty', vn                )    ! now j-velocity
+      CALL iom_rstput( 0, 0, inum, 'vovecrtz', wn                )    ! now k-velocity
       IF( ALLOCATED(ahtu) ) THEN
-         CALL histdef( id_i, "ahtu"    , "u-eddy diffusivity"    , "m2/s"    ,   &   ! zonal current
-            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
-         CALL histdef( id_i, "ahtv"    , "v-eddy diffusivity"    , "m2/s"    ,   &   ! meridonal current
-            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
-      ENDIF
-      IF( ALLOCATED(ahmt) ) THEN 
-         CALL histdef( id_i, "ahmt"    , "t-eddy viscosity"      , "m2/s"    ,   &   ! zonal current
-            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
-         CALL histdef( id_i, "ahmf"    , "f-eddy viscosity"      , "m2/s"    ,   &   ! meridonal current
-            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
-      ENDIF
-         !
-      CALL histdef( id_i, "sowaflup", "Net Upward Water Flux" , "Kg/m2/S",   &   ! net freshwater 
-         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
-      CALL histdef( id_i, "sohefldo", "Net Downward Heat Flux", "W/m2"   ,   &   ! net heat flux
-         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
-      CALL histdef( id_i, "soshfldo", "Shortwave Radiation"   , "W/m2"   ,   &   ! solar flux
-         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
-      CALL histdef( id_i, "soicecov", "Ice fraction"          , "[0,1]"  ,   &   ! fr_i
-         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
-      CALL histdef( id_i, "sozotaux", "Zonal Wind Stress"     , "N/m2"   ,   &   ! i-wind stress
-         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
-      CALL histdef( id_i, "sometauy", "Meridional Wind Stress", "N/m2"   ,   &   ! j-wind stress
-         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
-      IF( .NOT.ln_linssh ) THEN
-         CALL histdef( id_i, "vovvldep", "T point depth"         , "m"      , &   ! t-point depth
-            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
-         CALL histdef( id_i, "vovvle3t", "T point thickness"     , "m"      , &   ! t-point depth
-            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
-      ENDIF
-      !
+         CALL iom_rstput( 0, 0, inum,  'ahtu', ahtu              )    ! aht at u-point
+         CALL iom_rstput( 0, 0, inum,  'ahtv', ahtv              )    ! aht at v-point
+      ENDIF
+      IF( ALLOCATED(ahmt) ) THEN
+         CALL iom_rstput( 0, 0, inum,  'ahmt', ahmt              )    ! ahmt at u-point
+         CALL iom_rstput( 0, 0, inum,  'ahmf', ahmf              )    ! ahmf at v-point
+      ENDIF
+      CALL iom_rstput( 0, 0, inum, 'sowaflup', emp - rnf         )    ! freshwater budget
+      CALL iom_rstput( 0, 0, inum, 'sohefldo', qsr + qns         )    ! total heat flux
+      CALL iom_rstput( 0, 0, inum, 'soshfldo', qsr               )    ! solar heat flux
+      CALL iom_rstput( 0, 0, inum, 'soicecov', fr_i              )    ! ice fraction
+      CALL iom_rstput( 0, 0, inum, 'sozotaux', utau              )    ! i-wind stress
+      CALL iom_rstput( 0, 0, inum, 'sometauy', vtau              )    ! j-wind stress
+      IF(  .NOT.ln_linssh  ) THEN             
+         CALL iom_rstput( 0, 0, inum, 'vovvldep', gdept_n        )    !  T-cell depth 
+         CALL iom_rstput( 0, 0, inum, 'vovvle3t', e3t_n          )    !  T-cell thickness  
+      END IF
       IF( ln_wave .AND. ln_sdw ) THEN
-         CALL histdef( id_i, "sdzocrtx", "Stokes Drift Zonal"    , "m/s"    , &   ! StokesDrift zonal current
-            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
-         CALL histdef( id_i, "sdmecrty", "Stokes Drift Merid"    , "m/s"    , &   ! StokesDrift meridonal current
-            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
-         CALL histdef( id_i, "sdvecrtz", "Stokes Drift Vert"     , "m/s"    , &   ! StokesDrift vertical current
-            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
-      ENDIF
-
+         CALL iom_rstput( 0, 0, inum, 'sdzocrtx', usd            )    ! now StokesDrift i-velocity
+         CALL iom_rstput( 0, 0, inum, 'sdmecrty', vsd            )    ! now StokesDrift j-velocity
+         CALL iom_rstput( 0, 0, inum, 'sdvecrtz', wsd            )    ! now StokesDrift k-velocity
+      ENDIF
+ 
 #if defined key_si3
       IF( nn_ice == 2 ) THEN   ! condition needed in case agrif + ice-model but no-ice in child grid
-         CALL ice_wri_state( kt, id_i, nh_i )
-      ENDIF
-#else
-      CALL histend( id_i, snc4chunks=snc4set )
+         CALL ice_wri_state( inum )
+      ENDIF
 #endif
-
-      ! 2. Start writing data
-      ! ---------------------
-      ! idex(1) est utilise ssi l'avant dernier argument est diffferent de 
-      ! la taille du tableau en sortie. Dans ce cas , l'avant dernier argument
-      ! donne le nombre d'elements, et idex la liste des indices a sortir
-      idex(1) = 1   ! init to avoid compil warning
-
-      ! Write all fields on T grid
-      CALL histwrite( id_i, "votemper", kt, tsn(:,:,:,jp_tem), jpi*jpj*jpk, idex )    ! now temperature
-      CALL histwrite( id_i, "vosaline", kt, tsn(:,:,:,jp_sal), jpi*jpj*jpk, idex )    ! now salinity
-      CALL histwrite( id_i, "sossheig", kt, sshn             , jpi*jpj    , idex )    ! sea surface height
-      CALL histwrite( id_i, "vozocrtx", kt, un               , jpi*jpj*jpk, idex )    ! now i-velocity
-      CALL histwrite( id_i, "vomecrty", kt, vn               , jpi*jpj*jpk, idex )    ! now j-velocity
-      CALL histwrite( id_i, "vovecrtz", kt, wn               , jpi*jpj*jpk, idex )    ! now k-velocity
-      !
-      IF( ALLOCATED(ahtu) ) THEN
-         CALL histwrite( id_i, "ahtu"    , kt, ahtu             , jpi*jpj*jpk, idex )    ! aht at u-point
-         CALL histwrite( id_i, "ahtv"    , kt, ahtv             , jpi*jpj*jpk, idex )    !  -  at v-point
-      ENDIF
-      IF( ALLOCATED(ahmt) ) THEN
-         CALL histwrite( id_i, "ahmt"    , kt, ahmt             , jpi*jpj*jpk, idex )    ! ahm at t-point
-         CALL histwrite( id_i, "ahmf"    , kt, ahmf             , jpi*jpj*jpk, idex )    !  -  at f-point
-      ENDIF
-      !
-      CALL histwrite( id_i, "sowaflup", kt, emp - rnf        , jpi*jpj    , idex )    ! freshwater budget
-      CALL histwrite( id_i, "sohefldo", kt, qsr + qns        , jpi*jpj    , idex )    ! total heat flux
-      CALL histwrite( id_i, "soshfldo", kt, qsr              , jpi*jpj    , idex )    ! solar heat flux
-      CALL histwrite( id_i, "soicecov", kt, fr_i             , jpi*jpj    , idex )    ! ice fraction
-      CALL histwrite( id_i, "sozotaux", kt, utau             , jpi*jpj    , idex )    ! i-wind stress
-      CALL histwrite( id_i, "sometauy", kt, vtau             , jpi*jpj    , idex )    ! j-wind stress
-
-      IF(  .NOT.ln_linssh  ) THEN             
-         CALL histwrite( id_i, "vovvldep", kt, gdept_n(:,:,:), jpi*jpj*jpk, idex )!  T-cell depth 
-         CALL histwrite( id_i, "vovvle3t", kt, e3t_n (:,:,:) , jpi*jpj*jpk, idex )!  T-cell thickness  
-      END IF 
- 
-      IF( ln_wave .AND. ln_sdw ) THEN
-         CALL histwrite( id_i, "sdzocrtx", kt, usd           , jpi*jpj*jpk, idex)     ! now StokesDrift i-velocity
-         CALL histwrite( id_i, "sdmecrty", kt, vsd           , jpi*jpj*jpk, idex)     ! now StokesDrift j-velocity
-         CALL histwrite( id_i, "sdvecrtz", kt, wsd           , jpi*jpj*jpk, idex)     ! now StokesDrift k-velocity
-      ENDIF
-
-      ! 3. Close the file
-      ! -----------------
-      CALL histclo( id_i )
-#if ! defined key_iomput
-      IF( ninist /= 1  ) THEN
-         CALL histclo( nid_T )
-         CALL histclo( nid_U )
-         CALL histclo( nid_V )
-         CALL histclo( nid_W )
-      ENDIF
-#endif
+      !
+      CALL iom_close( inum )
       ! 
    END SUBROUTINE dia_wri_state
Index: NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/tests/CANAL/MY_SRC/diawri.F90
===================================================================
--- NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/tests/CANAL/MY_SRC/diawri.F90	(revision 10345)
+++ NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/tests/CANAL/MY_SRC/diawri.F90	(revision 10358)
@@ -52,4 +52,5 @@
 
 #if defined key_si3
+   USE ice 
    USE icewri 
 #endif
@@ -89,5 +90,4 @@
    !!   'key_iomput'                                        use IOM library
    !!----------------------------------------------------------------------
-
    INTEGER FUNCTION dia_wri_alloc()
       !
@@ -96,5 +96,5 @@
    END FUNCTION dia_wri_alloc
 
-
+   
    SUBROUTINE dia_wri( kt )
       !!---------------------------------------------------------------------
@@ -123,5 +123,5 @@
       ! Output the initial state and forcings
       IF( ninist == 1 ) THEN                       
-         CALL dia_wri_state( 'output.init', kt )
+         CALL dia_wri_state( 'output.init' )
          ninist = 0
       ENDIF
@@ -129,6 +129,6 @@
       ! Output of initial vertical scale factor
       CALL iom_put("e3t_0", e3t_0(:,:,:) )
-      CALL iom_put("e3u_0", e3t_0(:,:,:) )
-      CALL iom_put("e3v_0", e3t_0(:,:,:) )
+      CALL iom_put("e3u_0", e3u_0(:,:,:) )
+      CALL iom_put("e3v_0", e3v_0(:,:,:) )
       !
       CALL iom_put( "e3t" , e3t_n(:,:,:) )
@@ -505,5 +505,5 @@
    !!   Default option                                  use IOIPSL  library
    !!----------------------------------------------------------------------
-   
+
    INTEGER FUNCTION dia_wri_alloc()
       !!----------------------------------------------------------------------
@@ -519,6 +519,6 @@
       !
    END FUNCTION dia_wri_alloc
-  
-  
+
+   
    SUBROUTINE dia_wri( kt )
       !!---------------------------------------------------------------------
@@ -551,5 +551,5 @@
       !
       IF( ninist == 1 ) THEN     !==  Output the initial state and forcings  ==!
-         CALL dia_wri_state( 'output.init', kt )
+         CALL dia_wri_state( 'output.init' )
          ninist = 0
       ENDIF
@@ -974,5 +974,5 @@
 #endif
 
-   SUBROUTINE dia_wri_state( cdfile_name, kt )
+   SUBROUTINE dia_wri_state( cdfile_name )
       !!---------------------------------------------------------------------
       !!                 ***  ROUTINE dia_wri_state  ***
@@ -988,158 +988,56 @@
       !!----------------------------------------------------------------------
       CHARACTER (len=* ), INTENT( in ) ::   cdfile_name      ! name of the file created
-      INTEGER           , INTENT( in ) ::   kt               ! ocean time-step index
-      !! 
-      CHARACTER (len=32) :: clname
-      CHARACTER (len=40) :: clop
-      INTEGER  ::   id_i , nz_i, nh_i       
-      INTEGER, DIMENSION(1) ::   idex             ! local workspace
-      REAL(wp) ::   zsto, zout, zmax, zjulian
+      !!
+      INTEGER :: inum
       !!----------------------------------------------------------------------
       ! 
-      ! 0. Initialisation
-      ! -----------------
-
-      ! Define name, frequency of output and means
-      clname = cdfile_name
-      IF( .NOT. Agrif_Root() ) clname = TRIM(Agrif_CFixed())//'_'//TRIM(clname)
-      zsto = rdt
-      clop = "inst(x)"           ! no use of the mask value (require less cpu time)
-      zout = rdt
-      zmax = ( nitend - nit000 + 1 ) * rdt
-
       IF(lwp) WRITE(numout,*)
       IF(lwp) WRITE(numout,*) 'dia_wri_state : single instantaneous ocean state'
       IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~   and forcing fields file created '
-      IF(lwp) WRITE(numout,*) '                and named :', clname, '.nc'
-
-
-      ! 1. Define NETCDF files and fields at beginning of first time step
-      ! -----------------------------------------------------------------
-
-      ! Compute julian date from starting date of the run
-      CALL ymds2ju( nyear, nmonth, nday, rdt, zjulian )         ! time axis 
-      zjulian = zjulian - adatrj   !   set calendar origin to the beginning of the experiment
-      CALL histbeg( clname, jpi, glamt, jpj, gphit,   &
-          1, jpi, 1, jpj, nit000-1, zjulian, rdt, nh_i, id_i, domain_id=nidom, snc4chunks=snc4set ) ! Horizontal grid : glamt and gphit
-      CALL histvert( id_i, "deptht", "Vertical T levels",   &    ! Vertical grid : gdept
-          "m", jpk, gdept_1d, nz_i, "down")
-
-      ! Declare all the output fields as NetCDF variables
-
-      CALL histdef( id_i, "vosaline", "Salinity"              , "PSU"    ,   &   ! salinity
-         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
-      CALL histdef( id_i, "votemper", "Temperature"           , "C"      ,   &   ! temperature
-         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
-      CALL histdef( id_i, "sossheig", "Sea Surface Height"    , "m"      ,   &  ! ssh
-         &          jpi, jpj, nh_i, 1  , 1, 1  , nz_i, 32, clop, zsto, zout )
-      CALL histdef( id_i, "vozocrtx", "Zonal Current"         , "m/s"    ,   &   ! zonal current
-         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
-      CALL histdef( id_i, "vomecrty", "Meridional Current"    , "m/s"    ,   &   ! meridonal current
-         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
-      CALL histdef( id_i, "vovecrtz", "Vertical Velocity"     , "m/s"    ,   &   ! vertical current
-         &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
-         !
+      IF(lwp) WRITE(numout,*) '                and named :', cdfile_name, '...nc'
+
+#if defined key_si3
+     CALL iom_open( TRIM(cdfile_name), inum, ldwrt = .TRUE., kiolib = jprstlib, kdlev = jpl )
+#else
+     CALL iom_open( TRIM(cdfile_name), inum, ldwrt = .TRUE., kiolib = jprstlib )
+#endif
+
+      CALL iom_rstput( 0, 0, inum, 'votemper', tsn(:,:,:,jp_tem) )    ! now temperature
+      CALL iom_rstput( 0, 0, inum, 'vosaline', tsn(:,:,:,jp_sal) )    ! now salinity
+      CALL iom_rstput( 0, 0, inum, 'sossheig', sshn              )    ! sea surface height
+      CALL iom_rstput( 0, 0, inum, 'vozocrtx', un                )    ! now i-velocity
+      CALL iom_rstput( 0, 0, inum, 'vomecrty', vn                )    ! now j-velocity
+      CALL iom_rstput( 0, 0, inum, 'vovecrtz', wn                )    ! now k-velocity
       IF( ALLOCATED(ahtu) ) THEN
-         CALL histdef( id_i, "ahtu"    , "u-eddy diffusivity"    , "m2/s"    ,   &   ! zonal current
-            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
-         CALL histdef( id_i, "ahtv"    , "v-eddy diffusivity"    , "m2/s"    ,   &   ! meridonal current
-            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
-      ENDIF
-      IF( ALLOCATED(ahmt) ) THEN 
-         CALL histdef( id_i, "ahmt"    , "t-eddy viscosity"      , "m2/s"    ,   &   ! zonal current
-            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
-         CALL histdef( id_i, "ahmf"    , "f-eddy viscosity"      , "m2/s"    ,   &   ! meridonal current
-            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 
-      ENDIF
-         !
-      CALL histdef( id_i, "sowaflup", "Net Upward Water Flux" , "Kg/m2/S",   &   ! net freshwater 
-         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
-      CALL histdef( id_i, "sohefldo", "Net Downward Heat Flux", "W/m2"   ,   &   ! net heat flux
-         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
-      CALL histdef( id_i, "soshfldo", "Shortwave Radiation"   , "W/m2"   ,   &   ! solar flux
-         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
-      CALL histdef( id_i, "soicecov", "Ice fraction"          , "[0,1]"  ,   &   ! fr_i
-         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
-      CALL histdef( id_i, "sozotaux", "Zonal Wind Stress"     , "N/m2"   ,   &   ! i-wind stress
-         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
-      CALL histdef( id_i, "sometauy", "Meridional Wind Stress", "N/m2"   ,   &   ! j-wind stress
-         &          jpi, jpj, nh_i, 1  , 1, 1  , -99 , 32, clop, zsto, zout )
-      IF( .NOT.ln_linssh ) THEN
-         CALL histdef( id_i, "vovvldep", "T point depth"         , "m"      , &   ! t-point depth
-            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
-         CALL histdef( id_i, "vovvle3t", "T point thickness"     , "m"      , &   ! t-point depth
-            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
-      ENDIF
-      !
+         CALL iom_rstput( 0, 0, inum,  'ahtu', ahtu              )    ! aht at u-point
+         CALL iom_rstput( 0, 0, inum,  'ahtv', ahtv              )    ! aht at v-point
+      ENDIF
+      IF( ALLOCATED(ahmt) ) THEN
+         CALL iom_rstput( 0, 0, inum,  'ahmt', ahmt              )    ! ahmt at u-point
+         CALL iom_rstput( 0, 0, inum,  'ahmf', ahmf              )    ! ahmf at v-point
+      ENDIF
+      CALL iom_rstput( 0, 0, inum, 'sowaflup', emp - rnf         )    ! freshwater budget
+      CALL iom_rstput( 0, 0, inum, 'sohefldo', qsr + qns         )    ! total heat flux
+      CALL iom_rstput( 0, 0, inum, 'soshfldo', qsr               )    ! solar heat flux
+      CALL iom_rstput( 0, 0, inum, 'soicecov', fr_i              )    ! ice fraction
+      CALL iom_rstput( 0, 0, inum, 'sozotaux', utau              )    ! i-wind stress
+      CALL iom_rstput( 0, 0, inum, 'sometauy', vtau              )    ! j-wind stress
+      IF(  .NOT.ln_linssh  ) THEN             
+         CALL iom_rstput( 0, 0, inum, 'vovvldep', gdept_n        )    !  T-cell depth 
+         CALL iom_rstput( 0, 0, inum, 'vovvle3t', e3t_n          )    !  T-cell thickness  
+      END IF
       IF( ln_wave .AND. ln_sdw ) THEN
-         CALL histdef( id_i, "sdzocrtx", "Stokes Drift Zonal"    , "m/s"    , &   ! StokesDrift zonal current
-            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
-         CALL histdef( id_i, "sdmecrty", "Stokes Drift Merid"    , "m/s"    , &   ! StokesDrift meridonal current
-            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
-         CALL histdef( id_i, "sdvecrtz", "Stokes Drift Vert"     , "m/s"    , &   ! StokesDrift vertical current
-            &          jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout )
-      ENDIF
-
+         CALL iom_rstput( 0, 0, inum, 'sdzocrtx', usd            )    ! now StokesDrift i-velocity
+         CALL iom_rstput( 0, 0, inum, 'sdmecrty', vsd            )    ! now StokesDrift j-velocity
+         CALL iom_rstput( 0, 0, inum, 'sdvecrtz', wsd            )    ! now StokesDrift k-velocity
+      ENDIF
+ 
 #if defined key_si3
       IF( nn_ice == 2 ) THEN   ! condition needed in case agrif + ice-model but no-ice in child grid
-         CALL ice_wri_state( kt, id_i, nh_i )
-      ENDIF
-#else
-      CALL histend( id_i, snc4chunks=snc4set )
+         CALL ice_wri_state( inum )
+      ENDIF
 #endif
-
-      ! 2. Start writing data
-      ! ---------------------
-      ! idex(1) est utilise ssi l'avant dernier argument est diffferent de 
-      ! la taille du tableau en sortie. Dans ce cas , l'avant dernier argument
-      ! donne le nombre d'elements, et idex la liste des indices a sortir
-      idex(1) = 1   ! init to avoid compil warning
-
-      ! Write all fields on T grid
-      CALL histwrite( id_i, "votemper", kt, tsn(:,:,:,jp_tem), jpi*jpj*jpk, idex )    ! now temperature
-      CALL histwrite( id_i, "vosaline", kt, tsn(:,:,:,jp_sal), jpi*jpj*jpk, idex )    ! now salinity
-      CALL histwrite( id_i, "sossheig", kt, sshn             , jpi*jpj    , idex )    ! sea surface height
-      CALL histwrite( id_i, "vozocrtx", kt, un               , jpi*jpj*jpk, idex )    ! now i-velocity
-      CALL histwrite( id_i, "vomecrty", kt, vn               , jpi*jpj*jpk, idex )    ! now j-velocity
-      CALL histwrite( id_i, "vovecrtz", kt, wn               , jpi*jpj*jpk, idex )    ! now k-velocity
-      !
-      IF( ALLOCATED(ahtu) ) THEN
-         CALL histwrite( id_i, "ahtu"    , kt, ahtu             , jpi*jpj*jpk, idex )    ! aht at u-point
-         CALL histwrite( id_i, "ahtv"    , kt, ahtv             , jpi*jpj*jpk, idex )    !  -  at v-point
-      ENDIF
-      IF( ALLOCATED(ahmt) ) THEN
-         CALL histwrite( id_i, "ahmt"    , kt, ahmt             , jpi*jpj*jpk, idex )    ! ahm at t-point
-         CALL histwrite( id_i, "ahmf"    , kt, ahmf             , jpi*jpj*jpk, idex )    !  -  at f-point
-      ENDIF
-      !
-      CALL histwrite( id_i, "sowaflup", kt, emp - rnf        , jpi*jpj    , idex )    ! freshwater budget
-      CALL histwrite( id_i, "sohefldo", kt, qsr + qns        , jpi*jpj    , idex )    ! total heat flux
-      CALL histwrite( id_i, "soshfldo", kt, qsr              , jpi*jpj    , idex )    ! solar heat flux
-      CALL histwrite( id_i, "soicecov", kt, fr_i             , jpi*jpj    , idex )    ! ice fraction
-      CALL histwrite( id_i, "sozotaux", kt, utau             , jpi*jpj    , idex )    ! i-wind stress
-      CALL histwrite( id_i, "sometauy", kt, vtau             , jpi*jpj    , idex )    ! j-wind stress
-
-      IF(  .NOT.ln_linssh  ) THEN             
-         CALL histwrite( id_i, "vovvldep", kt, gdept_n(:,:,:), jpi*jpj*jpk, idex )!  T-cell depth 
-         CALL histwrite( id_i, "vovvle3t", kt, e3t_n (:,:,:) , jpi*jpj*jpk, idex )!  T-cell thickness  
-      END IF 
- 
-      IF( ln_wave .AND. ln_sdw ) THEN
-         CALL histwrite( id_i, "sdzocrtx", kt, usd           , jpi*jpj*jpk, idex)     ! now StokesDrift i-velocity
-         CALL histwrite( id_i, "sdmecrty", kt, vsd           , jpi*jpj*jpk, idex)     ! now StokesDrift j-velocity
-         CALL histwrite( id_i, "sdvecrtz", kt, wsd           , jpi*jpj*jpk, idex)     ! now StokesDrift k-velocity
-      ENDIF
-
-      ! 3. Close the file
-      ! -----------------
-      CALL histclo( id_i )
-#if ! defined key_iomput
-      IF( ninist /= 1  ) THEN
-         CALL histclo( nid_T )
-         CALL histclo( nid_U )
-         CALL histclo( nid_V )
-         CALL histclo( nid_W )
-      ENDIF
-#endif
+      !
+      CALL iom_close( inum )
       ! 
    END SUBROUTINE dia_wri_state
Index: NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/tests/CANAL/MY_SRC/stpctl.F90
===================================================================
--- NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/tests/CANAL/MY_SRC/stpctl.F90	(revision 10345)
+++ NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/tests/CANAL/MY_SRC/stpctl.F90	(revision 10358)
@@ -32,5 +32,6 @@
    PUBLIC stp_ctl           ! routine called by step.F90
 
-   INTEGER  ::   idrun, idtime, idssh, idu, ids, istatus
+   INTEGER  ::   idrun, idtime, idssh, idu, ids1, ids2, istatus
+   LOGICAL  ::   lsomeoce
    !!----------------------------------------------------------------------
    !! NEMO/OCE 4.0 , NEMO Consortium (2018)
@@ -60,12 +61,10 @@
       INTEGER, INTENT(inout) ::   kindic   ! error indicator
       !!
-      INTEGER  ::   ji, jj, jk             ! dummy loop indices
-      INTEGER  ::   iih, ijh               ! local integers
-      INTEGER  ::   iiu, iju, iku          !   -       -
-      INTEGER  ::   iis, ijs, iks          !   -       -
-      REAL(wp) ::   zzz                    ! local real 
-      INTEGER , DIMENSION(3) ::   ilocu, ilocs
-      INTEGER , DIMENSION(2) ::   iloch
-      REAL(wp), DIMENSION(4) ::   zmax
+      INTEGER                ::   ji, jj, jk          ! dummy loop indices
+      INTEGER, DIMENSION(2)  ::   ih                  ! min/max loc indices
+      INTEGER, DIMENSION(3)  ::   iu, is1, is2        ! min/max loc indices
+      REAL(wp)               ::   zzz                 ! local real 
+      REAL(wp), DIMENSION(5) ::   zmax
+      CHARACTER(len=20) :: clname
       !!----------------------------------------------------------------------
       !
@@ -75,19 +74,22 @@
          WRITE(numout,*) '~~~~~~~'
          !                                ! open time.step file
-         CALL ctl_opn( numstp, 'time.step', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea )
+         IF( lwm ) CALL ctl_opn( numstp, 'time.step', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea )
          !                                ! open run.stat file
-         CALL ctl_opn( numrun, 'run.stat', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea )
-
-         IF( lwm ) THEN
-            istatus = NF90_CREATE( 'run.stat.nc', NF90_CLOBBER, idrun )
+         IF( ln_ctl .AND. lwm ) THEN
+            CALL ctl_opn( numrun, 'run.stat', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea )
+            clname = 'run.stat.nc'
+            IF( .NOT. Agrif_Root() )   clname = TRIM(Agrif_CFixed())//"_"//TRIM(clname)
+            istatus = NF90_CREATE( TRIM(clname), NF90_CLOBBER, idrun )
             istatus = NF90_DEF_DIM( idrun, 'time', NF90_UNLIMITED, idtime )
             istatus = NF90_DEF_VAR( idrun, 'abs_ssh_max', NF90_DOUBLE, (/ idtime /), idssh )
-            istatus = NF90_DEF_VAR( idrun,   'abs_u_max', NF90_DOUBLE, (/ idtime /), idu )
+            istatus = NF90_DEF_VAR( idrun,   'abs_u_max', NF90_DOUBLE, (/ idtime /), idu   )
+            istatus = NF90_DEF_VAR( idrun,       's_min', NF90_DOUBLE, (/ idtime /), ids1  )
+            istatus = NF90_DEF_VAR( idrun,       's_max', NF90_DOUBLE, (/ idtime /), ids2  )
             istatus = NF90_ENDDEF(idrun)
          ENDIF
-         
       ENDIF
+      IF( kt == nit000 )   lsomeoce = COUNT( ssmask(:,:) == 1._wp ) > 0
       !
-      IF(lwp) THEN                        !==  current time step  ==!   ("time.step" file)
+      IF(lwm) THEN                        !==  current time step  ==!   ("time.step" file)
          WRITE ( numstp, '(1x, i8)' )   kt
          REWIND( numstp )
@@ -101,67 +103,72 @@
       ENDIF
       zmax(2) = MAXVAL(  ABS( un(:,:,:) )  )                                  ! velocity max (zonal only)
-!     zmax(3) = MAXVAL( -tsn(:,:,:,jp_sal) , mask = tmask(:,:,:) == 1._wp )   ! minus salinity max
-      zmax(3) = 0.0_wp
-      zmax(4) = REAL( nstop , wp )                                            ! stop indicator
+      zmax(3) = MAXVAL( -tsn(:,:,:,jp_sal) , mask = tmask(:,:,:) == 1._wp )   ! minus salinity max
+      zmax(4) = MAXVAL(  tsn(:,:,:,jp_sal) , mask = tmask(:,:,:) == 1._wp )   !       salinity max
+      zmax(5) = REAL( nstop , wp )                                            ! stop indicator
       !
-      IF( lk_mpp ) THEN
-         CALL mpp_max_multiple( zmax(:), 4 )    ! max over the global domain
-         !
-         nstop = INT( zmax(4) )                 ! nstop indicator sheared among all local domains
+      IF( lk_mpp .AND. ln_ctl ) THEN
+         CALL mpp_max( "stpctl", zmax )          ! max over the global domain
+         nstop = NINT( zmax(5) )                 ! nstop indicator sheared among all local domains
       ENDIF
-      !
-      IF( MOD( kt, nwrite ) == 1 .AND. lwp ) THEN
-         WRITE(numout,*) ' ==>> time-step= ', kt, ' |ssh| max: ',   zmax(1), ' |U| max: ', zmax(2)
-      ENDIF
-      !
-      IF (  zmax(1) >  50._wp .OR.   &                     ! too large sea surface height ( > 10 m)
-         &  zmax(2) >  20._wp .OR.   &                     ! too large velocity ( > 10 m/s)
-!!$         &  zmax(3) >=  0._wp .OR.   &                     ! negative or zero sea surface salinity
-         &  ISNAN( zmax(1) + zmax(2) + zmax(3) )  ) THEN   ! NaN encounter in the tests
-         IF( lk_mpp ) THEN
-            CALL mpp_maxloc( 'stpctl', ABS(sshn)        , ssmask(:,:)  , zzz, iih, ijh )
-            CALL mpp_maxloc( 'stpctl', ABS(un)          , umask (:,:,:), zzz, iiu, iju, iku )
-!           CALL mpp_minloc( 'stpctl', tsn(:,:,:,jp_sal), tmask (:,:,:), zzz, iis, ijs, iks )
-         ELSE
-            iloch = MINLOC( ABS( sshn(:,:)   )                               )
-            ilocu = MAXLOC( ABS( un  (:,:,:) )                               )
-!           ilocs = MINLOC( tsn(:,:,:,jp_sal) , mask = tmask(:,:,:) == 1._wp )
-            iih = iloch(1) + nimpp - 1   ;   ijh = iloch(2) + njmpp - 1
-            iiu = ilocu(1) + nimpp - 1   ;   iju = ilocu(2) + njmpp - 1   ;   iku = ilocu(3)
-!           iis = ilocs(1) + nimpp - 1   ;   ijs = ilocs(2) + njmpp - 1   ;   iks = ilocu(3)
-         ENDIF
-         IF(lwp) THEN
-            WRITE(numout,cform_err)
-            WRITE(numout,*) ' stp_ctl: |ssh| > 50 m   or   |U| > 20 m/s   or   NaN encounter in the tests'
-            WRITE(numout,*) ' ======= '
-            WRITE(numout,9100) kt,   zmax(1), iih, ijh
-            WRITE(numout,9200) kt,   zmax(2), iiu, iju, iku
-!!$            WRITE(numout,9300) kt, - zmax(3), iis, ijs, iks
-            WRITE(numout,*)
-            WRITE(numout,*) '          output of last computed fields in output.abort.nc file'
-         ENDIF
-         kindic = -3
-         !
-         nstop = nstop + 1                            ! increase nstop by 1 (on all local domains)
-         CALL dia_wri_state( 'output.abort', kt )     ! create an output.abort file
-         !
-      ENDIF
-9100  FORMAT (' kt=',i8,'   |ssh| max: ',1pg11.4,', at  i j  : ',2i5)
-9200  FORMAT (' kt=',i8,'   |U|   max: ',1pg11.4,', at  i j k: ',3i5)
-9300  FORMAT (' kt=',i8,'   S     min: ',1pg11.4,', at  i j  : ',2i5)
-      !
-      !                                            !==  run statistics  ==!   ("run.stat" file)
-!     IF(lwp) WRITE(numrun,9400) kt, zmax(1), zmax(2), - zmax(3)
-      IF(lwp) WRITE(numrun,9400) kt, zmax(1), zmax(2)
-      IF( lwm ) THEN
+      !                                   !==  run statistics  ==!   ("run.stat" files)
+      IF( ln_ctl .AND. lwm ) THEN
+         WRITE(numrun,9500) kt, zmax(1), zmax(2), -zmax(3), zmax(4)
          istatus = NF90_PUT_VAR( idrun, idssh, (/ zmax(1)/), (/kt/), (/1/) )
          istatus = NF90_PUT_VAR( idrun,   idu, (/ zmax(2)/), (/kt/), (/1/) )
-!        istatus = NF90_PUT_VAR( idrun,   ids, (/-zmax(3)/), (/kt/), (/1/) )
+         istatus = NF90_PUT_VAR( idrun,  ids1, (/-zmax(3)/), (/kt/), (/1/) )
+         istatus = NF90_PUT_VAR( idrun,  ids2, (/ zmax(4)/), (/kt/), (/1/) )
          IF( MOD( kt , 100 ) == 0 ) istatus = NF90_SYNC(idrun)
          IF( kt == nitend         ) istatus = NF90_CLOSE(idrun)
       END IF
+      !                                   !==  error handling  ==!
+      IF( ( ln_ctl .OR. lsomeoce ) .AND. (   &             ! have use mpp_max (because ln_ctl=.T.) or contains some ocean points
+         &  zmax(1) >   50._wp .OR.   &                    ! too large sea surface height ( > 50 m )
+         &  zmax(2) >   20._wp .OR.   &                    ! too large velocity ( > 20 m/s)
+         &  zmax(3) >= 100._wp .OR.   &                    ! too small sea surface salinity ( < -100 )
+         &  zmax(4) >= 100._wp .OR.   &                    ! too large sea surface salinity ( > 100 )
+         &  zmax(4) < -100._wp .OR.   &                    ! too large sea surface salinity (keep this line for sea-ice)
+         &  ISNAN( zmax(1) + zmax(2) + zmax(3) )  ) THEN   ! NaN encounter in the tests
+         IF( lk_mpp .AND. ln_ctl ) THEN
+            CALL mpp_maxloc( 'stpctl', ABS(sshn)        , ssmask(:,:)  , zzz, ih  )
+            CALL mpp_maxloc( 'stpctl', ABS(un)          , umask (:,:,:), zzz, iu  )
+            CALL mpp_minloc( 'stpctl', tsn(:,:,:,jp_sal), tmask (:,:,:), zzz, is1 )
+            CALL mpp_maxloc( 'stpctl', tsn(:,:,:,jp_sal), tmask (:,:,:), zzz, is2 )
+         ELSE
+            ih(:)  = MAXLOC( ABS( sshn(:,:)   )                              ) + (/ nimpp - 1, njmpp - 1    /)
+            iu(:)  = MAXLOC( ABS( un  (:,:,:) )                              ) + (/ nimpp - 1, njmpp - 1, 0 /)
+            is1(:) = MINLOC( tsn(:,:,:,jp_sal), mask = tmask(:,:,:) == 1._wp ) + (/ nimpp - 1, njmpp - 1, 0 /)
+            is2(:) = MAXLOC( tsn(:,:,:,jp_sal), mask = tmask(:,:,:) == 1._wp ) + (/ nimpp - 1, njmpp - 1, 0 /)
+         ENDIF
+         IF( numout == 6 )   &   ! force to open ocean.output file
+            CALL ctl_opn( numout, 'ocean.output', 'APPEND', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE. )
+
+         WRITE(numout,cform_err)
+         WRITE(numout,*) ' stp_ctl: |ssh| > 50 m  or  |U| > 20 m/s  or  S <= -100  or  S >= 100  or  NaN encounter in the tests'
+         WRITE(numout,*) ' ======= '
+         IF( lk_mpp .AND. .NOT. ln_ctl ) WRITE(numout,*) 'E R R O R message from sub-domain: ', narea
+         WRITE(numout,9100) kt,   zmax(1), ih(1) , ih(2)
+         WRITE(numout,9200) kt,   zmax(2), iu(1) , iu(2) , iu(3)
+         WRITE(numout,9300) kt, - zmax(3), is1(1), is1(2), is1(3)
+         WRITE(numout,9400) kt,   zmax(4), is2(1), is2(2), is2(3)
+         WRITE(numout,*)
+         WRITE(numout,*) '          output of last computed fields in output.abort.nc file'
+         
+         CALL dia_wri_state( 'output.abort' )     ! create an output.abort file
+         
+         IF( ln_ctl ) THEN
+            kindic = -3
+            nstop = nstop + 1                            ! increase nstop by 1 (on all local domains)
+         ELSE
+            CALL ctl_stop()
+            CALL mppstop(ld_force_abort = .true.)
+         ENDIF
+         !
+      ENDIF
       !
-!9400  FORMAT(' it :', i8, '    |ssh|_max: ', D23.16, ' |U|_max: ', D23.16,' S_min: ', D23.16)
-9400  FORMAT(' it :', i8, '    |ssh|_max: ', D23.16, ' |U|_max: ', D23.16)
+9100  FORMAT (' kt=',i8,'   |ssh| max: ',1pg11.4,', at  i j  : ',2i5)
+9200  FORMAT (' kt=',i8,'   |U|   max: ',1pg11.4,', at  i j k: ',3i5)
+9300  FORMAT (' kt=',i8,'   S     min: ',1pg11.4,', at  i j k: ',3i5)
+9400  FORMAT (' kt=',i8,'   S     max: ',1pg11.4,', at  i j k: ',3i5)
+9500  FORMAT(' it :', i8, '    |ssh|_max: ', D23.16, ' |U|_max: ', D23.16,' S_min: ', D23.16,' S_max: ', D23.16)
       !
    END SUBROUTINE stp_ctl
