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 branches/2017/dev_r7881_no_wrk_alloc/NEMOGCM/NEMO/OPA_SRC/USR – NEMO

source: branches/2017/dev_r7881_no_wrk_alloc/NEMOGCM/NEMO/OPA_SRC/USR/usrdef_zgr.F90 @ 8016

Last change on this file since 8016 was 8016, checked in by timgraham, 7 years ago

Delete some remaining "USE wrk_array" lines
File size: 12.8 KB
Line 
1MODULE usrdef_zgr
2   !!======================================================================
3   !!                       ***  MODULE  usrdef_zgr  ***
4   !!
5   !!                       ===  GYRE configuration  ===
6   !!
7   !! User defined : vertical coordinate system of a user configuration
8   !!======================================================================
9   !! History :  4.0  ! 2016-06  (G. Madec)  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   !
22   USE in_out_manager ! I/O manager
23   USE lbclnk         ! ocean lateral boundary conditions (or mpp link)
24   USE lib_mpp        ! distributed memory computing library
25   USE timing         ! Timing
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/OPA 4.0 , NEMO Consortium (2016)
36   !! $Id$
37   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
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  ::   inum   ! local logical unit
63      REAL(WP) ::   z_zco, z_zps, z_sco, z_cav
64      REAL(wp), DIMENSION(jpi,jpj) ::   z2d   ! 2D workspace
65      !!----------------------------------------------------------------------
66      !
67      IF(lwp) WRITE(numout,*)
68      IF(lwp) WRITE(numout,*) 'usr_def_zgr : GYRE configuration (z-coordinate closed flat box ocean without cavities)'
69      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
70      !
71      !
72      ! type of vertical coordinate
73      ! ---------------------------
74      ld_zco    = .TRUE.         ! GYRE case:  z-coordinate without ocean cavities
75      ld_zps    = .FALSE.
76      ld_sco    = .FALSE.
77      ld_isfcav = .FALSE.
78      !
79      !
80      ! Build the vertical coordinate system
81      ! ------------------------------------
82      CALL zgr_z( pdept_1d, pdepw_1d, pe3t_1d , pe3w_1d )   ! Reference z-coordinate system
83      !
84      CALL zgr_msk_top_bot( k_top , k_bot )                 ! masked top and bottom ocean t-level indices
85      !
86      !                                                     ! z-coordinate (3D arrays) from the 1D z-coord.
87      CALL zgr_zco( pdept_1d, pdepw_1d, pe3t_1d, pe3w_1d,   &   ! in  : 1D reference vertical coordinate
88         &          pdept   , pdepw   ,                     &   ! out : 3D t & w-points depth
89         &          pe3t    , pe3u    , pe3v   , pe3f   ,   &   !       vertical scale factors
90         &          pe3w    , pe3uw   , pe3vw             )     !           -      -      -
91      !
92   END SUBROUTINE usr_def_zgr
93
94
95   SUBROUTINE zgr_z( pdept_1d, pdepw_1d, pe3t_1d , pe3w_1d )   ! 1D reference vertical coordinate
96      !!----------------------------------------------------------------------
97      !!                   ***  ROUTINE zgr_z  ***
98      !!
99      !! ** Purpose :   set the 1D depth of model levels and the resulting
100      !!              vertical scale factors.
101      !!
102      !! ** Method  :   1D z-coordinate system (use in all type of coordinate)
103      !!       The depth of model levels is set from dep(k), an analytical function:
104      !!                   w-level: depw_1d  = dep(k)
105      !!                   t-level: dept_1d  = dep(k+0.5)
106      !!       The scale factors are the discrete derivative of the depth:
107      !!                   e3w_1d(jk) = dk[ dept_1d ]
108      !!                   e3t_1d(jk) = dk[ depw_1d ]
109      !!           with at top and bottom :
110      !!                   e3w_1d( 1 ) = 2 * ( dept_1d( 1 ) - depw_1d( 1 ) )
111      !!                   e3t_1d(jpk) = 2 * ( dept_1d(jpk) - depw_1d(jpk) )
112      !!       The depth are then re-computed from the sum of e3. This ensures
113      !!    that depths are identical when reading domain configuration file.
114      !!    Indeed, only e3. are saved in this file, depth are compute by a call
115      !!    to the e3_to_depth subroutine.
116      !!
117      !!       Here the Madec & Imbard (1996) function is used.
118      !!
119      !! ** Action  : - pdept_1d, pdepw_1d : depth of T- and W-point (m)
120      !!              - pe3t_1d , pe3w_1d  : scale factors at T- and W-levels (m)
121      !!
122      !! Reference : Marti, Madec & Delecluse, 1992, JGR, 97, No8, 12,763-12,766.
123      !!             Madec and Imbard, 1996, Clim. Dyn.
124      !!----------------------------------------------------------------------
125      REAL(wp), DIMENSION(:)    , INTENT(out) ::   pdept_1d, pdepw_1d   ! 1D grid-point depth        [m]
126      REAL(wp), DIMENSION(:)    , INTENT(out) ::   pe3t_1d , pe3w_1d    ! 1D vertical scale factors  [m]
127      !
128      INTEGER  ::   jk       ! dummy loop indices
129      REAL(wp) ::   zt, zw   ! local scalars
130      REAL(wp) ::   zsur, za0, za1, zkth, zacr   ! Values for the Madec & Imbard (1996) function 
131      !!----------------------------------------------------------------------
132      !
133      IF( nn_timing == 1 )  CALL timing_start('zgr_z')
134      !
135      ! Set parameters of z(k) function
136      ! -------------------------------
137      zsur = -2033.194295283385_wp       
138      za0  =   155.8325369664153_wp 
139      za1  =   146.3615918601890_wp
140      zkth =    17.28520372419791_wp   
141      zacr =     5.0_wp       
142      !
143      IF(lwp) THEN            ! Parameter print
144         WRITE(numout,*)
145         WRITE(numout,*) '    zgr_z   : Reference vertical z-coordinates '
146         WRITE(numout,*) '    ~~~~~~~'
147         WRITE(numout,*) '       GYRE case : MI96 function with the following coefficients :'
148         WRITE(numout,*) '                 zsur = ', zsur
149         WRITE(numout,*) '                 za0  = ', za0
150         WRITE(numout,*) '                 za1  = ', za1
151         WRITE(numout,*) '                 zkth = ', zkth
152         WRITE(numout,*) '                 zacr = ', zacr
153      ENDIF
154
155      !
156      ! 1D Reference z-coordinate    (using Madec & Imbard 1996 function)
157      ! -------------------------
158      !
159      DO jk = 1, jpk          ! depth at T and W-points
160         zw = REAL( jk , wp )
161         zt = REAL( jk , wp ) + 0.5_wp
162         pdepw_1d(jk) = ( zsur + za0 * zw + za1 * zacr *  LOG( COSH( (zw-zkth) / zacr ) )  )
163         pdept_1d(jk) = ( zsur + za0 * zt + za1 * zacr *  LOG( COSH( (zt-zkth) / zacr ) )  )
164      END DO
165      !
166      !                       ! e3t and e3w from depth
167      CALL depth_to_e3( pdept_1d, pdepw_1d, pe3t_1d, pe3w_1d ) 
168      !
169      !                       ! recompute depths from SUM(e3)  <== needed
170      CALL e3_to_depth( pe3t_1d, pe3w_1d, pdept_1d, pdepw_1d ) 
171      !
172      IF(lwp) THEN                        ! control print
173         WRITE(numout,*)
174         WRITE(numout,*) '              Reference 1D z-coordinate depth and scale factors:'
175         WRITE(numout, "(9x,' level  gdept_1d  gdepw_1d  e3t_1d   e3w_1d  ')" )
176         WRITE(numout, "(10x, i4, 4f9.2)" ) ( jk, pdept_1d(jk), pdepw_1d(jk), pe3t_1d(jk), pe3w_1d(jk), jk = 1, jpk )
177      ENDIF
178      !
179      IF( nn_timing == 1 )  CALL timing_stop('zgr_z')
180      !
181   END SUBROUTINE zgr_z
182
183
184   SUBROUTINE zgr_msk_top_bot( k_top , k_bot )
185      !!----------------------------------------------------------------------
186      !!                    ***  ROUTINE zgr_msk_top_bot  ***
187      !!
188      !! ** Purpose :   set the masked top and bottom ocean t-levels
189      !!
190      !! ** Method  :   GYRE case = closed flat box ocean without ocean cavities
191      !!                   k_top = 1     except along north, south, east and west boundaries
192      !!                   k_bot = jpk-1 except along north, south, east and west boundaries
193      !!
194      !! ** Action  : - k_top : first wet ocean level index
195      !!              - k_bot : last  wet ocean level index
196      !!----------------------------------------------------------------------
197      INTEGER , DIMENSION(:,:), INTENT(out) ::   k_top , k_bot   ! first & last wet ocean level
198      !
199      REAL(wp), DIMENSION(jpi,jpj) ::   z2d   ! 2D local workspace
200      !!----------------------------------------------------------------------
201      !
202      IF(lwp) WRITE(numout,*)
203      IF(lwp) WRITE(numout,*) '    zgr_top_bot : defines the top and bottom wet ocean levels.'
204      IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~'
205      IF(lwp) WRITE(numout,*) '       GYRE case : closed flat box ocean without ocean cavities'
206      !
207      z2d(:,:) = REAL( jpkm1 , wp )          ! flat bottom
208      !
209      CALL lbc_lnk( z2d, 'T', 1. )           ! set surrounding land to zero (here jperio=0 ==>> closed)
210      !
211      k_bot(:,:) = INT( z2d(:,:) )           ! =jpkm1 over the ocean point, =0 elsewhere
212      !
213      k_top(:,:) = MIN( 1 , k_bot(:,:) )     ! = 1    over the ocean point, =0 elsewhere
214      !
215   END SUBROUTINE zgr_msk_top_bot
216   
217
218   SUBROUTINE zgr_zco( pdept_1d, pdepw_1d, pe3t_1d, pe3w_1d,   &   ! in : 1D reference vertical coordinate
219      &                pdept   , pdepw   ,                     &   ! out: 3D t & w-points depth
220      &                pe3t    , pe3u    , pe3v   , pe3f   ,   &   !      vertical scale factors
221      &                pe3w    , pe3uw   , pe3vw             )     !          -      -      -
222      !!----------------------------------------------------------------------
223      !!                  ***  ROUTINE zgr_zco  ***
224      !!
225      !! ** Purpose :   define the reference z-coordinate system
226      !!
227      !! ** Method  :   set 3D coord. arrays to reference 1D array
228      !!----------------------------------------------------------------------
229      REAL(wp), DIMENSION(:)    , INTENT(in   ) ::   pdept_1d, pdepw_1d          ! 1D grid-point depth       [m]
230      REAL(wp), DIMENSION(:)    , INTENT(in   ) ::   pe3t_1d , pe3w_1d           ! 1D vertical scale factors [m]
231      REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pdept, pdepw                ! grid-point depth          [m]
232      REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pe3t , pe3u , pe3v , pe3f   ! vertical scale factors    [m]
233      REAL(wp), DIMENSION(:,:,:), INTENT(  out) ::   pe3w , pe3uw, pe3vw         !    -       -      -
234      !
235      INTEGER  ::   jk
236      !!----------------------------------------------------------------------
237      !
238      IF( nn_timing == 1 )  CALL timing_start('zgr_zco')
239      !
240      DO jk = 1, jpk
241         pdept(:,:,jk) = pdept_1d(jk)
242         pdepw(:,:,jk) = pdepw_1d(jk)
243         pe3t (:,:,jk) = pe3t_1d (jk)
244         pe3u (:,:,jk) = pe3t_1d (jk)
245         pe3v (:,:,jk) = pe3t_1d (jk)
246         pe3f (:,:,jk) = pe3t_1d (jk)
247         pe3w (:,:,jk) = pe3w_1d (jk)
248         pe3uw(:,:,jk) = pe3w_1d (jk)
249         pe3vw(:,:,jk) = pe3w_1d (jk)
250      END DO
251      !
252      IF( nn_timing == 1 )  CALL timing_stop('zgr_zco')
253      !
254   END SUBROUTINE zgr_zco
255
256   !!======================================================================
257END MODULE usrdef_zgr
Note: See TracBrowser for help on using the repository browser.