New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 13216 – NEMO

Changeset 13216


Ignore:
Timestamp:
2020-07-02T11:25:49+02:00 (4 years ago)
Author:
rblod
Message:

Merge dev_r12973_AGRIF_CMEMS

Location:
NEMO/trunk
Files:
35 edited

Legend:

Unmodified
Added
Removed
  • NEMO/trunk

    • Property svn:externals
      •  

        old new  
        22^/utils/build/makenemo@HEAD   makenemo 
        33^/utils/build/mk@HEAD         mk 
        4 ^/utils/tools@HEAD            tools 
        5 ^/vendors/AGRIF/dev@HEAD      ext/AGRIF 
         4^/utils/tools/@HEAD           tools 
         5^/vendors/AGRIF/dev_r12970_AGRIF_CMEMS      ext/AGRIF 
        66^/vendors/FCM@HEAD            ext/FCM 
        77^/vendors/IOIPSL@HEAD         ext/IOIPSL 
  • NEMO/trunk/cfgs/AGRIF_DEMO/EXPREF/2_namelist_cfg

    r13208 r13216  
    158158!----------------------------------------------------------------------- 
    159159   ln_spc_dyn    = .true.  !  use 0 as special value for dynamics 
    160    rn_sponge_tra = 1440.   !  coefficient for tracer   sponge layer [m2/s] 
    161    rn_sponge_dyn = 1440.   !  coefficient for dynamics sponge layer [m2/s] 
    162160   ln_chk_bathy  = .true.  !  =T  check the parent bathymetry 
    163161/ 
  • NEMO/trunk/cfgs/AGRIF_DEMO/EXPREF/3_namelist_cfg

    r12994 r13216  
    158158!----------------------------------------------------------------------- 
    159159   ln_spc_dyn    = .true.  !  use 0 as special value for dynamics 
    160    rn_sponge_tra =  480.   !  coefficient for tracer   sponge layer [m2/s] 
    161    rn_sponge_dyn =  480.   !  coefficient for dynamics sponge layer [m2/s] 
    162160   ln_chk_bathy  = .true.  !  =T  check the parent bathymetry 
    163161/ 
  • NEMO/trunk/cfgs/SHARED/namelist_ref

    r13214 r13216  
    643643&namagrif      !  AGRIF zoom                                            ("key_agrif") 
    644644!----------------------------------------------------------------------- 
    645    ln_agrif_2way = .true.  !  activate two way nesting 
    646    ln_spc_dyn    = .true.  !  use 0 as special value for dynamics 
    647    rn_sponge_tra = 2880.   !  coefficient for tracer   sponge layer [m2/s] 
    648    rn_sponge_dyn = 2880.   !  coefficient for dynamics sponge layer [m2/s] 
    649    rn_trelax_tra = 0.01    !  inverse of relaxation time (in steps) for tracers [] 
    650    rn_trelax_dyn = 0.01    !  inverse of relaxation time (in steps) for dynamics [] 
    651    ln_chk_bathy  = .false. !  =T  check the parent bathymetry 
     645   ln_agrif_2way   = .true.  !  activate two way nesting 
     646   ln_init_chfrpar = .false. !  initialize child grids from parent 
     647   ln_spc_dyn      = .true.  !  use 0 as special value for dynamics 
     648   rn_sponge_tra   = 0.002   !  coefficient for tracer   sponge layer [] 
     649   rn_sponge_dyn   = 0.002   !  coefficient for dynamics sponge layer [] 
     650   rn_trelax_tra   = 0.01    !  inverse of relaxation time (in steps) for tracers [] 
     651   rn_trelax_dyn   = 0.01    !  inverse of relaxation time (in steps) for dynamics [] 
     652   ln_chk_bathy    = .false. !  =T  check the parent bathymetry 
    652653/ 
    653654!----------------------------------------------------------------------- 
  • NEMO/trunk/src/ICE/iceistate.F90

    r12736 r13216  
    3232   USE lib_fortran    ! fortran utilities (glob_sum + no signed zero) 
    3333   USE fldread        ! read input fields 
     34 
     35# if defined key_agrif 
     36   USE agrif_oce 
     37   USE agrif_ice 
     38   USE agrif_ice_interp  
     39# endif    
    3440 
    3541   IMPLICIT NONE 
     
    168174      ! 2) overwrite some of the fields with namelist parameters or netcdf file 
    169175      !------------------------------------------------------------------------ 
     176 
     177 
    170178      IF( ln_iceini ) THEN 
    171179         !                             !---------------! 
    172          IF( ln_iceini_file )THEN      ! Read a file   ! 
    173             !                          !---------------! 
    174             WHERE( ff_t(:,:) >= 0._wp )   ;   zswitch(:,:) = 1._wp 
    175             ELSEWHERE                     ;   zswitch(:,:) = 0._wp 
     180          
     181         IF( Agrif_Root() ) THEN 
     182 
     183            IF( ln_iceini_file )THEN      ! Read a file   ! 
     184               !                          !---------------! 
     185               WHERE( ff_t(:,:) >= 0._wp )   ;   zswitch(:,:) = 1._wp 
     186               ELSEWHERE                     ;   zswitch(:,:) = 0._wp 
     187               END WHERE 
     188               ! 
     189               CALL fld_read( kt, 1, si ) ! input fields provided at the current time-step 
     190               ! 
     191               ! -- mandatory fields -- ! 
     192               zht_i_ini(:,:) = si(jp_hti)%fnow(:,:,1) * tmask(:,:,1) 
     193               zht_s_ini(:,:) = si(jp_hts)%fnow(:,:,1) * tmask(:,:,1) 
     194               zat_i_ini(:,:) = si(jp_ati)%fnow(:,:,1) * tmask(:,:,1) 
     195 
     196               ! -- optional fields -- ! 
     197               !    if fields do not exist then set them to the values present in the namelist (except for snow and surface temperature) 
     198               ! 
     199               ! ice salinity 
     200               IF( TRIM(si(jp_smi)%clrootname) == 'NOT USED' ) & 
     201                  &     si(jp_smi)%fnow(:,:,1) = ( rn_smi_ini_n * zswitch + rn_smi_ini_s * (1._wp - zswitch) ) * tmask(:,:,1) 
     202               ! 
     203               ! temperatures 
     204               IF    ( TRIM(si(jp_tmi)%clrootname) == 'NOT USED' .AND. TRIM(si(jp_tsu)%clrootname) == 'NOT USED' .AND. & 
     205                  &    TRIM(si(jp_tms)%clrootname) == 'NOT USED' ) THEN 
     206                  si(jp_tmi)%fnow(:,:,1) = ( rn_tmi_ini_n * zswitch + rn_tmi_ini_s * (1._wp - zswitch) ) * tmask(:,:,1) 
     207                  si(jp_tsu)%fnow(:,:,1) = ( rn_tsu_ini_n * zswitch + rn_tsu_ini_s * (1._wp - zswitch) ) * tmask(:,:,1) 
     208                  si(jp_tms)%fnow(:,:,1) = ( rn_tms_ini_n * zswitch + rn_tms_ini_s * (1._wp - zswitch) ) * tmask(:,:,1) 
     209               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 
     210                  si(jp_tmi)%fnow(:,:,1) = 0.5_wp * ( si(jp_tms)%fnow(:,:,1) + 271.15 ) 
     211               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 
     212                  si(jp_tmi)%fnow(:,:,1) = 0.5_wp * ( si(jp_tsu)%fnow(:,:,1) + 271.15 ) 
     213               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 
     214                  si(jp_tsu)%fnow(:,:,1) = si(jp_tms)%fnow(:,:,1) 
     215               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 
     216                  si(jp_tsu)%fnow(:,:,1) = si(jp_tmi)%fnow(:,:,1) 
     217               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 
     218                  si(jp_tms)%fnow(:,:,1) = si(jp_tsu)%fnow(:,:,1) 
     219               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 
     220                  si(jp_tms)%fnow(:,:,1) = si(jp_tmi)%fnow(:,:,1) 
     221               ENDIF 
     222               ! 
     223               ! pond concentration 
     224               IF( TRIM(si(jp_apd)%clrootname) == 'NOT USED' ) & 
     225                  &     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. 
     226                  &                              * si(jp_ati)%fnow(:,:,1)  
     227               ! 
     228               ! pond depth 
     229               IF( TRIM(si(jp_hpd)%clrootname) == 'NOT USED' ) & 
     230                  &     si(jp_hpd)%fnow(:,:,1) = ( rn_hpd_ini_n * zswitch + rn_hpd_ini_s * (1._wp - zswitch) ) * tmask(:,:,1) 
     231               ! 
     232               zsm_i_ini(:,:) = si(jp_smi)%fnow(:,:,1) * tmask(:,:,1) 
     233               ztm_i_ini(:,:) = si(jp_tmi)%fnow(:,:,1) * tmask(:,:,1) 
     234               zt_su_ini(:,:) = si(jp_tsu)%fnow(:,:,1) * tmask(:,:,1) 
     235               ztm_s_ini(:,:) = si(jp_tms)%fnow(:,:,1) * tmask(:,:,1) 
     236               zapnd_ini(:,:) = si(jp_apd)%fnow(:,:,1) * tmask(:,:,1) 
     237               zhpnd_ini(:,:) = si(jp_hpd)%fnow(:,:,1) * tmask(:,:,1) 
     238               ! 
     239               ! change the switch for the following 
     240               WHERE( zat_i_ini(:,:) > 0._wp )   ;   zswitch(:,:) = tmask(:,:,1)  
     241               ELSEWHERE                         ;   zswitch(:,:) = 0._wp 
     242               END WHERE 
     243 
     244               !                          !---------------! 
     245            ELSE                          ! Read namelist ! 
     246               !                          !---------------! 
     247               ! no ice if (sst - Tfreez) >= thresold 
     248               WHERE( ( sst_m(:,:) - (t_bo(:,:) - rt0) ) * tmask(:,:,1) >= rn_thres_sst )   ;   zswitch(:,:) = 0._wp  
     249               ELSEWHERE                                                                    ;   zswitch(:,:) = tmask(:,:,1) 
     250               END WHERE 
     251               ! 
     252               ! assign initial thickness, concentration, snow depth and salinity to an hemisphere-dependent array 
     253               WHERE( ff_t(:,:) >= 0._wp ) 
     254                  zht_i_ini(:,:) = rn_hti_ini_n * zswitch(:,:) 
     255                  zht_s_ini(:,:) = rn_hts_ini_n * zswitch(:,:) 
     256                  zat_i_ini(:,:) = rn_ati_ini_n * zswitch(:,:) 
     257                  zsm_i_ini(:,:) = rn_smi_ini_n * zswitch(:,:) 
     258                  ztm_i_ini(:,:) = rn_tmi_ini_n * zswitch(:,:) 
     259                  zt_su_ini(:,:) = rn_tsu_ini_n * zswitch(:,:) 
     260                  ztm_s_ini(:,:) = rn_tms_ini_n * zswitch(:,:) 
     261                  zapnd_ini(:,:) = rn_apd_ini_n * zswitch(:,:) * zat_i_ini(:,:) ! rn_apd = pond fraction => rn_apd * a_i = pond conc.  
     262                  zhpnd_ini(:,:) = rn_hpd_ini_n * zswitch(:,:) 
     263               ELSEWHERE 
     264                  zht_i_ini(:,:) = rn_hti_ini_s * zswitch(:,:) 
     265                  zht_s_ini(:,:) = rn_hts_ini_s * zswitch(:,:) 
     266                  zat_i_ini(:,:) = rn_ati_ini_s * zswitch(:,:) 
     267                  zsm_i_ini(:,:) = rn_smi_ini_s * zswitch(:,:) 
     268                  ztm_i_ini(:,:) = rn_tmi_ini_s * zswitch(:,:) 
     269                  zt_su_ini(:,:) = rn_tsu_ini_s * zswitch(:,:) 
     270                  ztm_s_ini(:,:) = rn_tms_ini_s * zswitch(:,:) 
     271                  zapnd_ini(:,:) = rn_apd_ini_s * zswitch(:,:) * zat_i_ini(:,:) ! rn_apd = pond fraction => rn_apd * a_i = pond conc. 
     272                  zhpnd_ini(:,:) = rn_hpd_ini_s * zswitch(:,:) 
     273               END WHERE 
     274               ! 
     275            ENDIF 
     276 
     277 
     278 
     279            ! make sure ponds = 0 if no ponds scheme 
     280            IF ( .NOT.ln_pnd ) THEN 
     281               zapnd_ini(:,:) = 0._wp 
     282               zhpnd_ini(:,:) = 0._wp 
     283            ENDIF 
     284             
     285            !-------------! 
     286            ! fill fields ! 
     287            !-------------! 
     288            ! select ice covered grid points 
     289            npti = 0 ; nptidx(:) = 0 
     290            DO_2D_11_11 
     291               IF ( zht_i_ini(ji,jj) > 0._wp ) THEN 
     292                  npti         = npti  + 1 
     293                  nptidx(npti) = (jj - 1) * jpi + ji 
     294               ENDIF 
     295            END_2D 
     296 
     297            ! move to 1D arrays: (jpi,jpj) -> (jpi*jpj) 
     298            CALL tab_2d_1d( npti, nptidx(1:npti), h_i_1d (1:npti)  , zht_i_ini ) 
     299            CALL tab_2d_1d( npti, nptidx(1:npti), h_s_1d (1:npti)  , zht_s_ini ) 
     300            CALL tab_2d_1d( npti, nptidx(1:npti), at_i_1d(1:npti)  , zat_i_ini ) 
     301            CALL tab_2d_1d( npti, nptidx(1:npti), t_i_1d (1:npti,1), ztm_i_ini ) 
     302            CALL tab_2d_1d( npti, nptidx(1:npti), t_s_1d (1:npti,1), ztm_s_ini ) 
     303            CALL tab_2d_1d( npti, nptidx(1:npti), t_su_1d(1:npti)  , zt_su_ini ) 
     304            CALL tab_2d_1d( npti, nptidx(1:npti), s_i_1d (1:npti)  , zsm_i_ini ) 
     305            CALL tab_2d_1d( npti, nptidx(1:npti), a_ip_1d(1:npti)  , zapnd_ini ) 
     306            CALL tab_2d_1d( npti, nptidx(1:npti), h_ip_1d(1:npti)  , zhpnd_ini ) 
     307 
     308            ! allocate temporary arrays 
     309            ALLOCATE( zhi_2d(npti,jpl), zhs_2d(npti,jpl), zai_2d (npti,jpl), & 
     310               &      zti_2d(npti,jpl), zts_2d(npti,jpl), ztsu_2d(npti,jpl), zsi_2d(npti,jpl), zaip_2d(npti,jpl), zhip_2d(npti,jpl) ) 
     311             
     312            ! distribute 1-cat into jpl-cat: (jpi*jpj) -> (jpi*jpj,jpl) 
     313            CALL ice_var_itd( h_i_1d(1:npti)  , h_s_1d(1:npti)  , at_i_1d(1:npti),                                                   & 
     314               &              zhi_2d          , zhs_2d          , zai_2d         ,                                                   & 
     315               &              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), & 
     316               &              zti_2d          , zts_2d          , ztsu_2d        , zsi_2d        , zaip_2d        , zhip_2d ) 
     317 
     318            ! move to 3D arrays: (jpi*jpj,jpl) -> (jpi,jpj,jpl) 
     319            DO jl = 1, jpl 
     320               zti_3d(:,:,jl) = rt0 * tmask(:,:,1) 
     321               zts_3d(:,:,jl) = rt0 * tmask(:,:,1) 
     322            END DO 
     323            CALL tab_2d_3d( npti, nptidx(1:npti), zhi_2d   , h_i    ) 
     324            CALL tab_2d_3d( npti, nptidx(1:npti), zhs_2d   , h_s    ) 
     325            CALL tab_2d_3d( npti, nptidx(1:npti), zai_2d   , a_i    ) 
     326            CALL tab_2d_3d( npti, nptidx(1:npti), zti_2d   , zti_3d ) 
     327            CALL tab_2d_3d( npti, nptidx(1:npti), zts_2d   , zts_3d ) 
     328            CALL tab_2d_3d( npti, nptidx(1:npti), ztsu_2d  , t_su   ) 
     329            CALL tab_2d_3d( npti, nptidx(1:npti), zsi_2d   , s_i    ) 
     330            CALL tab_2d_3d( npti, nptidx(1:npti), zaip_2d  , a_ip   ) 
     331            CALL tab_2d_3d( npti, nptidx(1:npti), zhip_2d  , h_ip   ) 
     332 
     333            ! deallocate temporary arrays 
     334            DEALLOCATE( zhi_2d, zhs_2d, zai_2d , & 
     335               &        zti_2d, zts_2d, ztsu_2d, zsi_2d, zaip_2d, zhip_2d ) 
     336 
     337            ! calculate extensive and intensive variables 
     338            CALL ice_var_salprof ! for sz_i 
     339            DO jl = 1, jpl 
     340               DO_2D_11_11 
     341                  v_i (ji,jj,jl) = h_i(ji,jj,jl) * a_i(ji,jj,jl) 
     342                  v_s (ji,jj,jl) = h_s(ji,jj,jl) * a_i(ji,jj,jl) 
     343                  sv_i(ji,jj,jl) = MIN( MAX( rn_simin , s_i(ji,jj,jl) ) , rn_simax ) * v_i(ji,jj,jl) 
     344               END_2D 
     345            END DO 
     346            ! 
     347            DO jl = 1, jpl 
     348               DO_3D_11_11( 1, nlay_s ) 
     349                  t_s(ji,jj,jk,jl) = zts_3d(ji,jj,jl) 
     350                  e_s(ji,jj,jk,jl) = zswitch(ji,jj) * v_s(ji,jj,jl) * r1_nlay_s * & 
     351                     &               rhos * ( rcpi * ( rt0 - t_s(ji,jj,jk,jl) ) + rLfus ) 
     352               END_3D 
     353            END DO 
     354            ! 
     355            DO jl = 1, jpl 
     356               DO_3D_11_11( 1, nlay_i ) 
     357                  t_i (ji,jj,jk,jl) = zti_3d(ji,jj,jl)  
     358                  ztmelts          = - rTmlt * sz_i(ji,jj,jk,jl) + rt0 ! melting temperature in K 
     359                  e_i(ji,jj,jk,jl) = zswitch(ji,jj) * v_i(ji,jj,jl) * r1_nlay_i * & 
     360                     &               rhoi * (  rcpi  * ( ztmelts - t_i(ji,jj,jk,jl) ) + & 
     361                     &                         rLfus * ( 1._wp - (ztmelts-rt0) / MIN( (t_i(ji,jj,jk,jl)-rt0), -epsi20 ) ) & 
     362                     &                       - rcp   * ( ztmelts - rt0 ) ) 
     363               END_3D 
     364            END DO 
     365 
     366            ! Melt ponds 
     367            WHERE( a_i > epsi10 ) 
     368               a_ip_frac(:,:,:) = a_ip(:,:,:) / a_i(:,:,:) 
     369            ELSEWHERE 
     370               a_ip_frac(:,:,:) = 0._wp 
    176371            END WHERE 
     372            v_ip(:,:,:) = h_ip(:,:,:) * a_ip(:,:,:) 
     373              
     374            ! specific temperatures for coupled runs 
     375            tn_ice(:,:,:) = t_su(:,:,:) 
     376            t1_ice(:,:,:) = t_i (:,:,1,:) 
    177377            ! 
    178             CALL fld_read( kt, 1, si ) ! input fields provided at the current time-step 
    179             ! 
    180             ! -- mandatory fields -- ! 
    181             zht_i_ini(:,:) = si(jp_hti)%fnow(:,:,1) * tmask(:,:,1) 
    182             zht_s_ini(:,:) = si(jp_hts)%fnow(:,:,1) * tmask(:,:,1) 
    183             zat_i_ini(:,:) = si(jp_ati)%fnow(:,:,1) * tmask(:,:,1) 
    184  
    185             ! -- optional fields -- ! 
    186             !    if fields do not exist then set them to the values present in the namelist (except for snow and surface temperature) 
    187             ! 
    188             ! ice salinity 
    189             IF( TRIM(si(jp_smi)%clrootname) == 'NOT USED' ) & 
    190                &     si(jp_smi)%fnow(:,:,1) = ( rn_smi_ini_n * zswitch + rn_smi_ini_s * (1._wp - zswitch) ) * tmask(:,:,1) 
    191             ! 
    192             ! temperatures 
    193             IF    ( TRIM(si(jp_tmi)%clrootname) == 'NOT USED' .AND. TRIM(si(jp_tsu)%clrootname) == 'NOT USED' .AND. & 
    194                &    TRIM(si(jp_tms)%clrootname) == 'NOT USED' ) THEN 
    195                si(jp_tmi)%fnow(:,:,1) = ( rn_tmi_ini_n * zswitch + rn_tmi_ini_s * (1._wp - zswitch) ) * tmask(:,:,1) 
    196                si(jp_tsu)%fnow(:,:,1) = ( rn_tsu_ini_n * zswitch + rn_tsu_ini_s * (1._wp - zswitch) ) * tmask(:,:,1) 
    197                si(jp_tms)%fnow(:,:,1) = ( rn_tms_ini_n * zswitch + rn_tms_ini_s * (1._wp - zswitch) ) * tmask(:,:,1) 
    198             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 
    199                si(jp_tmi)%fnow(:,:,1) = 0.5_wp * ( si(jp_tms)%fnow(:,:,1) + 271.15 ) 
    200             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 
    201                si(jp_tmi)%fnow(:,:,1) = 0.5_wp * ( si(jp_tsu)%fnow(:,:,1) + 271.15 ) 
    202             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 
    203                si(jp_tsu)%fnow(:,:,1) = si(jp_tms)%fnow(:,:,1) 
    204             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 
    205                si(jp_tsu)%fnow(:,:,1) = si(jp_tmi)%fnow(:,:,1) 
    206             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 
    207                si(jp_tms)%fnow(:,:,1) = si(jp_tsu)%fnow(:,:,1) 
    208             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 
    209                si(jp_tms)%fnow(:,:,1) = si(jp_tmi)%fnow(:,:,1) 
    210             ENDIF 
    211             ! 
    212             ! pond concentration 
    213             IF( TRIM(si(jp_apd)%clrootname) == 'NOT USED' ) & 
    214                &     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. 
    215                &                              * si(jp_ati)%fnow(:,:,1)  
    216             ! 
    217             ! pond depth 
    218             IF( TRIM(si(jp_hpd)%clrootname) == 'NOT USED' ) & 
    219                &     si(jp_hpd)%fnow(:,:,1) = ( rn_hpd_ini_n * zswitch + rn_hpd_ini_s * (1._wp - zswitch) ) * tmask(:,:,1) 
    220             ! 
    221             zsm_i_ini(:,:) = si(jp_smi)%fnow(:,:,1) * tmask(:,:,1) 
    222             ztm_i_ini(:,:) = si(jp_tmi)%fnow(:,:,1) * tmask(:,:,1) 
    223             zt_su_ini(:,:) = si(jp_tsu)%fnow(:,:,1) * tmask(:,:,1) 
    224             ztm_s_ini(:,:) = si(jp_tms)%fnow(:,:,1) * tmask(:,:,1) 
    225             zapnd_ini(:,:) = si(jp_apd)%fnow(:,:,1) * tmask(:,:,1) 
    226             zhpnd_ini(:,:) = si(jp_hpd)%fnow(:,:,1) * tmask(:,:,1) 
    227             ! 
    228             ! change the switch for the following 
    229             WHERE( zat_i_ini(:,:) > 0._wp )   ;   zswitch(:,:) = tmask(:,:,1)  
    230             ELSEWHERE                         ;   zswitch(:,:) = 0._wp 
     378          
     379#if  defined key_agrif 
     380         ELSE 
     381  
     382            Agrif_SpecialValue    = -9999. 
     383            Agrif_UseSpecialValue = .TRUE. 
     384            CALL Agrif_init_variable(tra_iceini_id,procname=interp_tra_ice) 
     385            use_sign_north = .TRUE. 
     386            sign_north = -1. 
     387            CALL Agrif_init_variable(u_iceini_id  ,procname=interp_u_ice) 
     388            CALL Agrif_init_variable(v_iceini_id  ,procname=interp_v_ice) 
     389            Agrif_SpecialValue    = 0._wp 
     390            use_sign_north = .FALSE. 
     391            Agrif_UseSpecialValue = .FALSE. 
     392        ! lbc ????  
     393   ! Here we know : a_i, v_i, v_s, sv_i, oa_i, a_ip, v_ip, t_su, e_s, e_i 
     394            CALL ice_var_glo2eqv 
     395            CALL ice_var_zapsmall 
     396            CALL ice_var_agg(2) 
     397 
     398            ! Melt ponds 
     399            WHERE( a_i > epsi10 ) 
     400               a_ip_frac(:,:,:) = a_ip(:,:,:) / a_i(:,:,:) 
     401            ELSEWHERE 
     402               a_ip_frac(:,:,:) = 0._wp 
    231403            END WHERE 
    232             !                          !---------------! 
    233          ELSE                          ! Read namelist ! 
    234             !                          !---------------! 
    235             ! no ice if (sst - Tfreez) >= thresold 
    236             WHERE( ( sst_m(:,:) - (t_bo(:,:) - rt0) ) * tmask(:,:,1) >= rn_thres_sst )   ;   zswitch(:,:) = 0._wp  
    237             ELSEWHERE                                                                    ;   zswitch(:,:) = tmask(:,:,1) 
    238             END WHERE 
    239             ! 
    240             ! assign initial thickness, concentration, snow depth and salinity to an hemisphere-dependent array 
    241             WHERE( ff_t(:,:) >= 0._wp ) 
    242                zht_i_ini(:,:) = rn_hti_ini_n * zswitch(:,:) 
    243                zht_s_ini(:,:) = rn_hts_ini_n * zswitch(:,:) 
    244                zat_i_ini(:,:) = rn_ati_ini_n * zswitch(:,:) 
    245                zsm_i_ini(:,:) = rn_smi_ini_n * zswitch(:,:) 
    246                ztm_i_ini(:,:) = rn_tmi_ini_n * zswitch(:,:) 
    247                zt_su_ini(:,:) = rn_tsu_ini_n * zswitch(:,:) 
    248                ztm_s_ini(:,:) = rn_tms_ini_n * zswitch(:,:) 
    249                zapnd_ini(:,:) = rn_apd_ini_n * zswitch(:,:) * zat_i_ini(:,:) ! rn_apd = pond fraction => rn_apd * a_i = pond conc.  
    250                zhpnd_ini(:,:) = rn_hpd_ini_n * zswitch(:,:) 
     404            WHERE( a_ip > 0._wp )       ! ???????     
     405               h_ip(:,:,:) = v_ip(:,:,:) / a_ip(:,:,:) 
    251406            ELSEWHERE 
    252                zht_i_ini(:,:) = rn_hti_ini_s * zswitch(:,:) 
    253                zht_s_ini(:,:) = rn_hts_ini_s * zswitch(:,:) 
    254                zat_i_ini(:,:) = rn_ati_ini_s * zswitch(:,:) 
    255                zsm_i_ini(:,:) = rn_smi_ini_s * zswitch(:,:) 
    256                ztm_i_ini(:,:) = rn_tmi_ini_s * zswitch(:,:) 
    257                zt_su_ini(:,:) = rn_tsu_ini_s * zswitch(:,:) 
    258                ztm_s_ini(:,:) = rn_tms_ini_s * zswitch(:,:) 
    259                zapnd_ini(:,:) = rn_apd_ini_s * zswitch(:,:) * zat_i_ini(:,:) ! rn_apd = pond fraction => rn_apd * a_i = pond conc. 
    260                zhpnd_ini(:,:) = rn_hpd_ini_s * zswitch(:,:) 
    261             END WHERE 
    262             ! 
    263          ENDIF 
    264  
    265          ! make sure ponds = 0 if no ponds scheme 
    266          IF ( .NOT.ln_pnd ) THEN 
    267             zapnd_ini(:,:) = 0._wp 
    268             zhpnd_ini(:,:) = 0._wp 
    269          ENDIF 
    270           
    271          !-------------! 
    272          ! fill fields ! 
    273          !-------------! 
    274          ! select ice covered grid points 
    275          npti = 0 ; nptidx(:) = 0 
    276          DO_2D_11_11 
    277             IF ( zht_i_ini(ji,jj) > 0._wp ) THEN 
    278                npti         = npti  + 1 
    279                nptidx(npti) = (jj - 1) * jpi + ji 
    280             ENDIF 
    281          END_2D 
    282  
    283          ! move to 1D arrays: (jpi,jpj) -> (jpi*jpj) 
    284          CALL tab_2d_1d( npti, nptidx(1:npti), h_i_1d (1:npti)  , zht_i_ini ) 
    285          CALL tab_2d_1d( npti, nptidx(1:npti), h_s_1d (1:npti)  , zht_s_ini ) 
    286          CALL tab_2d_1d( npti, nptidx(1:npti), at_i_1d(1:npti)  , zat_i_ini ) 
    287          CALL tab_2d_1d( npti, nptidx(1:npti), t_i_1d (1:npti,1), ztm_i_ini ) 
    288          CALL tab_2d_1d( npti, nptidx(1:npti), t_s_1d (1:npti,1), ztm_s_ini ) 
    289          CALL tab_2d_1d( npti, nptidx(1:npti), t_su_1d(1:npti)  , zt_su_ini ) 
    290          CALL tab_2d_1d( npti, nptidx(1:npti), s_i_1d (1:npti)  , zsm_i_ini ) 
    291          CALL tab_2d_1d( npti, nptidx(1:npti), a_ip_1d(1:npti)  , zapnd_ini ) 
    292          CALL tab_2d_1d( npti, nptidx(1:npti), h_ip_1d(1:npti)  , zhpnd_ini ) 
    293  
    294          ! allocate temporary arrays 
    295          ALLOCATE( zhi_2d(npti,jpl), zhs_2d(npti,jpl), zai_2d (npti,jpl), & 
    296             &      zti_2d(npti,jpl), zts_2d(npti,jpl), ztsu_2d(npti,jpl), zsi_2d(npti,jpl), zaip_2d(npti,jpl), zhip_2d(npti,jpl) ) 
    297           
    298          ! distribute 1-cat into jpl-cat: (jpi*jpj) -> (jpi*jpj,jpl) 
    299          CALL ice_var_itd( h_i_1d(1:npti)  , h_s_1d(1:npti)  , at_i_1d(1:npti),                                                   & 
    300             &              zhi_2d          , zhs_2d          , zai_2d         ,                                                   & 
    301             &              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), & 
    302             &              zti_2d          , zts_2d          , ztsu_2d        , zsi_2d        , zaip_2d        , zhip_2d ) 
    303  
    304          ! move to 3D arrays: (jpi*jpj,jpl) -> (jpi,jpj,jpl) 
    305          DO jl = 1, jpl 
    306             zti_3d(:,:,jl) = rt0 * tmask(:,:,1) 
    307             zts_3d(:,:,jl) = rt0 * tmask(:,:,1) 
    308          END DO 
    309          CALL tab_2d_3d( npti, nptidx(1:npti), zhi_2d   , h_i    ) 
    310          CALL tab_2d_3d( npti, nptidx(1:npti), zhs_2d   , h_s    ) 
    311          CALL tab_2d_3d( npti, nptidx(1:npti), zai_2d   , a_i    ) 
    312          CALL tab_2d_3d( npti, nptidx(1:npti), zti_2d   , zti_3d ) 
    313          CALL tab_2d_3d( npti, nptidx(1:npti), zts_2d   , zts_3d ) 
    314          CALL tab_2d_3d( npti, nptidx(1:npti), ztsu_2d  , t_su   ) 
    315          CALL tab_2d_3d( npti, nptidx(1:npti), zsi_2d   , s_i    ) 
    316          CALL tab_2d_3d( npti, nptidx(1:npti), zaip_2d  , a_ip   ) 
    317          CALL tab_2d_3d( npti, nptidx(1:npti), zhip_2d  , h_ip   ) 
    318  
    319          ! deallocate temporary arrays 
    320          DEALLOCATE( zhi_2d, zhs_2d, zai_2d , & 
    321             &        zti_2d, zts_2d, ztsu_2d, zsi_2d, zaip_2d, zhip_2d ) 
    322  
    323          ! calculate extensive and intensive variables 
    324          CALL ice_var_salprof ! for sz_i 
    325          DO jl = 1, jpl 
    326             DO_2D_11_11 
    327                v_i (ji,jj,jl) = h_i(ji,jj,jl) * a_i(ji,jj,jl) 
    328                v_s (ji,jj,jl) = h_s(ji,jj,jl) * a_i(ji,jj,jl) 
    329                sv_i(ji,jj,jl) = MIN( MAX( rn_simin , s_i(ji,jj,jl) ) , rn_simax ) * v_i(ji,jj,jl) 
    330             END_2D 
    331          END DO 
    332          ! 
    333          DO jl = 1, jpl 
    334             DO_3D_11_11( 1, nlay_s ) 
    335                t_s(ji,jj,jk,jl) = zts_3d(ji,jj,jl) 
    336                e_s(ji,jj,jk,jl) = zswitch(ji,jj) * v_s(ji,jj,jl) * r1_nlay_s * & 
    337                   &               rhos * ( rcpi * ( rt0 - t_s(ji,jj,jk,jl) ) + rLfus ) 
    338             END_3D 
    339          END DO 
    340          ! 
    341          DO jl = 1, jpl 
    342             DO_3D_11_11( 1, nlay_i ) 
    343                t_i (ji,jj,jk,jl) = zti_3d(ji,jj,jl)  
    344                ztmelts          = - rTmlt * sz_i(ji,jj,jk,jl) + rt0 ! melting temperature in K 
    345                e_i(ji,jj,jk,jl) = zswitch(ji,jj) * v_i(ji,jj,jl) * r1_nlay_i * & 
    346                   &               rhoi * (  rcpi  * ( ztmelts - t_i(ji,jj,jk,jl) ) + & 
    347                   &                         rLfus * ( 1._wp - (ztmelts-rt0) / MIN( (t_i(ji,jj,jk,jl)-rt0), -epsi20 ) ) & 
    348                   &                       - rcp   * ( ztmelts - rt0 ) ) 
    349             END_3D 
    350          END DO 
    351  
    352          ! Melt ponds 
    353          WHERE( a_i > epsi10 ) 
    354             a_ip_frac(:,:,:) = a_ip(:,:,:) / a_i(:,:,:) 
    355          ELSEWHERE 
    356             a_ip_frac(:,:,:) = 0._wp 
    357          END WHERE 
    358          v_ip(:,:,:) = h_ip(:,:,:) * a_ip(:,:,:) 
    359            
    360          ! specific temperatures for coupled runs 
    361          tn_ice(:,:,:) = t_su(:,:,:) 
    362          t1_ice(:,:,:) = t_i (:,:,1,:) 
    363          ! 
     407               h_ip(:,:,:) = 0._wp 
     408            END WHERE    
     409 
     410            tn_ice(:,:,:) = t_su(:,:,:) 
     411            t1_ice(:,:,:) = t_i (:,:,1,:) 
     412#endif 
     413          ENDIF ! Agrif_Root 
    364414      ENDIF ! ln_iceini 
    365415      ! 
  • NEMO/trunk/src/ICE/icestp.F90

    r12489 r13216  
    240240      CALL par_init                ! set some ice run parameters 
    241241      ! 
     242#if defined key_agrif 
     243      CALL Agrif_Declare_Var_ice  !  "      "   "   "      "  Sea ice 
     244#endif 
     245      ! 
    242246      !                                ! Allocate the ice arrays (sbc_ice already allocated in sbc_init) 
    243247      ierr =        ice_alloc        ()      ! ice variables 
  • NEMO/trunk/src/NST/agrif_ice.F90

    r10068 r13216  
    1616 
    1717   INTEGER, PUBLIC ::  u_ice_id, v_ice_id, tra_ice_id 
     18   INTEGER, PUBLIC ::  u_iceini_id, v_iceini_id, tra_iceini_id 
    1819   INTEGER, PUBLIC ::  nbstep_ice = 0    ! child time position in sea-ice model 
    1920 
  • NEMO/trunk/src/NST/agrif_ice_interp.F90

    r10069 r13216  
    1414   !!---------------------------------------------------------------------- 
    1515   !!  agrif_interp_ice    : interpolation of ice at "after" sea-ice time step 
    16    !!  agrif_interp_u_ice   : atomic routine to interpolate u_ice  
    17    !!  agrif_interp_v_ice   : atomic routine to interpolate v_ice  
    18    !!  agrif_interp_tra_ice : atomic routine to interpolate ice properties  
     16   !!  interp_u_ice   : atomic routine to interpolate u_ice  
     17   !!  interp_v_ice   : atomic routine to interpolate v_ice  
     18   !!  interp_tra_ice : atomic routine to interpolate ice properties  
    1919   !!---------------------------------------------------------------------- 
    2020   USE par_oce 
     
    2323   USE ice 
    2424   USE agrif_ice 
     25   USE agrif_oce 
    2526   USE phycst , ONLY: rt0 
    2627    
     
    2930 
    3031   PUBLIC   agrif_interp_ice   ! called by agrif_user.F90 
     32   PUBLIC   interp_tra_ice, interp_u_ice, interp_v_ice  ! called by iceistate.F90 
    3133 
    3234   !!---------------------------------------------------------------------- 
     
    6870      Agrif_SpecialValue    = -9999. 
    6971      Agrif_UseSpecialValue = .TRUE. 
     72 
     73      use_sign_north = .TRUE. 
     74      sign_north = -1. 
     75      if (cd_type == 'T') use_sign_north = .FALSE. 
     76 
    7077      SELECT CASE( cd_type ) 
    7178      CASE('U')   ;   CALL Agrif_Bc_variable( u_ice_id  , procname=interp_u_ice  , calledweight=zbeta ) 
     
    7582      Agrif_SpecialValue    = 0._wp 
    7683      Agrif_UseSpecialValue = .FALSE. 
     84       
     85      use_sign_north = .FALSE. 
    7786      ! 
    7887   END SUBROUTINE agrif_interp_ice 
     
    156165      ! and it is ok since we conserve tracers (same as in the ocean). 
    157166      ALLOCATE( ztab(SIZE(a_i,1),SIZE(a_i,2),SIZE(ptab,3)) ) 
    158       
     167 
    159168      IF( before ) THEN  ! parent grid 
    160169         jm = 1 
  • NEMO/trunk/src/NST/agrif_ice_update.F90

    r12377 r13216  
    6666      CALL Agrif_Update_Variable( tra_ice_id , locupdate=(/1,0/), procname = update_tra_ice  ) 
    6767#endif 
     68      use_sign_north = .TRUE. 
     69      sign_north = -1. 
     70 
    6871# if ! defined DECAL_FEEDBACK 
    6972      CALL Agrif_Update_Variable( u_ice_id   , procname = update_u_ice    ) 
     
    7376      CALL Agrif_Update_Variable( v_ice_id   , locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname=update_v_ice) 
    7477#endif 
     78      use_sign_north = .FALSE. 
    7579!      CALL Agrif_Update_Variable( tra_ice_id , locupdate=(/0,2/), procname = update_tra_ice  ) 
    7680!      CALL Agrif_Update_Variable( u_ice_id   , locupdate=(/0,1/), procname = update_u_ice    ) 
  • NEMO/trunk/src/NST/agrif_oce.F90

    r12377 r13216  
    1919   
    2020   !                                              !!* Namelist namagrif: AGRIF parameters 
     21   LOGICAL , PUBLIC ::   ln_init_chfrpar = .FALSE. !: set child grids initial state from parent 
    2122   LOGICAL , PUBLIC ::   ln_agrif_2way = .TRUE.    !: activate two way nesting  
    2223   LOGICAL , PUBLIC ::   ln_spc_dyn    = .FALSE.   !: use zeros (.false.) or not (.true.) in 
     
    2930   ! 
    3031   INTEGER , PUBLIC, PARAMETER ::   nn_sponge_len = 2  !: Sponge width (in number of parent grid points) 
     32 
    3133   LOGICAL , PUBLIC :: spongedoneT = .FALSE.       !: tracer   sponge layer indicator 
    3234   LOGICAL , PUBLIC :: spongedoneU = .FALSE.       !: dynamics sponge layer indicator 
     
    4951   INTEGER , PUBLIC,              SAVE                 ::  Kbb_a, Kmm_a, Krhs_a   !: AGRIF module-specific copies of time-level indices 
    5052 
    51 # if defined key_vertical 
    5253   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ht0_parent, hu0_parent, hv0_parent 
    5354   INTEGER,  PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: mbkt_parent, mbku_parent, mbkv_parent 
    54 # endif 
    5555 
    5656   INTEGER, PUBLIC :: tsn_id                                                  ! AGRIF profile for tracers interpolation and update 
     
    5858   INTEGER, PUBLIC :: un_update_id, vn_update_id                              ! AGRIF profiles for udpates 
    5959   INTEGER, PUBLIC :: tsn_sponge_id, un_sponge_id, vn_sponge_id               ! AGRIF profiles for sponge layers 
     60   INTEGER, PUBLIC :: tsini_id, uini_id, vini_id, sshini_id                   ! AGRIF profile for initialization 
    6061# if defined key_top 
    6162   INTEGER, PUBLIC :: trn_id, trn_sponge_id 
     
    6869   INTEGER, PUBLIC :: mbkt_id, ht0_id 
    6970   INTEGER, PUBLIC :: kindic_agr 
     71 
     72   ! North fold 
     73!$AGRIF_DO_NOT_TREAT 
     74   LOGICAL, PUBLIC :: use_sign_north 
     75   REAL, PUBLIC :: sign_north 
     76   LOGICAL, PUBLIC :: l_ini_child = .FALSE. 
     77# if defined key_vertical 
     78   LOGICAL, PUBLIC :: l_vremap    = .TRUE. 
     79# else 
     80   LOGICAL, PUBLIC :: l_vremap    = .FALSE. 
     81# endif 
     82!$AGRIF_END_DO_NOT_TREAT 
    7083    
    7184   !!---------------------------------------------------------------------- 
     
    91104         &      tabspongedone_trn(jpi,jpj),           & 
    92105# endif    
    93 # if defined key_vertical 
    94106         &      ht0_parent(jpi,jpj), mbkt_parent(jpi,jpj),  & 
    95107         &      hu0_parent(jpi,jpj), mbku_parent(jpi,jpj),  & 
    96108         &      hv0_parent(jpi,jpj), mbkv_parent(jpi,jpj),  & 
    97 # endif       
    98109         &      tabspongedone_u  (jpi,jpj),           & 
    99110         &      tabspongedone_v  (jpi,jpj), STAT = ierr(1) ) 
  • NEMO/trunk/src/NST/agrif_oce_interp.F90

    r12377 r13216  
    3434   USE lib_mpp 
    3535   USE vremap 
     36   USE lbclnk 
    3637  
    3738   IMPLICIT NONE 
     
    4445   PUBLIC   interpunb, interpvnb , interpub2b, interpvb2b 
    4546   PUBLIC   interpe3t 
    46 #if defined key_vertical 
    4747   PUBLIC   interpht0, interpmbkt 
    48 # endif 
     48   PUBLIC   agrif_initts, agrif_initssh 
     49 
    4950   INTEGER ::   bdy_tinterp = 0 
    5051 
     
    8990      Agrif_UseSpecialValue = ln_spc_dyn 
    9091      ! 
     92      use_sign_north = .TRUE. 
     93      sign_north = -1. 
    9194      CALL Agrif_Bc_variable( un_interp_id, procname=interpun ) 
    9295      CALL Agrif_Bc_variable( vn_interp_id, procname=interpvn ) 
     96      use_sign_north = .FALSE. 
    9397      ! 
    9498      Agrif_UseSpecialValue = .FALSE. 
    9599      ! 
    96100      ! --- West --- ! 
    97       ibdy1 = 2 
    98       ibdy2 = 1+nbghostcells  
    99       ! 
    100       IF( .NOT.ln_dynspg_ts ) THEN  ! Store transport 
     101      IF( lk_west ) THEN 
     102         ibdy1 = 2 
     103         ibdy2 = 1+nbghostcells  
     104         ! 
     105         IF( .NOT.ln_dynspg_ts ) THEN  ! Store transport 
     106            DO ji = mi0(ibdy1), mi1(ibdy2) 
     107               uu_b(ji,:,Krhs_a) = 0._wp 
     108 
     109               DO jk = 1, jpkm1 
     110                  DO jj = 1, jpj 
     111                     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) 
     112                  END DO 
     113               END DO 
     114 
     115               DO jj = 1, jpj 
     116                  uu_b(ji,jj,Krhs_a) = uu_b(ji,jj,Krhs_a) * r1_hu(ji,jj,Krhs_a) 
     117               END DO 
     118            END DO 
     119         ENDIF 
     120         ! 
    101121         DO ji = mi0(ibdy1), mi1(ibdy2) 
    102             uu_b(ji,:,Krhs_a) = 0._wp 
    103  
     122            zub(ji,:) = 0._wp    ! Correct transport 
    104123            DO jk = 1, jpkm1 
    105124               DO jj = 1, jpj 
    106                   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) 
    107                END DO 
    108             END DO 
    109  
    110             DO jj = 1, jpj 
    111                uu_b(ji,jj,Krhs_a) = uu_b(ji,jj,Krhs_a) * r1_hu(ji,jj,Krhs_a) 
    112             END DO 
    113          END DO 
    114       ENDIF 
    115       ! 
    116       DO ji = mi0(ibdy1), mi1(ibdy2) 
    117          zub(ji,:) = 0._wp    ! Correct transport 
    118          DO jk = 1, jpkm1 
    119             DO jj = 1, jpj 
    120                zub(ji,jj) = zub(ji,jj) &  
    121                   & + e3u(ji,jj,jk,Krhs_a)  * uu(ji,jj,jk,Krhs_a)*umask(ji,jj,jk) 
    122             END DO 
    123          END DO 
    124          DO jj=1,jpj 
    125             zub(ji,jj) = zub(ji,jj) * r1_hu(ji,jj,Krhs_a) 
    126          END DO 
    127              
    128          DO jk = 1, jpkm1 
    129             DO jj = 1, jpj 
    130                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) 
    131             END DO 
    132          END DO 
    133       END DO 
    134              
    135       IF( ln_dynspg_ts ) THEN       ! Set tangential velocities to time splitting estimate 
    136          DO ji = mi0(ibdy1), mi1(ibdy2) 
    137             zvb(ji,:) = 0._wp 
     125                  zub(ji,jj) = zub(ji,jj) &  
     126                     & + e3u(ji,jj,jk,Krhs_a)  * uu(ji,jj,jk,Krhs_a)*umask(ji,jj,jk) 
     127               END DO 
     128            END DO 
     129            DO jj=1,jpj 
     130               zub(ji,jj) = zub(ji,jj) * r1_hu(ji,jj,Krhs_a) 
     131            END DO 
     132                
    138133            DO jk = 1, jpkm1 
    139134               DO jj = 1, jpj 
    140                   zvb(ji,jj) = zvb(ji,jj) + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk) 
    141                END DO 
    142             END DO 
    143             DO jj = 1, jpj 
    144                zvb(ji,jj) = zvb(ji,jj) * r1_hv(ji,jj,Krhs_a) 
    145             END DO 
     135                  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) 
     136               END DO 
     137            END DO 
     138         END DO 
     139                
     140         IF( ln_dynspg_ts ) THEN       ! Set tangential velocities to time splitting estimate 
     141            DO ji = mi0(ibdy1), mi1(ibdy2) 
     142               zvb(ji,:) = 0._wp 
     143               DO jk = 1, jpkm1 
     144                  DO jj = 1, jpj 
     145                     zvb(ji,jj) = zvb(ji,jj) + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk) 
     146                  END DO 
     147               END DO 
     148               DO jj = 1, jpj 
     149                  zvb(ji,jj) = zvb(ji,jj) * r1_hv(ji,jj,Krhs_a) 
     150               END DO 
     151               DO jk = 1, jpkm1 
     152                  DO jj = 1, jpj 
     153                     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) 
     154                  END DO 
     155               END DO 
     156            END DO 
     157         ENDIF 
     158      ENDIF 
     159 
     160      ! --- East --- ! 
     161      IF( lk_east) THEN 
     162         ibdy1 = jpiglo-1-nbghostcells 
     163         ibdy2 = jpiglo-2  
     164         ! 
     165         IF( .NOT.ln_dynspg_ts ) THEN  ! Store transport 
     166            DO ji = mi0(ibdy1), mi1(ibdy2) 
     167               uu_b(ji,:,Krhs_a) = 0._wp 
     168               DO jk = 1, jpkm1 
     169                  DO jj = 1, jpj 
     170                     uu_b(ji,jj,Krhs_a) = uu_b(ji,jj,Krhs_a) &  
     171                         & + e3u(ji,jj,jk,Krhs_a) * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk) 
     172                  END DO 
     173               END DO 
     174               DO jj = 1, jpj 
     175                  uu_b(ji,jj,Krhs_a) = uu_b(ji,jj,Krhs_a) * r1_hu(ji,jj,Krhs_a) 
     176               END DO 
     177            END DO 
     178         ENDIF 
     179         ! 
     180         DO ji = mi0(ibdy1), mi1(ibdy2) 
     181            zub(ji,:) = 0._wp    ! Correct transport 
    146182            DO jk = 1, jpkm1 
    147183               DO jj = 1, jpj 
    148                   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) 
    149                END DO 
    150             END DO 
    151          END DO 
    152       ENDIF 
    153  
    154       ! --- East --- ! 
    155       ibdy1 = jpiglo-1-nbghostcells 
    156       ibdy2 = jpiglo-2  
    157       ! 
    158       IF( .NOT.ln_dynspg_ts ) THEN  ! Store transport 
    159          DO ji = mi0(ibdy1), mi1(ibdy2) 
    160             uu_b(ji,:,Krhs_a) = 0._wp 
     184                  zub(ji,jj) = zub(ji,jj) &  
     185                     & + e3u(ji,jj,jk,Krhs_a)  * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk) 
     186               END DO 
     187            END DO 
     188            DO jj=1,jpj 
     189               zub(ji,jj) = zub(ji,jj) * r1_hu(ji,jj,Krhs_a) 
     190            END DO 
     191                
    161192            DO jk = 1, jpkm1 
    162193               DO jj = 1, jpj 
    163                   uu_b(ji,jj,Krhs_a) = uu_b(ji,jj,Krhs_a) &  
    164                       & + e3u(ji,jj,jk,Krhs_a) * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk) 
    165                END DO 
    166             END DO 
    167             DO jj = 1, jpj 
    168                uu_b(ji,jj,Krhs_a) = uu_b(ji,jj,Krhs_a) * r1_hu(ji,jj,Krhs_a) 
    169             END DO 
    170          END DO 
    171       ENDIF 
    172       ! 
    173       DO ji = mi0(ibdy1), mi1(ibdy2) 
    174          zub(ji,:) = 0._wp    ! Correct transport 
    175          DO jk = 1, jpkm1 
    176             DO jj = 1, jpj 
    177                zub(ji,jj) = zub(ji,jj) &  
    178                   & + e3u(ji,jj,jk,Krhs_a)  * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk) 
    179             END DO 
    180          END DO 
    181          DO jj=1,jpj 
    182             zub(ji,jj) = zub(ji,jj) * r1_hu(ji,jj,Krhs_a) 
    183          END DO 
    184              
    185          DO jk = 1, jpkm1 
    186             DO jj = 1, jpj 
    187                uu(ji,jj,jk,Krhs_a) = ( uu(ji,jj,jk,Krhs_a) &  
    188                  & + uu_b(ji,jj,Krhs_a)-zub(ji,jj))*umask(ji,jj,jk) 
    189             END DO 
    190          END DO 
    191       END DO 
    192              
    193       IF( ln_dynspg_ts ) THEN       ! Set tangential velocities to time splitting estimate 
    194          ibdy1 = jpiglo-nbghostcells 
    195          ibdy2 = jpiglo-1  
    196          DO ji = mi0(ibdy1), mi1(ibdy2) 
    197             zvb(ji,:) = 0._wp 
    198             DO jk = 1, jpkm1 
     194                  uu(ji,jj,jk,Krhs_a) = ( uu(ji,jj,jk,Krhs_a) &  
     195                    & + uu_b(ji,jj,Krhs_a)-zub(ji,jj))*umask(ji,jj,jk) 
     196               END DO 
     197            END DO 
     198         END DO 
     199                
     200         IF( ln_dynspg_ts ) THEN       ! Set tangential velocities to time splitting estimate 
     201            ibdy1 = jpiglo-nbghostcells 
     202            ibdy2 = jpiglo-1  
     203            DO ji = mi0(ibdy1), mi1(ibdy2) 
     204               zvb(ji,:) = 0._wp 
     205               DO jk = 1, jpkm1 
     206                  DO jj = 1, jpj 
     207                     zvb(ji,jj) = zvb(ji,jj) & 
     208                        & + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk) 
     209                  END DO 
     210               END DO 
    199211               DO jj = 1, jpj 
    200                   zvb(ji,jj) = zvb(ji,jj) & 
     212                  zvb(ji,jj) = zvb(ji,jj) * r1_hv(ji,jj,Krhs_a) 
     213               END DO 
     214               DO jk = 1, jpkm1 
     215                  DO jj = 1, jpj 
     216                     vv(ji,jj,jk,Krhs_a) = ( vv(ji,jj,jk,Krhs_a) &  
     217                         & + vv_b(ji,jj,Krhs_a)-zvb(ji,jj)) * vmask(ji,jj,jk) 
     218                  END DO 
     219               END DO 
     220            END DO 
     221         ENDIF 
     222      ENDIF 
     223 
     224      ! --- South --- ! 
     225      IF( lk_south ) THEN 
     226         jbdy1 = 2 
     227         jbdy2 = 1+nbghostcells  
     228         ! 
     229         IF( .NOT.ln_dynspg_ts ) THEN  ! Store transport 
     230            DO jj = mj0(jbdy1), mj1(jbdy2) 
     231               vv_b(:,jj,Krhs_a) = 0._wp 
     232               DO jk = 1, jpkm1 
     233                  DO ji = 1, jpi 
     234                     vv_b(ji,jj,Krhs_a) = vv_b(ji,jj,Krhs_a) &  
     235                         & + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk) 
     236                  END DO 
     237               END DO 
     238               DO ji=1,jpi 
     239                  vv_b(ji,jj,Krhs_a) = vv_b(ji,jj,Krhs_a) * r1_hv(ji,jj,Krhs_a)      
     240               END DO 
     241            END DO 
     242         ENDIF 
     243         ! 
     244         DO jj = mj0(jbdy1), mj1(jbdy2) 
     245            zvb(:,jj) = 0._wp    ! Correct transport 
     246            DO jk=1,jpkm1 
     247               DO ji=1,jpi 
     248                  zvb(ji,jj) = zvb(ji,jj) &  
    201249                     & + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk) 
    202250               END DO 
    203251            END DO 
    204             DO jj = 1, jpj 
     252            DO ji = 1, jpi 
    205253               zvb(ji,jj) = zvb(ji,jj) * r1_hv(ji,jj,Krhs_a) 
    206254            END DO 
    207             DO jk = 1, jpkm1 
    208                DO jj = 1, jpj 
    209                   vv(ji,jj,jk,Krhs_a) = ( vv(ji,jj,jk,Krhs_a) &  
    210                       & + vv_b(ji,jj,Krhs_a)-zvb(ji,jj)) * vmask(ji,jj,jk) 
    211                END DO 
    212             END DO 
    213          END DO 
    214       ENDIF 
    215  
    216       ! --- South --- ! 
    217       jbdy1 = 2 
    218       jbdy2 = 1+nbghostcells  
    219       ! 
    220       IF( .NOT.ln_dynspg_ts ) THEN  ! Store transport 
    221          DO jj = mj0(jbdy1), mj1(jbdy2) 
    222             vv_b(:,jj,Krhs_a) = 0._wp 
     255 
    223256            DO jk = 1, jpkm1 
    224257               DO ji = 1, jpi 
    225                   vv_b(ji,jj,Krhs_a) = vv_b(ji,jj,Krhs_a) &  
    226                       & + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk) 
    227                END DO 
    228             END DO 
    229             DO ji=1,jpi 
    230                vv_b(ji,jj,Krhs_a) = vv_b(ji,jj,Krhs_a) * r1_hv(ji,jj,Krhs_a)      
    231             END DO 
    232          END DO 
    233       ENDIF 
    234       ! 
    235       DO jj = mj0(jbdy1), mj1(jbdy2) 
    236          zvb(:,jj) = 0._wp    ! Correct transport 
    237          DO jk=1,jpkm1 
    238             DO ji=1,jpi 
    239                zvb(ji,jj) = zvb(ji,jj) &  
    240                   & + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk) 
    241             END DO 
    242          END DO 
    243          DO ji = 1, jpi 
    244             zvb(ji,jj) = zvb(ji,jj) * r1_hv(ji,jj,Krhs_a) 
    245          END DO 
    246  
    247          DO jk = 1, jpkm1 
     258                  vv(ji,jj,jk,Krhs_a) = ( vv(ji,jj,jk,Krhs_a) &  
     259                    & + vv_b(ji,jj,Krhs_a) - zvb(ji,jj) ) * vmask(ji,jj,jk) 
     260               END DO 
     261            END DO 
     262         END DO 
     263                
     264         IF( ln_dynspg_ts ) THEN       ! Set tangential velocities to time splitting estimate 
     265            DO jj = mj0(jbdy1), mj1(jbdy2) 
     266               zub(:,jj) = 0._wp 
     267               DO jk = 1, jpkm1 
     268                  DO ji = 1, jpi 
     269                     zub(ji,jj) = zub(ji,jj) &  
     270                        & + e3u(ji,jj,jk,Krhs_a) * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk) 
     271                  END DO 
     272               END DO 
     273               DO ji = 1, jpi 
     274                  zub(ji,jj) = zub(ji,jj) * r1_hu(ji,jj,Krhs_a) 
     275               END DO 
     276                   
     277               DO jk = 1, jpkm1 
     278                  DO ji = 1, jpi 
     279                     uu(ji,jj,jk,Krhs_a) = ( uu(ji,jj,jk,Krhs_a) &  
     280                       & + uu_b(ji,jj,Krhs_a) - zub(ji,jj) ) * umask(ji,jj,jk) 
     281                  END DO 
     282               END DO 
     283            END DO 
     284         ENDIF 
     285      ENDIF 
     286 
     287      ! --- North --- ! 
     288      IF( lk_north ) THEN 
     289         jbdy1 = jpjglo-1-nbghostcells 
     290         jbdy2 = jpjglo-2  
     291         ! 
     292         IF( .NOT.ln_dynspg_ts ) THEN  ! Store transport 
     293            DO jj = mj0(jbdy1), mj1(jbdy2) 
     294               vv_b(:,jj,Krhs_a) = 0._wp 
     295               DO jk = 1, jpkm1 
     296                  DO ji = 1, jpi 
     297                     vv_b(ji,jj,Krhs_a) = vv_b(ji,jj,Krhs_a) &  
     298                         & + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk) 
     299                  END DO 
     300               END DO 
     301               DO ji=1,jpi 
     302                  vv_b(ji,jj,Krhs_a) = vv_b(ji,jj,Krhs_a) * r1_hv(ji,jj,Krhs_a) 
     303               END DO 
     304            END DO 
     305         ENDIF 
     306         ! 
     307         DO jj = mj0(jbdy1), mj1(jbdy2) 
     308            zvb(:,jj) = 0._wp    ! Correct transport 
     309            DO jk=1,jpkm1 
     310               DO ji=1,jpi 
     311                  zvb(ji,jj) = zvb(ji,jj) &  
     312                     & + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk) 
     313               END DO 
     314            END DO 
    248315            DO ji = 1, jpi 
    249                vv(ji,jj,jk,Krhs_a) = ( vv(ji,jj,jk,Krhs_a) &  
    250                  & + vv_b(ji,jj,Krhs_a) - zvb(ji,jj) ) * vmask(ji,jj,jk) 
    251             END DO 
    252          END DO 
    253       END DO 
    254              
    255       IF( ln_dynspg_ts ) THEN       ! Set tangential velocities to time splitting estimate 
    256          DO jj = mj0(jbdy1), mj1(jbdy2) 
    257             zub(:,jj) = 0._wp 
     316               zvb(ji,jj) = zvb(ji,jj) * r1_hv(ji,jj,Krhs_a) 
     317            END DO 
     318 
    258319            DO jk = 1, jpkm1 
    259320               DO ji = 1, jpi 
    260                   zub(ji,jj) = zub(ji,jj) &  
    261                      & + e3u(ji,jj,jk,Krhs_a) * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk) 
    262                END DO 
    263             END DO 
    264             DO ji = 1, jpi 
    265                zub(ji,jj) = zub(ji,jj) * r1_hu(ji,jj,Krhs_a) 
    266             END DO 
     321                  vv(ji,jj,jk,Krhs_a) = ( vv(ji,jj,jk,Krhs_a) &  
     322                    & + vv_b(ji,jj,Krhs_a) - zvb(ji,jj) ) * vmask(ji,jj,jk) 
     323               END DO 
     324            END DO 
     325         END DO 
    267326                
    268             DO jk = 1, jpkm1 
     327         IF( ln_dynspg_ts ) THEN       ! Set tangential velocities to time splitting estimate 
     328            jbdy1 = jpjglo-nbghostcells 
     329            jbdy2 = jpjglo-1 
     330            DO jj = mj0(jbdy1), mj1(jbdy2) 
     331               zub(:,jj) = 0._wp 
     332               DO jk = 1, jpkm1 
     333                  DO ji = 1, jpi 
     334                     zub(ji,jj) = zub(ji,jj) &  
     335                        & + e3u(ji,jj,jk,Krhs_a) * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk) 
     336                  END DO 
     337               END DO 
    269338               DO ji = 1, jpi 
    270                   uu(ji,jj,jk,Krhs_a) = ( uu(ji,jj,jk,Krhs_a) &  
    271                     & + uu_b(ji,jj,Krhs_a) - zub(ji,jj) ) * umask(ji,jj,jk) 
    272                END DO 
    273             END DO 
    274          END DO 
    275       ENDIF 
    276  
    277       ! --- North --- ! 
    278       jbdy1 = jpjglo-1-nbghostcells 
    279       jbdy2 = jpjglo-2  
    280       ! 
    281       IF( .NOT.ln_dynspg_ts ) THEN  ! Store transport 
    282          DO jj = mj0(jbdy1), mj1(jbdy2) 
    283             vv_b(:,jj,Krhs_a) = 0._wp 
    284             DO jk = 1, jpkm1 
    285                DO ji = 1, jpi 
    286                   vv_b(ji,jj,Krhs_a) = vv_b(ji,jj,Krhs_a) &  
    287                       & + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk) 
    288                END DO 
    289             END DO 
    290             DO ji=1,jpi 
    291                vv_b(ji,jj,Krhs_a) = vv_b(ji,jj,Krhs_a) * r1_hv(ji,jj,Krhs_a) 
    292             END DO 
    293          END DO 
    294       ENDIF 
    295       ! 
    296       DO jj = mj0(jbdy1), mj1(jbdy2) 
    297          zvb(:,jj) = 0._wp    ! Correct transport 
    298          DO jk=1,jpkm1 
    299             DO ji=1,jpi 
    300                zvb(ji,jj) = zvb(ji,jj) &  
    301                   & + e3v(ji,jj,jk,Krhs_a) * vv(ji,jj,jk,Krhs_a) * vmask(ji,jj,jk) 
    302             END DO 
    303          END DO 
    304          DO ji = 1, jpi 
    305             zvb(ji,jj) = zvb(ji,jj) * r1_hv(ji,jj,Krhs_a) 
    306          END DO 
    307  
    308          DO jk = 1, jpkm1 
    309             DO ji = 1, jpi 
    310                vv(ji,jj,jk,Krhs_a) = ( vv(ji,jj,jk,Krhs_a) &  
    311                  & + vv_b(ji,jj,Krhs_a) - zvb(ji,jj) ) * vmask(ji,jj,jk) 
    312             END DO 
    313          END DO 
    314       END DO 
    315              
    316       IF( ln_dynspg_ts ) THEN       ! Set tangential velocities to time splitting estimate 
    317          jbdy1 = jpjglo-nbghostcells 
    318          jbdy2 = jpjglo-1 
    319          DO jj = mj0(jbdy1), mj1(jbdy2) 
    320             zub(:,jj) = 0._wp 
    321             DO jk = 1, jpkm1 
    322                DO ji = 1, jpi 
    323                   zub(ji,jj) = zub(ji,jj) &  
    324                      & + e3u(ji,jj,jk,Krhs_a) * uu(ji,jj,jk,Krhs_a) * umask(ji,jj,jk) 
    325                END DO 
    326             END DO 
    327             DO ji = 1, jpi 
    328                zub(ji,jj) = zub(ji,jj) * r1_hu(ji,jj,Krhs_a) 
    329             END DO 
    330                 
    331             DO jk = 1, jpkm1 
    332                DO ji = 1, jpi 
    333                   uu(ji,jj,jk,Krhs_a) = ( uu(ji,jj,jk,Krhs_a) &  
    334                     & + uu_b(ji,jj,Krhs_a) - zub(ji,jj) ) * umask(ji,jj,jk) 
    335                END DO 
    336             END DO 
    337          END DO 
     339                  zub(ji,jj) = zub(ji,jj) * r1_hu(ji,jj,Krhs_a) 
     340               END DO 
     341                   
     342               DO jk = 1, jpkm1 
     343                  DO ji = 1, jpi 
     344                     uu(ji,jj,jk,Krhs_a) = ( uu(ji,jj,jk,Krhs_a) &  
     345                       & + uu_b(ji,jj,Krhs_a) - zub(ji,jj) ) * umask(ji,jj,jk) 
     346                  END DO 
     347               END DO 
     348            END DO 
     349         ENDIF 
    338350      ENDIF 
    339351      ! 
     
    354366      ! 
    355367      !--- West ---! 
    356       istart = 2 
    357       iend   = nbghostcells+1 
    358       DO ji = mi0(istart), mi1(iend) 
    359          DO jj=1,jpj 
    360             va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj) 
    361             ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj) 
    362          END DO 
    363       END DO 
     368      IF( lk_west ) THEN 
     369         istart = 2 
     370         iend   = nbghostcells+1 
     371         DO ji = mi0(istart), mi1(iend) 
     372            DO jj=1,jpj 
     373               va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj) 
     374               ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj) 
     375            END DO 
     376         END DO 
     377      ENDIF 
    364378      ! 
    365379      !--- East ---! 
    366       istart = jpiglo-nbghostcells 
    367       iend   = jpiglo-1 
    368       DO ji = mi0(istart), mi1(iend) 
    369          DO jj=1,jpj 
    370             va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj) 
    371          END DO 
    372       END DO 
    373       istart = jpiglo-nbghostcells-1 
    374       iend   = jpiglo-2 
    375       DO ji = mi0(istart), mi1(iend) 
    376          DO jj=1,jpj 
    377             ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj) 
    378          END DO 
    379       END DO 
     380      IF( lk_east ) THEN 
     381         istart = jpiglo-nbghostcells 
     382         iend   = jpiglo-1 
     383         DO ji = mi0(istart), mi1(iend) 
     384 
     385            DO jj=1,jpj 
     386               va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj) 
     387            END DO 
     388         END DO 
     389         istart = jpiglo-nbghostcells-1 
     390         iend   = jpiglo-2 
     391         DO ji = mi0(istart), mi1(iend) 
     392            DO jj=1,jpj 
     393               ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj) 
     394            END DO 
     395         END DO 
     396      ENDIF  
    380397      ! 
    381398      !--- South ---! 
    382       jstart = 2 
    383       jend   = nbghostcells+1 
    384       DO jj = mj0(jstart), mj1(jend) 
    385          DO ji=1,jpi 
    386             ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj) 
    387             va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj) 
    388          END DO 
    389       END DO 
     399      IF( lk_south ) THEN 
     400         jstart = 2 
     401         jend   = nbghostcells+1 
     402         DO jj = mj0(jstart), mj1(jend) 
     403 
     404            DO ji=1,jpi 
     405               ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj) 
     406               va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj) 
     407            END DO 
     408         END DO 
     409      ENDIF        
    390410      ! 
    391411      !--- North ---! 
    392       jstart = jpjglo-nbghostcells 
    393       jend   = jpjglo-1 
    394       DO jj = mj0(jstart), mj1(jend) 
    395          DO ji=1,jpi 
    396             ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj) 
    397          END DO 
    398       END DO 
    399       jstart = jpjglo-nbghostcells-1 
    400       jend   = jpjglo-2 
    401       DO jj = mj0(jstart), mj1(jend) 
    402          DO ji=1,jpi 
    403             va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj) 
    404          END DO 
    405       END DO 
     412      IF( lk_north ) THEN 
     413         jstart = jpjglo-nbghostcells 
     414         jend   = jpjglo-1 
     415         DO jj = mj0(jstart), mj1(jend) 
     416            DO ji=1,jpi 
     417               ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj) 
     418            END DO 
     419         END DO 
     420         jstart = jpjglo-nbghostcells-1 
     421         jend   = jpjglo-2 
     422         DO jj = mj0(jstart), mj1(jend) 
     423            DO ji=1,jpi 
     424               va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj) 
     425            END DO 
     426         END DO 
     427      ENDIF  
    406428      ! 
    407429   END SUBROUTINE Agrif_dyn_ts 
     
    421443      ! 
    422444      !--- West ---! 
    423       istart = 2 
    424       iend   = nbghostcells+1 
    425       DO ji = mi0(istart), mi1(iend) 
    426          DO jj=1,jpj 
    427             zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj) 
    428             zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj) 
    429          END DO 
    430       END DO 
     445      IF( lk_west ) THEN 
     446         istart = 2 
     447         iend   = nbghostcells+1 
     448         DO ji = mi0(istart), mi1(iend) 
     449            DO jj=1,jpj 
     450               zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj) 
     451               zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj) 
     452            END DO 
     453         END DO 
     454      ENDIF 
    431455      ! 
    432456      !--- East ---! 
    433       istart = jpiglo-nbghostcells 
    434       iend   = jpiglo-1 
    435       DO ji = mi0(istart), mi1(iend) 
    436          DO jj=1,jpj 
    437             zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj) 
    438          END DO 
    439       END DO 
    440       istart = jpiglo-nbghostcells-1 
    441       iend   = jpiglo-2 
    442       DO ji = mi0(istart), mi1(iend) 
    443          DO jj=1,jpj 
    444             zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj) 
    445          END DO 
    446       END DO 
     457      IF( lk_east ) THEN 
     458         istart = jpiglo-nbghostcells 
     459         iend   = jpiglo-1 
     460         DO ji = mi0(istart), mi1(iend) 
     461            DO jj=1,jpj 
     462               zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj) 
     463            END DO 
     464         END DO 
     465         istart = jpiglo-nbghostcells-1 
     466         iend   = jpiglo-2 
     467         DO ji = mi0(istart), mi1(iend) 
     468            DO jj=1,jpj 
     469               zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj) 
     470            END DO 
     471         END DO 
     472      ENDIF 
    447473      ! 
    448474      !--- South ---! 
    449       jstart = 2 
    450       jend   = nbghostcells+1 
    451       DO jj = mj0(jstart), mj1(jend) 
    452          DO ji=1,jpi 
    453             zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj) 
    454             zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj) 
    455          END DO 
    456       END DO 
     475      IF( lk_south ) THEN 
     476         jstart = 2 
     477         jend   = nbghostcells+1 
     478         DO jj = mj0(jstart), mj1(jend) 
     479            DO ji=1,jpi 
     480               zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj) 
     481               zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj) 
     482            END DO 
     483         END DO 
     484      ENDIF 
    457485      ! 
    458486      !--- North ---! 
    459       jstart = jpjglo-nbghostcells 
    460       jend   = jpjglo-1 
    461       DO jj = mj0(jstart), mj1(jend) 
    462          DO ji=1,jpi 
    463             zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj) 
    464          END DO 
    465       END DO 
    466       jstart = jpjglo-nbghostcells-1 
    467       jend   = jpjglo-2 
    468       DO jj = mj0(jstart), mj1(jend) 
    469          DO ji=1,jpi 
    470             zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj) 
    471          END DO 
    472       END DO 
     487      IF( lk_north ) THEN 
     488         jstart = jpjglo-nbghostcells 
     489         jend   = jpjglo-1 
     490         DO jj = mj0(jstart), mj1(jend) 
     491            DO ji=1,jpi 
     492               zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj) 
     493            END DO 
     494         END DO 
     495         jstart = jpjglo-nbghostcells-1 
     496         jend   = jpjglo-2 
     497         DO jj = mj0(jstart), mj1(jend) 
     498            DO ji=1,jpi 
     499               zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj) 
     500            END DO 
     501         END DO 
     502      ENDIF 
    473503      ! 
    474504   END SUBROUTINE Agrif_dyn_ts_flux 
     
    494524      Agrif_SpecialValue = 0._wp 
    495525      Agrif_UseSpecialValue = ln_spc_dyn 
     526 
     527      use_sign_north = .TRUE. 
     528      sign_north = -1. 
     529 
    496530      ! 
    497531      ! Set bdy time interpolation stage to 0 (latter incremented locally do deal with corners) 
     
    518552      ENDIF 
    519553      Agrif_UseSpecialValue = .FALSE. 
     554      use_sign_north = .FALSE. 
    520555      !  
    521556   END SUBROUTINE Agrif_dta_ts 
     
    542577      ! 
    543578      ! --- West --- ! 
    544       istart = 2 
    545       iend   = 1 + nbghostcells 
    546       DO ji = mi0(istart), mi1(iend) 
    547          DO jj = 1, jpj 
    548             ssh(ji,jj,Krhs_a) = hbdy(ji,jj) 
     579      IF(lk_west) THEN 
     580         istart = 2 
     581         iend   = 1 + nbghostcells 
     582         DO ji = mi0(istart), mi1(iend) 
     583            DO jj = 1, jpj 
     584               ssh(ji,jj,Krhs_a) = hbdy(ji,jj) 
     585            ENDDO 
    549586         ENDDO 
    550       ENDDO 
     587      ENDIF 
    551588      ! 
    552589      ! --- East --- ! 
    553       istart = jpiglo - nbghostcells 
    554       iend   = jpiglo - 1 
    555       DO ji = mi0(istart), mi1(iend) 
    556          DO jj = 1, jpj 
    557             ssh(ji,jj,Krhs_a) = hbdy(ji,jj) 
     590      IF(lk_east) THEN 
     591         istart = jpiglo - nbghostcells 
     592         iend   = jpiglo - 1 
     593         DO ji = mi0(istart), mi1(iend) 
     594            DO jj = 1, jpj 
     595               ssh(ji,jj,Krhs_a) = hbdy(ji,jj) 
     596            ENDDO 
    558597         ENDDO 
    559       ENDDO 
     598      ENDIF 
    560599      ! 
    561600      ! --- South --- ! 
    562       jstart = 2 
    563       jend   = 1 + nbghostcells 
    564       DO jj = mj0(jstart), mj1(jend) 
    565          DO ji = 1, jpi 
    566             ssh(ji,jj,Krhs_a) = hbdy(ji,jj) 
     601      IF(lk_south) THEN 
     602         jstart = 2 
     603         jend   = 1 + nbghostcells 
     604         DO jj = mj0(jstart), mj1(jend) 
     605            DO ji = 1, jpi 
     606               ssh(ji,jj,Krhs_a) = hbdy(ji,jj) 
     607            ENDDO 
    567608         ENDDO 
    568       ENDDO 
     609      ENDIF 
    569610      ! 
    570611      ! --- North --- ! 
    571       jstart = jpjglo - nbghostcells 
    572       jend   = jpjglo - 1 
    573       DO jj = mj0(jstart), mj1(jend) 
    574          DO ji = 1, jpi 
    575             ssh(ji,jj,Krhs_a) = hbdy(ji,jj) 
     612      IF(lk_north) THEN 
     613         jstart = jpjglo - nbghostcells 
     614         jend   = jpjglo - 1 
     615         DO jj = mj0(jstart), mj1(jend) 
     616            DO ji = 1, jpi 
     617               ssh(ji,jj,Krhs_a) = hbdy(ji,jj) 
     618            ENDDO 
    576619         ENDDO 
    577       ENDDO 
     620      ENDIF 
    578621      ! 
    579622   END SUBROUTINE Agrif_ssh 
     
    593636      ! 
    594637      ! --- West --- ! 
    595       istart = 2 
    596       iend   = 1+nbghostcells 
    597       DO ji = mi0(istart), mi1(iend) 
    598          DO jj = 1, jpj 
    599             ssha_e(ji,jj) = hbdy(ji,jj) 
     638      IF(lk_west) THEN 
     639         istart = 2 
     640         iend   = 1+nbghostcells 
     641         DO ji = mi0(istart), mi1(iend) 
     642            DO jj = 1, jpj 
     643               ssha_e(ji,jj) = hbdy(ji,jj) 
     644            ENDDO 
    600645         ENDDO 
    601       ENDDO 
     646      ENDIF 
    602647      ! 
    603648      ! --- East --- ! 
    604       istart = jpiglo - nbghostcells 
    605       iend   = jpiglo - 1 
    606       DO ji = mi0(istart), mi1(iend) 
    607          DO jj = 1, jpj 
    608             ssha_e(ji,jj) = hbdy(ji,jj) 
     649      IF(lk_east) THEN 
     650         istart = jpiglo - nbghostcells 
     651         iend   = jpiglo - 1 
     652         DO ji = mi0(istart), mi1(iend) 
     653            DO jj = 1, jpj 
     654               ssha_e(ji,jj) = hbdy(ji,jj) 
     655            ENDDO 
    609656         ENDDO 
    610       ENDDO 
     657      ENDIF 
    611658      ! 
    612659      ! --- South --- ! 
    613       jstart = 2 
    614       jend   = 1+nbghostcells 
    615       DO jj = mj0(jstart), mj1(jend) 
    616          DO ji = 1, jpi 
    617             ssha_e(ji,jj) = hbdy(ji,jj) 
     660      IF(lk_south) THEN 
     661         jstart = 2 
     662         jend   = 1+nbghostcells 
     663         DO jj = mj0(jstart), mj1(jend) 
     664            DO ji = 1, jpi 
     665               ssha_e(ji,jj) = hbdy(ji,jj) 
     666            ENDDO 
    618667         ENDDO 
    619       ENDDO 
     668      ENDIF 
    620669      ! 
    621670      ! --- North --- ! 
    622       jstart = jpjglo - nbghostcells 
    623       jend   = jpjglo - 1 
    624       DO jj = mj0(jstart), mj1(jend) 
    625          DO ji = 1, jpi 
    626             ssha_e(ji,jj) = hbdy(ji,jj) 
     671      IF(lk_north) THEN 
     672         jstart = jpjglo - nbghostcells 
     673         jend   = jpjglo - 1 
     674         DO jj = mj0(jstart), mj1(jend) 
     675            DO ji = 1, jpi 
     676               ssha_e(ji,jj) = hbdy(ji,jj) 
     677            ENDDO 
    627678         ENDDO 
    628       ENDDO 
     679      ENDIF 
    629680      ! 
    630681   END SUBROUTINE Agrif_ssh_ts 
     
    662713      INTEGER  ::   ji, jj, jk, jn  ! dummy loop indices 
    663714      INTEGER  ::   N_in, N_out 
     715      INTEGER  :: item 
    664716      ! vertical interpolation: 
    665717      REAL(wp) :: zhtot 
    666718      REAL(wp), DIMENSION(k1:k2,1:jpts) :: tabin 
    667       REAL(wp), DIMENSION(k1:k2) :: h_in 
    668       REAL(wp), DIMENSION(1:jpk) :: h_out 
    669       !!---------------------------------------------------------------------- 
    670  
    671       IF( before ) THEN          
     719      REAL(wp), DIMENSION(k1:k2) :: h_in, z_in 
     720      REAL(wp), DIMENSION(1:jpk) :: h_out, z_out 
     721      !!---------------------------------------------------------------------- 
     722 
     723      IF( before ) THEN 
     724 
     725         item = Kmm_a 
     726         IF( l_ini_child )   Kmm_a = Kbb_a   
     727 
    672728         DO jn = 1,jpts 
    673729            DO jk=k1,k2 
     
    678734              END DO 
    679735           END DO 
    680         END DO 
    681  
    682 # if defined key_vertical 
    683         ! Interpolate thicknesses 
    684         ! Warning: these are masked, hence extrapolated prior interpolation. 
    685         DO jk=k1,k2 
    686            DO jj=j1,j2 
    687               DO ji=i1,i2 
    688                   ptab(ji,jj,jk,jpts+1) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kmm_a) 
    689               END DO 
    690            END DO 
    691         END DO 
    692  
    693         ! Extrapolate thicknesses in partial bottom cells: 
    694         ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on 
    695         IF (ln_zps) THEN 
    696            DO jj=j1,j2 
    697               DO ji=i1,i2 
    698                   jk = mbkt(ji,jj) 
    699                   ptab(ji,jj,jk,jpts+1) = 0._wp 
    700               END DO 
    701            END DO            
    702         END IF 
    703       
    704         ! Save ssh at last level: 
    705         IF (.NOT.ln_linssh) THEN 
    706            ptab(i1:i2,j1:j2,k2,jpts+1) = ssh(i1:i2,j1:j2,Kmm_a)*tmask(i1:i2,j1:j2,1)  
    707         ELSE 
    708            ptab(i1:i2,j1:j2,k2,jpts+1) = 0._wp 
    709         END IF       
    710 # endif 
     736         END DO 
     737 
     738         IF( l_vremap .OR. l_ini_child) THEN 
     739            ! Interpolate thicknesses 
     740            ! Warning: these are masked, hence extrapolated prior interpolation. 
     741            DO jk=k1,k2 
     742               DO jj=j1,j2 
     743                  DO ji=i1,i2 
     744                      ptab(ji,jj,jk,jpts+1) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kmm_a) 
     745 
     746                  END DO 
     747               END DO 
     748            END DO 
     749 
     750            ! Extrapolate thicknesses in partial bottom cells: 
     751            ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on 
     752            IF (ln_zps) THEN 
     753               DO jj=j1,j2 
     754                  DO ji=i1,i2 
     755                      jk = mbkt(ji,jj) 
     756                      ptab(ji,jj,jk,jpts+1) = 0._wp 
     757                  END DO 
     758               END DO            
     759            END IF 
     760         
     761            ! Save ssh at last level: 
     762            IF (.NOT.ln_linssh) THEN 
     763               ptab(i1:i2,j1:j2,k2,jpts+1) = ssh(i1:i2,j1:j2,Kmm_a)*tmask(i1:i2,j1:j2,1)  
     764            ELSE 
     765               ptab(i1:i2,j1:j2,k2,jpts+1) = 0._wp 
     766            END IF       
     767         ENDIF 
     768         Kmm_a = item 
     769 
    711770      ELSE  
    712  
    713 # if defined key_vertical  
    714          IF (ln_linssh) ptab(i1:i2,j1:j2,k2,n2) = 0._wp  
    715              
    716          DO jj=j1,j2 
    717             DO ji=i1,i2 
    718                ts(ji,jj,:,:,Krhs_a) = 0._wp 
    719                N_in = mbkt_parent(ji,jj) 
    720                zhtot = 0._wp 
    721                DO jk=1,N_in !k2 = jpk of parent grid 
    722                   IF (jk==N_in) THEN 
    723                      h_in(jk) = ht0_parent(ji,jj) + ptab(ji,jj,k2,n2) - zhtot 
    724                   ELSE 
    725                      h_in(jk) = ptab(ji,jj,jk,n2) 
     771         item = Krhs_a 
     772         IF( l_ini_child )   Krhs_a = Kbb_a   
     773 
     774         IF( l_vremap .OR. l_ini_child ) THEN 
     775            IF (ln_linssh) ptab(i1:i2,j1:j2,k2,n2) = 0._wp  
     776                
     777            DO jj=j1,j2 
     778               DO ji=i1,i2 
     779                  ts(ji,jj,:,:,Krhs_a) = 0.                   
     780               !   IF( l_ini_child) ts(ji,jj,:,:,Krhs_a) = ptab(ji,jj,:,1:jpts) 
     781                  N_in = mbkt_parent(ji,jj) 
     782                  zhtot = 0._wp 
     783                  DO jk=1,N_in !k2 = jpk of parent grid 
     784                     IF (jk==N_in) THEN 
     785                        h_in(jk) = ht0_parent(ji,jj) + ptab(ji,jj,k2,n2) - zhtot 
     786                     ELSE 
     787                        h_in(jk) = ptab(ji,jj,jk,n2) 
     788                     ENDIF 
     789                     zhtot = zhtot + h_in(jk) 
     790                     tabin(jk,:) = ptab(ji,jj,jk,n1:n2-1) 
     791                  END DO 
     792                  z_in(1) = 0.5_wp * h_in(1) - zhtot + ht0_parent(ji,jj) 
     793                  DO jk=2,N_in 
     794                     z_in(jk) = z_in(jk-1) + 0.5_wp * h_in(jk) 
     795                  ENDDO 
     796 
     797                  N_out = 0 
     798                  DO jk=1,jpk ! jpk of child grid 
     799                     IF (tmask(ji,jj,jk) == 0._wp) EXIT  
     800                     N_out = N_out + 1 
     801                     h_out(jk) = e3t(ji,jj,jk,Krhs_a) 
     802                  ENDDO 
     803 
     804                  z_out(1) = 0.5_wp * h_out(1) - SUM(h_out(1:N_out)) + ht_0(ji,jj) 
     805                  DO jk=2,N_out 
     806                     z_out(jk) = z_out(jk-1) + 0.5_wp * h_out(jk) 
     807                  ENDDO 
     808 
     809                  IF (N_in*N_out > 0) THEN 
     810                     IF( l_ini_child ) THEN 
     811                        CALL remap_linear(tabin(1:N_in,1:jpts),z_in(1:N_in),ts(ji,jj,1:N_out,1:jpts,Krhs_a),          & 
     812                                      &   z_out(1:N_out),N_in,N_out,jpts)   
     813                     ELSE  
     814                        CALL reconstructandremap(tabin(1:N_in,1:jpts),h_in(1:N_in),ts(ji,jj,1:N_out,1:jpts,Krhs_a),   & 
     815                                      &   h_out(1:N_out),N_in,N_out,jpts)   
     816                     ENDIF 
    726817                  ENDIF 
    727                   zhtot = zhtot + h_in(jk) 
    728                   tabin(jk,:) = ptab(ji,jj,jk,n1:n2-1) 
    729                END DO 
    730                N_out = 0 
    731                DO jk=1,jpk ! jpk of child grid 
    732                   IF (tmask(ji,jj,jk) == 0._wp) EXIT  
    733                   N_out = N_out + 1 
    734                   h_out(jk) = e3t(ji,jj,jk,Krhs_a) 
    735818               ENDDO 
    736                IF (N_in*N_out > 0) THEN 
    737                   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) 
    738                ENDIF 
    739819            ENDDO 
    740          ENDDO 
    741 # else 
    742          ! 
    743          DO jn=1, jpts 
    744             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)  
    745          END DO 
    746 # endif 
     820            Krhs_a = item 
     821  
     822         ELSE 
     823          
     824            DO jn=1, jpts 
     825                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)  
     826            END DO 
     827         ENDIF 
    747828 
    748829      ENDIF 
     
    780861      REAL(wp) :: zrhoy, zhtot 
    781862      ! vertical interpolation: 
    782       REAL(wp), DIMENSION(k1:k2) :: tabin, h_in 
    783       REAL(wp), DIMENSION(1:jpk) :: h_out 
    784       INTEGER  :: N_in, N_out 
     863      REAL(wp), DIMENSION(k1:k2) :: tabin, h_in, z_in 
     864      REAL(wp), DIMENSION(1:jpk) :: h_out, z_out 
     865      INTEGER  :: N_in, N_out,item 
    785866      REAL(wp) :: h_diff 
    786867      !!---------------------------------------------     
    787868      ! 
    788869      IF (before) THEN  
     870 
     871         item = Kmm_a 
     872         IF( l_ini_child )   Kmm_a = Kbb_a      
     873 
    789874         DO jk=1,jpk 
    790875            DO jj=j1,j2 
    791876               DO ji=i1,i2 
    792877                  ptab(ji,jj,jk,1) = (e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a) * uu(ji,jj,jk,Kmm_a)*umask(ji,jj,jk))  
    793 # if defined key_vertical 
    794                   ! Interpolate thicknesses (masked for subsequent extrapolation) 
    795                   ptab(ji,jj,jk,2) = umask(ji,jj,jk) * e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a) 
    796 # endif 
    797                END DO 
    798             END DO 
    799          END DO 
    800 # if defined key_vertical 
     878                  IF( l_vremap .OR. l_ini_child) THEN 
     879                     ! Interpolate thicknesses (masked for subsequent extrapolation) 
     880                     ptab(ji,jj,jk,2) = umask(ji,jj,jk) * e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a) 
     881                  ENDIF 
     882               END DO 
     883            END DO 
     884         END DO 
     885 
     886        IF( l_vremap .OR. l_ini_child) THEN 
    801887         ! Extrapolate thicknesses in partial bottom cells: 
    802888         ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on 
    803          IF (ln_zps) THEN 
    804             DO jj=j1,j2 
    805                DO ji=i1,i2 
    806                   jk = mbku(ji,jj) 
    807                   ptab(ji,jj,jk,2) = 0._wp 
    808                END DO 
    809             END DO            
    810          END IF 
    811         ! Save ssh at last level: 
    812         ptab(i1:i2,j1:j2,k2,2) = 0._wp 
    813         IF (.NOT.ln_linssh) THEN 
    814            ! This vertical sum below should be replaced by the sea-level at U-points (optimization): 
    815            DO jk=1,jpk 
    816               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) 
    817            END DO 
    818            ptab(i1:i2,j1:j2,k2,2) = ptab(i1:i2,j1:j2,k2,2) - hu_0(i1:i2,j1:j2) 
    819         END IF  
    820 # endif 
     889            IF (ln_zps) THEN 
     890               DO jj=j1,j2 
     891                  DO ji=i1,i2 
     892                     jk = mbku(ji,jj) 
     893                     ptab(ji,jj,jk,2) = 0._wp 
     894                  END DO 
     895               END DO            
     896            END IF 
     897 
     898           ! Save ssh at last level: 
     899           ptab(i1:i2,j1:j2,k2,2) = 0._wp 
     900           IF (.NOT.ln_linssh) THEN 
     901              ! This vertical sum below should be replaced by the sea-level at U-points (optimization): 
     902              DO jk=1,jpk 
     903                 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) 
     904              END DO 
     905              ptab(i1:i2,j1:j2,k2,2) = ptab(i1:i2,j1:j2,k2,2) - hu_0(i1:i2,j1:j2) 
     906           END IF 
     907        ENDIF 
     908 
     909         Kmm_a = item 
    821910         ! 
    822911      ELSE 
    823912         zrhoy = Agrif_rhoy() 
    824 # if defined key_vertical 
     913 
     914        IF( l_vremap .OR. l_ini_child) THEN 
    825915! VERTICAL REFINEMENT BEGIN 
    826916 
    827          IF (ln_linssh) ptab(i1:i2,j1:j2,k2,2) = 0._wp  
    828  
    829          DO ji=i1,i2 
    830             DO jj=j1,j2 
    831                uu(ji,jj,:,Krhs_a) = 0._wp 
    832                N_in = mbku_parent(ji,jj) 
    833                zhtot = 0._wp 
    834                DO jk=1,N_in 
    835                   IF (jk==N_in) THEN 
    836                      h_in(jk) = hu0_parent(ji,jj) + ptab(ji,jj,k2,2) - zhtot 
    837                   ELSE 
    838                      h_in(jk) = ptab(ji,jj,jk,2)/(e2u(ji,jj)*zrhoy)  
    839                   ENDIF 
    840                   zhtot = zhtot + h_in(jk) 
    841                   tabin(jk) = ptab(ji,jj,jk,1)/(e2u(ji,jj)*zrhoy*h_in(jk)) 
    842               ENDDO 
    843                    
    844               N_out = 0 
    845               DO jk=1,jpk 
    846                  if (umask(ji,jj,jk) == 0) EXIT 
    847                  N_out = N_out + 1 
    848                  h_out(N_out) = e3u(ji,jj,jk,Krhs_a) 
    849               ENDDO 
    850               IF (N_in*N_out > 0) THEN 
    851                  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) 
    852               ENDIF 
     917            IF (ln_linssh) ptab(i1:i2,j1:j2,k2,2) = 0._wp  
     918 
     919            DO ji=i1,i2 
     920               DO jj=j1,j2 
     921                  uu(ji,jj,:,Krhs_a) = 0._wp 
     922                  N_in = mbku_parent(ji,jj) 
     923                  zhtot = 0._wp 
     924                  DO jk=1,N_in 
     925                     IF (jk==N_in) THEN 
     926                        h_in(jk) = hu0_parent(ji,jj) + ptab(ji,jj,k2,2) - zhtot 
     927                     ELSE 
     928                        h_in(jk) = ptab(ji,jj,jk,2)/(e2u(ji,jj)*zrhoy)  
     929                     ENDIF 
     930                     zhtot = zhtot + h_in(jk) 
     931                     IF( h_in(jk) .GT. 0. ) THEN 
     932                     tabin(jk) = ptab(ji,jj,jk,1)/(e2u(ji,jj)*zrhoy*h_in(jk)) 
     933                     ELSE 
     934                     tabin(jk) = 0. 
     935                     ENDIF 
     936                 ENDDO 
     937                 z_in(1) = 0.5_wp * h_in(1) - zhtot + hu0_parent(ji,jj)  
     938                 DO jk=2,N_in 
     939                    z_in(jk) = z_in(jk-1) + 0.5_wp * h_in(jk) 
     940                 ENDDO 
     941                      
     942                 N_out = 0 
     943                 DO jk=1,jpk 
     944                    IF (umask(ji,jj,jk) == 0) EXIT 
     945                    N_out = N_out + 1 
     946                    h_out(N_out) = e3u(ji,jj,jk,Krhs_a) 
     947                 ENDDO 
     948 
     949                 z_out(1) = 0.5_wp * h_out(1) - SUM(h_out(1:N_out)) + hu_0(ji,jj) 
     950                 DO jk=2,N_out 
     951                    z_out(jk) = z_out(jk-1) + 0.5_wp * h_out(jk)  
     952                 ENDDO   
     953 
     954                 IF (N_in*N_out > 0) THEN 
     955                     IF( l_ini_child ) THEN 
     956                        CALL remap_linear       (tabin(1:N_in),z_in(1:N_in),uu(ji,jj,1:N_out,Krhs_a),z_out(1:N_out),N_in,N_out,1) 
     957                     ELSE 
     958                        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) 
     959                     ENDIF    
     960                 ENDIF 
     961               ENDDO 
    853962            ENDDO 
    854          ENDDO 
    855  
    856 # else 
    857          DO jk = 1, jpkm1 
    858             DO jj=j1,j2 
    859                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) ) 
    860             END DO 
    861          END DO 
    862 # endif 
     963         ELSE 
     964            DO jk = 1, jpkm1 
     965               DO jj=j1,j2 
     966                  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) ) 
     967               END DO 
     968            END DO 
     969         ENDIF 
    863970 
    864971      ENDIF 
     
    878985      REAL(wp) :: zrhox 
    879986      ! vertical interpolation: 
    880       REAL(wp), DIMENSION(k1:k2) :: tabin, h_in 
    881       REAL(wp), DIMENSION(1:jpk) :: h_out 
    882       INTEGER  :: N_in, N_out 
     987      REAL(wp), DIMENSION(k1:k2) :: tabin, h_in, z_in 
     988      REAL(wp), DIMENSION(1:jpk) :: h_out, z_out 
     989      INTEGER  :: N_in, N_out, item 
    883990      REAL(wp) :: h_diff, zhtot 
    884991      !!---------------------------------------------     
    885992      !       
    886       IF (before) THEN           
     993      IF (before) THEN    
     994 
     995         item = Kmm_a 
     996         IF( l_ini_child )   Kmm_a = Kbb_a      
     997        
    887998         DO jk=k1,k2 
    888999            DO jj=j1,j2 
    8891000               DO ji=i1,i2 
    8901001                  ptab(ji,jj,jk,1) = (e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a) * vv(ji,jj,jk,Kmm_a)*vmask(ji,jj,jk)) 
    891 # if defined key_vertical 
    892                   ! Interpolate thicknesses (masked for subsequent extrapolation) 
    893                   ptab(ji,jj,jk,2) = vmask(ji,jj,jk) * e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a) 
    894 # endif 
    895                END DO 
    896             END DO 
    897          END DO 
    898 # if defined key_vertical 
     1002                  IF( l_vremap .OR. l_ini_child) THEN 
     1003                     ! Interpolate thicknesses (masked for subsequent extrapolation) 
     1004                     ptab(ji,jj,jk,2) = vmask(ji,jj,jk) * e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a) 
     1005                  ENDIF 
     1006               END DO 
     1007            END DO 
     1008         END DO 
     1009 
     1010         IF( l_vremap .OR. l_ini_child) THEN 
    8991011         ! Extrapolate thicknesses in partial bottom cells: 
    9001012         ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on 
    901          IF (ln_zps) THEN 
     1013            IF (ln_zps) THEN 
     1014               DO jj=j1,j2 
     1015                  DO ji=i1,i2 
     1016                     jk = mbkv(ji,jj) 
     1017                     ptab(ji,jj,jk,2) = 0._wp 
     1018                  END DO 
     1019               END DO            
     1020            END IF 
     1021            ! Save ssh at last level: 
     1022            ptab(i1:i2,j1:j2,k2,2) = 0._wp 
     1023            IF (.NOT.ln_linssh) THEN 
     1024               ! This vertical sum below should be replaced by the sea-level at V-points (optimization): 
     1025               DO jk=1,jpk 
     1026                  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) 
     1027               END DO 
     1028               ptab(i1:i2,j1:j2,k2,2) = ptab(i1:i2,j1:j2,k2,2) - hv_0(i1:i2,j1:j2) 
     1029            END IF  
     1030         ENDIF 
     1031         item = Kmm_a 
     1032 
     1033      ELSE        
     1034         zrhox = Agrif_rhox() 
     1035 
     1036         IF( l_vremap .OR. l_ini_child ) THEN 
     1037 
     1038            IF (ln_linssh) ptab(i1:i2,j1:j2,k2,2) = 0._wp  
     1039 
    9021040            DO jj=j1,j2 
    9031041               DO ji=i1,i2 
    904                   jk = mbkv(ji,jj) 
    905                   ptab(ji,jj,jk,2) = 0._wp 
    906                END DO 
    907             END DO            
    908          END IF 
    909         ! Save ssh at last level: 
    910         ptab(i1:i2,j1:j2,k2,2) = 0._wp 
    911         IF (.NOT.ln_linssh) THEN 
    912            ! This vertical sum below should be replaced by the sea-level at V-points (optimization): 
    913            DO jk=1,jpk 
    914               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) 
    915            END DO 
    916            ptab(i1:i2,j1:j2,k2,2) = ptab(i1:i2,j1:j2,k2,2) - hv_0(i1:i2,j1:j2) 
    917         END IF  
    918 # endif 
    919       ELSE        
    920          zrhox = Agrif_rhox() 
    921 # if defined key_vertical 
    922  
    923          IF (ln_linssh) ptab(i1:i2,j1:j2,k2,2) = 0._wp  
    924  
    925          DO jj=j1,j2 
    926             DO ji=i1,i2 
    927                vv(ji,jj,:,Krhs_a) = 0._wp 
    928                N_in = mbkv_parent(ji,jj) 
    929                zhtot = 0._wp 
    930                DO jk=1,N_in 
    931                   IF (jk==N_in) THEN 
    932                      h_in(jk) = hv0_parent(ji,jj) + ptab(ji,jj,k2,2) - zhtot 
    933                   ELSE 
    934                      h_in(jk) = ptab(ji,jj,jk,2)/(e1v(ji,jj)*zrhox)  
     1042                  vv(ji,jj,:,Krhs_a) = 0._wp 
     1043                  N_in = mbkv_parent(ji,jj) 
     1044                  zhtot = 0._wp 
     1045                  DO jk=1,N_in 
     1046                     IF (jk==N_in) THEN 
     1047                        h_in(jk) = hv0_parent(ji,jj) + ptab(ji,jj,k2,2) - zhtot 
     1048                     ELSE 
     1049                        h_in(jk) = ptab(ji,jj,jk,2)/(e1v(ji,jj)*zrhox)  
     1050                     ENDIF 
     1051                     zhtot = zhtot + h_in(jk) 
     1052                     IF( h_in(jk) .GT. 0. ) THEN 
     1053                       tabin(jk) = ptab(ji,jj,jk,1)/(e1v(ji,jj)*zrhox*h_in(jk)) 
     1054                     ELSE 
     1055                       tabin(jk)  = 0. 
     1056                     ENDIF  
     1057                  ENDDO 
     1058 
     1059                  z_in(1) = 0.5_wp * h_in(1) - zhtot + hv0_parent(ji,jj) 
     1060                  DO jk=2,N_in 
     1061                     z_in(jk) = z_in(jk-1) + 0.5_wp * h_in(jk) 
     1062                  ENDDO 
     1063 
     1064                  N_out = 0 
     1065                  DO jk=1,jpk 
     1066                     IF (vmask(ji,jj,jk) == 0) EXIT 
     1067                     N_out = N_out + 1 
     1068                     h_out(N_out) = e3v(ji,jj,jk,Krhs_a) 
     1069                  ENDDO 
     1070 
     1071                  z_out(1) = 0.5_wp * h_out(1) - SUM(h_out(1:N_out)) + hv_0(ji,jj) 
     1072                  DO jk=2,N_out 
     1073                     z_out(jk) = z_out(jk-1) + 0.5_wp * h_out(jk) 
     1074                  ENDDO 
     1075  
     1076                  IF (N_in*N_out > 0) THEN 
     1077                     IF( l_ini_child ) THEN 
     1078                        CALL remap_linear       (tabin(1:N_in),z_in(1:N_in),vv(ji,jj,1:N_out,Krhs_a),z_out(1:N_out),N_in,N_out,1) 
     1079                     ELSE 
     1080                        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) 
     1081                     ENDIF    
    9351082                  ENDIF 
    936                   zhtot = zhtot + h_in(jk) 
    937                   tabin(jk) = ptab(ji,jj,jk,1)/(e1v(ji,jj)*zrhox*h_in(jk)) 
    938               ENDDO 
    939           
    940                N_out = 0 
    941                DO jk=1,jpk 
    942                   if (vmask(ji,jj,jk) == 0) EXIT 
    943                   N_out = N_out + 1 
    944                   h_out(N_out) = e3v(ji,jj,jk,Krhs_a) 
    945                END DO 
    946                IF (N_in*N_out > 0) THEN 
    947                   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) 
    948                ENDIF 
    949             END DO 
    950          END DO 
    951 # else 
    952          DO jk = 1, jpkm1 
    953             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) ) 
    954          END DO 
    955 # endif 
     1083               END DO 
     1084            END DO 
     1085         ELSE 
     1086            DO jk = 1, jpkm1 
     1087               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) ) 
     1088            END DO 
     1089         ENDIF 
    9561090      ENDIF 
    9571091      !         
     
    11631297   END SUBROUTINE interpe3t 
    11641298 
    1165  
    11661299   SUBROUTINE interpavm( ptab, i1, i2, j1, j2, k1, k2, m1, m2, before ) 
    11671300      !!---------------------------------------------------------------------- 
     
    11851318              END DO 
    11861319           END DO 
    1187         END DO 
    1188  
    1189 # if defined key_vertical 
    1190         ! Interpolate thicknesses 
    1191         ! Warning: these are masked, hence extrapolated prior interpolation. 
    1192         DO jk=k1,k2 
    1193            DO jj=j1,j2 
    1194               DO ji=i1,i2 
    1195                   ptab(ji,jj,jk,2) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kmm_a) 
    1196               END DO 
    1197            END DO 
    1198         END DO 
    1199  
    1200         ! Extrapolate thicknesses in partial bottom cells: 
    1201         ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on 
    1202         IF (ln_zps) THEN 
    1203            DO jj=j1,j2 
    1204               DO ji=i1,i2 
    1205                   jk = mbkt(ji,jj) 
    1206                   ptab(ji,jj,jk,2) = 0._wp 
    1207               END DO 
    1208            END DO            
    1209         END IF 
    1210       
    1211         ! Save ssh at last level: 
    1212         IF (.NOT.ln_linssh) THEN 
    1213            ptab(i1:i2,j1:j2,k2,2) = ssh(i1:i2,j1:j2,Kmm_a)*tmask(i1:i2,j1:j2,1)  
    1214         ELSE 
    1215            ptab(i1:i2,j1:j2,k2,2) = 0._wp 
    1216         END IF       
    1217 # endif 
     1320         END DO 
     1321 
     1322         IF( l_vremap ) THEN 
     1323            ! Interpolate thicknesses 
     1324            ! Warning: these are masked, hence extrapolated prior interpolation. 
     1325            DO jk=k1,k2 
     1326               DO jj=j1,j2 
     1327                  DO ji=i1,i2 
     1328                      ptab(ji,jj,jk,2) = tmask(ji,jj,jk) * e3t(ji,jj,jk,Kmm_a) 
     1329                  END DO 
     1330               END DO 
     1331            END DO 
     1332 
     1333            ! Extrapolate thicknesses in partial bottom cells: 
     1334            ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on 
     1335            IF (ln_zps) THEN 
     1336               DO jj=j1,j2 
     1337                  DO ji=i1,i2 
     1338                      jk = mbkt(ji,jj) 
     1339                      ptab(ji,jj,jk,2) = 0._wp 
     1340                  END DO 
     1341               END DO            
     1342            END IF 
     1343         
     1344           ! Save ssh at last level: 
     1345            IF (.NOT.ln_linssh) THEN 
     1346               ptab(i1:i2,j1:j2,k2,2) = ssh(i1:i2,j1:j2,Kmm_a)*tmask(i1:i2,j1:j2,1)  
     1347            ELSE 
     1348               ptab(i1:i2,j1:j2,k2,2) = 0._wp 
     1349            END IF       
     1350          ENDIF 
     1351 
    12181352      ELSE  
    1219 #ifdef key_vertical          
    1220          IF (ln_linssh) ptab(i1:i2,j1:j2,k2,2) = 0._wp  
    1221          avm_k(i1:i2,j1:j2,k1:k2) = 0._wp 
    1222              
    1223          DO jj = j1, j2 
    1224             DO ji =i1, i2 
    1225                N_in = mbkt_parent(ji,jj) 
    1226                IF ( tmask(ji,jj,1) == 0._wp) N_in = 0 
    1227                z_in(N_in+1) = ht0_parent(ji,jj) + ptab(ji,jj,k2,2) 
    1228                DO jk = N_in, 1, -1  ! Parent vertical grid                
    1229                      z_in(jk) = z_in(jk+1) - ptab(ji,jj,jk,2) 
    1230                     tabin(jk) = ptab(ji,jj,jk,1) 
    1231                END DO 
    1232                N_out = mbkt(ji,jj)  
    1233                DO jk = 1, N_out        ! Child vertical grid 
    1234                   z_out(jk) = gdepw(ji,jj,jk,Kmm_a) 
     1353 
     1354         IF( l_vremap ) THEN 
     1355            IF (ln_linssh) ptab(i1:i2,j1:j2,k2,2) = 0._wp  
     1356            avm_k(i1:i2,j1:j2,k1:k2) = 0._wp 
     1357                
     1358            DO jj = j1, j2 
     1359               DO ji =i1, i2 
     1360                  N_in = mbkt_parent(ji,jj) 
     1361                  IF ( tmask(ji,jj,1) == 0._wp) N_in = 0 
     1362                  z_in(N_in+1) = ht0_parent(ji,jj) + ptab(ji,jj,k2,2) 
     1363                  DO jk = N_in, 1, -1  ! Parent vertical grid                
     1364                        z_in(jk) = z_in(jk+1) - ptab(ji,jj,jk,2) 
     1365                       tabin(jk) = ptab(ji,jj,jk,1) 
     1366                  END DO 
     1367                  N_out = mbkt(ji,jj)  
     1368                  DO jk = 1, N_out        ! Child vertical grid 
     1369                     z_out(jk) = gdepw(ji,jj,jk,Kmm_a) 
     1370                  ENDDO 
     1371                  IF (N_in*N_out > 0) THEN 
     1372                     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) 
     1373                  ENDIF 
    12351374               ENDDO 
    1236                IF (N_in*N_out > 0) THEN 
    1237                   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) 
    1238                ENDIF 
    12391375            ENDDO 
    1240          ENDDO 
    1241 #else 
    1242          avm_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2,1) 
    1243 #endif 
     1376         ELSE 
     1377            avm_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2,1) 
     1378         ENDIF 
    12441379      ENDIF 
    12451380      ! 
    12461381   END SUBROUTINE interpavm 
    12471382 
    1248 # if defined key_vertical 
    12491383   SUBROUTINE interpmbkt( ptab, i1, i2, j1, j2, before ) 
    12501384      !!---------------------------------------------------------------------- 
     
    12821416      ! 
    12831417   END SUBROUTINE interpht0 
    1284 #endif 
    1285  
     1418 
     1419   SUBROUTINE agrif_initts(tabres,i1,i2,j1,j2,k1,k2,m1,m2,before) 
     1420       INTEGER :: i1, i2, j1, j2, k1, k2, m1, m2 
     1421       REAL(wp):: tabres(i1:i2,j1:j2,k1:k2,m1:m2) 
     1422       LOGICAL :: before 
     1423 
     1424       INTEGER :: jm 
     1425 
     1426       IF (before) THEN 
     1427         DO jm=1,jpts 
     1428             tabres(i1:i2,j1:j2,k1:k2,jm) = ts(i1:i2,j1:j2,k1:k2,jm,Kbb_a) 
     1429         END DO 
     1430       ELSE 
     1431         DO jm=1,jpts 
     1432             ts(i1:i2,j1:j2,k1:k2,jm,Kbb_a)=tabres(i1:i2,j1:j2,k1:k2,jm) 
     1433         END DO 
     1434       ENDIF 
     1435   END SUBROUTINE agrif_initts  
     1436 
     1437   SUBROUTINE agrif_initssh( ptab, i1, i2, j1, j2, before ) 
     1438      !!---------------------------------------------------------------------- 
     1439      !!                  ***  ROUTINE interpsshn  *** 
     1440      !!----------------------------------------------------------------------   
     1441      INTEGER                         , INTENT(in   ) ::   i1, i2, j1, j2 
     1442      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) ::   ptab 
     1443      LOGICAL                         , INTENT(in   ) ::   before 
     1444      ! 
     1445      !!----------------------------------------------------------------------   
     1446      ! 
     1447      IF( before) THEN 
     1448         ptab(i1:i2,j1:j2) = ssh(i1:i2,j1:j2,Kbb_a) 
     1449      ELSE 
     1450         ssh(i1:i2,j1:j2,Kbb_a) = ptab(i1:i2,j1:j2)*tmask(i1:i2,j1:j2,1) 
     1451      ENDIF 
     1452      ! 
     1453   END SUBROUTINE agrif_initssh 
     1454    
    12861455#else 
    12871456   !!---------------------------------------------------------------------- 
  • NEMO/trunk/src/NST/agrif_oce_sponge.F90

    r12489 r13216  
    8080      Agrif_SpecialValue=0. 
    8181      Agrif_UseSpecialValue = ln_spc_dyn 
     82      use_sign_north = .TRUE. 
     83      sign_north = -1. 
    8284      ! 
    8385      tabspongedone_u = .FALSE. 
     
    9092      ! 
    9193      Agrif_UseSpecialValue = .FALSE. 
     94      use_sign_north = .FALSE. 
    9295#endif 
    9396      ! 
     
    127130 
    128131         ! --- West --- ! 
    129          ztabramp(:,:) = 0._wp 
    130          ind1 = 1+nbghostcells 
    131          DO ji = mi0(ind1), mi1(ind1)                 
    132             ztabramp(ji,:) = ssumask(ji,:) 
    133          END DO 
    134          ! 
    135          zmskwest(:) = 0._wp 
    136          zmskwest(1:jpj) = MAXVAL(ztabramp(:,:), dim=1) 
     132         IF( lk_west) THEN 
     133            ztabramp(:,:) = 0._wp 
     134            ind1 = 1+nbghostcells 
     135            DO ji = mi0(ind1), mi1(ind1)                 
     136               ztabramp(ji,:) = ssumask(ji,:) 
     137            END DO 
     138            ! 
     139            zmskwest(:) = 0._wp 
     140            zmskwest(1:jpj) = MAXVAL(ztabramp(:,:), dim=1) 
     141         ENDIF 
    137142 
    138143         ! --- East --- ! 
    139          ztabramp(:,:) = 0._wp 
    140          ind1 = jpiglo - nbghostcells - 1 
    141          DO ji = mi0(ind1), mi1(ind1)                  
    142             ztabramp(ji,:) = ssumask(ji,:) 
    143          END DO 
    144          ! 
    145          zmskeast(:) = 0._wp 
    146          zmskeast(1:jpj) = MAXVAL(ztabramp(:,:), dim=1) 
     144         IF( lk_east ) THEN 
     145            ztabramp(:,:) = 0._wp 
     146            ind1 = jpiglo - nbghostcells - 1 
     147            DO ji = mi0(ind1), mi1(ind1)                  
     148               ztabramp(ji,:) = ssumask(ji,:) 
     149            END DO 
     150            ! 
     151            zmskeast(:) = 0._wp 
     152            zmskeast(1:jpj) = MAXVAL(ztabramp(:,:), dim=1) 
     153         ENDIF 
    147154 
    148155         ! --- South --- ! 
    149          ztabramp(:,:) = 0._wp 
    150          ind1 = 1+nbghostcells 
    151          DO jj = mj0(ind1), mj1(ind1)                  
    152             ztabramp(:,jj) = ssvmask(:,jj) 
    153          END DO 
    154          ! 
    155          zmsksouth(:) = 0._wp 
    156          zmsksouth(1:jpi) = MAXVAL(ztabramp(:,:), dim=2) 
     156         IF( lk_south ) THEN 
     157            ztabramp(:,:) = 0._wp 
     158            ind1 = 1+nbghostcells 
     159            DO jj = mj0(ind1), mj1(ind1)                  
     160               ztabramp(:,jj) = ssvmask(:,jj) 
     161            END DO 
     162            ! 
     163            zmsksouth(:) = 0._wp 
     164            zmsksouth(1:jpi) = MAXVAL(ztabramp(:,:), dim=2) 
     165         ENDIF 
    157166 
    158167         ! --- North --- ! 
    159          ztabramp(:,:) = 0._wp 
    160          ind1 = jpjglo - nbghostcells - 1 
    161          DO jj = mj0(ind1), mj1(ind1)                  
    162             ztabramp(:,jj) = ssvmask(:,jj) 
    163          END DO 
    164          ! 
    165          zmsknorth(:) = 0._wp 
    166          zmsknorth(1:jpi) = MAXVAL(ztabramp(:,:), dim=2) 
     168         IF( lk_north) THEN 
     169            ztabramp(:,:) = 0._wp 
     170            ind1 = jpjglo - nbghostcells - 1 
     171            DO jj = mj0(ind1), mj1(ind1)                  
     172               ztabramp(:,jj) = ssvmask(:,jj) 
     173            END DO 
     174            ! 
     175            zmsknorth(:) = 0._wp 
     176            zmsknorth(1:jpi) = MAXVAL(ztabramp(:,:), dim=2) 
     177         ENDIF 
     178 
    167179         ! JC: SPONGE MASKING TO BE SORTED OUT: 
    168180         zmskwest(:)  = 1._wp 
     
    192204 
    193205         ! --- West --- ! 
    194          ind1 = 1+nbghostcells 
    195          ind2 = 1+nbghostcells + ispongearea  
    196          DO ji = mi0(ind1), mi1(ind2)    
    197             DO jj = 1, jpj                
    198                ztabramp(ji,jj) = REAL( ind2 - mig(ji) ) * z1_ispongearea * zmskwest(jj) 
    199             END DO 
    200          END DO 
    201  
    202          ! ghost cells: 
    203          ind1 = 1 
    204          ind2 = nbghostcells + 1 
    205          DO ji = mi0(ind1), mi1(ind2)    
    206             DO jj = 1, jpj                
    207                ztabramp(ji,jj) = zmskwest(jj) 
    208             END DO 
    209          END DO 
     206         IF(lk_west) THEN 
     207            ind1 = 1+nbghostcells 
     208            ind2 = 1+nbghostcells + ispongearea  
     209            DO ji = mi0(ind1), mi1(ind2)    
     210               DO jj = 1, jpj                
     211                  ztabramp(ji,jj) = REAL( ind2 - mig(ji) ) * z1_ispongearea * zmskwest(jj) 
     212               END DO 
     213            END DO          
     214 
     215            ! ghost cells: 
     216            ind1 = 1 
     217            ind2 = nbghostcells + 1 
     218            DO ji = mi0(ind1), mi1(ind2)    
     219               DO jj = 1, jpj                
     220                  ztabramp(ji,jj) = zmskwest(jj) 
     221               END DO 
     222            END DO 
     223         ENDIF 
    210224 
    211225         ! --- East --- ! 
    212          ind1 = jpiglo - nbghostcells - ispongearea 
    213          ind2 = jpiglo - nbghostcells 
    214          DO ji = mi0(ind1), mi1(ind2) 
    215             DO jj = 1, jpj 
    216                ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL( mig(ji) - ind1 ) * z1_ispongearea) * zmskeast(jj) 
    217             ENDDO 
    218          END DO 
    219  
    220          ! ghost cells: 
    221          ind1 = jpiglo - nbghostcells 
    222          ind2 = jpiglo 
    223          DO ji = mi0(ind1), mi1(ind2) 
    224             DO jj = 1, jpj 
    225                ztabramp(ji,jj) = zmskeast(jj) 
    226             ENDDO 
    227          END DO 
     226         IF(lk_east) THEN 
     227            ind1 = jpiglo - nbghostcells - ispongearea 
     228            ind2 = jpiglo - nbghostcells 
     229            DO ji = mi0(ind1), mi1(ind2) 
     230 
     231               DO jj = 1, jpj 
     232                  ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL( mig(ji) - ind1 ) * z1_ispongearea) * zmskeast(jj) 
     233               ENDDO 
     234            END DO 
     235 
     236            ! ghost cells: 
     237            ind1 = jpiglo - nbghostcells 
     238            ind2 = jpiglo 
     239            DO ji = mi0(ind1), mi1(ind2) 
     240 
     241               DO jj = 1, jpj 
     242                  ztabramp(ji,jj) = zmskeast(jj) 
     243               ENDDO 
     244            END DO 
     245         ENDIF 
    228246 
    229247         ! --- South --- ! 
    230          ind1 = 1+nbghostcells 
    231          ind2 = 1+nbghostcells + jspongearea 
    232          DO jj = mj0(ind1), mj1(ind2)  
    233             DO ji = 1, jpi 
    234                ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL( ind2 - mjg(jj) ) * z1_jspongearea) * zmsksouth(ji) 
    235             END DO 
    236          END DO 
    237  
    238          ! ghost cells: 
    239          ind1 = 1 
    240          ind2 = nbghostcells + 1 
    241          DO jj = mj0(ind1), mj1(ind2)  
    242             DO ji = 1, jpi 
    243                ztabramp(ji,jj) = zmsksouth(ji) 
    244             END DO 
    245          END DO 
     248         IF( lk_south ) THEN  
     249            ind1 = 1+nbghostcells 
     250            ind2 = 1+nbghostcells + jspongearea 
     251            DO jj = mj0(ind1), mj1(ind2)  
     252               DO ji = 1, jpi 
     253                  ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL( ind2 - mjg(jj) ) * z1_jspongearea) * zmsksouth(ji) 
     254               END DO 
     255            END DO 
     256 
     257            ! ghost cells: 
     258            ind1 = 1 
     259            ind2 = nbghostcells + 1 
     260            DO jj = mj0(ind1), mj1(ind2)  
     261               DO ji = 1, jpi 
     262                  ztabramp(ji,jj) = zmsksouth(ji) 
     263               END DO 
     264            END DO 
     265         ENDIF 
    246266 
    247267         ! --- North --- ! 
    248          ind1 = jpjglo - nbghostcells - jspongearea 
    249          ind2 = jpjglo - nbghostcells 
    250          DO jj = mj0(ind1), mj1(ind2) 
    251             DO ji = 1, jpi 
    252                ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL( mjg(jj) - ind1 ) * z1_jspongearea) * zmsknorth(ji) 
    253             END DO 
    254          END DO 
    255  
    256          ! ghost cells: 
    257          ind1 = jpjglo - nbghostcells 
    258          ind2 = jpjglo 
    259          DO jj = mj0(ind1), mj1(ind2) 
    260             DO ji = 1, jpi 
    261                ztabramp(ji,jj) = zmsknorth(ji) 
    262             END DO 
    263          END DO 
    264  
     268         IF( lk_north ) THEN   
     269            ind1 = jpjglo - nbghostcells - jspongearea 
     270            ind2 = jpjglo - nbghostcells 
     271            DO jj = mj0(ind1), mj1(ind2) 
     272               DO ji = 1, jpi 
     273                  ztabramp(ji,jj) = MAX( ztabramp(ji,jj), REAL( mjg(jj) - ind1 ) * z1_jspongearea) * zmsknorth(ji) 
     274               END DO 
     275            END DO 
     276 
     277            ! ghost cells: 
     278            ind1 = jpjglo - nbghostcells 
     279            ind2 = jpjglo 
     280            DO jj = mj0(ind1), mj1(ind2) 
     281               DO ji = 1, jpi 
     282                  ztabramp(ji,jj) = zmsknorth(ji) 
     283               END DO 
     284            END DO 
     285         ENDIF 
     286       
    265287      ENDIF 
    266288 
     
    334356      INTEGER  ::   ji, jj, jk, jn   ! dummy loop indices 
    335357      INTEGER  ::   iku, ikv 
    336       REAL(wp) :: ztsa, zabe1, zabe2, zbtr, zhtot, ztrelax 
     358      REAL(wp) :: ztsa, zabe1, zabe2, zbtr, zhtot 
    337359      REAL(wp), DIMENSION(i1:i2,j1:j2,jpk) :: ztu, ztv 
    338360      REAL(wp), DIMENSION(i1:i2,j1:j2,jpk,n1:n2) ::tsbdiff 
     
    438460         ENDDO 
    439461 
    440          !* set relaxation time scale 
    441          IF( l_1st_euler .AND. lk_agrif_fstep ) THEN   ;   ztrelax =   rn_trelax_tra  / (        rn_Dt ) 
    442          ELSE                                          ;   ztrelax =   rn_trelax_tra  / (2._wp * rn_Dt ) 
    443          ENDIF 
    444  
    445462         DO jn = 1, jpts             
    446463            DO jk = 1, jpkm1 
     
    448465               DO jj = j1,j2 
    449466                  DO ji = i1,i2-1 
    450                      zabe1 = rn_sponge_tra * fspu(ji,jj) * umask(ji,jj,jk) * e2_e1u(ji,jj) * e3u(ji,jj,jk,Kmm_a) 
     467                     zabe1 = rn_sponge_tra * r1_Dt * fspu(ji,jj) * umask(ji,jj,jk) * e1e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a) 
    451468                     ztu(ji,jj,jk) = zabe1 * ( tsbdiff(ji+1,jj  ,jk,jn) - tsbdiff(ji,jj,jk,jn) )  
    452469                  END DO 
     
    455472               DO ji = i1,i2 
    456473                  DO jj = j1,j2-1 
    457                      zabe2 = rn_sponge_tra * fspv(ji,jj) * vmask(ji,jj,jk) * e1_e2v(ji,jj) * e3v(ji,jj,jk,Kmm_a) 
     474                     zabe2 = rn_sponge_tra * r1_Dt * fspv(ji,jj) * vmask(ji,jj,jk) * e1e2v(ji,jj) * e3v(ji,jj,jk,Kmm_a) 
    458475                     ztv(ji,jj,jk) = zabe2 * ( tsbdiff(ji  ,jj+1,jk,jn) - tsbdiff(ji,jj,jk,jn) ) 
    459476                  END DO 
     
    480497                        ! horizontal diffusive trends 
    481498                        ztsa = zbtr * (  ztu(ji,jj,jk) - ztu(ji-1,jj,jk) + ztv(ji,jj,jk) - ztv(ji,jj-1,jk)  ) & 
    482                              &  - ztrelax * fspt(ji,jj) * tsbdiff(ji,jj,jk,jn)  
     499                             &  - rn_trelax_tra * r1_Dt * fspt(ji,jj) * tsbdiff(ji,jj,jk,jn)  
    483500                        ! add it to the general tracer trends 
    484501                        ts(ji,jj,jk,jn,Krhs_a) = ts(ji,jj,jk,jn,Krhs_a) + ztsa 
     
    507524 
    508525      ! sponge parameters  
    509       REAL(wp) :: ze2u, ze1v, zua, zva, zbtr, zhtot, ztrelax 
     526      REAL(wp) :: ze2u, ze1v, zua, zva, zbtr, zhtot 
    510527      REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk) :: ubdiff 
    511528      REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk) :: rotdiff, hdivdiff 
     
    595612         ubdiff(i1:i2,j1:j2,:) = (uu(i1:i2,j1:j2,:,Kbb_a) - tabres(i1:i2,j1:j2,:,1))*umask(i1:i2,j1:j2,:) 
    596613#endif 
    597          !* set relaxation time scale 
    598          IF( l_1st_euler .AND. lk_agrif_fstep ) THEN   ;   ztrelax =   rn_trelax_dyn  / (        rn_Dt ) 
    599          ELSE                                          ;   ztrelax =   rn_trelax_dyn  / (2._wp * rn_Dt ) 
    600          ENDIF 
    601614         ! 
    602615         DO jk = 1, jpkm1                                 ! Horizontal slab 
     
    608621            DO jj = j1,j2 
    609622               DO ji = i1+1,i2   ! vector opt. 
    610                   zbtr = r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kbb_a) * rn_sponge_dyn * fspt(ji,jj) 
     623                  zbtr = rn_sponge_dyn * r1_Dt * fspt(ji,jj) / e3t(ji,jj,jk,Kbb_a) 
    611624                  hdivdiff(ji,jj,jk) = (  e2u(ji  ,jj)*e3u(ji  ,jj,jk,Kbb_a) * ubdiff(ji  ,jj,jk) & 
    612625                                     &   -e2u(ji-1,jj)*e3u(ji-1,jj,jk,Kbb_a) * ubdiff(ji-1,jj,jk) ) * zbtr 
     
    616629            DO jj = j1,j2-1 
    617630               DO ji = i1,i2   ! vector opt. 
    618                   zbtr = r1_e1e2f(ji,jj) * e3f(ji,jj,jk) * rn_sponge_dyn * fspf(ji,jj) 
     631                  zbtr = rn_sponge_dyn * r1_Dt * fspf(ji,jj) * e3f(ji,jj,jk)  
    619632                  rotdiff(ji,jj,jk) = ( -e1u(ji,jj+1) * ubdiff(ji,jj+1,jk)   & 
    620633                                    &   +e1u(ji,jj  ) * ubdiff(ji,jj  ,jk) ) * fmask(ji,jj,jk) * zbtr  
     
    633646                     zua = - ( ze2u - rotdiff (ji,jj-1,jk) ) / ( e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a) )   & 
    634647                         & + ( hdivdiff(ji+1,jj,jk) - ze1v ) * r1_e1u(ji,jj) &  
    635                          & - ztrelax * fspu(ji,jj) * ubdiff(ji,jj,jk) 
     648                         & - rn_trelax_dyn * r1_Dt * fspu(ji,jj) * ubdiff(ji,jj,jk) 
    636649 
    637650                     ! add it to the general momentum trends 
     
    646659 
    647660         jmax = j2-1 
     661        ! IF (lk_north) jmax = MIN(jmax,nlcj-nbghostcells-2)   ! North 
    648662         IF ((nbondj == 1).OR.(nbondj == 2)) jmax = MIN(jmax,nlcj-nbghostcells-2)   ! North 
    649663 
     
    684698      ! 
    685699      INTEGER  ::   ji, jj, jk, imax 
    686       REAL(wp) ::   ze2u, ze1v, zua, zva, zbtr, zhtot, ztrelax 
     700      REAL(wp) ::   ze2u, ze1v, zua, zva, zbtr, zhtot 
    687701      REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk) :: vbdiff 
    688702      REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk) :: rotdiff, hdivdiff 
     
    771785         vbdiff(i1:i2,j1:j2,:) = (vv(i1:i2,j1:j2,:,Kbb_a) - tabres(i1:i2,j1:j2,:,1))*vmask(i1:i2,j1:j2,:) 
    772786# endif 
    773          !* set relaxation time scale 
    774          IF( l_1st_euler .AND. lk_agrif_fstep ) THEN   ;   ztrelax =   rn_trelax_dyn  / (        rn_Dt ) 
    775          ELSE                                          ;   ztrelax =   rn_trelax_dyn  / (2._wp * rn_Dt ) 
    776          ENDIF 
    777787         ! 
    778788         DO jk = 1, jpkm1                                 ! Horizontal slab 
     
    784794            DO jj = j1+1,j2 
    785795               DO ji = i1,i2   ! vector opt. 
    786                   zbtr = r1_e1e2t(ji,jj) / e3t(ji,jj,jk,Kbb_a) * rn_sponge_dyn * fspt(ji,jj) 
     796                  zbtr = rn_sponge_dyn * r1_Dt * fspt(ji,jj) / e3t(ji,jj,jk,Kbb_a) 
    787797                  hdivdiff(ji,jj,jk) = ( e1v(ji,jj  ) * e3v(ji,jj  ,jk,Kbb_a) * vbdiff(ji,jj  ,jk)  & 
    788798                                     &  -e1v(ji,jj-1) * e3v(ji,jj-1,jk,Kbb_a) * vbdiff(ji,jj-1,jk)  ) * zbtr 
     
    791801            DO jj = j1,j2 
    792802               DO ji = i1,i2-1   ! vector opt. 
    793                   zbtr = r1_e1e2f(ji,jj) * e3f(ji,jj,jk) * rn_sponge_dyn * fspf(ji,jj) 
     803                  zbtr = rn_sponge_dyn * r1_Dt * fspf(ji,jj) * e3f(ji,jj,jk)  
    794804                  rotdiff(ji,jj,jk) = ( e2v(ji+1,jj) * vbdiff(ji+1,jj,jk) &  
    795805                                    &  -e2v(ji  ,jj) * vbdiff(ji  ,jj,jk)  ) * fmask(ji,jj,jk) * zbtr 
     
    802812 
    803813         imax = i2 - 1 
     814      !   IF(lk_east) imax = MIN(imax,nlci-nbghostcells-2)   ! East 
    804815         IF ((nbondi == 1).OR.(nbondi == 2))   imax = MIN(imax,nlci-nbghostcells-2)   ! East 
    805816 
     
    808819               IF( .NOT. tabspongedone_u(ji,jj) ) THEN 
    809820                  DO jk = 1, jpkm1 
    810                      uu(ji,jj,jk,Krhs_a) = uu(ji,jj,jk,Krhs_a)                                                               & 
     821                     uu(ji,jj,jk,Krhs_a) = uu(ji,jj,jk,Krhs_a)                                                     & 
    811822                        & - ( rotdiff (ji  ,jj,jk) - rotdiff (ji,jj-1,jk)) / ( e2u(ji,jj) * e3u(ji,jj,jk,Kmm_a) )  & 
    812823                        & + ( hdivdiff(ji+1,jj,jk) - hdivdiff(ji,jj  ,jk)) * r1_e1u(ji,jj) 
     
    822833               IF( .NOT. tabspongedone_v(ji,jj) ) THEN 
    823834                  DO jk = 1, jpkm1 
    824                      vv(ji,jj,jk,Krhs_a) = vv(ji,jj,jk,Krhs_a)                                                                  & 
     835                     vv(ji,jj,jk,Krhs_a) = vv(ji,jj,jk,Krhs_a)                                                        & 
    825836                        &  + ( rotdiff (ji,jj  ,jk) - rotdiff (ji-1,jj,jk) ) / ( e1v(ji,jj) * e3v(ji,jj,jk,Kmm_a) )   & 
    826                         &  + ( hdivdiff(ji,jj+1,jk) - hdivdiff(ji  ,jj,jk) ) * r1_e2v(ji,jj)                      & 
    827                         &  - ztrelax * fspv(ji,jj) * vbdiff(ji,jj,jk) 
     837                        &  + ( hdivdiff(ji,jj+1,jk) - hdivdiff(ji  ,jj,jk) ) * r1_e2v(ji,jj)                          & 
     838                        &  - rn_trelax_dyn * r1_Dt * fspv(ji,jj) * vbdiff(ji,jj,jk) 
    828839                  END DO 
    829840               ENDIF 
  • NEMO/trunk/src/NST/agrif_oce_update.F90

    r12489 r13216  
    2626   USE domvvl         ! Need interpolation routines  
    2727   USE vremap         ! Vertical remapping 
     28   USE lbclnk  
    2829 
    2930   IMPLICIT NONE 
     
    8586      Agrif_UseSpecialValueInUpdate = .FALSE. 
    8687      Agrif_SpecialValueFineGrid = 0. 
     88 
     89      use_sign_north = .TRUE. 
     90      sign_north = -1. 
     91 
    8792      !      
    8893# if ! defined DECAL_FEEDBACK 
     
    127132      END IF 
    128133      ! 
     134      use_sign_north = .FALSE. 
     135      ! 
    129136   END SUBROUTINE Agrif_Update_Dyn 
    130137 
     
    148155#  if defined VOL_REFLUX 
    149156      IF ( ln_dynspg_ts.AND.ln_bt_fw ) THEN 
     157         use_sign_north = .TRUE. 
     158         sign_north = -1. 
    150159         ! Refluxing on ssh: 
    151160#  if defined DECAL_FEEDBACK_2D 
     
    156165         CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/ 0, 0/),locupdate2=(/-1,-1/),procname = reflux_sshv) 
    157166#  endif 
     167         use_sign_north = .FALSE. 
    158168      END IF 
    159169#  endif 
     
    826836   SUBROUTINE correct_v_bdy( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before, nb, ndir ) 
    827837      !!--------------------------------------------- 
    828       !!           *** ROUTINE correct_u_bdy *** 
     838      !!           *** ROUTINE correct_v_bdy *** 
    829839      !!--------------------------------------------- 
    830840      INTEGER                                     , INTENT(in   ) :: i1, i2, j1, j2, k1, k2, n1, n2 
  • NEMO/trunk/src/NST/agrif_top_interp.F90

    r12377 r13216  
    119119            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)  
    120120         END DO 
    121  
    122121      ENDIF 
    123122      ! 
  • NEMO/trunk/src/NST/agrif_user.F90

    r12489 r13216  
    2828      ! 
    2929      !                    !* Agrif initialization 
    30       CALL agrif_nemo_init 
    31       CALL Agrif_InitValues_cont_dom 
    3230      CALL Agrif_InitValues_cont 
    3331# if defined key_top 
     
    4038   END SUBROUTINE Agrif_initvalues 
    4139 
    42    SUBROUTINE Agrif_InitValues_cont_dom 
    43       !!---------------------------------------------------------------------- 
    44       !!                 *** ROUTINE Agrif_InitValues_cont_dom *** 
    45       !!---------------------------------------------------------------------- 
    46       ! 
    47       CALL agrif_declare_var_dom 
    48       ! 
    49    END SUBROUTINE Agrif_InitValues_cont_dom 
    50  
    51    SUBROUTINE agrif_declare_var_dom 
    52       !!---------------------------------------------------------------------- 
    53       !!                 *** ROUTINE agrif_declare_var_dom *** 
    54       !!---------------------------------------------------------------------- 
    55       USE par_oce, ONLY:  nbghostcells       
     40   SUBROUTINE agrif_istate( Kbb, Kmm, Kaa ) 
     41 
     42       USE domvvl 
     43       USE domain 
     44       USE par_oce 
     45       USE agrif_oce 
     46       USE agrif_oce_interp 
     47       USE oce 
     48       USE lib_mpp 
     49       USe lbclnk 
     50 
     51      INTEGER, INTENT(in)  :: Kbb, Kmm, Kaa 
     52      INTEGER :: jn 
     53 
     54      IF(lwp) WRITE(numout,*) ' ' 
     55      IF(lwp) WRITE(numout,*) 'AGRIF: interp child initial state from parent' 
     56      IF(lwp) WRITE(numout,*) ' ' 
     57 
     58      l_ini_child = .TRUE. 
     59      Agrif_SpecialValue    = 0._wp 
     60      Agrif_UseSpecialValue = .TRUE. 
     61      uu(:,:,:,:) = 0.  ;  vv(:,:,:,:) = 0.   ;  ts(:,:,:,:,:) = 0. 
     62        
     63      Krhs_a = Kbb ; Kmm_a = Kbb 
     64 
     65      ! Brutal fix to pas 1x1 refinment.  
     66  !    IF(Agrif_Irhox() == 1) THEN 
     67  !       CALL Agrif_Init_Variable(tsini_id, procname=agrif_initts) 
     68  !    ELSE 
     69      CALL Agrif_Init_Variable(tsini_id, procname=interptsn) 
     70 
     71  !    ENDIF 
     72! just for VORTEX because Parent velocities can actually be exactly zero 
     73!      Agrif_UseSpecialValue = .FALSE. 
     74      Agrif_UseSpecialValue = ln_spc_dyn 
     75      use_sign_north = .TRUE. 
     76      sign_north = -1. 
     77      CALL Agrif_Init_Variable(uini_id , procname=interpun ) 
     78      CALL Agrif_Init_Variable(vini_id , procname=interpvn ) 
     79      use_sign_north = .FALSE. 
     80 
     81      Agrif_UseSpecialValue = .FALSE.            ! 
     82      l_ini_child = .FALSE. 
     83 
     84      Krhs_a = Kaa ; Kmm_a = Kmm 
     85 
     86      DO jn = 1, jpts 
     87         ts(:,:,:,jn,Kbb) = ts(:,:,:,jn,Kbb)*tmask(:,:,:) 
     88      END DO 
     89      uu(:,:,:,Kbb) =  uu(:,:,:,Kbb) * umask(:,:,:)      
     90      vv(:,:,:,Kbb) =  vv(:,:,:,Kbb) * vmask(:,:,:)  
     91 
     92 
     93      CALL lbc_lnk_multi( 'agrif_istate', uu(:,:,:,Kbb), 'U', -1. , vv(:,:,:,Kbb), 'V', -1. ) 
     94      CALL lbc_lnk( 'agrif_istate', ts(:,:,:,:,Kbb), 'T', 1. ) 
     95 
     96   END SUBROUTINE agrif_istate    
     97 
     98   SUBROUTINE agrif_declare_var_ini 
     99      !!---------------------------------------------------------------------- 
     100      !!                 *** ROUTINE agrif_declare_var *** 
     101      !!---------------------------------------------------------------------- 
     102      USE agrif_util 
     103      USE agrif_oce 
     104      USE par_oce 
     105      USE zdf_oce  
     106      USE oce 
     107      USE dom_oce 
    56108      ! 
    57109      IMPLICIT NONE 
    58110      ! 
    59111      INTEGER :: ind1, ind2, ind3 
    60       !!---------------------------------------------------------------------- 
     112      External :: nemo_mapping 
     113      !!---------------------------------------------------------------------- 
     114 
     115! In case of East-West periodicity, prevent AGRIF interpolation at east and west boundaries 
     116! The procnames will not be called at these boundaries 
     117      IF (jperio == 1) THEN 
     118         CALL Agrif_Set_NearCommonBorderX(.TRUE.) 
     119         CALL Agrif_Set_DistantCommonBorderX(.TRUE.) 
     120      ENDIF 
     121 
     122      IF ( .NOT. lk_south ) THEN 
     123         CALL Agrif_Set_NearCommonBorderY(.TRUE.) 
     124      ENDIF 
    61125 
    62126      ! 1. Declaration of the type of variable which have to be interpolated 
    63127      !--------------------------------------------------------------------- 
    64128      ind1 =     nbghostcells 
    65       ind2 = 1 + nbghostcells 
    66       ind3 = 2 + nbghostcells 
    67       CALL agrif_declare_variable((/1,2/),(/ind2,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),e1u_id) 
    68       CALL agrif_declare_variable((/2,1/),(/ind3,ind2/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),e2v_id) 
    69  
     129      ind2 = 2 + nbghostcells_x 
     130      ind3 = 2 + nbghostcells_y_s 
     131 
     132      CALL agrif_declare_variable((/2,2,0/),(/ind2,ind3,0/),(/'x','y','N'/),(/1,1,1/),(/nlci,nlcj,jpk/),e3t_id) 
     133      CALL agrif_declare_variable((/2,2/)  ,(/ind2,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),mbkt_id) 
     134      CALL agrif_declare_variable((/2,2/)  ,(/ind2,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),ht0_id) 
     135 
     136      CALL agrif_declare_variable((/1,2/),(/ind2-1,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),e1u_id) 
     137      CALL agrif_declare_variable((/2,1/),(/ind2,ind3-1/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),e2v_id) 
     138 
     139    
     140      ! Initial or restart velues 
     141       
     142      CALL agrif_declare_variable((/2,2,0,0/),(/ind2  ,ind3  ,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,jpts+1/),tsini_id) 
     143      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 )  
     144      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 ) 
     145      CALL agrif_declare_variable((/2,2/)    ,(/ind2,ind3/)        ,(/'x','y'/),(/1,1/),(/nlci,nlcj/),sshini_id) 
     146      !  
     147      
    70148      ! 2. Type of interpolation 
    71149      !------------------------- 
     150      CALL Agrif_Set_bcinterp(e3t_id,interp=AGRIF_constant) 
     151 
     152      CALL Agrif_Set_bcinterp(mbkt_id,interp=AGRIF_constant) 
     153      CALL Agrif_Set_interp  (mbkt_id,interp=AGRIF_constant) 
     154      CALL Agrif_Set_bcinterp(ht0_id ,interp=AGRIF_constant) 
     155      CALL Agrif_Set_interp  (ht0_id ,interp=AGRIF_constant) 
     156 
    72157      CALL Agrif_Set_bcinterp( e1u_id, interp1=Agrif_linear, interp2=AGRIF_ppm    ) 
    73158      CALL Agrif_Set_bcinterp( e2v_id, interp1=AGRIF_ppm   , interp2=Agrif_linear ) 
    74159 
    75       ! 3. Location of interpolation 
     160      ! Initial fields 
     161      CALL Agrif_Set_bcinterp(tsini_id ,interp=AGRIF_linear) 
     162      CALL Agrif_Set_interp  (tsini_id ,interp=AGRIF_linear) 
     163      CALL Agrif_Set_bcinterp(uini_id  ,interp=AGRIF_linear) 
     164      CALL Agrif_Set_interp  (uini_id  ,interp=AGRIF_linear) 
     165      CALL Agrif_Set_bcinterp(vini_id  ,interp=AGRIF_linear) 
     166      CALL Agrif_Set_interp  (vini_id  ,interp=AGRIF_linear) 
     167      CALL Agrif_Set_bcinterp(sshini_id,interp=AGRIF_linear) 
     168      CALL Agrif_Set_interp  (sshini_id,interp=AGRIF_linear) 
     169 
     170       ! 3. Location of interpolation 
    76171      !----------------------------- 
     172!      CALL Agrif_Set_bc(  e3t_id, (/-nn_sponge_len*Agrif_irhox(),ind1-1/) )   
     173! JC: check near the boundary only until matching in sponge has been sorted out: 
     174      CALL Agrif_Set_bc(  e3t_id, (/0,ind1-1/) )   
     175 
     176      ! extend the interpolation zone by 1 more point than necessary: 
     177      ! RB check here 
     178      CALL Agrif_Set_bc(  mbkt_id, (/-nn_sponge_len*Agrif_irhox()-2,ind1/) ) 
     179      CALL Agrif_Set_bc(  ht0_id,  (/-nn_sponge_len*Agrif_irhox()-2,ind1/) ) 
     180       
    77181      CALL Agrif_Set_bc(e1u_id,(/0,ind1-1/)) 
    78       CALL Agrif_Set_bc(e2v_id,(/0,ind1-1/)) 
     182      CALL Agrif_Set_bc(e2v_id,(/0,ind1-1/))   
     183 
     184      CALL Agrif_Set_bc( tsini_id , (/0,ind1-1/) ) ! if west,  rhox=3 and nbghost=3: columns 2 to 4 
     185      CALL Agrif_Set_bc( uini_id  , (/0,ind1-1/) )  
     186      CALL Agrif_Set_bc( vini_id  , (/0,ind1-1/) ) 
     187      CALL Agrif_Set_bc( sshini_id, (/0,ind1-1/) ) 
    79188 
    80189      ! 4. Update type 
     
    87196      CALL Agrif_Set_Updatetype(e2v_id,update1 = Agrif_Update_Average, update2=Agrif_Update_Copy) 
    88197#endif 
    89  
    90    END SUBROUTINE agrif_declare_var_dom 
    91  
    92    SUBROUTINE Agrif_InitValues_cont 
    93       !!---------------------------------------------------------------------- 
    94       !!                 *** ROUTINE Agrif_InitValues_cont *** 
    95       !!---------------------------------------------------------------------- 
    96       USE agrif_oce 
     198       
     199   !   CALL Agrif_Set_ExternalMapping(nemo_mapping) 
     200      ! 
     201   END SUBROUTINE agrif_declare_var_ini 
     202 
     203 
     204   SUBROUTINE Agrif_Init_Domain( Kbb, Kmm, Kaa )  
     205      !!---------------------------------------------------------------------- 
     206      !!                 *** ROUTINE Agrif_InitValues_cont_dom *** 
     207      !!---------------------------------------------------------------------- 
     208   
     209         !!---------------------------------------------------------------------- 
     210         !!                 *** ROUTINE Agrif_InitValues_cont *** 
     211         !! 
     212         !! ** Purpose ::   Declaration of variables to be interpolated 
     213         !!---------------------------------------------------------------------- 
     214      USE agrif_oce_update 
    97215      USE agrif_oce_interp 
    98216      USE agrif_oce_sponge 
     217      USE Agrif_Util 
     218      USE oce  
    99219      USE dom_oce 
    100       USE oce 
     220      USE zdf_oce 
     221      USE nemogcm 
     222      USE agrif_oce 
     223      ! 
     224      USE lbclnk 
    101225      USE lib_mpp 
    102       USE lbclnk 
     226      USE in_out_manager 
    103227      ! 
    104228      IMPLICIT NONE 
    105229      ! 
    106       INTEGER :: ji, jj 
     230      INTEGER, INTENT(in) ::  Kbb, Kmm, Kaa 
     231      ! 
    107232      LOGICAL :: check_namelist 
    108233      CHARACTER(len=15) :: cl_check1, cl_check2, cl_check3, cl_check4  
    109 #if defined key_vertical 
    110234      REAL(wp), DIMENSION(jpi,jpj) ::   zk   ! workspace 
    111 #endif 
    112       !!---------------------------------------------------------------------- 
    113  
    114       ! 1. Declaration of the type of variable which have to be interpolated 
    115       !--------------------------------------------------------------------- 
    116       CALL agrif_declare_var 
    117  
    118       ! 2. First interpolations of potentially non zero fields 
    119       !------------------------------------------------------- 
    120  
    121 #if defined key_vertical 
     235      INTEGER :: ji, jj, jk 
     236      !!---------------------------------------------------------------------- 
     237     
     238     ! CALL Agrif_Declare_Var_ini 
     239 
     240      IF( agrif_oce_alloc()  > 0 )   CALL ctl_warn('agrif agrif_oce_alloc: allocation of arrays failed') 
     241 
    122242      ! Build consistent parent bathymetry and number of levels 
    123243      ! on the child grid  
     
    126246      mbkt_parent(:,:) = 0 
    127247      ! 
    128       CALL Agrif_Bc_variable(ht0_id ,calledweight=1.,procname=interpht0 ) 
    129       CALL Agrif_Bc_variable(mbkt_id,calledweight=1.,procname=interpmbkt) 
     248  !    CALL Agrif_Bc_variable(ht0_id ,calledweight=1.,procname=interpht0 ) 
     249  !    CALL Agrif_Bc_variable(mbkt_id,calledweight=1.,procname=interpmbkt) 
     250      CALL Agrif_Init_Variable(ht0_id , procname=interpht0 ) 
     251      CALL Agrif_Init_Variable(mbkt_id, procname=interpmbkt) 
    130252      ! 
    131253      ! Assume step wise change of bathymetry near interface 
     
    149271      ENDIF 
    150272      ! 
    151       CALL lbc_lnk( 'Agrif_InitValues_cont', hu0_parent, 'U', 1. ) 
    152       CALL lbc_lnk( 'Agrif_InitValues_cont', hv0_parent, 'V', 1. ) 
     273      CALL lbc_lnk( 'Agrif_Init_Domain', hu0_parent, 'U', 1. ) 
     274      CALL lbc_lnk( 'Agrif_Init_Domain', hv0_parent, 'V', 1. ) 
    153275      zk(:,:) = REAL( mbku_parent(:,:), wp )   ;   CALL lbc_lnk( 'Agrif_InitValues_cont', zk, 'U', 1. ) 
    154       mbku_parent(:,:) = MAX( NINT( zk(:,:) ), 1 ) 
     276      mbku_parent(:,:) = MAX( NINT( zk(:,:) ), 1 ) ; 
    155277      zk(:,:) = REAL( mbkv_parent(:,:), wp )   ;   CALL lbc_lnk( 'Agrif_InitValues_cont', zk, 'V', 1. ) 
    156278      mbkv_parent(:,:) = MAX( NINT( zk(:,:) ), 1 )    
    157 #endif 
    158  
     279 
     280      IF ( ln_init_chfrpar ) THEN  
     281         CALL Agrif_Init_Variable(sshini_id, procname=agrif_initssh) 
     282         CALL lbc_lnk( 'Agrif_Init_Domain', ssh(:,:,Kbb), 'T', 1. ) 
     283         DO jk = 1, jpk 
     284               e3t(:,:,jk,Kbb) =  e3t_0(:,:,jk) * ( ht_0(:,:) + ssh(:,:,Kbb)  ) & 
     285                        &             / ( ht_0(:,:) + 1._wp - ssmask(:,:) ) * tmask(:,:,jk)   & 
     286                        &              + e3t_0(:,:,jk) * ( 1._wp - tmask(:,:,jk) ) 
     287         END DO 
     288      ENDIF 
     289 
     290      ! check if masks and bathymetries match 
     291      IF(ln_chk_bathy) THEN 
     292         Agrif_UseSpecialValue = .FALSE. 
     293         ! 
     294         IF(lwp) WRITE(numout,*) ' ' 
     295         IF(lwp) WRITE(numout,*) 'AGRIF: Check Bathymetry and masks near bdys. Level: ', Agrif_Level() 
     296         ! 
     297         kindic_agr = 0 
     298         IF( .NOT. l_vremap ) THEN 
     299            ! 
     300            ! check if tmask and vertical scale factors agree with parent in sponge area: 
     301            CALL Agrif_Bc_variable(e3t_id,calledweight=1.,procname=interpe3t) 
     302            ! 
     303         ELSE 
     304            ! 
     305            ! In case of vertical interpolation, check only that total depths agree between child and parent: 
     306            DO ji = 1, jpi 
     307               DO jj = 1, jpj 
     308                  IF ((mbkt_parent(ji,jj)/=0).AND.(ABS(ht0_parent(ji,jj)-ht_0(ji,jj))>1.e-3)) kindic_agr = kindic_agr + 1 
     309                  IF ((mbku_parent(ji,jj)/=0).AND.(ABS(hu0_parent(ji,jj)-hu_0(ji,jj))>1.e-3)) kindic_agr = kindic_agr + 1 
     310                  IF ((mbkv_parent(ji,jj)/=0).AND.(ABS(hv0_parent(ji,jj)-hv_0(ji,jj))>1.e-3)) kindic_agr = kindic_agr + 1 
     311               END DO 
     312            END DO 
     313 
     314            CALL mpp_sum( 'agrif_user', kindic_agr ) 
     315            IF( kindic_agr /= 0 ) THEN 
     316               CALL ctl_stop('==> Child Bathymetry is NOT correct near boundaries.') 
     317            ELSE 
     318               IF(lwp) WRITE(numout,*) '==> Child Bathymetry is ok near boundaries.' 
     319               IF(lwp) WRITE(numout,*) ' ' 
     320            ENDIF   
     321         ENDIF 
     322      ENDIF 
     323 
     324      IF( l_vremap ) THEN 
     325      ! Additional constrain that should be removed someday: 
     326         IF ( Agrif_Parent(jpk).GT.jpk ) THEN 
     327            CALL ctl_stop( ' With l_vremap, child grids must have jpk greater or equal to the parent value' ) 
     328         ENDIF 
     329      ENDIF 
     330      ! 
     331   END SUBROUTINE Agrif_Init_Domain 
     332 
     333 
     334   SUBROUTINE Agrif_InitValues_cont 
     335         !!---------------------------------------------------------------------- 
     336         !!                 *** ROUTINE Agrif_InitValues_cont *** 
     337         !! 
     338         !! ** Purpose ::   Declaration of variables to be interpolated 
     339         !!---------------------------------------------------------------------- 
     340      USE agrif_oce_update 
     341      USE agrif_oce_interp 
     342      USE agrif_oce_sponge 
     343      USE Agrif_Util 
     344      USE oce  
     345      USE dom_oce 
     346      USE zdf_oce 
     347      USE nemogcm 
     348      USE agrif_oce 
     349      ! 
     350      USE lbclnk 
     351      USE lib_mpp 
     352      USE in_out_manager 
     353      ! 
     354      IMPLICIT NONE 
     355      ! 
     356      LOGICAL :: check_namelist 
     357      CHARACTER(len=15) :: cl_check1, cl_check2, cl_check3, cl_check4  
     358      REAL(wp), DIMENSION(jpi,jpj) ::   zk   ! workspace 
     359      INTEGER :: ji, jj 
     360 
     361      ! 1. Declaration of the type of variable which have to be interpolated 
     362      !--------------------------------------------------------------------- 
     363      CALL agrif_declare_var 
     364 
     365      ! 2. First interpolations of potentially non zero fields 
     366      !------------------------------------------------------- 
    159367      Agrif_SpecialValue    = 0._wp 
    160368      Agrif_UseSpecialValue = .TRUE. 
     
    163371      tabspongedone_tsn = .FALSE. 
    164372      CALL Agrif_Bc_variable(tsn_sponge_id,calledweight=1.,procname=interptsn_sponge) 
    165       ! reset ts(:,:,:,:,Krhs_a) to zero 
     373      ! reset tsa to zero 
    166374      ts(:,:,:,:,Krhs_a) = 0._wp 
    167375 
    168376      Agrif_UseSpecialValue = ln_spc_dyn 
     377      use_sign_north = .TRUE. 
     378      sign_north = -1. 
    169379      CALL Agrif_Bc_variable(un_interp_id,calledweight=1.,procname=interpun) 
    170380      CALL Agrif_Bc_variable(vn_interp_id,calledweight=1.,procname=interpvn) 
     
    175385      tabspongedone_v = .FALSE. 
    176386      CALL Agrif_Bc_variable(vn_sponge_id,calledweight=1.,procname=interpvn_sponge) 
     387      use_sign_north = .FALSE. 
    177388      uu(:,:,:,Krhs_a) = 0._wp 
    178389      vv(:,:,:,Krhs_a) = 0._wp 
     
    185396      IF ( ln_dynspg_ts ) THEN 
    186397         Agrif_UseSpecialValue = ln_spc_dyn 
     398         use_sign_north = .TRUE. 
     399         sign_north = -1. 
    187400         CALL Agrif_Bc_variable(unb_id,calledweight=1.,procname=interpunb) 
    188401         CALL Agrif_Bc_variable(vnb_id,calledweight=1.,procname=interpvnb) 
    189402         CALL Agrif_Bc_variable(ub2b_interp_id,calledweight=1.,procname=interpub2b) 
    190403         CALL Agrif_Bc_variable(vb2b_interp_id,calledweight=1.,procname=interpvb2b) 
     404         use_sign_north = .FALSE. 
    191405         ubdy(:,:) = 0._wp 
    192406         vbdy(:,:) = 0._wp 
    193407      ENDIF 
    194  
    195       Agrif_UseSpecialValue = .FALSE. 
    196  
    197       ! 3. Some controls 
     408      Agrif_UseSpecialValue = .FALSE.  
     409 
    198410      !----------------- 
    199411      check_namelist = .TRUE. 
    200412 
    201413      IF( check_namelist ) THEN  
    202  
    203          ! Check time steps            
    204          IF( NINT(Agrif_Rhot()) * NINT(rn_Dt) .NE. Agrif_Parent(rn_Dt) ) THEN 
    205             WRITE(cl_check1,*)  NINT(Agrif_Parent(rn_Dt)) 
    206             WRITE(cl_check2,*)  NINT(rn_Dt) 
    207             WRITE(cl_check3,*)  NINT(Agrif_Parent(rn_Dt)/Agrif_Rhot()) 
    208             CALL ctl_stop( 'Incompatible time step between ocean grids',   & 
    209                   &               'parent grid value : '//cl_check1    ,   &  
    210                   &               'child  grid value : '//cl_check2    ,   &  
    211                   &               'value on child grid should be changed to : '//cl_check3 ) 
    212          ENDIF 
    213  
    214          ! Check run length 
    215          IF( Agrif_IRhot() * (Agrif_Parent(nitend)- & 
    216                Agrif_Parent(nit000)+1) .NE. (nitend-nit000+1) ) THEN 
    217             WRITE(cl_check1,*)  (Agrif_Parent(nit000)-1)*Agrif_IRhot() + 1 
    218             WRITE(cl_check2,*)   Agrif_Parent(nitend)   *Agrif_IRhot() 
    219             CALL ctl_warn( 'Incompatible run length between grids'                      ,   & 
    220                   &               'nit000 on fine grid will be changed to : '//cl_check1,   & 
    221                   &               'nitend on fine grid will be changed to : '//cl_check2    ) 
    222             nit000 = (Agrif_Parent(nit000)-1)*Agrif_IRhot() + 1 
    223             nitend =  Agrif_Parent(nitend)   *Agrif_IRhot() 
    224          ENDIF 
    225  
    226414         ! Check free surface scheme 
    227415         IF ( ( Agrif_Parent(ln_dynspg_ts ).AND.ln_dynspg_exp ).OR.& 
     
    251439            STOP 
    252440         ENDIF 
    253  
    254       ENDIF 
    255  
    256       ! check if masks and bathymetries match 
    257       IF(ln_chk_bathy) THEN 
    258          Agrif_UseSpecialValue = .FALSE. 
    259          ! 
    260          IF(lwp) WRITE(numout,*) ' ' 
    261          IF(lwp) WRITE(numout,*) 'AGRIF: Check Bathymetry and masks near bdys. Level: ', Agrif_Level() 
    262          ! 
    263          kindic_agr = 0 
    264 # if ! defined key_vertical 
    265          ! 
    266          ! check if tmask and vertical scale factors agree with parent in sponge area: 
    267          CALL Agrif_Bc_variable(e3t_id,calledweight=1.,procname=interpe3t) 
    268          ! 
    269 # else 
    270          ! 
    271          ! In case of vertical interpolation, check only that total depths agree between child and parent: 
    272          DO ji = 1, jpi 
    273             DO jj = 1, jpj 
    274                IF ((mbkt_parent(ji,jj)/=0).AND.(ABS(ht0_parent(ji,jj)-ht_0(ji,jj))>1.e-3)) kindic_agr = kindic_agr + 1 
    275                IF ((mbku_parent(ji,jj)/=0).AND.(ABS(hu0_parent(ji,jj)-hu_0(ji,jj))>1.e-3)) kindic_agr = kindic_agr + 1 
    276                IF ((mbkv_parent(ji,jj)/=0).AND.(ABS(hv0_parent(ji,jj)-hv_0(ji,jj))>1.e-3)) kindic_agr = kindic_agr + 1 
    277             END DO 
    278          END DO 
    279 # endif 
    280          CALL mpp_sum( 'agrif_user', kindic_agr ) 
    281          IF( kindic_agr /= 0 ) THEN 
    282             CALL ctl_stop('==> Child Bathymetry is NOT correct near boundaries.') 
    283          ELSE 
    284             IF(lwp) WRITE(numout,*) '==> Child Bathymetry is ok near boundaries.' 
    285             IF(lwp) WRITE(numout,*) ' ' 
    286          END IF   
    287          !     
    288       ENDIF 
    289  
    290 # if defined key_vertical 
    291       ! Additional constrain that should be removed someday: 
    292       IF ( Agrif_Parent(jpk).GT.jpk ) THEN 
    293     CALL ctl_stop( ' With key_vertical, child grids must have jpk greater or equal to the parent value' ) 
    294       ENDIF 
    295 # endif 
    296       !  
     441      ENDIF 
     442 
    297443   END SUBROUTINE Agrif_InitValues_cont 
    298444 
     
    314460      ! 1. Declaration of the type of variable which have to be interpolated 
    315461      !--------------------------------------------------------------------- 
     462 
    316463      ind1 =     nbghostcells 
    317       ind2 = 1 + nbghostcells 
    318       ind3 = 2 + nbghostcells 
     464      ind2 = 2 + nbghostcells_x 
     465      ind3 = 2 + nbghostcells_y_s 
     466 
    319467# if defined key_vertical 
    320       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) 
    321       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) 
    322  
    323       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) 
    324       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) 
    325       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) 
    326       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) 
    327       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) 
    328       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) 
     468      CALL agrif_declare_variable((/2,2,0,0/),(/ind2,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,jpts+1/),tsn_id) 
     469      CALL agrif_declare_variable((/2,2,0,0/),(/ind2,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,jpts+1/),tsn_sponge_id) 
     470      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) 
     471      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) 
     472      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) 
     473      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) 
     474      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) 
     475      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) 
    329476# else 
    330       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) 
    331       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) 
    332  
    333       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) 
    334       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) 
    335       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) 
    336       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) 
    337       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) 
    338       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) 
     477      CALL agrif_declare_variable((/2,2,0,0/),(/ind2,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,jpts/),tsn_id) 
     478      CALL agrif_declare_variable((/2,2,0,0/),(/ind2,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,jpts/),tsn_sponge_id) 
     479      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) 
     480      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) 
     481      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) 
     482      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) 
     483      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) 
     484      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) 
    339485# endif 
    340  
    341       CALL agrif_declare_variable((/2,2,0/),(/ind3,ind3,0/),(/'x','y','N'/),(/1,1,1/),(/nlci,nlcj,jpk/),e3t_id) 
    342  
    343 # if defined key_vertical 
    344       CALL agrif_declare_variable((/2,2/),(/ind3,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),mbkt_id) 
    345       CALL agrif_declare_variable((/2,2/),(/ind3,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),ht0_id) 
    346 # endif 
    347  
    348       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) 
    349  
    350       CALL agrif_declare_variable((/1,2/),(/ind2,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),unb_id) 
    351       CALL agrif_declare_variable((/2,1/),(/ind3,ind2/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),vnb_id) 
    352       CALL agrif_declare_variable((/1,2/),(/ind2,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),ub2b_interp_id) 
    353       CALL agrif_declare_variable((/2,1/),(/ind3,ind2/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),vb2b_interp_id) 
    354       CALL agrif_declare_variable((/1,2/),(/ind2,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),ub2b_update_id) 
    355       CALL agrif_declare_variable((/2,1/),(/ind3,ind2/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),vb2b_update_id) 
    356  
    357       CALL agrif_declare_variable((/2,2/),(/ind3,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),sshn_id) 
    358  
    359       IF( ln_zdftke.OR.ln_zdfgls ) THEN 
     486      CALL agrif_declare_variable((/1,2/),(/ind2-1,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),unb_id) 
     487      CALL agrif_declare_variable((/2,1/),(/ind2,ind3-1/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),vnb_id) 
     488      CALL agrif_declare_variable((/1,2/),(/ind2-1,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),ub2b_interp_id) 
     489      CALL agrif_declare_variable((/2,1/),(/ind2,ind3-1/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),vb2b_interp_id) 
     490      CALL agrif_declare_variable((/1,2/),(/ind2-1,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),ub2b_update_id) 
     491      CALL agrif_declare_variable((/2,1/),(/ind2,ind3-1/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),vb2b_update_id) 
     492 
     493      CALL agrif_declare_variable((/2,2/),(/ind2,ind3/),(/'x','y'/),(/1,1/),(/nlci,nlcj/),sshn_id) 
     494 
     495 
     496      IF( ln_zdftke.OR.ln_zdfgls ) THEN  ! logical not known at this point 
    360497!         CALL agrif_declare_variable((/2,2,0/),(/ind3,ind3,0/),(/'x','y','N'/),(/1,1,1/),(/nlci,nlcj,jpk/), en_id) 
    361498!         CALL agrif_declare_variable((/2,2,0/),(/ind3,ind3,0/),(/'x','y','N'/),(/1,1,1/),(/nlci,nlcj,jpk/),avt_id) 
    362499# if defined key_vertical 
    363          CALL agrif_declare_variable((/2,2,0,0/),(/ind3,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,2/),avm_id) 
     500         CALL agrif_declare_variable((/2,2,0,0/),(/ind2,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,2/),avm_id) 
    364501# else 
    365          CALL agrif_declare_variable((/2,2,0,0/),(/ind3,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,1/),avm_id) 
     502         CALL agrif_declare_variable((/2,2,0,0/),(/ind2,ind3,0,0/),(/'x','y','N','N'/),(/1,1,1,1/),(/nlci,nlcj,jpk,1/),avm_id) 
    366503# endif 
    367504      ENDIF 
    368  
     505      
    369506      ! 2. Type of interpolation 
    370507      !------------------------- 
    371508      CALL Agrif_Set_bcinterp(tsn_id,interp=AGRIF_linear) 
    372  
    373509      CALL Agrif_Set_bcinterp(un_interp_id,interp1=Agrif_linear,interp2=AGRIF_ppm) 
    374510      CALL Agrif_Set_bcinterp(vn_interp_id,interp1=AGRIF_ppm,interp2=Agrif_linear) 
    375511 
    376512      CALL Agrif_Set_bcinterp(tsn_sponge_id,interp=AGRIF_linear) 
     513      CALL Agrif_Set_bcinterp(un_sponge_id,interp1=Agrif_linear,interp2=AGRIF_ppm) 
     514      CALL Agrif_Set_bcinterp(vn_sponge_id,interp1=AGRIF_ppm,interp2=Agrif_linear) 
    377515 
    378516      CALL Agrif_Set_bcinterp(sshn_id,interp=AGRIF_linear) 
     
    390528!< 
    391529 
    392       CALL Agrif_Set_bcinterp(un_sponge_id,interp1=Agrif_linear,interp2=AGRIF_ppm) 
    393       CALL Agrif_Set_bcinterp(vn_sponge_id,interp1=AGRIF_ppm,interp2=Agrif_linear) 
    394  
    395       CALL Agrif_Set_bcinterp(e3t_id,interp=AGRIF_constant) 
    396  
    397 # if defined key_vertical 
    398       CALL Agrif_Set_bcinterp(mbkt_id,interp=AGRIF_constant) 
    399       CALL Agrif_Set_bcinterp(ht0_id ,interp=AGRIF_constant) 
    400 # endif 
    401  
    402       IF( ln_zdftke.OR.ln_zdfgls )   CALL Agrif_Set_bcinterp( avm_id, interp=AGRIF_linear ) 
    403  
    404       ! 3. Location of interpolation 
     530      IF( ln_zdftke.OR.ln_zdfgls )  CALL Agrif_Set_bcinterp( avm_id, interp=AGRIF_linear ) 
     531     
     532 
     533       ! 3. Location of interpolation 
    405534      !----------------------------- 
    406535      CALL Agrif_Set_bc(       tsn_id, (/0,ind1-1/) ) ! if west,  rhox=3 and nbghost=3: columns 2 to 4 
     
    418547      CALL Agrif_Set_bc( vb2b_interp_id, (/0,ind1-1/) ) 
    419548 
    420 !      CALL Agrif_Set_bc(  e3t_id, (/-nn_sponge_len*Agrif_irhox(),ind1-1/) )   
    421 ! JC: check near the boundary only until matching in sponge has been sorted out: 
    422       CALL Agrif_Set_bc(  e3t_id, (/0,ind1-1/) )   
    423  
    424 # if defined key_vertical  
    425       ! extend the interpolation zone by 1 more point than necessary: 
    426       CALL Agrif_Set_bc(  mbkt_id, (/-nn_sponge_len*Agrif_irhox()-2,ind1/) ) 
    427       CALL Agrif_Set_bc(  ht0_id,  (/-nn_sponge_len*Agrif_irhox()-2,ind1/) ) 
    428 # endif 
    429  
    430       IF( ln_zdftke.OR.ln_zdfgls )   CALL Agrif_Set_bc( avm_id, (/0,ind1/) ) 
     549      IF( ln_zdftke.OR.ln_zdfgls ) CALL Agrif_Set_bc( avm_id, (/0,ind1/) ) 
    431550 
    432551      ! 4. Update type 
    433552      !---------------  
    434       CALL Agrif_Set_Updatetype(scales_t_id, update = AGRIF_Update_Average) 
    435553 
    436554# if defined UPD_HIGH 
     
    444562      CALL Agrif_Set_Updatetype(e3t_id, update = Agrif_Update_Full_Weighting) 
    445563 
    446       IF( ln_zdftke.OR.ln_zdfgls ) THEN 
     564  !    IF( ln_zdftke.OR.ln_zdfgls ) THEN 
    447565!         CALL Agrif_Set_Updatetype( en_id, update = AGRIF_Update_Full_Weighting) 
    448566!         CALL Agrif_Set_Updatetype(avt_id, update = AGRIF_Update_Full_Weighting) 
    449567!         CALL Agrif_Set_Updatetype(avm_id, update = AGRIF_Update_Full_Weighting) 
    450       ENDIF 
     568   !   ENDIF 
    451569 
    452570#else 
     
    460578      CALL Agrif_Set_Updatetype(e3t_id, update = AGRIF_Update_Average) 
    461579 
    462       IF( ln_zdftke.OR.ln_zdfgls ) THEN 
     580 !     IF( ln_zdftke.OR.ln_zdfgls ) THEN 
    463581!         CALL Agrif_Set_Updatetype( en_id, update = AGRIF_Update_Average) 
    464582!         CALL Agrif_Set_Updatetype(avt_id, update = AGRIF_Update_Average) 
    465583!         CALL Agrif_Set_Updatetype(avm_id, update = AGRIF_Update_Average) 
    466       ENDIF 
     584 !     ENDIF 
    467585 
    468586#endif 
     
    472590#if defined key_si3 
    473591SUBROUTINE Agrif_InitValues_cont_ice 
    474       !!---------------------------------------------------------------------- 
    475       !!                 *** ROUTINE Agrif_InitValues_cont_ice *** 
    476       !!---------------------------------------------------------------------- 
    477592      USE Agrif_Util 
    478593      USE sbc_oce, ONLY : nn_fsbc  ! clem: necessary otherwise Agrif_Parent(nn_fsbc) = nn_fsbc 
     
    482597      USE agrif_ice_interp 
    483598      USE lib_mpp 
    484       ! 
    485       IMPLICIT NONE 
    486       !!---------------------------------------------------------------------- 
    487       ! 
    488       ! Declaration of the type of variable which have to be interpolated (parent=>child) 
    489       !---------------------------------------------------------------------------------- 
    490       CALL agrif_declare_var_ice 
     599      !!---------------------------------------------------------------------- 
     600      !!                 *** ROUTINE Agrif_InitValues_cont_ice *** 
     601      !!---------------------------------------------------------------------- 
    491602 
    492603      ! Controls 
     
    495606      !          the run must satisfy CFL=Uice/(dx/dt) < 0.6/nn_fsbc(child) 
    496607      !          therefore, if nn_fsbc(child)>1 one must reduce the time-step in proportion to nn_fsbc(child), which is not acceptable 
    497       !       If a solution is found, the following stop could be removed because the rest of the code take nn_fsbc(child) into account 
     608      !       If a solution is found, the following stop could be removed because the rest of the code take nn_fsbc(child) into account      
    498609      IF( nn_fsbc > 1 )  CALL ctl_stop('nn_fsbc(child) must be set to 1 otherwise agrif and sea-ice may not work properly') 
    499610 
     
    516627      !!                 *** ROUTINE agrif_declare_var_ice *** 
    517628      !!---------------------------------------------------------------------- 
     629 
    518630      USE Agrif_Util 
    519631      USE ice 
    520       USE par_oce, ONLY : nbghostcells 
     632      USE par_oce, ONLY : nbghostcells, nbghostcells_x, nbghostcells_y_s 
    521633      ! 
    522634      IMPLICIT NONE 
    523635      ! 
    524636      INTEGER :: ind1, ind2, ind3 
    525       !!---------------------------------------------------------------------- 
     637         !!---------------------------------------------------------------------- 
    526638      ! 
    527639      ! 1. Declaration of the type of variable which have to be interpolated (parent=>child) 
     
    532644      !                            2,2 = two ghost lines 
    533645      !------------------------------------------------------------------------------------- 
     646 
    534647      ind1 =     nbghostcells 
    535       ind2 = 1 + nbghostcells 
    536       ind3 = 2 + nbghostcells 
    537       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) 
    538       CALL agrif_declare_variable((/1,2/)  ,(/ind2,ind3/)  ,(/'x','y'/)    ,(/1,1/)  ,(/nlci,nlcj/)                      ,u_ice_id  ) 
    539       CALL agrif_declare_variable((/2,1/)  ,(/ind3,ind2/)  ,(/'x','y'/)    ,(/1,1/)  ,(/nlci,nlcj/)                      ,v_ice_id  ) 
     648      ind2 = 2 + nbghostcells_x 
     649      ind3 = 2 + nbghostcells_y_s 
     650      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) 
     651      CALL agrif_declare_variable((/1,2/)  ,(/ind2-1,ind3/)  ,(/'x','y'/)    ,(/1,1/)  ,(/nlci,nlcj/)                      ,u_ice_id  ) 
     652      CALL agrif_declare_variable((/2,1/)  ,(/ind2,ind3-1/)  ,(/'x','y'/)    ,(/1,1/)  ,(/nlci,nlcj/)                      ,v_ice_id  ) 
     653 
     654      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) 
     655      CALL agrif_declare_variable((/1,2/)  ,(/ind2-1,ind3/)  ,(/'x','y'/)    ,(/1,1/)  ,(/nlci,nlcj/)                      ,u_iceini_id  ) 
     656      CALL agrif_declare_variable((/2,1/)  ,(/ind2,ind3-1/)  ,(/'x','y'/)    ,(/1,1/)  ,(/nlci,nlcj/)                      ,v_iceini_id  ) 
    540657 
    541658      ! 2. Set interpolations (normal & tangent to the grid cell for velocities) 
     
    545662      CALL Agrif_Set_bcinterp(v_ice_id  , interp1 = AGRIF_ppm   ,interp2 = Agrif_linear) 
    546663 
     664      CALL Agrif_Set_bcinterp(tra_iceini_id, interp  = AGRIF_linear) 
     665      CALL Agrif_Set_interp  (tra_iceini_id, interp  = AGRIF_linear) 
     666      CALL Agrif_Set_bcinterp(u_iceini_id  , interp  = AGRIF_linear  ) 
     667      CALL Agrif_Set_interp  (u_iceini_id  , interp  = AGRIF_linear   ) 
     668      CALL Agrif_Set_bcinterp(v_iceini_id  , interp  = AGRIF_linear) 
     669      CALL Agrif_Set_interp  (v_iceini_id  , interp  = AGRIF_linear) 
     670 
    547671      ! 3. Set location of interpolations 
    548672      !---------------------------------- 
     
    550674      CALL Agrif_Set_bc(u_ice_id  ,(/0,ind1/)) 
    551675      CALL Agrif_Set_bc(v_ice_id  ,(/0,ind1/)) 
     676 
     677      CALL Agrif_Set_bc(tra_iceini_id,(/0,ind1/)) 
     678      CALL Agrif_Set_bc(u_iceini_id  ,(/0,ind1/)) 
     679      CALL Agrif_Set_bc(v_iceini_id  ,(/0,ind1/)) 
    552680 
    553681      ! 4. Set update type in case 2 ways (child=>parent) (normal & tangent to the grid cell for velocities) 
     
    557685      CALL Agrif_Set_Updatetype(u_ice_id  , update1 = Agrif_Update_Average       , update2 = Agrif_Update_Full_Weighting) 
    558686      CALL Agrif_Set_Updatetype(v_ice_id  , update1 = Agrif_Update_Full_Weighting, update2 = Agrif_Update_Average       ) 
    559 #else 
     687# else 
    560688      CALL Agrif_Set_Updatetype(tra_ice_id, update  = AGRIF_Update_Average) 
    561689      CALL Agrif_Set_Updatetype(u_ice_id  , update1 = Agrif_Update_Copy   , update2 = Agrif_Update_Average) 
    562690      CALL Agrif_Set_Updatetype(v_ice_id  , update1 = Agrif_Update_Average, update2 = Agrif_Update_Copy   ) 
    563 #endif 
     691# endif 
    564692 
    565693   END SUBROUTINE agrif_declare_var_ice 
     
    585713      USE agrif_top_sponge 
    586714      !! 
    587       IMPLICIT NONE 
    588       ! 
    589       CHARACTER(len=10) :: cl_check1, cl_check2, cl_check3 
    590       LOGICAL :: check_namelist 
    591       !!---------------------------------------------------------------------- 
    592  
    593       ! 1. Declaration of the type of variable which have to be interpolated 
    594       !--------------------------------------------------------------------- 
    595       CALL agrif_declare_var_top 
    596  
    597       ! 2. First interpolations of potentially non zero fields 
    598       !------------------------------------------------------- 
    599       Agrif_SpecialValue=0._wp 
    600       Agrif_UseSpecialValue = .TRUE. 
    601       CALL Agrif_Bc_variable(trn_id,calledweight=1.,procname=interptrn) 
    602       Agrif_UseSpecialValue = .FALSE. 
    603       CALL Agrif_Sponge 
    604       tabspongedone_trn = .FALSE. 
    605       CALL Agrif_Bc_variable(trn_sponge_id,calledweight=1.,procname=interptrn_sponge) 
    606       ! reset ts(:,:,:,:,Krhs_a) to zero 
    607       tr(:,:,:,:,Krhs_a) = 0._wp 
    608  
    609       ! 3. Some controls 
    610       !----------------- 
    611       check_namelist = .TRUE. 
    612  
    613       IF( check_namelist ) THEN 
    614          ! Check time steps 
    615       IF( NINT(Agrif_Rhot()) * NINT(rn_Dt) .NE. Agrif_Parent(rn_Dt) ) THEN 
    616          WRITE(cl_check1,*)  Agrif_Parent(rn_Dt) 
    617          WRITE(cl_check2,*)  rn_Dt 
    618          WRITE(cl_check3,*)  rn_Dt*Agrif_Rhot() 
     715   
     716   !! 
     717   IMPLICIT NONE 
     718   ! 
     719   CHARACTER(len=10) :: cl_check1, cl_check2, cl_check3 
     720   LOGICAL :: check_namelist 
     721      !!---------------------------------------------------------------------- 
     722 
     723 
     724   ! 1. Declaration of the type of variable which have to be interpolated 
     725   !--------------------------------------------------------------------- 
     726   CALL agrif_declare_var_top 
     727 
     728   ! 2. First interpolations of potentially non zero fields 
     729   !------------------------------------------------------- 
     730   Agrif_SpecialValue=0. 
     731   Agrif_UseSpecialValue = .TRUE. 
     732   CALL Agrif_Bc_variable(trn_id,calledweight=1.,procname=interptrn) 
     733   Agrif_UseSpecialValue = .FALSE. 
     734   CALL Agrif_Sponge 
     735   tabspongedone_trn = .FALSE. 
     736   CALL Agrif_Bc_variable(trn_sponge_id,calledweight=1.,procname=interptrn_sponge) 
     737   ! reset tsa to zero 
     738   tra(:,:,:,:) = 0. 
     739 
     740   ! 3. Some controls 
     741   !----------------- 
     742   check_namelist = .TRUE. 
     743 
     744   IF( check_namelist ) THEN 
     745      ! Check time steps 
     746      IF( NINT(Agrif_Rhot()) * NINT(rdt) .NE. Agrif_Parent(rdt) ) THEN 
     747         WRITE(cl_check1,*)  Agrif_Parent(rdt) 
     748         WRITE(cl_check2,*)  rdt 
     749         WRITE(cl_check3,*)  rdt*Agrif_Rhot() 
    619750         CALL ctl_stop( 'incompatible time step between grids',   & 
    620751               &               'parent grid value : '//cl_check1    ,   &  
     
    635766         nitend =  Agrif_Parent(nitend)   *Agrif_IRhot() 
    636767      ENDIF 
    637  
    638768   ENDIF 
    639769   ! 
     
    655785      !!---------------------------------------------------------------------- 
    656786 
     787 
     788 
     789!RB_CMEMS : declare here init for top       
    657790      ! 1. Declaration of the type of variable which have to be interpolated 
    658791      !--------------------------------------------------------------------- 
    659792      ind1 =     nbghostcells 
    660       ind2 = 1 + nbghostcells 
    661       ind3 = 2 + nbghostcells 
     793      ind2 = 2 + nbghostcells_x 
     794      ind3 = 2 + nbghostcells_y_s 
    662795# if defined key_vertical 
    663       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) 
    664       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) 
     796      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) 
     797      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) 
    665798# else 
    666       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) 
    667       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) 
     799! LAURENT: STRANGE why (3,3) here ? 
     800      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) 
     801      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) 
    668802# endif 
    669803 
     
    705839      !!                     *** ROUTINE agrif_init *** 
    706840      !!---------------------------------------------------------------------- 
    707       USE agrif_oce  
    708       USE agrif_ice 
    709       USE in_out_manager 
    710       USE lib_mpp 
     841   USE agrif_oce  
     842   USE agrif_ice 
     843   USE dom_oce 
     844   USE in_out_manager 
     845   USE lib_mpp 
    711846      !! 
    712847      IMPLICIT NONE 
    713848      ! 
    714849      INTEGER  ::   ios                 ! Local integer output status for namelist read 
    715       NAMELIST/namagrif/ ln_agrif_2way, rn_sponge_tra, rn_sponge_dyn, rn_trelax_tra, rn_trelax_dyn, & 
     850      NAMELIST/namagrif/ ln_agrif_2way, ln_init_chfrpar, rn_sponge_tra, rn_sponge_dyn, rn_trelax_tra, rn_trelax_dyn, & 
    716851                       & ln_spc_dyn, ln_chk_bathy 
    717852      !!-------------------------------------------------------------------------------------- 
     
    729864         WRITE(numout,*) '   Namelist namagrif : set AGRIF parameters' 
    730865         WRITE(numout,*) '      Two way nesting activated ln_agrif_2way         = ', ln_agrif_2way 
    731          WRITE(numout,*) '      sponge coefficient for tracers    rn_sponge_tra = ', rn_sponge_tra, ' m^2/s' 
    732          WRITE(numout,*) '      sponge coefficient for dynamics   rn_sponge_tra = ', rn_sponge_dyn, ' m^2/s' 
    733          WRITE(numout,*) '      time relaxation for tracers       rn_trelax_tra = ', rn_trelax_tra, ' ad.' 
    734          WRITE(numout,*) '      time relaxation for dynamics      rn_trelax_dyn = ', rn_trelax_dyn, ' ad.' 
     866         WRITE(numout,*) '      child initial state from parent ln_init_chfrpar = ', ln_init_chfrpar 
     867         WRITE(numout,*) '      ad. sponge coeft for tracers      rn_sponge_tra = ', rn_sponge_tra 
     868         WRITE(numout,*) '      ad. sponge coeft for dynamics     rn_sponge_tra = ', rn_sponge_dyn 
     869         WRITE(numout,*) '      ad. time relaxation for tracers   rn_trelax_tra = ', rn_trelax_tra 
     870         WRITE(numout,*) '      ad. time relaxation for dynamics  rn_trelax_dyn = ', rn_trelax_dyn 
    735871         WRITE(numout,*) '      use special values for dynamics   ln_spc_dyn    = ', ln_spc_dyn 
    736872         WRITE(numout,*) '      check bathymetry                  ln_chk_bathy  = ', ln_chk_bathy 
    737873      ENDIF 
    738       ! 
    739       ! 
    740       IF( agrif_oce_alloc()  > 0 )   CALL ctl_warn('agrif agrif_oce_alloc: allocation of arrays failed') 
     874 
     875      lk_west  = .NOT. ( Agrif_Ix() == 1 ) 
     876      lk_east  = .NOT. ( Agrif_Ix() + nbcellsx/AGRIF_Irhox() == Agrif_Parent(jpiglo) -1 ) 
     877      lk_south = .NOT. ( Agrif_Iy() == 1 ) 
     878      lk_north = .NOT. ( Agrif_Iy() + nbcellsy/AGRIF_Irhoy() == Agrif_Parent(jpjglo) -1 ) 
     879 
     880      ! 
     881      ! Set the number of ghost cells according to periodicity 
     882      nbghostcells_x = nbghostcells 
     883      nbghostcells_y_s = nbghostcells 
     884      nbghostcells_y_n = nbghostcells 
     885      ! 
     886      IF ( jperio == 1 ) nbghostcells_x = 0 
     887      IF ( .NOT. lk_south ) nbghostcells_y_s = 0 
     888 
     889      ! Some checks 
     890      IF( jpiglo /= nbcellsx + 2 + 2*nbghostcells_x )   & 
     891          CALL ctl_stop( 'STOP', 'agrif_nemo_init: Agrif children requires jpiglo == nbcellsx + 2 + 2*nbghostcells_x' ) 
     892      IF( jpjglo /= nbcellsy + 2 + nbghostcells_y_s + nbghostcells_y_n )   & 
     893          CALL ctl_stop( 'STOP', 'agrif_nemo_init: Agrif children requires jpjglo == nbcellsy + 2 + nbghostcells_y_s + nbghostcells_y_n' ) 
     894      IF( ln_use_jattr )   CALL ctl_stop( 'STOP', 'agrif_nemo_init:Agrif children requires ln_use_jattr = .false. ' ) 
    741895      ! 
    742896   END SUBROUTINE agrif_nemo_init 
    743897 
    744898# if defined key_mpp_mpi 
    745  
    746899   SUBROUTINE Agrif_InvLoc( indloc, nprocloc, i, indglob ) 
    747900      !!---------------------------------------------------------------------- 
     
    803956# endif 
    804957 
     958   SUBROUTINE nemo_mapping(ndim,ptx,pty,bounds,bounds_chunks,correction_required,nb_chunks) 
     959      !!---------------------------------------------------------------------- 
     960      !!                   *** ROUTINE Nemo_mapping *** 
     961      !!---------------------------------------------------------------------- 
     962      USE dom_oce 
     963      !! 
     964      IMPLICIT NONE 
     965      ! 
     966      INTEGER :: ndim 
     967      INTEGER :: ptx, pty 
     968      INTEGER, DIMENSION(ndim,2,2) :: bounds 
     969      INTEGER, DIMENSION(:,:,:,:), ALLOCATABLE :: bounds_chunks 
     970      LOGICAL, DIMENSION(:), ALLOCATABLE :: correction_required 
     971      INTEGER :: nb_chunks 
     972      ! 
     973      INTEGER :: i 
     974 
     975      IF (agrif_debug_interp) THEN 
     976         DO i=1,ndim 
     977            WRITE(*,*) 'direction = ',i,bounds(i,1,2),bounds(i,2,2) 
     978         ENDDO 
     979      ENDIF 
     980 
     981      IF( bounds(2,2,2) > jpjglo) THEN 
     982         IF( bounds(2,1,2) <=jpjglo) THEN 
     983            nb_chunks = 2 
     984            ALLOCATE(bounds_chunks(nb_chunks,ndim,2,2)) 
     985            ALLOCATE(correction_required(nb_chunks)) 
     986            DO i = 1,nb_chunks 
     987               bounds_chunks(i,:,:,:) = bounds 
     988            END DO 
     989         
     990      ! FIRST CHUNCK (for j<=jpjglo)    
     991 
     992      ! Original indices 
     993            bounds_chunks(1,1,1,1) = bounds(1,1,2) 
     994            bounds_chunks(1,1,2,1) = bounds(1,2,2) 
     995            bounds_chunks(1,2,1,1) = bounds(2,1,2) 
     996            bounds_chunks(1,2,2,1) = jpjglo 
     997 
     998            bounds_chunks(1,1,1,2) = bounds(1,1,2) 
     999            bounds_chunks(1,1,2,2) = bounds(1,2,2) 
     1000            bounds_chunks(1,2,1,2) = bounds(2,1,2) 
     1001            bounds_chunks(1,2,2,2) = jpjglo 
     1002 
     1003      ! Correction required or not 
     1004            correction_required(1)=.FALSE. 
     1005        
     1006      ! SECOND CHUNCK (for j>jpjglo) 
     1007 
     1008      ! Original indices 
     1009            bounds_chunks(2,1,1,1) = bounds(1,1,2) 
     1010            bounds_chunks(2,1,2,1) = bounds(1,2,2) 
     1011            bounds_chunks(2,2,1,1) = jpjglo-2 
     1012            bounds_chunks(2,2,2,1) = bounds(2,2,2) 
     1013 
     1014      ! Where to find them 
     1015      ! We use the relation TAB(ji,jj)=TAB(jpiglo-ji+2,jpjglo-2-(jj-jpjglo)) 
     1016 
     1017            IF( ptx == 2) THEN ! T, V points 
     1018               bounds_chunks(2,1,1,2) = jpiglo-bounds(1,2,2)+2 
     1019               bounds_chunks(2,1,2,2) = jpiglo-bounds(1,1,2)+2 
     1020            ELSE ! U, F points 
     1021               bounds_chunks(2,1,1,2) = jpiglo-bounds(1,2,2)+1 
     1022               bounds_chunks(2,1,2,2) = jpiglo-bounds(1,1,2)+1        
     1023            ENDIF 
     1024 
     1025            IF( pty == 2) THEN ! T, U points 
     1026               bounds_chunks(2,2,1,2) = jpjglo-2-(bounds(2,2,2) -jpjglo) 
     1027               bounds_chunks(2,2,2,2) = jpjglo-2-(jpjglo-2      -jpjglo) 
     1028            ELSE ! V, F points 
     1029               bounds_chunks(2,2,1,2) = jpjglo-3-(bounds(2,2,2) -jpjglo) 
     1030               bounds_chunks(2,2,2,2) = jpjglo-3-(jpjglo-2      -jpjglo) 
     1031            ENDIF 
     1032      ! Correction required or not 
     1033            correction_required(2)=.TRUE. 
     1034 
     1035         ELSE 
     1036            nb_chunks = 1 
     1037            ALLOCATE(bounds_chunks(nb_chunks,ndim,2,2)) 
     1038            ALLOCATE(correction_required(nb_chunks)) 
     1039            DO i=1,nb_chunks 
     1040               bounds_chunks(i,:,:,:) = bounds 
     1041            END DO 
     1042 
     1043            bounds_chunks(1,1,1,1) = bounds(1,1,2) 
     1044            bounds_chunks(1,1,2,1) = bounds(1,2,2) 
     1045            bounds_chunks(1,2,1,1) = bounds(2,1,2) 
     1046            bounds_chunks(1,2,2,1) = bounds(2,2,2) 
     1047 
     1048            bounds_chunks(1,1,1,2) = jpiglo-bounds(1,2,2)+2 
     1049            bounds_chunks(1,1,2,2) = jpiglo-bounds(1,1,2)+2 
     1050 
     1051            bounds_chunks(1,2,1,2) = jpjglo-2-(bounds(2,2,2)-jpjglo) 
     1052            bounds_chunks(1,2,2,2) = jpjglo-2-(bounds(2,1,2)-jpjglo) 
     1053 
     1054            IF( ptx == 2) THEN ! T, V points 
     1055               bounds_chunks(1,1,1,2) = jpiglo-bounds(1,2,2)+2 
     1056               bounds_chunks(1,1,2,2) = jpiglo-bounds(1,1,2)+2 
     1057            ELSE ! U, F points 
     1058               bounds_chunks(1,1,1,2) = jpiglo-bounds(1,2,2)+1 
     1059               bounds_chunks(1,1,2,2) = jpiglo-bounds(1,1,2)+1        
     1060            ENDIF 
     1061 
     1062            IF (pty == 2) THEN ! T, U points 
     1063               bounds_chunks(1,2,1,2) = jpjglo-2-(bounds(2,2,2) -jpjglo) 
     1064               bounds_chunks(1,2,2,2) = jpjglo-2-(bounds(2,1,2) -jpjglo) 
     1065            ELSE ! V, F points 
     1066               bounds_chunks(1,2,1,2) = jpjglo-3-(bounds(2,2,2) -jpjglo) 
     1067               bounds_chunks(1,2,2,2) = jpjglo-3-(bounds(2,1,2) -jpjglo) 
     1068            ENDIF 
     1069 
     1070            correction_required(1)=.TRUE.           
     1071         ENDIF 
     1072 
     1073      ELSE IF (bounds(1,1,2) < 1) THEN 
     1074         IF (bounds(1,2,2) > 0) THEN 
     1075            nb_chunks = 2 
     1076            ALLOCATE(correction_required(nb_chunks)) 
     1077            correction_required=.FALSE. 
     1078            ALLOCATE(bounds_chunks(nb_chunks,ndim,2,2)) 
     1079            DO i=1,nb_chunks 
     1080               bounds_chunks(i,:,:,:) = bounds 
     1081            END DO 
     1082               
     1083            bounds_chunks(1,1,1,2) = bounds(1,1,2)+jpiglo-2 
     1084            bounds_chunks(1,1,2,2) = 1+jpiglo-2 
     1085           
     1086            bounds_chunks(1,1,1,1) = bounds(1,1,2) 
     1087            bounds_chunks(1,1,2,1) = 1 
     1088        
     1089            bounds_chunks(2,1,1,2) = 2 
     1090            bounds_chunks(2,1,2,2) = bounds(1,2,2) 
     1091           
     1092            bounds_chunks(2,1,1,1) = 2 
     1093            bounds_chunks(2,1,2,1) = bounds(1,2,2) 
     1094 
     1095         ELSE 
     1096            nb_chunks = 1 
     1097            ALLOCATE(correction_required(nb_chunks)) 
     1098            correction_required=.FALSE. 
     1099            ALLOCATE(bounds_chunks(nb_chunks,ndim,2,2)) 
     1100            DO i=1,nb_chunks 
     1101               bounds_chunks(i,:,:,:) = bounds 
     1102            END DO     
     1103            bounds_chunks(1,1,1,2) = bounds(1,1,2)+jpiglo-2 
     1104            bounds_chunks(1,1,2,2) = bounds(1,2,2)+jpiglo-2 
     1105           
     1106            bounds_chunks(1,1,1,1) = bounds(1,1,2) 
     1107           bounds_chunks(1,1,2,1) = bounds(1,2,2) 
     1108         ENDIF 
     1109      ELSE 
     1110         nb_chunks=1   
     1111         ALLOCATE(correction_required(nb_chunks)) 
     1112         correction_required=.FALSE. 
     1113         ALLOCATE(bounds_chunks(nb_chunks,ndim,2,2)) 
     1114         DO i=1,nb_chunks 
     1115            bounds_chunks(i,:,:,:) = bounds 
     1116         END DO 
     1117         bounds_chunks(1,1,1,2) = bounds(1,1,2) 
     1118         bounds_chunks(1,1,2,2) = bounds(1,2,2) 
     1119         bounds_chunks(1,2,1,2) = bounds(2,1,2) 
     1120         bounds_chunks(1,2,2,2) = bounds(2,2,2) 
     1121           
     1122         bounds_chunks(1,1,1,1) = bounds(1,1,2) 
     1123         bounds_chunks(1,1,2,1) = bounds(1,2,2) 
     1124         bounds_chunks(1,2,1,1) = bounds(2,1,2) 
     1125         bounds_chunks(1,2,2,1) = bounds(2,2,2)               
     1126      ENDIF 
     1127         
     1128   END SUBROUTINE nemo_mapping 
     1129 
     1130   FUNCTION agrif_external_switch_index(ptx,pty,i1,isens) 
     1131 
     1132   USE dom_oce 
     1133 
     1134   INTEGER :: ptx, pty, i1, isens 
     1135   INTEGER :: agrif_external_switch_index 
     1136 
     1137   IF( isens == 1 ) THEN 
     1138      IF( ptx == 2 ) THEN ! T, V points 
     1139         agrif_external_switch_index = jpiglo-i1+2 
     1140      ELSE ! U, F points 
     1141         agrif_external_switch_index = jpiglo-i1+1       
     1142      ENDIF 
     1143   ELSE IF( isens ==2 ) THEN 
     1144      IF ( pty == 2 ) THEN ! T, U points 
     1145         agrif_external_switch_index = jpjglo-2-(i1 -jpjglo) 
     1146      ELSE ! V, F points 
     1147         agrif_external_switch_index = jpjglo-3-(i1 -jpjglo) 
     1148      ENDIF 
     1149   ENDIF 
     1150 
     1151   END FUNCTION agrif_external_switch_index 
     1152 
     1153   SUBROUTINE Correct_field(tab2d,i1,i2,j1,j2) 
     1154      !!---------------------------------------------------------------------- 
     1155      !!                   *** ROUTINE Correct_field *** 
     1156      !!---------------------------------------------------------------------- 
     1157    
     1158   USE dom_oce 
     1159   USE agrif_oce 
     1160 
     1161   INTEGER :: i1,i2,j1,j2 
     1162   REAL(wp), DIMENSION(i1:i2,j1:j2) :: tab2d 
     1163 
     1164   INTEGER :: i,j 
     1165   REAL(wp), DIMENSION(i1:i2,j1:j2) :: tab2dtemp 
     1166 
     1167   tab2dtemp = tab2d 
     1168 
     1169   IF( .NOT. use_sign_north ) THEN 
     1170      DO j=j1,j2 
     1171         DO i=i1,i2 
     1172            tab2d(i,j)=tab2dtemp(i2-(i-i1),j2-(j-j1)) 
     1173         END DO 
     1174      END DO 
     1175   ELSE 
     1176      DO j=j1,j2 
     1177         DO i=i1,i2 
     1178            tab2d(i,j)=sign_north * tab2dtemp(i2-(i-i1),j2-(j-j1)) 
     1179         END DO 
     1180      END DO 
     1181   ENDIF 
     1182 
     1183   END SUBROUTINE Correct_field 
     1184 
    8051185#else 
    8061186   SUBROUTINE Subcalledbyagrif 
  • NEMO/trunk/src/OCE/DOM/dom_oce.F90

    r12933 r13216  
    215215#if defined key_agrif 
    216216   LOGICAL, PUBLIC, PARAMETER ::   lk_agrif = .TRUE.    !: agrif flag 
     217   LOGICAL, PUBLIC            ::   lk_south, lk_north, lk_west, lk_east !: Child grid boundaries (interpolation or not) 
    217218#else 
    218219   LOGICAL, PUBLIC, PARAMETER ::   lk_agrif = .FALSE.   !: agrif flag 
  • NEMO/trunk/src/OCE/DOM/domain.F90

    r12489 r13216  
    187187      ENDIF 
    188188      ! 
     189 
    189190      IF( lk_c1d         )   CALL cor_c1d       ! 1D configuration: Coriolis set at T-point 
    190191      ! 
     192 
     193#if defined key_agrif 
     194      IF( .NOT. Agrif_Root() ) CALL Agrif_Init_Domain( Kbb, Kmm, Kaa ) 
     195#endif 
    191196      IF( ln_meshmask    )   CALL dom_wri       ! Create a domain file 
    192197      IF( .NOT.ln_rstart )   CALL dom_ctl       ! Domain control 
     
    307312902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namrun in configuration namelist' ) 
    308313      IF(lwm) WRITE ( numond, namrun ) 
     314 
     315#if defined key_agrif 
     316      IF( .NOT. Agrif_Root() ) THEN 
     317            nn_it000 = (Agrif_Parent(nn_it000)-1)*Agrif_IRhot() + 1 
     318            nn_itend =  Agrif_Parent(nn_itend)   *Agrif_IRhot() 
     319      ENDIF 
     320#endif 
    309321      ! 
    310322      IF(lwp) THEN                  ! control print 
     
    403415      IF(lwm) WRITE( numond, namdom ) 
    404416      ! 
     417#if defined key_agrif 
     418      IF( .NOT. Agrif_Root() ) THEN 
     419            rn_Dt = Agrif_Parent(rn_Dt) / Agrif_Rhot() 
     420      ENDIF 
     421#endif 
     422      ! 
    405423      IF(lwp) THEN 
    406424         WRITE(numout,*) 
  • NEMO/trunk/src/OCE/DOM/istate.F90

    r13100 r13216  
    3434   USE lib_mpp         ! MPP library 
    3535   USE restart         ! restart 
     36#if defined key_agrif 
     37   USE agrif_oce_interp 
     38   USE agrif_oce 
     39#endif    
    3640 
    3741   IMPLICIT NONE 
     
    6973!!gm  Why not include in the first call of dta_tsd ?   
    7074!!gm  probably associated with the use of internal damping... 
    71                      CALL dta_tsd_init        ! Initialisation of T & S input data 
     75       CALL dta_tsd_init        ! Initialisation of T & S input data 
    7276!!gm to be moved in usrdef of C1D case 
    7377!      IF( lk_c1d )   CALL dta_uvd_init        ! Initialization of U & V input data 
     
    8387#endif 
    8488 
     89#if defined key_agrif 
     90      IF ( (.NOT.Agrif_root()).AND.ln_init_chfrpar ) THEN 
     91         numror = 0                           ! define numror = 0 -> no restart file to read 
     92         ln_1st_euler = .true.                ! Set time-step indicator at nit000 (euler forward) 
     93         CALL day_init  
     94         CALL agrif_istate( Kbb, Kmm, Kaa )   ! Interp from parent 
     95         ! 
     96         ts  (:,:,:,:,Kmm) = ts (:,:,:,:,Kbb)  
     97         ssh (:,:,Kmm)     = ssh(:,:,Kbb) 
     98         uu   (:,:,:,Kmm)   = uu  (:,:,:,Kbb) 
     99         vv   (:,:,:,Kmm)   = vv  (:,:,:,Kbb) 
     100      ELSE 
     101#endif 
    85102      IF( ln_rstart ) THEN                    ! Restart from a file 
    86103         !                                    ! ------------------- 
     
    99116            ! 
    100117            ssh(:,:,Kbb)   = 0._wp               ! set the ocean at rest 
     118            uu  (:,:,:,Kbb) = 0._wp 
     119            vv  (:,:,:,Kbb) = 0._wp   
     120            ! 
    101121            IF( ll_wd ) THEN 
    102122               ssh(:,:,Kbb) =  -ssh_ref  ! Added in 30 here for bathy that adds 30 as Iterative test CEOD  
     
    110130               END_2D 
    111131            ENDIF  
    112             uu  (:,:,:,Kbb) = 0._wp 
    113             vv  (:,:,:,Kbb) = 0._wp   
    114             ! 
     132             ! 
    115133         ELSE                                 ! user defined initial T and S 
    116134            CALL usr_def_istate( gdept(:,:,:,Kbb), tmask, ts(:,:,:,:,Kbb), uu(:,:,:,Kbb), vv(:,:,:,Kbb), ssh(:,:,Kbb)  )          
     
    147165         !  
    148166      ENDIF  
     167#if defined key_agrif 
     168      ENDIF 
     169#endif 
    149170      !  
    150171      ! Initialize "now" and "before" barotropic velocities: 
  • NEMO/trunk/src/OCE/DYN/dynspg_ts.F90

    r12489 r13216  
    513513         END_2D 
    514514         ! 
     515         ! Duplicate sea level across open boundaries (this is only cosmetic if linssh=T) 
     516         IF( ln_bdy )   CALL bdy_ssh( ssha_e ) 
     517#if defined key_agrif 
     518         IF( .NOT.Agrif_Root() )   CALL agrif_ssh_ts( jn ) 
     519#endif 
     520 
    515521         CALL lbc_lnk_multi( 'dynspg_ts', ssha_e, 'T', 1._wp,  zhU, 'U', -1._wp,  zhV, 'V', -1._wp ) 
    516522         ! 
     
    525531         END IF 
    526532         ! 
    527          ! Duplicate sea level across open boundaries (this is only cosmetic if linssh=T) 
    528          IF( ln_bdy )   CALL bdy_ssh( ssha_e ) 
    529 #if defined key_agrif 
    530          IF( .NOT.Agrif_Root() )   CALL agrif_ssh_ts( jn ) 
    531 #endif 
    532533         !   
    533534         ! Sea Surface Height at u-,v-points (vvl case only) 
     
    643644         ENDIF 
    644645        
    645          IF( .NOT.ln_linssh ) THEN   !* Update ocean depth (variable volume case only) 
     646         IF( .NOT.ln_linssh ) THEN !* Update ocean depth (variable volume case only) 
    646647            hu_e (2:jpim1,2:jpjm1) = hu_0(2:jpim1,2:jpjm1) + zsshu_a(2:jpim1,2:jpjm1) 
    647648            hur_e(2:jpim1,2:jpjm1) = ssumask(2:jpim1,2:jpjm1) / ( hu_e(2:jpim1,2:jpjm1) + 1._wp - ssumask(2:jpim1,2:jpjm1) ) 
    648649            hv_e (2:jpim1,2:jpjm1) = hv_0(2:jpim1,2:jpjm1) + zsshv_a(2:jpim1,2:jpjm1) 
    649650            hvr_e(2:jpim1,2:jpjm1) = ssvmask(2:jpim1,2:jpjm1) / ( hv_e(2:jpim1,2:jpjm1) + 1._wp - ssvmask(2:jpim1,2:jpjm1) ) 
     651         ENDIF 
     652         !                                                 ! open boundaries 
     653         IF( ln_bdy )   CALL bdy_dyn2d( jn, ua_e, va_e, un_e, vn_e, hur_e, hvr_e, ssha_e ) 
     654#if defined key_agrif                                                            
     655         IF( .NOT.Agrif_Root() )  CALL agrif_dyn_ts( jn )  ! Agrif 
     656#endif 
     657         ! 
     658         IF( .NOT.ln_linssh ) THEN   !* Update ocean depth (variable volume case only) 
    650659            CALL lbc_lnk_multi( 'dynspg_ts', ua_e , 'U', -1._wp, va_e , 'V', -1._wp  & 
    651660                 &                         , hu_e , 'U',  1._wp, hv_e , 'V',  1._wp  & 
     
    655664         ENDIF 
    656665         ! 
    657          ! 
    658          !                                                 ! open boundaries 
    659          IF( ln_bdy )   CALL bdy_dyn2d( jn, ua_e, va_e, un_e, vn_e, hur_e, hvr_e, ssha_e ) 
    660 #if defined key_agrif                                                            
    661          IF( .NOT.Agrif_Root() )  CALL agrif_dyn_ts( jn )  ! Agrif 
    662 #endif 
    663666         !                                             !* Swap 
    664667         !                                             !  ---- 
  • NEMO/trunk/src/OCE/DYN/sshwzv.F90

    r12965 r13216  
    204204         ! Mask vertical velocity at first/last columns/row  
    205205         ! inside computational domain (cosmetic)  
    206          ! --- West --- ! 
    207          DO ji = mi0(2), mi1(2) 
    208             DO jj = 1, jpj 
    209                pww(ji,jj,:) = 0._wp  
     206         ! --- West --- !          
     207         IF( lk_west) THEN 
     208            DO ji = mi0(2), mi1(2) 
     209               DO jj = 1, jpj 
     210                  pww(ji,jj,:) = 0._wp  
     211               ENDDO 
    210212            ENDDO 
    211          ENDDO 
     213         ENDIF 
    212214         ! 
    213215         ! --- East --- ! 
    214          DO ji = mi0(jpiglo-1), mi1(jpiglo-1) 
    215             DO jj = 1, jpj 
    216                pww(ji,jj,:) = 0._wp 
     216         IF( lk_east) THEN 
     217            DO ji = mi0(jpiglo-1), mi1(jpiglo-1) 
     218               DO jj = 1, jpj 
     219                  pww(ji,jj,:) = 0._wp 
     220               ENDDO 
    217221            ENDDO 
    218          ENDDO 
     222         ENDIF 
    219223         ! 
    220224         ! --- South --- ! 
    221          DO jj = mj0(2), mj1(2) 
    222             DO ji = 1, jpi 
    223                pww(ji,jj,:) = 0._wp 
     225         IF( lk_south) THEN 
     226            DO jj = mj0(2), mj1(2) 
     227               DO ji = 1, jpi 
     228                  pww(ji,jj,:) = 0._wp 
     229               ENDDO 
    224230            ENDDO 
    225          ENDDO 
     231         ENDIF 
    226232         ! 
    227233         ! --- North --- ! 
    228          DO jj = mj0(jpjglo-1), mj1(jpjglo-1) 
    229             DO ji = 1, jpi 
    230                pww(ji,jj,:) = 0._wp 
     234         IF( lk_north) THEN 
     235            DO jj = mj0(jpjglo-1), mj1(jpjglo-1) 
     236               DO ji = 1, jpi 
     237                  pww(ji,jj,:) = 0._wp 
     238               ENDDO 
    231239            ENDDO 
    232          ENDDO 
     240         ENDIF 
     241         ! 
    233242      ENDIF  
    234243#endif  
  • NEMO/trunk/src/OCE/FLO/floblk.F90

    r12649 r13216  
    4141      INTEGER, INTENT( in  ) ::   Kbb, Kmm ! ocean time level indices 
    4242      !! 
     43#ifndef key_agrif 
     44 
     45!RB super quick fix to compile with agrif 
     46 
    4347      INTEGER :: jfl              ! dummy loop arguments 
    4448      INTEGER :: ind, ifin, iloop 
     
    364368         GO TO 222 
    365369      ENDIF 
     370#endif 
    366371      ! 
    367372      ! 
  • NEMO/trunk/src/OCE/LBC/lib_mpp.F90

    r13062 r13216  
    137137 
    138138   ! Communications summary report 
    139    CHARACTER(len=128), DIMENSION(:), ALLOCATABLE ::   crname_lbc                   !: names of lbc_lnk calling routines 
    140    CHARACTER(len=128), DIMENSION(:), ALLOCATABLE ::   crname_glb                   !: names of global comm calling routines 
    141    CHARACTER(len=128), DIMENSION(:), ALLOCATABLE ::   crname_dlg                   !: names of delayed global comm calling routines 
     139   CHARACTER(len=lca), DIMENSION(:), ALLOCATABLE ::   crname_lbc                   !: names of lbc_lnk calling routines 
     140   CHARACTER(len=lca), DIMENSION(:), ALLOCATABLE ::   crname_glb                   !: names of global comm calling routines 
     141   CHARACTER(len=lca), DIMENSION(:), ALLOCATABLE ::   crname_dlg                   !: names of delayed global comm calling routines 
    142142   INTEGER, PUBLIC                               ::   ncom_stp = 0                 !: copy of time step # istp 
    143143   INTEGER, PUBLIC                               ::   ncom_fsbc = 1                !: copy of sbc time step # nn_fsbc 
  • NEMO/trunk/src/OCE/LBC/mppini.F90

    r12377 r13216  
    102102            &           'the domain is lay out for distributed memory computing!' ) 
    103103         ! 
     104#if defined key_agrif 
     105    IF (.NOT.agrif_root()) THEN 
     106      call agrif_nemo_init() 
     107    ENDIF 
     108#endif 
    104109   END SUBROUTINE mpp_init 
    105110 
     
    333338#if defined key_agrif 
    334339      IF( .NOT. Agrif_Root() ) THEN       ! AGRIF children: specific setting (cf. agrif_user.F90) 
    335          IF( jpiglo /= nbcellsx + 2 + 2*nbghostcells )   & 
    336             CALL ctl_stop( 'STOP', 'mpp_init: Agrif children requires jpiglo == nbcellsx + 2 + 2*nbghostcells' ) 
    337          IF( jpjglo /= nbcellsy + 2 + 2*nbghostcells )   & 
    338             CALL ctl_stop( 'STOP', 'mpp_init: Agrif children requires jpjglo == nbcellsy + 2 + 2*nbghostcells' ) 
    339          IF( ln_use_jattr )   CALL ctl_stop( 'STOP', 'mpp_init:Agrif children requires ln_use_jattr = .false. ' ) 
     340         CALL agrif_nemo_init() 
    340341      ENDIF 
    341342#endif 
  • NEMO/trunk/src/OCE/OBS/diaobs.F90

    r12489 r13216  
    9494   TYPE(obs_prof), PUBLIC, POINTER, DIMENSION(:) ::   profdataqc   !: Profile data after quality control 
    9595 
    96    CHARACTER(len=6), PUBLIC, DIMENSION(:), ALLOCATABLE ::   cobstypesprof, cobstypessurf   !: Profile & surface obs types 
     96   CHARACTER(len=lca), PUBLIC, DIMENSION(:), ALLOCATABLE ::   cobstypesprof, cobstypessurf   !: Profile & surface obs types 
    9797 
    9898   !!---------------------------------------------------------------------- 
  • NEMO/trunk/src/OCE/SBC/sbcmod.F90

    r12994 r13216  
    121121#endif 
    122122      ! 
     123      ! 
    123124      IF(lwp) THEN                  !* Control print 
    124125         WRITE(numout,*) '   Namelist namsbc (partly overwritten with CPP key setting)' 
  • NEMO/trunk/src/OCE/STO/stopar.F90

    r12933 r13216  
    5656   INTEGER,          DIMENSION(:),       ALLOCATABLE :: sto3d_ord  ! order of autoregressive process 
    5757 
    58    CHARACTER(len=1), DIMENSION(:),       ALLOCATABLE :: sto2d_typ  ! nature of grid point (T, U, V, W, F, I) 
    59    CHARACTER(len=1), DIMENSION(:),       ALLOCATABLE :: sto3d_typ  ! nature of grid point (T, U, V, W, F, I) 
     58   CHARACTER(len=lca), DIMENSION(:),       ALLOCATABLE :: sto2d_typ  ! nature of grid point (T, U, V, W, F, I) 
     59   CHARACTER(len=lca), DIMENSION(:),       ALLOCATABLE :: sto3d_typ  ! nature of grid point (T, U, V, W, F, I) 
    6060   REAL(wp),         DIMENSION(:),       ALLOCATABLE :: sto2d_sgn  ! control of the sign accross the north fold 
    6161   REAL(wp),         DIMENSION(:),       ALLOCATABLE :: sto3d_sgn  ! control of the sign accross the north fold 
  • NEMO/trunk/src/OCE/USR/usrdef_hgr.F90

    r12489 r13216  
    9595#if defined key_agrif 
    9696      ! ! Upper left longitude and latitude from parent: 
     97      ! Laurent: Should be modify in case of an east-west cyclic parent grid 
    9798      IF (.NOT.Agrif_root()) THEN 
    9899         zlam0 = zlam1 + Agrif_irhox() * REAL(Agrif_Parent(jpjglo)-2 , wp) * ze1deg * zcos_alpha  & 
  • NEMO/trunk/src/OCE/USR/usrdef_nam.F90

    r12377 r13216  
    7474#if defined key_agrif 
    7575      IF( .NOT. Agrif_Root() ) THEN 
    76          kpi  = nbcellsx + 2 + 2*nbghostcells 
    77          kpj  = nbcellsy + 2 + 2*nbghostcells 
     76         kpi  = nbcellsx + 2 + 2*nbghostcells_x 
     77         kpj  = nbcellsy + 2 + 2*nbghostcells_y_s 
    7878      ENDIF 
    7979#endif 
  • NEMO/trunk/src/OCE/nemogcm.F90

    r13011 r13216  
    143143#if defined key_agrif 
    144144      Kbb_a = Nbb; Kmm_a = Nnn; Krhs_a = Nrhs   ! agrif_oce module copies of time level indices 
    145       CALL Agrif_Declare_Var_dom   ! AGRIF: set the meshes for DOM 
    146       CALL Agrif_Declare_Var       !  "      "   "   "      "  DYN/TRA  
     145      CALL Agrif_Declare_Var       !  "      "   "   "      "  DYN/TRA 
    147146# if defined key_top 
    148147      CALL Agrif_Declare_Var_top   !  "      "   "   "      "  TOP 
    149 # endif 
    150 # if defined key_si3 
    151       CALL Agrif_Declare_Var_ice   !  "      "   "   "      "  Sea ice 
    152148# endif 
    153149#endif 
     
    400396      ! Initialise time level indices 
    401397      Nbb = 1; Nnn = 2; Naa = 3; Nrhs = Naa 
    402  
     398#if defined key_agrif 
     399      Kbb_a = Nbb; Kmm_a = Nnn; Krhs_a = Nrhs   ! agrif_oce module copies of time level indices 
     400#endif  
    403401      !                             !-------------------------------! 
    404402      !                             !  NEMO general initialization  ! 
     
    415413      IF( lk_c1d       )   CALL     c1d_init        ! 1D column configuration 
    416414                           CALL     wad_init        ! Wetting and drying options 
     415 
     416#if defined key_agrif 
     417     CALL Agrif_Declare_Var_ini   !  "      "   "   "      "  DOM 
     418#endif 
    417419                           CALL     dom_init( Nbb, Nnn, Naa, "OPA") ! Domain 
     420 
     421 
     422 
    418423      IF( ln_crs       )   CALL     crs_init(      Nnn )       ! coarsened grid: domain initialization  
    419424      IF( sn_cfctl%l_prtctl )   & 
     
    436441      ENDIF 
    437442      ! 
    438        
     443 
    439444                           CALL  istate_init( Nbb, Nnn, Naa )    ! ocean initial state (Dynamics and tracers) 
    440445 
  • NEMO/trunk/src/OCE/par_kind.F90

    r10068 r13216  
    3131    
    3232   !                                                                !!** Integer ** 
    33    INTEGER, PUBLIC, PARAMETER ::   lc = 256                          !: Lenght of Character strings 
     33   INTEGER, PUBLIC, PARAMETER ::   lc  = 256                          !: Lenght of Character strings 
     34   INTEGER, PUBLIC, PARAMETER ::   lca = 400                          !: Lenght of Character arrays 
    3435 
    3536   !!---------------------------------------------------------------------- 
  • NEMO/trunk/src/OCE/par_oce.F90

    r12377 r13216  
    4747   ! global domain size for AGRIF     !!! * total AGRIF computational domain * 
    4848   INTEGER, PUBLIC            ::   nbug_in_agrif_conv_do_not_remove_or_modify = 1 - 1 
    49    INTEGER, PUBLIC, PARAMETER ::   nbghostcells = 3                             !: number of ghost cells 
    50    INTEGER, PUBLIC            ::   nbcellsx   ! = jpiglo - 2 - 2*nbghostcells   !: number of cells in i-direction 
    51    INTEGER, PUBLIC            ::   nbcellsy   ! = jpjglo - 2 - 2*nbghostcells   !: number of cells in j-direction 
     49   INTEGER, PUBLIC, PARAMETER ::   nbghostcells = 3 !: number of ghost cells: default value 
     50   INTEGER, PUBLIC            ::   nbghostcells_x   !: number of ghost cells in i-direction 
     51   INTEGER, PUBLIC            ::   nbghostcells_y_s   !: number of ghost cells in j-direction at south 
     52   INTEGER, PUBLIC            ::   nbghostcells_y_n   !: number of ghost cells in j-direction at north                       !: number of ghost cells 
     53   INTEGER, PUBLIC            ::   nbcellsx   ! = jpiglo - 2 - 2*nbghostcells_x   !: number of cells in i-direction 
     54   INTEGER, PUBLIC            ::   nbcellsy   ! = jpjglo - 2 - 2*nbghostcells-y   !: number of cells in j-direction 
    5255 
    5356   ! local domain size                !!! * local computational domain * 
  • NEMO/trunk/src/OCE/stpctl.F90

    r13136 r13216  
    119119      !                                   !==            test of local extrema           ==! 
    120120      !                                   !==  done by all processes at every time step  ==! 
    121       ! 
    122       ! define zmax default value. needed for land processors 
    123       IF( ll_colruns ) THEN    ! default value: must not be kept when calling mpp_max -> must be as small as possible 
    124          zmax(:) = -HUGE(1._wp) 
    125       ELSE                     ! default value: must not give true for any of the tests bellow (-> avoid manipulating HUGE...) 
    126          zmax(:) =  0._wp 
    127          zmax(3) = -1._wp      ! avoid salinity minimum at 0. 
    128       ENDIF 
    129       ! 
    130121      llmsk(:,:,1) = ssmask(:,:) == 1._wp 
    131       IF( COUNT( llmsk(:,:,1) ) > 0 ) THEN   ! avoid huge values sent back for land processors... 
    132          IF( ll_wd ) THEN 
    133             zmax(1) = MAXVAL( ABS( ssh(:,:,Kmm) + ssh_ref ), mask = llmsk(:,:,1) )   ! ssh max 
     122      IF( ll_wd ) THEN 
     123         zmax(1) = MAXVAL( ABS( ssh(:,:,Kmm) + ssh_ref ), mask = llmsk(:,:,1) )   ! ssh max 
     124      ELSE 
     125         zmax(1) = MAXVAL( ABS( ssh(:,:,Kmm)           ), mask = llmsk(:,:,1) )   ! ssh max 
     126      ENDIF 
     127      llmsk(:,:,:) = umask(:,:,:) == 1._wp 
     128      zmax(2) = MAXVAL(  ABS( uu(:,:,:,Kmm) ), mask = llmsk )                     ! velocity max (zonal only) 
     129      llmsk(:,:,:) = tmask(:,:,:) == 1._wp 
     130      zmax(3) = MAXVAL( -ts(:,:,:,jp_sal,Kmm), mask = llmsk )                     ! minus salinity max 
     131      zmax(4) = MAXVAL(  ts(:,:,:,jp_sal,Kmm), mask = llmsk )                     !       salinity max 
     132      IF( ll_colruns .OR. jpnij == 1 ) THEN     ! following variables are used only in the netcdf file 
     133         zmax(5) = MAXVAL( -ts(:,:,:,jp_tem,Kmm), mask = llmsk )                  ! minus temperature max 
     134         zmax(6) = MAXVAL(  ts(:,:,:,jp_tem,Kmm), mask = llmsk )                  !       temperature max 
     135         IF( ln_zad_Aimp ) THEN 
     136            zmax(7) = MAXVAL(   Cu_adv(:,:,:)   , mask = llmsk )                  ! partitioning coeff. max 
     137            llmsk(:,:,:) = wmask(:,:,:) == 1._wp 
     138            zmax(8) = MAXVAL(  ABS( wi(:,:,:) ) , mask = llmsk )                  ! implicit vertical vel. max 
    134139         ELSE 
    135             zmax(1) = MAXVAL( ABS( ssh(:,:,Kmm)           ), mask = llmsk(:,:,1) )   ! ssh max 
    136          ENDIF 
    137       ENDIF 
    138       zmax(2) = MAXVAL( ABS( uu(:,:,:,Kmm) ) )                                       ! velocity max (zonal only) 
    139       llmsk(:,:,:) = tmask(:,:,:) == 1._wp 
    140       IF( COUNT( llmsk(:,:,:) ) > 0 ) THEN   ! avoid huge values sent back for land processors... 
    141          zmax(3) = MAXVAL( -ts(:,:,:,jp_sal,Kmm), mask = llmsk )                     ! minus salinity max 
    142          zmax(4) = MAXVAL(  ts(:,:,:,jp_sal,Kmm), mask = llmsk )                     !       salinity max 
    143          IF( ll_colruns .OR. jpnij == 1 ) THEN     ! following variables are used only in the netcdf file 
    144             zmax(5) = MAXVAL( -ts(:,:,:,jp_tem,Kmm), mask = llmsk )                  ! minus temperature max 
    145             zmax(6) = MAXVAL(  ts(:,:,:,jp_tem,Kmm), mask = llmsk )                  !       temperature max 
    146             IF( ln_zad_Aimp ) THEN 
    147                zmax(7) = MAXVAL(   Cu_adv(:,:,:)   , mask = llmsk )                  ! partitioning coeff. max 
    148                llmsk(:,:,:) = wmask(:,:,:) == 1._wp 
    149                IF( COUNT( llmsk(:,:,:) ) > 0 ) THEN   ! avoid huge values sent back for land processors... 
    150                   zmax(8) = MAXVAL(ABS( wi(:,:,:) ), mask = llmsk )                  ! implicit vertical vel. max 
    151                ENDIF 
    152             ENDIF 
    153          ENDIF 
     140            zmax(7:8) = 0._wp 
     141         ENDIF 
     142      ELSE 
     143         zmax(5:8) = 0._wp 
    154144      ENDIF 
    155145      zmax(9) = REAL( nstop, wp )                                              ! stop indicator 
  • NEMO/trunk/src/SAS/nemogcm.F90

    r13011 r13216  
    9191#if defined key_agrif 
    9292      Kbb_a = Nbb; Kmm_a = Nnn; Krhs_a = Nrhs   ! agrif_oce module copies of time level indices 
    93       CALL Agrif_Declare_Var_dom   ! AGRIF: set the meshes for DOM 
    9493      CALL Agrif_Declare_Var       !  "      "   "   "      "  DYN/TRA  
    9594# if defined key_top 
    9695      CALL Agrif_Declare_Var_top   !  "      "   "   "      "  TOP 
    97 # endif 
    98 # if defined key_si3 
    99       CALL Agrif_Declare_Var_ice   !  "      "   "   "      "  Sea ice 
    10096# endif 
    10197#endif 
     
    345341      ! Initialise time level indices 
    346342      Nbb = 1; Nnn = 2; Naa = 3; Nrhs = Naa 
     343#if defined key_agrif 
     344      Kbb_a = Nbb; Kmm_a = Nnn; Krhs_a = Nrhs   ! agrif_oce module copies of time level indices 
     345#endif  
    347346 
    348347      !                             !-------------------------------! 
     
    358357                           CALL phy_cst         ! Physical constants 
    359358                           CALL eos_init        ! Equation of seawater 
     359#if defined key_agrif 
     360     CALL Agrif_Declare_Var_ini   !  "      "   "   "      "  DOM 
     361#endif 
    360362                           CALL dom_init( Nbb, Nnn, Naa, 'SAS') ! Domain 
    361363      IF( sn_cfctl%l_prtctl )   & 
  • NEMO/trunk/tests/CANAL/MY_SRC/stpctl.F90

    r13214 r13216  
    1919   USE dom_oce         ! ocean space and time domain variables  
    2020   USE c1d             ! 1D vertical configuration 
     21   USE zdf_oce ,  ONLY : ln_zad_Aimp       ! ocean vertical physics variables 
     22   USE wet_dry,   ONLY : ll_wd, ssh_ref    ! reference depth for negative bathy 
     23   !   
    2124   USE diawri          ! Standard run outputs       (dia_wri_state routine) 
    22    ! 
    2325   USE in_out_manager  ! I/O manager 
    2426   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
    2527   USE lib_mpp         ! distributed memory computing 
    26    USE zdf_oce ,  ONLY : ln_zad_Aimp       ! ocean vertical physics variables 
    27    USE wet_dry,   ONLY : ll_wd, ssh_ref    ! reference depth for negative bathy 
    28  
     28   ! 
    2929   USE netcdf          ! NetCDF library 
    3030   IMPLICIT NONE 
     
    3333   PUBLIC stp_ctl           ! routine called by step.F90 
    3434 
    35    INTEGER  ::   idrun, idtime, idssh, idu, ids1, ids2, idt1, idt2, idc1, idw1, istatus 
    36    LOGICAL  ::   lsomeoce 
     35   INTEGER                ::   nrunid   ! netcdf file id 
     36   INTEGER, DIMENSION(8)  ::   nvarid   ! netcdf variable id 
    3737   !!---------------------------------------------------------------------- 
    3838   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    4242CONTAINS 
    4343 
    44    SUBROUTINE stp_ctl( kt, Kbb, Kmm, kindic ) 
     44   SUBROUTINE stp_ctl( kt, Kmm ) 
    4545      !!---------------------------------------------------------------------- 
    4646      !!                    ***  ROUTINE stp_ctl  *** 
     
    5050      !! ** Method  : - Save the time step in numstp 
    5151      !!              - Print it each 50 time steps 
    52       !!              - Stop the run IF problem encountered by setting indic=-3 
     52      !!              - Stop the run IF problem encountered by setting nstop > 0 
    5353      !!                Problems checked: |ssh| maximum larger than 10 m 
    5454      !!                                  |U|   maximum larger than 10 m/s  
     
    5757      !! ** Actions :   "time.step" file = last ocean time-step 
    5858      !!                "run.stat"  file = run statistics 
    59       !!                nstop indicator sheared among all local domain (lk_mpp=T) 
     59      !!                 nstop indicator sheared among all local domain 
    6060      !!---------------------------------------------------------------------- 
    6161      INTEGER, INTENT(in   ) ::   kt       ! ocean time-step index 
    62       INTEGER, INTENT(in   ) ::   Kbb, Kmm      ! ocean time level index 
    63       INTEGER, INTENT(inout) ::   kindic   ! error indicator 
    64       !! 
    65       INTEGER                ::   ji, jj, jk          ! dummy loop indices 
    66       INTEGER, DIMENSION(2)  ::   ih                  ! min/max loc indices 
    67       INTEGER, DIMENSION(3)  ::   iu, is1, is2        ! min/max loc indices 
    68       REAL(wp)               ::   zzz                 ! local real  
    69       REAL(wp), DIMENSION(9) ::   zmax 
    70       LOGICAL                ::   ll_wrtstp, ll_colruns, ll_wrtruns 
    71       CHARACTER(len=20) :: clname 
    72       !!---------------------------------------------------------------------- 
    73       ! 
    74       ll_wrtstp  = ( MOD( kt, sn_cfctl%ptimincr ) == 0 ) .OR. ( kt == nitend ) 
    75       ll_colruns = ll_wrtstp .AND. ( sn_cfctl%l_runstat ) 
    76       ll_wrtruns = ll_colruns .AND. lwm 
    77       IF( kt == nit000 .AND. lwp ) THEN 
    78          WRITE(numout,*) 
    79          WRITE(numout,*) 'stp_ctl : time-stepping control' 
    80          WRITE(numout,*) '~~~~~~~' 
    81          !                                ! open time.step file 
    82          IF( lwm ) CALL ctl_opn( numstp, 'time.step', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea ) 
    83          !                                ! open run.stat file(s) at start whatever 
    84          !                                ! the value of sn_cfctl%ptimincr 
    85          IF( lwm .AND. ( sn_cfctl%l_runstat ) ) THEN 
     62      INTEGER, INTENT(in   ) ::   Kmm      ! ocean time level index 
     63      !! 
     64      INTEGER                         ::   ji                                    ! dummy loop indices 
     65      INTEGER                         ::   idtime, istatus 
     66      INTEGER , DIMENSION(9)          ::   iareasum, iareamin, iareamax 
     67      INTEGER , DIMENSION(3,4)        ::   iloc                                  ! min/max loc indices 
     68      REAL(wp)                        ::   zzz                                   ! local real  
     69      REAL(wp), DIMENSION(9)          ::   zmax, zmaxlocal 
     70      LOGICAL                         ::   ll_wrtstp, ll_colruns, ll_wrtruns 
     71      LOGICAL, DIMENSION(jpi,jpj,jpk) ::   llmsk 
     72      CHARACTER(len=20)               ::   clname 
     73      !!---------------------------------------------------------------------- 
     74      IF( nstop > 0 .AND. ngrdstop > -1 )   RETURN   !   stpctl was already called by a child grid 
     75      ! 
     76      ll_wrtstp  = ( MOD( kt-nit000, sn_cfctl%ptimincr ) == 0 ) .OR. ( kt == nitend ) 
     77      ll_colruns = ll_wrtstp .AND. sn_cfctl%l_runstat .AND. jpnij > 1  
     78      ll_wrtruns = ( ll_colruns .OR. jpnij == 1 ) .AND. lwm 
     79      ! 
     80      IF( kt == nit000 ) THEN 
     81         ! 
     82         IF( lwp ) THEN 
     83            WRITE(numout,*) 
     84            WRITE(numout,*) 'stp_ctl : time-stepping control' 
     85            WRITE(numout,*) '~~~~~~~' 
     86         ENDIF 
     87         !                                ! open time.step    ascii file, done only by 1st subdomain 
     88         IF( lwm )   CALL ctl_opn( numstp, 'time.step', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea ) 
     89         ! 
     90         IF( ll_wrtruns ) THEN 
     91            !                             ! open run.stat     ascii file, done only by 1st subdomain 
    8692            CALL ctl_opn( numrun, 'run.stat', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, narea ) 
     93            !                             ! open run.stat.nc netcdf file, done only by 1st subdomain 
    8794            clname = 'run.stat.nc' 
    8895            IF( .NOT. Agrif_Root() )   clname = TRIM(Agrif_CFixed())//"_"//TRIM(clname) 
    89             istatus = NF90_CREATE( TRIM(clname), NF90_CLOBBER, idrun ) 
    90             istatus = NF90_DEF_DIM( idrun, 'time', NF90_UNLIMITED, idtime ) 
    91             istatus = NF90_DEF_VAR( idrun, 'abs_ssh_max', NF90_DOUBLE, (/ idtime /), idssh ) 
    92             istatus = NF90_DEF_VAR( idrun,   'abs_u_max', NF90_DOUBLE, (/ idtime /), idu  ) 
    93             istatus = NF90_DEF_VAR( idrun,       's_min', NF90_DOUBLE, (/ idtime /), ids1 ) 
    94             istatus = NF90_DEF_VAR( idrun,       's_max', NF90_DOUBLE, (/ idtime /), ids2 ) 
    95             istatus = NF90_DEF_VAR( idrun,       't_min', NF90_DOUBLE, (/ idtime /), idt1 ) 
    96             istatus = NF90_DEF_VAR( idrun,       't_max', NF90_DOUBLE, (/ idtime /), idt2 ) 
     96            istatus = NF90_CREATE( TRIM(clname), NF90_CLOBBER, nrunid ) 
     97            istatus = NF90_DEF_DIM( nrunid, 'time', NF90_UNLIMITED, idtime ) 
     98            istatus = NF90_DEF_VAR( nrunid, 'abs_ssh_max', NF90_DOUBLE, (/ idtime /), nvarid(1) ) 
     99            istatus = NF90_DEF_VAR( nrunid,   'abs_u_max', NF90_DOUBLE, (/ idtime /), nvarid(2) ) 
     100            istatus = NF90_DEF_VAR( nrunid,       's_min', NF90_DOUBLE, (/ idtime /), nvarid(3) ) 
     101            istatus = NF90_DEF_VAR( nrunid,       's_max', NF90_DOUBLE, (/ idtime /), nvarid(4) ) 
     102            istatus = NF90_DEF_VAR( nrunid,       't_min', NF90_DOUBLE, (/ idtime /), nvarid(5) ) 
     103            istatus = NF90_DEF_VAR( nrunid,       't_max', NF90_DOUBLE, (/ idtime /), nvarid(6) ) 
    97104            IF( ln_zad_Aimp ) THEN 
    98                istatus = NF90_DEF_VAR( idrun,   'abs_wi_max', NF90_DOUBLE, (/ idtime /), idw1 ) 
    99                istatus = NF90_DEF_VAR( idrun,       'Cf_max', NF90_DOUBLE, (/ idtime /), idc1 ) 
     105               istatus = NF90_DEF_VAR( nrunid,   'Cf_max', NF90_DOUBLE, (/ idtime /), nvarid(7) ) 
     106               istatus = NF90_DEF_VAR( nrunid,'abs_wi_max',NF90_DOUBLE, (/ idtime /), nvarid(8) ) 
    100107            ENDIF 
    101             istatus = NF90_ENDDEF(idrun) 
    102             zmax(8:9) = 0._wp    ! initialise to zero in case ln_zad_Aimp option is not in use 
    103          ENDIF 
    104       ENDIF 
    105       IF( kt == nit000 )   lsomeoce = COUNT( ssmask(:,:) == 1._wp ) > 0 
    106       ! 
    107       IF(lwm .AND. ll_wrtstp) THEN        !==  current time step  ==!   ("time.step" file) 
     108            istatus = NF90_ENDDEF(nrunid) 
     109         ENDIF 
     110         !     
     111      ENDIF 
     112      ! 
     113      !                                   !==              write current time step              ==! 
     114      !                                   !==  done only by 1st subdomain at writting timestep  ==! 
     115      IF( lwm .AND. ll_wrtstp ) THEN 
    108116         WRITE ( numstp, '(1x, i8)' )   kt 
    109117         REWIND( numstp ) 
    110118      ENDIF 
    111       ! 
    112       !                                   !==  test of extrema  ==! 
     119      !                                   !==            test of local extrema           ==! 
     120      !                                   !==  done by all processes at every time step  ==! 
     121      llmsk(:,:,1) = ssmask(:,:) == 1._wp 
    113122      IF( ll_wd ) THEN 
    114          zmax(1) = MAXVAL(  ABS( ssh(:,:,Kmm) + ssh_ref*tmask(:,:,1) )  )        ! ssh max  
     123         zmax(1) = MAXVAL( ABS( ssh(:,:,Kmm) + ssh_ref ), mask = llmsk(:,:,1) )   ! ssh max 
    115124      ELSE 
    116          zmax(1) = MAXVAL(  ABS( ssh(:,:,Kmm) )  )                               ! ssh max 
    117       ENDIF 
    118       zmax(2) = MAXVAL(  ABS( uu(:,:,:,Kmm) )  )                                  ! velocity max (zonal only) 
    119       zmax(3) = MAXVAL( -ts(:,:,:,jp_sal,Kmm) , mask = tmask(:,:,:) == 1._wp )   ! minus salinity max 
    120       zmax(4) = MAXVAL(  ts(:,:,:,jp_sal,Kmm) , mask = tmask(:,:,:) == 1._wp )   !       salinity max 
    121       zmax(5) = MAXVAL( -ts(:,:,:,jp_tem,Kmm) , mask = tmask(:,:,:) == 1._wp )   ! minus temperature max 
    122       zmax(6) = MAXVAL(  ts(:,:,:,jp_tem,Kmm) , mask = tmask(:,:,:) == 1._wp )   !       temperature max 
    123       zmax(7) = REAL( nstop , wp )                                            ! stop indicator 
    124       IF( ln_zad_Aimp ) THEN 
    125          zmax(8) = MAXVAL(  ABS( wi(:,:,:) ) , mask = wmask(:,:,:) == 1._wp ) ! implicit vertical vel. max 
    126          zmax(9) = MAXVAL(   Cu_adv(:,:,:)   , mask = tmask(:,:,:) == 1._wp ) ! partitioning coeff. max 
    127       ENDIF 
    128       ! 
     125         zmax(1) = MAXVAL( ABS( ssh(:,:,Kmm)           ), mask = llmsk(:,:,1) )   ! ssh max 
     126      ENDIF 
     127      llmsk(:,:,:) = umask(:,:,:) == 1._wp 
     128      zmax(2) = MAXVAL(  ABS( uu(:,:,:,Kmm) ), mask = llmsk )                     ! velocity max (zonal only) 
     129      llmsk(:,:,:) = tmask(:,:,:) == 1._wp 
     130      zmax(3) = MAXVAL( -ts(:,:,:,jp_sal,Kmm), mask = llmsk )                     ! minus salinity max 
     131      zmax(4) = MAXVAL(  ts(:,:,:,jp_sal,Kmm), mask = llmsk )                     !       salinity max 
     132      IF( ll_colruns .OR. jpnij == 1 ) THEN     ! following variables are used only in the netcdf file 
     133         zmax(5) = MAXVAL( -ts(:,:,:,jp_tem,Kmm), mask = llmsk )                  ! minus temperature max 
     134         zmax(6) = MAXVAL(  ts(:,:,:,jp_tem,Kmm), mask = llmsk )                  !       temperature max 
     135         IF( ln_zad_Aimp ) THEN 
     136            zmax(7) = MAXVAL(   Cu_adv(:,:,:)   , mask = llmsk )                  ! partitioning coeff. max 
     137            llmsk(:,:,:) = wmask(:,:,:) == 1._wp 
     138            zmax(8) = MAXVAL(  ABS( wi(:,:,:) ) , mask = llmsk )                  ! implicit vertical vel. max 
     139         ELSE 
     140            zmax(7:8) = 0._wp 
     141         ENDIF 
     142      ELSE 
     143         zmax(5:8) = 0._wp 
     144      ENDIF 
     145      zmax(9) = REAL( nstop, wp )                                              ! stop indicator 
     146      !                                   !==               get global extrema             ==! 
     147      !                                   !==  done by all processes if writting run.stat  ==! 
    129148      IF( ll_colruns ) THEN 
     149         zmaxlocal(:) = zmax(:) 
    130150         CALL mpp_max( "stpctl", zmax )          ! max over the global domain 
    131          nstop = NINT( zmax(7) )                 ! nstop indicator sheared among all local domains 
    132       ENDIF 
    133       !                                   !==  run statistics  ==!   ("run.stat" files) 
     151         nstop = NINT( zmax(9) )                 ! update nstop indicator (now sheared among all local domains) 
     152      ENDIF 
     153      !                                   !==              write "run.stat" files              ==! 
     154      !                                   !==  done only by 1st subdomain at writting timestep  ==! 
    134155      IF( ll_wrtruns ) THEN 
    135156         WRITE(numrun,9500) kt, zmax(1), zmax(2), -zmax(3), zmax(4) 
    136          istatus = NF90_PUT_VAR( idrun, idssh, (/ zmax(1)/), (/kt/), (/1/) ) 
    137          istatus = NF90_PUT_VAR( idrun,   idu, (/ zmax(2)/), (/kt/), (/1/) ) 
    138          istatus = NF90_PUT_VAR( idrun,  ids1, (/-zmax(3)/), (/kt/), (/1/) ) 
    139          istatus = NF90_PUT_VAR( idrun,  ids2, (/ zmax(4)/), (/kt/), (/1/) ) 
    140          istatus = NF90_PUT_VAR( idrun,  idt1, (/-zmax(5)/), (/kt/), (/1/) ) 
    141          istatus = NF90_PUT_VAR( idrun,  idt2, (/ zmax(6)/), (/kt/), (/1/) ) 
     157         istatus = NF90_PUT_VAR( nrunid, nvarid(1), (/ zmax(1)/), (/kt/), (/1/) ) 
     158         istatus = NF90_PUT_VAR( nrunid, nvarid(2), (/ zmax(2)/), (/kt/), (/1/) ) 
     159         istatus = NF90_PUT_VAR( nrunid, nvarid(3), (/-zmax(3)/), (/kt/), (/1/) ) 
     160         istatus = NF90_PUT_VAR( nrunid, nvarid(4), (/ zmax(4)/), (/kt/), (/1/) ) 
     161         istatus = NF90_PUT_VAR( nrunid, nvarid(5), (/-zmax(5)/), (/kt/), (/1/) ) 
     162         istatus = NF90_PUT_VAR( nrunid, nvarid(6), (/ zmax(6)/), (/kt/), (/1/) ) 
    142163         IF( ln_zad_Aimp ) THEN 
    143             istatus = NF90_PUT_VAR( idrun,  idw1, (/ zmax(8)/), (/kt/), (/1/) ) 
    144             istatus = NF90_PUT_VAR( idrun,  idc1, (/ zmax(9)/), (/kt/), (/1/) ) 
    145          ENDIF 
    146          IF( MOD( kt , 100 ) == 0 ) istatus = NF90_SYNC(idrun) 
    147          IF( kt == nitend         ) istatus = NF90_CLOSE(idrun) 
     164            istatus = NF90_PUT_VAR( nrunid, nvarid(7), (/ zmax(7)/), (/kt/), (/1/) ) 
     165            istatus = NF90_PUT_VAR( nrunid, nvarid(8), (/ zmax(8)/), (/kt/), (/1/) ) 
     166         ENDIF 
     167         IF( kt == nitend )   istatus = NF90_CLOSE(nrunid) 
    148168      END IF 
    149       !                                   !==  error handling  ==! 
    150       IF( ( sn_cfctl%l_glochk .OR. lsomeoce ) .AND. (   &  ! domain contains some ocean points, check for sensible ranges 
    151          &  zmax(1) >   20._wp .OR.   &                    ! too large sea surface height ( > 20 m ) 
    152          &  zmax(2) >   10._wp .OR.   &                    ! too large velocity ( > 10 m/s) 
    153 !!$         &  zmax(3) >=   0._wp .OR.   &                    ! negative or zero sea surface salinity 
    154 !!$         &  zmax(4) >= 100._wp .OR.   &                    ! too large sea surface salinity ( > 100 ) 
    155 !!$         &  zmax(4) <    0._wp .OR.   &                    ! too large sea surface salinity (keep this line for sea-ice) 
    156          &  ISNAN( zmax(1) + zmax(2) + zmax(3) ) ) ) THEN   ! NaN encounter in the tests 
    157          IF( lk_mpp .AND. sn_cfctl%l_glochk ) THEN 
    158             ! have use mpp_max (because sn_cfctl%l_glochk=.T. and distributed) 
    159             CALL mpp_maxloc( 'stpctl', ABS(ssh(:,:,Kmm))        , ssmask(:,:)  , zzz, ih  ) 
    160             CALL mpp_maxloc( 'stpctl', ABS(uu(:,:,:,Kmm))          , umask (:,:,:), zzz, iu  ) 
    161             CALL mpp_minloc( 'stpctl', ts(:,:,:,jp_sal,Kmm), tmask (:,:,:), zzz, is1 ) 
    162             CALL mpp_maxloc( 'stpctl', ts(:,:,:,jp_sal,Kmm), tmask (:,:,:), zzz, is2 ) 
     169      !                                   !==               error handling               ==! 
     170      !                                   !==  done by all processes at every time step  ==! 
     171      ! 
     172      IF(   zmax(1) >   20._wp .OR.   &                   ! too large sea surface height ( > 20 m ) 
     173         &  zmax(2) >   10._wp .OR.   &                   ! too large velocity ( > 10 m/s) 
     174!!$         &  zmax(3) >=   0._wp .OR.   &                   ! negative or zero sea surface salinity 
     175!!$         &  zmax(4) >= 100._wp .OR.   &                   ! too large sea surface salinity ( > 100 ) 
     176!!$         &  zmax(4) <    0._wp .OR.   &                   ! too large sea surface salinity (keep this line for sea-ice) 
     177         &  ISNAN( zmax(1) + zmax(2) + zmax(3) ) .OR.   &               ! NaN encounter in the tests 
     178         &  ABS(   zmax(1) + zmax(2) + zmax(3) ) > HUGE(1._wp) ) THEN   ! Infinity encounter in the tests 
     179         ! 
     180         iloc(:,:) = 0 
     181         IF( ll_colruns ) THEN   ! zmax is global, so it is the same on all subdomains -> no dead lock with mpp_maxloc 
     182            ! first: close the netcdf file, so we can read it 
     183            IF( lwm .AND. kt /= nitend )   istatus = NF90_CLOSE(nrunid) 
     184            ! get global loc on the min/max 
     185            CALL mpp_maxloc( 'stpctl', ABS(ssh(:,:,         Kmm)), ssmask(:,:  ), zzz, iloc(1:2,1) )   ! mpp_maxloc ok if mask = F  
     186            CALL mpp_maxloc( 'stpctl', ABS( uu(:,:,:,       Kmm)),  umask(:,:,:), zzz, iloc(1:3,2) ) 
     187            CALL mpp_minloc( 'stpctl',      ts(:,:,:,jp_sal,Kmm) ,  tmask(:,:,:), zzz, iloc(1:3,3) ) 
     188            CALL mpp_maxloc( 'stpctl',      ts(:,:,:,jp_sal,Kmm) ,  tmask(:,:,:), zzz, iloc(1:3,4) ) 
     189            ! find which subdomain has the max. 
     190            iareamin(:) = jpnij+1   ;   iareamax(:) = 0   ;   iareasum(:) = 0 
     191            DO ji = 1, 9 
     192               IF( zmaxlocal(ji) == zmax(ji) ) THEN 
     193                  iareamin(ji) = narea   ;   iareamax(ji) = narea   ;   iareasum(ji) = 1 
     194               ENDIF 
     195            END DO 
     196            CALL mpp_min( "stpctl", iareamin )         ! min over the global domain 
     197            CALL mpp_max( "stpctl", iareamax )         ! max over the global domain 
     198            CALL mpp_sum( "stpctl", iareasum )         ! sum over the global domain 
     199         ELSE                    ! find local min and max locations: 
     200            ! if we are here, this means that the subdomain contains some oce points -> no need to test the mask used in maxloc 
     201            iloc(1:2,1) = MAXLOC( ABS( ssh(:,:,         Kmm)), mask = ssmask(:,:  ) == 1._wp ) + (/ nimpp - 1, njmpp - 1    /) 
     202            iloc(1:3,2) = MAXLOC( ABS(  uu(:,:,:,       Kmm)), mask =  umask(:,:,:) == 1._wp ) + (/ nimpp - 1, njmpp - 1, 0 /) 
     203            iloc(1:3,3) = MINLOC(       ts(:,:,:,jp_sal,Kmm) , mask =  tmask(:,:,:) == 1._wp ) + (/ nimpp - 1, njmpp - 1, 0 /) 
     204            iloc(1:3,4) = MAXLOC(       ts(:,:,:,jp_sal,Kmm) , mask =  tmask(:,:,:) == 1._wp ) + (/ nimpp - 1, njmpp - 1, 0 /) 
     205            iareamin(:) = narea   ;   iareamax(:) = narea   ;   iareasum(:) = 1         ! this is local information 
     206         ENDIF 
     207         ! 
     208         WRITE(ctmp1,*) ' stp_ctl: |ssh| > 20 m  or  |U| > 10 m/s  or  S <= 0  or  S >= 100  or  NaN encounter in the tests' 
     209         CALL wrt_line( ctmp2, kt, '|ssh| max',  zmax(1), iloc(:,1), iareasum(1), iareamin(1), iareamax(1) ) 
     210         CALL wrt_line( ctmp3, kt, '|U|   max',  zmax(2), iloc(:,2), iareasum(2), iareamin(2), iareamax(2) ) 
     211         CALL wrt_line( ctmp4, kt, 'Sal   min', -zmax(3), iloc(:,3), iareasum(3), iareamin(3), iareamax(3) ) 
     212         CALL wrt_line( ctmp5, kt, 'Sal   max',  zmax(4), iloc(:,4), iareasum(4), iareamin(4), iareamax(4) ) 
     213         IF( Agrif_Root() ) THEN 
     214            WRITE(ctmp6,*) '      ===> output of last computed fields in output.abort* files' 
    163215         ELSE 
    164             ! find local min and max locations 
    165             ih(:)  = MAXLOC( ABS( ssh(:,:,Kmm)   )                              ) + (/ nimpp - 1, njmpp - 1    /) 
    166             iu(:)  = MAXLOC( ABS( uu  (:,:,:,Kmm) )                              ) + (/ nimpp - 1, njmpp - 1, 0 /) 
    167             is1(:) = MINLOC( ts(:,:,:,jp_sal,Kmm), mask = tmask(:,:,:) == 1._wp ) + (/ nimpp - 1, njmpp - 1, 0 /) 
    168             is2(:) = MAXLOC( ts(:,:,:,jp_sal,Kmm), mask = tmask(:,:,:) == 1._wp ) + (/ nimpp - 1, njmpp - 1, 0 /) 
    169          ENDIF 
    170           
    171          WRITE(ctmp1,*) ' stp_ctl: |ssh| > 20 m  or  |U| > 10 m/s  or  S <= 0  or  S >= 100  or  NaN encounter in the tests' 
    172          WRITE(ctmp2,9100) kt,   zmax(1), ih(1) , ih(2) 
    173          WRITE(ctmp3,9200) kt,   zmax(2), iu(1) , iu(2) , iu(3) 
    174          WRITE(ctmp4,9300) kt, - zmax(3), is1(1), is1(2), is1(3) 
    175          WRITE(ctmp5,9400) kt,   zmax(4), is2(1), is2(2), is2(3) 
    176          WRITE(ctmp6,*) '      ===> output of last computed fields in output.abort.nc file' 
    177           
     216            WRITE(ctmp6,*) '      ===> output of last computed fields in '//TRIM(Agrif_CFixed())//'_output.abort* files' 
     217         ENDIF 
     218         ! 
    178219         CALL dia_wri_state( Kmm, 'output.abort' )     ! create an output.abort file 
    179           
    180          IF( .NOT. sn_cfctl%l_glochk ) THEN 
    181             WRITE(ctmp8,*) 'E R R O R message from sub-domain: ', narea 
    182             CALL ctl_stop( 'STOP', ctmp1, ' ', ctmp8, ' ', ctmp2, ctmp3, ctmp4, ctmp5, ctmp6 ) 
    183          ELSE 
    184             CALL ctl_stop( ctmp1, ' ', ctmp2, ctmp3, ctmp4, ctmp5, ' ', ctmp6, ' ' ) 
    185          ENDIF 
    186  
    187          kindic = -3 
    188          ! 
    189       ENDIF 
    190       ! 
    191 9100  FORMAT (' kt=',i8,'   |ssh| max: ',1pg11.4,', at  i j  : ',2i5) 
    192 9200  FORMAT (' kt=',i8,'   |U|   max: ',1pg11.4,', at  i j k: ',3i5) 
    193 9300  FORMAT (' kt=',i8,'   S     min: ',1pg11.4,', at  i j k: ',3i5) 
    194 9400  FORMAT (' kt=',i8,'   S     max: ',1pg11.4,', at  i j k: ',3i5) 
     220         ! 
     221         IF( ll_colruns .or. jpnij == 1 ) THEN   ! all processes synchronized -> use lwp to print in opened ocean.output files 
     222            IF(lwp) THEN   ;   CALL ctl_stop( ctmp1, ' ', ctmp2, ctmp3, ctmp4, ctmp5, ' ', ctmp6 ) 
     223            ELSE           ;   nstop = MAX(1, nstop)   ! make sure nstop > 0 (automatically done when calling ctl_stop) 
     224            ENDIF 
     225         ELSE                                    ! only mpi subdomains with errors are here -> STOP now 
     226            CALL ctl_stop( 'STOP', ctmp1, ' ', ctmp2, ctmp3, ctmp4, ctmp5, ' ', ctmp6 ) 
     227         ENDIF 
     228         ! 
     229      ENDIF 
     230      ! 
     231      IF( nstop > 0 ) THEN                                                  ! an error was detected and we did not abort yet... 
     232         ngrdstop = Agrif_Fixed()                                           ! store which grid got this error 
     233         IF( .NOT. ll_colruns .AND. jpnij > 1 )   CALL ctl_stop( 'STOP' )   ! we must abort here to avoid MPI deadlock 
     234      ENDIF 
     235      ! 
    1952369500  FORMAT(' it :', i8, '    |ssh|_max: ', D23.16, ' |U|_max: ', D23.16,' S_min: ', D23.16,' S_max: ', D23.16) 
    196237      ! 
    197238   END SUBROUTINE stp_ctl 
     239 
     240 
     241   SUBROUTINE wrt_line( cdline, kt, cdprefix, pval, kloc, ksum, kmin, kmax ) 
     242      !!---------------------------------------------------------------------- 
     243      !!                     ***  ROUTINE wrt_line  *** 
     244      !! 
     245      !! ** Purpose :   write information line 
     246      !! 
     247      !!---------------------------------------------------------------------- 
     248      CHARACTER(len=*),      INTENT(  out) ::   cdline 
     249      CHARACTER(len=*),      INTENT(in   ) ::   cdprefix 
     250      REAL(wp),              INTENT(in   ) ::   pval 
     251      INTEGER, DIMENSION(3), INTENT(in   ) ::   kloc 
     252      INTEGER,               INTENT(in   ) ::   kt, ksum, kmin, kmax 
     253      ! 
     254      CHARACTER(len=80) ::   clsuff 
     255      CHARACTER(len=9 ) ::   clkt, clsum, clmin, clmax 
     256      CHARACTER(len=9 ) ::   cli, clj, clk 
     257      CHARACTER(len=1 ) ::   clfmt 
     258      CHARACTER(len=4 ) ::   cl4   ! needed to be able to compile with Agrif, I don't know why 
     259      INTEGER           ::   ifmtk 
     260      !!---------------------------------------------------------------------- 
     261      WRITE(clkt , '(i9)') kt 
     262       
     263      WRITE(clfmt, '(i1)') INT(LOG10(REAL(jpnij  ,wp))) + 1     ! how many digits to we need to write ? (we decide max = 9) 
     264      !!! WRITE(clsum, '(i'//clfmt//')') ksum                   ! this is creating a compilation error with AGRIF 
     265      cl4 = '(i'//clfmt//')'   ;   WRITE(clsum, cl4) ksum 
     266      WRITE(clfmt, '(i1)') INT(LOG10(REAL(MAX(1,jpnij-1),wp))) + 1    ! how many digits to we need to write ? (we decide max = 9) 
     267      cl4 = '(i'//clfmt//')'   ;   WRITE(clmin, cl4) kmin-1 
     268                                   WRITE(clmax, cl4) kmax-1 
     269      ! 
     270      WRITE(clfmt, '(i1)') INT(LOG10(REAL(jpiglo,wp))) + 1      ! how many digits to we need to write jpiglo? (we decide max = 9) 
     271      cl4 = '(i'//clfmt//')'   ;   WRITE(cli, cl4) kloc(1)      ! this is ok with AGRIF 
     272      WRITE(clfmt, '(i1)') INT(LOG10(REAL(jpjglo,wp))) + 1      ! how many digits to we need to write jpjglo? (we decide max = 9) 
     273      cl4 = '(i'//clfmt//')'   ;   WRITE(clj, cl4) kloc(2)      ! this is ok with AGRIF 
     274      ! 
     275      IF( ksum == 1 ) THEN   ;   WRITE(clsuff,9100) TRIM(clmin) 
     276      ELSE                   ;   WRITE(clsuff,9200) TRIM(clsum), TRIM(clmin), TRIM(clmax) 
     277      ENDIF 
     278      IF(kloc(3) == 0) THEN 
     279         ifmtk = INT(LOG10(REAL(jpk,wp))) + 1                   ! how many digits to we need to write jpk? (we decide max = 9) 
     280         clk = REPEAT(' ', ifmtk)                               ! create the equivalent in blank string 
     281         WRITE(cdline,9300) TRIM(ADJUSTL(clkt)), TRIM(ADJUSTL(cdprefix)), pval, TRIM(cli), TRIM(clj), clk(1:ifmtk), TRIM(clsuff) 
     282      ELSE 
     283         WRITE(clfmt, '(i1)') INT(LOG10(REAL(jpk,wp))) + 1      ! how many digits to we need to write jpk? (we decide max = 9) 
     284         !!! WRITE(clk, '(i'//clfmt//')') kloc(3)               ! this is creating a compilation error with AGRIF 
     285         cl4 = '(i'//clfmt//')'   ;   WRITE(clk, cl4) kloc(3)   ! this is ok with AGRIF 
     286         WRITE(cdline,9400) TRIM(ADJUSTL(clkt)), TRIM(ADJUSTL(cdprefix)), pval, TRIM(cli), TRIM(clj),    TRIM(clk), TRIM(clsuff) 
     287      ENDIF 
     288      ! 
     2899100  FORMAT('MPI rank ', a) 
     2909200  FORMAT('found in ', a, ' MPI tasks, spread out among ranks ', a, ' to ', a) 
     2919300  FORMAT('kt ', a, ' ', a, ' ', 1pg11.4, ' at i j   ', a, ' ', a, ' ', a, ' ', a) 
     2929400  FORMAT('kt ', a, ' ', a, ' ', 1pg11.4, ' at i j k ', a, ' ', a, ' ', a, ' ', a) 
     293      ! 
     294   END SUBROUTINE wrt_line 
     295 
    198296 
    199297   !!====================================================================== 
  • NEMO/trunk/tests/VORTEX/EXPREF/1_namelist_cfg

    r12489 r13216  
    9898&namagrif      !  AGRIF zoom                                            ("key_agrif") 
    9999!----------------------------------------------------------------------- 
    100    ln_spc_dyn    = .true.  !  use 0 as special value for dynamics 
    101    rn_sponge_tra =  800.   !  coefficient for tracer   sponge layer [m2/s] 
    102    rn_sponge_dyn =  800.   !  coefficient for dynamics sponge layer [m2/s] 
    103    ln_chk_bathy  = .FALSE. ! 
     100   rn_sponge_tra =  0.00768   !  coefficient for tracer   sponge layer [] 
     101   rn_sponge_dyn =  0.00768   !  coefficient for dynamics sponge layer [] 
    104102/ 
    105103!!====================================================================== 
Note: See TracChangeset for help on using the changeset viewer.