Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/ICE/iceistate.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/ICE/iceistate.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/ICE/iceistate.F90	(revision 13026)
@@ -32,4 +32,10 @@
    USE lib_fortran    ! fortran utilities (glob_sum + no signed zero)
    USE fldread        ! read input fields
+
+# if defined key_agrif
+   USE agrif_oce
+   USE agrif_ice
+   USE agrif_ice_interp 
+# endif   
 
    IMPLICIT NONE
@@ -168,198 +174,242 @@
       ! 2) overwrite some of the fields with namelist parameters or netcdf file
       !------------------------------------------------------------------------
+
+
       IF( ln_iceini ) THEN
          !                             !---------------!
-         IF( ln_iceini_file )THEN      ! Read a file   !
-            !                          !---------------!
-            WHERE( ff_t(:,:) >= 0._wp )   ;   zswitch(:,:) = 1._wp
-            ELSEWHERE                     ;   zswitch(:,:) = 0._wp
+         
+         IF( Agrif_Root() ) THEN
+
+            IF( ln_iceini_file )THEN      ! Read a file   !
+               !                          !---------------!
+               WHERE( ff_t(:,:) >= 0._wp )   ;   zswitch(:,:) = 1._wp
+               ELSEWHERE                     ;   zswitch(:,:) = 0._wp
+               END WHERE
+               !
+               CALL fld_read( kt, 1, si ) ! input fields provided at the current time-step
+               !
+               ! -- mandatory fields -- !
+               zht_i_ini(:,:) = si(jp_hti)%fnow(:,:,1) * tmask(:,:,1)
+               zht_s_ini(:,:) = si(jp_hts)%fnow(:,:,1) * tmask(:,:,1)
+               zat_i_ini(:,:) = si(jp_ati)%fnow(:,:,1) * tmask(:,:,1)
+
+               ! -- optional fields -- !
+               !    if fields do not exist then set them to the values present in the namelist (except for snow and surface temperature)
+               !
+               ! ice salinity
+               IF( TRIM(si(jp_smi)%clrootname) == 'NOT USED' ) &
+                  &     si(jp_smi)%fnow(:,:,1) = ( rn_smi_ini_n * zswitch + rn_smi_ini_s * (1._wp - zswitch) ) * tmask(:,:,1)
+               !
+               ! temperatures
+               IF    ( TRIM(si(jp_tmi)%clrootname) == 'NOT USED' .AND. TRIM(si(jp_tsu)%clrootname) == 'NOT USED' .AND. &
+                  &    TRIM(si(jp_tms)%clrootname) == 'NOT USED' ) THEN
+                  si(jp_tmi)%fnow(:,:,1) = ( rn_tmi_ini_n * zswitch + rn_tmi_ini_s * (1._wp - zswitch) ) * tmask(:,:,1)
+                  si(jp_tsu)%fnow(:,:,1) = ( rn_tsu_ini_n * zswitch + rn_tsu_ini_s * (1._wp - zswitch) ) * tmask(:,:,1)
+                  si(jp_tms)%fnow(:,:,1) = ( rn_tms_ini_n * zswitch + rn_tms_ini_s * (1._wp - zswitch) ) * tmask(:,:,1)
+               ELSEIF( TRIM(si(jp_tmi)%clrootname) == 'NOT USED' .AND. TRIM(si(jp_tms)%clrootname) /= 'NOT USED' ) THEN ! if T_s is read and not T_i, set T_i = (T_s + T_freeze)/2
+                  si(jp_tmi)%fnow(:,:,1) = 0.5_wp * ( si(jp_tms)%fnow(:,:,1) + 271.15 )
+               ELSEIF( TRIM(si(jp_tmi)%clrootname) == 'NOT USED' .AND. TRIM(si(jp_tsu)%clrootname) /= 'NOT USED' ) THEN ! if T_su is read and not T_i, set T_i = (T_su + T_freeze)/2
+                  si(jp_tmi)%fnow(:,:,1) = 0.5_wp * ( si(jp_tsu)%fnow(:,:,1) + 271.15 )
+               ELSEIF( TRIM(si(jp_tsu)%clrootname) == 'NOT USED' .AND. TRIM(si(jp_tms)%clrootname) /= 'NOT USED' ) THEN ! if T_s is read and not T_su, set T_su = T_s
+                  si(jp_tsu)%fnow(:,:,1) = si(jp_tms)%fnow(:,:,1)
+               ELSEIF( TRIM(si(jp_tsu)%clrootname) == 'NOT USED' .AND. TRIM(si(jp_tmi)%clrootname) /= 'NOT USED' ) THEN ! if T_i is read and not T_su, set T_su = T_i
+                  si(jp_tsu)%fnow(:,:,1) = si(jp_tmi)%fnow(:,:,1)
+               ELSEIF( TRIM(si(jp_tms)%clrootname) == 'NOT USED' .AND. TRIM(si(jp_tsu)%clrootname) /= 'NOT USED' ) THEN ! if T_su is read and not T_s, set T_s = T_su
+                  si(jp_tms)%fnow(:,:,1) = si(jp_tsu)%fnow(:,:,1)
+               ELSEIF( TRIM(si(jp_tms)%clrootname) == 'NOT USED' .AND. TRIM(si(jp_tmi)%clrootname) /= 'NOT USED' ) THEN ! if T_i is read and not T_s, set T_s = T_i
+                  si(jp_tms)%fnow(:,:,1) = si(jp_tmi)%fnow(:,:,1)
+               ENDIF
+               !
+               ! pond concentration
+               IF( TRIM(si(jp_apd)%clrootname) == 'NOT USED' ) &
+                  &     si(jp_apd)%fnow(:,:,1) = ( rn_apd_ini_n * zswitch + rn_apd_ini_s * (1._wp - zswitch) ) * tmask(:,:,1) & ! rn_apd = pond fraction => rn_apnd * a_i = pond conc.
+                  &                              * si(jp_ati)%fnow(:,:,1) 
+               !
+               ! pond depth
+               IF( TRIM(si(jp_hpd)%clrootname) == 'NOT USED' ) &
+                  &     si(jp_hpd)%fnow(:,:,1) = ( rn_hpd_ini_n * zswitch + rn_hpd_ini_s * (1._wp - zswitch) ) * tmask(:,:,1)
+               !
+               zsm_i_ini(:,:) = si(jp_smi)%fnow(:,:,1) * tmask(:,:,1)
+               ztm_i_ini(:,:) = si(jp_tmi)%fnow(:,:,1) * tmask(:,:,1)
+               zt_su_ini(:,:) = si(jp_tsu)%fnow(:,:,1) * tmask(:,:,1)
+               ztm_s_ini(:,:) = si(jp_tms)%fnow(:,:,1) * tmask(:,:,1)
+               zapnd_ini(:,:) = si(jp_apd)%fnow(:,:,1) * tmask(:,:,1)
+               zhpnd_ini(:,:) = si(jp_hpd)%fnow(:,:,1) * tmask(:,:,1)
+               !
+               ! change the switch for the following
+               WHERE( zat_i_ini(:,:) > 0._wp )   ;   zswitch(:,:) = tmask(:,:,1) 
+               ELSEWHERE                         ;   zswitch(:,:) = 0._wp
+               END WHERE
+
+               !                          !---------------!
+            ELSE                          ! Read namelist !
+               !                          !---------------!
+               ! no ice if (sst - Tfreez) >= thresold
+               WHERE( ( sst_m(:,:) - (t_bo(:,:) - rt0) ) * tmask(:,:,1) >= rn_thres_sst )   ;   zswitch(:,:) = 0._wp 
+               ELSEWHERE                                                                    ;   zswitch(:,:) = tmask(:,:,1)
+               END WHERE
+               !
+               ! assign initial thickness, concentration, snow depth and salinity to an hemisphere-dependent array
+               WHERE( ff_t(:,:) >= 0._wp )
+                  zht_i_ini(:,:) = rn_hti_ini_n * zswitch(:,:)
+                  zht_s_ini(:,:) = rn_hts_ini_n * zswitch(:,:)
+                  zat_i_ini(:,:) = rn_ati_ini_n * zswitch(:,:)
+                  zsm_i_ini(:,:) = rn_smi_ini_n * zswitch(:,:)
+                  ztm_i_ini(:,:) = rn_tmi_ini_n * zswitch(:,:)
+                  zt_su_ini(:,:) = rn_tsu_ini_n * zswitch(:,:)
+                  ztm_s_ini(:,:) = rn_tms_ini_n * zswitch(:,:)
+                  zapnd_ini(:,:) = rn_apd_ini_n * zswitch(:,:) * zat_i_ini(:,:) ! rn_apd = pond fraction => rn_apd * a_i = pond conc. 
+                  zhpnd_ini(:,:) = rn_hpd_ini_n * zswitch(:,:)
+               ELSEWHERE
+                  zht_i_ini(:,:) = rn_hti_ini_s * zswitch(:,:)
+                  zht_s_ini(:,:) = rn_hts_ini_s * zswitch(:,:)
+                  zat_i_ini(:,:) = rn_ati_ini_s * zswitch(:,:)
+                  zsm_i_ini(:,:) = rn_smi_ini_s * zswitch(:,:)
+                  ztm_i_ini(:,:) = rn_tmi_ini_s * zswitch(:,:)
+                  zt_su_ini(:,:) = rn_tsu_ini_s * zswitch(:,:)
+                  ztm_s_ini(:,:) = rn_tms_ini_s * zswitch(:,:)
+                  zapnd_ini(:,:) = rn_apd_ini_s * zswitch(:,:) * zat_i_ini(:,:) ! rn_apd = pond fraction => rn_apd * a_i = pond conc.
+                  zhpnd_ini(:,:) = rn_hpd_ini_s * zswitch(:,:)
+               END WHERE
+               !
+            ENDIF
+
+
+
+            ! make sure ponds = 0 if no ponds scheme
+            IF ( .NOT.ln_pnd ) THEN
+               zapnd_ini(:,:) = 0._wp
+               zhpnd_ini(:,:) = 0._wp
+            ENDIF
+            
+            !-------------!
+            ! fill fields !
+            !-------------!
+            ! select ice covered grid points
+            npti = 0 ; nptidx(:) = 0
+            DO_2D_11_11
+               IF ( zht_i_ini(ji,jj) > 0._wp ) THEN
+                  npti         = npti  + 1
+                  nptidx(npti) = (jj - 1) * jpi + ji
+               ENDIF
+            END_2D
+
+            ! move to 1D arrays: (jpi,jpj) -> (jpi*jpj)
+            CALL tab_2d_1d( npti, nptidx(1:npti), h_i_1d (1:npti)  , zht_i_ini )
+            CALL tab_2d_1d( npti, nptidx(1:npti), h_s_1d (1:npti)  , zht_s_ini )
+            CALL tab_2d_1d( npti, nptidx(1:npti), at_i_1d(1:npti)  , zat_i_ini )
+            CALL tab_2d_1d( npti, nptidx(1:npti), t_i_1d (1:npti,1), ztm_i_ini )
+            CALL tab_2d_1d( npti, nptidx(1:npti), t_s_1d (1:npti,1), ztm_s_ini )
+            CALL tab_2d_1d( npti, nptidx(1:npti), t_su_1d(1:npti)  , zt_su_ini )
+            CALL tab_2d_1d( npti, nptidx(1:npti), s_i_1d (1:npti)  , zsm_i_ini )
+            CALL tab_2d_1d( npti, nptidx(1:npti), a_ip_1d(1:npti)  , zapnd_ini )
+            CALL tab_2d_1d( npti, nptidx(1:npti), h_ip_1d(1:npti)  , zhpnd_ini )
+
+            ! allocate temporary arrays
+            ALLOCATE( zhi_2d(npti,jpl), zhs_2d(npti,jpl), zai_2d (npti,jpl), &
+               &      zti_2d(npti,jpl), zts_2d(npti,jpl), ztsu_2d(npti,jpl), zsi_2d(npti,jpl), zaip_2d(npti,jpl), zhip_2d(npti,jpl) )
+            
+            ! distribute 1-cat into jpl-cat: (jpi*jpj) -> (jpi*jpj,jpl)
+            CALL ice_var_itd( h_i_1d(1:npti)  , h_s_1d(1:npti)  , at_i_1d(1:npti),                                                   &
+               &              zhi_2d          , zhs_2d          , zai_2d         ,                                                   &
+               &              t_i_1d(1:npti,1), t_s_1d(1:npti,1), t_su_1d(1:npti), s_i_1d(1:npti), a_ip_1d(1:npti), h_ip_1d(1:npti), &
+               &              zti_2d          , zts_2d          , ztsu_2d        , zsi_2d        , zaip_2d        , zhip_2d )
+
+            ! move to 3D arrays: (jpi*jpj,jpl) -> (jpi,jpj,jpl)
+            DO jl = 1, jpl
+               zti_3d(:,:,jl) = rt0 * tmask(:,:,1)
+               zts_3d(:,:,jl) = rt0 * tmask(:,:,1)
+            END DO
+            CALL tab_2d_3d( npti, nptidx(1:npti), zhi_2d   , h_i    )
+            CALL tab_2d_3d( npti, nptidx(1:npti), zhs_2d   , h_s    )
+            CALL tab_2d_3d( npti, nptidx(1:npti), zai_2d   , a_i    )
+            CALL tab_2d_3d( npti, nptidx(1:npti), zti_2d   , zti_3d )
+            CALL tab_2d_3d( npti, nptidx(1:npti), zts_2d   , zts_3d )
+            CALL tab_2d_3d( npti, nptidx(1:npti), ztsu_2d  , t_su   )
+            CALL tab_2d_3d( npti, nptidx(1:npti), zsi_2d   , s_i    )
+            CALL tab_2d_3d( npti, nptidx(1:npti), zaip_2d  , a_ip   )
+            CALL tab_2d_3d( npti, nptidx(1:npti), zhip_2d  , h_ip   )
+
+            ! deallocate temporary arrays
+            DEALLOCATE( zhi_2d, zhs_2d, zai_2d , &
+               &        zti_2d, zts_2d, ztsu_2d, zsi_2d, zaip_2d, zhip_2d )
+
+            ! calculate extensive and intensive variables
+            CALL ice_var_salprof ! for sz_i
+            DO jl = 1, jpl
+               DO_2D_11_11
+                  v_i (ji,jj,jl) = h_i(ji,jj,jl) * a_i(ji,jj,jl)
+                  v_s (ji,jj,jl) = h_s(ji,jj,jl) * a_i(ji,jj,jl)
+                  sv_i(ji,jj,jl) = MIN( MAX( rn_simin , s_i(ji,jj,jl) ) , rn_simax ) * v_i(ji,jj,jl)
+               END_2D
+            END DO
+            !
+            DO jl = 1, jpl
+               DO_3D_11_11( 1, nlay_s )
+                  t_s(ji,jj,jk,jl) = zts_3d(ji,jj,jl)
+                  e_s(ji,jj,jk,jl) = zswitch(ji,jj) * v_s(ji,jj,jl) * r1_nlay_s * &
+                     &               rhos * ( rcpi * ( rt0 - t_s(ji,jj,jk,jl) ) + rLfus )
+               END_3D
+            END DO
+            !
+            DO jl = 1, jpl
+               DO_3D_11_11( 1, nlay_i )
+                  t_i (ji,jj,jk,jl) = zti_3d(ji,jj,jl) 
+                  ztmelts          = - rTmlt * sz_i(ji,jj,jk,jl) + rt0 ! melting temperature in K
+                  e_i(ji,jj,jk,jl) = zswitch(ji,jj) * v_i(ji,jj,jl) * r1_nlay_i * &
+                     &               rhoi * (  rcpi  * ( ztmelts - t_i(ji,jj,jk,jl) ) + &
+                     &                         rLfus * ( 1._wp - (ztmelts-rt0) / MIN( (t_i(ji,jj,jk,jl)-rt0), -epsi20 ) ) &
+                     &                       - rcp   * ( ztmelts - rt0 ) )
+               END_3D
+            END DO
+
+            ! Melt ponds
+            WHERE( a_i > epsi10 )
+               a_ip_frac(:,:,:) = a_ip(:,:,:) / a_i(:,:,:)
+            ELSEWHERE
+               a_ip_frac(:,:,:) = 0._wp
             END WHERE
+            v_ip(:,:,:) = h_ip(:,:,:) * a_ip(:,:,:)
+             
+            ! specific temperatures for coupled runs
+            tn_ice(:,:,:) = t_su(:,:,:)
+            t1_ice(:,:,:) = t_i (:,:,1,:)
             !
-            CALL fld_read( kt, 1, si ) ! input fields provided at the current time-step
-            !
-            ! -- mandatory fields -- !
-            zht_i_ini(:,:) = si(jp_hti)%fnow(:,:,1) * tmask(:,:,1)
-            zht_s_ini(:,:) = si(jp_hts)%fnow(:,:,1) * tmask(:,:,1)
-            zat_i_ini(:,:) = si(jp_ati)%fnow(:,:,1) * tmask(:,:,1)
-
-            ! -- optional fields -- !
-            !    if fields do not exist then set them to the values present in the namelist (except for snow and surface temperature)
-            !
-            ! ice salinity
-            IF( TRIM(si(jp_smi)%clrootname) == 'NOT USED' ) &
-               &     si(jp_smi)%fnow(:,:,1) = ( rn_smi_ini_n * zswitch + rn_smi_ini_s * (1._wp - zswitch) ) * tmask(:,:,1)
-            !
-            ! temperatures
-            IF    ( TRIM(si(jp_tmi)%clrootname) == 'NOT USED' .AND. TRIM(si(jp_tsu)%clrootname) == 'NOT USED' .AND. &
-               &    TRIM(si(jp_tms)%clrootname) == 'NOT USED' ) THEN
-               si(jp_tmi)%fnow(:,:,1) = ( rn_tmi_ini_n * zswitch + rn_tmi_ini_s * (1._wp - zswitch) ) * tmask(:,:,1)
-               si(jp_tsu)%fnow(:,:,1) = ( rn_tsu_ini_n * zswitch + rn_tsu_ini_s * (1._wp - zswitch) ) * tmask(:,:,1)
-               si(jp_tms)%fnow(:,:,1) = ( rn_tms_ini_n * zswitch + rn_tms_ini_s * (1._wp - zswitch) ) * tmask(:,:,1)
-            ELSEIF( TRIM(si(jp_tmi)%clrootname) == 'NOT USED' .AND. TRIM(si(jp_tms)%clrootname) /= 'NOT USED' ) THEN ! if T_s is read and not T_i, set T_i = (T_s + T_freeze)/2
-               si(jp_tmi)%fnow(:,:,1) = 0.5_wp * ( si(jp_tms)%fnow(:,:,1) + 271.15 )
-            ELSEIF( TRIM(si(jp_tmi)%clrootname) == 'NOT USED' .AND. TRIM(si(jp_tsu)%clrootname) /= 'NOT USED' ) THEN ! if T_su is read and not T_i, set T_i = (T_su + T_freeze)/2
-               si(jp_tmi)%fnow(:,:,1) = 0.5_wp * ( si(jp_tsu)%fnow(:,:,1) + 271.15 )
-            ELSEIF( TRIM(si(jp_tsu)%clrootname) == 'NOT USED' .AND. TRIM(si(jp_tms)%clrootname) /= 'NOT USED' ) THEN ! if T_s is read and not T_su, set T_su = T_s
-               si(jp_tsu)%fnow(:,:,1) = si(jp_tms)%fnow(:,:,1)
-            ELSEIF( TRIM(si(jp_tsu)%clrootname) == 'NOT USED' .AND. TRIM(si(jp_tmi)%clrootname) /= 'NOT USED' ) THEN ! if T_i is read and not T_su, set T_su = T_i
-               si(jp_tsu)%fnow(:,:,1) = si(jp_tmi)%fnow(:,:,1)
-            ELSEIF( TRIM(si(jp_tms)%clrootname) == 'NOT USED' .AND. TRIM(si(jp_tsu)%clrootname) /= 'NOT USED' ) THEN ! if T_su is read and not T_s, set T_s = T_su
-               si(jp_tms)%fnow(:,:,1) = si(jp_tsu)%fnow(:,:,1)
-            ELSEIF( TRIM(si(jp_tms)%clrootname) == 'NOT USED' .AND. TRIM(si(jp_tmi)%clrootname) /= 'NOT USED' ) THEN ! if T_i is read and not T_s, set T_s = T_i
-               si(jp_tms)%fnow(:,:,1) = si(jp_tmi)%fnow(:,:,1)
-            ENDIF
-            !
-            ! pond concentration
-            IF( TRIM(si(jp_apd)%clrootname) == 'NOT USED' ) &
-               &     si(jp_apd)%fnow(:,:,1) = ( rn_apd_ini_n * zswitch + rn_apd_ini_s * (1._wp - zswitch) ) * tmask(:,:,1) & ! rn_apd = pond fraction => rn_apnd * a_i = pond conc.
-               &                              * si(jp_ati)%fnow(:,:,1) 
-            !
-            ! pond depth
-            IF( TRIM(si(jp_hpd)%clrootname) == 'NOT USED' ) &
-               &     si(jp_hpd)%fnow(:,:,1) = ( rn_hpd_ini_n * zswitch + rn_hpd_ini_s * (1._wp - zswitch) ) * tmask(:,:,1)
-            !
-            zsm_i_ini(:,:) = si(jp_smi)%fnow(:,:,1) * tmask(:,:,1)
-            ztm_i_ini(:,:) = si(jp_tmi)%fnow(:,:,1) * tmask(:,:,1)
-            zt_su_ini(:,:) = si(jp_tsu)%fnow(:,:,1) * tmask(:,:,1)
-            ztm_s_ini(:,:) = si(jp_tms)%fnow(:,:,1) * tmask(:,:,1)
-            zapnd_ini(:,:) = si(jp_apd)%fnow(:,:,1) * tmask(:,:,1)
-            zhpnd_ini(:,:) = si(jp_hpd)%fnow(:,:,1) * tmask(:,:,1)
-            !
-            ! change the switch for the following
-            WHERE( zat_i_ini(:,:) > 0._wp )   ;   zswitch(:,:) = tmask(:,:,1) 
-            ELSEWHERE                         ;   zswitch(:,:) = 0._wp
+         
+#if  defined key_agrif
+         ELSE
+ 
+            Agrif_SpecialValue    = -9999.
+            Agrif_UseSpecialValue = .TRUE.
+            CALL Agrif_init_variable(tra_iceini_id,procname=interp_tra_ice)
+            use_sign_north = .TRUE.
+            sign_north = -1.
+            CALL Agrif_init_variable(u_iceini_id  ,procname=interp_u_ice)
+            CALL Agrif_init_variable(v_iceini_id  ,procname=interp_v_ice)
+            Agrif_SpecialValue    = 0._wp
+            use_sign_north = .FALSE.
+            Agrif_UseSpecialValue = .FALSE.
+        ! lbc ???? 
+   ! Here we know : a_i, v_i, v_s, sv_i, oa_i, a_ip, v_ip, t_su, e_s, e_i
+            CALL ice_var_glo2eqv
+            CALL ice_var_zapsmall
+            CALL ice_var_agg(2)
+
+            ! Melt ponds
+            WHERE( a_i > epsi10 )
+               a_ip_frac(:,:,:) = a_ip(:,:,:) / a_i(:,:,:)
+            ELSEWHERE
+               a_ip_frac(:,:,:) = 0._wp
             END WHERE
-            !                          !---------------!
-         ELSE                          ! Read namelist !
-            !                          !---------------!
-            ! no ice if (sst - Tfreez) >= thresold
-            WHERE( ( sst_m(:,:) - (t_bo(:,:) - rt0) ) * tmask(:,:,1) >= rn_thres_sst )   ;   zswitch(:,:) = 0._wp 
-            ELSEWHERE                                                                    ;   zswitch(:,:) = tmask(:,:,1)
-            END WHERE
-            !
-            ! assign initial thickness, concentration, snow depth and salinity to an hemisphere-dependent array
-            WHERE( ff_t(:,:) >= 0._wp )
-               zht_i_ini(:,:) = rn_hti_ini_n * zswitch(:,:)
-               zht_s_ini(:,:) = rn_hts_ini_n * zswitch(:,:)
-               zat_i_ini(:,:) = rn_ati_ini_n * zswitch(:,:)
-               zsm_i_ini(:,:) = rn_smi_ini_n * zswitch(:,:)
-               ztm_i_ini(:,:) = rn_tmi_ini_n * zswitch(:,:)
-               zt_su_ini(:,:) = rn_tsu_ini_n * zswitch(:,:)
-               ztm_s_ini(:,:) = rn_tms_ini_n * zswitch(:,:)
-               zapnd_ini(:,:) = rn_apd_ini_n * zswitch(:,:) * zat_i_ini(:,:) ! rn_apd = pond fraction => rn_apd * a_i = pond conc. 
-               zhpnd_ini(:,:) = rn_hpd_ini_n * zswitch(:,:)
+            WHERE( a_ip > 0._wp )       ! ???????    
+               h_ip(:,:,:) = v_ip(:,:,:) / a_ip(:,:,:)
             ELSEWHERE
-               zht_i_ini(:,:) = rn_hti_ini_s * zswitch(:,:)
-               zht_s_ini(:,:) = rn_hts_ini_s * zswitch(:,:)
-               zat_i_ini(:,:) = rn_ati_ini_s * zswitch(:,:)
-               zsm_i_ini(:,:) = rn_smi_ini_s * zswitch(:,:)
-               ztm_i_ini(:,:) = rn_tmi_ini_s * zswitch(:,:)
-               zt_su_ini(:,:) = rn_tsu_ini_s * zswitch(:,:)
-               ztm_s_ini(:,:) = rn_tms_ini_s * zswitch(:,:)
-               zapnd_ini(:,:) = rn_apd_ini_s * zswitch(:,:) * zat_i_ini(:,:) ! rn_apd = pond fraction => rn_apd * a_i = pond conc.
-               zhpnd_ini(:,:) = rn_hpd_ini_s * zswitch(:,:)
-            END WHERE
-            !
-         ENDIF
-
-         ! make sure ponds = 0 if no ponds scheme
-         IF ( .NOT.ln_pnd ) THEN
-            zapnd_ini(:,:) = 0._wp
-            zhpnd_ini(:,:) = 0._wp
-         ENDIF
-         
-         !-------------!
-         ! fill fields !
-         !-------------!
-         ! select ice covered grid points
-         npti = 0 ; nptidx(:) = 0
-         DO_2D_11_11
-            IF ( zht_i_ini(ji,jj) > 0._wp ) THEN
-               npti         = npti  + 1
-               nptidx(npti) = (jj - 1) * jpi + ji
-            ENDIF
-         END_2D
-
-         ! move to 1D arrays: (jpi,jpj) -> (jpi*jpj)
-         CALL tab_2d_1d( npti, nptidx(1:npti), h_i_1d (1:npti)  , zht_i_ini )
-         CALL tab_2d_1d( npti, nptidx(1:npti), h_s_1d (1:npti)  , zht_s_ini )
-         CALL tab_2d_1d( npti, nptidx(1:npti), at_i_1d(1:npti)  , zat_i_ini )
-         CALL tab_2d_1d( npti, nptidx(1:npti), t_i_1d (1:npti,1), ztm_i_ini )
-         CALL tab_2d_1d( npti, nptidx(1:npti), t_s_1d (1:npti,1), ztm_s_ini )
-         CALL tab_2d_1d( npti, nptidx(1:npti), t_su_1d(1:npti)  , zt_su_ini )
-         CALL tab_2d_1d( npti, nptidx(1:npti), s_i_1d (1:npti)  , zsm_i_ini )
-         CALL tab_2d_1d( npti, nptidx(1:npti), a_ip_1d(1:npti)  , zapnd_ini )
-         CALL tab_2d_1d( npti, nptidx(1:npti), h_ip_1d(1:npti)  , zhpnd_ini )
-
-         ! allocate temporary arrays
-         ALLOCATE( zhi_2d(npti,jpl), zhs_2d(npti,jpl), zai_2d (npti,jpl), &
-            &      zti_2d(npti,jpl), zts_2d(npti,jpl), ztsu_2d(npti,jpl), zsi_2d(npti,jpl), zaip_2d(npti,jpl), zhip_2d(npti,jpl) )
-         
-         ! distribute 1-cat into jpl-cat: (jpi*jpj) -> (jpi*jpj,jpl)
-         CALL ice_var_itd( h_i_1d(1:npti)  , h_s_1d(1:npti)  , at_i_1d(1:npti),                                                   &
-            &              zhi_2d          , zhs_2d          , zai_2d         ,                                                   &
-            &              t_i_1d(1:npti,1), t_s_1d(1:npti,1), t_su_1d(1:npti), s_i_1d(1:npti), a_ip_1d(1:npti), h_ip_1d(1:npti), &
-            &              zti_2d          , zts_2d          , ztsu_2d        , zsi_2d        , zaip_2d        , zhip_2d )
-
-         ! move to 3D arrays: (jpi*jpj,jpl) -> (jpi,jpj,jpl)
-         DO jl = 1, jpl
-            zti_3d(:,:,jl) = rt0 * tmask(:,:,1)
-            zts_3d(:,:,jl) = rt0 * tmask(:,:,1)
-         END DO
-         CALL tab_2d_3d( npti, nptidx(1:npti), zhi_2d   , h_i    )
-         CALL tab_2d_3d( npti, nptidx(1:npti), zhs_2d   , h_s    )
-         CALL tab_2d_3d( npti, nptidx(1:npti), zai_2d   , a_i    )
-         CALL tab_2d_3d( npti, nptidx(1:npti), zti_2d   , zti_3d )
-         CALL tab_2d_3d( npti, nptidx(1:npti), zts_2d   , zts_3d )
-         CALL tab_2d_3d( npti, nptidx(1:npti), ztsu_2d  , t_su   )
-         CALL tab_2d_3d( npti, nptidx(1:npti), zsi_2d   , s_i    )
-         CALL tab_2d_3d( npti, nptidx(1:npti), zaip_2d  , a_ip   )
-         CALL tab_2d_3d( npti, nptidx(1:npti), zhip_2d  , h_ip   )
-
-         ! deallocate temporary arrays
-         DEALLOCATE( zhi_2d, zhs_2d, zai_2d , &
-            &        zti_2d, zts_2d, ztsu_2d, zsi_2d, zaip_2d, zhip_2d )
-
-         ! calculate extensive and intensive variables
-         CALL ice_var_salprof ! for sz_i
-         DO jl = 1, jpl
-            DO_2D_11_11
-               v_i (ji,jj,jl) = h_i(ji,jj,jl) * a_i(ji,jj,jl)
-               v_s (ji,jj,jl) = h_s(ji,jj,jl) * a_i(ji,jj,jl)
-               sv_i(ji,jj,jl) = MIN( MAX( rn_simin , s_i(ji,jj,jl) ) , rn_simax ) * v_i(ji,jj,jl)
-            END_2D
-         END DO
-         !
-         DO jl = 1, jpl
-            DO_3D_11_11( 1, nlay_s )
-               t_s(ji,jj,jk,jl) = zts_3d(ji,jj,jl)
-               e_s(ji,jj,jk,jl) = zswitch(ji,jj) * v_s(ji,jj,jl) * r1_nlay_s * &
-                  &               rhos * ( rcpi * ( rt0 - t_s(ji,jj,jk,jl) ) + rLfus )
-            END_3D
-         END DO
-         !
-         DO jl = 1, jpl
-            DO_3D_11_11( 1, nlay_i )
-               t_i (ji,jj,jk,jl) = zti_3d(ji,jj,jl) 
-               ztmelts          = - rTmlt * sz_i(ji,jj,jk,jl) + rt0 ! melting temperature in K
-               e_i(ji,jj,jk,jl) = zswitch(ji,jj) * v_i(ji,jj,jl) * r1_nlay_i * &
-                  &               rhoi * (  rcpi  * ( ztmelts - t_i(ji,jj,jk,jl) ) + &
-                  &                         rLfus * ( 1._wp - (ztmelts-rt0) / MIN( (t_i(ji,jj,jk,jl)-rt0), -epsi20 ) ) &
-                  &                       - rcp   * ( ztmelts - rt0 ) )
-            END_3D
-         END DO
-
-         ! Melt ponds
-         WHERE( a_i > epsi10 )
-            a_ip_frac(:,:,:) = a_ip(:,:,:) / a_i(:,:,:)
-         ELSEWHERE
-            a_ip_frac(:,:,:) = 0._wp
-         END WHERE
-         v_ip(:,:,:) = h_ip(:,:,:) * a_ip(:,:,:)
-          
-         ! specific temperatures for coupled runs
-         tn_ice(:,:,:) = t_su(:,:,:)
-         t1_ice(:,:,:) = t_i (:,:,1,:)
-         !
+               h_ip(:,:,:) = 0._wp
+            END WHERE   
+
+            tn_ice(:,:,:) = t_su(:,:,:)
+            t1_ice(:,:,:) = t_i (:,:,1,:)
+#endif
+          ENDIF ! Agrif_Root
       ENDIF ! ln_iceini
       !
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/ICE/icestp.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/ICE/icestp.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/ICE/icestp.F90	(revision 13026)
@@ -240,4 +240,8 @@
       CALL par_init                ! set some ice run parameters
       !
+#if defined key_agrif
+      CALL Agrif_Declare_Var_ice  !  "      "   "   "      "  Sea ice
+#endif
+      !
       !                                ! Allocate the ice arrays (sbc_ice already allocated in sbc_init)
       ierr =        ice_alloc        ()      ! ice variables
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_ice.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_ice.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_ice.F90	(revision 13026)
@@ -16,4 +16,5 @@
 
    INTEGER, PUBLIC ::  u_ice_id, v_ice_id, tra_ice_id
+   INTEGER, PUBLIC ::  u_iceini_id, v_iceini_id, tra_iceini_id
    INTEGER, PUBLIC ::  nbstep_ice = 0    ! child time position in sea-ice model
 
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_ice_interp.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_ice_interp.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_ice_interp.F90	(revision 13026)
@@ -14,7 +14,7 @@
    !!----------------------------------------------------------------------
    !!  agrif_interp_ice    : interpolation of ice at "after" sea-ice time step
-   !!  agrif_interp_u_ice   : atomic routine to interpolate u_ice 
-   !!  agrif_interp_v_ice   : atomic routine to interpolate v_ice 
-   !!  agrif_interp_tra_ice : atomic routine to interpolate ice properties 
+   !!  interp_u_ice   : atomic routine to interpolate u_ice 
+   !!  interp_v_ice   : atomic routine to interpolate v_ice 
+   !!  interp_tra_ice : atomic routine to interpolate ice properties 
    !!----------------------------------------------------------------------
    USE par_oce
@@ -23,4 +23,5 @@
    USE ice
    USE agrif_ice
+   USE agrif_oce
    USE phycst , ONLY: rt0
    
@@ -29,4 +30,5 @@
 
    PUBLIC   agrif_interp_ice   ! called by agrif_user.F90
+   PUBLIC   interp_tra_ice, interp_u_ice, interp_v_ice  ! called by iceistate.F90
 
    !!----------------------------------------------------------------------
@@ -68,4 +70,9 @@
       Agrif_SpecialValue    = -9999.
       Agrif_UseSpecialValue = .TRUE.
+
+      use_sign_north = .TRUE.
+      sign_north = -1.
+      if (cd_type == 'T') use_sign_north = .FALSE.
+
       SELECT CASE( cd_type )
       CASE('U')   ;   CALL Agrif_Bc_variable( u_ice_id  , procname=interp_u_ice  , calledweight=zbeta )
@@ -75,4 +82,6 @@
       Agrif_SpecialValue    = 0._wp
       Agrif_UseSpecialValue = .FALSE.
+      
+      use_sign_north = .FALSE.
       !
    END SUBROUTINE agrif_interp_ice
@@ -156,5 +165,5 @@
       ! and it is ok since we conserve tracers (same as in the ocean).
       ALLOCATE( ztab(SIZE(a_i,1),SIZE(a_i,2),SIZE(ptab,3)) )
-     
+
       IF( before ) THEN  ! parent grid
          jm = 1
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_ice_update.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_ice_update.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_ice_update.F90	(revision 13026)
@@ -1,2 +1,3 @@
+#define TWO_WAY
 #undef DECAL_FEEDBACK  /* SEPARATION of INTERFACES*/
 
@@ -66,4 +67,7 @@
       CALL Agrif_Update_Variable( tra_ice_id , locupdate=(/1,0/), procname = update_tra_ice  )
 #endif
+      use_sign_north = .TRUE.
+      sign_north = -1.
+
 # if ! defined DECAL_FEEDBACK
       CALL Agrif_Update_Variable( u_ice_id   , procname = update_u_ice    )
@@ -73,4 +77,5 @@
       CALL Agrif_Update_Variable( v_ice_id   , locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname=update_v_ice)
 #endif
+      use_sign_north = .FALSE.
 !      CALL Agrif_Update_Variable( tra_ice_id , locupdate=(/0,2/), procname = update_tra_ice  )
 !      CALL Agrif_Update_Variable( u_ice_id   , locupdate=(/0,1/), procname = update_u_ice    )
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_oce.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_oce.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_oce.F90	(revision 13026)
@@ -29,4 +29,6 @@
    !
    INTEGER , PUBLIC, PARAMETER ::   nn_sponge_len = 2  !: Sponge width (in number of parent grid points)
+   LOGICAL , PUBLIC ::   ln_bry_south  = .TRUE. !: Is the South boundary open ?
+
    LOGICAL , PUBLIC :: spongedoneT = .FALSE.       !: tracer   sponge layer indicator
    LOGICAL , PUBLIC :: spongedoneU = .FALSE.       !: dynamics sponge layer indicator
@@ -49,8 +51,6 @@
    INTEGER , PUBLIC,              SAVE                 ::  Kbb_a, Kmm_a, Krhs_a   !: AGRIF module-specific copies of time-level indices
 
-# if defined key_vertical
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ht0_parent, hu0_parent, hv0_parent
    INTEGER,  PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: mbkt_parent, mbku_parent, mbkv_parent
-# endif
 
    INTEGER, PUBLIC :: tsn_id                                                  ! AGRIF profile for tracers interpolation and update
@@ -58,4 +58,5 @@
    INTEGER, PUBLIC :: un_update_id, vn_update_id                              ! AGRIF profiles for udpates
    INTEGER, PUBLIC :: tsn_sponge_id, un_sponge_id, vn_sponge_id               ! AGRIF profiles for sponge layers
+   INTEGER, PUBLIC :: tsini_id, uini_id, vini_id, sshini_id                   ! AGRIF profile for initialization
 # if defined key_top
    INTEGER, PUBLIC :: trn_id, trn_sponge_id
@@ -68,4 +69,16 @@
    INTEGER, PUBLIC :: mbkt_id, ht0_id
    INTEGER, PUBLIC :: kindic_agr
+
+   ! North fold
+!$AGRIF_DO_NOT_TREAT
+   LOGICAL, PUBLIC :: use_sign_north
+   REAL, PUBLIC :: sign_north
+   LOGICAL, PUBLIC :: l_ini_child = .FALSE.
+# if defined key_vertical
+   LOGICAL, PUBLIC :: l_vremap    = .TRUE.
+# else
+   LOGICAL, PUBLIC :: l_vremap    = .FALSE.
+# endif
+!$AGRIF_END_DO_NOT_TREAT
    
    !!----------------------------------------------------------------------
@@ -91,9 +104,7 @@
          &      tabspongedone_trn(jpi,jpj),           &
 # endif   
-# if defined key_vertical
          &      ht0_parent(jpi,jpj), mbkt_parent(jpi,jpj),  &
          &      hu0_parent(jpi,jpj), mbku_parent(jpi,jpj),  &
          &      hv0_parent(jpi,jpj), mbkv_parent(jpi,jpj),  &
-# endif      
          &      tabspongedone_u  (jpi,jpj),           &
          &      tabspongedone_v  (jpi,jpj), STAT = ierr(1) )
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_oce_interp.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_oce_interp.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_oce_interp.F90	(revision 13026)
@@ -34,4 +34,5 @@
    USE lib_mpp
    USE vremap
+   USE lbclnk
  
    IMPLICIT NONE
@@ -44,7 +45,7 @@
    PUBLIC   interpunb, interpvnb , interpub2b, interpvb2b
    PUBLIC   interpe3t
-#if defined key_vertical
    PUBLIC   interpht0, interpmbkt
-# endif
+   PUBLIC   agrif_initts, agrif_initssh
+
    INTEGER ::   bdy_tinterp = 0
 
@@ -89,251 +90,262 @@
       Agrif_UseSpecialValue = ln_spc_dyn
       !
+      use_sign_north = .TRUE.
+      sign_north = -1.
       CALL Agrif_Bc_variable( un_interp_id, procname=interpun )
       CALL Agrif_Bc_variable( vn_interp_id, procname=interpvn )
+      use_sign_north = .FALSE.
       !
       Agrif_UseSpecialValue = .FALSE.
       !
       ! --- West --- !
-      ibdy1 = 2
-      ibdy2 = 1+nbghostcells 
-      !
-      IF( .NOT.ln_dynspg_ts ) THEN  ! Store transport
+      IF( lk_west ) THEN
+         ibdy1 = 2
+         ibdy2 = 1+nbghostcells 
+         !
+         IF( .NOT.ln_dynspg_ts ) THEN  ! Store transport
+            DO ji = mi0(ibdy1), mi1(ibdy2)
+               uu_b(ji,:,Krhs_a) = 0._wp
+
+               DO jk = 1, jpkm1
+                  DO jj = 1, jpj
+                     uu_b(ji,jj,Krhs_a) = uu_b(ji,jj,Krhs_a) + e3u(ji,jj,jk,Krhs_a) * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk)
+                  END DO
+               END DO
+
+               DO jj = 1, jpj
+                  uu_b(ji,jj,Krhs_a) = uu_b(ji,jj,Krhs_a) * r1_hu(ji,jj,Krhs_a)
+               END DO
+            END DO
+         ENDIF
+         !
          DO ji = mi0(ibdy1), mi1(ibdy2)
-            uu_b(ji,:,Krhs_a) = 0._wp
-
+            zub(ji,:) = 0._wp    ! Correct transport
             DO jk = 1, jpkm1
                DO jj = 1, jpj
-                  uu_b(ji,jj,Krhs_a) = uu_b(ji,jj,Krhs_a) + e3u(ji,jj,jk,Krhs_a) * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk)
-               END DO
-            END DO
-
-            DO jj = 1, jpj
-               uu_b(ji,jj,Krhs_a) = uu_b(ji,jj,Krhs_a) * r1_hu(ji,jj,Krhs_a)
-            END DO
-         END DO
-      ENDIF
-      !
-      DO ji = mi0(ibdy1), mi1(ibdy2)
-         zub(ji,:) = 0._wp    ! Correct transport
-         DO jk = 1, jpkm1
-            DO jj = 1, jpj
-               zub(ji,jj) = zub(ji,jj) & 
-                  & + e3u(ji,jj,jk,Krhs_a)  * uu(ji,jj,jk,Krhs_a)*umask(ji,jj,jk)
-            END DO
-         END DO
-         DO jj=1,jpj
-            zub(ji,jj) = zub(ji,jj) * r1_hu(ji,jj,Krhs_a)
-         END DO
-            
-         DO jk = 1, jpkm1
-            DO jj = 1, jpj
-               uu(ji,jj,jk,Krhs_a) = ( uu(ji,jj,jk,Krhs_a) + uu_b(ji,jj,Krhs_a)-zub(ji,jj)) * umask(ji,jj,jk)
-            END DO
-         END DO
-      END DO
-            
-      IF( ln_dynspg_ts ) THEN       ! Set tangential velocities to time splitting estimate
-         DO ji = mi0(ibdy1), mi1(ibdy2)
-            zvb(ji,:) = 0._wp
+                  zub(ji,jj) = zub(ji,jj) & 
+                     & + e3u(ji,jj,jk,Krhs_a)  * uu(ji,jj,jk,Krhs_a)*umask(ji,jj,jk)
+               END DO
+            END DO
+            DO jj=1,jpj
+               zub(ji,jj) = zub(ji,jj) * r1_hu(ji,jj,Krhs_a)
+            END DO
+               
             DO jk = 1, jpkm1
                DO jj = 1, jpj
-                  zvb(ji,jj) = zvb(ji,jj) + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk)
-               END DO
-            END DO
-            DO jj = 1, jpj
-               zvb(ji,jj) = zvb(ji,jj) * r1_hv(ji,jj,Krhs_a)
-            END DO
+                  uu(ji,jj,jk,Krhs_a) = ( uu(ji,jj,jk,Krhs_a) + uu_b(ji,jj,Krhs_a)-zub(ji,jj)) * umask(ji,jj,jk)
+               END DO
+            END DO
+         END DO
+               
+         IF( ln_dynspg_ts ) THEN       ! Set tangential velocities to time splitting estimate
+            DO ji = mi0(ibdy1), mi1(ibdy2)
+               zvb(ji,:) = 0._wp
+               DO jk = 1, jpkm1
+                  DO jj = 1, jpj
+                     zvb(ji,jj) = zvb(ji,jj) + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk)
+                  END DO
+               END DO
+               DO jj = 1, jpj
+                  zvb(ji,jj) = zvb(ji,jj) * r1_hv(ji,jj,Krhs_a)
+               END DO
+               DO jk = 1, jpkm1
+                  DO jj = 1, jpj
+                     vv(ji,jj,jk,Krhs_a) = ( vv(ji,jj,jk,Krhs_a) + vv_b(ji,jj,Krhs_a)-zvb(ji,jj))*vmask(ji,jj,jk)
+                  END DO
+               END DO
+            END DO
+         ENDIF
+      ENDIF
+
+      ! --- East --- !
+      IF( lk_east) THEN
+         ibdy1 = jpiglo-1-nbghostcells
+         ibdy2 = jpiglo-2 
+         !
+         IF( .NOT.ln_dynspg_ts ) THEN  ! Store transport
+            DO ji = mi0(ibdy1), mi1(ibdy2)
+               uu_b(ji,:,Krhs_a) = 0._wp
+               DO jk = 1, jpkm1
+                  DO jj = 1, jpj
+                     uu_b(ji,jj,Krhs_a) = uu_b(ji,jj,Krhs_a) & 
+                         & + e3u(ji,jj,jk,Krhs_a) * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk)
+                  END DO
+               END DO
+               DO jj = 1, jpj
+                  uu_b(ji,jj,Krhs_a) = uu_b(ji,jj,Krhs_a) * r1_hu(ji,jj,Krhs_a)
+               END DO
+            END DO
+         ENDIF
+         !
+         DO ji = mi0(ibdy1), mi1(ibdy2)
+            zub(ji,:) = 0._wp    ! Correct transport
             DO jk = 1, jpkm1
                DO jj = 1, jpj
-                  vv(ji,jj,jk,Krhs_a) = ( vv(ji,jj,jk,Krhs_a) + vv_b(ji,jj,Krhs_a)-zvb(ji,jj))*vmask(ji,jj,jk)
-               END DO
-            END DO
-         END DO
-      ENDIF
-
-      ! --- East --- !
-      ibdy1 = jpiglo-1-nbghostcells
-      ibdy2 = jpiglo-2 
-      !
-      IF( .NOT.ln_dynspg_ts ) THEN  ! Store transport
-         DO ji = mi0(ibdy1), mi1(ibdy2)
-            uu_b(ji,:,Krhs_a) = 0._wp
+                  zub(ji,jj) = zub(ji,jj) & 
+                     & + e3u(ji,jj,jk,Krhs_a)  * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk)
+               END DO
+            END DO
+            DO jj=1,jpj
+               zub(ji,jj) = zub(ji,jj) * r1_hu(ji,jj,Krhs_a)
+            END DO
+               
             DO jk = 1, jpkm1
                DO jj = 1, jpj
-                  uu_b(ji,jj,Krhs_a) = uu_b(ji,jj,Krhs_a) & 
-                      & + e3u(ji,jj,jk,Krhs_a) * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk)
-               END DO
-            END DO
-            DO jj = 1, jpj
-               uu_b(ji,jj,Krhs_a) = uu_b(ji,jj,Krhs_a) * r1_hu(ji,jj,Krhs_a)
-            END DO
-         END DO
-      ENDIF
-      !
-      DO ji = mi0(ibdy1), mi1(ibdy2)
-         zub(ji,:) = 0._wp    ! Correct transport
-         DO jk = 1, jpkm1
-            DO jj = 1, jpj
-               zub(ji,jj) = zub(ji,jj) & 
-                  & + e3u(ji,jj,jk,Krhs_a)  * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk)
-            END DO
-         END DO
-         DO jj=1,jpj
-            zub(ji,jj) = zub(ji,jj) * r1_hu(ji,jj,Krhs_a)
-         END DO
-            
-         DO jk = 1, jpkm1
-            DO jj = 1, jpj
-               uu(ji,jj,jk,Krhs_a) = ( uu(ji,jj,jk,Krhs_a) & 
-                 & + uu_b(ji,jj,Krhs_a)-zub(ji,jj))*umask(ji,jj,jk)
-            END DO
-         END DO
-      END DO
-            
-      IF( ln_dynspg_ts ) THEN       ! Set tangential velocities to time splitting estimate
-         ibdy1 = jpiglo-nbghostcells
-         ibdy2 = jpiglo-1 
-         DO ji = mi0(ibdy1), mi1(ibdy2)
-            zvb(ji,:) = 0._wp
-            DO jk = 1, jpkm1
+                  uu(ji,jj,jk,Krhs_a) = ( uu(ji,jj,jk,Krhs_a) & 
+                    & + uu_b(ji,jj,Krhs_a)-zub(ji,jj))*umask(ji,jj,jk)
+               END DO
+            END DO
+         END DO
+               
+         IF( ln_dynspg_ts ) THEN       ! Set tangential velocities to time splitting estimate
+            ibdy1 = jpiglo-nbghostcells
+            ibdy2 = jpiglo-1 
+            DO ji = mi0(ibdy1), mi1(ibdy2)
+               zvb(ji,:) = 0._wp
+               DO jk = 1, jpkm1
+                  DO jj = 1, jpj
+                     zvb(ji,jj) = zvb(ji,jj) &
+                        & + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk)
+                  END DO
+               END DO
                DO jj = 1, jpj
-                  zvb(ji,jj) = zvb(ji,jj) &
+                  zvb(ji,jj) = zvb(ji,jj) * r1_hv(ji,jj,Krhs_a)
+               END DO
+               DO jk = 1, jpkm1
+                  DO jj = 1, jpj
+                     vv(ji,jj,jk,Krhs_a) = ( vv(ji,jj,jk,Krhs_a) & 
+                         & + vv_b(ji,jj,Krhs_a)-zvb(ji,jj)) * vmask(ji,jj,jk)
+                  END DO
+               END DO
+            END DO
+         ENDIF
+      ENDIF
+
+      ! --- South --- !
+      IF( lk_south ) THEN
+         jbdy1 = 2
+         jbdy2 = 1+nbghostcells 
+         !
+         IF( .NOT.ln_dynspg_ts ) THEN  ! Store transport
+            DO jj = mj0(jbdy1), mj1(jbdy2)
+               vv_b(:,jj,Krhs_a) = 0._wp
+               DO jk = 1, jpkm1
+                  DO ji = 1, jpi
+                     vv_b(ji,jj,Krhs_a) = vv_b(ji,jj,Krhs_a) & 
+                         & + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk)
+                  END DO
+               END DO
+               DO ji=1,jpi
+                  vv_b(ji,jj,Krhs_a) = vv_b(ji,jj,Krhs_a) * r1_hv(ji,jj,Krhs_a)     
+               END DO
+            END DO
+         ENDIF
+         !
+         DO jj = mj0(jbdy1), mj1(jbdy2)
+            zvb(:,jj) = 0._wp    ! Correct transport
+            DO jk=1,jpkm1
+               DO ji=1,jpi
+                  zvb(ji,jj) = zvb(ji,jj) & 
                      & + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk)
                END DO
             END DO
-            DO jj = 1, jpj
+            DO ji = 1, jpi
                zvb(ji,jj) = zvb(ji,jj) * r1_hv(ji,jj,Krhs_a)
             END DO
-            DO jk = 1, jpkm1
-               DO jj = 1, jpj
-                  vv(ji,jj,jk,Krhs_a) = ( vv(ji,jj,jk,Krhs_a) & 
-                      & + vv_b(ji,jj,Krhs_a)-zvb(ji,jj)) * vmask(ji,jj,jk)
-               END DO
-            END DO
-         END DO
-      ENDIF
-
-      ! --- South --- !
-      jbdy1 = 2
-      jbdy2 = 1+nbghostcells 
-      !
-      IF( .NOT.ln_dynspg_ts ) THEN  ! Store transport
-         DO jj = mj0(jbdy1), mj1(jbdy2)
-            vv_b(:,jj,Krhs_a) = 0._wp
+
             DO jk = 1, jpkm1
                DO ji = 1, jpi
-                  vv_b(ji,jj,Krhs_a) = vv_b(ji,jj,Krhs_a) & 
-                      & + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk)
-               END DO
-            END DO
-            DO ji=1,jpi
-               vv_b(ji,jj,Krhs_a) = vv_b(ji,jj,Krhs_a) * r1_hv(ji,jj,Krhs_a)     
-            END DO
-         END DO
-      ENDIF
-      !
-      DO jj = mj0(jbdy1), mj1(jbdy2)
-         zvb(:,jj) = 0._wp    ! Correct transport
-         DO jk=1,jpkm1
-            DO ji=1,jpi
-               zvb(ji,jj) = zvb(ji,jj) & 
-                  & + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk)
-            END DO
-         END DO
-         DO ji = 1, jpi
-            zvb(ji,jj) = zvb(ji,jj) * r1_hv(ji,jj,Krhs_a)
-         END DO
-
-         DO jk = 1, jpkm1
+                  vv(ji,jj,jk,Krhs_a) = ( vv(ji,jj,jk,Krhs_a) & 
+                    & + vv_b(ji,jj,Krhs_a) - zvb(ji,jj) ) * vmask(ji,jj,jk)
+               END DO
+            END DO
+         END DO
+               
+         IF( ln_dynspg_ts ) THEN       ! Set tangential velocities to time splitting estimate
+            DO jj = mj0(jbdy1), mj1(jbdy2)
+               zub(:,jj) = 0._wp
+               DO jk = 1, jpkm1
+                  DO ji = 1, jpi
+                     zub(ji,jj) = zub(ji,jj) & 
+                        & + e3u(ji,jj,jk,Krhs_a) * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk)
+                  END DO
+               END DO
+               DO ji = 1, jpi
+                  zub(ji,jj) = zub(ji,jj) * r1_hu(ji,jj,Krhs_a)
+               END DO
+                  
+               DO jk = 1, jpkm1
+                  DO ji = 1, jpi
+                     uu(ji,jj,jk,Krhs_a) = ( uu(ji,jj,jk,Krhs_a) & 
+                       & + uu_b(ji,jj,Krhs_a) - zub(ji,jj) ) * umask(ji,jj,jk)
+                  END DO
+               END DO
+            END DO
+         ENDIF
+      ENDIF
+
+      ! --- North --- !
+      IF( lk_north ) THEN
+         jbdy1 = jpjglo-1-nbghostcells
+         jbdy2 = jpjglo-2 
+         !
+         IF( .NOT.ln_dynspg_ts ) THEN  ! Store transport
+            DO jj = mj0(jbdy1), mj1(jbdy2)
+               vv_b(:,jj,Krhs_a) = 0._wp
+               DO jk = 1, jpkm1
+                  DO ji = 1, jpi
+                     vv_b(ji,jj,Krhs_a) = vv_b(ji,jj,Krhs_a) & 
+                         & + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk)
+                  END DO
+               END DO
+               DO ji=1,jpi
+                  vv_b(ji,jj,Krhs_a) = vv_b(ji,jj,Krhs_a) * r1_hv(ji,jj,Krhs_a)
+               END DO
+            END DO
+         ENDIF
+         !
+         DO jj = mj0(jbdy1), mj1(jbdy2)
+            zvb(:,jj) = 0._wp    ! Correct transport
+            DO jk=1,jpkm1
+               DO ji=1,jpi
+                  zvb(ji,jj) = zvb(ji,jj) & 
+                     & + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk)
+               END DO
+            END DO
             DO ji = 1, jpi
-               vv(ji,jj,jk,Krhs_a) = ( vv(ji,jj,jk,Krhs_a) & 
-                 & + vv_b(ji,jj,Krhs_a) - zvb(ji,jj) ) * vmask(ji,jj,jk)
-            END DO
-         END DO
-      END DO
-            
-      IF( ln_dynspg_ts ) THEN       ! Set tangential velocities to time splitting estimate
-         DO jj = mj0(jbdy1), mj1(jbdy2)
-            zub(:,jj) = 0._wp
+               zvb(ji,jj) = zvb(ji,jj) * r1_hv(ji,jj,Krhs_a)
+            END DO
+
             DO jk = 1, jpkm1
                DO ji = 1, jpi
-                  zub(ji,jj) = zub(ji,jj) & 
-                     & + e3u(ji,jj,jk,Krhs_a) * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk)
-               END DO
-            END DO
-            DO ji = 1, jpi
-               zub(ji,jj) = zub(ji,jj) * r1_hu(ji,jj,Krhs_a)
-            END DO
+                  vv(ji,jj,jk,Krhs_a) = ( vv(ji,jj,jk,Krhs_a) & 
+                    & + vv_b(ji,jj,Krhs_a) - zvb(ji,jj) ) * vmask(ji,jj,jk)
+               END DO
+            END DO
+         END DO
                
-            DO jk = 1, jpkm1
+         IF( ln_dynspg_ts ) THEN       ! Set tangential velocities to time splitting estimate
+            jbdy1 = jpjglo-nbghostcells
+            jbdy2 = jpjglo-1
+            DO jj = mj0(jbdy1), mj1(jbdy2)
+               zub(:,jj) = 0._wp
+               DO jk = 1, jpkm1
+                  DO ji = 1, jpi
+                     zub(ji,jj) = zub(ji,jj) & 
+                        & + e3u(ji,jj,jk,Krhs_a) * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk)
+                  END DO
+               END DO
                DO ji = 1, jpi
-                  uu(ji,jj,jk,Krhs_a) = ( uu(ji,jj,jk,Krhs_a) & 
-                    & + uu_b(ji,jj,Krhs_a) - zub(ji,jj) ) * umask(ji,jj,jk)
-               END DO
-            END DO
-         END DO
-      ENDIF
-
-      ! --- North --- !
-      jbdy1 = jpjglo-1-nbghostcells
-      jbdy2 = jpjglo-2 
-      !
-      IF( .NOT.ln_dynspg_ts ) THEN  ! Store transport
-         DO jj = mj0(jbdy1), mj1(jbdy2)
-            vv_b(:,jj,Krhs_a) = 0._wp
-            DO jk = 1, jpkm1
-               DO ji = 1, jpi
-                  vv_b(ji,jj,Krhs_a) = vv_b(ji,jj,Krhs_a) & 
-                      & + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk)
-               END DO
-            END DO
-            DO ji=1,jpi
-               vv_b(ji,jj,Krhs_a) = vv_b(ji,jj,Krhs_a) * r1_hv(ji,jj,Krhs_a)
-            END DO
-         END DO
-      ENDIF
-      !
-      DO jj = mj0(jbdy1), mj1(jbdy2)
-         zvb(:,jj) = 0._wp    ! Correct transport
-         DO jk=1,jpkm1
-            DO ji=1,jpi
-               zvb(ji,jj) = zvb(ji,jj) & 
-                  & + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk)
-            END DO
-         END DO
-         DO ji = 1, jpi
-            zvb(ji,jj) = zvb(ji,jj) * r1_hv(ji,jj,Krhs_a)
-         END DO
-
-         DO jk = 1, jpkm1
-            DO ji = 1, jpi
-               vv(ji,jj,jk,Krhs_a) = ( vv(ji,jj,jk,Krhs_a) & 
-                 & + vv_b(ji,jj,Krhs_a) - zvb(ji,jj) ) * vmask(ji,jj,jk)
-            END DO
-         END DO
-      END DO
-            
-      IF( ln_dynspg_ts ) THEN       ! Set tangential velocities to time splitting estimate
-         jbdy1 = jpjglo-nbghostcells
-         jbdy2 = jpjglo-1
-         DO jj = mj0(jbdy1), mj1(jbdy2)
-            zub(:,jj) = 0._wp
-            DO jk = 1, jpkm1
-               DO ji = 1, jpi
-                  zub(ji,jj) = zub(ji,jj) & 
-                     & + e3u(ji,jj,jk,Krhs_a) * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk)
-               END DO
-            END DO
-            DO ji = 1, jpi
-               zub(ji,jj) = zub(ji,jj) * r1_hu(ji,jj,Krhs_a)
-            END DO
-               
-            DO jk = 1, jpkm1
-               DO ji = 1, jpi
-                  uu(ji,jj,jk,Krhs_a) = ( uu(ji,jj,jk,Krhs_a) & 
-                    & + uu_b(ji,jj,Krhs_a) - zub(ji,jj) ) * umask(ji,jj,jk)
-               END DO
-            END DO
-         END DO
+                  zub(ji,jj) = zub(ji,jj) * r1_hu(ji,jj,Krhs_a)
+               END DO
+                  
+               DO jk = 1, jpkm1
+                  DO ji = 1, jpi
+                     uu(ji,jj,jk,Krhs_a) = ( uu(ji,jj,jk,Krhs_a) & 
+                       & + uu_b(ji,jj,Krhs_a) - zub(ji,jj) ) * umask(ji,jj,jk)
+                  END DO
+               END DO
+            END DO
+         ENDIF
       ENDIF
       !
@@ -354,54 +366,64 @@
       !
       !--- West ---!
-      istart = 2
-      iend   = nbghostcells+1
-      DO ji = mi0(istart), mi1(iend)
-         DO jj=1,jpj
-            va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj)
-            ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj)
-         END DO
-      END DO
+      IF( lk_west ) THEN
+         istart = 2
+         iend   = nbghostcells+1
+         DO ji = mi0(istart), mi1(iend)
+            DO jj=1,jpj
+               va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj)
+               ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj)
+            END DO
+         END DO
+      ENDIF
       !
       !--- East ---!
-      istart = jpiglo-nbghostcells
-      iend   = jpiglo-1
-      DO ji = mi0(istart), mi1(iend)
-         DO jj=1,jpj
-            va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj)
-         END DO
-      END DO
-      istart = jpiglo-nbghostcells-1
-      iend   = jpiglo-2
-      DO ji = mi0(istart), mi1(iend)
-         DO jj=1,jpj
-            ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj)
-         END DO
-      END DO
+      IF( lk_east ) THEN
+         istart = jpiglo-nbghostcells
+         iend   = jpiglo-1
+         DO ji = mi0(istart), mi1(iend)
+
+            DO jj=1,jpj
+               va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj)
+            END DO
+         END DO
+         istart = jpiglo-nbghostcells-1
+         iend   = jpiglo-2
+         DO ji = mi0(istart), mi1(iend)
+            DO jj=1,jpj
+               ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj)
+            END DO
+         END DO
+      ENDIF 
       !
       !--- South ---!
-      jstart = 2
-      jend   = nbghostcells+1
-      DO jj = mj0(jstart), mj1(jend)
-         DO ji=1,jpi
-            ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj)
-            va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj)
-         END DO
-      END DO
+      IF( lk_south ) THEN
+         jstart = 2
+         jend   = nbghostcells+1
+         DO jj = mj0(jstart), mj1(jend)
+
+            DO ji=1,jpi
+               ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj)
+               va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj)
+            END DO
+         END DO
+      ENDIF       
       !
       !--- North ---!
-      jstart = jpjglo-nbghostcells
-      jend   = jpjglo-1
-      DO jj = mj0(jstart), mj1(jend)
-         DO ji=1,jpi
-            ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj)
-         END DO
-      END DO
-      jstart = jpjglo-nbghostcells-1
-      jend   = jpjglo-2
-      DO jj = mj0(jstart), mj1(jend)
-         DO ji=1,jpi
-            va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj)
-         END DO
-      END DO
+      IF( lk_north ) THEN
+         jstart = jpjglo-nbghostcells
+         jend   = jpjglo-1
+         DO jj = mj0(jstart), mj1(jend)
+            DO ji=1,jpi
+               ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj)
+            END DO
+         END DO
+         jstart = jpjglo-nbghostcells-1
+         jend   = jpjglo-2
+         DO jj = mj0(jstart), mj1(jend)
+            DO ji=1,jpi
+               va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj)
+            END DO
+         END DO
+      ENDIF 
       !
    END SUBROUTINE Agrif_dyn_ts
@@ -421,54 +443,62 @@
       !
       !--- West ---!
-      istart = 2
-      iend   = nbghostcells+1
-      DO ji = mi0(istart), mi1(iend)
-         DO jj=1,jpj
-            zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj)
-            zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj)
-         END DO
-      END DO
+      IF( lk_west ) THEN
+         istart = 2
+         iend   = nbghostcells+1
+         DO ji = mi0(istart), mi1(iend)
+            DO jj=1,jpj
+               zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj)
+               zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj)
+            END DO
+         END DO
+      ENDIF
       !
       !--- East ---!
-      istart = jpiglo-nbghostcells
-      iend   = jpiglo-1
-      DO ji = mi0(istart), mi1(iend)
-         DO jj=1,jpj
-            zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj)
-         END DO
-      END DO
-      istart = jpiglo-nbghostcells-1
-      iend   = jpiglo-2
-      DO ji = mi0(istart), mi1(iend)
-         DO jj=1,jpj
-            zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj)
-         END DO
-      END DO
+      IF( lk_east ) THEN
+         istart = jpiglo-nbghostcells
+         iend   = jpiglo-1
+         DO ji = mi0(istart), mi1(iend)
+            DO jj=1,jpj
+               zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj)
+            END DO
+         END DO
+         istart = jpiglo-nbghostcells-1
+         iend   = jpiglo-2
+         DO ji = mi0(istart), mi1(iend)
+            DO jj=1,jpj
+               zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj)
+            END DO
+         END DO
+      ENDIF
       !
       !--- South ---!
-      jstart = 2
-      jend   = nbghostcells+1
-      DO jj = mj0(jstart), mj1(jend)
-         DO ji=1,jpi
-            zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj)
-            zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj)
-         END DO
-      END DO
+      IF( lk_south ) THEN
+         jstart = 2
+         jend   = nbghostcells+1
+         DO jj = mj0(jstart), mj1(jend)
+            DO ji=1,jpi
+               zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj)
+               zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj)
+            END DO
+         END DO
+      ENDIF
       !
       !--- North ---!
-      jstart = jpjglo-nbghostcells
-      jend   = jpjglo-1
-      DO jj = mj0(jstart), mj1(jend)
-         DO ji=1,jpi
-            zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj)
-         END DO
-      END DO
-      jstart = jpjglo-nbghostcells-1
-      jend   = jpjglo-2
-      DO jj = mj0(jstart), mj1(jend)
-         DO ji=1,jpi
-            zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj)
-         END DO
-      END DO
+      IF( lk_north ) THEN
+         jstart = jpjglo-nbghostcells
+         jend   = jpjglo-1
+         DO jj = mj0(jstart), mj1(jend)
+            DO ji=1,jpi
+               zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj)
+            END DO
+         END DO
+         jstart = jpjglo-nbghostcells-1
+         jend   = jpjglo-2
+         DO jj = mj0(jstart), mj1(jend)
+            DO ji=1,jpi
+               zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj)
+            END DO
+         END DO
+      ENDIF
       !
    END SUBROUTINE Agrif_dyn_ts_flux
@@ -494,4 +524,8 @@
       Agrif_SpecialValue = 0._wp
       Agrif_UseSpecialValue = ln_spc_dyn
+
+      use_sign_north = .TRUE.
+      sign_north = -1.
+
       !
       ! Set bdy time interpolation stage to 0 (latter incremented locally do deal with corners)
@@ -518,4 +552,5 @@
       ENDIF
       Agrif_UseSpecialValue = .FALSE.
+      use_sign_north = .FALSE.
       ! 
    END SUBROUTINE Agrif_dta_ts
@@ -542,38 +577,46 @@
       !
       ! --- West --- !
-      istart = 2
-      iend   = 1 + nbghostcells
-      DO ji = mi0(istart), mi1(iend)
-         DO jj = 1, jpj
-            ssh(ji,jj,Krhs_a) = hbdy(ji,jj)
+      IF(lk_west) THEN
+         istart = 2
+         iend   = 1 + nbghostcells
+         DO ji = mi0(istart), mi1(iend)
+            DO jj = 1, jpj
+               ssh(ji,jj,Krhs_a) = hbdy(ji,jj)
+            ENDDO
          ENDDO
-      ENDDO
+      ENDIF
       !
       ! --- East --- !
-      istart = jpiglo - nbghostcells
-      iend   = jpiglo - 1
-      DO ji = mi0(istart), mi1(iend)
-         DO jj = 1, jpj
-            ssh(ji,jj,Krhs_a) = hbdy(ji,jj)
+      IF(lk_east) THEN
+         istart = jpiglo - nbghostcells
+         iend   = jpiglo - 1
+         DO ji = mi0(istart), mi1(iend)
+            DO jj = 1, jpj
+               ssh(ji,jj,Krhs_a) = hbdy(ji,jj)
+            ENDDO
          ENDDO
-      ENDDO
+      ENDIF
       !
       ! --- South --- !
-      jstart = 2
-      jend   = 1 + nbghostcells
-      DO jj = mj0(jstart), mj1(jend)
-         DO ji = 1, jpi
-            ssh(ji,jj,Krhs_a) = hbdy(ji,jj)
+      IF(lk_south) THEN
+         jstart = 2
+         jend   = 1 + nbghostcells
+         DO jj = mj0(jstart), mj1(jend)
+            DO ji = 1, jpi
+               ssh(ji,jj,Krhs_a) = hbdy(ji,jj)
+            ENDDO
          ENDDO
-      ENDDO
+      ENDIF
       !
       ! --- North --- !
-      jstart = jpjglo - nbghostcells
-      jend   = jpjglo - 1
-      DO jj = mj0(jstart), mj1(jend)
-         DO ji = 1, jpi
-            ssh(ji,jj,Krhs_a) = hbdy(ji,jj)
+      IF(lk_north) THEN
+         jstart = jpjglo - nbghostcells
+         jend   = jpjglo - 1
+         DO jj = mj0(jstart), mj1(jend)
+            DO ji = 1, jpi
+               ssh(ji,jj,Krhs_a) = hbdy(ji,jj)
+            ENDDO
          ENDDO
-      ENDDO
+      ENDIF
       !
    END SUBROUTINE Agrif_ssh
@@ -593,38 +636,46 @@
       !
       ! --- West --- !
-      istart = 2
-      iend   = 1+nbghostcells
-      DO ji = mi0(istart), mi1(iend)
-         DO jj = 1, jpj
-            ssha_e(ji,jj) = hbdy(ji,jj)
+      IF(lk_west) THEN
+         istart = 2
+         iend   = 1+nbghostcells
+         DO ji = mi0(istart), mi1(iend)
+            DO jj = 1, jpj
+               ssha_e(ji,jj) = hbdy(ji,jj)
+            ENDDO
          ENDDO
-      ENDDO
+      ENDIF
       !
       ! --- East --- !
-      istart = jpiglo - nbghostcells
-      iend   = jpiglo - 1
-      DO ji = mi0(istart), mi1(iend)
-         DO jj = 1, jpj
-            ssha_e(ji,jj) = hbdy(ji,jj)
+      IF(lk_east) THEN
+         istart = jpiglo - nbghostcells
+         iend   = jpiglo - 1
+         DO ji = mi0(istart), mi1(iend)
+            DO jj = 1, jpj
+               ssha_e(ji,jj) = hbdy(ji,jj)
+            ENDDO
          ENDDO
-      ENDDO
+      ENDIF
       !
       ! --- South --- !
-      jstart = 2
-      jend   = 1+nbghostcells
-      DO jj = mj0(jstart), mj1(jend)
-         DO ji = 1, jpi
-            ssha_e(ji,jj) = hbdy(ji,jj)
+      IF(lk_south) THEN
+         jstart = 2
+         jend   = 1+nbghostcells
+         DO jj = mj0(jstart), mj1(jend)
+            DO ji = 1, jpi
+               ssha_e(ji,jj) = hbdy(ji,jj)
+            ENDDO
          ENDDO
-      ENDDO
+      ENDIF
       !
       ! --- North --- !
-      jstart = jpjglo - nbghostcells
-      jend   = jpjglo - 1
-      DO jj = mj0(jstart), mj1(jend)
-         DO ji = 1, jpi
-            ssha_e(ji,jj) = hbdy(ji,jj)
+      IF(lk_north) THEN
+         jstart = jpjglo - nbghostcells
+         jend   = jpjglo - 1
+         DO jj = mj0(jstart), mj1(jend)
+            DO ji = 1, jpi
+               ssha_e(ji,jj) = hbdy(ji,jj)
+            ENDDO
          ENDDO
-      ENDDO
+      ENDIF
       !
    END SUBROUTINE Agrif_ssh_ts
@@ -662,4 +713,5 @@
       INTEGER  ::   ji, jj, jk, jn  ! dummy loop indices
       INTEGER  ::   N_in, N_out
+      INTEGER  :: item
       ! vertical interpolation:
       REAL(wp) :: zhtot
@@ -669,5 +721,9 @@
       !!----------------------------------------------------------------------
 
-      IF( before ) THEN         
+      IF( before ) THEN
+
+         item = Kmm_a
+         IF( l_ini_child )   Kmm_a = Kbb_a  
+
          DO jn = 1,jpts
             DO jk=k1,k2
@@ -678,71 +734,85 @@
               END DO
            END DO
-        END DO
-
-# if defined key_vertical
-        ! Interpolate thicknesses
-        ! Warning: these are masked, hence extrapolated prior interpolation.
-        DO jk=k1,k2
-           DO jj=j1,j2
-              DO ji=i1,i2
-                  ptab(ji,jj,jk,jpts+1) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kmm_a)
-              END DO
-           END DO
-        END DO
-
-        ! Extrapolate thicknesses in partial bottom cells:
-        ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on
-        IF (ln_zps) THEN
-           DO jj=j1,j2
-              DO ji=i1,i2
-                  jk = mbkt(ji,jj)
-                  ptab(ji,jj,jk,jpts+1) = 0._wp
-              END DO
-           END DO           
-        END IF
-     
-        ! Save ssh at last level:
-        IF (.NOT.ln_linssh) THEN
-           ptab(i1:i2,j1:j2,k2,jpts+1) = ssh(i1:i2,j1:j2,Kmm_a)*tmask(i1:i2,j1:j2,1) 
-        ELSE
-           ptab(i1:i2,j1:j2,k2,jpts+1) = 0._wp
-        END IF      
-# endif
+         END DO
+
+         IF( l_vremap .OR. l_ini_child) THEN
+            ! Interpolate thicknesses
+            ! Warning: these are masked, hence extrapolated prior interpolation.
+            DO jk=k1,k2
+               DO jj=j1,j2
+                  DO ji=i1,i2
+                      ptab(ji,jj,jk,jpts+1) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kmm_a)
+
+                  END DO
+               END DO
+            END DO
+
+            ! Extrapolate thicknesses in partial bottom cells:
+            ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on
+            IF (ln_zps) THEN
+               DO jj=j1,j2
+                  DO ji=i1,i2
+                      jk = mbkt(ji,jj)
+                      ptab(ji,jj,jk,jpts+1) = 0._wp
+                  END DO
+               END DO           
+            END IF
+        
+            ! Save ssh at last level:
+            IF (.NOT.ln_linssh) THEN
+               ptab(i1:i2,j1:j2,k2,jpts+1) = ssh(i1:i2,j1:j2,Kmm_a)*tmask(i1:i2,j1:j2,1) 
+            ELSE
+               ptab(i1:i2,j1:j2,k2,jpts+1) = 0._wp
+            END IF      
+         ENDIF
+         Kmm_a = item
+
       ELSE 
-
-# if defined key_vertical 
-         IF (ln_linssh) ptab(i1:i2,j1:j2,k2,n2) = 0._wp 
-            
-         DO jj=j1,j2
-            DO ji=i1,i2
-               ts(ji,jj,:,:,Krhs_a) = 0._wp
-               N_in = mbkt_parent(ji,jj)
-               zhtot = 0._wp
-               DO jk=1,N_in !k2 = jpk of parent grid
-                  IF (jk==N_in) THEN
-                     h_in(jk) = ht0_parent(ji,jj) + ptab(ji,jj,k2,n2) - zhtot
-                  ELSE
-                     h_in(jk) = ptab(ji,jj,jk,n2)
+         item = Krhs_a
+         IF( l_ini_child )   Krhs_a = Kbb_a  
+
+         IF( l_vremap .OR. l_ini_child ) THEN
+            IF (ln_linssh) ptab(i1:i2,j1:j2,k2,n2) = 0._wp 
+               
+            DO jj=j1,j2
+               DO ji=i1,i2
+                  ts(ji,jj,:,:,Krhs_a) = 0.                  
+               !   IF( l_ini_child) ts(ji,jj,:,:,Krhs_a) = ptab(ji,jj,:,1:jpts)
+                  N_in = mbkt_parent(ji,jj)
+                  zhtot = 0._wp
+                  DO jk=1,N_in !k2 = jpk of parent grid
+                     IF (jk==N_in) THEN
+                        h_in(jk) = ht0_parent(ji,jj) + ptab(ji,jj,k2,n2) - zhtot
+                     ELSE
+                        h_in(jk) = ptab(ji,jj,jk,n2)
+                     ENDIF
+                     zhtot = zhtot + h_in(jk)
+                     tabin(jk,:) = ptab(ji,jj,jk,n1:n2-1)
+                  END DO
+                  N_out = 0
+                  DO jk=1,jpk ! jpk of child grid
+                     IF (tmask(ji,jj,jk) == 0._wp) EXIT 
+                     N_out = N_out + 1
+                     h_out(jk) = e3t(ji,jj,jk,Krhs_a)
+                  ENDDO
+                  IF (N_in*N_out > 0) THEN
+                     IF( l_ini_child ) THEN
+                        CALL remap_linear(tabin(1:N_in,1:jpts),h_in(1:N_in),ts(ji,jj,1:N_out,1:jpts,Krhs_a),          &
+                                      &   h_out(1:N_out),N_in,N_out,jpts)  
+                     ELSE 
+                        CALL reconstructandremap(tabin(1:N_in,1:jpts),h_in(1:N_in),ts(ji,jj,1:N_out,1:jpts,Krhs_a),   &
+                                      &   h_out(1:N_out),N_in,N_out,jpts)  
+                     ENDIF
                   ENDIF
-                  zhtot = zhtot + h_in(jk)
-                  tabin(jk,:) = ptab(ji,jj,jk,n1:n2-1)
-               END DO
-               N_out = 0
-               DO jk=1,jpk ! jpk of child grid
-                  IF (tmask(ji,jj,jk) == 0._wp) EXIT 
-                  N_out = N_out + 1
-                  h_out(jk) = e3t(ji,jj,jk,Krhs_a)
                ENDDO
-               IF (N_in*N_out > 0) THEN
-                  CALL reconstructandremap(tabin(1:N_in,1:jpts),h_in(1:N_in),ts(ji,jj,1:N_out,1:jpts,Krhs_a),h_out(1:N_out),N_in,N_out,jpts)
-               ENDIF
             ENDDO
-         ENDDO
-# else
-         !
-         DO jn=1, jpts
-            ts(i1:i2,j1:j2,1:jpk,jn,Krhs_a)=ptab(i1:i2,j1:j2,1:jpk,jn)*tmask(i1:i2,j1:j2,1:jpk) 
-         END DO
-# endif
+            Krhs_a = item
+ 
+         ELSE
+         
+            DO jn=1, jpts
+                ts(i1:i2,j1:j2,1:jpk,jn,Krhs_a)=ptab(i1:i2,j1:j2,1:jpk,jn)*tmask(i1:i2,j1:j2,1:jpk) 
+            END DO
+         ENDIF
 
       ENDIF
@@ -782,83 +852,99 @@
       REAL(wp), DIMENSION(k1:k2) :: tabin, h_in
       REAL(wp), DIMENSION(1:jpk) :: h_out
-      INTEGER  :: N_in, N_out
+      INTEGER  :: N_in, N_out,item
       REAL(wp) :: h_diff
       !!---------------------------------------------    
       !
       IF (before) THEN 
+
+         item = Kmm_a
+         IF( l_ini_child )   Kmm_a = Kbb_a     
+
          DO jk=1,jpk
             DO jj=j1,j2
                DO ji=i1,i2
                   ptab(ji,jj,jk,1) = (e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a) * uu(ji,jj,jk,Kmm_a)*umask(ji,jj,jk)) 
-# if defined key_vertical
-                  ! Interpolate thicknesses (masked for subsequent extrapolation)
-                  ptab(ji,jj,jk,2) = umask(ji,jj,jk) * e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a)
-# endif
-               END DO
-            END DO
-         END DO
-# if defined key_vertical
+                  IF( l_vremap .OR. l_ini_child) THEN
+                     ! Interpolate thicknesses (masked for subsequent extrapolation)
+                     ptab(ji,jj,jk,2) = umask(ji,jj,jk) * e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a)
+                  ENDIF
+               END DO
+            END DO
+         END DO
+
+        IF( l_vremap .OR. l_ini_child) THEN
          ! Extrapolate thicknesses in partial bottom cells:
          ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on
-         IF (ln_zps) THEN
-            DO jj=j1,j2
-               DO ji=i1,i2
-                  jk = mbku(ji,jj)
-                  ptab(ji,jj,jk,2) = 0._wp
-               END DO
-            END DO           
-         END IF
-        ! Save ssh at last level:
-        ptab(i1:i2,j1:j2,k2,2) = 0._wp
-        IF (.NOT.ln_linssh) THEN
-           ! This vertical sum below should be replaced by the sea-level at U-points (optimization):
-           DO jk=1,jpk
-              ptab(i1:i2,j1:j2,k2,2) = ptab(i1:i2,j1:j2,k2,2) + e3u(i1:i2,j1:j2,jk,Kmm_a) * umask(i1:i2,j1:j2,jk)
-           END DO
-           ptab(i1:i2,j1:j2,k2,2) = ptab(i1:i2,j1:j2,k2,2) - hu_0(i1:i2,j1:j2)
-        END IF 
-# endif
+            IF (ln_zps) THEN
+               DO jj=j1,j2
+                  DO ji=i1,i2
+                     jk = mbku(ji,jj)
+                     ptab(ji,jj,jk,2) = 0._wp
+                  END DO
+               END DO           
+            END IF
+
+           ! Save ssh at last level:
+           ptab(i1:i2,j1:j2,k2,2) = 0._wp
+           IF (.NOT.ln_linssh) THEN
+              ! This vertical sum below should be replaced by the sea-level at U-points (optimization):
+              DO jk=1,jpk
+                 ptab(i1:i2,j1:j2,k2,2) = ptab(i1:i2,j1:j2,k2,2) + e3u(i1:i2,j1:j2,jk,Kmm_a) * umask(i1:i2,j1:j2,jk)
+              END DO
+              ptab(i1:i2,j1:j2,k2,2) = ptab(i1:i2,j1:j2,k2,2) - hu_0(i1:i2,j1:j2)
+           END IF
+        ENDIF
+
+         Kmm_a = item
          !
       ELSE
          zrhoy = Agrif_rhoy()
-# if defined key_vertical
+
+        IF( l_vremap .OR. l_ini_child) THEN
 ! VERTICAL REFINEMENT BEGIN
 
-         IF (ln_linssh) ptab(i1:i2,j1:j2,k2,2) = 0._wp 
-
-         DO ji=i1,i2
-            DO jj=j1,j2
-               uu(ji,jj,:,Krhs_a) = 0._wp
-               N_in = mbku_parent(ji,jj)
-               zhtot = 0._wp
-               DO jk=1,N_in
-                  IF (jk==N_in) THEN
-                     h_in(jk) = hu0_parent(ji,jj) + ptab(ji,jj,k2,2) - zhtot
-                  ELSE
-                     h_in(jk) = ptab(ji,jj,jk,2)/(e2u(ji,jj)*zrhoy) 
-                  ENDIF
-                  zhtot = zhtot + h_in(jk)
-                  tabin(jk) = ptab(ji,jj,jk,1)/(e2u(ji,jj)*zrhoy*h_in(jk))
-              ENDDO
-                  
-              N_out = 0
-              DO jk=1,jpk
-                 if (umask(ji,jj,jk) == 0) EXIT
-                 N_out = N_out + 1
-                 h_out(N_out) = e3u(ji,jj,jk,Krhs_a)
-              ENDDO
-              IF (N_in*N_out > 0) THEN
-                 CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),uu(ji,jj,1:N_out,Krhs_a),h_out(1:N_out),N_in,N_out,1)
-              ENDIF
+            IF (ln_linssh) ptab(i1:i2,j1:j2,k2,2) = 0._wp 
+
+            DO ji=i1,i2
+               DO jj=j1,j2
+                  uu(ji,jj,:,Krhs_a) = 0._wp
+                  N_in = mbku_parent(ji,jj)
+                  zhtot = 0._wp
+                  DO jk=1,N_in
+                     IF (jk==N_in) THEN
+                        h_in(jk) = hu0_parent(ji,jj) + ptab(ji,jj,k2,2) - zhtot
+                     ELSE
+                        h_in(jk) = ptab(ji,jj,jk,2)/(e2u(ji,jj)*zrhoy) 
+                     ENDIF
+                     zhtot = zhtot + h_in(jk)
+                     IF( h_in(jk) .GT. 0. ) THEN
+                     tabin(jk) = ptab(ji,jj,jk,1)/(e2u(ji,jj)*zrhoy*h_in(jk))
+                     ELSE
+                     tabin(jk) = 0.
+                     ENDIF
+                 ENDDO
+                     
+                 N_out = 0
+                 DO jk=1,jpk
+                    IF (umask(ji,jj,jk) == 0) EXIT
+                    N_out = N_out + 1
+                    h_out(N_out) = e3u(ji,jj,jk,Krhs_a)
+                 ENDDO
+                 IF (N_in*N_out > 0) THEN
+                     IF( l_ini_child ) THEN
+                        CALL remap_linear       (tabin(1:N_in),h_in(1:N_in),uu(ji,jj,1:N_out,Krhs_a),h_out(1:N_out),N_in,N_out,1)
+                     ELSE
+                        CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),uu(ji,jj,1:N_out,Krhs_a),h_out(1:N_out),N_in,N_out,1)
+                     ENDIF   
+                 ENDIF
+               ENDDO
             ENDDO
-         ENDDO
-
-# else
-         DO jk = 1, jpkm1
-            DO jj=j1,j2
-               uu(i1:i2,jj,jk,Krhs_a) = ptab(i1:i2,jj,jk,1) / ( zrhoy * e2u(i1:i2,jj) * e3u(i1:i2,jj,jk,Krhs_a) )
-            END DO
-         END DO
-# endif
+         ELSE
+            DO jk = 1, jpkm1
+               DO jj=j1,j2
+                  uu(i1:i2,jj,jk,Krhs_a) = ptab(i1:i2,jj,jk,1) / ( zrhoy * e2u(i1:i2,jj) * e3u(i1:i2,jj,jk,Krhs_a) )
+               END DO
+            END DO
+         ENDIF
 
       ENDIF
@@ -880,78 +966,94 @@
       REAL(wp), DIMENSION(k1:k2) :: tabin, h_in
       REAL(wp), DIMENSION(1:jpk) :: h_out
-      INTEGER  :: N_in, N_out
+      INTEGER  :: N_in, N_out, item
       REAL(wp) :: h_diff, zhtot
       !!---------------------------------------------    
       !      
-      IF (before) THEN          
+      IF (before) THEN   
+
+         item = Kmm_a
+         IF( l_ini_child )   Kmm_a = Kbb_a     
+       
          DO jk=k1,k2
             DO jj=j1,j2
                DO ji=i1,i2
                   ptab(ji,jj,jk,1) = (e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a) * vv(ji,jj,jk,Kmm_a)*vmask(ji,jj,jk))
-# if defined key_vertical
-                  ! Interpolate thicknesses (masked for subsequent extrapolation)
-                  ptab(ji,jj,jk,2) = vmask(ji,jj,jk) * e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a)
-# endif
-               END DO
-            END DO
-         END DO
-# if defined key_vertical
+                  IF( l_vremap .OR. l_ini_child) THEN
+                     ! Interpolate thicknesses (masked for subsequent extrapolation)
+                     ptab(ji,jj,jk,2) = vmask(ji,jj,jk) * e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a)
+                  ENDIF
+               END DO
+            END DO
+         END DO
+
+         IF( l_vremap .OR. l_ini_child) THEN
          ! Extrapolate thicknesses in partial bottom cells:
          ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on
-         IF (ln_zps) THEN
+            IF (ln_zps) THEN
+               DO jj=j1,j2
+                  DO ji=i1,i2
+                     jk = mbkv(ji,jj)
+                     ptab(ji,jj,jk,2) = 0._wp
+                  END DO
+               END DO           
+            END IF
+            ! Save ssh at last level:
+            ptab(i1:i2,j1:j2,k2,2) = 0._wp
+            IF (.NOT.ln_linssh) THEN
+               ! This vertical sum below should be replaced by the sea-level at V-points (optimization):
+               DO jk=1,jpk
+                  ptab(i1:i2,j1:j2,k2,2) = ptab(i1:i2,j1:j2,k2,2) + e3v(i1:i2,j1:j2,jk,Kmm_a) * vmask(i1:i2,j1:j2,jk)
+               END DO
+               ptab(i1:i2,j1:j2,k2,2) = ptab(i1:i2,j1:j2,k2,2) - hv_0(i1:i2,j1:j2)
+            END IF 
+         ENDIF
+         item = Kmm_a
+
+      ELSE       
+         zrhox = Agrif_rhox()
+
+         IF( l_vremap .OR. l_ini_child ) THEN
+
+            IF (ln_linssh) ptab(i1:i2,j1:j2,k2,2) = 0._wp 
+
             DO jj=j1,j2
                DO ji=i1,i2
-                  jk = mbkv(ji,jj)
-                  ptab(ji,jj,jk,2) = 0._wp
-               END DO
-            END DO           
-         END IF
-        ! Save ssh at last level:
-        ptab(i1:i2,j1:j2,k2,2) = 0._wp
-        IF (.NOT.ln_linssh) THEN
-           ! This vertical sum below should be replaced by the sea-level at V-points (optimization):
-           DO jk=1,jpk
-              ptab(i1:i2,j1:j2,k2,2) = ptab(i1:i2,j1:j2,k2,2) + e3v(i1:i2,j1:j2,jk,Kmm_a) * vmask(i1:i2,j1:j2,jk)
-           END DO
-           ptab(i1:i2,j1:j2,k2,2) = ptab(i1:i2,j1:j2,k2,2) - hv_0(i1:i2,j1:j2)
-        END IF 
-# endif
-      ELSE       
-         zrhox = Agrif_rhox()
-# if defined key_vertical
-
-         IF (ln_linssh) ptab(i1:i2,j1:j2,k2,2) = 0._wp 
-
-         DO jj=j1,j2
-            DO ji=i1,i2
-               vv(ji,jj,:,Krhs_a) = 0._wp
-               N_in = mbkv_parent(ji,jj)
-               zhtot = 0._wp
-               DO jk=1,N_in
-                  IF (jk==N_in) THEN
-                     h_in(jk) = hv0_parent(ji,jj) + ptab(ji,jj,k2,2) - zhtot
-                  ELSE
-                     h_in(jk) = ptab(ji,jj,jk,2)/(e1v(ji,jj)*zrhox) 
+                  vv(ji,jj,:,Krhs_a) = 0._wp
+                  N_in = mbkv_parent(ji,jj)
+                  zhtot = 0._wp
+                  DO jk=1,N_in
+                     IF (jk==N_in) THEN
+                        h_in(jk) = hv0_parent(ji,jj) + ptab(ji,jj,k2,2) - zhtot
+                     ELSE
+                        h_in(jk) = ptab(ji,jj,jk,2)/(e1v(ji,jj)*zrhox) 
+                     ENDIF
+                     zhtot = zhtot + h_in(jk)
+                    IF( h_in(jk) .GT. 0. ) THEN
+                     tabin(jk) = ptab(ji,jj,jk,1)/(e1v(ji,jj)*zrhox*h_in(jk))
+                    ELSE
+                     tabin(jk)  = 0.
+                    ENDIF 
+                 ENDDO
+            
+                  N_out = 0
+                  DO jk=1,jpk
+                     if (vmask(ji,jj,jk) == 0) EXIT
+                     N_out = N_out + 1
+                     h_out(N_out) = e3v(ji,jj,jk,Krhs_a)
+                  END DO
+                  IF (N_in*N_out > 0) THEN
+                     IF( l_ini_child ) THEN
+                        CALL remap_linear       (tabin(1:N_in),h_in(1:N_in),vv(ji,jj,1:N_out,Krhs_a),h_out(1:N_out),N_in,N_out,1)
+                     ELSE
+                        CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),vv(ji,jj,1:N_out,Krhs_a),h_out(1:N_out),N_in,N_out,1)
+                     ENDIF   
                   ENDIF
-                  zhtot = zhtot + h_in(jk)
-                  tabin(jk) = ptab(ji,jj,jk,1)/(e1v(ji,jj)*zrhox*h_in(jk))
-              ENDDO
-         
-               N_out = 0
-               DO jk=1,jpk
-                  if (vmask(ji,jj,jk) == 0) EXIT
-                  N_out = N_out + 1
-                  h_out(N_out) = e3v(ji,jj,jk,Krhs_a)
-               END DO
-               IF (N_in*N_out > 0) THEN
-                  call reconstructandremap(tabin(1:N_in),h_in(1:N_in),vv(ji,jj,1:N_out,Krhs_a),h_out(1:N_out),N_in,N_out,1)
-               ENDIF
-            END DO
-         END DO
-# else
-         DO jk = 1, jpkm1
-            vv(i1:i2,j1:j2,jk,Krhs_a) = ptab(i1:i2,j1:j2,jk,1) / ( zrhox * e1v(i1:i2,j1:j2) * e3v(i1:i2,j1:j2,jk,Krhs_a) )
-         END DO
-# endif
+               END DO
+            END DO
+         ELSE
+            DO jk = 1, jpkm1
+               vv(i1:i2,j1:j2,jk,Krhs_a) = ptab(i1:i2,j1:j2,jk,1) / ( zrhox * e1v(i1:i2,j1:j2) * e3v(i1:i2,j1:j2,jk,Krhs_a) )
+            END DO
+         ENDIF
       ENDIF
       !        
@@ -1163,5 +1265,4 @@
    END SUBROUTINE interpe3t
 
-
    SUBROUTINE interpavm( ptab, i1, i2, j1, j2, k1, k2, m1, m2, before )
       !!----------------------------------------------------------------------
@@ -1185,66 +1286,67 @@
               END DO
            END DO
-        END DO
-
-# if defined key_vertical
-        ! Interpolate thicknesses
-        ! Warning: these are masked, hence extrapolated prior interpolation.
-        DO jk=k1,k2
-           DO jj=j1,j2
-              DO ji=i1,i2
-                  ptab(ji,jj,jk,2) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kmm_a)
-              END DO
-           END DO
-        END DO
-
-        ! Extrapolate thicknesses in partial bottom cells:
-        ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on
-        IF (ln_zps) THEN
-           DO jj=j1,j2
-              DO ji=i1,i2
-                  jk = mbkt(ji,jj)
-                  ptab(ji,jj,jk,2) = 0._wp
-              END DO
-           END DO           
-        END IF
-     
-        ! Save ssh at last level:
-        IF (.NOT.ln_linssh) THEN
-           ptab(i1:i2,j1:j2,k2,2) = ssh(i1:i2,j1:j2,Kmm_a)*tmask(i1:i2,j1:j2,1) 
-        ELSE
-           ptab(i1:i2,j1:j2,k2,2) = 0._wp
-        END IF      
-# endif
+         END DO
+
+         IF( l_vremap ) THEN
+            ! Interpolate thicknesses
+            ! Warning: these are masked, hence extrapolated prior interpolation.
+            DO jk=k1,k2
+               DO jj=j1,j2
+                  DO ji=i1,i2
+                      ptab(ji,jj,jk,2) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kmm_a)
+                  END DO
+               END DO
+            END DO
+
+            ! Extrapolate thicknesses in partial bottom cells:
+            ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on
+            IF (ln_zps) THEN
+               DO jj=j1,j2
+                  DO ji=i1,i2
+                      jk = mbkt(ji,jj)
+                      ptab(ji,jj,jk,2) = 0._wp
+                  END DO
+               END DO           
+            END IF
+        
+           ! Save ssh at last level:
+            IF (.NOT.ln_linssh) THEN
+               ptab(i1:i2,j1:j2,k2,2) = ssh(i1:i2,j1:j2,Kmm_a)*tmask(i1:i2,j1:j2,1) 
+            ELSE
+               ptab(i1:i2,j1:j2,k2,2) = 0._wp
+            END IF      
+          ENDIF
+
       ELSE 
-#ifdef key_vertical         
-         IF (ln_linssh) ptab(i1:i2,j1:j2,k2,2) = 0._wp 
-         avm_k(i1:i2,j1:j2,k1:k2) = 0._wp
-            
-         DO jj = j1, j2
-            DO ji =i1, i2
-               N_in = mbkt_parent(ji,jj)
-               IF ( tmask(ji,jj,1) == 0._wp) N_in = 0
-               z_in(N_in+1) = ht0_parent(ji,jj) + ptab(ji,jj,k2,2)
-               DO jk = N_in, 1, -1  ! Parent vertical grid               
-                     z_in(jk) = z_in(jk+1) - ptab(ji,jj,jk,2)
-                    tabin(jk) = ptab(ji,jj,jk,1)
-               END DO
-               N_out = mbkt(ji,jj) 
-               DO jk = 1, N_out        ! Child vertical grid
-                  z_out(jk) = gdepw(ji,jj,jk,Kmm_a)
+
+         IF( l_vremap ) THEN
+            IF (ln_linssh) ptab(i1:i2,j1:j2,k2,2) = 0._wp 
+            avm_k(i1:i2,j1:j2,k1:k2) = 0._wp
+               
+            DO jj = j1, j2
+               DO ji =i1, i2
+                  N_in = mbkt_parent(ji,jj)
+                  IF ( tmask(ji,jj,1) == 0._wp) N_in = 0
+                  z_in(N_in+1) = ht0_parent(ji,jj) + ptab(ji,jj,k2,2)
+                  DO jk = N_in, 1, -1  ! Parent vertical grid               
+                        z_in(jk) = z_in(jk+1) - ptab(ji,jj,jk,2)
+                       tabin(jk) = ptab(ji,jj,jk,1)
+                  END DO
+                  N_out = mbkt(ji,jj) 
+                  DO jk = 1, N_out        ! Child vertical grid
+                     z_out(jk) = gdepw(ji,jj,jk,Kmm_a)
+                  ENDDO
+                  IF (N_in*N_out > 0) THEN
+                     CALL remap_linear(tabin(1:N_in),z_in(1:N_in),avm_k(ji,jj,1:N_out),z_out(1:N_out),N_in,N_out,1)
+                  ENDIF
                ENDDO
-               IF (N_in*N_out > 0) THEN
-                  CALL remap_linear(tabin(1:N_in),z_in(1:N_in),avm_k(ji,jj,1:N_out),z_out(1:N_out),N_in,N_out,1)
-               ENDIF
             ENDDO
-         ENDDO
-#else
-         avm_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2,1)
-#endif
+         ELSE
+            avm_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2,1)
+         ENDIF
       ENDIF
       !
    END SUBROUTINE interpavm
 
-# if defined key_vertical
    SUBROUTINE interpmbkt( ptab, i1, i2, j1, j2, before )
       !!----------------------------------------------------------------------
@@ -1282,6 +1384,41 @@
       !
    END SUBROUTINE interpht0
-#endif
-
+
+   SUBROUTINE agrif_initts(tabres,i1,i2,j1,j2,k1,k2,m1,m2,before)
+       INTEGER :: i1, i2, j1, j2, k1, k2, m1, m2
+       REAL(wp):: tabres(i1:i2,j1:j2,k1:k2,m1:m2)
+       LOGICAL :: before
+
+       INTEGER :: jm
+
+       IF (before) THEN
+         DO jm=1,jpts
+             tabres(i1:i2,j1:j2,k1:k2,jm) = ts(i1:i2,j1:j2,k1:k2,jm,Kbb_a)
+         END DO
+       ELSE
+         DO jm=1,jpts
+             ts(i1:i2,j1:j2,k1:k2,jm,Kbb_a)=tabres(i1:i2,j1:j2,k1:k2,jm)
+         END DO
+       ENDIF
+   END SUBROUTINE agrif_initts 
+
+   SUBROUTINE agrif_initssh( ptab, i1, i2, j1, j2, before )
+      !!----------------------------------------------------------------------
+      !!                  ***  ROUTINE interpsshn  ***
+      !!----------------------------------------------------------------------  
+      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2
+      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   ptab
+      LOGICAL                         , INTENT(in   ) ::   before
+      !
+      !!----------------------------------------------------------------------  
+      !
+      IF( before) THEN
+         ptab(i1:i2,j1:j2) = ssh(i1:i2,j1:j2,Kbb_a)
+      ELSE
+         ssh(i1:i2,j1:j2,Kbb_a) = ptab(i1:i2,j1:j2)*tmask(i1:i2,j1:j2,1)
+      ENDIF
+      !
+   END SUBROUTINE agrif_initssh
+   
 #else
    !!----------------------------------------------------------------------
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_oce_sponge.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_oce_sponge.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_oce_sponge.F90	(revision 13026)
@@ -80,4 +80,6 @@
       Agrif_SpecialValue=0.
       Agrif_UseSpecialValue = ln_spc_dyn
+      use_sign_north = .TRUE.
+      sign_north = -1.
       !
       tabspongedone_u = .FALSE.
@@ -90,4 +92,5 @@
       !
       Agrif_UseSpecialValue = .FALSE.
+      use_sign_north = .FALSE.
 #endif
       !
@@ -127,42 +130,51 @@
 
          ! --- West --- !
-         ztabramp(:,:) = 0._wp
-         ind1 = 1+nbghostcells
-         DO ji = mi0(ind1), mi1(ind1)                
-            ztabramp(ji,:) = ssumask(ji,:)
-         END DO
-         !
-         zmskwest(:) = 0._wp
-         zmskwest(1:jpj) = MAXVAL(ztabramp(:,:), dim=1)
+         IF( lk_west) THEN
+            ztabramp(:,:) = 0._wp
+            ind1 = 1+nbghostcells
+            DO ji = mi0(ind1), mi1(ind1)                
+               ztabramp(ji,:) = ssumask(ji,:)
+            END DO
+            !
+            zmskwest(:) = 0._wp
+            zmskwest(1:jpj) = MAXVAL(ztabramp(:,:), dim=1)
+         ENDIF
 
          ! --- East --- !
-         ztabramp(:,:) = 0._wp
-         ind1 = jpiglo - nbghostcells - 1
-         DO ji = mi0(ind1), mi1(ind1)                 
-            ztabramp(ji,:) = ssumask(ji,:)
-         END DO
-         !
-         zmskeast(:) = 0._wp
-         zmskeast(1:jpj) = MAXVAL(ztabramp(:,:), dim=1)
+         IF( lk_east ) THEN
+            ztabramp(:,:) = 0._wp
+            ind1 = jpiglo - nbghostcells - 1
+            DO ji = mi0(ind1), mi1(ind1)                 
+               ztabramp(ji,:) = ssumask(ji,:)
+            END DO
+            !
+            zmskeast(:) = 0._wp
+            zmskeast(1:jpj) = MAXVAL(ztabramp(:,:), dim=1)
+         ENDIF
 
          ! --- South --- !
-         ztabramp(:,:) = 0._wp
-         ind1 = 1+nbghostcells
-         DO jj = mj0(ind1), mj1(ind1)                 
-            ztabramp(:,jj) = ssvmask(:,jj)
-         END DO
-         !
-         zmsksouth(:) = 0._wp
-         zmsksouth(1:jpi) = MAXVAL(ztabramp(:,:), dim=2)
+         IF( lk_south ) THEN
+            ztabramp(:,:) = 0._wp
+            ind1 = 1+nbghostcells
+            DO jj = mj0(ind1), mj1(ind1)                 
+               ztabramp(:,jj) = ssvmask(:,jj)
+            END DO
+            !
+            zmsksouth(:) = 0._wp
+            zmsksouth(1:jpi) = MAXVAL(ztabramp(:,:), dim=2)
+         ENDIF
 
          ! --- North --- !
-         ztabramp(:,:) = 0._wp
-         ind1 = jpjglo - nbghostcells - 1
-         DO jj = mj0(ind1), mj1(ind1)                 
-            ztabramp(:,jj) = ssvmask(:,jj)
-         END DO
-         !
-         zmsknorth(:) = 0._wp
-         zmsknorth(1:jpi) = MAXVAL(ztabramp(:,:), dim=2)
+         IF( lk_north) THEN
+            ztabramp(:,:) = 0._wp
+            ind1 = jpjglo - nbghostcells - 1
+            DO jj = mj0(ind1), mj1(ind1)                 
+               ztabramp(:,jj) = ssvmask(:,jj)
+            END DO
+            !
+            zmsknorth(:) = 0._wp
+            zmsknorth(1:jpi) = MAXVAL(ztabramp(:,:), dim=2)
+         ENDIF
+
          ! JC: SPONGE MASKING TO BE SORTED OUT:
          zmskwest(:)  = 1._wp
@@ -192,75 +204,85 @@
 
          ! --- West --- !
-         ind1 = 1+nbghostcells
-         ind2 = 1+nbghostcells + ispongearea 
-         DO ji = mi0(ind1), mi1(ind2)   
-            DO jj = 1, jpj               
-               ztabramp(ji,jj) = REAL( ind2 - mig(ji) ) * z1_ispongearea * zmskwest(jj)
-            END DO
-         END DO
-
-         ! ghost cells:
-         ind1 = 1
-         ind2 = nbghostcells + 1
-         DO ji = mi0(ind1), mi1(ind2)   
-            DO jj = 1, jpj               
-               ztabramp(ji,jj) = zmskwest(jj)
-            END DO
-         END DO
+         IF(lk_west) THEN
+            ind1 = 1+nbghostcells
+            ind2 = 1+nbghostcells + ispongearea 
+            DO ji = mi0(ind1), mi1(ind2)   
+               DO jj = 1, jpj               
+                  ztabramp(ji,jj) = REAL( ind2 - mig(ji) ) * z1_ispongearea * zmskwest(jj)
+               END DO
+            END DO         
+
+            ! ghost cells:
+            ind1 = 1
+            ind2 = nbghostcells + 1
+            DO ji = mi0(ind1), mi1(ind2)   
+               DO jj = 1, jpj               
+                  ztabramp(ji,jj) = zmskwest(jj)
+               END DO
+            END DO
+         ENDIF
 
          ! --- East --- !
-         ind1 = jpiglo - nbghostcells - ispongearea
-         ind2 = jpiglo - nbghostcells
-         DO ji = mi0(ind1), mi1(ind2)
-            DO jj = 1, jpj
-               ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL( mig(ji) - ind1 ) * z1_ispongearea) * zmskeast(jj)
-            ENDDO
-         END DO
-
-         ! ghost cells:
-         ind1 = jpiglo - nbghostcells
-         ind2 = jpiglo
-         DO ji = mi0(ind1), mi1(ind2)
-            DO jj = 1, jpj
-               ztabramp(ji,jj) = zmskeast(jj)
-            ENDDO
-         END DO
+         IF(lk_east) THEN
+            ind1 = jpiglo - nbghostcells - ispongearea
+            ind2 = jpiglo - nbghostcells
+            DO ji = mi0(ind1), mi1(ind2)
+
+               DO jj = 1, jpj
+                  ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL( mig(ji) - ind1 ) * z1_ispongearea) * zmskeast(jj)
+               ENDDO
+            END DO
+
+            ! ghost cells:
+            ind1 = jpiglo - nbghostcells
+            ind2 = jpiglo
+            DO ji = mi0(ind1), mi1(ind2)
+
+               DO jj = 1, jpj
+                  ztabramp(ji,jj) = zmskeast(jj)
+               ENDDO
+            END DO
+         ENDIF
 
          ! --- South --- !
-         ind1 = 1+nbghostcells
-         ind2 = 1+nbghostcells + jspongearea
-         DO jj = mj0(ind1), mj1(ind2) 
-            DO ji = 1, jpi
-               ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL( ind2 - mjg(jj) ) * z1_jspongearea) * zmsksouth(ji)
-            END DO
-         END DO
-
-         ! ghost cells:
-         ind1 = 1
-         ind2 = nbghostcells + 1
-         DO jj = mj0(ind1), mj1(ind2) 
-            DO ji = 1, jpi
-               ztabramp(ji,jj) = zmsksouth(ji)
-            END DO
-         END DO
+         IF( lk_south ) THEN 
+            ind1 = 1+nbghostcells
+            ind2 = 1+nbghostcells + jspongearea
+            DO jj = mj0(ind1), mj1(ind2) 
+               DO ji = 1, jpi
+                  ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL( ind2 - mjg(jj) ) * z1_jspongearea) * zmsksouth(ji)
+               END DO
+            END DO
+
+            ! ghost cells:
+            ind1 = 1
+            ind2 = nbghostcells + 1
+            DO jj = mj0(ind1), mj1(ind2) 
+               DO ji = 1, jpi
+                  ztabramp(ji,jj) = zmsksouth(ji)
+               END DO
+            END DO
+         ENDIF
 
          ! --- North --- !
-         ind1 = jpjglo - nbghostcells - jspongearea
-         ind2 = jpjglo - nbghostcells
-         DO jj = mj0(ind1), mj1(ind2)
-            DO ji = 1, jpi
-               ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL( mjg(jj) - ind1 ) * z1_jspongearea) * zmsknorth(ji)
-            END DO
-         END DO
-
-         ! ghost cells:
-         ind1 = jpjglo - nbghostcells
-         ind2 = jpjglo
-         DO jj = mj0(ind1), mj1(ind2)
-            DO ji = 1, jpi
-               ztabramp(ji,jj) = zmsknorth(ji)
-            END DO
-         END DO
-
+         IF( lk_north ) THEN  
+            ind1 = jpjglo - nbghostcells - jspongearea
+            ind2 = jpjglo - nbghostcells
+            DO jj = mj0(ind1), mj1(ind2)
+               DO ji = 1, jpi
+                  ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL( mjg(jj) - ind1 ) * z1_jspongearea) * zmsknorth(ji)
+               END DO
+            END DO
+
+            ! ghost cells:
+            ind1 = jpjglo - nbghostcells
+            ind2 = jpjglo
+            DO jj = mj0(ind1), mj1(ind2)
+               DO ji = 1, jpi
+                  ztabramp(ji,jj) = zmsknorth(ji)
+               END DO
+            END DO
+         ENDIF
+      
       ENDIF
 
@@ -646,4 +668,5 @@
 
          jmax = j2-1
+        ! IF (lk_north) jmax = MIN(jmax,nlcj-nbghostcells-2)   ! North
          IF ((nbondj == 1).OR.(nbondj == 2)) jmax = MIN(jmax,nlcj-nbghostcells-2)   ! North
 
@@ -802,4 +825,5 @@
 
          imax = i2 - 1
+      !   IF(lk_east) imax = MIN(imax,nlci-nbghostcells-2)   ! East
          IF ((nbondi == 1).OR.(nbondi == 2))   imax = MIN(imax,nlci-nbghostcells-2)   ! East
 
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_oce_update.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_oce_update.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_oce_update.F90	(revision 13026)
@@ -26,4 +26,5 @@
    USE domvvl         ! Need interpolation routines 
    USE vremap         ! Vertical remapping
+   USE lbclnk 
 
    IMPLICIT NONE
@@ -81,8 +82,13 @@
       IF (Agrif_Root()) RETURN
       !
+#if defined TWO_WAY
       IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Update momentum from grid Number',Agrif_Fixed()
 
       Agrif_UseSpecialValueInUpdate = .FALSE.
       Agrif_SpecialValueFineGrid = 0.
+
+      use_sign_north = .TRUE.
+      sign_north = -1.
+
       !     
 # if ! defined DECAL_FEEDBACK
@@ -126,4 +132,6 @@
 #  endif
       END IF
+#endif
+      use_sign_north = .FALSE.
       !
    END SUBROUTINE Agrif_Update_Dyn
@@ -148,4 +156,6 @@
 #  if defined VOL_REFLUX
       IF ( ln_dynspg_ts.AND.ln_bt_fw ) THEN
+         use_sign_north = .TRUE.
+         sign_north = -1.
          ! Refluxing on ssh:
 #  if defined DECAL_FEEDBACK_2D
@@ -156,4 +166,5 @@
          CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/ 0, 0/),locupdate2=(/-1,-1/),procname = reflux_sshv)
 #  endif
+         use_sign_north = .FALSE.
       END IF
 #  endif
@@ -826,5 +837,5 @@
    SUBROUTINE correct_v_bdy( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before, nb, ndir )
       !!---------------------------------------------
-      !!           *** ROUTINE correct_u_bdy ***
+      !!           *** ROUTINE correct_v_bdy ***
       !!---------------------------------------------
       INTEGER                                     , INTENT(in   ) :: i1, i2, j1, j2, k1, k2, n1, n2
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_top_interp.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_top_interp.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_top_interp.F90	(revision 13026)
@@ -119,5 +119,4 @@
             tr(i1:i2,j1:j2,1:jpk,jn,Krhs_a)=ptab_child(i1:i2,j1:j2,1:jpk,jn)*tmask(i1:i2,j1:j2,1:jpk) 
          END DO
-
       ENDIF
       !
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_user.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_user.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/NST/agrif_user.F90	(revision 13026)
@@ -28,6 +28,4 @@
       !
       !                    !* Agrif initialization
-      CALL agrif_nemo_init
-      CALL Agrif_InitValues_cont_dom
       CALL Agrif_InitValues_cont
 # if defined key_top
@@ -40,41 +38,147 @@
    END SUBROUTINE Agrif_initvalues
 
-   SUBROUTINE Agrif_InitValues_cont_dom
-      !!----------------------------------------------------------------------
-      !!                 *** ROUTINE Agrif_InitValues_cont_dom ***
-      !!----------------------------------------------------------------------
-      !
-      CALL agrif_declare_var_dom
-      !
-   END SUBROUTINE Agrif_InitValues_cont_dom
-
-   SUBROUTINE agrif_declare_var_dom
-      !!----------------------------------------------------------------------
-      !!                 *** ROUTINE agrif_declare_var_dom ***
-      !!----------------------------------------------------------------------
-      USE par_oce, ONLY:  nbghostcells      
+   SUBROUTINE agrif_istate( Kbb, Kmm, Kaa )
+
+       USE domvvl
+       USE domain
+       USE par_oce
+       USE agrif_oce
+       USE agrif_oce_interp
+       USE oce
+       USE lib_mpp
+       USe lbclnk
+
+      INTEGER, INTENT(in)  :: Kbb, Kmm, Kaa
+      INTEGER :: jn
+
+      l_ini_child = .TRUE.
+      Agrif_SpecialValue    = 0._wp
+      Agrif_UseSpecialValue = .TRUE.
+      uu(:,:,:,:) = 0.  ;  vv(:,:,:,:) = 0.   ;  ts(:,:,:,:,:) = 0.
+       
+      Krhs_a = Kbb ; Kmm_a = Kbb
+
+      ! Brutal fix to pas 1x1 refinment. 
+  !    IF(Agrif_Irhox() == 1) THEN
+         CALL Agrif_Init_Variable(tsini_id, procname=agrif_initts)
+  !    ELSE
+   !      CALL Agrif_Init_Variable(tsini_id, procname=interptsn)
+
+  !    ENDIF
+      Agrif_UseSpecialValue = ln_spc_dyn
+      use_sign_north = .TRUE.
+      sign_north = -1.
+  !    CALL Agrif_Init_Variable(uini_id , procname=interpun )
+  !    CALL Agrif_Init_Variable(vini_id , procname=interpvn )
+       use_sign_north = .FALSE.
+
+      Agrif_UseSpecialValue = .FALSE.            !
+      l_ini_child = .FALSE.
+      Krhs_a = Kaa ; Kmm_a = Kmm
+
+      DO jn = 1, jpts
+         ts(:,:,:,jn,Kbb) = ts(:,:,:,jn,Kbb)*tmask(:,:,:)
+      END DO
+      uu(:,:,:,Kbb) =  uu(:,:,:,Kbb) * umask(:,:,:)     
+      vv(:,:,:,Kbb) =  vv(:,:,:,Kbb) * vmask(:,:,:) 
+
+
+      CALL lbc_lnk_multi( 'agrif_istate', uu(:,:,:,Kbb), 'U', -1. , vv(:,:,:,Kbb), 'V', -1. )
+      CALL lbc_lnk( 'agrif_istate', ts(:,:,:,:,Kbb), 'T', 1. )
+
+   END SUBROUTINE agrif_istate   
+
+   SUBROUTINE agrif_declare_var_ini
+      !!----------------------------------------------------------------------
+      !!                 *** ROUTINE agrif_declare_var ***
+      !!----------------------------------------------------------------------
+      USE agrif_util
+      USE agrif_oce
+      USE par_oce
+      USE zdf_oce 
+      USE oce
+      USE dom_oce
       !
       IMPLICIT NONE
       !
       INTEGER :: ind1, ind2, ind3
-      !!----------------------------------------------------------------------
+      External :: nemo_mapping
+      !!----------------------------------------------------------------------
+
+! In case of East-West periodicity, prevent AGRIF interpolation at east and west boundaries
+! The procnames will not be called at these boundaries
+      IF (jperio == 1) THEN
+         CALL Agrif_Set_NearCommonBorderX(.TRUE.)
+         CALL Agrif_Set_DistantCommonBorderX(.TRUE.)
+      ENDIF
+
+      IF ( .NOT. ln_bry_south) THEN
+         CALL Agrif_Set_NearCommonBorderY(.TRUE.)
+      ENDIF
 
       ! 1. Declaration of the type of variable which have to be interpolated
       !---------------------------------------------------------------------
       ind1 =     nbghostcells
-      ind2 = 1 + nbghostcells
-      ind3 = 2 + nbghostcells
-      CALL agrif_declare_variable((/1,2/),(/ind2,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),e1u_id)
-      CALL agrif_declare_variable((/2,1/),(/ind3,ind2/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),e2v_id)
-
+      ind2 = 2 + nbghostcells_x
+      ind3 = 2 + nbghostcells_y_s
+
+      CALL agrif_declare_variable((/2,2,0/),(/ind3,ind3,0/),(/'x','y','N'/),(/1,1,1/),(/nlci,nlcj,jpk/),e3t_id)
+      CALL agrif_declare_variable((/2,2/),(/ind3,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),mbkt_id)
+      CALL agrif_declare_variable((/2,2/),(/ind3,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),ht0_id)
+
+      CALL agrif_declare_variable((/1,2/),(/ind2-1,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),e1u_id)
+      CALL agrif_declare_variable((/2,1/),(/ind2,ind3-1/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),e2v_id)
+
+   
+      ! Initial or restart velues
+      CALL Agrif_Set_MaskMaxSearch(25)
+      !
+      CALL agrif_declare_variable((/2,2,0,0/),(/ind3,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,jpts+1/),tsini_id)
+      CALL agrif_declare_variable((/1,2,0,0/),(/ind2-1,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,2/)     ,uini_id ) 
+      CALL agrif_declare_variable((/2,1,0,0/),(/ind2,ind3-1,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,2/)     ,vini_id )
+      CALL agrif_declare_variable((/2,2/),(/ind3,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),sshini_id)
+      ! 
+      CALL Agrif_Set_MaskMaxSearch(5)
+     
       ! 2. Type of interpolation
       !-------------------------
+      CALL Agrif_Set_bcinterp(e3t_id,interp=AGRIF_constant)
+
+      CALL Agrif_Set_bcinterp(mbkt_id,interp=AGRIF_constant)
+      CALL Agrif_Set_interp  (mbkt_id,interp=AGRIF_constant)
+      CALL Agrif_Set_bcinterp(ht0_id ,interp=AGRIF_constant)
+      CALL Agrif_Set_interp  (ht0_id ,interp=AGRIF_constant)
+
       CALL Agrif_Set_bcinterp( e1u_id, interp1=Agrif_linear, interp2=AGRIF_ppm    )
       CALL Agrif_Set_bcinterp( e2v_id, interp1=AGRIF_ppm   , interp2=Agrif_linear )
 
-      ! 3. Location of interpolation
+      ! Initial fields
+      CALL Agrif_Set_bcinterp(tsini_id ,interp=AGRIF_linear)
+      CALL Agrif_Set_interp  (tsini_id ,interp=AGRIF_linear)
+      CALL Agrif_Set_bcinterp(uini_id  ,interp=AGRIF_linear)
+      CALL Agrif_Set_interp  (uini_id  ,interp=AGRIF_linear)
+      CALL Agrif_Set_bcinterp(vini_id  ,interp=AGRIF_linear)
+      CALL Agrif_Set_interp  (vini_id  ,interp=AGRIF_linear)
+      CALL Agrif_Set_bcinterp(sshini_id,interp=AGRIF_linear)
+      CALL Agrif_Set_interp  (sshini_id,interp=AGRIF_linear)
+
+       ! 3. Location of interpolation
       !-----------------------------
+!      CALL Agrif_Set_bc(  e3t_id, (/-nn_sponge_len*Agrif_irhox(),ind1-1/) )  
+! JC: check near the boundary only until matching in sponge has been sorted out:
+      CALL Agrif_Set_bc(  e3t_id, (/0,ind1-1/) )  
+
+      ! extend the interpolation zone by 1 more point than necessary:
+      ! RB check here
+      CALL Agrif_Set_bc(  mbkt_id, (/-nn_sponge_len*Agrif_irhox()-2,ind1/) )
+      CALL Agrif_Set_bc(  ht0_id,  (/-nn_sponge_len*Agrif_irhox()-2,ind1/) )
+      
       CALL Agrif_Set_bc(e1u_id,(/0,ind1-1/))
-      CALL Agrif_Set_bc(e2v_id,(/0,ind1-1/))
+      CALL Agrif_Set_bc(e2v_id,(/0,ind1-1/))  
+
+      CALL Agrif_Set_bc( tsini_id , (/0,ind1-1/) ) ! if west,  rhox=3 and nbghost=3: columns 2 to 4
+      CALL Agrif_Set_bc( uini_id  , (/0,ind1-1/) ) 
+      CALL Agrif_Set_bc( vini_id  , (/0,ind1-1/) )
+      CALL Agrif_Set_bc( sshini_id, (/0,ind1-1/) )
 
       ! 4. Update type
@@ -87,37 +191,63 @@
       CALL Agrif_Set_Updatetype(e2v_id,update1 = Agrif_Update_Average, update2=Agrif_Update_Copy)
 #endif
-
-   END SUBROUTINE agrif_declare_var_dom
-
-   SUBROUTINE Agrif_InitValues_cont
-      !!----------------------------------------------------------------------
-      !!                 *** ROUTINE Agrif_InitValues_cont ***
-      !!----------------------------------------------------------------------
-      USE agrif_oce
+      
+      CALL Agrif_Set_ExternalMapping(nemo_mapping)
+      !
+   END SUBROUTINE agrif_declare_var_ini
+
+
+   SUBROUTINE Agrif_Init_Domain( Kbb, Kmm, Kaa ) 
+      !!----------------------------------------------------------------------
+      !!                 *** ROUTINE Agrif_InitValues_cont_dom ***
+      !!----------------------------------------------------------------------
+  
+         !!----------------------------------------------------------------------
+         !!                 *** ROUTINE Agrif_InitValues_cont ***
+         !!
+         !! ** Purpose ::   Declaration of variables to be interpolated
+         !!----------------------------------------------------------------------
+      USE agrif_oce_update
       USE agrif_oce_interp
       USE agrif_oce_sponge
+      USE Agrif_Util
+      USE oce 
       USE dom_oce
-      USE oce
+      USE zdf_oce
+      USE nemogcm
+      USE agrif_oce
+      !
+      USE lbclnk
       USE lib_mpp
-      USE lbclnk
+      USE in_out_manager
       !
       IMPLICIT NONE
       !
-      INTEGER :: ji, jj
+      INTEGER, INTENT(in) ::  Kbb, Kmm, Kaa
+      !
       LOGICAL :: check_namelist
       CHARACTER(len=15) :: cl_check1, cl_check2, cl_check3, cl_check4 
-#if defined key_vertical
       REAL(wp), DIMENSION(jpi,jpj) ::   zk   ! workspace
-#endif
-      !!----------------------------------------------------------------------
-
-      ! 1. Declaration of the type of variable which have to be interpolated
-      !---------------------------------------------------------------------
-      CALL agrif_declare_var
-
-      ! 2. First interpolations of potentially non zero fields
-      !-------------------------------------------------------
-
-#if defined key_vertical
+      INTEGER :: ji, jj, jk, iminspon
+      !!----------------------------------------------------------------------
+    
+     ! CALL Agrif_Declare_Var_ini
+
+      IF( agrif_oce_alloc()  > 0 )   CALL ctl_warn('agrif agrif_oce_alloc: allocation of arrays failed')
+
+    !  lk_west  = ( ((nbondi == -1) .OR. (nbondi == 2) ).AND. .NOT. (jperio == 1 .OR. jperio == 4 .OR. jperio == 6))
+    !  lk_east  = ( ((nbondi ==  1) .OR. (nbondi == 2) ).AND. .NOT. (jperio == 1 .OR. jperio == 4 .OR. jperio == 6))
+    !  lk_south = ( ((nbondj == -1) .OR. (nbondj == 2) ).AND. ln_bry_south)
+    !  lk_north = ( ((nbondj ==  1) .OR. (nbondj == 2) ))
+    
+      lk_west  = ( .NOT. (jperio == 1 .OR. jperio == 4 .OR. jperio == 6) )
+      lk_east  = ( .NOT. (jperio == 1 .OR. jperio == 4 .OR. jperio == 6) )
+      lk_south = ln_bry_south
+      lk_north = .true.
+
+      ! Check sponge length:
+      iminspon = MIN(FLOOR(REAL(jpiglo-4)/REAL(2*Agrif_irhox())), FLOOR(REAL(jpjglo-4)/REAL(2*Agrif_irhox())) )
+      IF (lk_mpp) iminspon = MIN(iminspon,FLOOR(REAL(jpi-2)/REAL(Agrif_irhox())), FLOOR(REAL(jpj-2)/REAL(Agrif_irhox())) )
+      IF (nn_sponge_len > iminspon)  CALL ctl_stop('agrif sponge length is too large') 
+      
       ! Build consistent parent bathymetry and number of levels
       ! on the child grid 
@@ -126,6 +256,8 @@
       mbkt_parent(:,:) = 0
       !
-      CALL Agrif_Bc_variable(ht0_id ,calledweight=1.,procname=interpht0 )
-      CALL Agrif_Bc_variable(mbkt_id,calledweight=1.,procname=interpmbkt)
+  !    CALL Agrif_Bc_variable(ht0_id ,calledweight=1.,procname=interpht0 )
+  !    CALL Agrif_Bc_variable(mbkt_id,calledweight=1.,procname=interpmbkt)
+      CALL Agrif_Init_Variable(ht0_id , procname=interpht0 )
+      CALL Agrif_Init_Variable(mbkt_id, procname=interpmbkt)
       !
       ! Assume step wise change of bathymetry near interface
@@ -149,12 +281,97 @@
       ENDIF
       !
-      CALL lbc_lnk( 'Agrif_InitValues_cont', hu0_parent, 'U', 1. )
-      CALL lbc_lnk( 'Agrif_InitValues_cont', hv0_parent, 'V', 1. )
+      CALL lbc_lnk( 'Agrif_Init_Domain', hu0_parent, 'U', 1. )
+      CALL lbc_lnk( 'Agrif_Init_Domain', hv0_parent, 'V', 1. )
       zk(:,:) = REAL( mbku_parent(:,:), wp )   ;   CALL lbc_lnk( 'Agrif_InitValues_cont', zk, 'U', 1. )
-      mbku_parent(:,:) = MAX( NINT( zk(:,:) ), 1 )
+      mbku_parent(:,:) = MAX( NINT( zk(:,:) ), 1 ) ;
       zk(:,:) = REAL( mbkv_parent(:,:), wp )   ;   CALL lbc_lnk( 'Agrif_InitValues_cont', zk, 'V', 1. )
       mbkv_parent(:,:) = MAX( NINT( zk(:,:) ), 1 )   
-#endif
-
+
+
+      CALL Agrif_Init_Variable(sshini_id, procname=agrif_initssh)
+      CALL lbc_lnk( 'Agrif_Init_Domain', ssh(:,:,Kbb), 'T', 1. )
+      DO jk = 1, jpk
+            e3t(:,:,jk,Kbb) =  e3t_0(:,:,jk) * ( ht_0(:,:) + ssh(:,:,Kbb)  ) &
+      &                            / ( ht_0(:,:) + 1._wp - ssmask(:,:) ) * tmask(:,:,jk)   &
+                     &              + e3t_0(:,:,jk) * ( 1._wp - tmask(:,:,jk) )
+      END DO
+
+      ! check if masks and bathymetries match
+      IF(ln_chk_bathy) THEN
+         Agrif_UseSpecialValue = .FALSE.
+         !
+         IF(lwp) WRITE(numout,*) ' '
+         IF(lwp) WRITE(numout,*) 'AGRIF: Check Bathymetry and masks near bdys. Level: ', Agrif_Level()
+         !
+         kindic_agr = 0
+         IF( .NOT. l_vremap ) THEN
+            !
+            ! check if tmask and vertical scale factors agree with parent in sponge area:
+            CALL Agrif_Bc_variable(e3t_id,calledweight=1.,procname=interpe3t)
+            !
+         ELSE
+            !
+            ! In case of vertical interpolation, check only that total depths agree between child and parent:
+            DO ji = 1, jpi
+               DO jj = 1, jpj
+                  IF ((mbkt_parent(ji,jj)/=0).AND.(ABS(ht0_parent(ji,jj)-ht_0(ji,jj))>1.e-3)) kindic_agr = kindic_agr + 1
+                  IF ((mbku_parent(ji,jj)/=0).AND.(ABS(hu0_parent(ji,jj)-hu_0(ji,jj))>1.e-3)) kindic_agr = kindic_agr + 1
+                  IF ((mbkv_parent(ji,jj)/=0).AND.(ABS(hv0_parent(ji,jj)-hv_0(ji,jj))>1.e-3)) kindic_agr = kindic_agr + 1
+               END DO
+            END DO
+
+            CALL mpp_sum( 'agrif_user', kindic_agr )
+            IF( kindic_agr /= 0 ) THEN
+               CALL ctl_stop('==> Child Bathymetry is NOT correct near boundaries.')
+            ELSE
+               IF(lwp) WRITE(numout,*) '==> Child Bathymetry is ok near boundaries.'
+               IF(lwp) WRITE(numout,*) ' '
+            ENDIF  
+         ENDIF
+      ENDIF
+
+      IF( l_vremap ) THEN
+      ! Additional constrain that should be removed someday:
+         IF ( Agrif_Parent(jpk).GT.jpk ) THEN
+            CALL ctl_stop( ' With l_vremap, child grids must have jpk greater or equal to the parent value' )
+         ENDIF
+      ENDIF
+      !
+   END SUBROUTINE Agrif_Init_Domain
+
+
+   SUBROUTINE Agrif_InitValues_cont
+         !!----------------------------------------------------------------------
+         !!                 *** ROUTINE Agrif_InitValues_cont ***
+         !!
+         !! ** Purpose ::   Declaration of variables to be interpolated
+         !!----------------------------------------------------------------------
+      USE agrif_oce_update
+      USE agrif_oce_interp
+      USE agrif_oce_sponge
+      USE Agrif_Util
+      USE oce 
+      USE dom_oce
+      USE zdf_oce
+      USE nemogcm
+      USE agrif_oce
+      !
+      USE lbclnk
+      USE lib_mpp
+      USE in_out_manager
+      !
+      IMPLICIT NONE
+      !
+      LOGICAL :: check_namelist
+      CHARACTER(len=15) :: cl_check1, cl_check2, cl_check3, cl_check4 
+      REAL(wp), DIMENSION(jpi,jpj) ::   zk   ! workspace
+      INTEGER :: ji, jj
+
+      ! 1. Declaration of the type of variable which have to be interpolated
+      !---------------------------------------------------------------------
+      CALL agrif_declare_var
+
+      ! 2. First interpolations of potentially non zero fields
+      !-------------------------------------------------------
       Agrif_SpecialValue    = 0._wp
       Agrif_UseSpecialValue = .TRUE.
@@ -163,8 +380,10 @@
       tabspongedone_tsn = .FALSE.
       CALL Agrif_Bc_variable(tsn_sponge_id,calledweight=1.,procname=interptsn_sponge)
-      ! reset ts(:,:,:,:,Krhs_a) to zero
+      ! reset tsa to zero
       ts(:,:,:,:,Krhs_a) = 0._wp
 
       Agrif_UseSpecialValue = ln_spc_dyn
+      use_sign_north = .TRUE.
+      sign_north = -1.
       CALL Agrif_Bc_variable(un_interp_id,calledweight=1.,procname=interpun)
       CALL Agrif_Bc_variable(vn_interp_id,calledweight=1.,procname=interpvn)
@@ -175,4 +394,5 @@
       tabspongedone_v = .FALSE.
       CALL Agrif_Bc_variable(vn_sponge_id,calledweight=1.,procname=interpvn_sponge)
+      use_sign_north = .FALSE.
       uu(:,:,:,Krhs_a) = 0._wp
       vv(:,:,:,Krhs_a) = 0._wp
@@ -185,43 +405,20 @@
       IF ( ln_dynspg_ts ) THEN
          Agrif_UseSpecialValue = ln_spc_dyn
+         use_sign_north = .TRUE.
+         sign_north = -1.
          CALL Agrif_Bc_variable(unb_id,calledweight=1.,procname=interpunb)
          CALL Agrif_Bc_variable(vnb_id,calledweight=1.,procname=interpvnb)
          CALL Agrif_Bc_variable(ub2b_interp_id,calledweight=1.,procname=interpub2b)
          CALL Agrif_Bc_variable(vb2b_interp_id,calledweight=1.,procname=interpvb2b)
+         use_sign_north = .FALSE.
          ubdy(:,:) = 0._wp
          vbdy(:,:) = 0._wp
       ENDIF
-
-      Agrif_UseSpecialValue = .FALSE.
-
-      ! 3. Some controls
+      Agrif_UseSpecialValue = .FALSE. 
+
       !-----------------
       check_namelist = .TRUE.
 
       IF( check_namelist ) THEN 
-
-         ! Check time steps           
-         IF( NINT(Agrif_Rhot()) * NINT(rn_Dt) .NE. Agrif_Parent(rn_Dt) ) THEN
-            WRITE(cl_check1,*)  NINT(Agrif_Parent(rn_Dt))
-            WRITE(cl_check2,*)  NINT(rn_Dt)
-            WRITE(cl_check3,*)  NINT(Agrif_Parent(rn_Dt)/Agrif_Rhot())
-            CALL ctl_stop( 'Incompatible time step between ocean grids',   &
-                  &               'parent grid value : '//cl_check1    ,   & 
-                  &               'child  grid value : '//cl_check2    ,   & 
-                  &               'value on child grid should be changed to : '//cl_check3 )
-         ENDIF
-
-         ! Check run length
-         IF( Agrif_IRhot() * (Agrif_Parent(nitend)- &
-               Agrif_Parent(nit000)+1) .NE. (nitend-nit000+1) ) THEN
-            WRITE(cl_check1,*)  (Agrif_Parent(nit000)-1)*Agrif_IRhot() + 1
-            WRITE(cl_check2,*)   Agrif_Parent(nitend)   *Agrif_IRhot()
-            CALL ctl_warn( 'Incompatible run length between grids'                      ,   &
-                  &               'nit000 on fine grid will be changed to : '//cl_check1,   &
-                  &               'nitend on fine grid will be changed to : '//cl_check2    )
-            nit000 = (Agrif_Parent(nit000)-1)*Agrif_IRhot() + 1
-            nitend =  Agrif_Parent(nitend)   *Agrif_IRhot()
-         ENDIF
-
          ! Check free surface scheme
          IF ( ( Agrif_Parent(ln_dynspg_ts ).AND.ln_dynspg_exp ).OR.&
@@ -251,48 +448,6 @@
             STOP
          ENDIF
-
-      ENDIF
-
-      ! check if masks and bathymetries match
-      IF(ln_chk_bathy) THEN
-         Agrif_UseSpecialValue = .FALSE.
-         !
-         IF(lwp) WRITE(numout,*) ' '
-         IF(lwp) WRITE(numout,*) 'AGRIF: Check Bathymetry and masks near bdys. Level: ', Agrif_Level()
-         !
-         kindic_agr = 0
-# if ! defined key_vertical
-         !
-         ! check if tmask and vertical scale factors agree with parent in sponge area:
-         CALL Agrif_Bc_variable(e3t_id,calledweight=1.,procname=interpe3t)
-         !
-# else
-         !
-         ! In case of vertical interpolation, check only that total depths agree between child and parent:
-         DO ji = 1, jpi
-            DO jj = 1, jpj
-               IF ((mbkt_parent(ji,jj)/=0).AND.(ABS(ht0_parent(ji,jj)-ht_0(ji,jj))>1.e-3)) kindic_agr = kindic_agr + 1
-               IF ((mbku_parent(ji,jj)/=0).AND.(ABS(hu0_parent(ji,jj)-hu_0(ji,jj))>1.e-3)) kindic_agr = kindic_agr + 1
-               IF ((mbkv_parent(ji,jj)/=0).AND.(ABS(hv0_parent(ji,jj)-hv_0(ji,jj))>1.e-3)) kindic_agr = kindic_agr + 1
-            END DO
-         END DO
-# endif
-         CALL mpp_sum( 'agrif_user', kindic_agr )
-         IF( kindic_agr /= 0 ) THEN
-            CALL ctl_stop('==> Child Bathymetry is NOT correct near boundaries.')
-         ELSE
-            IF(lwp) WRITE(numout,*) '==> Child Bathymetry is ok near boundaries.'
-            IF(lwp) WRITE(numout,*) ' '
-         END IF  
-         !    
-      ENDIF
-
-# if defined key_vertical
-      ! Additional constrain that should be removed someday:
-      IF ( Agrif_Parent(jpk).GT.jpk ) THEN
-	 CALL ctl_stop( ' With key_vertical, child grids must have jpk greater or equal to the parent value' )
-      ENDIF
-# endif
-      ! 
+      ENDIF
+
    END SUBROUTINE Agrif_InitValues_cont
 
@@ -314,37 +469,30 @@
       ! 1. Declaration of the type of variable which have to be interpolated
       !---------------------------------------------------------------------
+
       ind1 =     nbghostcells
-      ind2 = 1 + nbghostcells
-      ind3 = 2 + nbghostcells
+      ind2 = 2 + nbghostcells_x
+      ind3 = 2 + nbghostcells_y_s
+
 # if defined key_vertical
       CALL agrif_declare_variable((/2,2,0,0/),(/ind3,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,jpts+1/),tsn_id)
       CALL agrif_declare_variable((/2,2,0,0/),(/ind3,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,jpts+1/),tsn_sponge_id)
 
-      CALL agrif_declare_variable((/1,2,0,0/),(/ind2,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,2/),un_interp_id)
-      CALL agrif_declare_variable((/2,1,0,0/),(/ind3,ind2,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,2/),vn_interp_id)
-      CALL agrif_declare_variable((/1,2,0,0/),(/ind2,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,2/),un_update_id)
-      CALL agrif_declare_variable((/2,1,0,0/),(/ind3,ind2,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,2/),vn_update_id)
-      CALL agrif_declare_variable((/1,2,0,0/),(/ind2,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,2/),un_sponge_id)
-      CALL agrif_declare_variable((/2,1,0,0/),(/ind3,ind2,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,2/),vn_sponge_id)
+      CALL agrif_declare_variable((/1,2,0,0/),(/ind2-1,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,2/),un_interp_id) !
+      CALL agrif_declare_variable((/2,1,0,0/),(/ind2,ind3-1,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,2/),vn_interp_id)
+      CALL agrif_declare_variable((/1,2,0,0/),(/ind2-1,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,2/),un_update_id)
+      CALL agrif_declare_variable((/2,1,0,0/),(/ind2,ind3-1,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,2/),vn_update_id)
+      CALL agrif_declare_variable((/1,2,0,0/),(/ind2-1,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,2/),un_sponge_id)
+      CALL agrif_declare_variable((/2,1,0,0/),(/ind2,ind3-1,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,2/),vn_sponge_id)
 # else
       CALL agrif_declare_variable((/2,2,0,0/),(/ind3,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,jpts/),tsn_id)
       CALL agrif_declare_variable((/2,2,0,0/),(/ind3,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,jpts/),tsn_sponge_id)
 
-      CALL agrif_declare_variable((/1,2,0,0/),(/ind2,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,1/),un_interp_id)
-      CALL agrif_declare_variable((/2,1,0,0/),(/ind3,ind2,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,1/),vn_interp_id)
-      CALL agrif_declare_variable((/1,2,0,0/),(/ind2,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,1/),un_update_id)
-      CALL agrif_declare_variable((/2,1,0,0/),(/ind3,ind2,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,1/),vn_update_id)
-      CALL agrif_declare_variable((/1,2,0,0/),(/ind2,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,1/),un_sponge_id)
-      CALL agrif_declare_variable((/2,1,0,0/),(/ind3,ind2,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,1/),vn_sponge_id)
+      CALL agrif_declare_variable((/1,2,0,0/),(/ind2-1,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,1/),un_interp_id)
+      CALL agrif_declare_variable((/2,1,0,0/),(/ind2,ind3-1,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,1/),vn_interp_id)
+      CALL agrif_declare_variable((/1,2,0,0/),(/ind2-1,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,1/),un_update_id)
+      CALL agrif_declare_variable((/2,1,0,0/),(/ind2,ind3-1,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,1/),vn_update_id)
+      CALL agrif_declare_variable((/1,2,0,0/),(/ind2-1,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,1/),un_sponge_id)
+      CALL agrif_declare_variable((/2,1,0,0/),(/ind2,ind3-1,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,1/),vn_sponge_id)
 # endif
-
-      CALL agrif_declare_variable((/2,2,0/),(/ind3,ind3,0/),(/'x','y','N'/),(/1,1,1/),(/nlci,nlcj,jpk/),e3t_id)
-
-# if defined key_vertical
-      CALL agrif_declare_variable((/2,2/),(/ind3,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),mbkt_id)
-      CALL agrif_declare_variable((/2,2/),(/ind3,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),ht0_id)
-# endif
-
-      CALL agrif_declare_variable((/2,2,0,0/),(/ind3,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,3/),scales_t_id)
 
       CALL agrif_declare_variable((/1,2/),(/ind2,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),unb_id)
@@ -357,5 +505,6 @@
       CALL agrif_declare_variable((/2,2/),(/ind3,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),sshn_id)
 
-      IF( ln_zdftke.OR.ln_zdfgls ) THEN
+
+      IF( ln_zdftke.OR.ln_zdfgls ) THEN  ! logical not known at this point
 !         CALL agrif_declare_variable((/2,2,0/),(/ind3,ind3,0/),(/'x','y','N'/),(/1,1,1/),(/nlci,nlcj,jpk/), en_id)
 !         CALL agrif_declare_variable((/2,2,0/),(/ind3,ind3,0/),(/'x','y','N'/),(/1,1,1/),(/nlci,nlcj,jpk/),avt_id)
@@ -366,13 +515,14 @@
 # endif
       ENDIF
-
+     
       ! 2. Type of interpolation
       !-------------------------
       CALL Agrif_Set_bcinterp(tsn_id,interp=AGRIF_linear)
-
       CALL Agrif_Set_bcinterp(un_interp_id,interp1=Agrif_linear,interp2=AGRIF_ppm)
       CALL Agrif_Set_bcinterp(vn_interp_id,interp1=AGRIF_ppm,interp2=Agrif_linear)
 
       CALL Agrif_Set_bcinterp(tsn_sponge_id,interp=AGRIF_linear)
+      CALL Agrif_Set_bcinterp(un_sponge_id,interp1=Agrif_linear,interp2=AGRIF_ppm)
+      CALL Agrif_Set_bcinterp(vn_sponge_id,interp1=AGRIF_ppm,interp2=Agrif_linear)
 
       CALL Agrif_Set_bcinterp(sshn_id,interp=AGRIF_linear)
@@ -390,17 +540,8 @@
 !<
 
-      CALL Agrif_Set_bcinterp(un_sponge_id,interp1=Agrif_linear,interp2=AGRIF_ppm)
-      CALL Agrif_Set_bcinterp(vn_sponge_id,interp1=AGRIF_ppm,interp2=Agrif_linear)
-
-      CALL Agrif_Set_bcinterp(e3t_id,interp=AGRIF_constant)
-
-# if defined key_vertical
-      CALL Agrif_Set_bcinterp(mbkt_id,interp=AGRIF_constant)
-      CALL Agrif_Set_bcinterp(ht0_id ,interp=AGRIF_constant)
-# endif
-
-      IF( ln_zdftke.OR.ln_zdfgls )   CALL Agrif_Set_bcinterp( avm_id, interp=AGRIF_linear )
-
-      ! 3. Location of interpolation
+      IF( ln_zdftke.OR.ln_zdfgls )  CALL Agrif_Set_bcinterp( avm_id, interp=AGRIF_linear )
+    
+
+       ! 3. Location of interpolation
       !-----------------------------
       CALL Agrif_Set_bc(       tsn_id, (/0,ind1-1/) ) ! if west,  rhox=3 and nbghost=3: columns 2 to 4
@@ -418,19 +559,8 @@
       CALL Agrif_Set_bc( vb2b_interp_id, (/0,ind1-1/) )
 
-!      CALL Agrif_Set_bc(  e3t_id, (/-nn_sponge_len*Agrif_irhox(),ind1-1/) )  
-! JC: check near the boundary only until matching in sponge has been sorted out:
-      CALL Agrif_Set_bc(  e3t_id, (/0,ind1-1/) )  
-
-# if defined key_vertical 
-      ! extend the interpolation zone by 1 more point than necessary:
-      CALL Agrif_Set_bc(  mbkt_id, (/-nn_sponge_len*Agrif_irhox()-2,ind1/) )
-      CALL Agrif_Set_bc(  ht0_id,  (/-nn_sponge_len*Agrif_irhox()-2,ind1/) )
-# endif
-
-      IF( ln_zdftke.OR.ln_zdfgls )   CALL Agrif_Set_bc( avm_id, (/0,ind1/) )
+      IF( ln_zdftke.OR.ln_zdfgls ) CALL Agrif_Set_bc( avm_id, (/0,ind1/) )
 
       ! 4. Update type
       !--------------- 
-      CALL Agrif_Set_Updatetype(scales_t_id, update = AGRIF_Update_Average)
 
 # if defined UPD_HIGH
@@ -444,9 +574,9 @@
       CALL Agrif_Set_Updatetype(e3t_id, update = Agrif_Update_Full_Weighting)
 
-      IF( ln_zdftke.OR.ln_zdfgls ) THEN
+  !    IF( ln_zdftke.OR.ln_zdfgls ) THEN
 !         CALL Agrif_Set_Updatetype( en_id, update = AGRIF_Update_Full_Weighting)
 !         CALL Agrif_Set_Updatetype(avt_id, update = AGRIF_Update_Full_Weighting)
 !         CALL Agrif_Set_Updatetype(avm_id, update = AGRIF_Update_Full_Weighting)
-      ENDIF
+   !   ENDIF
 
 #else
@@ -460,9 +590,9 @@
       CALL Agrif_Set_Updatetype(e3t_id, update = AGRIF_Update_Average)
 
-      IF( ln_zdftke.OR.ln_zdfgls ) THEN
+ !     IF( ln_zdftke.OR.ln_zdfgls ) THEN
 !         CALL Agrif_Set_Updatetype( en_id, update = AGRIF_Update_Average)
 !         CALL Agrif_Set_Updatetype(avt_id, update = AGRIF_Update_Average)
 !         CALL Agrif_Set_Updatetype(avm_id, update = AGRIF_Update_Average)
-      ENDIF
+ !     ENDIF
 
 #endif
@@ -472,7 +602,4 @@
 #if defined key_si3
 SUBROUTINE Agrif_InitValues_cont_ice
-      !!----------------------------------------------------------------------
-      !!                 *** ROUTINE Agrif_InitValues_cont_ice ***
-      !!----------------------------------------------------------------------
       USE Agrif_Util
       USE sbc_oce, ONLY : nn_fsbc  ! clem: necessary otherwise Agrif_Parent(nn_fsbc) = nn_fsbc
@@ -482,11 +609,7 @@
       USE agrif_ice_interp
       USE lib_mpp
-      !
-      IMPLICIT NONE
-      !!----------------------------------------------------------------------
-      !
-      ! Declaration of the type of variable which have to be interpolated (parent=>child)
-      !----------------------------------------------------------------------------------
-      CALL agrif_declare_var_ice
+      !!----------------------------------------------------------------------
+      !!                 *** ROUTINE Agrif_InitValues_cont_ice ***
+      !!----------------------------------------------------------------------
 
       ! Controls
@@ -495,5 +618,5 @@
       !          the run must satisfy CFL=Uice/(dx/dt) < 0.6/nn_fsbc(child)
       !          therefore, if nn_fsbc(child)>1 one must reduce the time-step in proportion to nn_fsbc(child), which is not acceptable
-      !       If a solution is found, the following stop could be removed because the rest of the code take nn_fsbc(child) into account
+      !       If a solution is found, the following stop could be removed because the rest of the code take nn_fsbc(child) into account     
       IF( nn_fsbc > 1 )  CALL ctl_stop('nn_fsbc(child) must be set to 1 otherwise agrif and sea-ice may not work properly')
 
@@ -516,12 +639,13 @@
       !!                 *** ROUTINE agrif_declare_var_ice ***
       !!----------------------------------------------------------------------
+
       USE Agrif_Util
       USE ice
-      USE par_oce, ONLY : nbghostcells
+      USE par_oce, ONLY : nbghostcells, nbghostcells_x, nbghostcells_y_s
       !
       IMPLICIT NONE
       !
       INTEGER :: ind1, ind2, ind3
-      !!----------------------------------------------------------------------
+         !!----------------------------------------------------------------------
       !
       ! 1. Declaration of the type of variable which have to be interpolated (parent=>child)
@@ -532,10 +656,17 @@
       !                            2,2 = two ghost lines
       !-------------------------------------------------------------------------------------
+
       ind1 =     nbghostcells
-      ind2 = 1 + nbghostcells
-      ind3 = 2 + nbghostcells
-      CALL agrif_declare_variable((/2,2,0/),(/ind3,ind3,0/),(/'x','y','N'/),(/1,1,1/),(/nlci,nlcj,jpl*(8+nlay_s+nlay_i)/),tra_ice_id)
-      CALL agrif_declare_variable((/1,2/)  ,(/ind2,ind3/)  ,(/'x','y'/)    ,(/1,1/)  ,(/nlci,nlcj/)                      ,u_ice_id  )
-      CALL agrif_declare_variable((/2,1/)  ,(/ind3,ind2/)  ,(/'x','y'/)    ,(/1,1/)  ,(/nlci,nlcj/)                      ,v_ice_id  )
+      ind2 = 2 + nbghostcells_x
+      ind3 = 2 + nbghostcells_y_s
+      CALL agrif_declare_variable((/2,2,0/),(/ind2,ind3,0/),(/'x','y','N'/),(/1,1,1/),(/nlci,nlcj,jpl*(8+nlay_s+nlay_i)/),tra_ice_id)
+      CALL agrif_declare_variable((/1,2/)  ,(/ind2-1,ind3/)  ,(/'x','y'/)    ,(/1,1/)  ,(/nlci,nlcj/)                      ,u_ice_id  )
+      CALL agrif_declare_variable((/2,1/)  ,(/ind2,ind3-1/)  ,(/'x','y'/)    ,(/1,1/)  ,(/nlci,nlcj/)                      ,v_ice_id  )
+
+      CALL Agrif_Set_MaskMaxSearch(25)
+      CALL agrif_declare_variable((/2,2,0/),(/ind3,ind3,0/),(/'x','y','N'/),(/1,1,1/),(/nlci,nlcj,jpl*(8+nlay_s+nlay_i)/),tra_iceini_id)
+      CALL agrif_declare_variable((/1,2/)  ,(/ind2-1,ind3/)  ,(/'x','y'/)    ,(/1,1/)  ,(/nlci,nlcj/)                      ,u_iceini_id  )
+      CALL agrif_declare_variable((/2,1/)  ,(/ind2,ind3-1/)  ,(/'x','y'/)    ,(/1,1/)  ,(/nlci,nlcj/)                      ,v_iceini_id  )
+      CALL Agrif_Set_MaskMaxSearch(5)
 
       ! 2. Set interpolations (normal & tangent to the grid cell for velocities)
@@ -545,4 +676,11 @@
       CALL Agrif_Set_bcinterp(v_ice_id  , interp1 = AGRIF_ppm   ,interp2 = Agrif_linear)
 
+      CALL Agrif_Set_bcinterp(tra_iceini_id, interp  = AGRIF_linear)
+      CALL Agrif_Set_interp  (tra_iceini_id, interp  = AGRIF_linear)
+      CALL Agrif_Set_bcinterp(u_iceini_id  , interp  = AGRIF_linear  )
+      CALL Agrif_Set_interp  (u_iceini_id  , interp  = AGRIF_linear   )
+      CALL Agrif_Set_bcinterp(v_iceini_id  , interp  = AGRIF_linear)
+      CALL Agrif_Set_interp  (v_iceini_id  , interp  = AGRIF_linear)
+
       ! 3. Set location of interpolations
       !----------------------------------
@@ -550,4 +688,8 @@
       CALL Agrif_Set_bc(u_ice_id  ,(/0,ind1/))
       CALL Agrif_Set_bc(v_ice_id  ,(/0,ind1/))
+
+      CALL Agrif_Set_bc(tra_iceini_id,(/0,ind1/))
+      CALL Agrif_Set_bc(u_iceini_id  ,(/0,ind1/))
+      CALL Agrif_Set_bc(v_iceini_id  ,(/0,ind1/))
 
       ! 4. Set update type in case 2 ways (child=>parent) (normal & tangent to the grid cell for velocities)
@@ -557,9 +699,9 @@
       CALL Agrif_Set_Updatetype(u_ice_id  , update1 = Agrif_Update_Average       , update2 = Agrif_Update_Full_Weighting)
       CALL Agrif_Set_Updatetype(v_ice_id  , update1 = Agrif_Update_Full_Weighting, update2 = Agrif_Update_Average       )
-#else
+# else
       CALL Agrif_Set_Updatetype(tra_ice_id, update  = AGRIF_Update_Average)
       CALL Agrif_Set_Updatetype(u_ice_id  , update1 = Agrif_Update_Copy   , update2 = Agrif_Update_Average)
       CALL Agrif_Set_Updatetype(v_ice_id  , update1 = Agrif_Update_Average, update2 = Agrif_Update_Copy   )
-#endif
+# endif
 
    END SUBROUTINE agrif_declare_var_ice
@@ -585,36 +727,39 @@
       USE agrif_top_sponge
       !!
-      IMPLICIT NONE
-      !
-      CHARACTER(len=10) :: cl_check1, cl_check2, cl_check3
-      LOGICAL :: check_namelist
-      !!----------------------------------------------------------------------
-
-      ! 1. Declaration of the type of variable which have to be interpolated
-      !---------------------------------------------------------------------
-      CALL agrif_declare_var_top
-
-      ! 2. First interpolations of potentially non zero fields
-      !-------------------------------------------------------
-      Agrif_SpecialValue=0._wp
-      Agrif_UseSpecialValue = .TRUE.
-      CALL Agrif_Bc_variable(trn_id,calledweight=1.,procname=interptrn)
-      Agrif_UseSpecialValue = .FALSE.
-      CALL Agrif_Sponge
-      tabspongedone_trn = .FALSE.
-      CALL Agrif_Bc_variable(trn_sponge_id,calledweight=1.,procname=interptrn_sponge)
-      ! reset ts(:,:,:,:,Krhs_a) to zero
-      tr(:,:,:,:,Krhs_a) = 0._wp
-
-      ! 3. Some controls
-      !-----------------
-      check_namelist = .TRUE.
-
-      IF( check_namelist ) THEN
-         ! Check time steps
-      IF( NINT(Agrif_Rhot()) * NINT(rn_Dt) .NE. Agrif_Parent(rn_Dt) ) THEN
-         WRITE(cl_check1,*)  Agrif_Parent(rn_Dt)
-         WRITE(cl_check2,*)  rn_Dt
-         WRITE(cl_check3,*)  rn_Dt*Agrif_Rhot()
+  
+   !!
+   IMPLICIT NONE
+   !
+   CHARACTER(len=10) :: cl_check1, cl_check2, cl_check3
+   LOGICAL :: check_namelist
+      !!----------------------------------------------------------------------
+
+
+   ! 1. Declaration of the type of variable which have to be interpolated
+   !---------------------------------------------------------------------
+   CALL agrif_declare_var_top
+
+   ! 2. First interpolations of potentially non zero fields
+   !-------------------------------------------------------
+   Agrif_SpecialValue=0.
+   Agrif_UseSpecialValue = .TRUE.
+   CALL Agrif_Bc_variable(trn_id,calledweight=1.,procname=interptrn)
+   Agrif_UseSpecialValue = .FALSE.
+   CALL Agrif_Sponge
+   tabspongedone_trn = .FALSE.
+   CALL Agrif_Bc_variable(trn_sponge_id,calledweight=1.,procname=interptrn_sponge)
+   ! reset tsa to zero
+   tra(:,:,:,:) = 0.
+
+   ! 3. Some controls
+   !-----------------
+   check_namelist = .TRUE.
+
+   IF( check_namelist ) THEN
+      ! Check time steps
+      IF( NINT(Agrif_Rhot()) * NINT(rdt) .NE. Agrif_Parent(rdt) ) THEN
+         WRITE(cl_check1,*)  Agrif_Parent(rdt)
+         WRITE(cl_check2,*)  rdt
+         WRITE(cl_check3,*)  rdt*Agrif_Rhot()
          CALL ctl_stop( 'incompatible time step between grids',   &
                &               'parent grid value : '//cl_check1    ,   & 
@@ -635,5 +780,4 @@
          nitend =  Agrif_Parent(nitend)   *Agrif_IRhot()
       ENDIF
-
    ENDIF
    !
@@ -655,15 +799,19 @@
       !!----------------------------------------------------------------------
 
+
+
+!RB_CMEMS : declare here init for top      
       ! 1. Declaration of the type of variable which have to be interpolated
       !---------------------------------------------------------------------
       ind1 =     nbghostcells
-      ind2 = 1 + nbghostcells
-      ind3 = 2 + nbghostcells
+      ind2 = 2 + nbghostcells_x
+      ind3 = 2 + nbghostcells_y_s
 # if defined key_vertical
-      CALL agrif_declare_variable((/2,2,0,0/),(/ind3,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,jptra+1/),trn_id)
-      CALL agrif_declare_variable((/2,2,0,0/),(/ind3,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,jptra+1/),trn_sponge_id)
+      CALL agrif_declare_variable((/2,2,0,0/),(/ind2,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,jptra+1/),trn_id)
+      CALL agrif_declare_variable((/2,2,0,0/),(/ind2,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,jptra+1/),trn_sponge_id)
 # else
-      CALL agrif_declare_variable((/2,2,0,0/),(/ind3,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,jptra/),trn_id)
-      CALL agrif_declare_variable((/2,2,0,0/),(/ind3,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,jptra/),trn_sponge_id)
+! LAURENT: STRANGE why (3,3) here ?
+      CALL agrif_declare_variable((/2,2,0,0/),(/3,3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,jptra/),trn_id)
+      CALL agrif_declare_variable((/2,2,0,0/),(/3,3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,jptra/),trn_sponge_id)
 # endif
 
@@ -705,8 +853,9 @@
       !!                     *** ROUTINE agrif_init ***
       !!----------------------------------------------------------------------
-      USE agrif_oce 
-      USE agrif_ice
-      USE in_out_manager
-      USE lib_mpp
+   USE agrif_oce 
+   USE agrif_ice
+   USE dom_oce
+   USE in_out_manager
+   USE lib_mpp
       !!
       IMPLICIT NONE
@@ -714,5 +863,5 @@
       INTEGER  ::   ios                 ! Local integer output status for namelist read
       NAMELIST/namagrif/ ln_agrif_2way, rn_sponge_tra, rn_sponge_dyn, rn_trelax_tra, rn_trelax_dyn, &
-                       & ln_spc_dyn, ln_chk_bathy
+                       & ln_spc_dyn, ln_chk_bathy, ln_bry_south
       !!--------------------------------------------------------------------------------------
       !
@@ -735,13 +884,25 @@
          WRITE(numout,*) '      use special values for dynamics   ln_spc_dyn    = ', ln_spc_dyn
          WRITE(numout,*) '      check bathymetry                  ln_chk_bathy  = ', ln_chk_bathy
-      ENDIF
-      !
-      !
-      IF( agrif_oce_alloc()  > 0 )   CALL ctl_warn('agrif agrif_oce_alloc: allocation of arrays failed')
+         WRITE(numout,*) '      south boundary                    ln_bry_south  = ', ln_bry_south
+      ENDIF
+      !
+      ! Set the number of ghost cells according to periodicity
+      nbghostcells_x = nbghostcells
+      nbghostcells_y_s = nbghostcells
+      nbghostcells_y_n = nbghostcells
+      !
+      IF ( jperio == 1 ) nbghostcells_x = 0
+      IF ( .NOT. ln_bry_south ) nbghostcells_y_s = 0
+
+      ! Some checks
+      IF( jpiglo /= nbcellsx + 2 + 2*nbghostcells_x )   &
+          CALL ctl_stop( 'STOP', 'agrif_nemo_init: Agrif children requires jpiglo == nbcellsx + 2 + 2*nbghostcells_x' )
+      IF( jpjglo /= nbcellsy + 2 + nbghostcells_y_s + nbghostcells_y_n )   &
+          CALL ctl_stop( 'STOP', 'agrif_nemo_init: Agrif children requires jpjglo == nbcellsy + 2 + nbghostcells_y_s + nbghostcells_y_n' )
+      IF( ln_use_jattr )   CALL ctl_stop( 'STOP', 'agrif_nemo_init:Agrif children requires ln_use_jattr = .false. ' )
       !
    END SUBROUTINE agrif_nemo_init
 
 # if defined key_mpp_mpi
-
    SUBROUTINE Agrif_InvLoc( indloc, nprocloc, i, indglob )
       !!----------------------------------------------------------------------
@@ -803,4 +964,231 @@
 # endif
 
+   SUBROUTINE nemo_mapping(ndim,ptx,pty,bounds,bounds_chunks,correction_required,nb_chunks)
+      !!----------------------------------------------------------------------
+      !!                   *** ROUTINE Nemo_mapping ***
+      !!----------------------------------------------------------------------
+      USE dom_oce
+      !!
+      IMPLICIT NONE
+      !
+      INTEGER :: ndim
+      INTEGER :: ptx, pty
+      INTEGER, DIMENSION(ndim,2,2) :: bounds
+      INTEGER, DIMENSION(:,:,:,:), ALLOCATABLE :: bounds_chunks
+      LOGICAL, DIMENSION(:), ALLOCATABLE :: correction_required
+      INTEGER :: nb_chunks
+      !
+      INTEGER :: i
+
+      IF (agrif_debug_interp) THEN
+         DO i=1,ndim
+            WRITE(*,*) 'direction = ',i,bounds(i,1,2),bounds(i,2,2)
+         ENDDO
+      ENDIF
+
+      IF( bounds(2,2,2) > jpjglo) THEN
+         IF( bounds(2,1,2) <=jpjglo) THEN
+            nb_chunks = 2
+            ALLOCATE(bounds_chunks(nb_chunks,ndim,2,2))
+            ALLOCATE(correction_required(nb_chunks))
+            DO i = 1,nb_chunks
+               bounds_chunks(i,:,:,:) = bounds
+            END DO
+        
+      ! FIRST CHUNCK (for j<=jpjglo)   
+
+      ! Original indices
+            bounds_chunks(1,1,1,1) = bounds(1,1,2)
+            bounds_chunks(1,1,2,1) = bounds(1,2,2)
+            bounds_chunks(1,2,1,1) = bounds(2,1,2)
+            bounds_chunks(1,2,2,1) = jpjglo
+
+            bounds_chunks(1,1,1,2) = bounds(1,1,2)
+            bounds_chunks(1,1,2,2) = bounds(1,2,2)
+            bounds_chunks(1,2,1,2) = bounds(2,1,2)
+            bounds_chunks(1,2,2,2) = jpjglo
+
+      ! Correction required or not
+            correction_required(1)=.FALSE.
+       
+      ! SECOND CHUNCK (for j>jpjglo)
+
+      ! Original indices
+            bounds_chunks(2,1,1,1) = bounds(1,1,2)
+            bounds_chunks(2,1,2,1) = bounds(1,2,2)
+            bounds_chunks(2,2,1,1) = jpjglo-2
+            bounds_chunks(2,2,2,1) = bounds(2,2,2)
+
+      ! Where to find them
+      ! We use the relation TAB(ji,jj)=TAB(jpiglo-ji+2,jpjglo-2-(jj-jpjglo))
+
+            IF( ptx == 2) THEN ! T, V points
+               bounds_chunks(2,1,1,2) = jpiglo-bounds(1,2,2)+2
+               bounds_chunks(2,1,2,2) = jpiglo-bounds(1,1,2)+2
+            ELSE ! U, F points
+               bounds_chunks(2,1,1,2) = jpiglo-bounds(1,2,2)+1
+               bounds_chunks(2,1,2,2) = jpiglo-bounds(1,1,2)+1       
+            ENDIF
+
+            IF( pty == 2) THEN ! T, U points
+               bounds_chunks(2,2,1,2) = jpjglo-2-(bounds(2,2,2) -jpjglo)
+               bounds_chunks(2,2,2,2) = jpjglo-2-(jpjglo-2      -jpjglo)
+            ELSE ! V, F points
+               bounds_chunks(2,2,1,2) = jpjglo-3-(bounds(2,2,2) -jpjglo)
+               bounds_chunks(2,2,2,2) = jpjglo-3-(jpjglo-2      -jpjglo)
+            ENDIF
+      ! Correction required or not
+            correction_required(2)=.TRUE.
+
+         ELSE
+            nb_chunks = 1
+            ALLOCATE(bounds_chunks(nb_chunks,ndim,2,2))
+            ALLOCATE(correction_required(nb_chunks))
+            DO i=1,nb_chunks
+               bounds_chunks(i,:,:,:) = bounds
+            END DO
+
+            bounds_chunks(1,1,1,1) = bounds(1,1,2)
+            bounds_chunks(1,1,2,1) = bounds(1,2,2)
+            bounds_chunks(1,2,1,1) = bounds(2,1,2)
+            bounds_chunks(1,2,2,1) = bounds(2,2,2)
+
+            bounds_chunks(1,1,1,2) = jpiglo-bounds(1,2,2)+2
+            bounds_chunks(1,1,2,2) = jpiglo-bounds(1,1,2)+2
+
+            bounds_chunks(1,2,1,2) = jpjglo-2-(bounds(2,2,2)-jpjglo)
+            bounds_chunks(1,2,2,2) = jpjglo-2-(bounds(2,1,2)-jpjglo)
+
+            IF( ptx == 2) THEN ! T, V points
+               bounds_chunks(1,1,1,2) = jpiglo-bounds(1,2,2)+2
+               bounds_chunks(1,1,2,2) = jpiglo-bounds(1,1,2)+2
+            ELSE ! U, F points
+               bounds_chunks(1,1,1,2) = jpiglo-bounds(1,2,2)+1
+               bounds_chunks(1,1,2,2) = jpiglo-bounds(1,1,2)+1    	
+            ENDIF
+
+            IF (pty == 2) THEN ! T, U points
+               bounds_chunks(1,2,1,2) = jpjglo-2-(bounds(2,2,2) -jpjglo)
+               bounds_chunks(1,2,2,2) = jpjglo-2-(bounds(2,1,2) -jpjglo)
+            ELSE ! V, F points
+               bounds_chunks(1,2,1,2) = jpjglo-3-(bounds(2,2,2) -jpjglo)
+               bounds_chunks(1,2,2,2) = jpjglo-3-(bounds(2,1,2) -jpjglo)
+            ENDIF
+
+            correction_required(1)=.TRUE.          
+         ENDIF
+
+      ELSE IF (bounds(1,1,2) < 1) THEN
+         IF (bounds(1,2,2) > 0) THEN
+            nb_chunks = 2
+            ALLOCATE(correction_required(nb_chunks))
+            correction_required=.FALSE.
+            ALLOCATE(bounds_chunks(nb_chunks,ndim,2,2))
+            DO i=1,nb_chunks
+               bounds_chunks(i,:,:,:) = bounds
+            END DO
+              
+            bounds_chunks(1,1,1,2) = bounds(1,1,2)+jpiglo-2
+            bounds_chunks(1,1,2,2) = 1+jpiglo-2
+          
+            bounds_chunks(1,1,1,1) = bounds(1,1,2)
+            bounds_chunks(1,1,2,1) = 1
+       
+            bounds_chunks(2,1,1,2) = 2
+            bounds_chunks(2,1,2,2) = bounds(1,2,2)
+          
+            bounds_chunks(2,1,1,1) = 2
+            bounds_chunks(2,1,2,1) = bounds(1,2,2)
+
+         ELSE
+            nb_chunks = 1
+            ALLOCATE(correction_required(nb_chunks))
+            correction_required=.FALSE.
+            ALLOCATE(bounds_chunks(nb_chunks,ndim,2,2))
+            DO i=1,nb_chunks
+               bounds_chunks(i,:,:,:) = bounds
+            END DO    
+            bounds_chunks(1,1,1,2) = bounds(1,1,2)+jpiglo-2
+            bounds_chunks(1,1,2,2) = bounds(1,2,2)+jpiglo-2
+          
+            bounds_chunks(1,1,1,1) = bounds(1,1,2)
+           bounds_chunks(1,1,2,1) = bounds(1,2,2)
+         ENDIF
+      ELSE
+         nb_chunks=1  
+         ALLOCATE(correction_required(nb_chunks))
+         correction_required=.FALSE.
+         ALLOCATE(bounds_chunks(nb_chunks,ndim,2,2))
+         DO i=1,nb_chunks
+            bounds_chunks(i,:,:,:) = bounds
+         END DO
+         bounds_chunks(1,1,1,2) = bounds(1,1,2)
+         bounds_chunks(1,1,2,2) = bounds(1,2,2)
+         bounds_chunks(1,2,1,2) = bounds(2,1,2)
+         bounds_chunks(1,2,2,2) = bounds(2,2,2)
+          
+         bounds_chunks(1,1,1,1) = bounds(1,1,2)
+         bounds_chunks(1,1,2,1) = bounds(1,2,2)
+         bounds_chunks(1,2,1,1) = bounds(2,1,2)
+         bounds_chunks(1,2,2,1) = bounds(2,2,2)              
+      ENDIF
+        
+   END SUBROUTINE nemo_mapping
+
+   FUNCTION agrif_external_switch_index(ptx,pty,i1,isens)
+
+   USE dom_oce
+
+   INTEGER :: ptx, pty, i1, isens
+   INTEGER :: agrif_external_switch_index
+
+   IF( isens == 1 ) THEN
+      IF( ptx == 2 ) THEN ! T, V points
+         agrif_external_switch_index = jpiglo-i1+2
+      ELSE ! U, F points
+         agrif_external_switch_index = jpiglo-i1+1      
+      ENDIF
+   ELSE IF( isens ==2 ) THEN
+      IF ( pty == 2 ) THEN ! T, U points
+         agrif_external_switch_index = jpjglo-2-(i1 -jpjglo)
+      ELSE ! V, F points
+         agrif_external_switch_index = jpjglo-3-(i1 -jpjglo)
+      ENDIF
+   ENDIF
+
+   END function agrif_external_switch_index
+
+   SUBROUTINE Correct_field(tab2d,i1,i2,j1,j2)
+      !!----------------------------------------------------------------------
+      !!                   *** ROUTINE Correct_field ***
+      !!----------------------------------------------------------------------
+   
+   USE dom_oce
+   USE agrif_oce
+
+   INTEGER :: i1,i2,j1,j2
+   REAL(wp), DIMENSION(i1:i2,j1:j2) :: tab2d
+
+   INTEGER :: i,j
+   REAL(wp), DIMENSION(i1:i2,j1:j2) :: tab2dtemp
+
+   tab2dtemp = tab2d
+
+   IF( .NOT. use_sign_north ) THEN
+      DO j=j1,j2
+         DO i=i1,i2
+            tab2d(i,j)=tab2dtemp(i2-(i-i1),j2-(j-j1))
+         END DO
+      END DO
+   ELSE
+      DO j=j1,j2
+         DO i=i1,i2
+            tab2d(i,j)=sign_north * tab2dtemp(i2-(i-i1),j2-(j-j1))
+         END DO
+      END DO
+   ENDIF
+
+   END SUBROUTINE Correct_field
+
 #else
    SUBROUTINE Subcalledbyagrif
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/DOM/dom_oce.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/DOM/dom_oce.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/DOM/dom_oce.F90	(revision 13026)
@@ -215,4 +215,5 @@
 #if defined key_agrif
    LOGICAL, PUBLIC, PARAMETER ::   lk_agrif = .TRUE.    !: agrif flag
+   LOGICAL, PUBLIC            ::   lk_south, lk_north, lk_west, lk_east !: Child grid boundaries (interpolation or not)
 #else
    LOGICAL, PUBLIC, PARAMETER ::   lk_agrif = .FALSE.   !: agrif flag
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/DOM/domain.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/DOM/domain.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/DOM/domain.F90	(revision 13026)
@@ -187,6 +187,11 @@
       ENDIF
       !
+
       IF( lk_c1d         )   CALL cor_c1d       ! 1D configuration: Coriolis set at T-point
       !
+
+#if defined key_agrif
+      IF( .NOT. Agrif_Root() ) CALL Agrif_Init_Domain( Kbb, Kmm, Kaa )
+#endif
       IF( ln_meshmask    )   CALL dom_wri       ! Create a domain file
       IF( .NOT.ln_rstart )   CALL dom_ctl       ! Domain control
@@ -307,4 +312,11 @@
 902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namrun in configuration namelist' )
       IF(lwm) WRITE ( numond, namrun )
+
+#if defined key_agrif
+      IF( .NOT. Agrif_Root() ) THEN
+            nn_it000 = (Agrif_Parent(nn_it000)-1)*Agrif_IRhot() + 1
+            nn_itend =  Agrif_Parent(nn_itend)   *Agrif_IRhot()
+      ENDIF
+#endif
       !
       IF(lwp) THEN                  ! control print
@@ -403,4 +415,10 @@
       IF(lwm) WRITE( numond, namdom )
       !
+#if defined key_agrif
+      IF( .NOT. Agrif_Root() ) THEN
+            rn_Dt = Agrif_Parent(rn_Dt) / Agrif_Rhot()
+      ENDIF
+#endif
+      !
       IF(lwp) THEN
          WRITE(numout,*)
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/DOM/istate.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/DOM/istate.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/DOM/istate.F90	(revision 13026)
@@ -35,4 +35,8 @@
    USE lib_mpp         ! MPP library
    USE restart         ! restart
+#if defined key_agrif
+   USE agrif_oce_interp
+   USE agrif_oce
+#endif   
 
    IMPLICIT NONE
@@ -70,5 +74,5 @@
 !!gm  Why not include in the first call of dta_tsd ?  
 !!gm  probably associated with the use of internal damping...
-                     CALL dta_tsd_init        ! Initialisation of T & S input data
+       CALL dta_tsd_init        ! Initialisation of T & S input data
 !!gm to be moved in usrdef of C1D case
 !      IF( lk_c1d )   CALL dta_uvd_init        ! Initialization of U & V input data
@@ -97,7 +101,15 @@
          !
          IF( ln_tsd_init ) THEN               
-            CALL dta_tsd( nit000, ts(:,:,:,:,Kbb) )       ! read 3D T and S data at nit000
+            IF( Agrif_root() ) THEN
+               CALL dta_tsd( nit000, ts(:,:,:,:,Kbb) )       ! read 3D T and S data at nit000
+               ssh(:,:,Kbb)   = 0._wp               ! set the ocean at rest
+               uu  (:,:,:,Kbb) = 0._wp
+               vv  (:,:,:,Kbb) = 0._wp  
+#if defined key_agrif
+            ELSE
+               CALL agrif_istate( Kbb, Kmm, Kaa ) 
+#endif
+            ENDIF
             !
-            ssh(:,:,Kbb)   = 0._wp               ! set the ocean at rest
             IF( ll_wd ) THEN
                ssh(:,:,Kbb) =  -ssh_ref  ! Added in 30 here for bathy that adds 30 as Iterative test CEOD 
@@ -111,7 +123,5 @@
                END_2D
             ENDIF 
-            uu  (:,:,:,Kbb) = 0._wp
-            vv  (:,:,:,Kbb) = 0._wp  
-            !
+             !
          ELSE                                 ! user defined initial T and S
             CALL usr_def_istate( gdept(:,:,:,Kbb), tmask, ts(:,:,:,:,Kbb), uu(:,:,:,Kbb), vv(:,:,:,Kbb), ssh(:,:,Kbb)  )         
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/DYN/sshwzv.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/DYN/sshwzv.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/DYN/sshwzv.F90	(revision 13026)
@@ -204,31 +204,40 @@
          ! Mask vertical velocity at first/last columns/row 
          ! inside computational domain (cosmetic) 
-         ! --- West --- !
-         DO ji = mi0(2), mi1(2)
-            DO jj = 1, jpj
-               pww(ji,jj,:) = 0._wp 
+         ! --- West --- !         
+         IF( lk_west) THEN
+            DO ji = mi0(2), mi1(2)
+               DO jj = 1, jpj
+                  pww(ji,jj,:) = 0._wp 
+               ENDDO
             ENDDO
-         ENDDO
+         ENDIF
          !
          ! --- East --- !
-         DO ji = mi0(jpiglo-1), mi1(jpiglo-1)
-            DO jj = 1, jpj
-               pww(ji,jj,:) = 0._wp
+         IF( lk_east) THEN
+            DO ji = mi0(jpiglo-1), mi1(jpiglo-1)
+               DO jj = 1, jpj
+                  pww(ji,jj,:) = 0._wp
+               ENDDO
             ENDDO
-         ENDDO
+         ENDIF
          !
          ! --- South --- !
-         DO jj = mj0(2), mj1(2)
-            DO ji = 1, jpi
-               pww(ji,jj,:) = 0._wp
+         IF( lk_south) THEN
+            DO jj = mj0(2), mj1(2)
+               DO ji = 1, jpi
+                  pww(ji,jj,:) = 0._wp
+               ENDDO
             ENDDO
-         ENDDO
+         ENDIF
          !
          ! --- North --- !
-         DO jj = mj0(jpjglo-1), mj1(jpjglo-1)
-            DO ji = 1, jpi
-               pww(ji,jj,:) = 0._wp
+         IF( lk_north) THEN
+            DO jj = mj0(jpjglo-1), mj1(jpjglo-1)
+               DO ji = 1, jpi
+                  pww(ji,jj,:) = 0._wp
+               ENDDO
             ENDDO
-         ENDDO
+         ENDIF
+         !
       ENDIF 
 #endif 
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/FLO/floblk.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/FLO/floblk.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/FLO/floblk.F90	(revision 13026)
@@ -41,4 +41,8 @@
       INTEGER, INTENT( in  ) ::   Kbb, Kmm ! ocean time level indices
       !!
+#ifndef key_agrif
+
+!RB super quick fix to compile with agrif
+
       INTEGER :: jfl              ! dummy loop arguments
       INTEGER :: ind, ifin, iloop
@@ -364,4 +368,5 @@
          GO TO 222
       ENDIF
+#endif
       !
       !
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/LBC/lib_mpp.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/LBC/lib_mpp.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/LBC/lib_mpp.F90	(revision 13026)
@@ -137,7 +137,7 @@
 
    ! Communications summary report
-   CHARACTER(len=128), DIMENSION(:), ALLOCATABLE ::   crname_lbc                   !: names of lbc_lnk calling routines
-   CHARACTER(len=128), DIMENSION(:), ALLOCATABLE ::   crname_glb                   !: names of global comm calling routines
-   CHARACTER(len=128), DIMENSION(:), ALLOCATABLE ::   crname_dlg                   !: names of delayed global comm calling routines
+   CHARACTER(len=lca), DIMENSION(:), ALLOCATABLE ::   crname_lbc                   !: names of lbc_lnk calling routines
+   CHARACTER(len=lca), DIMENSION(:), ALLOCATABLE ::   crname_glb                   !: names of global comm calling routines
+   CHARACTER(len=lca), DIMENSION(:), ALLOCATABLE ::   crname_dlg                   !: names of delayed global comm calling routines
    INTEGER, PUBLIC                               ::   ncom_stp = 0                 !: copy of time step # istp
    INTEGER, PUBLIC                               ::   ncom_fsbc = 1                !: copy of sbc time step # nn_fsbc
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/LBC/mppini.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/LBC/mppini.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/LBC/mppini.F90	(revision 13026)
@@ -102,4 +102,9 @@
             &           'the domain is lay out for distributed memory computing!' )
          !
+#if defined key_agrif
+    IF (.NOT.agrif_root()) THEN
+      call agrif_nemo_init()
+    ENDIF
+#endif
    END SUBROUTINE mpp_init
 
@@ -333,9 +338,5 @@
 #if defined key_agrif
       IF( .NOT. Agrif_Root() ) THEN       ! AGRIF children: specific setting (cf. agrif_user.F90)
-         IF( jpiglo /= nbcellsx + 2 + 2*nbghostcells )   &
-            CALL ctl_stop( 'STOP', 'mpp_init: Agrif children requires jpiglo == nbcellsx + 2 + 2*nbghostcells' )
-         IF( jpjglo /= nbcellsy + 2 + 2*nbghostcells )   &
-            CALL ctl_stop( 'STOP', 'mpp_init: Agrif children requires jpjglo == nbcellsy + 2 + 2*nbghostcells' )
-         IF( ln_use_jattr )   CALL ctl_stop( 'STOP', 'mpp_init:Agrif children requires ln_use_jattr = .false. ' )
+         CALL agrif_nemo_init()
       ENDIF
 #endif
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/OBS/diaobs.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/OBS/diaobs.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/OBS/diaobs.F90	(revision 13026)
@@ -94,5 +94,5 @@
    TYPE(obs_prof), PUBLIC, POINTER, DIMENSION(:) ::   profdataqc   !: Profile data after quality control
 
-   CHARACTER(len=6), PUBLIC, DIMENSION(:), ALLOCATABLE ::   cobstypesprof, cobstypessurf   !: Profile & surface obs types
+   CHARACTER(len=lca), PUBLIC, DIMENSION(:), ALLOCATABLE ::   cobstypesprof, cobstypessurf   !: Profile & surface obs types
 
    !!----------------------------------------------------------------------
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/SBC/sbcmod.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/SBC/sbcmod.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/SBC/sbcmod.F90	(revision 13026)
@@ -122,8 +122,8 @@
       !                             !* overwrite namelist parameter using CPP key information
 #if defined key_agrif
-      IF( Agrif_Root() ) THEN                ! AGRIF zoom (cf r1242: possibility to run without ice in fine grid)
-         IF( lk_si3  )   nn_ice      = 2
-         IF( lk_cice )   nn_ice      = 3
-      ENDIF
+ !     IF( Agrif_Root() ) THEN                ! AGRIF zoom (cf r1242: possibility to run without ice in fine grid)
+ !        IF( lk_si3  )   nn_ice      = 2
+ !        IF( lk_cice )   nn_ice      = 3
+ !     ENDIF
 !!GS: TBD
 !#else
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/STO/stopar.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/STO/stopar.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/STO/stopar.F90	(revision 13026)
@@ -56,6 +56,6 @@
    INTEGER,          DIMENSION(:),       ALLOCATABLE :: sto3d_ord  ! order of autoregressive process
 
-   CHARACTER(len=1), DIMENSION(:),       ALLOCATABLE :: sto2d_typ  ! nature of grid point (T, U, V, W, F, I)
-   CHARACTER(len=1), DIMENSION(:),       ALLOCATABLE :: sto3d_typ  ! nature of grid point (T, U, V, W, F, I)
+   CHARACTER(len=lca), DIMENSION(:),       ALLOCATABLE :: sto2d_typ  ! nature of grid point (T, U, V, W, F, I)
+   CHARACTER(len=lca), DIMENSION(:),       ALLOCATABLE :: sto3d_typ  ! nature of grid point (T, U, V, W, F, I)
    REAL(wp),         DIMENSION(:),       ALLOCATABLE :: sto2d_sgn  ! control of the sign accross the north fold
    REAL(wp),         DIMENSION(:),       ALLOCATABLE :: sto3d_sgn  ! control of the sign accross the north fold
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/USR/usrdef_hgr.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/USR/usrdef_hgr.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/USR/usrdef_hgr.F90	(revision 13026)
@@ -95,4 +95,5 @@
 #if defined key_agrif
       ! ! Upper left longitude and latitude from parent:
+      ! Laurent: Should be modify in case of an east-west cyclic parent grid
       IF (.NOT.Agrif_root()) THEN
          zlam0 = zlam1 + Agrif_irhox() * REAL(Agrif_Parent(jpjglo)-2 , wp) * ze1deg * zcos_alpha  &
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/USR/usrdef_nam.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/USR/usrdef_nam.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/USR/usrdef_nam.F90	(revision 13026)
@@ -74,6 +74,6 @@
 #if defined key_agrif
       IF( .NOT. Agrif_Root() ) THEN
-         kpi  = nbcellsx + 2 + 2*nbghostcells
-         kpj  = nbcellsy + 2 + 2*nbghostcells
+         kpi  = nbcellsx + 2 + 2*nbghostcells_x
+         kpj  = nbcellsy + 2 + 2*nbghostcells_y_s
       ENDIF
 #endif
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/nemogcm.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/nemogcm.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/nemogcm.F90	(revision 13026)
@@ -143,11 +143,7 @@
 #if defined key_agrif
       Kbb_a = Nbb; Kmm_a = Nnn; Krhs_a = Nrhs   ! agrif_oce module copies of time level indices
-      CALL Agrif_Declare_Var_dom   ! AGRIF: set the meshes for DOM
-      CALL Agrif_Declare_Var       !  "      "   "   "      "  DYN/TRA 
+      CALL Agrif_Declare_Var       !  "      "   "   "      "  DYN/TRA
 # if defined key_top
       CALL Agrif_Declare_Var_top   !  "      "   "   "      "  TOP
-# endif
-# if defined key_si3
-      CALL Agrif_Declare_Var_ice   !  "      "   "   "      "  Sea ice
 # endif
 #endif
@@ -399,5 +395,7 @@
       ! Initialise time level indices
       Nbb = 1; Nnn = 2; Naa = 3; Nrhs = Naa
-
+#if defined key_agrif
+      Kbb_a = Nbb; Kmm_a = Nnn; Krhs_a = Nrhs   ! agrif_oce module copies of time level indices
+#endif 
       !                             !-------------------------------!
       !                             !  NEMO general initialization  !
@@ -414,5 +412,12 @@
       IF( lk_c1d       )   CALL     c1d_init        ! 1D column configuration
                            CALL     wad_init        ! Wetting and drying options
+
+#if defined key_agrif
+     CALL Agrif_Declare_Var_ini   !  "      "   "   "      "  DOM
+#endif
                            CALL     dom_init( Nbb, Nnn, Naa, "OPA") ! Domain
+
+
+
       IF( ln_crs       )   CALL     crs_init(      Nnn )       ! coarsened grid: domain initialization 
       IF( sn_cfctl%l_prtctl )   &
@@ -435,5 +440,5 @@
       ENDIF
       !
-      
+
                            CALL  istate_init( Nbb, Nnn, Naa )    ! ocean initial state (Dynamics and tracers)
 
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/par_kind.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/par_kind.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/par_kind.F90	(revision 13026)
@@ -31,5 +31,6 @@
    
    !                                                                !!** Integer **
-   INTEGER, PUBLIC, PARAMETER ::   lc = 256                          !: Lenght of Character strings
+   INTEGER, PUBLIC, PARAMETER ::   lc  = 256                          !: Lenght of Character strings
+   INTEGER, PUBLIC, PARAMETER ::   lca = 400                          !: Lenght of Character arrays
 
    !!----------------------------------------------------------------------
Index: NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/par_oce.F90
===================================================================
--- NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/par_oce.F90	(revision 13023)
+++ NEMO/branches/2020/dev_r12973_AGRIF_CMEMS/src/OCE/par_oce.F90	(revision 13026)
@@ -47,7 +47,10 @@
    ! global domain size for AGRIF     !!! * total AGRIF computational domain *
    INTEGER, PUBLIC            ::   nbug_in_agrif_conv_do_not_remove_or_modify = 1 - 1
-   INTEGER, PUBLIC, PARAMETER ::   nbghostcells = 3                             !: number of ghost cells
-   INTEGER, PUBLIC            ::   nbcellsx   ! = jpiglo - 2 - 2*nbghostcells   !: number of cells in i-direction
-   INTEGER, PUBLIC            ::   nbcellsy   ! = jpjglo - 2 - 2*nbghostcells   !: number of cells in j-direction
+   INTEGER, PUBLIC, PARAMETER ::   nbghostcells = 3 !: number of ghost cells: default value
+   INTEGER, PUBLIC            ::   nbghostcells_x   !: number of ghost cells in i-direction
+   INTEGER, PUBLIC            ::   nbghostcells_y_s   !: number of ghost cells in j-direction at south
+   INTEGER, PUBLIC            ::   nbghostcells_y_n   !: number of ghost cells in j-direction at north                       !: number of ghost cells
+   INTEGER, PUBLIC            ::   nbcellsx   ! = jpiglo - 2 - 2*nbghostcells_x   !: number of cells in i-direction
+   INTEGER, PUBLIC            ::   nbcellsy   ! = jpjglo - 2 - 2*nbghostcells-y   !: number of cells in j-direction
 
    ! local domain size                !!! * local computational domain *
