Index: NEMO/branches/2020/dev_r12512_HPC-04_mcastril_Mixed_Precision_implementation/src/ICE/iceitd.F90
===================================================================
--- NEMO/branches/2020/dev_r12512_HPC-04_mcastril_Mixed_Precision_implementation/src/ICE/iceitd.F90	(revision 12632)
+++ NEMO/branches/2020/dev_r12512_HPC-04_mcastril_Mixed_Precision_implementation/src/ICE/iceitd.F90	(revision 12701)
@@ -148,6 +148,11 @@
                !    Note: hn(t+1) must not be too close to either HR or HL otherwise a division by nearly 0 is possible 
                !          in itd_glinear in the case (HR-HL) = 3(Hice - HL) or = 3(HR - Hice)
+# if defined key_single
+               IF( a_i_2d(ji,jl  ) > epsi10 .AND. h_i_2d(ji,jl  ) > ( zhbnew(ji,jl) - epsi06 ) )   nptidx(ji) = 0
+               IF( a_i_2d(ji,jl+1) > epsi10 .AND. h_i_2d(ji,jl+1) < ( zhbnew(ji,jl) + epsi06 ) )   nptidx(ji) = 0
+# else
                IF( a_i_2d(ji,jl  ) > epsi10 .AND. h_i_2d(ji,jl  ) > ( zhbnew(ji,jl) - epsi10 ) )   nptidx(ji) = 0
                IF( a_i_2d(ji,jl+1) > epsi10 .AND. h_i_2d(ji,jl+1) < ( zhbnew(ji,jl) + epsi10 ) )   nptidx(ji) = 0
+# endif
                !
                ! 2) Hn-1 < Hn* < Hn+1  
@@ -170,6 +175,11 @@
             !    h1(t) must not be too close to either HR or HL otherwise a division by nearly 0 is possible 
             !    in itd_glinear in the case (HR-HL) = 3(Hice - HL) or = 3(HR - Hice)
+# if defined key_single
+            IF( h_ib_2d(ji,1) < ( hi_max(0) + epsi06 ) )   nptidx(ji) = 0
+            IF( h_ib_2d(ji,1) > ( hi_max(1) - epsi06 ) )   nptidx(ji) = 0
+# else
             IF( h_ib_2d(ji,1) < ( hi_max(0) + epsi10 ) )   nptidx(ji) = 0
             IF( h_ib_2d(ji,1) > ( hi_max(1) - epsi10 ) )   nptidx(ji) = 0
+# endif
          END DO
          !
@@ -538,5 +548,9 @@
       ! 4) Update ice thickness and temperature
       !-------------------------------------------------------------------------------
+# if defined SINGLE_PRECISION
+      WHERE( a_i_2d(1:npti,:) >= epsi06 )
+# else
       WHERE( a_i_2d(1:npti,:) >= epsi20 )
+# endif
          h_i_2d (1:npti,:)  =  v_i_2d(1:npti,:) / a_i_2d(1:npti,:) 
          t_su_2d(1:npti,:)  =  zaTsfn(1:npti,:) / a_i_2d(1:npti,:) 
Index: NEMO/branches/2020/dev_r12512_HPC-04_mcastril_Mixed_Precision_implementation/src/ICE/icethd.F90
===================================================================
--- NEMO/branches/2020/dev_r12512_HPC-04_mcastril_Mixed_Precision_implementation/src/ICE/icethd.F90	(revision 12632)
+++ NEMO/branches/2020/dev_r12512_HPC-04_mcastril_Mixed_Precision_implementation/src/ICE/icethd.F90	(revision 12701)
@@ -218,5 +218,5 @@
                               CALL ice_thd_dh                           ! Ice-Snow thickness   
                               CALL ice_thd_pnd                          ! Melt ponds formation
-                              CALL ice_thd_ent( e_i_1d(1:npti,:) )      ! Ice enthalpy remapping
+                              CALL ice_thd_ent( e_i_1d(1:npti,:), .true. )      ! Ice enthalpy remapping
             ENDIF
                               CALL ice_thd_sal( ln_icedS )          ! --- Ice salinity --- !    
Index: NEMO/branches/2020/dev_r12512_HPC-04_mcastril_Mixed_Precision_implementation/src/ICE/icethd_do.F90
===================================================================
--- NEMO/branches/2020/dev_r12512_HPC-04_mcastril_Mixed_Precision_implementation/src/ICE/icethd_do.F90	(revision 12632)
+++ NEMO/branches/2020/dev_r12512_HPC-04_mcastril_Mixed_Precision_implementation/src/ICE/icethd_do.F90	(revision 12701)
@@ -385,5 +385,5 @@
             END DO
             ! --- Ice enthalpy remapping --- !
-            CALL ice_thd_ent( ze_i_2d(1:npti,:,jl) ) 
+            CALL ice_thd_ent( ze_i_2d(1:npti,:,jl), .false. ) 
          END DO
 
Index: NEMO/branches/2020/dev_r12512_HPC-04_mcastril_Mixed_Precision_implementation/src/ICE/icethd_ent.F90
===================================================================
--- NEMO/branches/2020/dev_r12512_HPC-04_mcastril_Mixed_Precision_implementation/src/ICE/icethd_ent.F90	(revision 12632)
+++ NEMO/branches/2020/dev_r12512_HPC-04_mcastril_Mixed_Precision_implementation/src/ICE/icethd_ent.F90	(revision 12701)
@@ -38,5 +38,5 @@
 CONTAINS
  
-   SUBROUTINE ice_thd_ent( qnew )
+   SUBROUTINE ice_thd_ent( qnew, compute_hfx_err )
       !!-------------------------------------------------------------------
       !!               ***   ROUTINE ice_thd_ent  ***
@@ -64,4 +64,6 @@
       !!-------------------------------------------------------------------
       REAL(wp), DIMENSION(:,:), INTENT(inout) ::   qnew             ! new enthlapies (J.m-3, remapped)
+      LOGICAL, INTENT(in)                     ::   compute_hfx_err  ! determines whether to compute diag.
+                                                                    ! error or not
       !
       INTEGER  :: ji         !  dummy loop indices
@@ -128,9 +130,11 @@
       ! comment: if input h_i_old and eh_i_old are already multiplied by a_i (as in icethd_do), 
       ! then we should not (* a_i) again but not important since this is just to check that remap error is ~0
-      DO ji = 1, npti
-         hfx_err_rem_1d(ji) = hfx_err_rem_1d(ji) + a_i_1d(ji) * r1_Dt_ice *  &
-            &               ( SUM( qnew(ji,1:nlay_i) ) * zhnew(ji) - SUM( eh_i_old(ji,0:nlay_i+1) ) ) 
-      END DO
-      
+      IF( compute_hfx_err ) THEN
+         DO ji = 1, npti
+            hfx_err_rem_1d(ji) = hfx_err_rem_1d(ji) + a_i_1d(ji) * r1_rdtice *  &
+               &               ( SUM( qnew(ji,1:nlay_i) ) * zhnew(ji) - SUM( eh_i_old(ji,0:nlay_i+1) ) )
+         END DO
+      END IF
+ 
    END SUBROUTINE ice_thd_ent
 
Index: NEMO/branches/2020/dev_r12512_HPC-04_mcastril_Mixed_Precision_implementation/src/OCE/ZDF/zdftke.F90
===================================================================
--- NEMO/branches/2020/dev_r12512_HPC-04_mcastril_Mixed_Precision_implementation/src/OCE/ZDF/zdftke.F90	(revision 12632)
+++ NEMO/branches/2020/dev_r12512_HPC-04_mcastril_Mixed_Precision_implementation/src/OCE/ZDF/zdftke.F90	(revision 12701)
@@ -305,5 +305,9 @@
          DO_3D_00_00( 2, jpkm1 )
             !                             ! local Richardson number
-            zri = MAX( rn2b(ji,jj,jk), 0._wp ) * p_avm(ji,jj,jk) / ( p_sh2(ji,jj,jk) + rn_bshear )
+            IF (rn2b(ji,jj,jk) <= 0.0_wp) then
+                zri = 0.0_wp
+            ELSE
+                zri = rn2b(ji,jj,jk) * p_avm(ji,jj,jk) / ( p_sh2(ji,jj,jk) + rn_bshear )
+            ENDIF
             !                             ! inverse of Prandtl number
             apdlr(ji,jj,jk) = MAX(  0.1_wp,  ri_cri / MAX( ri_cri , zri )  )
