Index: NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/cfgs/SHARED/field_def_nemo-ice.xml
===================================================================
--- NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/cfgs/SHARED/field_def_nemo-ice.xml	(revision 11370)
+++ NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/cfgs/SHARED/field_def_nemo-ice.xml	(revision 11371)
@@ -32,17 +32,17 @@
  
 	  <!-- general fields -->
-          <field id="icemass"      long_name="Sea-ice mass per area"                                   standard_name="sea_ice_amount"                            unit="kg/m2"/>
+          <field id="icemass"      long_name="Sea-ice mass per gridcell area"                          standard_name="sea_ice_amount"                            unit="kg/m2"/>
           <field id="snwmass"      long_name="Snow mass per area"                                      standard_name="liquid_water_content_of_surface_snow"      unit="kg/m2"/>
-          <field id="iceconc"      long_name="Sea-ice area fraction"                                   standard_name="sea_ice_area_fraction"                     unit=""    />
-          <field id="icevolu"      long_name="Sea-ice volume per area"                                 standard_name="sea_ice_thickness"                         unit="m"   />
-          <field id="icethic"      long_name="Sea-ice thickness per area"                              standard_name="sea_ice_thickness"                         unit="m"   detect_missing_value="true" />
-          <field id="snwthic"      long_name="Snow thickness per area"                                 standard_name="snow_thickness"                            unit="m"   detect_missing_value="true" />
-          <field id="icebrv"       long_name="brine volume"                                                                                                      unit="%"   />
+          <field id="iceconc"      long_name="Sea-ice area fraction"                                   standard_name="sea_ice_area_fraction"                     unit=""     />
+          <field id="icevolu"      long_name="Sea-ice volume per area"                                 standard_name="sea_ice_thickness"                         unit="m"    />
+          <field id="icethic"      long_name="Sea-ice thickness"                                       standard_name="sea_ice_thickness"                         unit="m"    />
+          <field id="snwthic"      long_name="Snow thickness"                                          standard_name="snow_thickness"                            unit="m"    />
+          <field id="icebrv"       long_name="brine volume"                                                                                                      unit="%"    />
           <field id="iceage"       long_name="ice age"                                                                                                           unit="days" detect_missing_value="true"/>
-          <field id="icehnew"      long_name="frazil ice collection thickness"                                                                                   unit="m"   />
-          <field id="snwvolu"      long_name="snow volume"                                                                                                       unit="m"   />
-          <field id="icefrb"       long_name="Sea-ice freeboard"                                       standard_name="sea_ice_freeboard"                         unit="m"   detect_missing_value="true" />
+          <field id="icehnew"      long_name="frazil ice collection thickness"                                                                                   unit="m"    />
+          <field id="snwvolu"      long_name="snow volume"                                                                                                       unit="m"    />
+          <field id="icefrb"       long_name="Sea-ice freeboard"                                       standard_name="sea_ice_freeboard"                         unit="m"    />
           <field id="icealb"       long_name="Sea-ice or snow albedo"                                  standard_name="sea_ice_albedo"                            unit=""    detect_missing_value="true" />
-          <field id="tau_icebfr"   long_name="ice friction on ocean bottom for landfast ice"                                                                     unit="N/2" detect_missing_value="true" />
+          <field id="tau_icebfr"   long_name="ice friction on ocean bottom for landfast ice"                                                                     unit="N/2"  />
 	  
 	  <!-- melt ponds -->
@@ -61,14 +61,14 @@
 	  <!-- salt -->
           <field id="icesalt"      long_name="Sea ice salinity"                                                                                                  unit="g/kg"  detect_missing_value="true" />
-          <field id="icesalm"      long_name="Mass of salt in sea ice per area"                        standard_name="sea_ice_salt_mass"                         unit="kg/m2" detect_missing_value="true" />
+          <field id="icesalm"      long_name="Mass of salt in sea ice per area"                        standard_name="sea_ice_salt_mass"                         unit="kg/m2" />
 	  
 	  <!-- momentum (rheology) -->
           <field id="uice"         long_name="X-component of sea ice velocity"                         standard_name="sea_ice_x_velocity"                        unit="m/s"  /> 
           <field id="vice"         long_name="Y-component of sea ice velocity"                         standard_name="sea_ice_y_velocity"                        unit="m/s"  />      
-          <field id="icevel"       long_name="Sea-ice speed"                                           standard_name="sea_ice_speed"                             unit="m/s"  detect_missing_value="true" />
-          <field id="utau_ai"      long_name="X-component of atmospheric stress on sea ice"            standard_name="surface_downward_x_stress"                 unit="N/m2" detect_missing_value="true" />
-          <field id="vtau_ai"      long_name="Y-component of atmospheric stress on sea ice"            standard_name="surface_downward_y_stress"                 unit="N/m2" detect_missing_value="true" />
-          <field id="utau_oi"      long_name="X-component of ocean stress on sea ice"                  standard_name="sea_ice_base_upward_x_stress"              unit="N/m2" detect_missing_value="true" />
-          <field id="vtau_oi"      long_name="Y-component of ocean stress on sea ice"                  standard_name="sea_ice_base_upward_y_stress"              unit="N/m2" detect_missing_value="true" />
+          <field id="icevel"       long_name="Sea-ice speed"                                           standard_name="sea_ice_speed"                             unit="m/s"  />
+          <field id="utau_ai"      long_name="X-component of atmospheric stress on sea ice"            standard_name="surface_downward_x_stress"                 unit="N/m2" />
+          <field id="vtau_ai"      long_name="Y-component of atmospheric stress on sea ice"            standard_name="surface_downward_y_stress"                 unit="N/m2" />
+          <field id="utau_oi"      long_name="X-component of ocean stress on sea ice"                  standard_name="sea_ice_base_upward_x_stress"              unit="N/m2" />
+          <field id="vtau_oi"      long_name="Y-component of ocean stress on sea ice"                  standard_name="sea_ice_base_upward_y_stress"              unit="N/m2" />
           <field id="isig1"        long_name="1st principal stress component for EVP rhg"                                                                        unit=""     />
           <field id="isig2"        long_name="2nd principal stress component for EVP rhg"                                                                        unit=""     />
@@ -76,7 +76,7 @@
           <field id="normstr"      long_name="Average normal stress in sea ice"                        standard_name="average_normal_stress"                     unit="N/m"  />
           <field id="sheastr"      long_name="Maximum shear stress in sea ice"                         standard_name="maximum_shear_stress"                      unit="N/m"  />
-          <field id="icestr"       long_name="Compressive sea ice strength"                            standard_name="compressive_strength_of_sea_ice"           unit="N/m"  detect_missing_value="true" />
-          <field id="icediv"       long_name="Divergence of the sea-ice velocity field"                standard_name="divergence_of_sea_ice_velocity"            unit="s-1"  detect_missing_value="true" />
-          <field id="iceshe"       long_name="Maximum shear of sea-ice velocity field"                 standard_name="maximum_shear_of_sea_ice_velocity"         unit="s-1"  detect_missing_value="true" />
+          <field id="icestr"       long_name="Compressive sea ice strength"                            standard_name="compressive_strength_of_sea_ice"           unit="N/m"  />
+          <field id="icediv"       long_name="Divergence of the sea-ice velocity field"                standard_name="divergence_of_sea_ice_velocity"            unit="s-1"  />
+          <field id="iceshe"       long_name="Maximum shear of sea-ice velocity field"                 standard_name="maximum_shear_of_sea_ice_velocity"         unit="s-1"  />
 	  
 	  <!-- surface heat fluxes -->
@@ -92,7 +92,7 @@
 	  
 	  <!-- trends -->
-          <field id="afxthd"       long_name="sea-ice area fraction change from thermodynamics"  standard_name="tendency_of_sea_ice_area_fraction_due_to_dynamics" unit="s-1" />
-          <field id="afxdyn"       long_name="sea-ice area fraction change from dynamics"        standard_name="tendency_of_sea_ice_area_fraction_due_to_dynamics" unit="s-1" />
-	  <field id="afxtot"       long_name="area tendency (total)"                                                                                               unit="s-1" />
+          <field id="afxthd"       long_name="sea-ice area fraction change from thermodynamics" standard_name="tendency_of_sea_ice_area_fraction_due_to_dynamics" unit="s-1" />
+          <field id="afxdyn"       long_name="sea-ice area fraction change from dynamics"       standard_name="tendency_of_sea_ice_area_fraction_due_to_dynamics" unit="s-1" />
+	  <field id="afxtot"       long_name="area tendency (total)"                                                                                              unit="s-1" />
 	  
 	  <!-- momentum (advection) -->
@@ -248,10 +248,10 @@
           <field id="snwtemp_cat"  long_name="Snow temperature per category"                     unit="degC"    detect_missing_value="true" />
           <field id="icettop_cat"  long_name="Ice/snow surface temperature per category"         unit="degC"    detect_missing_value="true" />
-          <field id="iceapnd_cat"  long_name="Ice melt pond concentration per category"          unit="%"       /> 
+          <field id="iceapnd_cat"  long_name="Ice melt pond concentration per category"          unit=""        /> 
           <field id="icehpnd_cat"  long_name="Ice melt pond thickness per category"              unit="m"       detect_missing_value="true" /> 
-          <field id="iceafpnd_cat" long_name="Ice melt pond fraction per category"               unit="m"       /> 
+          <field id="iceafpnd_cat" long_name="Ice melt pond fraction per category"               unit=""        /> 
           <field id="icemask_cat"  long_name="Fraction of time step with sea ice (per category)" unit=""        />
           <field id="iceage_cat"   long_name="Ice age per category"                              unit="days"    detect_missing_value="true" />
-          <field id="icebrv_cat"   long_name="Brine volume per category"                         unit="%"       />
+          <field id="icebrv_cat"   long_name="Brine volume per category"                         unit="%"       detect_missing_value="true" />
 
           <!-- Add-ons for SIMIP -->
Index: NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/ice.F90
===================================================================
--- NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/ice.F90	(revision 11370)
+++ NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/ice.F90	(revision 11371)
@@ -253,6 +253,4 @@
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   wfx_err_sub !: mass flux error after sublimation                        [kg.m-2.s-1]
 
-   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   afx_tot     !: ice concentration tendency (total)        [s-1]
-
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   sfx_bog     !: salt flux due to ice bottom growth                   [pss.kg.m-2.s-1 => g.m-2.s-1]
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   sfx_bom     !: salt flux due to ice bottom melt                     [pss.kg.m-2.s-1 => g.m-2.s-1]
@@ -411,5 +409,5 @@
          &      wfx_bog    (jpi,jpj) , wfx_dyn   (jpi,jpj) , wfx_bom(jpi,jpj) , wfx_sum(jpi,jpj) ,           &
          &      wfx_res    (jpi,jpj) , wfx_sni   (jpi,jpj) , wfx_opw(jpi,jpj) , wfx_spr(jpi,jpj) ,           &
-         &      afx_tot    (jpi,jpj) , rn_amax_2d(jpi,jpj),                                                  &
+         &      rn_amax_2d (jpi,jpj) ,                                                                       &
          &      qsb_ice_bot(jpi,jpj) , qlead     (jpi,jpj) ,                                                 &
          &      sfx_res    (jpi,jpj) , sfx_bri   (jpi,jpj) , sfx_dyn(jpi,jpj) , sfx_sub(jpi,jpj) , sfx_lam(jpi,jpj) ,  &
Index: NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icecor.F90
===================================================================
--- NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icecor.F90	(revision 11370)
+++ NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icecor.F90	(revision 11371)
@@ -84,5 +84,5 @@
       !                             !-----------------------------------------------------
       IF ( nn_icesal == 2 ) THEN    !  salinity must stay in bounds [Simin,Simax]        !
-      !                             !-----------------------------------------------------
+         !                          !-----------------------------------------------------
          zzc = rhoi * r1_rdtice
          DO jl = 1, jpl
@@ -117,12 +117,7 @@
             END DO
          END DO
-         CALL lbc_lnk_multi( 'icecor', u_ice, 'U', -1., v_ice, 'V', -1. )            ! lateral boundary conditions
+         CALL lbc_lnk_multi( 'icecor', u_ice, 'U', -1., v_ice, 'V', -1. )
       ENDIF
 
-!!gm I guess the trends are only out on demand 
-!!   So please, only do this is it exite an iom_use of on a these variables
-!!   furthermore, only allocate the diag_ arrays in this case 
-!!   and do the iom_put here so that it is only a local allocation
-!!gm 
       !                             !-----------------------------------------------------
       SELECT CASE( kn )             !  Diagnostics                                       !
@@ -143,7 +138,8 @@
          END DO
          !                       ! concentration tendency (dynamics)
-         zafx   (:,:) = SUM( a_i(:,:,:) - a_i_b(:,:,:), dim=3 ) * r1_rdtice 
-         afx_tot(:,:) = zafx(:,:)
-         IF( iom_use('afxdyn') )   CALL iom_put( 'afxdyn' , zafx(:,:) )
+         IF( iom_use('afxdyn') .OR. iom_use('afxthd') .OR. iom_use('afxtot') ) THEN 
+            zafx(:,:) = SUM( a_i(:,:,:) - a_i_b(:,:,:), dim=3 ) * r1_rdtice 
+            CALL iom_put( 'afxdyn' , zafx )
+         ENDIF
          !
       CASE( 2 )                        !--- thermo trend diagnostics & ice aging
@@ -164,8 +160,9 @@
          END DO
          !                       ! concentration tendency (total + thermo)
-         zafx   (:,:) = SUM( a_i(:,:,:) - a_i_b(:,:,:), dim=3 ) * r1_rdtice
-         afx_tot(:,:) = afx_tot(:,:) + zafx(:,:)
-         IF( iom_use('afxthd') )   CALL iom_put( 'afxthd' , zafx(:,:) )
-         IF( iom_use('afxtot') )   CALL iom_put( 'afxtot' , afx_tot(:,:) )
+         IF( iom_use('afxdyn') .OR. iom_use('afxthd') .OR. iom_use('afxtot') ) THEN 
+            zafx(:,:) = zafx(:,:) + SUM( a_i(:,:,:) - a_i_b(:,:,:), dim=3 ) * r1_rdtice
+            CALL iom_put( 'afxthd' , SUM( a_i(:,:,:) - a_i_b(:,:,:), dim=3 ) * r1_rdtice )
+            CALL iom_put( 'afxtot' , zafx )
+         ENDIF
          !
       END SELECT
@@ -174,5 +171,6 @@
       IF( ln_icediachk   )   CALL ice_cons_hsm(1, 'icecor', rdiag_v, rdiag_s, rdiag_t, rdiag_fv, rdiag_fs, rdiag_ft) ! conservation
       IF( ln_ctl         )   CALL ice_prt3D   ('icecor')                                                             ! prints
-      IF( ln_icectl .AND. kn == 2 )   CALL ice_prt( kt, iiceprt, jiceprt, 2, ' - Final state - ' )                   ! prints
+      IF( ln_icectl .AND. kn == 2 ) &
+         &                   CALL ice_prt     ( kt, iiceprt, jiceprt, 2, ' - Final state - ' )                       ! prints
       IF( ln_timing      )   CALL timing_stop ('icecor')                                                             ! timing
       !
Index: NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icedia.F90
===================================================================
--- NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icedia.F90	(revision 11370)
+++ NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icedia.F90	(revision 11371)
@@ -148,7 +148,7 @@
          zdiff_tem = zdiff_tem - ( frc_tembot - frc_temtop )
          
-         CALL iom_put( 'ibgvolume' , zdiff_vol     )   ! ice/snow volume  drift            (km3 equivalent ocean water)         
-         CALL iom_put( 'ibgsaltco' , zdiff_sal     )   ! ice salt content drift            (psu*km3 equivalent ocean water)
-         CALL iom_put( 'ibgheatco' , zdiff_tem     )   ! ice/snow heat content drift       (1.e20 J)
+         CALL iom_put( 'ibgvolume' , zdiff_vol )   ! ice/snow volume  drift            (km3 equivalent ocean water)         
+         CALL iom_put( 'ibgsaltco' , zdiff_sal )   ! ice salt content drift            (psu*km3 equivalent ocean water)
+         CALL iom_put( 'ibgheatco' , zdiff_tem )   ! ice/snow heat content drift       (1.e20 J)
          !
       ENDIF
Index: NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icedyn.F90
===================================================================
--- NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icedyn.F90	(revision 11370)
+++ NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icedyn.F90	(revision 11371)
@@ -75,7 +75,6 @@
       !!
       INTEGER  ::   ji, jj        ! dummy loop indices
-      REAL(wp) ::   zmiss_val     ! missing value retrieved from xios 
       REAL(wp) ::   zcoefu, zcoefv
-      REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   zdivu_i, zmsk00
+      REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   zdivu_i
       !!--------------------------------------------------------------------
       !
@@ -89,7 +88,4 @@
       ENDIF
       !                      
-      ! get missing value from xml
-      CALL iom_miss_val( "icethic", zmiss_val )
-      !
       ! retrieve thickness from volume for landfast param. and UMx advection scheme
       WHERE( a_i(:,:,:) >= epsi20 )
@@ -159,11 +155,5 @@
          CASE ( np_dynADV1D , np_dynADV2D )
 
-            ALLOCATE( zdivu_i(jpi,jpj), zmsk00(jpi,jpj) )
-
-            DO jj = 1, jpj
-               DO ji = 1, jpi
-                  zmsk00(ji,jj) = MAX( 0._wp , SIGN( 1._wp , at_i(ji,jj) - epsi06  ) ) ! 1 if ice    , 0 if no ice
-               END DO
-            END DO
+            ALLOCATE( zdivu_i(jpi,jpj) )
             DO jj = 2, jpjm1
                DO ji = 2, jpim1
@@ -174,7 +164,7 @@
             CALL lbc_lnk( 'icedyn', zdivu_i, 'T', 1. )
             ! output
-            CALL iom_put( "icediv" , zdivu_i * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )
-
-            DEALLOCATE( zdivu_i, zmsk00 )
+            CALL iom_put( 'icediv' , zdivu_i )
+
+            DEALLOCATE( zdivu_i )
 
          END SELECT
Index: NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icedyn_adv.F90
===================================================================
--- NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icedyn_adv.F90	(revision 11370)
+++ NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icedyn_adv.F90	(revision 11371)
@@ -100,9 +100,9 @@
       diag_trp_vi(:,:) = SUM(     v_i (:,:,:)          - v_i_b (:,:,:)                  , dim=3 ) * r1_rdtice
       diag_trp_vs(:,:) = SUM(     v_s (:,:,:)          - v_s_b (:,:,:)                  , dim=3 ) * r1_rdtice
-      IF( iom_use('icemtrp') )   CALL iom_put( "icemtrp" , diag_trp_vi * rhoi          )   ! ice mass transport
-      IF( iom_use('snwmtrp') )   CALL iom_put( "snwmtrp" , diag_trp_vs * rhos          )   ! snw mass transport
-      IF( iom_use('salmtrp') )   CALL iom_put( "salmtrp" , diag_trp_sv * rhoi * 1.e-03 )   ! salt mass transport (kg/m2/s)
-      IF( iom_use('dihctrp') )   CALL iom_put( "dihctrp" , -diag_trp_ei                )   ! advected ice heat content (W/m2)
-      IF( iom_use('dshctrp') )   CALL iom_put( "dshctrp" , -diag_trp_es                )   ! advected snw heat content (W/m2)
+      IF( iom_use('icemtrp') )   CALL iom_put( 'icemtrp' ,  diag_trp_vi * rhoi          )   ! ice mass transport
+      IF( iom_use('snwmtrp') )   CALL iom_put( 'snwmtrp' ,  diag_trp_vs * rhos          )   ! snw mass transport
+      IF( iom_use('salmtrp') )   CALL iom_put( 'salmtrp' ,  diag_trp_sv * rhoi * 1.e-03 )   ! salt mass transport (kg/m2/s)
+      IF( iom_use('dihctrp') )   CALL iom_put( 'dihctrp' , -diag_trp_ei                 )   ! advected ice heat content (W/m2)
+      IF( iom_use('dshctrp') )   CALL iom_put( 'dshctrp' , -diag_trp_es                 )   ! advected snw heat content (W/m2)
 
       ! controls
Index: NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icedyn_rhg_evp.F90
===================================================================
--- NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icedyn_rhg_evp.F90	(revision 11370)
+++ NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icedyn_rhg_evp.F90	(revision 11371)
@@ -117,5 +117,4 @@
       INTEGER ::   jter         ! local integers
       !
-      REAL(wp) ::   zmiss_val                                           ! missing value retrieved from xios 
       REAL(wp) ::   zrhoco                                              ! rau0 * rn_cio
       REAL(wp) ::   zdtevp, z1_dtevp                                    ! time step for subcycling
@@ -167,12 +166,4 @@
       REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   zsig1, zsig2, zsig3
       !! --- SIMIP diags
-      REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   zdiag_dssh_dx   ! X-direction sea-surface tilt term (N/m2)
-      REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   zdiag_dssh_dy   ! X-direction sea-surface tilt term (N/m2)
-      REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   zdiag_corstrx   ! X-direction coriolis stress (N/m2)
-      REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   zdiag_corstry   ! Y-direction coriolis stress (N/m2)
-      REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   zdiag_intstrx   ! X-direction internal stress (N/m2)
-      REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   zdiag_intstry   ! Y-direction internal stress (N/m2)
-      REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   zdiag_utau_oi   ! X-direction ocean-ice stress
-      REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   zdiag_vtau_oi   ! Y-direction ocean-ice stress  
       REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   zdiag_xmtrp_ice ! X-component of ice mass transport (kg/s)
       REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   zdiag_ymtrp_ice ! Y-component of ice mass transport (kg/s)
@@ -760,7 +751,4 @@
       ! 5) diagnostics
       !------------------------------------------------------------------------------!
-      ! get missing value from xml
-      CALL iom_miss_val( "icethic", zmiss_val )
-
       DO jj = 1, jpj
          DO ji = 1, jpi
@@ -769,27 +757,21 @@
       END DO
 
-      ! --- ice-ocean stress --- !
-      IF( iom_use('utau_oi') .OR. iom_use('vtau_oi') ) THEN
-         ALLOCATE( zdiag_utau_oi(jpi,jpj) , zdiag_vtau_oi(jpi,jpj) )
-         zdiag_utau_oi(:,:) = zrhoco * SQRT( ( u_ice (:,:) - u_oce (:,:) ) * ( u_ice (:,:) - u_oce (:,:) )   &
-            &                              + ( v_iceU(:,:) - v_oceU(:,:) ) * ( v_iceU(:,:) - v_oceU(:,:) ) ) &
-            &                        * ( u_oce(:,:) - u_ice(:,:) )
-         !
-         zdiag_vtau_oi(:,:) = zrhoco * SQRT( ( v_ice (:,:) - v_oce (:,:) ) * ( v_ice (:,:) - v_oce (:,:) )   &
-            &                              + ( u_iceV(:,:) - u_oceV(:,:) ) * ( u_iceV(:,:) - u_oceV(:,:) ) ) &
-            &                        * ( v_oce(:,:) - v_ice(:,:) )
-         !
-         CALL iom_put( 'utau_oi' , zdiag_utau_oi * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )
-         CALL iom_put( 'vtau_oi' , zdiag_vtau_oi * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )
-         DEALLOCATE( zdiag_utau_oi , zdiag_vtau_oi )
+      ! --- ice-ocean and ice-atm. stresses --- !
+      IF( iom_use('utau_oi') .OR. iom_use('vtau_oi') .OR. iom_use('utau_ai') .OR. iom_use('vtau_ai') ) THEN
+         CALL lbc_lnk_multi( 'icedyn_rhg_evp', ztaux_oi, 'U', -1., ztauy_oi, 'V', -1., &
+            &                                  zTauU_ia, 'U', -1., zTauV_ia, 'V', -1. )
+         CALL iom_put( 'utau_oi' , ztaux_oi * zmsk00 )
+         CALL iom_put( 'vtau_oi' , ztauy_oi * zmsk00 )
+         CALL iom_put( 'utau_ai' , zTauU_ia * zmsk00 )
+         CALL iom_put( 'vtau_ai' , zTauV_ia * zmsk00 )
       ENDIF
  
-      ! --- ice-bathy stress in case of landfast --- !
-      IF( iom_use('tau_icebfr') )   CALL iom_put( 'tau_icebfr', tau_icebfr * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )
+      ! --- ice-oceanbottom stress in case of landfast --- !
+      IF( iom_use('tau_icebfr') )   CALL iom_put( 'tau_icebfr', tau_icebfr * zmsk00 )
       
       ! --- divergence, shear and strength --- !
-      IF( iom_use('icediv') )   CALL iom_put( "icediv" , pdivu_i  * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )   ! divergence
-      IF( iom_use('iceshe') )   CALL iom_put( "iceshe" , pshear_i * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )   ! shear
-      IF( iom_use('icestr') )   CALL iom_put( "icestr" , strength * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )   ! Ice strength
+      IF( iom_use('icediv') )   CALL iom_put( 'icediv' , pdivu_i  * zmsk00 )   ! divergence
+      IF( iom_use('iceshe') )   CALL iom_put( 'iceshe' , pshear_i * zmsk00 )   ! shear
+      IF( iom_use('icestr') )   CALL iom_put( 'icestr' , strength * zmsk00 )   ! strength
 
       ! --- stress tensor --- !
@@ -819,7 +801,7 @@
          CALL lbc_lnk_multi( 'icedyn_rhg_evp', zsig1, 'T', 1., zsig2, 'T', 1., zsig3, 'T', 1. )
          !
-         CALL iom_put( "isig1" , zsig1 )
-         CALL iom_put( "isig2" , zsig2 )
-         CALL iom_put( "isig3" , zsig3 )
+         CALL iom_put( 'isig1' , zsig1 )
+         CALL iom_put( 'isig2' , zsig2 )
+         CALL iom_put( 'isig3' , zsig3 )
          !
          ! Stress tensor invariants (normal and shear stress N/m)
@@ -831,66 +813,55 @@
       
       ! --- SIMIP --- !
-      IF ( iom_use( 'dssh_dx'  ) .OR. iom_use( 'dssh_dy'  ) .OR. iom_use( 'corstrx'  ) .OR. iom_use( 'corstry'  ) .OR. &
-         & iom_use( 'intstrx'  ) .OR. iom_use( 'intstry'  ) .OR. iom_use( 'xmtrpice' ) .OR. iom_use( 'ymtrpice' ) .OR. &
-         & iom_use( 'xmtrpsnw' ) .OR. iom_use( 'ymtrpsnw' ) .OR. iom_use( 'xatrp'    ) .OR. iom_use( 'yatrp'    ) ) THEN
-         !
-         ALLOCATE( zdiag_dssh_dx  (jpi,jpj) , zdiag_dssh_dy  (jpi,jpj) , zdiag_corstrx  (jpi,jpj) , zdiag_corstry  (jpi,jpj) , &
-            &      zdiag_intstrx  (jpi,jpj) , zdiag_intstry  (jpi,jpj) , zdiag_xmtrp_ice(jpi,jpj) , zdiag_ymtrp_ice(jpi,jpj) , &
-            &      zdiag_xmtrp_snw(jpi,jpj) , zdiag_ymtrp_snw(jpi,jpj) , zdiag_xatrp    (jpi,jpj) , zdiag_yatrp    (jpi,jpj) )
+      IF(  iom_use('dssh_dx') .OR. iom_use('dssh_dy') .OR. &
+         & iom_use('corstrx') .OR. iom_use('corstry') .OR. iom_use('intstrx') .OR. iom_use('intstry') ) THEN
+         !
+         CALL lbc_lnk_multi( 'icedyn_rhg_evp', zspgU, 'U', -1., zspgV, 'V', -1., &
+            &                                  zCorx, 'U', -1., zCory, 'V', -1., zfU, 'U', -1., zfV, 'V', -1. )
+
+         CALL iom_put( 'dssh_dx' , zspgU * zmsk00 )   ! Sea-surface tilt term in force balance (x)
+         CALL iom_put( 'dssh_dy' , zspgV * zmsk00 )   ! Sea-surface tilt term in force balance (y)
+         CALL iom_put( 'corstrx' , zCorx * zmsk00 )   ! Coriolis force term in force balance (x)
+         CALL iom_put( 'corstry' , zCory * zmsk00 )   ! Coriolis force term in force balance (y)
+         CALL iom_put( 'intstrx' , zfU   * zmsk00 )   ! Internal force term in force balance (x)
+         CALL iom_put( 'intstry' , zfV   * zmsk00 )   ! Internal force term in force balance (y)
+      ENDIF
+
+      IF(  iom_use('xmtrpice') .OR. iom_use('ymtrpice') .OR. &
+         & iom_use('xmtrpsnw') .OR. iom_use('ymtrpsnw') .OR. iom_use('xatrp') .OR. iom_use('yatrp') ) THEN
+         !
+         ALLOCATE( zdiag_xmtrp_ice(jpi,jpj) , zdiag_ymtrp_ice(jpi,jpj) , &
+            &      zdiag_xmtrp_snw(jpi,jpj) , zdiag_ymtrp_snw(jpi,jpj) , zdiag_xatrp(jpi,jpj) , zdiag_yatrp(jpi,jpj) )
          !
          DO jj = 2, jpjm1
             DO ji = 2, jpim1
-               ! Stress terms of the momentum equation (N/m2)
-               zdiag_dssh_dx(ji,jj) = zspgU(ji,jj) * zmsk00(ji,jj)     ! sea surface slope stress term
-               zdiag_dssh_dy(ji,jj) = zspgV(ji,jj) * zmsk00(ji,jj)
-
-               zdiag_corstrx(ji,jj) = zCorx(ji,jj) * zmsk00(ji,jj)     ! Coriolis stress term
-               zdiag_corstry(ji,jj) = zCory(ji,jj) * zmsk00(ji,jj)
-
-               zdiag_intstrx(ji,jj) = zfU(ji,jj)   * zmsk00(ji,jj)     ! internal stress term
-               zdiag_intstry(ji,jj) = zfV(ji,jj)   * zmsk00(ji,jj)
-
                ! 2D ice mass, snow mass, area transport arrays (X, Y)
                zfac_x = 0.5 * u_ice(ji,jj) * e2u(ji,jj) * zmsk00(ji,jj)
                zfac_y = 0.5 * v_ice(ji,jj) * e1v(ji,jj) * zmsk00(ji,jj)
-                              
+
                zdiag_xmtrp_ice(ji,jj) = rhoi * zfac_x * ( vt_i(ji+1,jj) + vt_i(ji,jj) ) ! ice mass transport, X-component
                zdiag_ymtrp_ice(ji,jj) = rhoi * zfac_y * ( vt_i(ji,jj+1) + vt_i(ji,jj) ) !        ''           Y-   ''
-               
+
                zdiag_xmtrp_snw(ji,jj) = rhos * zfac_x * ( vt_s(ji+1,jj) + vt_s(ji,jj) ) ! snow mass transport, X-component
                zdiag_ymtrp_snw(ji,jj) = rhos * zfac_y * ( vt_s(ji,jj+1) + vt_s(ji,jj) ) !          ''          Y-   ''
-               
+
                zdiag_xatrp(ji,jj)     = zfac_x * ( at_i(ji+1,jj) + at_i(ji,jj) )        ! area transport,      X-component
                zdiag_yatrp(ji,jj)     = zfac_y * ( at_i(ji,jj+1) + at_i(ji,jj) )        !        ''            Y-   ''
-               
-            END DO
-         END DO
-         
-         CALL lbc_lnk_multi( 'icedyn_rhg_evp', &
-            &                zdiag_dssh_dx, 'U', -1., zdiag_dssh_dy, 'V', -1.,   &
-            &                zdiag_corstrx, 'U', -1., zdiag_corstry, 'V', -1.,   & 
-            &                zdiag_intstrx, 'U', -1., zdiag_intstry, 'V', -1.    )
-                  
-         CALL lbc_lnk_multi( 'icedyn_rhg_evp', &
-            &                zdiag_xmtrp_ice, 'U', -1., zdiag_xmtrp_snw, 'U', -1.,   &
-            &                zdiag_xatrp    , 'U', -1., zdiag_ymtrp_ice, 'V', -1.,   &
-            &                zdiag_ymtrp_snw, 'V', -1., zdiag_yatrp    , 'V', -1.    )
-         
-         CALL iom_put( 'dssh_dx'  ,  zdiag_dssh_dx   )   ! Sea-surface tilt term in force balance (x)
-         CALL iom_put( 'dssh_dy'  ,  zdiag_dssh_dy   )   ! Sea-surface tilt term in force balance (y)
-         CALL iom_put( 'corstrx'  ,  zdiag_corstrx   )   ! Coriolis force term in force balance (x)
-         CALL iom_put( 'corstry'  ,  zdiag_corstry   )   ! Coriolis force term in force balance (y)
-         CALL iom_put( 'intstrx'  ,  zdiag_intstrx   )   ! Internal force term in force balance (x)
-         CALL iom_put( 'intstry'  ,  zdiag_intstry   )   ! Internal force term in force balance (y)
-         CALL iom_put( 'xmtrpice' ,  zdiag_xmtrp_ice )   ! X-component of sea-ice mass transport (kg/s)
-         CALL iom_put( 'ymtrpice' ,  zdiag_ymtrp_ice )   ! Y-component of sea-ice mass transport 
-         CALL iom_put( 'xmtrpsnw' ,  zdiag_xmtrp_snw )   ! X-component of snow mass transport (kg/s)
-         CALL iom_put( 'ymtrpsnw' ,  zdiag_ymtrp_snw )   ! Y-component of snow mass transport
-         CALL iom_put( 'xatrp'    ,  zdiag_xatrp     )   ! X-component of ice area transport
-         CALL iom_put( 'yatrp'    ,  zdiag_yatrp     )   ! Y-component of ice area transport
-
-         DEALLOCATE( zdiag_dssh_dx   , zdiag_dssh_dy   , zdiag_corstrx   , zdiag_corstry   , &
-            &        zdiag_intstrx   , zdiag_intstry   , zdiag_xmtrp_ice , zdiag_ymtrp_ice , &
-            &        zdiag_xmtrp_snw , zdiag_ymtrp_snw , zdiag_xatrp     , zdiag_yatrp     )
+
+            END DO
+         END DO
+
+         CALL lbc_lnk_multi( 'icedyn_rhg_evp', zdiag_xmtrp_ice, 'U', -1., zdiag_ymtrp_ice, 'V', -1., &
+            &                                  zdiag_xmtrp_snw, 'U', -1., zdiag_ymtrp_snw, 'V', -1., &
+            &                                  zdiag_xatrp    , 'U', -1., zdiag_yatrp    , 'V', -1. )
+
+         CALL iom_put( 'xmtrpice' , zdiag_xmtrp_ice )   ! X-component of sea-ice mass transport (kg/s)
+         CALL iom_put( 'ymtrpice' , zdiag_ymtrp_ice )   ! Y-component of sea-ice mass transport 
+         CALL iom_put( 'xmtrpsnw' , zdiag_xmtrp_snw )   ! X-component of snow mass transport (kg/s)
+         CALL iom_put( 'ymtrpsnw' , zdiag_ymtrp_snw )   ! Y-component of snow mass transport
+         CALL iom_put( 'xatrp'    , zdiag_xatrp     )   ! X-component of ice area transport
+         CALL iom_put( 'yatrp'    , zdiag_yatrp     )   ! Y-component of ice area transport
+
+         DEALLOCATE( zdiag_xmtrp_ice , zdiag_ymtrp_ice , &
+            &        zdiag_xmtrp_snw , zdiag_ymtrp_snw , zdiag_xatrp , zdiag_yatrp )
 
       ENDIF
Index: NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icesbc.F90
===================================================================
--- NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icesbc.F90	(revision 11370)
+++ NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icesbc.F90	(revision 11371)
@@ -129,5 +129,5 @@
 
       ! get missing value from xml
-      CALL iom_miss_val( "icethic", zmiss_val )
+      CALL iom_miss_val( "icetemp", zmiss_val )
 
       ! --- cloud-sky and overcast-sky ice albedos --- !
@@ -168,8 +168,8 @@
          END WHERE
          ! ice albedo
-         CALL iom_put( "icealb" , zalb * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )
+         CALL iom_put( 'icealb' , zalb * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )
          ! ice+ocean albedo
          zalb(:,:) = SUM( alb_ice * a_i_b, dim=3 ) + rn_alb_oce * ( 1._wp - at_i_b )
-         CALL iom_put( "albedo" , zalb )
+         CALL iom_put( 'albedo' , zalb )
 
          DEALLOCATE( zalb, zmsk00 )
Index: NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icestp.F90
===================================================================
--- NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icestp.F90	(revision 11370)
+++ NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icestp.F90	(revision 11371)
@@ -425,6 +425,4 @@
       wfx_err_sub(:,:) = 0._wp
       !
-      afx_tot(:,:) = 0._wp   ;
-      !
       diag_heat(:,:) = 0._wp ;   diag_sice(:,:) = 0._wp
       diag_vice(:,:) = 0._wp ;   diag_vsnw(:,:) = 0._wp
Index: NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/iceupdate.F90
===================================================================
--- NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/iceupdate.F90	(revision 11370)
+++ NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/iceupdate.F90	(revision 11371)
@@ -198,88 +198,88 @@
       ! --- salt fluxes [kg/m2/s] --- !
       !                           ! sfxice =  sfxbog + sfxbom + sfxsum + sfxsni + sfxopw + sfxres + sfxdyn + sfxbri + sfxsub + sfxlam
-      IF( iom_use('sfxice'  ) )   CALL iom_put( "sfxice", sfx     * 1.e-03 )   ! salt flux from total ice growth/melt
-      IF( iom_use('sfxbog'  ) )   CALL iom_put( "sfxbog", sfx_bog * 1.e-03 )   ! salt flux from bottom growth
-      IF( iom_use('sfxbom'  ) )   CALL iom_put( "sfxbom", sfx_bom * 1.e-03 )   ! salt flux from bottom melting
-      IF( iom_use('sfxsum'  ) )   CALL iom_put( "sfxsum", sfx_sum * 1.e-03 )   ! salt flux from surface melting
-      IF( iom_use('sfxlam'  ) )   CALL iom_put( "sfxlam", sfx_lam * 1.e-03 )   ! salt flux from lateral melting
-      IF( iom_use('sfxsni'  ) )   CALL iom_put( "sfxsni", sfx_sni * 1.e-03 )   ! salt flux from snow ice formation
-      IF( iom_use('sfxopw'  ) )   CALL iom_put( "sfxopw", sfx_opw * 1.e-03 )   ! salt flux from open water formation
-      IF( iom_use('sfxdyn'  ) )   CALL iom_put( "sfxdyn", sfx_dyn * 1.e-03 )   ! salt flux from ridging rafting
-      IF( iom_use('sfxbri'  ) )   CALL iom_put( "sfxbri", sfx_bri * 1.e-03 )   ! salt flux from brines
-      IF( iom_use('sfxres'  ) )   CALL iom_put( "sfxres", sfx_res * 1.e-03 )   ! salt flux from undiagnosed processes
-      IF( iom_use('sfxsub'  ) )   CALL iom_put( "sfxsub", sfx_sub * 1.e-03 )   ! salt flux from sublimation
+      IF( iom_use('sfxice'  ) )   CALL iom_put( 'sfxice', sfx     * 1.e-03 )   ! salt flux from total ice growth/melt
+      IF( iom_use('sfxbog'  ) )   CALL iom_put( 'sfxbog', sfx_bog * 1.e-03 )   ! salt flux from bottom growth
+      IF( iom_use('sfxbom'  ) )   CALL iom_put( 'sfxbom', sfx_bom * 1.e-03 )   ! salt flux from bottom melting
+      IF( iom_use('sfxsum'  ) )   CALL iom_put( 'sfxsum', sfx_sum * 1.e-03 )   ! salt flux from surface melting
+      IF( iom_use('sfxlam'  ) )   CALL iom_put( 'sfxlam', sfx_lam * 1.e-03 )   ! salt flux from lateral melting
+      IF( iom_use('sfxsni'  ) )   CALL iom_put( 'sfxsni', sfx_sni * 1.e-03 )   ! salt flux from snow ice formation
+      IF( iom_use('sfxopw'  ) )   CALL iom_put( 'sfxopw', sfx_opw * 1.e-03 )   ! salt flux from open water formation
+      IF( iom_use('sfxdyn'  ) )   CALL iom_put( 'sfxdyn', sfx_dyn * 1.e-03 )   ! salt flux from ridging rafting
+      IF( iom_use('sfxbri'  ) )   CALL iom_put( 'sfxbri', sfx_bri * 1.e-03 )   ! salt flux from brines
+      IF( iom_use('sfxres'  ) )   CALL iom_put( 'sfxres', sfx_res * 1.e-03 )   ! salt flux from undiagnosed processes
+      IF( iom_use('sfxsub'  ) )   CALL iom_put( 'sfxsub', sfx_sub * 1.e-03 )   ! salt flux from sublimation
 
       ! --- mass fluxes [kg/m2/s] --- !
-      IF( iom_use('emp_oce' ) )   CALL iom_put( "emp_oce", emp_oce )   ! emp over ocean (taking into account the snow blown away from the ice)
-      IF( iom_use('emp_ice' ) )   CALL iom_put( "emp_ice", emp_ice )   ! emp over ice   (taking into account the snow blown away from the ice)
+      CALL iom_put( 'emp_oce', emp_oce )   ! emp over ocean (taking into account the snow blown away from the ice)
+      CALL iom_put( 'emp_ice', emp_ice )   ! emp over ice   (taking into account the snow blown away from the ice)
 
       !                           ! vfxice = vfxbog + vfxbom + vfxsum + vfxsni + vfxopw + vfxdyn + vfxres + vfxlam + vfxpnd
-      IF( iom_use('vfxice'  ) )   CALL iom_put( "vfxice" , wfx_ice )   ! mass flux from total ice growth/melt
-      IF( iom_use('vfxbog'  ) )   CALL iom_put( "vfxbog" , wfx_bog )   ! mass flux from bottom growth
-      IF( iom_use('vfxbom'  ) )   CALL iom_put( "vfxbom" , wfx_bom )   ! mass flux from bottom melt 
-      IF( iom_use('vfxsum'  ) )   CALL iom_put( "vfxsum" , wfx_sum )   ! mass flux from surface melt 
-      IF( iom_use('vfxlam'  ) )   CALL iom_put( "vfxlam" , wfx_lam )   ! mass flux from lateral melt 
-      IF( iom_use('vfxsni'  ) )   CALL iom_put( "vfxsni" , wfx_sni )   ! mass flux from snow-ice formation
-      IF( iom_use('vfxopw'  ) )   CALL iom_put( "vfxopw" , wfx_opw )   ! mass flux from growth in open water
-      IF( iom_use('vfxdyn'  ) )   CALL iom_put( "vfxdyn" , wfx_dyn )   ! mass flux from dynamics (ridging)
-      IF( iom_use('vfxres'  ) )   CALL iom_put( "vfxres" , wfx_res )   ! mass flux from undiagnosed processes 
-      IF( iom_use('vfxpnd'  ) )   CALL iom_put( "vfxpnd" , wfx_pnd )   ! mass flux from melt ponds
-      IF( iom_use('vfxsub'  ) )   CALL iom_put( "vfxsub" , wfx_ice_sub )   ! mass flux from ice sublimation (ice-atm.)
-      IF( iom_use('vfxsub_err') ) CALL iom_put( "vfxsub_err", wfx_err_sub )   ! "excess" of sublimation sent to ocean      
-
-      IF ( iom_use( "vfxthin" ) ) THEN   ! mass flux from ice growth in open water + thin ice (<20cm) => comparable to observations  
+      CALL iom_put( 'vfxice'    , wfx_ice     )   ! mass flux from total ice growth/melt
+      CALL iom_put( 'vfxbog'    , wfx_bog     )   ! mass flux from bottom growth
+      CALL iom_put( 'vfxbom'    , wfx_bom     )   ! mass flux from bottom melt 
+      CALL iom_put( 'vfxsum'    , wfx_sum     )   ! mass flux from surface melt 
+      CALL iom_put( 'vfxlam'    , wfx_lam     )   ! mass flux from lateral melt 
+      CALL iom_put( 'vfxsni'    , wfx_sni     )   ! mass flux from snow-ice formation
+      CALL iom_put( 'vfxopw'    , wfx_opw     )   ! mass flux from growth in open water
+      CALL iom_put( 'vfxdyn'    , wfx_dyn     )   ! mass flux from dynamics (ridging)
+      CALL iom_put( 'vfxres'    , wfx_res     )   ! mass flux from undiagnosed processes 
+      CALL iom_put( 'vfxpnd'    , wfx_pnd     )   ! mass flux from melt ponds
+      CALL iom_put( 'vfxsub'    , wfx_ice_sub )   ! mass flux from ice sublimation (ice-atm.)
+      CALL iom_put( 'vfxsub_err', wfx_err_sub )   ! "excess" of sublimation sent to ocean      
+
+      IF ( iom_use( 'vfxthin' ) ) THEN   ! mass flux from ice growth in open water + thin ice (<20cm) => comparable to observations  
          WHERE( hm_i(:,:) < 0.2 .AND. hm_i(:,:) > 0. ) ; z2d = wfx_bog
          ELSEWHERE                                     ; z2d = 0._wp
          END WHERE
-         CALL iom_put( "vfxthin", wfx_opw + z2d )
-      ENDIF
-
-      !                              ! vfxsnw = vfxsnw_sni + vfxsnw_dyn + vfxsnw_sum
-      IF( iom_use('vfxsnw'     ) )   CALL iom_put( "vfxsnw"     , wfx_snw     )   ! mass flux from total snow growth/melt
-      IF( iom_use('vfxsnw_sum' ) )   CALL iom_put( "vfxsnw_sum" , wfx_snw_sum )   ! mass flux from snow melt at the surface
-      IF( iom_use('vfxsnw_sni' ) )   CALL iom_put( "vfxsnw_sni" , wfx_snw_sni )   ! mass flux from snow melt during snow-ice formation 
-      IF( iom_use('vfxsnw_dyn' ) )   CALL iom_put( "vfxsnw_dyn" , wfx_snw_dyn )   ! mass flux from dynamics (ridging) 
-      IF( iom_use('vfxsnw_sub' ) )   CALL iom_put( "vfxsnw_sub" , wfx_snw_sub )   ! mass flux from snow sublimation (ice-atm.) 
-      IF( iom_use('vfxsnw_pre' ) )   CALL iom_put( "vfxsnw_pre" , wfx_spr     )   ! snow precip
+         CALL iom_put( 'vfxthin', wfx_opw + z2d )
+      ENDIF
+
+      !                            ! vfxsnw = vfxsnw_sni + vfxsnw_dyn + vfxsnw_sum
+      CALL iom_put( 'vfxsnw'     , wfx_snw     )   ! mass flux from total snow growth/melt
+      CALL iom_put( 'vfxsnw_sum' , wfx_snw_sum )   ! mass flux from snow melt at the surface
+      CALL iom_put( 'vfxsnw_sni' , wfx_snw_sni )   ! mass flux from snow melt during snow-ice formation 
+      CALL iom_put( 'vfxsnw_dyn' , wfx_snw_dyn )   ! mass flux from dynamics (ridging) 
+      CALL iom_put( 'vfxsnw_sub' , wfx_snw_sub )   ! mass flux from snow sublimation (ice-atm.) 
+      CALL iom_put( 'vfxsnw_pre' , wfx_spr     )   ! snow precip
 
       ! --- heat fluxes [W/m2] --- !
       !                              ! qt_atm_oi - qt_oce_ai = hfxdhc - ( dihctrp + dshctrp )
-      IF( iom_use('qsr_oce'    ) )   CALL iom_put( "qsr_oce"    , qsr_oce * ( 1._wp - at_i_b )                               )   !     solar flux at ocean surface
-      IF( iom_use('qns_oce'    ) )   CALL iom_put( "qns_oce"    , qns_oce * ( 1._wp - at_i_b ) + qemp_oce                    )   ! non-solar flux at ocean surface
-      IF( iom_use('qsr_ice'    ) )   CALL iom_put( "qsr_ice"    , SUM( qsr_ice * a_i_b, dim=3 )                              )   !     solar flux at ice surface
-      IF( iom_use('qns_ice'    ) )   CALL iom_put( "qns_ice"    , SUM( qns_ice * a_i_b, dim=3 ) + qemp_ice                   )   ! non-solar flux at ice surface
-      IF( iom_use('qtr_ice_bot') )   CALL iom_put( "qtr_ice_bot", SUM( qtr_ice_bot * a_i_b, dim=3 )                          )   !     solar flux transmitted thru ice
-      IF( iom_use('qtr_ice_top') )   CALL iom_put( "qtr_ice_top", SUM( qtr_ice_top * a_i_b, dim=3 )                          )   !     solar flux transmitted thru ice surface
-      IF( iom_use('qt_oce'     ) )   CALL iom_put( "qt_oce"     ,      ( qsr_oce + qns_oce ) * ( 1._wp - at_i_b ) + qemp_oce )
-      IF( iom_use('qt_ice'     ) )   CALL iom_put( "qt_ice"     , SUM( ( qns_ice + qsr_ice ) * a_i_b, dim=3 )     + qemp_ice )
-      IF( iom_use('qt_oce_ai'  ) )   CALL iom_put( "qt_oce_ai"  , qt_oce_ai * tmask(:,:,1)                                   )   ! total heat flux at the ocean   surface: interface oce-(ice+atm) 
-      IF( iom_use('qt_atm_oi'  ) )   CALL iom_put( "qt_atm_oi"  , qt_atm_oi * tmask(:,:,1)                                   )   ! total heat flux at the oce-ice surface: interface atm-(ice+oce) 
-      IF( iom_use('qemp_oce'   ) )   CALL iom_put( "qemp_oce"   , qemp_oce                                                   )   ! Downward Heat Flux from E-P over ocean
-      IF( iom_use('qemp_ice'   ) )   CALL iom_put( "qemp_ice"   , qemp_ice                                                   )   ! Downward Heat Flux from E-P over ice
+      IF( iom_use('qsr_oce'    ) )   CALL iom_put( 'qsr_oce'    , qsr_oce * ( 1._wp - at_i_b )                               )   !     solar flux at ocean surface
+      IF( iom_use('qns_oce'    ) )   CALL iom_put( 'qns_oce'    , qns_oce * ( 1._wp - at_i_b ) + qemp_oce                    )   ! non-solar flux at ocean surface
+      IF( iom_use('qsr_ice'    ) )   CALL iom_put( 'qsr_ice'    , SUM( qsr_ice * a_i_b, dim=3 )                              )   !     solar flux at ice surface
+      IF( iom_use('qns_ice'    ) )   CALL iom_put( 'qns_ice'    , SUM( qns_ice * a_i_b, dim=3 ) + qemp_ice                   )   ! non-solar flux at ice surface
+      IF( iom_use('qtr_ice_bot') )   CALL iom_put( 'qtr_ice_bot', SUM( qtr_ice_bot * a_i_b, dim=3 )                          )   !     solar flux transmitted thru ice
+      IF( iom_use('qtr_ice_top') )   CALL iom_put( 'qtr_ice_top', SUM( qtr_ice_top * a_i_b, dim=3 )                          )   !     solar flux transmitted thru ice surface
+      IF( iom_use('qt_oce'     ) )   CALL iom_put( 'qt_oce'     ,      ( qsr_oce + qns_oce ) * ( 1._wp - at_i_b ) + qemp_oce )
+      IF( iom_use('qt_ice'     ) )   CALL iom_put( 'qt_ice'     , SUM( ( qns_ice + qsr_ice ) * a_i_b, dim=3 )     + qemp_ice )
+      IF( iom_use('qt_oce_ai'  ) )   CALL iom_put( 'qt_oce_ai'  , qt_oce_ai * tmask(:,:,1)                                   )   ! total heat flux at the ocean   surface: interface oce-(ice+atm) 
+      IF( iom_use('qt_atm_oi'  ) )   CALL iom_put( 'qt_atm_oi'  , qt_atm_oi * tmask(:,:,1)                                   )   ! total heat flux at the oce-ice surface: interface atm-(ice+oce) 
+      IF( iom_use('qemp_oce'   ) )   CALL iom_put( 'qemp_oce'   , qemp_oce                                                   )   ! Downward Heat Flux from E-P over ocean
+      IF( iom_use('qemp_ice'   ) )   CALL iom_put( 'qemp_ice'   , qemp_ice                                                   )   ! Downward Heat Flux from E-P over ice
 
       ! heat fluxes from ice transformations
-      !                              ! hfxdhc = hfxbog + hfxbom + hfxsum + hfxopw + hfxdif + hfxsnw - ( hfxthd + hfxdyn + hfxres + hfxsub + hfxspr )
-      IF( iom_use('hfxbog'     ) )   CALL iom_put ("hfxbog"     , hfx_bog             )   ! heat flux used for ice bottom growth 
-      IF( iom_use('hfxbom'     ) )   CALL iom_put ("hfxbom"     , hfx_bom             )   ! heat flux used for ice bottom melt
-      IF( iom_use('hfxsum'     ) )   CALL iom_put ("hfxsum"     , hfx_sum             )   ! heat flux used for ice surface melt
-      IF( iom_use('hfxopw'     ) )   CALL iom_put ("hfxopw"     , hfx_opw             )   ! heat flux used for ice formation in open water
-      IF( iom_use('hfxdif'     ) )   CALL iom_put ("hfxdif"     , hfx_dif             )   ! heat flux used for ice temperature change
-      IF( iom_use('hfxsnw'     ) )   CALL iom_put ("hfxsnw"     , hfx_snw             )   ! heat flux used for snow melt 
-      IF( iom_use('hfxerr'     ) )   CALL iom_put ("hfxerr"     , hfx_err_dif         )   ! heat flux error after heat diffusion (included in qt_oce_ai)
+      !                            ! hfxdhc = hfxbog + hfxbom + hfxsum + hfxopw + hfxdif + hfxsnw - ( hfxthd + hfxdyn + hfxres + hfxsub + hfxspr )
+      CALL iom_put ('hfxbog'     , hfx_bog     )   ! heat flux used for ice bottom growth 
+      CALL iom_put ('hfxbom'     , hfx_bom     )   ! heat flux used for ice bottom melt
+      CALL iom_put ('hfxsum'     , hfx_sum     )   ! heat flux used for ice surface melt
+      CALL iom_put ('hfxopw'     , hfx_opw     )   ! heat flux used for ice formation in open water
+      CALL iom_put ('hfxdif'     , hfx_dif     )   ! heat flux used for ice temperature change
+      CALL iom_put ('hfxsnw'     , hfx_snw     )   ! heat flux used for snow melt 
+      CALL iom_put ('hfxerr'     , hfx_err_dif )   ! heat flux error after heat diffusion (included in qt_oce_ai)
 
       ! heat fluxes associated with mass exchange (freeze/melt/precip...)
-      IF( iom_use('hfxthd'     ) )   CALL iom_put ("hfxthd"     , hfx_thd             )   !  
-      IF( iom_use('hfxdyn'     ) )   CALL iom_put ("hfxdyn"     , hfx_dyn             )   !  
-      IF( iom_use('hfxres'     ) )   CALL iom_put ("hfxres"     , hfx_res             )   !  
-      IF( iom_use('hfxsub'     ) )   CALL iom_put ("hfxsub"     , hfx_sub             )   !  
-      IF( iom_use('hfxspr'     ) )   CALL iom_put ("hfxspr"     , hfx_spr             )   ! Heat flux from snow precip heat content 
+      CALL iom_put ('hfxthd'     , hfx_thd     )   !  
+      CALL iom_put ('hfxdyn'     , hfx_dyn     )   !  
+      CALL iom_put ('hfxres'     , hfx_res     )   !  
+      CALL iom_put ('hfxsub'     , hfx_sub     )   !  
+      CALL iom_put ('hfxspr'     , hfx_spr     )   ! Heat flux from snow precip heat content 
 
       ! other heat fluxes
-      IF( iom_use('hfxsensib'  ) )   CALL iom_put( "hfxsensib"  ,     -qsb_ice_bot * at_i_b         )   ! Sensible oceanic heat flux
-      IF( iom_use('hfxcndbot'  ) )   CALL iom_put( "hfxcndbot"  , SUM( qcn_ice_bot * a_i_b, dim=3 ) )   ! Bottom conduction flux
-      IF( iom_use('hfxcndtop'  ) )   CALL iom_put( "hfxcndtop"  , SUM( qcn_ice_top * a_i_b, dim=3 ) )   ! Surface conduction flux
+      IF( iom_use('hfxsensib'  ) )   CALL iom_put( 'hfxsensib'  ,     -qsb_ice_bot * at_i_b         )   ! Sensible oceanic heat flux
+      IF( iom_use('hfxcndbot'  ) )   CALL iom_put( 'hfxcndbot'  , SUM( qcn_ice_bot * a_i_b, dim=3 ) )   ! Bottom conduction flux
+      IF( iom_use('hfxcndtop'  ) )   CALL iom_put( 'hfxcndtop'  , SUM( qcn_ice_top * a_i_b, dim=3 ) )   ! Surface conduction flux
 
       ! diags
-      IF( iom_use('hfxdhc'     ) )   CALL iom_put ("hfxdhc"     , diag_heat           )   ! Heat content variation in snow and ice 
+      CALL iom_put ('hfxdhc'     , diag_heat   )   ! Heat content variation in snow and ice 
       !
       ! controls
@@ -413,5 +413,5 @@
       !! ** Method  :   use of IOM library
       !!----------------------------------------------------------------------
-      CHARACTER(len=*) , INTENT(in) ::   cdrw   ! "READ"/"WRITE" flag
+      CHARACTER(len=*) , INTENT(in) ::   cdrw   ! 'READ'/'WRITE' flag
       INTEGER, OPTIONAL, INTENT(in) ::   kt     ! ice time-step
       !
Index: NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icewri.F90
===================================================================
--- NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icewri.F90	(revision 11370)
+++ NEMO/branches/2019/dev_r10984_HPC-13_IRRMANN_BDY_optimization/src/ICE/icewri.F90	(revision 11371)
@@ -63,5 +63,5 @@
 
       ! get missing value from xml
-      CALL iom_miss_val( "icethic", zmiss_val )
+      CALL iom_miss_val( 'icetemp', zmiss_val )
 
       ! brine volume
@@ -91,51 +91,44 @@
       zrho1 = ( rau0 - rhoi ) * r1_rau0 ; zrho2 = rhos * r1_rau0
       ! masks
-      CALL iom_put( "icemask"  , zmsk00 )   ! ice mask 0%
-      CALL iom_put( "icemask05", zmsk05 )   ! ice mask 5%
-      CALL iom_put( "icemask15", zmsk15 )   ! ice mask 15%
-      CALL iom_put( "icepres"  , zmsk00 )   ! Ice presence (1 or 0) 
+      CALL iom_put( 'icemask'  , zmsk00 )   ! ice mask 0%
+      CALL iom_put( 'icemask05', zmsk05 )   ! ice mask 5%
+      CALL iom_put( 'icemask15', zmsk15 )   ! ice mask 15%
+      CALL iom_put( 'icepres'  , zmsk00 )   ! Ice presence (1 or 0) 
       !
       ! general fields
-      IF( iom_use('icemass' ) )   CALL iom_put( "icemass", rhoi * vt_i * zmsk00  )   ! Ice mass per cell area 
-      IF( iom_use('snwmass' ) )   CALL iom_put( "snwmass", rhos * vt_s * zmsksn  )   ! Snow mass per cell area
-      IF( iom_use('iceconc' ) )   CALL iom_put( "iceconc", at_i  * zmsk00        )   ! ice concentration
-      IF( iom_use('icevolu' ) )   CALL iom_put( "icevolu", vt_i  * zmsk00        )   ! ice volume = mean ice thickness over the cell
-      IF( iom_use('icethic' ) )   CALL iom_put( "icethic", hm_i  * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )   ! ice thickness
-      IF( iom_use('snwthic' ) )   CALL iom_put( "snwthic", hm_s  * zmsk00 + zmiss_val * ( 1._wp - zmsk00 )       )   ! snw thickness
-      IF( iom_use('icebrv'  ) )   CALL iom_put( "icebrv" , bvm_i * zmsk00 * 100. )   ! brine volume
-      IF( iom_use('iceage'  ) )   CALL iom_put( "iceage" , om_i * zmsk15 / rday + zmiss_val * ( 1._wp - zmsk15 ) )   ! ice age
-      IF( iom_use('icehnew' ) )   CALL iom_put( "icehnew", ht_i_new              )   ! new ice thickness formed in the leads
-      IF( iom_use('snwvolu' ) )   CALL iom_put( "snwvolu", vt_s  * zmsksn        )   ! snow volume
-      IF( iom_use('icefrb'  ) ) THEN
+      IF( iom_use('icemass' ) )   CALL iom_put( 'icemass', vt_i * rhoi * zmsk00 )                                           ! Ice mass per cell area 
+      IF( iom_use('snwmass' ) )   CALL iom_put( 'snwmass', vt_s * rhos * zmsksn )                                           ! Snow mass per cell area
+      IF( iom_use('iceconc' ) )   CALL iom_put( 'iceconc', at_i        * zmsk00 )                                           ! ice concentration
+      IF( iom_use('icevolu' ) )   CALL iom_put( 'icevolu', vt_i        * zmsk00 )                                           ! ice volume = mean ice thickness over the cell
+      IF( iom_use('icethic' ) )   CALL iom_put( 'icethic', hm_i        * zmsk00 )                                           ! ice thickness
+      IF( iom_use('snwthic' ) )   CALL iom_put( 'snwthic', hm_s        * zmsk00 )                                           ! snw thickness
+      IF( iom_use('icebrv'  ) )   CALL iom_put( 'icebrv' , bvm_i* 100. * zmsk00 )                                           ! brine volume
+      IF( iom_use('iceage'  ) )   CALL iom_put( 'iceage' , om_i / rday * zmsk15 + zmiss_val * ( 1._wp - zmsk15 ) )          ! ice age
+      IF( iom_use('icehnew' ) )   CALL iom_put( 'icehnew', ht_i_new             )                                           ! new ice thickness formed in the leads
+      IF( iom_use('snwvolu' ) )   CALL iom_put( 'snwvolu', vt_s        * zmsksn )                                           ! snow volume
+      IF( iom_use('icefrb'  ) ) THEN                                                                                        ! Ice freeboard
          z2d(:,:) = ( zrho1 * hm_i(:,:) - zrho2 * hm_s(:,:) )                                         
          WHERE( z2d < 0._wp )   z2d = 0._wp
-                                   CALL iom_put( "icefrb" , z2d * zmsk00 + zmiss_val * ( 1._wp - zmsk00 )          )   ! Ice freeboard
+                                  CALL iom_put( 'icefrb' , z2d * zmsk00         )
       ENDIF
-      !
       ! melt ponds
-      IF( iom_use('iceapnd' ) )   CALL iom_put( "iceapnd", at_ip  * zmsk00       )   ! melt pond total fraction
-      IF( iom_use('icevpnd' ) )   CALL iom_put( "icevpnd", vt_ip  * zmsk00       )   ! melt pond total volume per unit area
-      !
+      IF( iom_use('iceapnd' ) )   CALL iom_put( 'iceapnd', at_ip  * zmsk00      )                                           ! melt pond total fraction
+      IF( iom_use('icevpnd' ) )   CALL iom_put( 'icevpnd', vt_ip  * zmsk00      )                                           ! melt pond total volume per unit area
       ! salt
-      IF( iom_use('icesalt' ) )   CALL iom_put( "icesalt", sm_i                 * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )   ! mean ice salinity
-      IF( iom_use('icesalm' ) )   CALL iom_put( "icesalm", st_i * rhoi * 1.0e-3 * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )   ! Mass of salt in sea ice per cell area
-
+      IF( iom_use('icesalt' ) )   CALL iom_put( 'icesalt', sm_i                 * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) ) ! mean ice salinity
+      IF( iom_use('icesalm' ) )   CALL iom_put( 'icesalm', st_i * rhoi * 1.0e-3 * zmsk00 )                                  ! Mass of salt in sea ice per cell area
       ! heat
-      IF( iom_use('icetemp' ) )   CALL iom_put( "icetemp", ( tm_i  - rt0 ) * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )   ! ice mean temperature
-      IF( iom_use('snwtemp' ) )   CALL iom_put( "snwtemp", ( tm_s  - rt0 ) * zmsksn + zmiss_val * ( 1._wp - zmsksn ) )   ! snw mean temperature
-      IF( iom_use('icettop' ) )   CALL iom_put( "icettop", ( tm_su - rt0 ) * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )   ! temperature at the ice surface
-      IF( iom_use('icetbot' ) )   CALL iom_put( "icetbot", ( t_bo  - rt0 ) * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )   ! temperature at the ice bottom
-      IF( iom_use('icetsni' ) )   CALL iom_put( "icetsni", ( tm_si - rt0 ) * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )   ! temperature at the snow-ice interface
-      IF( iom_use('icehc'   ) )   CALL iom_put( "icehc"  ,  -et_i          * zmsk00 )   ! ice heat content
-      IF( iom_use('snwhc'   ) )   CALL iom_put( "snwhc"  ,  -et_s          * zmsksn )   ! snow heat content
-
+      IF( iom_use('icetemp' ) )   CALL iom_put( 'icetemp', ( tm_i  - rt0 ) * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )      ! ice mean temperature
+      IF( iom_use('snwtemp' ) )   CALL iom_put( 'snwtemp', ( tm_s  - rt0 ) * zmsksn + zmiss_val * ( 1._wp - zmsksn ) )      ! snw mean temperature
+      IF( iom_use('icettop' ) )   CALL iom_put( 'icettop', ( tm_su - rt0 ) * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )      ! temperature at the ice surface
+      IF( iom_use('icetbot' ) )   CALL iom_put( 'icetbot', ( t_bo  - rt0 ) * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )      ! temperature at the ice bottom
+      IF( iom_use('icetsni' ) )   CALL iom_put( 'icetsni', ( tm_si - rt0 ) * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )      ! temperature at the snow-ice interface
+      IF( iom_use('icehc'   ) )   CALL iom_put( 'icehc'  ,  -et_i          * zmsk00 )                                       ! ice heat content
+      IF( iom_use('snwhc'   ) )   CALL iom_put( 'snwhc'  ,  -et_s          * zmsksn )                                       ! snow heat content
       ! momentum
-      IF( iom_use('uice'    ) )   CALL iom_put( "uice"   , u_ice     )   ! ice velocity u
-      IF( iom_use('vice'    ) )   CALL iom_put( "vice"   , v_ice     )   ! ice velocity v
-      IF( iom_use('utau_ai' ) )   CALL iom_put( "utau_ai", utau_ice * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )   ! Wind stress u
-      IF( iom_use('vtau_ai' ) )   CALL iom_put( "vtau_ai", vtau_ice * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )   ! Wind stress v
-
-      IF( iom_use('icevel') .OR. iom_use('fasticepres') ) THEN 
-        ! module of ice velocity
+      IF( iom_use('uice'    ) )   CALL iom_put( 'uice'   , u_ice    )                                                       ! ice velocity u
+      IF( iom_use('vice'    ) )   CALL iom_put( 'vice'   , v_ice    )                                                       ! ice velocity v
+      !
+      IF( iom_use('icevel') .OR. iom_use('fasticepres') ) THEN                                                              ! module of ice velocity
          DO jj = 2 , jpjm1
             DO ji = 2 , jpim1
@@ -146,29 +139,28 @@
          END DO
          CALL lbc_lnk( 'icewri', z2d, 'T', 1. )
-         CALL iom_put( "icevel", z2d * zmsk00 + zmiss_val * ( 1._wp - zmsk00 ) )
-
-        ! record presence of fast ice
-         WHERE( z2d(:,:) < 5.e-04_wp .AND. zmsk15(:,:) == 1._wp ) ; zfast(:,:) = 1._wp
+         CALL iom_put( 'icevel', z2d )
+
+         WHERE( z2d(:,:) < 5.e-04_wp .AND. zmsk15(:,:) == 1._wp ) ; zfast(:,:) = 1._wp                                      ! record presence of fast ice
          ELSEWHERE                                                ; zfast(:,:) = 0._wp
          END WHERE
-         CALL iom_put( "fasticepres", zfast )
+         CALL iom_put( 'fasticepres', zfast )
       ENDIF
 
       ! --- category-dependent fields --- !
-      IF( iom_use('icemask_cat' ) )   CALL iom_put( "icemask_cat" , zmsk00l                                                    )   ! ice mask 0%
-      IF( iom_use('iceconc_cat' ) )   CALL iom_put( "iceconc_cat" , a_i * zmsk00l                                              )   ! area for categories
-      IF( iom_use('icethic_cat' ) )   CALL iom_put( "icethic_cat" , h_i * zmsk00l + zmiss_val * ( 1._wp - zmsk00l )            )   ! thickness for categories
-      IF( iom_use('snwthic_cat' ) )   CALL iom_put( "snwthic_cat" , h_s * zmsksnl + zmiss_val * ( 1._wp - zmsksnl )            )   ! snow depth for categories
-      IF( iom_use('icesalt_cat' ) )   CALL iom_put( "icesalt_cat" , s_i * zmsk00l + zmiss_val * ( 1._wp - zmsk00l )            )   ! salinity for categories
-      IF( iom_use('iceage_cat'  ) )   CALL iom_put( "iceage_cat"  , o_i * zmsk00l / rday + zmiss_val * ( 1._wp - zmsk00l )     )   ! ice age
-      IF( iom_use('icetemp_cat' ) )   CALL iom_put( "icetemp_cat" , ( SUM( t_i(:,:,:,:), dim=3 ) * r1_nlay_i - rt0 ) * zmsk00l &
-         &                                                          + zmiss_val * ( 1._wp - zmsk00l ) )                            ! ice temperature
-      IF( iom_use('snwtemp_cat' ) )   CALL iom_put( "snwtemp_cat" , ( SUM( t_s(:,:,:,:), dim=3 ) * r1_nlay_s - rt0 ) * zmsksnl &
-         &                                                          + zmiss_val * ( 1._wp - zmsksnl ) )                            ! snow temperature
-      IF( iom_use('icettop_cat' ) )   CALL iom_put( "icettop_cat" , ( t_su - rt0 ) * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) )   ! surface temperature
-      IF( iom_use('icebrv_cat'  ) )   CALL iom_put( "icebrv_cat"  ,   bv_i * 100.  * zmsk00l                                   )   ! brine volume
-      IF( iom_use('iceapnd_cat' ) )   CALL iom_put( "iceapnd_cat" ,   a_ip         * zmsk00l                                   )   ! melt pond frac for categories
-      IF( iom_use('icehpnd_cat' ) )   CALL iom_put( "icehpnd_cat" ,   h_ip         * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) )   ! melt pond frac for categories
-      IF( iom_use('iceafpnd_cat') )   CALL iom_put( "iceafpnd_cat",   a_ip_frac    * zmsk00l                                   )   ! melt pond frac for categories
+      IF( iom_use('icemask_cat' ) )   CALL iom_put( 'icemask_cat' ,                  zmsk00l                                   ) ! ice mask 0%
+      IF( iom_use('iceconc_cat' ) )   CALL iom_put( 'iceconc_cat' , a_i            * zmsk00l                                   ) ! area for categories
+      IF( iom_use('icethic_cat' ) )   CALL iom_put( 'icethic_cat' , h_i            * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) ) ! thickness for categories
+      IF( iom_use('snwthic_cat' ) )   CALL iom_put( 'snwthic_cat' , h_s            * zmsksnl + zmiss_val * ( 1._wp - zmsksnl ) ) ! snow depth for categories
+      IF( iom_use('icesalt_cat' ) )   CALL iom_put( 'icesalt_cat' , s_i            * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) ) ! salinity for categories
+      IF( iom_use('iceage_cat'  ) )   CALL iom_put( 'iceage_cat'  , o_i / rday     * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) ) ! ice age
+      IF( iom_use('icetemp_cat' ) )   CALL iom_put( 'icetemp_cat' , ( SUM( t_i, dim=3 ) * r1_nlay_i - rt0 ) &
+         &                                                                         * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) ) ! ice temperature
+      IF( iom_use('snwtemp_cat' ) )   CALL iom_put( 'snwtemp_cat' , ( SUM( t_s, dim=3 ) * r1_nlay_s - rt0 ) &
+         &                                                                         * zmsksnl + zmiss_val * ( 1._wp - zmsksnl ) ) ! snow temperature
+      IF( iom_use('icettop_cat' ) )   CALL iom_put( 'icettop_cat' , ( t_su - rt0 ) * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) ) ! surface temperature
+      IF( iom_use('icebrv_cat'  ) )   CALL iom_put( 'icebrv_cat'  ,   bv_i * 100.  * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) ) ! brine volume
+      IF( iom_use('iceapnd_cat' ) )   CALL iom_put( 'iceapnd_cat' ,   a_ip         * zmsk00l                                   ) ! melt pond frac for categories
+      IF( iom_use('icehpnd_cat' ) )   CALL iom_put( 'icehpnd_cat' ,   h_ip         * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) ) ! melt pond frac for categories
+      IF( iom_use('iceafpnd_cat') )   CALL iom_put( 'iceafpnd_cat',   a_ip_frac    * zmsk00l                                   ) ! melt pond frac for categories
 
       !------------------
@@ -176,18 +168,18 @@
       !------------------
       ! trends
-      IF( iom_use('dmithd') )   CALL iom_put( "dmithd", - wfx_bog - wfx_bom - wfx_sum - wfx_sni - wfx_opw - wfx_lam - wfx_res ) ! Sea-ice mass change from thermodynamics
-      IF( iom_use('dmidyn') )   CALL iom_put( "dmidyn", - wfx_dyn + rhoi * diag_trp_vi      )   ! Sea-ice mass change from dynamics(kg/m2/s)
-      IF( iom_use('dmiopw') )   CALL iom_put( "dmiopw", - wfx_opw                           )   ! Sea-ice mass change through growth in open water
-      IF( iom_use('dmibog') )   CALL iom_put( "dmibog", - wfx_bog                           )   ! Sea-ice mass change through basal growth
-      IF( iom_use('dmisni') )   CALL iom_put( "dmisni", - wfx_sni                           )   ! Sea-ice mass change through snow-to-ice conversion
-      IF( iom_use('dmisum') )   CALL iom_put( "dmisum", - wfx_sum                           )   ! Sea-ice mass change through surface melting
-      IF( iom_use('dmibom') )   CALL iom_put( "dmibom", - wfx_bom                           )   ! Sea-ice mass change through bottom melting
-      IF( iom_use('dmtsub') )   CALL iom_put( "dmtsub", - wfx_sub                           )   ! Sea-ice mass change through evaporation and sublimation
-      IF( iom_use('dmssub') )   CALL iom_put( "dmssub", - wfx_snw_sub                       )   ! Snow mass change through sublimation
-      IF( iom_use('dmisub') )   CALL iom_put( "dmisub", - wfx_ice_sub                       )   ! Sea-ice mass change through sublimation
-      IF( iom_use('dmsspr') )   CALL iom_put( "dmsspr", - wfx_spr                           )   ! Snow mass change through snow fall
-      IF( iom_use('dmsssi') )   CALL iom_put( "dmsssi",   wfx_sni*rhos*r1_rhoi              )   ! Snow mass change through snow-to-ice conversion
-      IF( iom_use('dmsmel') )   CALL iom_put( "dmsmel", - wfx_snw_sum                       )   ! Snow mass change through melt
-      IF( iom_use('dmsdyn') )   CALL iom_put( "dmsdyn", - wfx_snw_dyn + rhos * diag_trp_vs  )   ! Snow mass change through dynamics(kg/m2/s)
+      IF( iom_use('dmithd') )   CALL iom_put( 'dmithd', - wfx_bog - wfx_bom - wfx_sum - wfx_sni - wfx_opw - wfx_lam - wfx_res ) ! Sea-ice mass change from thermodynamics
+      IF( iom_use('dmidyn') )   CALL iom_put( 'dmidyn', - wfx_dyn + rhoi * diag_trp_vi                                        ) ! Sea-ice mass change from dynamics(kg/m2/s)
+      IF( iom_use('dmiopw') )   CALL iom_put( 'dmiopw', - wfx_opw                                                             ) ! Sea-ice mass change through growth in open water
+      IF( iom_use('dmibog') )   CALL iom_put( 'dmibog', - wfx_bog                                                             ) ! Sea-ice mass change through basal growth
+      IF( iom_use('dmisni') )   CALL iom_put( 'dmisni', - wfx_sni                                                             ) ! Sea-ice mass change through snow-to-ice conversion
+      IF( iom_use('dmisum') )   CALL iom_put( 'dmisum', - wfx_sum                                                             ) ! Sea-ice mass change through surface melting
+      IF( iom_use('dmibom') )   CALL iom_put( 'dmibom', - wfx_bom                                                             ) ! Sea-ice mass change through bottom melting
+      IF( iom_use('dmtsub') )   CALL iom_put( 'dmtsub', - wfx_sub                                                             ) ! Sea-ice mass change through evaporation and sublimation
+      IF( iom_use('dmssub') )   CALL iom_put( 'dmssub', - wfx_snw_sub                                                         ) ! Snow mass change through sublimation
+      IF( iom_use('dmisub') )   CALL iom_put( 'dmisub', - wfx_ice_sub                                                         ) ! Sea-ice mass change through sublimation
+      IF( iom_use('dmsspr') )   CALL iom_put( 'dmsspr', - wfx_spr                                                             ) ! Snow mass change through snow fall
+      IF( iom_use('dmsssi') )   CALL iom_put( 'dmsssi',   wfx_sni*rhos*r1_rhoi                                                ) ! Snow mass change through snow-to-ice conversion
+      IF( iom_use('dmsmel') )   CALL iom_put( 'dmsmel', - wfx_snw_sum                                                         ) ! Snow mass change through melt
+      IF( iom_use('dmsdyn') )   CALL iom_put( 'dmsdyn', - wfx_snw_dyn + rhos * diag_trp_vs                                    ) ! Snow mass change through dynamics(kg/m2/s)
       
       ! Global ice diagnostics
@@ -207,10 +199,10 @@
          IF( iom_use('SH_iceextt') )   zdiag_extt_sh = glob_sum( 'icewri',        ( 1._wp - z2d ) * e1e2t * 1.e-12 * zmsk15 )
          !
-         CALL iom_put( "NH_icearea" , zdiag_area_nh )
-         CALL iom_put( "NH_icevolu" , zdiag_volu_nh )
-         CALL iom_put( "NH_iceextt" , zdiag_extt_nh )
-         CALL iom_put( "SH_icearea" , zdiag_area_sh )
-         CALL iom_put( "SH_icevolu" , zdiag_volu_sh )
-         CALL iom_put( "SH_iceextt" , zdiag_extt_sh )
+         CALL iom_put( 'NH_icearea' , zdiag_area_nh )
+         CALL iom_put( 'NH_icevolu' , zdiag_volu_nh )
+         CALL iom_put( 'NH_iceextt' , zdiag_extt_nh )
+         CALL iom_put( 'SH_icearea' , zdiag_area_sh )
+         CALL iom_put( 'SH_icevolu' , zdiag_volu_sh )
+         CALL iom_put( 'SH_iceextt' , zdiag_extt_sh )
          !
       ENDIF
