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

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
usrdef_zgr.F90 in NEMO/trunk/tests/BENCH/MY_SRC – NEMO

source: NEMO/trunk/tests/BENCH/MY_SRC/usrdef_zgr.F90 @ 14433

Last change on this file since 14433 was 14433, checked in by smasson, 3 years ago

trunk: merge dev_r14312_MPI_Interface into the trunk, #2598

File size: 13.0 KB
Line 
1MODULE usrdef_zgr
2   !!======================================================================
3   !!                       ***  MODULE  usrdef_zgr  ***
4   !!
5   !!                      ===  BENCH configuration  ===
6   !!
7   !! User defined : vertical coordinate system of a user configuration
8   !!======================================================================
9   !! History :  4.0  !
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 phycst         ! physical constants
21   USE depth_e3       ! depth <=> e3
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   !!----------------------------------------------------------------------
33   !! NEMO/OPA 4.0 , NEMO Consortium (2016)
34   !! $Id$
35   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
36   !!----------------------------------------------------------------------
37CONTAINS             
38
39   SUBROUTINE usr_def_zgr( ld_zco  , ld_zps  , ld_sco  , ld_isfcav,    &   ! type of vertical coordinate
40      &                    pdept_1d, pdepw_1d, pe3t_1d , pe3w_1d  ,    &   ! 1D reference vertical coordinate
41      &                    pdept , pdepw ,                             &   ! 3D t & w-points depth
42      &                    pe3t  , pe3u  , pe3v   , pe3f ,             &   ! vertical scale factors
43      &                    pe3w  , pe3uw , pe3vw         ,             &   !     -      -      -
44      &                    k_top  , k_bot    )                             ! top & bottom ocean level
45      !!---------------------------------------------------------------------
46      !!              ***  ROUTINE usr_def_zgr  ***
47      !!
48      !! ** Purpose :   User defined the vertical coordinates
49      !!
50      !!----------------------------------------------------------------------
51      LOGICAL                   , INTENT(out) ::   ld_zco, ld_zps, ld_sco      ! vertical coordinate flags
52      LOGICAL                   , INTENT(out) ::   ld_isfcav                   ! under iceshelf cavity flag
53      REAL(wp), DIMENSION(:)    , INTENT(out) ::   pdept_1d, pdepw_1d          ! 1D grid-point depth     [m]
54      REAL(wp), DIMENSION(:)    , INTENT(out) ::   pe3t_1d , pe3w_1d           ! 1D grid-point depth     [m]
55      REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   pdept, pdepw                ! grid-point depth        [m]
56      REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   pe3t , pe3u , pe3v , pe3f   ! vertical scale factors  [m]
57      REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   pe3w , pe3uw, pe3vw         ! i-scale factors
58      INTEGER , DIMENSION(:,:)  , INTENT(out) ::   k_top, k_bot                ! first & last ocean level
59      !
60      INTEGER  ::   inum   ! local logical unit
61      REAL(WP) ::   z_zco, z_zps, z_sco, z_cav
62      REAL(wp), DIMENSION(jpi,jpj) ::   z2d   ! 2D workspace
63      !!----------------------------------------------------------------------
64      !
65      IF(lwp) WRITE(numout,*)
66      IF(lwp) WRITE(numout,*) 'usr_def_zgr : BENCH configuration (z-coordinate closed flat box ocean)'
67      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
68      !
69      !
70      ! type of vertical coordinate
71      ! ---------------------------
72      ld_zco    = .TRUE.         ! BENCH case:  z-coordinate without ocean cavities
73      ld_zps    = .FALSE.
74      ld_sco    = .FALSE.
75      ld_isfcav = .FALSE.
76      !
77      !
78      ! Build the vertical coordinate system
79      ! ------------------------------------
80      CALL zgr_z( pdept_1d, pdepw_1d, pe3t_1d , pe3w_1d )   ! Reference z-coordinate system
81      !
82      CALL zgr_msk_top_bot( k_top , k_bot )                 ! masked top and bottom ocean t-level indices
83      !
84      !                                                     ! z-coordinate (3D arrays) from the 1D z-coord.
85      CALL zgr_zco( pdept_1d, pdepw_1d, pe3t_1d, pe3w_1d,   &   ! in  : 1D reference vertical coordinate
86         &          pdept   , pdepw   ,                     &   ! out : 3D t & w-points depth
87         &          pe3t    , pe3u    , pe3v   , pe3f   ,   &   !       vertical scale factors
88         &          pe3w    , pe3uw   , pe3vw             )     !           -      -      -
89      !
90   END SUBROUTINE usr_def_zgr
91
92
93   SUBROUTINE zgr_z( pdept_1d, pdepw_1d, pe3t_1d , pe3w_1d )   ! 1D reference vertical coordinate
94      !!----------------------------------------------------------------------
95      !!                   ***  ROUTINE zgr_z  ***
96      !!
97      !! ** Purpose :   set the 1D depth of model levels and the resulting
98      !!              vertical scale factors.
99      !!
100      !! ** Method  :   1D z-coordinate system (use in all type of coordinate)
101      !!       The depth of model levels is set from dep(k), an analytical function:
102      !!                   w-level: depw_1d  = dep(k)
103      !!                   t-level: dept_1d  = dep(k+0.5)
104      !!       The scale factors are the discrete derivative of the depth:
105      !!                   e3w_1d(jk) = dk[ dept_1d ]
106      !!                   e3t_1d(jk) = dk[ depw_1d ]
107      !!           with at top and bottom :
108      !!                   e3w_1d( 1 ) = 2 * ( dept_1d( 1 ) - depw_1d( 1 ) )
109      !!                   e3t_1d(jpk) = 2 * ( dept_1d(jpk) - depw_1d(jpk) )
110      !!       The depth are then re-computed from the sum of e3. This ensures
111      !!    that depths are identical when reading domain configuration file.
112      !!    Indeed, only e3. are saved in this file, depth are compute by a call
113      !!    to the e3_to_depth subroutine.
114      !!
115      !!       Here the Madec & Imbard (1996) function is used.
116      !!
117      !! ** Action  : - pdept_1d, pdepw_1d : depth of T- and W-point (m)
118      !!              - pe3t_1d , pe3w_1d  : scale factors at T- and W-levels (m)
119      !!
120      !! Reference : Marti, Madec & Delecluse, 1992, JGR, 97, No8, 12,763-12,766.
121      !!             Madec and Imbard, 1996, Clim. Dyn.
122      !!----------------------------------------------------------------------
123      REAL(wp), DIMENSION(:)    , INTENT(out) ::   pdept_1d, pdepw_1d   ! 1D grid-point depth        [m]
124      REAL(wp), DIMENSION(:)    , INTENT(out) ::   pe3t_1d , pe3w_1d    ! 1D vertical scale factors  [m]
125      !
126      INTEGER  ::   jk       ! dummy loop indices
127      REAL(wp) ::   zd       ! local scalar
128      !!----------------------------------------------------------------------
129      !
130      zd = 5000./FLOAT(jpkm1)
131      !
132      IF(lwp) THEN            ! Parameter print
133         WRITE(numout,*)
134         WRITE(numout,*) '    zgr_z   : Reference vertical z-coordinates '
135         WRITE(numout,*) '    ~~~~~~~'
136         WRITE(numout,*) '       BENCH case : uniform vertical grid :'
137         WRITE(numout,*) '                     with thickness = ', zd
138      ENDIF
139
140      !
141      ! 1D Reference z-coordinate    (using Madec & Imbard 1996 function)
142      ! -------------------------
143      !
144      pdepw_1d(1) = 0._wp
145      pdept_1d(1) = 0.5_wp * zd
146      !
147      DO jk = 2, jpk          ! depth at T and W-points
148         pdepw_1d(jk) = pdepw_1d(jk-1) + zd 
149         pdept_1d(jk) = pdept_1d(jk-1) + zd 
150      END DO
151      !
152      !                       ! e3t and e3w from depth
153      CALL depth_to_e3( pdept_1d, pdepw_1d, pe3t_1d, pe3w_1d ) 
154      !
155      !                       ! recompute depths from SUM(e3)  <== needed
156      CALL e3_to_depth( pe3t_1d, pe3w_1d, pdept_1d, pdepw_1d ) 
157      !
158      IF(lwp) THEN                        ! control print
159         WRITE(numout,*)
160         WRITE(numout,*) '              Reference 1D z-coordinate depth and scale factors:'
161         WRITE(numout, "(9x,' level  gdept_1d  gdepw_1d  e3t_1d   e3w_1d  ')" )
162         WRITE(numout, "(10x, i4, 4f9.2)" ) ( jk, pdept_1d(jk), pdepw_1d(jk), pe3t_1d(jk), pe3w_1d(jk), jk = 1, jpk )
163      ENDIF
164      !
165   END SUBROUTINE zgr_z
166
167
168   SUBROUTINE zgr_msk_top_bot( k_top , k_bot )
169      !!----------------------------------------------------------------------
170      !!                    ***  ROUTINE zgr_msk_top_bot  ***
171      !!
172      !! ** Purpose :   set the masked top and bottom ocean t-levels
173      !!
174      !! ** Method  :   BENCH case = closed flat box ocean without ocean cavities
175      !!                   k_top = 1     except along north, south, east and west boundaries
176      !!                   k_bot = jpk-1 except along north, south, east and west boundaries
177      !!
178      !! ** Action  : - k_top : first wet ocean level index
179      !!              - k_bot : last  wet ocean level index
180      !!----------------------------------------------------------------------
181      INTEGER , DIMENSION(:,:), INTENT(out) ::   k_top , k_bot   ! first & last wet ocean level
182      !
183      REAL(wp), DIMENSION(jpi,jpj) ::   z2d   ! 2D local workspace
184      REAL(wp)                     ::   zmaxlam, zscl
185      !!----------------------------------------------------------------------
186      !
187      IF(lwp) WRITE(numout,*)
188      IF(lwp) WRITE(numout,*) '    zgr_top_bot : defines the top and bottom wet ocean levels.'
189      IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~'
190      IF(lwp) WRITE(numout,*) '       BENCH case : closed flat box ocean without ocean cavities'
191      !
192      z2d(:,:) = REAL( jpkm1 , wp )          ! flat bottom
193      !
194      !
195      ! BENCH should work without these 2 small islands on the 2 poles of the folding...
196      !   -> Comment out these lines if instabilities are too large...
197      !
198     
199!!$      IF( c_NFtype == 'T' ) THEN   ! add a small island in the upper corners to avoid model instabilities...
200!!$         z2d(mi0(       nn_hls):mi1(                  nn_hls+2 ),mj0(jpjglo-nn_hls-1):mj1(jpjglo-nn_hls+1)) = 0.
201!!$         z2d(mi0(jpiglo-nn_hls):mi1(MIN(jpiglo,jpiglo-nn_hls+2)),mj0(jpjglo-nn_hls-1):mj1(jpjglo-nn_hls+1)) = 0.
202!!$         z2d(mi0(jpiglo/2     ):mi1(           jpiglo/2     +2 ),mj0(jpjglo-nn_hls-1):mj1(jpjglo-nn_hls+1)) = 0.
203!!$      ENDIF
204!!$      !
205!!$      IF( c_NFtype == 'F' ) THEN   ! add a small island in the upper corners to avoid model instabilities...
206!!$         z2d(mi0(       nn_hls):mi1(       nn_hls+1),mj0(jpjglo-nn_hls):mj1(jpjglo-nn_hls+1)) = 0.
207!!$         z2d(mi0(jpiglo-nn_hls):mi1(jpiglo-nn_hls+1),mj0(jpjglo-nn_hls):mj1(jpjglo-nn_hls+1)) = 0.
208!!$         z2d(mi0(jpiglo/2     ):mi1(jpiglo/2     +1),mj0(jpjglo-nn_hls):mj1(jpjglo-nn_hls+1)) = 0.
209!!$      ENDIF
210
211      !
212      CALL lbc_lnk( 'usrdef_zgr', z2d, 'T', 1. )           ! set surrounding land to zero (closed boundaries)
213      !
214      k_bot(:,:) = INT( z2d(:,:) )           ! =jpkm1 over the ocean point, =0 elsewhere
215      !
216      k_top(:,:) = MIN( 1 , k_bot(:,:) )     ! = 1    over the ocean point, =0 elsewhere
217      !
218   END SUBROUTINE zgr_msk_top_bot
219   
220
221   SUBROUTINE zgr_zco( pdept_1d, pdepw_1d, pe3t_1d, pe3w_1d,   &   ! in : 1D reference vertical coordinate
222      &                pdept   , pdepw   ,                     &   ! out: 3D t & w-points depth
223      &                pe3t    , pe3u    , pe3v   , pe3f   ,   &   !      vertical scale factors
224      &                pe3w    , pe3uw   , pe3vw             )     !          -      -      -
225      !!----------------------------------------------------------------------
226      !!                  ***  ROUTINE zgr_zco  ***
227      !!
228      !! ** Purpose :   define the reference z-coordinate system
229      !!
230      !! ** Method  :   set 3D coord. arrays to reference 1D array
231      !!----------------------------------------------------------------------
232      REAL(wp), DIMENSION(:)    , INTENT(in   ) ::   pdept_1d, pdepw_1d          ! 1D grid-point depth       [m]
233      REAL(wp), DIMENSION(:)    , INTENT(in   ) ::   pe3t_1d , pe3w_1d           ! 1D vertical scale factors [m]
234      REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pdept, pdepw                ! grid-point depth          [m]
235      REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pe3t , pe3u , pe3v , pe3f   ! vertical scale factors    [m]
236      REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pe3w , pe3uw, pe3vw         !    -       -      -
237      !
238      INTEGER  ::   jk
239      !!----------------------------------------------------------------------
240      !
241      DO jk = 1, jpk
242         pdept(:,:,jk) = pdept_1d(jk)
243         pdepw(:,:,jk) = pdepw_1d(jk)
244         pe3t (:,:,jk) = pe3t_1d (jk)
245         pe3u (:,:,jk) = pe3t_1d (jk)
246         pe3v (:,:,jk) = pe3t_1d (jk)
247         pe3f (:,:,jk) = pe3t_1d (jk)
248         pe3w (:,:,jk) = pe3w_1d (jk)
249         pe3uw(:,:,jk) = pe3w_1d (jk)
250         pe3vw(:,:,jk) = pe3w_1d (jk)
251      END DO
252      !
253   END SUBROUTINE zgr_zco
254
255   !!======================================================================
256END MODULE usrdef_zgr
Note: See TracBrowser for help on using the repository browser.