Ignore:
Timestamp:
2020-02-12T15:39:06+01:00 (17 months ago)
Author:
acc
Message:

The big one. Merging all 2019 developments from the option 1 branch back onto the trunk.

This changeset reproduces 2019/dev_r11943_MERGE_2019 on the trunk using a 2-URL merge
onto a working copy of the trunk. I.e.:

svn merge —ignore-ancestry \

svn+ssh://acc@forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/NEMO/trunk \
svn+ssh://acc@forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/NEMO/branches/2019/dev_r11943_MERGE_2019 ./

The —ignore-ancestry flag avoids problems that may otherwise arise from the fact that
the merge history been trunk and branch may have been applied in a different order but
care has been taken before this step to ensure that all applicable fixes and updates
are present in the merge branch.

The trunk state just before this step has been branched to releases/release-4.0-HEAD
and that branch has been immediately tagged as releases/release-4.0.2. Any fixes
or additions in response to tickets on 4.0, 4.0.1 or 4.0.2 should be done on
releases/release-4.0-HEAD. From now on future 'point' releases (e.g. 4.0.2) will
remain unchanged with periodic releases as needs demand. Note release-4.0-HEAD is a
transitional naming convention. Future full releases, say 4.2, will have a release-4.2
branch which fulfills this role and the first point release (e.g. 4.2.0) will be made
immediately following the release branch creation.

2020 developments can be started from any trunk revision later than this one.

Location:
NEMO/trunk
Files:
5 edited

Legend:

Unmodified
Added
Removed
  • NEMO/trunk

    • Property svn:externals
      •  

        old new  
        33^/utils/build/mk@HEAD         mk 
        44^/utils/tools@HEAD            tools 
        5 ^/vendors/AGRIF/dev@HEAD      ext/AGRIF 
         5^/vendors/AGRIF/dev_r11615_ENHANCE-04_namelists_as_internalfiles_agrif@HEAD      ext/AGRIF 
        66^/vendors/FCM@HEAD            ext/FCM 
        77^/vendors/IOIPSL@HEAD         ext/IOIPSL 
  • NEMO/trunk/tests/VORTEX/MY_SRC/domvvl.F90

    r11536 r12377  
    88   !!            3.3  !  2011-10  (M. Leclair) totally rewrote domvvl: vvl option includes z_star and z_tilde coordinates 
    99   !!            3.6  !  2014-11  (P. Mathiot) add ice shelf capability 
     10   !!            4.1  !  2019-08  (A. Coward, D. Storkey) rename dom_vvl_sf_swp -> dom_vvl_sf_update for new timestepping 
    1011   !!---------------------------------------------------------------------- 
    1112 
     
    1314   !!   dom_vvl_init     : define initial vertical scale factors, depths and column thickness 
    1415   !!   dom_vvl_sf_nxt   : Compute next vertical scale factors 
    15    !!   dom_vvl_sf_swp   : Swap vertical scale factors and update the vertical grid 
     16   !!   dom_vvl_sf_update   : Swap vertical scale factors and update the vertical grid 
    1617   !!   dom_vvl_interpol : Interpolate vertical scale factors from one grid point to another 
    1718   !!   dom_vvl_rst      : read/write restart file 
     
    3637 
    3738   PUBLIC  dom_vvl_init       ! called by domain.F90 
     39   PUBLIC  dom_vvl_zgr        ! called by isfcpl.F90 
    3840   PUBLIC  dom_vvl_sf_nxt     ! called by step.F90 
    39    PUBLIC  dom_vvl_sf_swp     ! called by step.F90 
     41   PUBLIC  dom_vvl_sf_update  ! called by step.F90 
    4042   PUBLIC  dom_vvl_interpol   ! called by dynnxt.F90 
    4143 
     
    6163   REAL(wp)        , ALLOCATABLE, SAVE, DIMENSION(:,:)   :: frq_rst_hdv                 ! retoring period for low freq. divergence 
    6264 
    63    !! * Substitutions 
    64 #  include "vectopt_loop_substitute.h90" 
    6565   !!---------------------------------------------------------------------- 
    6666   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    9393 
    9494 
    95    SUBROUTINE dom_vvl_init 
     95   SUBROUTINE dom_vvl_init( Kbb, Kmm, Kaa ) 
    9696      !!---------------------------------------------------------------------- 
    9797      !!                ***  ROUTINE dom_vvl_init  *** 
     
    102102      !! ** Method  :  - use restart file and/or initialize 
    103103      !!               - interpolate scale factors 
     104      !! 
     105      !! ** Action  : - e3t_(n/b) and tilde_e3t_(n/b) 
     106      !!              - Regrid: e3[u/v](:,:,:,Kmm) 
     107      !!                        e3[u/v](:,:,:,Kmm)        
     108      !!                        e3w(:,:,:,Kmm)            
     109      !!                        e3[u/v]w_b 
     110      !!                        e3[u/v]w_n       
     111      !!                        gdept(:,:,:,Kmm), gdepw(:,:,:,Kmm) and gde3w 
     112      !!              - h(t/u/v)_0 
     113      !!              - frq_rst_e3t and frq_rst_hdv 
     114      !! 
     115      !! Reference  : Leclair, M., and G. Madec, 2011, Ocean Modelling. 
     116      !!---------------------------------------------------------------------- 
     117      INTEGER, INTENT(in) :: Kbb, Kmm, Kaa 
     118      ! 
     119      IF(lwp) WRITE(numout,*) 
     120      IF(lwp) WRITE(numout,*) 'dom_vvl_init : Variable volume activated' 
     121      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~' 
     122      ! 
     123      CALL dom_vvl_ctl     ! choose vertical coordinate (z_star, z_tilde or layer) 
     124      ! 
     125      !                    ! Allocate module arrays 
     126      IF( dom_vvl_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'dom_vvl_init : unable to allocate arrays' ) 
     127      ! 
     128      !                    ! Read or initialize e3t_(b/n), tilde_e3t_(b/n) and hdiv_lf 
     129      CALL dom_vvl_rst( nit000, Kbb, Kmm, 'READ' ) 
     130      e3t(:,:,jpk,Kaa) = e3t_0(:,:,jpk)  ! last level always inside the sea floor set one for all 
     131      ! 
     132      CALL dom_vvl_zgr(Kbb, Kmm, Kaa) ! interpolation scale factor, depth and water column 
     133      ! 
     134   END SUBROUTINE dom_vvl_init 
     135   ! 
     136   SUBROUTINE dom_vvl_zgr(Kbb, Kmm, Kaa) 
     137      !!---------------------------------------------------------------------- 
     138      !!                ***  ROUTINE dom_vvl_init  *** 
     139      !!                    
     140      !! ** Purpose :  Interpolation of all scale factors,  
     141      !!               depths and water column heights 
     142      !! 
     143      !! ** Method  :  - interpolate scale factors 
    104144      !! 
    105145      !! ** Action  : - e3t_(n/b) and tilde_e3t_(n/b) 
     
    115155      !! Reference  : Leclair, M., and G. Madec, 2011, Ocean Modelling. 
    116156      !!---------------------------------------------------------------------- 
     157      INTEGER, INTENT(in) :: Kbb, Kmm, Kaa 
     158      !!---------------------------------------------------------------------- 
    117159      INTEGER ::   ji, jj, jk 
    118160      INTEGER ::   ii0, ii1, ij0, ij1 
     
    120162      !!---------------------------------------------------------------------- 
    121163      ! 
    122       IF(lwp) WRITE(numout,*) 
    123       IF(lwp) WRITE(numout,*) 'dom_vvl_init : Variable volume activated' 
    124       IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~' 
    125       ! 
    126       CALL dom_vvl_ctl     ! choose vertical coordinate (z_star, z_tilde or layer) 
    127       ! 
    128       !                    ! Allocate module arrays 
    129       IF( dom_vvl_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'dom_vvl_init : unable to allocate arrays' ) 
    130       ! 
    131       !                    ! Read or initialize e3t_(b/n), tilde_e3t_(b/n) and hdiv_lf 
    132       CALL dom_vvl_rst( nit000, 'READ' ) 
    133       e3t_a(:,:,jpk) = e3t_0(:,:,jpk)  ! last level always inside the sea floor set one for all 
    134       ! 
    135164      !                    !== Set of all other vertical scale factors  ==!  (now and before) 
    136165      !                                ! Horizontal interpolation of e3t 
    137       CALL dom_vvl_interpol( e3t_b(:,:,:), e3u_b(:,:,:), 'U' )    ! from T to U 
    138       CALL dom_vvl_interpol( e3t_n(:,:,:), e3u_n(:,:,:), 'U' ) 
    139       CALL dom_vvl_interpol( e3t_b(:,:,:), e3v_b(:,:,:), 'V' )    ! from T to V  
    140       CALL dom_vvl_interpol( e3t_n(:,:,:), e3v_n(:,:,:), 'V' ) 
    141       CALL dom_vvl_interpol( e3u_n(:,:,:), e3f_n(:,:,:), 'F' )    ! from U to F 
     166      CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3u(:,:,:,Kbb), 'U' )    ! from T to U 
     167      CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3u(:,:,:,Kmm), 'U' ) 
     168      CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3v(:,:,:,Kbb), 'V' )    ! from T to V  
     169      CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3v(:,:,:,Kmm), 'V' ) 
     170      CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3f(:,:,:), 'F' )    ! from U to F 
    142171      !                                ! Vertical interpolation of e3t,u,v  
    143       CALL dom_vvl_interpol( e3t_n(:,:,:), e3w_n (:,:,:), 'W'  )  ! from T to W 
    144       CALL dom_vvl_interpol( e3t_b(:,:,:), e3w_b (:,:,:), 'W'  ) 
    145       CALL dom_vvl_interpol( e3u_n(:,:,:), e3uw_n(:,:,:), 'UW' )  ! from U to UW 
    146       CALL dom_vvl_interpol( e3u_b(:,:,:), e3uw_b(:,:,:), 'UW' ) 
    147       CALL dom_vvl_interpol( e3v_n(:,:,:), e3vw_n(:,:,:), 'VW' )  ! from V to UW 
    148       CALL dom_vvl_interpol( e3v_b(:,:,:), e3vw_b(:,:,:), 'VW' ) 
     172      CALL dom_vvl_interpol( e3t(:,:,:,Kmm), e3w (:,:,:,Kmm), 'W'  )  ! from T to W 
     173      CALL dom_vvl_interpol( e3t(:,:,:,Kbb), e3w (:,:,:,Kbb), 'W'  ) 
     174      CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3uw(:,:,:,Kmm), 'UW' )  ! from U to UW 
     175      CALL dom_vvl_interpol( e3u(:,:,:,Kbb), e3uw(:,:,:,Kbb), 'UW' ) 
     176      CALL dom_vvl_interpol( e3v(:,:,:,Kmm), e3vw(:,:,:,Kmm), 'VW' )  ! from V to UW 
     177      CALL dom_vvl_interpol( e3v(:,:,:,Kbb), e3vw(:,:,:,Kbb), 'VW' ) 
    149178 
    150179      ! We need to define e3[tuv]_a for AGRIF initialisation (should not be a problem for the restartability...) 
    151       e3t_a(:,:,:) = e3t_n(:,:,:) 
    152       e3u_a(:,:,:) = e3u_n(:,:,:) 
    153       e3v_a(:,:,:) = e3v_n(:,:,:) 
     180      e3t(:,:,:,Kaa) = e3t(:,:,:,Kmm) 
     181      e3u(:,:,:,Kaa) = e3u(:,:,:,Kmm) 
     182      e3v(:,:,:,Kaa) = e3v(:,:,:,Kmm) 
    154183      ! 
    155184      !                    !==  depth of t and w-point  ==!   (set the isf depth as it is in the initial timestep) 
    156       gdept_n(:,:,1) = 0.5_wp * e3w_n(:,:,1)       ! reference to the ocean surface (used for MLD and light penetration) 
    157       gdepw_n(:,:,1) = 0.0_wp 
    158       gde3w_n(:,:,1) = gdept_n(:,:,1) - sshn(:,:)  ! reference to a common level z=0 for hpg 
    159       gdept_b(:,:,1) = 0.5_wp * e3w_b(:,:,1) 
    160       gdepw_b(:,:,1) = 0.0_wp 
     185      gdept(:,:,1,Kmm) = 0.5_wp * e3w(:,:,1,Kmm)       ! reference to the ocean surface (used for MLD and light penetration) 
     186      gdepw(:,:,1,Kmm) = 0.0_wp 
     187      gde3w(:,:,1) = gdept(:,:,1,Kmm) - ssh(:,:,Kmm)  ! reference to a common level z=0 for hpg 
     188      gdept(:,:,1,Kbb) = 0.5_wp * e3w(:,:,1,Kbb) 
     189      gdepw(:,:,1,Kbb) = 0.0_wp 
    161190      DO jk = 2, jpk                               ! vertical sum 
    162191         DO jj = 1,jpj 
     
    165194               !                             ! 1 everywhere from mbkt to mikt + 1 or 1 (if no isf) 
    166195               !                             ! 0.5 where jk = mikt      
    167 !!gm ???????   BUG ?  gdept_n as well as gde3w_n  does not include the thickness of ISF ?? 
     196!!gm ???????   BUG ?  gdept(:,:,:,Kmm) as well as gde3w  does not include the thickness of ISF ?? 
    168197               zcoef = ( tmask(ji,jj,jk) - wmask(ji,jj,jk) ) 
    169                gdepw_n(ji,jj,jk) = gdepw_n(ji,jj,jk-1) + e3t_n(ji,jj,jk-1) 
    170                gdept_n(ji,jj,jk) =      zcoef  * ( gdepw_n(ji,jj,jk  ) + 0.5 * e3w_n(ji,jj,jk))  & 
    171                   &                + (1-zcoef) * ( gdept_n(ji,jj,jk-1) +       e3w_n(ji,jj,jk))  
    172                gde3w_n(ji,jj,jk) = gdept_n(ji,jj,jk) - sshn(ji,jj) 
    173                gdepw_b(ji,jj,jk) = gdepw_b(ji,jj,jk-1) + e3t_b(ji,jj,jk-1) 
    174                gdept_b(ji,jj,jk) =      zcoef  * ( gdepw_b(ji,jj,jk  ) + 0.5 * e3w_b(ji,jj,jk))  & 
    175                   &                + (1-zcoef) * ( gdept_b(ji,jj,jk-1) +       e3w_b(ji,jj,jk))  
     198               gdepw(ji,jj,jk,Kmm) = gdepw(ji,jj,jk-1,Kmm) + e3t(ji,jj,jk-1,Kmm) 
     199               gdept(ji,jj,jk,Kmm) =      zcoef  * ( gdepw(ji,jj,jk  ,Kmm) + 0.5 * e3w(ji,jj,jk,Kmm))  & 
     200                  &                + (1-zcoef) * ( gdept(ji,jj,jk-1,Kmm) +       e3w(ji,jj,jk,Kmm))  
     201               gde3w(ji,jj,jk) = gdept(ji,jj,jk,Kmm) - ssh(ji,jj,Kmm) 
     202               gdepw(ji,jj,jk,Kbb) = gdepw(ji,jj,jk-1,Kbb) + e3t(ji,jj,jk-1,Kbb) 
     203               gdept(ji,jj,jk,Kbb) =      zcoef  * ( gdepw(ji,jj,jk  ,Kbb) + 0.5 * e3w(ji,jj,jk,Kbb))  & 
     204                  &                + (1-zcoef) * ( gdept(ji,jj,jk-1,Kbb) +       e3w(ji,jj,jk,Kbb))  
    176205            END DO 
    177206         END DO 
     
    179208      ! 
    180209      !                    !==  thickness of the water column  !!   (ocean portion only) 
    181       ht_n(:,:) = e3t_n(:,:,1) * tmask(:,:,1)   !!gm  BUG  :  this should be 1/2 * e3w(k=1) .... 
    182       hu_b(:,:) = e3u_b(:,:,1) * umask(:,:,1) 
    183       hu_n(:,:) = e3u_n(:,:,1) * umask(:,:,1) 
    184       hv_b(:,:) = e3v_b(:,:,1) * vmask(:,:,1) 
    185       hv_n(:,:) = e3v_n(:,:,1) * vmask(:,:,1) 
     210      ht(:,:) = e3t(:,:,1,Kmm) * tmask(:,:,1)   !!gm  BUG  :  this should be 1/2 * e3w(k=1) .... 
     211      hu(:,:,Kbb) = e3u(:,:,1,Kbb) * umask(:,:,1) 
     212      hu(:,:,Kmm) = e3u(:,:,1,Kmm) * umask(:,:,1) 
     213      hv(:,:,Kbb) = e3v(:,:,1,Kbb) * vmask(:,:,1) 
     214      hv(:,:,Kmm) = e3v(:,:,1,Kmm) * vmask(:,:,1) 
    186215      DO jk = 2, jpkm1 
    187          ht_n(:,:) = ht_n(:,:) + e3t_n(:,:,jk) * tmask(:,:,jk) 
    188          hu_b(:,:) = hu_b(:,:) + e3u_b(:,:,jk) * umask(:,:,jk) 
    189          hu_n(:,:) = hu_n(:,:) + e3u_n(:,:,jk) * umask(:,:,jk) 
    190          hv_b(:,:) = hv_b(:,:) + e3v_b(:,:,jk) * vmask(:,:,jk) 
    191          hv_n(:,:) = hv_n(:,:) + e3v_n(:,:,jk) * vmask(:,:,jk) 
     216         ht(:,:) = ht(:,:) + e3t(:,:,jk,Kmm) * tmask(:,:,jk) 
     217         hu(:,:,Kbb) = hu(:,:,Kbb) + e3u(:,:,jk,Kbb) * umask(:,:,jk) 
     218         hu(:,:,Kmm) = hu(:,:,Kmm) + e3u(:,:,jk,Kmm) * umask(:,:,jk) 
     219         hv(:,:,Kbb) = hv(:,:,Kbb) + e3v(:,:,jk,Kbb) * vmask(:,:,jk) 
     220         hv(:,:,Kmm) = hv(:,:,Kmm) + e3v(:,:,jk,Kmm) * vmask(:,:,jk) 
    192221      END DO 
    193222      ! 
    194223      !                    !==  inverse of water column thickness   ==!   (u- and v- points) 
    195       r1_hu_b(:,:) = ssumask(:,:) / ( hu_b(:,:) + 1._wp - ssumask(:,:) )    ! _i mask due to ISF 
    196       r1_hu_n(:,:) = ssumask(:,:) / ( hu_n(:,:) + 1._wp - ssumask(:,:) ) 
    197       r1_hv_b(:,:) = ssvmask(:,:) / ( hv_b(:,:) + 1._wp - ssvmask(:,:) ) 
    198       r1_hv_n(:,:) = ssvmask(:,:) / ( hv_n(:,:) + 1._wp - ssvmask(:,:) ) 
     224      r1_hu(:,:,Kbb) = ssumask(:,:) / ( hu(:,:,Kbb) + 1._wp - ssumask(:,:) )    ! _i mask due to ISF 
     225      r1_hu(:,:,Kmm) = ssumask(:,:) / ( hu(:,:,Kmm) + 1._wp - ssumask(:,:) ) 
     226      r1_hv(:,:,Kbb) = ssvmask(:,:) / ( hv(:,:,Kbb) + 1._wp - ssvmask(:,:) ) 
     227      r1_hv(:,:,Kmm) = ssvmask(:,:) / ( hv(:,:,Kmm) + 1._wp - ssvmask(:,:) ) 
    199228 
    200229      !                    !==   z_tilde coordinate case  ==!   (Restoring frequencies) 
     
    263292      ENDIF 
    264293      ! 
    265    END SUBROUTINE dom_vvl_init 
    266  
    267  
    268    SUBROUTINE dom_vvl_sf_nxt( kt, kcall )  
     294   END SUBROUTINE dom_vvl_zgr 
     295 
     296 
     297   SUBROUTINE dom_vvl_sf_nxt( kt, Kbb, Kmm, Kaa, kcall )  
    269298      !!---------------------------------------------------------------------- 
    270299      !!                ***  ROUTINE dom_vvl_sf_nxt  *** 
     
    288317      !! Reference  : Leclair, M., and Madec, G. 2011, Ocean Modelling. 
    289318      !!---------------------------------------------------------------------- 
    290       INTEGER, INTENT( in )           ::   kt      ! time step 
    291       INTEGER, INTENT( in ), OPTIONAL ::   kcall   ! optional argument indicating call sequence 
     319      INTEGER, INTENT( in )           ::   kt             ! time step 
     320      INTEGER, INTENT( in )           ::   Kbb, Kmm, Kaa  ! time step 
     321      INTEGER, INTENT( in ), OPTIONAL ::   kcall          ! optional argument indicating call sequence 
    292322      ! 
    293323      INTEGER                ::   ji, jj, jk            ! dummy loop indices 
     
    321351      !                                                ! --------------------------------------------- ! 
    322352      ! 
    323       z_scale(:,:) = ( ssha(:,:) - sshb(:,:) ) * ssmask(:,:) / ( ht_0(:,:) + sshn(:,:) + 1. - ssmask(:,:) ) 
     353      z_scale(:,:) = ( ssh(:,:,Kaa) - ssh(:,:,Kbb) ) * ssmask(:,:) / ( ht_0(:,:) + ssh(:,:,Kmm) + 1. - ssmask(:,:) ) 
    324354      DO jk = 1, jpkm1 
    325          ! formally this is the same as e3t_a = e3t_0*(1+ssha/ht_0) 
    326          e3t_a(:,:,jk) = e3t_b(:,:,jk) + e3t_n(:,:,jk) * z_scale(:,:) * tmask(:,:,jk) 
     355         ! formally this is the same as e3t(:,:,:,Kaa) = e3t_0*(1+ssha/ht_0) 
     356         e3t(:,:,jk,Kaa) = e3t(:,:,jk,Kbb) + e3t(:,:,jk,Kmm) * z_scale(:,:) * tmask(:,:,jk) 
    327357      END DO 
    328358      ! 
     
    337367         zht(:,:)   = 0._wp 
    338368         DO jk = 1, jpkm1 
    339             zhdiv(:,:) = zhdiv(:,:) + e3t_n(:,:,jk) * hdivn(:,:,jk) 
    340             zht  (:,:) = zht  (:,:) + e3t_n(:,:,jk) * tmask(:,:,jk) 
     369            zhdiv(:,:) = zhdiv(:,:) + e3t(:,:,jk,Kmm) * hdiv(:,:,jk) 
     370            zht  (:,:) = zht  (:,:) + e3t(:,:,jk,Kmm) * tmask(:,:,jk) 
    341371         END DO 
    342372         zhdiv(:,:) = zhdiv(:,:) / ( zht(:,:) + 1. - tmask_i(:,:) ) 
     
    348378               DO jk = 1, jpkm1 
    349379                  hdiv_lf(:,:,jk) = hdiv_lf(:,:,jk) - rdt * frq_rst_hdv(:,:)   & 
    350                      &          * ( hdiv_lf(:,:,jk) - e3t_n(:,:,jk) * ( hdivn(:,:,jk) - zhdiv(:,:) ) ) 
     380                     &          * ( hdiv_lf(:,:,jk) - e3t(:,:,jk,Kmm) * ( hdiv(:,:,jk) - zhdiv(:,:) ) ) 
    351381               END DO 
    352382            ENDIF 
     
    361391         IF( ln_vvl_ztilde ) THEN     ! z_tilde case 
    362392            DO jk = 1, jpkm1 
    363                tilde_e3t_a(:,:,jk) = tilde_e3t_a(:,:,jk) - ( e3t_n(:,:,jk) * ( hdivn(:,:,jk) - zhdiv(:,:) ) - hdiv_lf(:,:,jk) ) 
     393               tilde_e3t_a(:,:,jk) = tilde_e3t_a(:,:,jk) - ( e3t(:,:,jk,Kmm) * ( hdiv(:,:,jk) - zhdiv(:,:) ) - hdiv_lf(:,:,jk) ) 
    364394            END DO 
    365395         ELSE                         ! layer case 
    366396            DO jk = 1, jpkm1 
    367                tilde_e3t_a(:,:,jk) = tilde_e3t_a(:,:,jk) -   e3t_n(:,:,jk) * ( hdivn(:,:,jk) - zhdiv(:,:) ) * tmask(:,:,jk) 
     397               tilde_e3t_a(:,:,jk) = tilde_e3t_a(:,:,jk) -   e3t(:,:,jk,Kmm) * ( hdiv(:,:,jk) - zhdiv(:,:) ) * tmask(:,:,jk) 
    368398            END DO 
    369399         ENDIF 
     
    383413         DO jk = 1, jpkm1        ! a - first derivative: diffusive fluxes 
    384414            DO jj = 1, jpjm1 
    385                DO ji = 1, fs_jpim1   ! vector opt. 
     415               DO ji = 1, jpim1   ! vector opt. 
    386416                  un_td(ji,jj,jk) = rn_ahe3 * umask(ji,jj,jk) * e2_e1u(ji,jj)           & 
    387417                     &            * ( tilde_e3t_b(ji,jj,jk) - tilde_e3t_b(ji+1,jj  ,jk) ) 
     
    401431         DO jk = 1, jpkm1        ! c - second derivative: divergence of diffusive fluxes 
    402432            DO jj = 2, jpjm1 
    403                DO ji = fs_2, fs_jpim1   ! vector opt. 
     433               DO ji = 2, jpim1   ! vector opt. 
    404434                  tilde_e3t_a(ji,jj,jk) = tilde_e3t_a(ji,jj,jk) + (   un_td(ji-1,jj  ,jk) - un_td(ji,jj,jk)    & 
    405435                     &                                          +     vn_td(ji  ,jj-1,jk) - vn_td(ji,jj,jk)    & 
     
    476506            zht(:,:)  = zht(:,:) + tilde_e3t_a(:,:,jk) * tmask(:,:,jk) 
    477507         END DO 
    478          z_scale(:,:) =  - zht(:,:) / ( ht_0(:,:) + sshn(:,:) + 1. - ssmask(:,:) ) 
     508         z_scale(:,:) =  - zht(:,:) / ( ht_0(:,:) + ssh(:,:,Kmm) + 1. - ssmask(:,:) ) 
    479509         DO jk = 1, jpkm1 
    480             dtilde_e3t_a(:,:,jk) = dtilde_e3t_a(:,:,jk) + e3t_n(:,:,jk) * z_scale(:,:) * tmask(:,:,jk) 
     510            dtilde_e3t_a(:,:,jk) = dtilde_e3t_a(:,:,jk) + e3t(:,:,jk,Kmm) * z_scale(:,:) * tmask(:,:,jk) 
    481511         END DO 
    482512 
     
    486516      !                                           ! ---baroclinic part--------- ! 
    487517         DO jk = 1, jpkm1 
    488             e3t_a(:,:,jk) = e3t_a(:,:,jk) + dtilde_e3t_a(:,:,jk) * tmask(:,:,jk) 
     518            e3t(:,:,jk,Kaa) = e3t(:,:,jk,Kaa) + dtilde_e3t_a(:,:,jk) * tmask(:,:,jk) 
    489519         END DO 
    490520      ENDIF 
     
    501531         zht(:,:) = 0.0_wp 
    502532         DO jk = 1, jpkm1 
    503             zht(:,:) = zht(:,:) + e3t_n(:,:,jk) * tmask(:,:,jk) 
    504          END DO 
    505          z_tmax = MAXVAL( tmask(:,:,1) * tmask_i(:,:) * ABS( ht_0(:,:) + sshn(:,:) - zht(:,:) ) ) 
     533            zht(:,:) = zht(:,:) + e3t(:,:,jk,Kmm) * tmask(:,:,jk) 
     534         END DO 
     535         z_tmax = MAXVAL( tmask(:,:,1) * tmask_i(:,:) * ABS( ht_0(:,:) + ssh(:,:,Kmm) - zht(:,:) ) ) 
    506536         CALL mpp_max( 'domvvl', z_tmax )                                ! max over the global domain 
    507          IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(ht_0+sshn-SUM(e3t_n))) =', z_tmax 
     537         IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(ht_0+sshn-SUM(e3t(:,:,:,Kmm)))) =', z_tmax 
    508538         ! 
    509539         zht(:,:) = 0.0_wp 
    510540         DO jk = 1, jpkm1 
    511             zht(:,:) = zht(:,:) + e3t_a(:,:,jk) * tmask(:,:,jk) 
    512          END DO 
    513          z_tmax = MAXVAL( tmask(:,:,1) * tmask_i(:,:) * ABS( ht_0(:,:) + ssha(:,:) - zht(:,:) ) ) 
     541            zht(:,:) = zht(:,:) + e3t(:,:,jk,Kaa) * tmask(:,:,jk) 
     542         END DO 
     543         z_tmax = MAXVAL( tmask(:,:,1) * tmask_i(:,:) * ABS( ht_0(:,:) + ssh(:,:,Kaa) - zht(:,:) ) ) 
    514544         CALL mpp_max( 'domvvl', z_tmax )                                ! max over the global domain 
    515          IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(ht_0+ssha-SUM(e3t_a))) =', z_tmax 
     545         IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(ht_0+ssha-SUM(e3t(:,:,:,Kaa)))) =', z_tmax 
    516546         ! 
    517547         zht(:,:) = 0.0_wp 
    518548         DO jk = 1, jpkm1 
    519             zht(:,:) = zht(:,:) + e3t_b(:,:,jk) * tmask(:,:,jk) 
    520          END DO 
    521          z_tmax = MAXVAL( tmask(:,:,1) * tmask_i(:,:) * ABS( ht_0(:,:) + sshb(:,:) - zht(:,:) ) ) 
     549            zht(:,:) = zht(:,:) + e3t(:,:,jk,Kbb) * tmask(:,:,jk) 
     550         END DO 
     551         z_tmax = MAXVAL( tmask(:,:,1) * tmask_i(:,:) * ABS( ht_0(:,:) + ssh(:,:,Kbb) - zht(:,:) ) ) 
    522552         CALL mpp_max( 'domvvl', z_tmax )                                ! max over the global domain 
    523          IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(ht_0+sshb-SUM(e3t_b))) =', z_tmax 
    524          ! 
    525          z_tmax = MAXVAL( tmask(:,:,1) *  ABS( sshb(:,:) ) ) 
     553         IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(ht_0+sshb-SUM(e3t(:,:,:,Kbb)))) =', z_tmax 
     554         ! 
     555         z_tmax = MAXVAL( tmask(:,:,1) *  ABS( ssh(:,:,Kbb) ) ) 
    526556         CALL mpp_max( 'domvvl', z_tmax )                                ! max over the global domain 
    527          IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(sshb))) =', z_tmax 
    528          ! 
    529          z_tmax = MAXVAL( tmask(:,:,1) *  ABS( sshn(:,:) ) ) 
     557         IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(ssh(:,:,Kbb)))) =', z_tmax 
     558         ! 
     559         z_tmax = MAXVAL( tmask(:,:,1) *  ABS( ssh(:,:,Kmm) ) ) 
    530560         CALL mpp_max( 'domvvl', z_tmax )                                ! max over the global domain 
    531          IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(sshn))) =', z_tmax 
    532          ! 
    533          z_tmax = MAXVAL( tmask(:,:,1) *  ABS( ssha(:,:) ) ) 
     561         IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(ssh(:,:,Kmm)))) =', z_tmax 
     562         ! 
     563         z_tmax = MAXVAL( tmask(:,:,1) *  ABS( ssh(:,:,Kaa) ) ) 
    534564         CALL mpp_max( 'domvvl', z_tmax )                                ! max over the global domain 
    535          IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(ssha))) =', z_tmax 
     565         IF( lwp    ) WRITE(numout, *) kt,' MAXVAL(abs(ssh(:,:,Kaa)))) =', z_tmax 
    536566      END IF 
    537567 
     
    540570      ! *********************************** ! 
    541571 
    542       CALL dom_vvl_interpol( e3t_a(:,:,:), e3u_a(:,:,:), 'U' ) 
    543       CALL dom_vvl_interpol( e3t_a(:,:,:), e3v_a(:,:,:), 'V' ) 
     572      CALL dom_vvl_interpol( e3t(:,:,:,Kaa), e3u(:,:,:,Kaa), 'U' ) 
     573      CALL dom_vvl_interpol( e3t(:,:,:,Kaa), e3v(:,:,:,Kaa), 'V' ) 
    544574 
    545575      ! *********************************** ! 
     
    547577      ! *********************************** ! 
    548578 
    549       hu_a(:,:) = e3u_a(:,:,1) * umask(:,:,1) 
    550       hv_a(:,:) = e3v_a(:,:,1) * vmask(:,:,1) 
     579      hu(:,:,Kaa) = e3u(:,:,1,Kaa) * umask(:,:,1) 
     580      hv(:,:,Kaa) = e3v(:,:,1,Kaa) * vmask(:,:,1) 
    551581      DO jk = 2, jpkm1 
    552          hu_a(:,:) = hu_a(:,:) + e3u_a(:,:,jk) * umask(:,:,jk) 
    553          hv_a(:,:) = hv_a(:,:) + e3v_a(:,:,jk) * vmask(:,:,jk) 
     582         hu(:,:,Kaa) = hu(:,:,Kaa) + e3u(:,:,jk,Kaa) * umask(:,:,jk) 
     583         hv(:,:,Kaa) = hv(:,:,Kaa) + e3v(:,:,jk,Kaa) * vmask(:,:,jk) 
    554584      END DO 
    555585      !                                        ! Inverse of the local depth 
    556586!!gm BUG ?  don't understand the use of umask_i here ..... 
    557       r1_hu_a(:,:) = ssumask(:,:) / ( hu_a(:,:) + 1._wp - ssumask(:,:) ) 
    558       r1_hv_a(:,:) = ssvmask(:,:) / ( hv_a(:,:) + 1._wp - ssvmask(:,:) ) 
     587      r1_hu(:,:,Kaa) = ssumask(:,:) / ( hu(:,:,Kaa) + 1._wp - ssumask(:,:) ) 
     588      r1_hv(:,:,Kaa) = ssvmask(:,:) / ( hv(:,:,Kaa) + 1._wp - ssvmask(:,:) ) 
    559589      ! 
    560590      IF( ln_timing )   CALL timing_stop('dom_vvl_sf_nxt') 
     
    563593 
    564594 
    565    SUBROUTINE dom_vvl_sf_swp( kt ) 
    566       !!---------------------------------------------------------------------- 
    567       !!                ***  ROUTINE dom_vvl_sf_swp  *** 
     595   SUBROUTINE dom_vvl_sf_update( kt, Kbb, Kmm, Kaa ) 
     596      !!---------------------------------------------------------------------- 
     597      !!                ***  ROUTINE dom_vvl_sf_update  *** 
    568598      !!                    
    569       !! ** Purpose :  compute time filter and swap of scale factors  
     599      !! ** Purpose :  for z tilde case: compute time filter and swap of scale factors  
    570600      !!               compute all depths and related variables for next time step 
    571601      !!               write outputs and restart file 
    572602      !! 
    573       !! ** Method  :  - swap of e3t with trick for volume/tracer conservation 
     603      !! ** Method  :  - swap of e3t with trick for volume/tracer conservation (ONLY FOR Z TILDE CASE) 
    574604      !!               - reconstruct scale factor at other grid points (interpolate) 
    575605      !!               - recompute depths and water height fields 
    576606      !! 
    577       !! ** Action  :  - e3t_(b/n), tilde_e3t_(b/n) and e3(u/v)_n ready for next time step 
     607      !! ** Action  :  - tilde_e3t_(b/n) ready for next time step 
    578608      !!               - Recompute: 
    579609      !!                    e3(u/v)_b        
    580       !!                    e3w_n            
     610      !!                    e3w(:,:,:,Kmm)            
    581611      !!                    e3(u/v)w_b       
    582612      !!                    e3(u/v)w_n       
    583       !!                    gdept_n, gdepw_n  and gde3w_n 
     613      !!                    gdept(:,:,:,Kmm), gdepw(:,:,:,Kmm)  and gde3w 
    584614      !!                    h(u/v) and h(u/v)r 
    585615      !! 
     
    587617      !!              Leclair, M., and G. Madec, 2011, Ocean Modelling. 
    588618      !!---------------------------------------------------------------------- 
    589       INTEGER, INTENT( in ) ::   kt   ! time step 
     619      INTEGER, INTENT( in ) ::   kt              ! time step 
     620      INTEGER, INTENT( in ) ::   Kbb, Kmm, Kaa   ! time level indices 
    590621      ! 
    591622      INTEGER  ::   ji, jj, jk   ! dummy loop indices 
     
    595626      IF( ln_linssh )   RETURN      ! No calculation in linear free surface 
    596627      ! 
    597       IF( ln_timing )   CALL timing_start('dom_vvl_sf_swp') 
     628      IF( ln_timing )   CALL timing_start('dom_vvl_sf_update') 
    598629      ! 
    599630      IF( kt == nit000 )   THEN 
    600631         IF(lwp) WRITE(numout,*) 
    601          IF(lwp) WRITE(numout,*) 'dom_vvl_sf_swp : - time filter and swap of scale factors' 
    602          IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~   - interpolate scale factors and compute depths for next time step' 
     632         IF(lwp) WRITE(numout,*) 'dom_vvl_sf_update : - interpolate scale factors and compute depths for next time step' 
     633         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~~~~' 
    603634      ENDIF 
    604635      ! 
     
    615646         tilde_e3t_n(:,:,:) = tilde_e3t_a(:,:,:) 
    616647      ENDIF 
    617       gdept_b(:,:,:) = gdept_n(:,:,:) 
    618       gdepw_b(:,:,:) = gdepw_n(:,:,:) 
    619  
    620       e3t_n(:,:,:) = e3t_a(:,:,:) 
    621       e3u_n(:,:,:) = e3u_a(:,:,:) 
    622       e3v_n(:,:,:) = e3v_a(:,:,:) 
    623648 
    624649      ! Compute all missing vertical scale factor and depths 
     
    626651      ! Horizontal scale factor interpolations 
    627652      ! -------------------------------------- 
    628       ! - ML - e3u_b and e3v_b are allready computed in dynnxt 
    629       ! - JC - hu_b, hv_b, hur_b, hvr_b also 
     653      ! - ML - e3u(:,:,:,Kbb) and e3v(:,:,:,Kbb) are allready computed in dynnxt 
     654      ! - JC - hu(:,:,:,Kbb), hv(:,:,:,:,Kbb), hur_b, hvr_b also 
    630655       
    631       CALL dom_vvl_interpol( e3u_n(:,:,:), e3f_n(:,:,:), 'F'  ) 
     656      CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3f(:,:,:), 'F'  ) 
    632657       
    633658      ! Vertical scale factor interpolations 
    634       CALL dom_vvl_interpol( e3t_n(:,:,:),  e3w_n(:,:,:), 'W'  ) 
    635       CALL dom_vvl_interpol( e3u_n(:,:,:), e3uw_n(:,:,:), 'UW' ) 
    636       CALL dom_vvl_interpol( e3v_n(:,:,:), e3vw_n(:,:,:), 'VW' ) 
    637       CALL dom_vvl_interpol( e3t_b(:,:,:),  e3w_b(:,:,:), 'W'  ) 
    638       CALL dom_vvl_interpol( e3u_b(:,:,:), e3uw_b(:,:,:), 'UW' ) 
    639       CALL dom_vvl_interpol( e3v_b(:,:,:), e3vw_b(:,:,:), 'VW' ) 
     659      CALL dom_vvl_interpol( e3t(:,:,:,Kmm),  e3w(:,:,:,Kmm), 'W'  ) 
     660      CALL dom_vvl_interpol( e3u(:,:,:,Kmm), e3uw(:,:,:,Kmm), 'UW' ) 
     661      CALL dom_vvl_interpol( e3v(:,:,:,Kmm), e3vw(:,:,:,Kmm), 'VW' ) 
     662      CALL dom_vvl_interpol( e3t(:,:,:,Kbb),  e3w(:,:,:,Kbb), 'W'  ) 
     663      CALL dom_vvl_interpol( e3u(:,:,:,Kbb), e3uw(:,:,:,Kbb), 'UW' ) 
     664      CALL dom_vvl_interpol( e3v(:,:,:,Kbb), e3vw(:,:,:,Kbb), 'VW' ) 
    640665 
    641666      ! t- and w- points depth (set the isf depth as it is in the initial step) 
    642       gdept_n(:,:,1) = 0.5_wp * e3w_n(:,:,1) 
    643       gdepw_n(:,:,1) = 0.0_wp 
    644       gde3w_n(:,:,1) = gdept_n(:,:,1) - sshn(:,:) 
     667      gdept(:,:,1,Kmm) = 0.5_wp * e3w(:,:,1,Kmm) 
     668      gdepw(:,:,1,Kmm) = 0.0_wp 
     669      gde3w(:,:,1) = gdept(:,:,1,Kmm) - ssh(:,:,Kmm) 
    645670      DO jk = 2, jpk 
    646671         DO jj = 1,jpj 
     
    649674                                                                 ! 1 for jk = mikt 
    650675               zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk)) 
    651                gdepw_n(ji,jj,jk) = gdepw_n(ji,jj,jk-1) + e3t_n(ji,jj,jk-1) 
    652                gdept_n(ji,jj,jk) =    zcoef  * ( gdepw_n(ji,jj,jk  ) + 0.5 * e3w_n(ji,jj,jk) )  & 
    653                    &             + (1-zcoef) * ( gdept_n(ji,jj,jk-1) +       e3w_n(ji,jj,jk) )  
    654                gde3w_n(ji,jj,jk) = gdept_n(ji,jj,jk) - sshn(ji,jj) 
     676               gdepw(ji,jj,jk,Kmm) = gdepw(ji,jj,jk-1,Kmm) + e3t(ji,jj,jk-1,Kmm) 
     677               gdept(ji,jj,jk,Kmm) =    zcoef  * ( gdepw(ji,jj,jk  ,Kmm) + 0.5 * e3w(ji,jj,jk,Kmm) )  & 
     678                   &             + (1-zcoef) * ( gdept(ji,jj,jk-1,Kmm) +       e3w(ji,jj,jk,Kmm) )  
     679               gde3w(ji,jj,jk) = gdept(ji,jj,jk,Kmm) - ssh(ji,jj,Kmm) 
    655680            END DO 
    656681         END DO 
     
    659684      ! Local depth and Inverse of the local depth of the water 
    660685      ! ------------------------------------------------------- 
    661       hu_n(:,:) = hu_a(:,:)   ;   r1_hu_n(:,:) = r1_hu_a(:,:) 
    662       hv_n(:,:) = hv_a(:,:)   ;   r1_hv_n(:,:) = r1_hv_a(:,:) 
    663       ! 
    664       ht_n(:,:) = e3t_n(:,:,1) * tmask(:,:,1) 
     686      ! 
     687      ht(:,:) = e3t(:,:,1,Kmm) * tmask(:,:,1) 
    665688      DO jk = 2, jpkm1 
    666          ht_n(:,:) = ht_n(:,:) + e3t_n(:,:,jk) * tmask(:,:,jk) 
     689         ht(:,:) = ht(:,:) + e3t(:,:,jk,Kmm) * tmask(:,:,jk) 
    667690      END DO 
    668691 
    669692      ! write restart file 
    670693      ! ================== 
    671       IF( lrst_oce  )   CALL dom_vvl_rst( kt, 'WRITE' ) 
    672       ! 
    673       IF( ln_timing )   CALL timing_stop('dom_vvl_sf_swp') 
    674       ! 
    675    END SUBROUTINE dom_vvl_sf_swp 
     694      IF( lrst_oce  )   CALL dom_vvl_rst( kt, Kbb, Kmm, 'WRITE' ) 
     695      ! 
     696      IF( ln_timing )   CALL timing_stop('dom_vvl_sf_update') 
     697      ! 
     698   END SUBROUTINE dom_vvl_sf_update 
    676699 
    677700 
     
    706729         DO jk = 1, jpk 
    707730            DO jj = 1, jpjm1 
    708                DO ji = 1, fs_jpim1   ! vector opt. 
     731               DO ji = 1, jpim1   ! vector opt. 
    709732                  pe3_out(ji,jj,jk) = 0.5_wp * (  umask(ji,jj,jk) * (1.0_wp - zlnwd) + zlnwd ) * r1_e1e2u(ji,jj)   & 
    710733                     &                       * (   e1e2t(ji  ,jj) * ( pe3_in(ji  ,jj,jk) - e3t_0(ji  ,jj,jk) )     & 
     
    719742         DO jk = 1, jpk 
    720743            DO jj = 1, jpjm1 
    721                DO ji = 1, fs_jpim1   ! vector opt. 
     744               DO ji = 1, jpim1   ! vector opt. 
    722745                  pe3_out(ji,jj,jk) = 0.5_wp * ( vmask(ji,jj,jk)  * (1.0_wp - zlnwd) + zlnwd ) * r1_e1e2v(ji,jj)   & 
    723746                     &                       * (   e1e2t(ji,jj  ) * ( pe3_in(ji,jj  ,jk) - e3t_0(ji,jj  ,jk) )     & 
     
    732755         DO jk = 1, jpk 
    733756            DO jj = 1, jpjm1 
    734                DO ji = 1, fs_jpim1   ! vector opt. 
     757               DO ji = 1, jpim1   ! vector opt. 
    735758                  pe3_out(ji,jj,jk) = 0.5_wp * (  umask(ji,jj,jk) * umask(ji,jj+1,jk) * (1.0_wp - zlnwd) + zlnwd ) & 
    736759                     &                       *    r1_e1e2f(ji,jj)                                                  & 
     
    783806 
    784807 
    785    SUBROUTINE dom_vvl_rst( kt, cdrw ) 
     808   SUBROUTINE dom_vvl_rst( kt, Kbb, Kmm, cdrw ) 
    786809      !!--------------------------------------------------------------------- 
    787810      !!                   ***  ROUTINE dom_vvl_rst  *** 
     
    795818      !!                they are set to 0. 
    796819      !!---------------------------------------------------------------------- 
    797       INTEGER         , INTENT(in) ::   kt     ! ocean time-step 
    798       CHARACTER(len=*), INTENT(in) ::   cdrw   ! "READ"/"WRITE" flag 
     820      INTEGER         , INTENT(in) ::   kt        ! ocean time-step 
     821      INTEGER         , INTENT(in) ::   Kbb, Kmm  ! ocean time level indices 
     822      CHARACTER(len=*), INTENT(in) ::   cdrw      ! "READ"/"WRITE" flag 
    799823      ! 
    800824      INTEGER ::   ji, jj, jk 
     
    806830         IF( ln_rstart ) THEN                   !* Read the restart file 
    807831            CALL rst_read_open                  !  open the restart file if necessary 
    808             CALL iom_get( numror, jpdom_autoglo, 'sshn'   , sshn, ldxios = lrxios    ) 
     832            CALL iom_get( numror, jpdom_autoglo, 'sshn'   , ssh(:,:,Kmm), ldxios = lrxios    ) 
    809833            ! 
    810834            id1 = iom_varid( numror, 'e3t_b', ldstop = .FALSE. ) 
     
    817841            !                             ! --------- ! 
    818842            IF( MIN( id1, id2 ) > 0 ) THEN       ! all required arrays exist 
    819                CALL iom_get( numror, jpdom_autoglo, 'e3t_b', e3t_b(:,:,:), ldxios = lrxios ) 
    820                CALL iom_get( numror, jpdom_autoglo, 'e3t_n', e3t_n(:,:,:), ldxios = lrxios ) 
     843               CALL iom_get( numror, jpdom_autoglo, 'e3t_b', e3t(:,:,:,Kbb), ldxios = lrxios ) 
     844               CALL iom_get( numror, jpdom_autoglo, 'e3t_n', e3t(:,:,:,Kmm), ldxios = lrxios ) 
    821845               ! needed to restart if land processor not computed  
    822                IF(lwp) write(numout,*) 'dom_vvl_rst : e3t_b and e3t_n found in restart files' 
     846               IF(lwp) write(numout,*) 'dom_vvl_rst : e3t(:,:,:,Kbb) and e3t(:,:,:,Kmm) found in restart files' 
    823847               WHERE ( tmask(:,:,:) == 0.0_wp )  
    824                   e3t_n(:,:,:) = e3t_0(:,:,:) 
    825                   e3t_b(:,:,:) = e3t_0(:,:,:) 
     848                  e3t(:,:,:,Kmm) = e3t_0(:,:,:) 
     849                  e3t(:,:,:,Kbb) = e3t_0(:,:,:) 
    826850               END WHERE 
    827851               IF( neuler == 0 ) THEN 
    828                   e3t_b(:,:,:) = e3t_n(:,:,:) 
     852                  e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm) 
    829853               ENDIF 
    830854            ELSE IF( id1 > 0 ) THEN 
    831                IF(lwp) write(numout,*) 'dom_vvl_rst WARNING : e3t_n not found in restart files' 
     855               IF(lwp) write(numout,*) 'dom_vvl_rst WARNING : e3t(:,:,:,Kmm) not found in restart files' 
    832856               IF(lwp) write(numout,*) 'e3t_n set equal to e3t_b.' 
    833857               IF(lwp) write(numout,*) 'neuler is forced to 0' 
    834                CALL iom_get( numror, jpdom_autoglo, 'e3t_b', e3t_b(:,:,:), ldxios = lrxios ) 
    835                e3t_n(:,:,:) = e3t_b(:,:,:) 
     858               CALL iom_get( numror, jpdom_autoglo, 'e3t_b', e3t(:,:,:,Kbb), ldxios = lrxios ) 
     859               e3t(:,:,:,Kmm) = e3t(:,:,:,Kbb) 
    836860               neuler = 0 
    837861            ELSE IF( id2 > 0 ) THEN 
    838                IF(lwp) write(numout,*) 'dom_vvl_rst WARNING : e3t_b not found in restart files' 
     862               IF(lwp) write(numout,*) 'dom_vvl_rst WARNING : e3t(:,:,:,Kbb) not found in restart files' 
    839863               IF(lwp) write(numout,*) 'e3t_b set equal to e3t_n.' 
    840864               IF(lwp) write(numout,*) 'neuler is forced to 0' 
    841                CALL iom_get( numror, jpdom_autoglo, 'e3t_n', e3t_n(:,:,:), ldxios = lrxios ) 
    842                e3t_b(:,:,:) = e3t_n(:,:,:) 
     865               CALL iom_get( numror, jpdom_autoglo, 'e3t_n', e3t(:,:,:,Kmm), ldxios = lrxios ) 
     866               e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm) 
    843867               neuler = 0 
    844868            ELSE 
    845                IF(lwp) write(numout,*) 'dom_vvl_rst WARNING : e3t_n not found in restart file' 
     869               IF(lwp) write(numout,*) 'dom_vvl_rst WARNING : e3t(:,:,:,Kmm) not found in restart file' 
    846870               IF(lwp) write(numout,*) 'Compute scale factor from sshn' 
    847871               IF(lwp) write(numout,*) 'neuler is forced to 0' 
    848872               DO jk = 1, jpk 
    849                   e3t_n(:,:,jk) =  e3t_0(:,:,jk) * ( ht_0(:,:) + sshn(:,:) ) & 
     873                  e3t(:,:,jk,Kmm) =  e3t_0(:,:,jk) * ( ht_0(:,:) + ssh(:,:,Kmm) ) & 
    850874                      &                          / ( ht_0(:,:) + 1._wp - ssmask(:,:) ) * tmask(:,:,jk)   & 
    851875                      &          + e3t_0(:,:,jk)                               * (1._wp -tmask(:,:,jk)) 
    852876               END DO 
    853                e3t_b(:,:,:) = e3t_n(:,:,:) 
     877               e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm) 
    854878               neuler = 0 
    855879            ENDIF 
     
    888912               IF( cn_cfg == 'wad' ) THEN 
    889913                  ! Wetting and drying test case 
    890                   CALL usr_def_istate( gdept_b, tmask, tsb, ub, vb, sshb  ) 
    891                   tsn  (:,:,:,:) = tsb (:,:,:,:)       ! set now values from to before ones 
    892                   sshn (:,:)     = sshb(:,:) 
    893                   un   (:,:,:)   = ub  (:,:,:) 
    894                   vn   (:,:,:)   = vb  (:,:,:) 
     914                  CALL usr_def_istate( gdept(:,:,:,Kbb), tmask, ts(:,:,:,:,Kbb), uu(:,:,:,Kbb), vv(:,:,:,Kbb), ssh(:,:,Kbb)  ) 
     915                  ts  (:,:,:,:,Kmm) = ts (:,:,:,:,Kbb)       ! set now values from to before ones 
     916                  ssh (:,:,Kmm)     = ssh(:,:,Kbb) 
     917                  uu   (:,:,:,Kmm)   = uu  (:,:,:,Kbb) 
     918                  vv   (:,:,:,Kmm)   = vv  (:,:,:,Kbb) 
    895919               ELSE 
    896920                  ! if not test case 
    897                   sshn(:,:) = -ssh_ref 
    898                   sshb(:,:) = -ssh_ref 
     921                  ssh(:,:,Kmm) = -ssh_ref 
     922                  ssh(:,:,Kbb) = -ssh_ref 
    899923 
    900924                  DO jj = 1, jpj 
    901925                     DO ji = 1, jpi 
    902926                        IF( ht_0(ji,jj)-ssh_ref <  rn_wdmin1 ) THEN ! if total depth is less than min depth 
    903  
    904                            sshb(ji,jj) = rn_wdmin1 - (ht_0(ji,jj) ) 
    905                            sshn(ji,jj) = rn_wdmin1 - (ht_0(ji,jj) ) 
    906                            ssha(ji,jj) = rn_wdmin1 - (ht_0(ji,jj) ) 
     927                           ssh(ji,jj,Kbb) = rn_wdmin1 - (ht_0(ji,jj) ) 
     928                           ssh(ji,jj,Kmm) = rn_wdmin1 - (ht_0(ji,jj) ) 
    907929                        ENDIF 
    908930                     ENDDO 
     
    912934               ! Adjust vertical metrics for all wad 
    913935               DO jk = 1, jpk 
    914                   e3t_n(:,:,jk) =  e3t_0(:,:,jk) * ( ht_0(:,:) + sshn(:,:)  ) & 
     936                  e3t(:,:,jk,Kmm) =  e3t_0(:,:,jk) * ( ht_0(:,:) + ssh(:,:,Kmm)  ) & 
    915937                    &                            / ( ht_0(:,:) + 1._wp - ssmask(:,:) ) * tmask(:,:,jk)   & 
    916938                    &            + e3t_0(:,:,jk) * ( 1._wp - tmask(:,:,jk) ) 
    917939               END DO 
    918                e3t_b(:,:,:) = e3t_n(:,:,:) 
     940               e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm) 
    919941 
    920942               DO ji = 1, jpi 
     
    928950            ELSE 
    929951               ! 
    930                ! usr_def_istate called here only to get sshb, that is needed to initialize e3t_b and e3t_n 
    931                CALL usr_def_istate( gdept_0, tmask, tsb, ub, vb, sshb  )   
    932                ! usr_def_istate will be called again in istate_init to initialize ts(bn), ssh(bn), u(bn) and v(bn) 
     952               ! usr_def_istate called here only to get ssh(Kbb) needed to initialize e3t(Kbb) and e3t(Kmm) 
     953               ! 
     954               CALL usr_def_istate( gdept_0, tmask, ts(:,:,:,:,Kbb), uu(:,:,:,Kbb), vv(:,:,:,Kbb), ssh(:,:,Kbb)  )   
     955               ! 
     956               ! usr_def_istate will be called again in istate_init to initialize ts, ssh, u and v 
    933957               ! 
    934958               DO jk=1,jpk 
    935                   e3t_b(:,:,jk) =  e3t_0(:,:,jk) * ( ht_0(:,:) + sshb(:,:) ) & 
     959                  e3t(:,:,jk,Kbb) =  e3t_0(:,:,jk) * ( ht_0(:,:) + ssh(:,:,Kbb) ) & 
    936960                    &                            / ( ht_0(:,:) + 1._wp - ssmask(:,:) ) * tmask(:,:,jk)   & 
    937                     &            + e3t_0(:,:,jk) * ( 1._wp - tmask(:,:,jk) )   ! make sure e3t_b != 0 on land points 
     961                    &            + e3t_0(:,:,jk) * ( 1._wp - tmask(:,:,jk) )   ! make sure e3t(:,:,:,Kbb) != 0 on land points 
    938962               END DO 
    939                e3t_n(:,:,:) = e3t_b(:,:,:) 
    940                sshn(:,:) = sshb(:,:)   ! needed later for gde3w 
    941 !!$                e3t_n(:,:,:)=e3t_0(:,:,:) 
    942 !!$                e3t_b(:,:,:)=e3t_0(:,:,:) 
     963               e3t(:,:,:,Kmm) = e3t(:,:,:,Kbb) 
     964               ssh(:,:,Kmm) = ssh(:,:,Kbb)                                     ! needed later for gde3w 
    943965               ! 
    944966            END IF           ! end of ll_wd edits 
     
    958980         !                                           ! all cases ! 
    959981         !                                           ! --------- ! 
    960          CALL iom_rstput( kt, nitrst, numrow, 'e3t_b', e3t_b(:,:,:), ldxios = lwxios ) 
    961          CALL iom_rstput( kt, nitrst, numrow, 'e3t_n', e3t_n(:,:,:), ldxios = lwxios ) 
     982         CALL iom_rstput( kt, nitrst, numrow, 'e3t_b', e3t(:,:,:,Kbb), ldxios = lwxios ) 
     983         CALL iom_rstput( kt, nitrst, numrow, 'e3t_n', e3t(:,:,:,Kmm), ldxios = lwxios ) 
    962984         !                                           ! ----------------------- ! 
    963985         IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN  ! z_tilde and layer cases ! 
     
    9921014      !!----------------------------------------------------------------------  
    9931015      ! 
    994       REWIND( numnam_ref )              ! Namelist nam_vvl in reference namelist :  
    9951016      READ  ( numnam_ref, nam_vvl, IOSTAT = ios, ERR = 901) 
    9961017901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nam_vvl in reference namelist' ) 
    997       REWIND( numnam_cfg )              ! Namelist nam_vvl in configuration namelist : Parameters of the run 
    9981018      READ  ( numnam_cfg, nam_vvl, IOSTAT = ios, ERR = 902 ) 
    9991019902   IF( ios >  0 ) CALL ctl_nam ( ios , 'nam_vvl in configuration namelist' ) 
     
    10341054      ! 
    10351055      IF( ioptio /= 1 )   CALL ctl_stop( 'Choose ONE vertical coordinate in namelist nam_vvl' ) 
    1036       IF( .NOT. ln_vvl_zstar .AND. ln_isf ) CALL ctl_stop( 'Only vvl_zstar has been tested with ice shelf cavity' ) 
    10371056      ! 
    10381057      IF(lwp) THEN                   ! Print the choice 
  • NEMO/trunk/tests/VORTEX/MY_SRC/usrdef_nam.F90

    r11536 r12377  
    6363      !!---------------------------------------------------------------------- 
    6464      ! 
    65       REWIND( numnam_cfg )          ! Namelist namusr_def (exist in namelist_cfg only) 
    6665      READ  ( numnam_cfg, namusr_def, IOSTAT = ios, ERR = 902 ) 
    6766902   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namusr_def in configuration namelist' ) 
  • NEMO/trunk/tests/VORTEX/MY_SRC/usrdef_sbc.F90

    r10074 r12377  
    3030   PUBLIC   usrdef_sbc_ice_flx  ! routine called by icestp.F90 for ice thermo 
    3131 
    32    !! * Substitutions 
    33 #  include "vectopt_loop_substitute.h90" 
    3432   !!---------------------------------------------------------------------- 
    3533   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
     
    3937CONTAINS 
    4038 
    41    SUBROUTINE usrdef_sbc_oce( kt ) 
     39   SUBROUTINE usrdef_sbc_oce( kt, Kbb ) 
    4240      !!--------------------------------------------------------------------- 
    4341      !!                    ***  ROUTINE usr_def_sbc  *** 
     
    5452      !!---------------------------------------------------------------------- 
    5553      INTEGER, INTENT(in) ::   kt   ! ocean time step 
     54      INTEGER, INTENT(in) ::   Kbb  ! ocean time index 
    5655      !!--------------------------------------------------------------------- 
    5756      ! 
     
    7978   END SUBROUTINE usrdef_sbc_ice_tau 
    8079 
    81    SUBROUTINE usrdef_sbc_ice_flx( kt ) 
     80 
     81   SUBROUTINE usrdef_sbc_ice_flx( kt, phs, phi ) 
    8282      INTEGER, INTENT(in) ::   kt   ! ocean time step 
     83      REAL(wp), DIMENSION(:,:,:), INTENT(in)  ::   phs    ! snow thickness 
     84      REAL(wp), DIMENSION(:,:,:), INTENT(in)  ::   phi    ! ice thickness 
    8385   END SUBROUTINE usrdef_sbc_ice_flx 
    8486 
  • NEMO/trunk/tests/VORTEX/MY_SRC/usrdef_zgr.F90

    r10425 r12377  
    2929   PUBLIC   usr_def_zgr        ! called by domzgr.F90 
    3030 
    31   !! * Substitutions 
    32 #  include "vectopt_loop_substitute.h90" 
    3331   !!---------------------------------------------------------------------- 
    3432   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
Note: See TracChangeset for help on using the changeset viewer.