Index: NEMO/trunk/src/OCE/DOM/domvvl.F90
===================================================================
--- NEMO/trunk/src/OCE/DOM/domvvl.F90	(revision 13982)
+++ NEMO/trunk/src/OCE/DOM/domvvl.F90	(revision 14053)
@@ -9,5 +9,5 @@
    !!            3.6  !  2014-11  (P. Mathiot) add ice shelf capability
    !!            4.1  !  2019-08  (A. Coward, D. Storkey) rename dom_vvl_sf_swp -> dom_vvl_sf_update for new timestepping
-   !!            4.x  ! 2020-02  (G. Madec, S. Techene) introduce ssh to h0 ratio
+   !!             -   !  2020-02  (G. Madec, S. Techene) introduce ssh to h0 ratio
    !!----------------------------------------------------------------------
 
@@ -766,9 +766,18 @@
       !! ** Purpose :   Read or write VVL file in restart file
       !!
-      !! ** Method  :   use of IOM library
-      !!                if the restart does not contain vertical scale factors,
-      !!                they are set to the _0 values
-      !!                if the restart does not contain vertical scale factors increments (z_tilde),
-      !!                they are set to 0.
+      !! ** Method  : * restart comes from a linear ssh simulation :
+      !!                   an attempt to read e3t_n stops simulation
+      !!              * restart comes from a z-star, z-tilde, or layer :
+      !!                   read e3t_n and e3t_b
+      !!              * restart comes from a z-star :
+      !!                   set tilde_e3t_n, tilde_e3t_n, and hdiv_lf to 0
+      !!              * restart comes from layer :
+      !!                   read tilde_e3t_n and tilde_e3t_b
+      !!                   set hdiv_lf to 0
+      !!              * restart comes from a z-tilde:
+      !!                   read tilde_e3t_n, tilde_e3t_b, and hdiv_lf
+      !!
+      !!              NB: if l_1st_euler = T (ln_1st_euler or ssh_b not found)
+      !!                   Kbb fields set to Kmm ones
       !!----------------------------------------------------------------------
       INTEGER         , INTENT(in) ::   kt        ! ocean time-step
@@ -776,154 +785,91 @@
       CHARACTER(len=*), INTENT(in) ::   cdrw      ! "READ"/"WRITE" flag
       !
-      INTEGER ::   ji, jj, jk
-      INTEGER ::   id1, id2, id3, id4, id5     ! local integers
-      !!----------------------------------------------------------------------
-      !
-      IF( TRIM(cdrw) == 'READ' ) THEN        ! Read/initialise 
-         !                                   ! ===============
-         IF( ln_rstart ) THEN                   !* Read the restart file
-            CALL rst_read_open                  !  open the restart file if necessary
-            CALL iom_get( numror, jpdom_auto, 'sshn'   , ssh(:,:,Kmm)    )
+      INTEGER ::   ji, jj, jk      ! dummy loop indices
+      INTEGER ::   id3, id4, id5   ! local integers
+      !!----------------------------------------------------------------------
+      !
+      !                                      !=====================!
+      IF( TRIM(cdrw) == 'READ' ) THEN        !  Read / initialise  !
+         !                                   !=====================!
+         !
+         IF( ln_rstart ) THEN                   !==  Read the restart file  ==!
             !
-            id1 = iom_varid( numror, 'e3t_b', ldstop = .FALSE. )
-            id2 = iom_varid( numror, 'e3t_n', ldstop = .FALSE. )
-            id3 = iom_varid( numror, 'tilde_e3t_b', ldstop = .FALSE. )
+            CALL rst_read_open                                          !*  open the restart file if necessary
+            !                                         ! --------- !
+            !                                         ! all cases !
+            !                                         ! --------- !
+            !
+            id3 = iom_varid( numror, 'tilde_e3t_b', ldstop = .FALSE. )  !*  check presence
             id4 = iom_varid( numror, 'tilde_e3t_n', ldstop = .FALSE. )
-            id5 = iom_varid( numror, 'hdiv_lf', ldstop = .FALSE. )
+            id5 = iom_varid( numror, 'hdiv_lf'    , ldstop = .FALSE. )
             !
-            !                             ! --------- !
-            !                             ! all cases !
-            !                             ! --------- !
-            !
-            IF( MIN( id1, id2 ) > 0 ) THEN       ! all required arrays exist
+            !                                                           !*  scale factors
+            IF(lwp) WRITE(numout,*)    '          Kmm scale factor read in the restart file'
+            CALL iom_get( numror, jpdom_auto, 'e3t_n', e3t(:,:,:,Kmm) )
+            WHERE ( tmask(:,:,:) == 0.0_wp ) 
+               e3t(:,:,:,Kmm) = e3t_0(:,:,:)
+            END WHERE
+            IF( l_1st_euler ) THEN                       ! euler
+               IF(lwp) WRITE(numout,*) '          Euler first time step : e3t(Kbb) = e3t(Kmm)'
+               e3t(:,:,:,Kmm) = e3t(:,:,:,Kbb)
+            ELSE                                         ! leap frog
+               IF(lwp) WRITE(numout,*) '          Kbb scale factor read in the restart file'
                CALL iom_get( numror, jpdom_auto, 'e3t_b', e3t(:,:,:,Kbb) )
-               CALL iom_get( numror, jpdom_auto, 'e3t_n', e3t(:,:,:,Kmm) )
-               ! needed to restart if land processor not computed 
-               IF(lwp) write(numout,*) 'dom_vvl_rst : e3t(:,:,:,Kbb) and e3t(:,:,:,Kmm) found in restart files'
                WHERE ( tmask(:,:,:) == 0.0_wp ) 
-                  e3t(:,:,:,Kmm) = e3t_0(:,:,:)
                   e3t(:,:,:,Kbb) = e3t_0(:,:,:)
                END WHERE
-               IF( l_1st_euler ) THEN
-                  e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm)
-               ENDIF
-            ELSE IF( id1 > 0 ) THEN
-               IF(lwp) write(numout,*) 'dom_vvl_rst WARNING : e3t(:,:,:,Kmm) not found in restart files'
-               IF(lwp) write(numout,*) 'e3t_n set equal to e3t_b.'
-               IF(lwp) write(numout,*) 'l_1st_euler is forced to true'
-               CALL iom_get( numror, jpdom_auto, 'e3t_b', e3t(:,:,:,Kbb) )
-               e3t(:,:,:,Kmm) = e3t(:,:,:,Kbb)
-               l_1st_euler = .true.
-            ELSE IF( id2 > 0 ) THEN
-               IF(lwp) write(numout,*) 'dom_vvl_rst WARNING : e3t(:,:,:,Kbb) not found in restart files'
-               IF(lwp) write(numout,*) 'e3t_b set equal to e3t_n.'
-               IF(lwp) write(numout,*) 'l_1st_euler is forced to true'
-               CALL iom_get( numror, jpdom_auto, 'e3t_n', e3t(:,:,:,Kmm) )
-               e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm)
-               l_1st_euler = .true.
-            ELSE
-               IF(lwp) write(numout,*) 'dom_vvl_rst WARNING : e3t(:,:,:,Kmm) not found in restart file'
-               IF(lwp) write(numout,*) 'Compute scale factor from sshn'
-               IF(lwp) write(numout,*) 'l_1st_euler is forced to true'
-               DO jk = 1, jpk
-                  e3t(:,:,jk,Kmm) =  e3t_0(:,:,jk) * ( ht_0(:,:) + ssh(:,:,Kmm) ) &
-                      &                          / ( ht_0(:,:) + 1._wp - ssmask(:,:) ) * tmask(:,:,jk)   &
-                      &          + e3t_0(:,:,jk)                               * (1._wp -tmask(:,:,jk))
-               END DO
-               e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm)
-               l_1st_euler = .true.
             ENDIF
-            !                             ! ----------- !
-            IF( ln_vvl_zstar ) THEN       ! z_star case !
-               !                          ! ----------- !
+            !                                         ! ------------ !
+            IF( ln_vvl_zstar ) THEN                   !  z_star case !
+               !                                      ! ------------ !
                IF( MIN( id3, id4 ) > 0 ) THEN
                   CALL ctl_stop( 'dom_vvl_rst: z_star cannot restart from a z_tilde or layer run' )
                ENDIF
-               !                          ! ----------------------- !
-            ELSE                          ! z_tilde and layer cases !
-               !                          ! ----------------------- !
-               IF( MIN( id3, id4 ) > 0 ) THEN  ! all required arrays exist
-                  CALL iom_get( numror, jpdom_auto, 'tilde_e3t_b', tilde_e3t_b(:,:,:) )
+               !                                      ! ------------------------ !
+            ELSE                                      !  z_tilde and layer cases !
+               !                                      ! ------------------------ !
+               !
+               IF( id4 > 0 ) THEN                                       !*  scale factor increments
+                  IF(lwp) WRITE(numout,*)    '          Kmm scale factor increments read in the restart file'
                   CALL iom_get( numror, jpdom_auto, 'tilde_e3t_n', tilde_e3t_n(:,:,:) )
-               ELSE                            ! one at least array is missing
+                  IF( l_1st_euler ) THEN                 ! euler
+                     IF(lwp) WRITE(numout,*) '          Euler first time step : tilde_e3t(Kbb) = tilde_e3t(Kmm)'
+                     tilde_e3t_b(:,:,:) = tilde_e3t_n(:,:,:)
+                  ELSE                                   ! leap frog
+                     IF(lwp) WRITE(numout,*) '          Kbb scale factor increments read in the restart file'
+                     CALL iom_get( numror, jpdom_auto, 'tilde_e3t_b', tilde_e3t_b(:,:,:) )
+                  ENDIF
+               ELSE 
                   tilde_e3t_b(:,:,:) = 0.0_wp
                   tilde_e3t_n(:,:,:) = 0.0_wp
                ENDIF
-               !                          ! ------------ !
-               IF( ln_vvl_ztilde ) THEN   ! z_tilde case !
-                  !                       ! ------------ !
+               !                                      ! ------------ !
+               IF( ln_vvl_ztilde ) THEN               ! z_tilde case !
+                  !                                   ! ------------ !
                   IF( id5 > 0 ) THEN  ! required array exists
                      CALL iom_get( numror, jpdom_auto, 'hdiv_lf', hdiv_lf(:,:,:) )
                   ELSE                ! array is missing
-                     hdiv_lf(:,:,:) = 0.0_wp
+                     hdiv_lf(:,:,:) = 0.0_wp 
                   ENDIF
                ENDIF
             ENDIF
             !
-         ELSE                                   !* Initialize at "rest"
+         ELSE                                   !==  Initialize at "rest" with ssh  ==!
             !
-
-            IF( ll_wd ) THEN   ! MJB ll_wd edits start here - these are essential 
-               !
-               IF( cn_cfg == 'wad' ) THEN
-                  ! Wetting and drying test case
-                  CALL usr_def_istate( gdept(:,:,:,Kbb), tmask, ts(:,:,:,:,Kbb), uu(:,:,:,Kbb), vv(:,:,:,Kbb), ssh(:,:,Kbb)  )
-                  ts  (:,:,:,:,Kmm) = ts (:,:,:,:,Kbb)       ! set now values from to before ones
-                  ssh (:,:,Kmm)     = ssh(:,:,Kbb)
-                  uu   (:,:,:,Kmm)   = uu  (:,:,:,Kbb)
-                  vv   (:,:,:,Kmm)   = vv  (:,:,:,Kbb)
-               ELSE
-                  ! if not test case
-                  ssh(:,:,Kmm) = -ssh_ref
-                  ssh(:,:,Kbb) = -ssh_ref
-
-                  DO_2D( 1, 1, 1, 1 )
-                     IF( ht_0(ji,jj)-ssh_ref <  rn_wdmin1 ) THEN ! if total depth is less than min depth
-                        ssh(ji,jj,Kbb) = rn_wdmin1 - (ht_0(ji,jj) )
-                        ssh(ji,jj,Kmm) = rn_wdmin1 - (ht_0(ji,jj) )
-                     ENDIF
-                  END_2D
-               ENDIF !If test case else
-
-               ! Adjust vertical metrics for all wad
-               DO jk = 1, jpk
-                  e3t(:,:,jk,Kmm) =  e3t_0(:,:,jk) * ( ht_0(:,:) + ssh(:,:,Kmm)  ) &
-                    &                            / ( ht_0(:,:) + 1._wp - ssmask(:,:) ) * tmask(:,:,jk)   &
-                    &            + e3t_0(:,:,jk) * ( 1._wp - tmask(:,:,jk) )
-               END DO
-               e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm)
-
-               DO_2D( 1, 1, 1, 1 )
-                  IF ( ht_0(ji,jj) .LE. 0.0 .AND. NINT( ssmask(ji,jj) ) .EQ. 1) THEN
-                     CALL ctl_stop( 'dom_vvl_rst: ht_0 must be positive at potentially wet points' )
-                  ENDIF
-               END_2D
-               !
-            ELSE
-               !
-               ! Just to read set ssh in fact, called latter once vertical grid
-               ! is set up:
-!               CALL usr_def_istate( gdept_0, tmask, ts(:,:,:,:,Kbb), uu(:,:,:,Kbb), vv(:,:,:,Kbb), ssh(:,:,Kbb)  )
-!               !
-!               DO jk=1,jpk
-!                  e3t(:,:,jk,Kbb) =  e3t_0(:,:,jk) * ( ht_0(:,:) + ssh(:,:,Kbb) ) &
-!                     &            / ( ht_0(:,:) + 1._wp -ssmask(:,:) ) * tmask(:,:,jk)
-!               END DO
-!               e3t(:,:,:,Kmm) = e3t(:,:,:,Kbb)
-                ssh(:,:,Kmm)=0._wp
-                e3t(:,:,:,Kmm)=e3t_0(:,:,:)
-                e3t(:,:,:,Kbb)=e3t_0(:,:,:)
-               !
-            END IF           ! end of ll_wd edits
-
+            DO jk = 1, jpk
+               e3t(:,:,jk,Kmm) =  e3t_0(:,:,jk) * ( 1._wp + ssh(:,:,Kmm) * r1_ht_0(:,:) * tmask(:,:,jk) )
+            END DO
+            e3t(:,:,:,Kbb) = e3t(:,:,:,Kmm)
+            !
             IF( ln_vvl_ztilde .OR. ln_vvl_layer) THEN
                tilde_e3t_b(:,:,:) = 0._wp
                tilde_e3t_n(:,:,:) = 0._wp
                IF( ln_vvl_ztilde ) hdiv_lf(:,:,:) = 0._wp
-            END IF
+            ENDIF
          ENDIF
-         !
-      ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN   ! Create restart file
-         !                                   ! ===================
+         !                                       !=======================!
+      ELSEIF( TRIM(cdrw) == 'WRITE' ) THEN       !  Create restart file  !
+         !                                       !=======================!
+         !
          IF(lwp) WRITE(numout,*) '---- dom_vvl_rst ----'
          !                                           ! --------- !
