From cd89b1dd07e30149eaa025a0a5456c986abf3163 Mon Sep 17 00:00:00 2001
From: Sibylle Techene <sibylle.techene@locean-ipsl.upmc.fr>
Date: Tue, 24 Jan 2023 13:56:48 +0000
Subject: [PATCH] Resolve "Fix issue introduced with vertical keys"

---
 src/OCE/DOM/dom_oce.F90           | 20 ++++---
 src/OCE/DOM/domzgr.F90            | 45 ++++++++--------
 src/OCE/USR/usrdef_zgr.F90        | 86 ++++++++++++++++++++++---------
 tests/BENCH/MY_SRC/usrdef_zgr.F90 | 80 +++++++++++++++++++++-------
 4 files changed, 153 insertions(+), 78 deletions(-)

diff --git a/src/OCE/DOM/dom_oce.F90 b/src/OCE/DOM/dom_oce.F90
index 75a49ead..05742adc 100644
--- a/src/OCE/DOM/dom_oce.F90
+++ b/src/OCE/DOM/dom_oce.F90
@@ -136,19 +136,19 @@ MODULE dom_oce
    LOGICAL, PUBLIC, PARAMETER ::   lk_ALE    = .FALSE.  !: ALE key flag
 #endif
 #if defined key_vco_1d
-   LOGICAL, PUBLIC, PARAMETER ::   lk_vco_1d   = .TRUE.   !: zco key flag
+   LOGICAL, PUBLIC, PARAMETER ::   lk_vco_1d   = .TRUE.   !: 1d key flag
 #else
-   LOGICAL, PUBLIC, PARAMETER ::   lk_vco_1d   = .FALSE.  !: zco key flag
+   LOGICAL, PUBLIC, PARAMETER ::   lk_vco_1d   = .FALSE.  !: 1d key flag
 #endif
 #if defined key_vco_1d3d
-   LOGICAL, PUBLIC, PARAMETER ::   lk_vco_1d3d = .TRUE.   !: zps key flag
+   LOGICAL, PUBLIC, PARAMETER ::   lk_vco_1d3d = .TRUE.   !: 1d3d key flag
 #else
-   LOGICAL, PUBLIC, PARAMETER ::   lk_vco_1d3d = .FALSE.  !: zps key flag
+   LOGICAL, PUBLIC, PARAMETER ::   lk_vco_1d3d = .FALSE.  !: 1d3d key flag
 #endif
 #if defined key_vco_3d
-   LOGICAL, PUBLIC, PARAMETER ::   lk_vco_3d   = .TRUE.   !: sco key flag
+   LOGICAL, PUBLIC, PARAMETER ::   lk_vco_3d   = .TRUE.   !: 3d key flag
 #else
-   LOGICAL, PUBLIC, PARAMETER ::   lk_vco_3d   = .FALSE.  !: sco key flag
+   LOGICAL, PUBLIC, PARAMETER ::   lk_vco_3d   = .FALSE.  !: 3d key flag
 #endif
 
 !!gm obsolescent feature replaced by key_xxx ==>>>  to be removed when z-tilde and or ALE key added (and domvvl removed)
@@ -362,15 +362,14 @@ CONTAINS
             !
             ALLOCATE(        gdept_1d(jpk) ,        gdepw_1d(jpk) ,       &
                &               e3t_1d(jpk) ,          e3w_1d(jpk) ,   STAT=ierr(ii) )
-               !
-         ELSEIF( lk_vco_1d3d ) THEN       !* zps :   allocate 1d vertical arrays, except t-level e3 !!st WHT not ??? 
-            !                 
+         ELSEIF( lk_vco_1d3d ) THEN
+            !                             !* zps :   allocate 1d vertical arrays for gdep and w-level e3 fields and t-level e3 fields
             ALLOCATE(        gdept_1d(jpk) ,        gdepw_1d(jpk) ,       &
                &               e3t_1d(jpk) ,          e3w_1d(jpk) ,       &
                &       e3t_3d(jpi,jpj,jpk) ,  e3u_3d(jpi,jpj,jpk) ,       &
                &       e3v_3d(jpi,jpj,jpk) ,  e3f_3d(jpi,jpj,jpk) ,   STAT=ierr(ii) )
          ELSEIF( lk_vco_3d ) THEN
-            !                             !* sco :   allocate 3d vertical arrays for all gdep and e3 fields (no more _1d) !!st WHT not ???
+            !                             !* sco :   allocate 3d vertical arrays for all gdep and e3 fields (no more _1d)
             ALLOCATE(        gdept_1d(jpk) ,        gdepw_1d(jpk) ,       &
                &               e3t_1d(jpk) ,          e3w_1d(jpk) ,       &
                &     gdept_3d(jpi,jpj,jpk) ,gdepw_3d(jpi,jpj,jpk) ,       &
@@ -382,7 +381,6 @@ CONTAINS
          !                                !-------------------------------------!
       ELSEIF( lk_ALE ) THEN               !-  combine time & space variations  -!   (vertical ALE coordinate)
          !                                !-------------------------------------!      NOT yet implemented
-         !!st ca ne devrait pas etre des *_0 qui varient dans le temps ? 
          ii = ii+1
          ALLOCATE(    ht(jpi,jpj,jpt) ,    hu(jpi,jpj,jpt) ,      hv(jpi,jpj,jpt) ,       &
             &                           r1_hu(jpi,jpj,jpt) , r1_hv  (jpi,jpj,jpt) ,   STAT=ierr(ii)  )
diff --git a/src/OCE/DOM/domzgr.F90 b/src/OCE/DOM/domzgr.F90
index 9029c88e..ca1139c2 100644
--- a/src/OCE/DOM/domzgr.F90
+++ b/src/OCE/DOM/domzgr.F90
@@ -87,7 +87,6 @@ CONTAINS
          WRITE(numout,*) '~~~~~~~'
          IF( ln_linssh ) WRITE(numout,*) '          linear free surface: the vertical mesh does not change in time'
       ENDIF
-CALL FLUSH(numout)
       !                             !==============================!
       IF( ln_read_cfg ) THEN        !==  read in domcfg.nc file  ==!
          !                          !==============================!
@@ -229,46 +228,44 @@ CALL FLUSH(numout)
             !                                !--------------------!
             IF( lk_vco_1d ) THEN             !--  z-coordinate  --!   use only 1D arrays for all gdep and e3 fields
                !                             !--------------------!
-               l_zco = .TRUE.                     ! old logical   ==> to be removed
-               l_zps = .FALSE.
-               l_sco = .FALSE.
                !
+               
                CALL usr_def_zgr( l_zco  , l_zps  , l_sco, ln_isfcav,   &
                   &              k_top   , k_bot                      ,   &    ! 1st & last ocean level
                   &              gdept_1d, gdepw_1d, e3t_1d, e3w_1d   )        ! 1D gridpoints depth
+               !
+               IF( l_sco )   CALL ctl_stop( 'STOP','domzgr: key_vco_1d and l_sco=T are incompatible. Fix usrdef_zgr !' )
+               IF( l_zps )   CALL ctl_stop( 'STOP','domzgr: key_vco_1d and l_zps=T are incompatible. Fix usrdef_zgr !' )
+               !
                !                             !-----------------------!
             ELSEIF( lk_vco_1d3d ) THEN       !--  z-partial cells  --!   use 3D t-level e3
                !                             !-----------------------!
-               l_zco = .FALSE.                     ! old logical   ==> to be removed
-               l_zps = .TRUE.
-               l_sco = .FALSE.
                !
-               CALL usr_def_zgr( l_zco  , l_zps  , l_sco, ln_isfcav,   &
+               CALL usr_def_zgr( l_zco  , l_zps  , l_sco, ln_isfcav   ,   &
                   &              k_top   , k_bot                      ,   &    ! 1st & last ocean level
                   &              gdept_1d, gdepw_1d, e3t_1d, e3w_1d   ,   &    ! 1D gridpoints depth
-                  &              e3t_3d  , e3u_3d  , e3v_3d, e3f_3d    )        ! vertical scale factors
+                  &              e3t_3d  , e3u_3d  , e3v_3d, e3f_3d       )    ! vertical scale factors
+               !
+               IF( l_sco )   CALL ctl_stop( 'STOP','domzgr: key_vco_1d3d and l_sco=T are incompatible. Fix usrdef_zgr !' )
                !
                ! make sure that periodicities are properly applied 
                CALL lbc_lnk( 'dom_zgr', e3t_3d, 'T', 1._wp,   e3u_3d, 'U', 1._wp,  e3v_3d, 'V', 1._wp, e3f_3d, 'F', 1._wp,   &
                     &                     kfillmode = jpfillcopy )   ! do not put 0 over closed boundaries
-
+               !
                !                             !--------------------!               
             ELSEIF( lk_vco_3d ) THEN         !--  s-coordinate  --!   use 3D for all gdep and e3 fields
                !                             !--------------------!
-               l_zco = .FALSE.                     ! old logical   ==> to be removed
-               l_zps = .FALSE.
-               l_sco = .TRUE.
                !
-               CALL usr_def_zgr( l_zco  , l_zps  , l_sco, ln_isfcav,   &
-               &              k_top   , k_bot                      ,   &    ! 1st & last ocean level
-               &              gdept_1d, gdepw_1d, e3t_1d, e3w_1d   ,   &    ! 1D gridpoints depth
-               &              e3t_3d  , e3u_3d  , e3v_3d, e3f_3d    ,   &    ! vertical scale factors
-               &              gdept_3d, gdepw_3d                    ,   &    ! gridpoints depth 
-               &              e3w_3d  , e3uw_3d , e3vw_3d           )        ! vertical scale factors
-               CALL lbc_lnk( 'dom_zgr', gdept_3d, 'T', 1._wp, gdepw_3d, 'W', 1._wp,                                          &
-               &                          e3t_3d, 'T', 1._wp,   e3u_3d, 'U', 1._wp,  e3v_3d, 'V', 1._wp, e3f_3d, 'F', 1._wp,   &
-               &                          e3w_3d, 'W', 1._wp,  e3uw_3d, 'U', 1._wp, e3vw_3d, 'V', 1._wp,   &   
-               &                     kfillmode = jpfillcopy )   ! do not put 0 over closed boundaries
+               CALL usr_def_zgr( l_zco  , l_zps  , l_sco, ln_isfcav   ,   &
+                  &              k_top   , k_bot                      ,   &    ! 1st & last ocean level
+                  &              gdept_1d, gdepw_1d, e3t_1d, e3w_1d   ,   &    ! 1D gridpoints depth
+                  &              e3t_3d  , e3u_3d  , e3v_3d, e3f_3d   ,   &    ! vertical scale factors
+                  &              gdept_3d, gdepw_3d                   ,   &    ! gridpoints depth 
+                  &              e3w_3d  , e3uw_3d , e3vw_3d               )   ! vertical scale factors
+               CALL lbc_lnk( 'dom_zgr', gdept_3d, 'T', 1._wp, gdepw_3d, 'W', 1._wp,                                            &
+                  &                       e3t_3d, 'T', 1._wp,   e3u_3d, 'U', 1._wp,  e3v_3d, 'V', 1._wp, e3f_3d, 'F', 1._wp,   &
+                  &                       e3w_3d, 'W', 1._wp,  e3uw_3d, 'U', 1._wp, e3vw_3d, 'V', 1._wp,                       &   
+                  &                     kfillmode = jpfillcopy )   ! do not put 0 over closed boundaries
             ENDIF
          ENDIF
          ztopbot(:,:,1) = REAL(k_top, wp)
@@ -311,7 +308,7 @@ CALL FLUSH(numout)
       !
       IF(lwp) THEN                     ! Control print
          WRITE(numout,*)
-         WRITE(numout,*) '   Type of vertical coordinate (read in ', TRIM( cn_domcfg ), ' file or set in userdef_nam) :'
+         WRITE(numout,*) '   Type of vertical coordinate (read in ', TRIM( cn_domcfg ), ' file or set in userdef_zgr) :'
          WRITE(numout,*) '      z-coordinate - full steps      l_zco    = ', l_zco
          WRITE(numout,*) '      z-coordinate - partial steps   l_zps    = ', l_zps
          WRITE(numout,*) '      s- or hybrid z-s-coordinate    l_sco    = ', l_sco
diff --git a/src/OCE/USR/usrdef_zgr.F90 b/src/OCE/USR/usrdef_zgr.F90
index 4523a2ff..47dfff1e 100644
--- a/src/OCE/USR/usrdef_zgr.F90
+++ b/src/OCE/USR/usrdef_zgr.F90
@@ -38,23 +38,23 @@ CONTAINS
    SUBROUTINE usr_def_zgr( ld_zco  , ld_zps  , ld_sco  , ld_isfcav,    &   ! type of vertical coordinate
       &                    k_top   , k_bot                        ,    &   ! top & bottom ocean level
       &                    pdept_1d, pdepw_1d, pe3t_1d , pe3w_1d  ,    &   ! 1D reference vertical coordinate
+      &                    pe3t  , pe3u  , pe3v   , pe3f ,             &   ! 3D t-level vertical scale factors
       &                    pdept , pdepw ,                             &   ! 3D t & w-points depth
-      &                    pe3t  , pe3u  , pe3v   , pe3f ,             &   ! vertical scale factors
-      &                    pe3w  , pe3uw , pe3vw                       )   !     -      -      -
+      &                    pe3w  , pe3uw , pe3vw                       )   ! 3D w-level vertical scale factors
       !!---------------------------------------------------------------------
       !!              ***  ROUTINE usr_def_zgr  ***
       !!
       !! ** Purpose :   User defined the vertical coordinates
       !!
       !!----------------------------------------------------------------------
-      LOGICAL                   , INTENT(out) ::   ld_zco, ld_zps, ld_sco      ! vertical coordinate flags
-      LOGICAL                   , INTENT(out) ::   ld_isfcav                   ! under iceshelf cavity flag
-      INTEGER , DIMENSION(:,:)  , INTENT(out) ::   k_top, k_bot                ! first & last ocean level
-      REAL(wp), DIMENSION(:)    , INTENT(out) ::   pdept_1d, pdepw_1d          ! 1D grid-point depth     [m]
-      REAL(wp), DIMENSION(:)    , INTENT(out) ::   pe3t_1d , pe3w_1d           ! 1D grid-point depth     [m]
-      REAL(wp), DIMENSION(:,:,:), OPTIONAL, INTENT(out) ::   pdept, pdepw                ! grid-point depth        [m]
-      REAL(wp), DIMENSION(:,:,:), OPTIONAL, INTENT(out) ::   pe3t , pe3u , pe3v , pe3f   ! vertical scale factors  [m]
-      REAL(wp), DIMENSION(:,:,:), OPTIONAL, INTENT(out) ::   pe3w , pe3uw, pe3vw         ! i-scale factors 
+      LOGICAL                             , INTENT(out) ::   ld_zco, ld_zps, ld_sco      ! vertical coordinate flags
+      LOGICAL                             , INTENT(out) ::   ld_isfcav                   ! under iceshelf cavity flag
+      INTEGER , DIMENSION(:,:)            , INTENT(out) ::   k_top, k_bot                ! first & last ocean level
+      REAL(wp), DIMENSION(:)              , INTENT(out) ::   pdept_1d, pdepw_1d          ! 1D grid-point depth             [m]
+      REAL(wp), DIMENSION(:)              , INTENT(out) ::   pe3t_1d , pe3w_1d           ! 1D grid-point depth             [m]
+      REAL(wp), DIMENSION(:,:,:), OPTIONAL, INTENT(out) ::   pe3t , pe3u , pe3v , pe3f   ! t-level vertical scale factors  [m]
+      REAL(wp), DIMENSION(:,:,:), OPTIONAL, INTENT(out) ::   pdept, pdepw                ! grid-point depth                [m]
+      REAL(wp), DIMENSION(:,:,:), OPTIONAL, INTENT(out) ::   pe3w , pe3uw, pe3vw         ! w-level vertical scale factors  [m]
       !!----------------------------------------------------------------------
       !
       IF(lwp) WRITE(numout,*)
@@ -76,11 +76,16 @@ CONTAINS
       !
       CALL zgr_msk_top_bot( k_top , k_bot )                 ! masked top and bottom ocean t-level indices
       !
-      IF( PRESENT( pe3t ) ) THEN                            ! z-coordinate (3D arrays) from the 1D z-coord.
-         CALL zgr_zco( pdept_1d, pdepw_1d, pe3t_1d, pe3w_1d,   &   ! in  : 1D reference vertical coordinate
-            &          pdept   , pdepw   ,                     &   ! out : 3D t & w-points depth
-            &          pe3t    , pe3u    , pe3v   , pe3f   ,   &   !       vertical scale factors
-            &          pe3w    , pe3uw   , pe3vw             )     !           -      -      -
+      IF( lk_vco_1d3d ) THEN                                ! z-coordinate (3D arrays) from the 1D z-coord.
+         CALL zgr_zco_1d3d( pe3t_1d,                   &               ! in  : 1D reference vertical coordinate
+            &               pe3t , pe3u , pe3v , pe3f  )               ! out : 3D vertical scale factors at t-level
+      ENDIF
+      !
+      IF( lk_vco_3d ) THEN                                  ! z-coordinate (3D arrays) from the 1D z-coord.
+         CALL zgr_zco_3d(   pdept_1d, pdepw_1d, pe3t_1d, pe3w_1d,   &   ! in  : 1D reference vertical coordinate
+            &               pe3t    , pe3u    , pe3v   , pe3f   ,   &   ! out : 3D vertical scale factors at t-level
+            &               pdept   , pdepw   ,                     &   !       3D t & w-points depth
+            &               pe3w    , pe3uw   , pe3vw               )   !       3D vertical scale factors at w-level
       ENDIF
       !
    END SUBROUTINE usr_def_zgr
@@ -203,10 +208,10 @@ CONTAINS
    END SUBROUTINE zgr_msk_top_bot
    
 
-   SUBROUTINE zgr_zco( pdept_1d, pdepw_1d, pe3t_1d, pe3w_1d,   &   ! in : 1D reference vertical coordinate
-      &                pdept   , pdepw   ,                     &   ! out: 3D t & w-points depth
-      &                pe3t    , pe3u    , pe3v   , pe3f   ,   &   !      vertical scale factors
-      &                pe3w    , pe3uw   , pe3vw             )     !          -      -      -
+   SUBROUTINE zgr_zco_3d( pdept_1d, pdepw_1d, pe3t_1d, pe3w_1d,   &   ! in : 1D reference vertical coordinate
+      &                   pe3t    , pe3u    , pe3v   , pe3f   ,   &   ! out: 3D vertical scale factors at t-level
+      &                   pdept   , pdepw   ,                     &   !      3D t & w-points depth
+      &                   pe3w    , pe3uw   , pe3vw             )     !      3D vertical scale factors at w-level
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE zgr_zco  ***
       !!
@@ -216,13 +221,18 @@ CONTAINS
       !!----------------------------------------------------------------------
       REAL(wp), DIMENSION(:)    , INTENT(in   ) ::   pdept_1d, pdepw_1d          ! 1D grid-point depth       [m]
       REAL(wp), DIMENSION(:)    , INTENT(in   ) ::   pe3t_1d , pe3w_1d           ! 1D vertical scale factors [m]
-      REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pdept, pdepw                ! grid-point depth          [m]
-      REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pe3t , pe3u , pe3v , pe3f   ! vertical scale factors    [m]
-      REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pe3w , pe3uw, pe3vw         !    -       -      -
+      REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pe3t , pe3u , pe3v , pe3f   ! 3D vertical scale factors [m]
+      REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pdept, pdepw                ! 3D grid-point depth       [m]
+      REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pe3w , pe3uw, pe3vw         ! 3D vertical scale factors [m]
       !
       INTEGER  ::   jk
       !!----------------------------------------------------------------------
       !
+      IF(lwp) WRITE(numout,*)
+      IF(lwp) WRITE(numout,*) '    zgr_zco_3d : defines depths and scale-factors.'
+      IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~'
+      IF(lwp) WRITE(numout,*) '       GYRE case : uniform'
+      !
       DO jk = 1, jpk
          pdept(:,:,jk) = pdept_1d(jk)
          pdepw(:,:,jk) = pdepw_1d(jk)
@@ -235,7 +245,37 @@ CONTAINS
          pe3vw(:,:,jk) = pe3w_1d (jk)
       END DO
       !
-   END SUBROUTINE zgr_zco
+   END SUBROUTINE zgr_zco_3d
+
+   
+   SUBROUTINE zgr_zco_1d3d( pe3t_1d,                          &   ! in : 1D reference vertical coordinate
+      &                     pe3t   , pe3u   , pe3v  , pe3f  )     ! out: vertical scale factors
+      !!----------------------------------------------------------------------
+      !!                  ***  ROUTINE zgr_zco  ***
+      !!
+      !! ** Purpose :   define the reference z-coordinate system
+      !!
+      !! ** Method  :   set 3D coord. arrays to reference 1D array 
+      !!----------------------------------------------------------------------
+      REAL(wp), DIMENSION(:)    , INTENT(in   ) ::   pe3t_1d                     ! 1D vertical scale factors [m]
+      REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pe3t , pe3u , pe3v , pe3f   ! 3D vertical scale factors [m]
+      !
+      INTEGER  ::   jk
+      !!----------------------------------------------------------------------
+      !
+      IF(lwp) WRITE(numout,*)
+      IF(lwp) WRITE(numout,*) '    zgr_zco_1d3d : defines t-level scale-factors.'
+      IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~~~'
+      IF(lwp) WRITE(numout,*) '       GYRE case : uniform'
+      !
+      DO jk = 1, jpk
+         pe3t(:,:,jk) = pe3t_1d (jk)
+         pe3u(:,:,jk) = pe3t_1d (jk)
+         pe3v(:,:,jk) = pe3t_1d (jk)
+         pe3f(:,:,jk) = pe3t_1d (jk)
+      END DO
+      !
+   END SUBROUTINE zgr_zco_1d3d
 
    !!======================================================================
 END MODULE usrdef_zgr
diff --git a/tests/BENCH/MY_SRC/usrdef_zgr.F90 b/tests/BENCH/MY_SRC/usrdef_zgr.F90
index 1f601c21..c95eb6df 100644
--- a/tests/BENCH/MY_SRC/usrdef_zgr.F90
+++ b/tests/BENCH/MY_SRC/usrdef_zgr.F90
@@ -39,20 +39,20 @@ CONTAINS
    SUBROUTINE usr_def_zgr( ld_zco  , ld_zps  , ld_sco  , ld_isfcav,    &   ! type of vertical coordinate
       &                    k_top   , k_bot                        ,    &   ! top & bottom ocean level
       &                    pdept_1d, pdepw_1d, pe3t_1d , pe3w_1d  ,    &   ! 1D reference vertical coordinate
-      &                    pe3t  , pe3u  , pe3v   , pe3f ,             &   ! vertical scale factors
+      &                    pe3t  , pe3u  , pe3v   , pe3f ,             &   ! 3D vertical scale factors at t-level
       &                    pdept , pdepw ,                             &   ! 3D t & w-points depth
-      &                    pe3w  , pe3uw , pe3vw                       )   ! vertical scale factors
+      &                    pe3w  , pe3uw , pe3vw                       )   ! 3D vertical scale factors at w-level
       !!---------------------------------------------------------------------
       !!              ***  ROUTINE usr_def_zgr  ***
       !!
       !! ** Purpose :   User defined the vertical coordinates
       !!
       !!----------------------------------------------------------------------
-      LOGICAL                   , INTENT(out) ::   ld_zco, ld_zps, ld_sco      ! vertical coordinate flags
-      LOGICAL                   , INTENT(out) ::   ld_isfcav                   ! under iceshelf cavity flag
-      INTEGER , DIMENSION(:,:)  , INTENT(out) ::   k_top, k_bot                ! first & last ocean level
-      REAL(wp), DIMENSION(:)    , INTENT(out) ::   pdept_1d, pdepw_1d          ! 1D grid-point depth     [m]
-      REAL(wp), DIMENSION(:)    , INTENT(out) ::   pe3t_1d , pe3w_1d           ! 1D grid-point depth     [m]
+      LOGICAL                             , INTENT(out) ::   ld_zco, ld_zps, ld_sco      ! vertical coordinate flags
+      LOGICAL                             , INTENT(out) ::   ld_isfcav                   ! under iceshelf cavity flag
+      INTEGER , DIMENSION(:,:)            , INTENT(out) ::   k_top, k_bot                ! first & last ocean level
+      REAL(wp), DIMENSION(:)              , INTENT(out) ::   pdept_1d, pdepw_1d          ! 1D grid-point depth     [m]
+      REAL(wp), DIMENSION(:)              , INTENT(out) ::   pe3t_1d , pe3w_1d           ! 1D grid-point depth     [m]
       REAL(wp), DIMENSION(:,:,:), OPTIONAL, INTENT(out) ::   pdept, pdepw                ! grid-point depth        [m]
       REAL(wp), DIMENSION(:,:,:), OPTIONAL, INTENT(out) ::   pe3t , pe3u , pe3v , pe3f   ! vertical scale factors  [m]
       REAL(wp), DIMENSION(:,:,:), OPTIONAL, INTENT(out) ::   pe3w , pe3uw, pe3vw         ! i-scale factors 
@@ -77,11 +77,16 @@ CONTAINS
       !
       CALL zgr_msk_top_bot( k_top , k_bot )                 ! masked top and bottom ocean t-level indices
       !
-      IF( PRESENT( pe3t ) ) THEN                            ! z-coordinate (3D arrays) from the 1D z-coord.
-         CALL zgr_zco( pdept_1d, pdepw_1d, pe3t_1d, pe3w_1d,   &   ! in  : 1D reference vertical coordinate
-            &          pdept   , pdepw   ,                     &   ! out : 3D t & w-points depth
-            &          pe3t    , pe3u    , pe3v   , pe3f   ,   &   !       vertical scale factors
-            &          pe3w    , pe3uw   , pe3vw             )     !           -      -      -
+      IF( lk_vco_1d3d ) THEN                                ! z-coordinate (3D arrays) from the 1D z-coord.
+         CALL zgr_zco_1d3d( pe3t_1d,                   &               ! in  : 1D reference vertical coordinate
+            &               pe3t , pe3u , pe3v , pe3f  )               ! out : 3D vertical scale factors at t-level
+      ENDIF
+      !
+      IF( lk_vco_3d ) THEN                                  ! z-coordinate (3D arrays) from the 1D z-coord.
+         CALL zgr_zco_3d(   pdept_1d, pdepw_1d, pe3t_1d, pe3w_1d,   &   ! in  : 1D reference vertical coordinate
+            &               pe3t    , pe3u    , pe3v   , pe3f   ,   &   ! out : 3D vertical scale factors at t-level
+            &               pdept   , pdepw   ,                     &   !       3D t & w-points depth
+            &               pe3w    , pe3uw   , pe3vw               )   !       3D vertical scale factors at w-level
       ENDIF
       !
    END SUBROUTINE usr_def_zgr
@@ -213,10 +218,10 @@ CONTAINS
    END SUBROUTINE zgr_msk_top_bot
    
 
-   SUBROUTINE zgr_zco( pdept_1d, pdepw_1d, pe3t_1d, pe3w_1d,   &   ! in : 1D reference vertical coordinate
-      &                pdept   , pdepw   ,                     &   ! out: 3D t & w-points depth
-      &                pe3t    , pe3u    , pe3v   , pe3f   ,   &   !      vertical scale factors
-      &                pe3w    , pe3uw   , pe3vw             )     !          -      -      -
+   SUBROUTINE zgr_zco_3d( pdept_1d, pdepw_1d, pe3t_1d, pe3w_1d,   &   ! in : 1D reference vertical coordinate
+      &                   pe3t    , pe3u    , pe3v   , pe3f   ,   &   ! out: 3D vertical scale factors at t-level
+      &                   pdept   , pdepw   ,                     &   !      3D t & w-points depth
+      &                   pe3w    , pe3uw   , pe3vw             )     !      3D vertical scale factors at w-level
       !!----------------------------------------------------------------------
       !!                  ***  ROUTINE zgr_zco  ***
       !!
@@ -226,13 +231,18 @@ CONTAINS
       !!----------------------------------------------------------------------
       REAL(wp), DIMENSION(:)    , INTENT(in   ) ::   pdept_1d, pdepw_1d          ! 1D grid-point depth       [m]
       REAL(wp), DIMENSION(:)    , INTENT(in   ) ::   pe3t_1d , pe3w_1d           ! 1D vertical scale factors [m]
-      REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pdept, pdepw                ! grid-point depth          [m]
-      REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pe3t , pe3u , pe3v , pe3f   ! vertical scale factors    [m]
-      REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pe3w , pe3uw, pe3vw         !    -       -      -
+      REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pe3t , pe3u , pe3v , pe3f   ! 3D vertical scale factors [m]
+      REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pdept, pdepw                ! 3D grid-point depth       [m]
+      REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pe3w , pe3uw, pe3vw         ! 3D vertical scale factors [m]
       !
       INTEGER  ::   jk
       !!----------------------------------------------------------------------
       !
+      IF(lwp) WRITE(numout,*)
+      IF(lwp) WRITE(numout,*) '    zgr_zco_3d : defines depths and scale-factors.'
+      IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~'
+      IF(lwp) WRITE(numout,*) '       BENCH case : uniform'
+      !
       DO jk = 1, jpk
          pdept(:,:,jk) = pdept_1d(jk)
          pdepw(:,:,jk) = pdepw_1d(jk)
@@ -245,7 +255,37 @@ CONTAINS
          pe3vw(:,:,jk) = pe3w_1d (jk)
       END DO
       !
-   END SUBROUTINE zgr_zco
+   END SUBROUTINE zgr_zco_3d
+
+   
+   SUBROUTINE zgr_zco_1d3d( pe3t_1d,                          &   ! in : 1D reference vertical coordinate
+      &                     pe3t   , pe3u   , pe3v  , pe3f  )     ! out: vertical scale factors
+      !!----------------------------------------------------------------------
+      !!                  ***  ROUTINE zgr_zco  ***
+      !!
+      !! ** Purpose :   define the reference z-coordinate system
+      !!
+      !! ** Method  :   set 3D coord. arrays to reference 1D array 
+      !!----------------------------------------------------------------------
+      REAL(wp), DIMENSION(:)    , INTENT(in   ) ::   pe3t_1d                     ! 1D vertical scale factors [m]
+      REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pe3t , pe3u , pe3v , pe3f   ! 3D vertical scale factors [m]
+      !
+      INTEGER  ::   jk
+      !!----------------------------------------------------------------------
+      !
+      IF(lwp) WRITE(numout,*)
+      IF(lwp) WRITE(numout,*) '    zgr_zco_1d3d : defines t-level scale-factors.'
+      IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~~~'
+      IF(lwp) WRITE(numout,*) '       BENCH case : uniform'
+      !
+      DO jk = 1, jpk
+         pe3t(:,:,jk) = pe3t_1d (jk)
+         pe3u(:,:,jk) = pe3t_1d (jk)
+         pe3v(:,:,jk) = pe3t_1d (jk)
+         pe3f(:,:,jk) = pe3t_1d (jk)
+      END DO
+      !
+   END SUBROUTINE zgr_zco_1d3d
 
    !!======================================================================
 END MODULE usrdef_zgr
-- 
GitLab