source: NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/tests/CANAL/MY_SRC/usrdef_zgr.F90 @ 10297

Last change on this file since 10297 was 10297, checked in by smasson, 2 years ago

dev_r10164_HPC09_ESIWACE_PREP_MERGE: action 2a: add report calls of mppmin/max/sum, see #2133

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