source: NEMO/trunk/cfgs/C1D_PAPA/MY_SRC/usrdef_zgr.F90 @ 9799

Last change on this file since 9799 was 9799, checked in by cbricaud, 2 years ago

Restore C1D_PAPA configuration, #2061 and #2087

File size: 9.1 KB
Line 
1MODULE usrdef_zgr
2   !!======================================================================
3   !!                       ***  MODULE  usrdef_zgr  ***
4   !!
5   !!                       ===  C1D_PAPA configuration  ===
6   !!
7   !! User defined : vertical coordinate system of a user configuration
8   !!======================================================================
9   !! History :  4.0  ! 2016-06  (R. Bourdalle-Badie)  Original code
10   !!----------------------------------------------------------------------
11
12   !!----------------------------------------------------------------------
13   !!   usr_def_zgr   : user defined vertical coordinate system
14   !!      zgr_z      : reference 1D z-coordinate
15   !!      zgr_top_bot: ocean top and bottom level indices
16   !!      zgr_zco    : 3D verticl coordinate in pure z-coordinate case
17   !!---------------------------------------------------------------------
18   USE oce            ! ocean variables
19   USE dom_oce        ! ocean domain
20   USE depth_e3       ! depth <=> e3
21   USE usrdef_nam     ! User defined : namelist variables
22   !
23   USE in_out_manager ! I/O manager
24   USE lbclnk         ! ocean lateral boundary conditions (or mpp link)
25   USE lib_mpp        ! distributed memory computing library
26
27   IMPLICIT NONE
28   PRIVATE
29
30   PUBLIC   usr_def_zgr        ! called by domzgr.F90
31
32   !! * Substitutions
33#  include "vectopt_loop_substitute.h90"
34   !!----------------------------------------------------------------------
35   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
36   !! $Id:$
37   !! Software governed by the CeCILL licence     (./LICENSE)
38   !!----------------------------------------------------------------------
39CONTAINS             
40
41   SUBROUTINE usr_def_zgr( ld_zco  , ld_zps  , ld_sco  , ld_isfcav,    &   ! type of vertical coordinate
42      &                    pdept_1d, pdepw_1d, pe3t_1d , pe3w_1d  ,    &   ! 1D reference vertical coordinate
43      &                    pdept , pdepw ,                             &   ! 3D t & w-points depth
44      &                    pe3t  , pe3u  , pe3v   , pe3f ,             &   ! vertical scale factors
45      &                    pe3w  , pe3uw , pe3vw         ,             &   !     -      -      -
46      &                    k_top  , k_bot    )                             ! top & bottom ocean level
47      !!---------------------------------------------------------------------
48      !!              ***  ROUTINE usr_def_zgr  ***
49      !!
50      !! ** Purpose :   User defined the vertical coordinates
51      !!
52      !!----------------------------------------------------------------------
53      LOGICAL                   , INTENT(out) ::   ld_zco, ld_zps, ld_sco      ! vertical coordinate flags
54      LOGICAL                   , INTENT(out) ::   ld_isfcav                   ! under iceshelf cavity flag
55      REAL(wp), DIMENSION(:)    , INTENT(out) ::   pdept_1d, pdepw_1d          ! 1D grid-point depth     [m]
56      REAL(wp), DIMENSION(:)    , INTENT(out) ::   pe3t_1d , pe3w_1d           ! 1D grid-point depth     [m]
57      REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   pdept, pdepw                ! grid-point depth        [m]
58      REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   pe3t , pe3u , pe3v , pe3f   ! vertical scale factors  [m]
59      REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   pe3w , pe3uw, pe3vw         ! i-scale factors
60      INTEGER , DIMENSION(:,:)  , INTENT(out) ::   k_top, k_bot                ! first & last ocean level
61      !
62      INTEGER  ::   ji, jj, jk        ! dummy indices
63      INTEGER  ::   ik                ! local integers
64      REAL(wp) ::   zfact, z1_jpkm1   ! local scalar
65      REAL(wp) ::   ze3min            ! local scalar
66      REAL(wp) ::   zt, zw            ! local scalars
67      REAL(wp) ::   zsur, za0, za1, zkth, zacr        ! Values for the Madec & Imbard (1996) function
68      REAL(wp) ::   za2, zkth2, zacr2                 ! Values for optional double tanh function set from parameters
69      REAL(wp), DIMENSION(jpi,jpj) ::   zht, zhu, z2d ! 2D workspace
70
71      !!----------------------------------------------------------------------
72      !
73      IF(lwp) WRITE(numout,*)
74      IF(lwp) WRITE(numout,*) 'usr_def_zgr : C1D configuration (zps-coordinate closed box ocean without cavities)'
75      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
76      !
77      ! type of vertical coordinate
78      ! ---------------------------
79      ld_zco    = .FALSE.         ! C1D case:  z-coordinate without ocean cavities
80      ld_zps    = .TRUE.
81      ld_sco    = .FALSE.
82      ld_isfcav = .FALSE.
83      !
84      ! Build the vertical coordinate system
85      ! ------------------------------------
86      !
87      ! Set parameters of z(k) function
88      ! -------------------------------
89      zsur =   -3958.95137127683
90      za0  =    103.953009600000
91      za1  =    2.41595126900000
92      zkth =    15.3510137000000
93      zacr =    7.00000000000000
94      za2  =    100.760928500000
95      zkth2=    48.0298937200000
96      zacr2=    13.0000000000000
97      !
98      IF(lwp) THEN            ! Parameter print
99         WRITE(numout,*)
100         WRITE(numout,*) '     zgr_z75L   : Reference vertical z-coordinates '
101         WRITE(numout,*) '     ~~~~~~~'
102         WRITE(numout,*) '       C1D case : L75 function with the following coefficients :'
103         WRITE(numout,*) '                 zsur = ', zsur
104         WRITE(numout,*) '                 za0  = ', za0
105         WRITE(numout,*) '                 za1  = ', za1
106         WRITE(numout,*) '                 zkth = ', zkth
107         WRITE(numout,*) '                 zacr = ', zacr
108         WRITE(numout,*) '                 za2  = ', za2
109         WRITE(numout,*) '                 zkth2= ', zkth2
110         WRITE(numout,*) '                 zacr2= ', zacr2
111      ENDIF
112
113      !                       !==  UNmasked meter bathymetry  ==!
114      !
115      zht(:,:) = rn_bathy
116      !
117      DO jk = 1, jpk          ! depth at T and W-points
118         zw = REAL( jk , wp )
119         zt = REAL( jk , wp ) + 0.5_wp
120         pdepw_1d(jk) = ( zsur + za0 * zw + za1 * zacr * LOG ( COSH( (zw-zkth ) / zacr  ) )    &
121                  &                    + za2 * zacr2* LOG ( COSH( (zw-zkth2) / zacr2 ) )  )
122         pdept_1d(jk) = ( zsur + za0 * zt + za1 * zacr * LOG ( COSH( (zt-zkth ) / zacr  ) )    &
123                  &                    + za2 * zacr2* LOG ( COSH( (zt-zkth2) / zacr2 ) )  )
124      END DO
125      !
126      !                       ! e3t and e3w from depth
127      CALL depth_to_e3( pdept_1d, pdepw_1d, pe3t_1d, pe3w_1d )
128      !
129      !                       ! recompute depths from SUM(e3)  <== needed
130      CALL e3_to_depth( pe3t_1d, pe3w_1d, pdept_1d, pdepw_1d )
131      !
132      IF(lwp) THEN                        ! control print
133         WRITE(numout,*)
134         WRITE(numout,*) '              Reference 1D z-coordinate depth and scale factors:'
135         WRITE(numout, "(9x,' level  gdept_1d  gdepw_1d  e3t_1d   e3w_1d  ')" )
136         WRITE(numout, "(10x, i4, 4f9.2)" ) ( jk, pdept_1d(jk), pdepw_1d(jk), pe3t_1d(jk), pe3w_1d(jk), jk = 1, jpk )
137      ENDIF
138      !
139      !                       !==  top masked level bathymetry  ==!  (all coordinates)
140      !
141      k_top(:,:) = 1 
142      !                                   !* bottom ocean compute from the depth of grid-points
143      k_bot(:,:) = jpkm1
144      DO jk = jpkm1, 1, -1
145        ze3min = 0.1_wp * pe3t_1d (jk)
146         WHERE( zht(:,:) < pdepw_1d(jk) + ze3min )   k_bot(:,:) = jk-1
147      END DO
148      !
149      !                                !* vertical coordinate system
150      DO jk = 1, jpk                      ! initialization to the reference z-coordinate
151         pdept(:,:,jk) = pdept_1d(jk)
152         pdepw(:,:,jk) = pdepw_1d(jk)
153         pe3t (:,:,jk) = pe3t_1d (jk)
154         pe3u (:,:,jk) = pe3t_1d (jk)
155         pe3v (:,:,jk) = pe3t_1d (jk)
156         pe3f (:,:,jk) = pe3t_1d (jk)
157         pe3w (:,:,jk) = pe3w_1d (jk)
158         pe3uw(:,:,jk) = pe3w_1d (jk)
159         pe3vw(:,:,jk) = pe3w_1d (jk)
160      END DO
161      DO jj = 1, jpj                      ! bottom scale factors and depth at T- and W-points
162         DO ji = 1, jpi
163            ik = k_bot(ji,jj)
164            pdepw(ji,jj,ik+1) = MIN( zht(ji,jj) , pdepw_1d(ik+1) )
165            pe3t (ji,jj,ik  ) = pdepw(ji,jj,ik+1) - pdepw(ji,jj,ik)
166            pe3t (ji,jj,ik+1) = pe3t (ji,jj,ik  ) 
167            !
168            pdept(ji,jj,ik  ) = pdepw(ji,jj,ik  ) + pe3t (ji,jj,ik  ) * 0.5_wp
169            pdept(ji,jj,ik+1) = pdepw(ji,jj,ik+1) + pe3t (ji,jj,ik+1) * 0.5_wp
170            pe3w (ji,jj,ik+1) = pdept(ji,jj,ik+1) - pdept(ji,jj,ik)              ! = pe3t (ji,jj,ik  )
171         END DO
172      END DO         
173      !                                   ! bottom scale factors and depth at  U-, V-, UW and VW-points
174      !                                   ! usually Computed as the minimum of neighbooring scale factors
175      pe3u (:,:,:) = pe3t(:,:,:)          ! HERE C1D configuration :
176      pe3v (:,:,:) = pe3t(:,:,:)          !    e3 increases with k-index
177      pe3f (:,:,:) = pe3t(:,:,:)          !    so e3 minimum of (i,i+1) points is (i) point
178      pe3uw(:,:,:) = pe3w(:,:,:)          !    in j-direction e3v=e3t and e3f=e3v
179      pe3vw(:,:,:) = pe3w(:,:,:)          !    ==>>  no need of lbc_lnk calls
180      !     
181      !
182   END SUBROUTINE usr_def_zgr
183
184   !!======================================================================
185END MODULE usrdef_zgr
Note: See TracBrowser for help on using the repository browser.