diff --git a/src/ICE/icedyn_adv.F90 b/src/ICE/icedyn_adv.F90
index a6d2e957735e1c0663619f52f663e237f61c3215..6afabac00bcf4adc66c8cbd98cfb577df3151843 100644
--- a/src/ICE/icedyn_adv.F90
+++ b/src/ICE/icedyn_adv.F90
@@ -25,7 +25,6 @@ MODULE icedyn_adv
    USE lib_mpp        ! MPP library
    USE lib_fortran    ! fortran utilities (glob_sum + no signed zero)
    USE timing         ! Timing
-   USE prtctl         ! Print control
 
    IMPLICIT NONE
    PRIVATE
@@ -108,6 +107,8 @@ CONTAINS
       ! controls
       IF( ln_icediachk )   CALL ice_cons_hsm(1, 'icedyn_adv', rdiag_v, rdiag_s, rdiag_t, rdiag_fv, rdiag_fs, rdiag_ft) ! conservation
       IF( ln_icectl    )   CALL ice_prt     (kt, iiceprt, jiceprt,-1, ' - ice dyn & trp - ')                           ! prints
+      IF( sn_cfctl%l_prtctl )   &
+         &                 CALL ice_prt3D('icedyn_adv')                                                                ! prints
       IF( ln_timing    )   CALL timing_stop ('icedyn_adv')                                                             ! timing
       !
    END SUBROUTINE ice_dyn_adv
diff --git a/src/ICE/icedyn_adv_umx.F90 b/src/ICE/icedyn_adv_umx.F90
index de133b1e9f26fa896d94f4ca9f2b1be3d133f57a..f15c2349f2a5f7fa8df0ab4d6493a3f8b2ffc86a 100644
--- a/src/ICE/icedyn_adv_umx.F90
+++ b/src/ICE/icedyn_adv_umx.F90
@@ -375,7 +375,8 @@ CONTAINS
          zati2(:,:) = SUM( pa_i(:,:,:), dim=3 )
          DO_2D( 0, 0, 0, 0 )
             pato_i(ji,jj) = pato_i(ji,jj) - ( zati2(ji,jj) - zati1(ji,jj) ) &
-               &                          - ( zudy(ji,jj) - zudy(ji-1,jj) + zvdx(ji,jj) - zvdx(ji,jj-1) ) * r1_e1e2t(ji,jj) * zdt
+               &                          - (   ( zudy(ji,jj) - zudy(ji-1,jj) ) &   ! ad () for NP repro
+               &                              + ( zvdx(ji,jj) - zvdx(ji,jj-1) )   ) * r1_e1e2t(ji,jj) * zdt
          END_2D
          CALL lbc_lnk( 'icedyn_adv_umx', pato_i, 'T',  1._wp )
          !
diff --git a/src/ICE/icedyn_rdgrft.F90 b/src/ICE/icedyn_rdgrft.F90
index 0e21ffeffd18c7ab4e112365e4f025721223b7f0..395fa5a7ef941f0407a3b72df3b3f3dcd0680b9b 100644
--- a/src/ICE/icedyn_rdgrft.F90
+++ b/src/ICE/icedyn_rdgrft.F90
@@ -990,8 +990,8 @@ CONTAINS
          DO_2D( 0, 0, 0, 0 )
             IF ( SUM( a_i(ji,jj,:) ) > 0._wp ) THEN
                zworka(ji,jj) = ( 4._wp * strength(ji,jj)              &
-                  &                    + ( strength(ji-1,jj) * tmask(ji-1,jj,1) + strength(ji+1,jj) * tmask(ji+1,jj,1) ) &
-                  &                    + ( strength(ji,jj-1) * tmask(ji,jj-1,1) + strength(ji,jj+1) * tmask(ji,jj+1,1) ) &
+                  &                    + ( ( strength(ji-1,jj) * tmask(ji-1,jj,1) + strength(ji+1,jj) * tmask(ji+1,jj,1) ) &
+                  &                      + ( strength(ji,jj-1) * tmask(ji,jj-1,1) + strength(ji,jj+1) * tmask(ji,jj+1,1) ) ) &
                   &            ) / ( 4._wp + tmask(ji-1,jj,1) + tmask(ji+1,jj,1) + tmask(ji,jj-1,1) + tmask(ji,jj+1,1) )
             ELSE
                zworka(ji,jj) = 0._wp
diff --git a/src/ICE/icedyn_rhg_evp.F90 b/src/ICE/icedyn_rhg_evp.F90
index 6c1fd8aa3e57ffe33670697983378691838a19d8..136f13d122b5d03a3416906f729bae42f3befb41 100644
--- a/src/ICE/icedyn_rhg_evp.F90
+++ b/src/ICE/icedyn_rhg_evp.F90
@@ -470,17 +470,17 @@ CONTAINS
          ! (brackets added to fix the order of floating point operations for halo 1 - halo 2 compatibility)
          DO_2D( nn_hls-1, nn_hls-1, nn_hls-1, nn_hls-1 )
             !                   !--- U points
-            zfU(ji,jj) = 0.5_wp * ( (( zs1(ji+1,jj) - zs1(ji,jj) ) * e2u(ji,jj)                                             &
-               &                  + ( zs2(ji+1,jj) * e2t(ji+1,jj) * e2t(ji+1,jj) - zs2(ji,jj) * e2t(ji,jj) * e2t(ji,jj)    &
-               &                    ) * r1_e2u(ji,jj))                                                                      &
+            zfU(ji,jj) = 0.5_wp * ( ( ( zs1(ji+1,jj) - zs1(ji,jj) ) * e2u(ji,jj)                                             &
+               &                    + ( zs2(ji+1,jj) * e2t(ji+1,jj) * e2t(ji+1,jj) - zs2(ji,jj) * e2t(ji,jj) * e2t(ji,jj)    &
+               &                      ) * r1_e2u(ji,jj) )                                                                      &
                &                  + ( zs12(ji,jj) * e1f(ji,jj) * e1f(ji,jj) - zs12(ji,jj-1) * e1f(ji,jj-1) * e1f(ji,jj-1)  &
                &                    ) * 2._wp * r1_e1u(ji,jj)                                                              &
                &                  ) * r1_e1e2u(ji,jj)
             !
             !                !--- V points
-            zfV(ji,jj) = 0.5_wp * ( (( zs1(ji,jj+1) - zs1(ji,jj) ) * e1v(ji,jj)                                             &
-               &                  - ( zs2(ji,jj+1) * e1t(ji,jj+1) * e1t(ji,jj+1) - zs2(ji,jj) * e1t(ji,jj) * e1t(ji,jj)    &
-               &                    ) * r1_e1v(ji,jj))                                                                      &
+            zfV(ji,jj) = 0.5_wp * ( ( ( zs1(ji,jj+1) - zs1(ji,jj) ) * e1v(ji,jj)                                             &
+               &                    - ( zs2(ji,jj+1) * e1t(ji,jj+1) * e1t(ji,jj+1) - zs2(ji,jj) * e1t(ji,jj) * e1t(ji,jj)    &
+               &                      ) * r1_e1v(ji,jj) )                                                                      &
                &                  + ( zs12(ji,jj) * e2f(ji,jj) * e2f(ji,jj) - zs12(ji-1,jj) * e2f(ji-1,jj) * e2f(ji-1,jj)  &
                &                    ) * 2._wp * r1_e2v(ji,jj)                                                              &
                &                  ) * r1_e1e2v(ji,jj)
@@ -754,8 +754,8 @@ CONTAINS
          zten_i(ji,jj) = zdt
 
          ! shear**2 at T points (doc eq. A16)
-         zds2 = ( zds(ji,jj  ) * zds(ji,jj  ) * e1e2f(ji,jj  ) + zds(ji-1,jj  ) * zds(ji-1,jj  ) * e1e2f(ji-1,jj  )  &
-            &   + zds(ji,jj-1) * zds(ji,jj-1) * e1e2f(ji,jj-1) + zds(ji-1,jj-1) * zds(ji-1,jj-1) * e1e2f(ji-1,jj-1)  &
+         zds2 = ( ( zds(ji,jj  ) * zds(ji,jj  ) * e1e2f(ji,jj  ) + zds(ji-1,jj  ) * zds(ji-1,jj  ) * e1e2f(ji-1,jj  ) )  & ! add () 
+            &   + ( zds(ji,jj-1) * zds(ji,jj-1) * e1e2f(ji,jj-1) + zds(ji-1,jj-1) * zds(ji-1,jj-1) * e1e2f(ji-1,jj-1) )  & ! NP rep
             &   ) * 0.25_wp * r1_e1e2t(ji,jj)
 
          ! maximum shear rate at T points (includes tension, output only)
@@ -765,8 +765,8 @@ CONTAINS
          zshear(ji,jj)   = SQRT( zds2 ) * zmsk(ji,jj)
 
          ! divergence at T points
-         pdivu_i(ji,jj) = ( e2u(ji,jj) * u_ice(ji,jj) - e2u(ji-1,jj) * u_ice(ji-1,jj)   &
-            &             + e1v(ji,jj) * v_ice(ji,jj) - e1v(ji,jj-1) * v_ice(ji,jj-1)   &
+         pdivu_i(ji,jj) = ( ( e2u(ji,jj) * u_ice(ji,jj) - e2u(ji-1,jj) * u_ice(ji-1,jj) )  & ! add () for NP repro
+            &             + ( e1v(ji,jj) * v_ice(ji,jj) - e1v(ji,jj-1) * v_ice(ji,jj-1) )  &
             &             ) * r1_e1e2t(ji,jj) * zmsk(ji,jj)
 
          ! delta at T points
diff --git a/src/OCE/DYN/dynspg_ts.F90 b/src/OCE/DYN/dynspg_ts.F90
index bd40bf5409e9d9b034118816df9ad6f59fc5cf35..ddb60a5068b70a06066150a1e3ef41f2e2f7edc0 100644
--- a/src/OCE/DYN/dynspg_ts.F90
+++ b/src/OCE/DYN/dynspg_ts.F90
@@ -559,7 +559,7 @@ CONTAINS
          !--        ssh    =  ssh   - delta_t' * [ frc + div( flux      ) ]      --!
          !-------------------------------------------------------------------------!
          DO_2D( 0, 0, 0, 0 )
-            zhdiv = (   zhU(ji,jj) - zhU(ji-1,jj) + zhV(ji,jj) - zhV(ji,jj-1)   ) * r1_e1e2t(ji,jj)
+            zhdiv = (   ( zhU(ji,jj) - zhU(ji-1,jj) ) + ( zhV(ji,jj) - zhV(ji,jj-1) )  ) * r1_e1e2t(ji,jj)
             ssha_e(ji,jj) = (  sshn_e(ji,jj) - rDt_e * ( ssh_frc(ji,jj) + zhdiv )  ) * ssmask(ji,jj)
          END_2D
          !
diff --git a/src/OCE/IOM/prtctl.F90 b/src/OCE/IOM/prtctl.F90
index 8601142c5050693521e02d703e79360cab20b573..bb00edc5b9521d1caa08311ffaafe90d4571cc5a 100644
--- a/src/OCE/IOM/prtctl.F90
+++ b/src/OCE/IOM/prtctl.F90
@@ -132,8 +132,9 @@ CONTAINS
       CHARACTER(len=*)                    , INTENT(in), OPTIONAL ::   clinfo3
       INTEGER                             , INTENT(in), OPTIONAL ::   kdim
       !
-      CHARACTER(len=30) :: cl1, cl2
+      CHARACTER(len=30) :: cl1, cl2, cl3
       CHARACTER(len=6) :: clfmt
+      CHARACTER(len=1) :: cli1
       INTEGER ::  jn, jl, kdir
       INTEGER ::  ipi1, ipi2, ipim1, ipim2
       INTEGER ::  isht1, isht2, ishtm1, ishtm2
@@ -299,12 +300,28 @@ CONTAINS
                ENDIF
                ! replace .false. by .true. to switch on theses prints of the last inner line
                IF( .FALSE. .AND. l_IdoNFold .AND. jj1e == Nje0 - isht1 ) THEN
-                  IF( PRESENT(tab2d_1) ) WRITE(inum, *) 'Last line '//TRIM(cl1)//' ', tab2d_1(ii1s:ii1e,jj1e)
-                  IF( PRESENT(tab3d_1) ) WRITE(inum, *) 'Last line '//TRIM(cl1)//' ', tab3d_1(ii1s:ii1e,jj1e,1:kdir)
+                  IF( PRESENT(tab2d_1) ) THEN
+                     WRITE(cli1, '(i1)') INT(LOG10(REAL(ii1e-ii1s+1,wp))) + 1     ! how many digits to we need to write ?
+                     WRITE(cl3, "(i"//cli1//")") ii1e-ii1s+1
+                     WRITE(inum, "(a,"//TRIM(cl3)//clfmt//")") 'Last line '//TRIM(cl1)//' ', tab2d_1(ii1s:ii1e,jj1e)
+                  ENDIF
+                  IF( PRESENT(tab3d_1) ) THEN
+                     WRITE(cli1, '(i1)') INT(LOG10(REAL((ii1e-ii1s+1)*kdir,wp))) + 1     ! how many digits to we need to write ?
+                     WRITE(cl3, "(i"//cli1//")") (ii1e-ii1s+1)*kdir
+                     WRITE(inum, "(a,"//TRIM(cl3)//clfmt//")") 'Last line '//TRIM(cl1)//' ', tab3d_1(ii1s:ii1e,jj1e,1:kdir)
+                  ENDIF
                ENDIF
                IF( .FALSE. .AND. l_IdoNFold .AND. jj2e == Nje0 - isht2 ) THEN
-                  IF( PRESENT(tab2d_2) ) WRITE(inum, *) 'Last line '//TRIM(cl2)//' ', tab2d_2(ii2s:ii2e,jj2e)
-                  IF( PRESENT(tab3d_2) ) WRITE(inum, *) 'Last line '//TRIM(cl2)//' ', tab3d_2(ii2s:ii2e,jj2e,1:kdir)
+                  IF( PRESENT(tab2d_2) ) THEN
+                     WRITE(cli1, '(i1)') INT(LOG10(REAL(ii2e-ii2s+1,wp))) + 1     ! how many digits to we need to write ?
+                     WRITE(cl3, "(i"//cli1//")") ii2e-ii2s+1
+                     WRITE(inum, "(a,"//TRIM(cl3)//clfmt//")") 'Last line '//TRIM(cl2)//' ', tab2d_2(ii2s:ii2e,jj2e)
+                  ENDIF
+                  IF( PRESENT(tab3d_2) ) THEN
+                     WRITE(cli1, '(i1)') INT(LOG10(REAL((ii2e-ii2s+1)*kdir,wp))) + 1     ! how many digits to we need to write ?
+                     WRITE(cl3, "(i"//cli1//")") (ii2e-ii2s+1)*kdir
+                     WRITE(inum, "(a,"//TRIM(cl3)//clfmt//")") 'Last line '//TRIM(cl2)//' ', tab3d_2(ii2s:ii2e,jj2e,1:kdir)
+                  ENDIF
                ENDIF
             END DO
          ENDIF