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/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/CONFIG/OVERFLOW/MY_SRC – NEMO

source: branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/CONFIG/OVERFLOW/MY_SRC/usrdef_zgr.F90 @ 6906

Last change on this file since 6906 was 6906, checked in by flavoni, 8 years ago

#1692 and ROBUST-3 : commit zgr and some update for LOCK EXCHANGE case

File size: 12.8 KB
Line 
1MODULE usrdef_zgr
2   !!==============================================================================
3   !!                       ***  MODULE usrdef_zgr  ***
4   !! Ocean domain : user defined vertical coordinate system
5   !!
6   !!                       ===      OVERFLOW case      ===
7   !!
8   !!==============================================================================
9   !! History :  4.0  ! 2016-08  (G. Madec)  Original code
10   !!----------------------------------------------------------------------
11
12   !!----------------------------------------------------------------------
13   !!   usr_def_zgr      : user defined vertical coordinate system (required)
14   !!       zgr_z1d      : reference 1D z-coordinate
15   !!       zgr_zps      : 3D vertical coordinate in z-partial cell coordinate
16   !!---------------------------------------------------------------------
17   USE oce               ! ocean variables
18   USE dom_oce  ,  ONLY: ln_zco, ln_zps, ln_sco   ! ocean space and time domain
19   USE dom_oce  ,  ONLY: nimpp, njmpp             ! ocean space and time domain
20   USE dom_oce  ,  ONLY: glamt                    ! ocean space and time domain
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   USE wrk_nemo          ! Memory allocation
27   USE timing            ! Timing
28
29   IMPLICIT NONE
30   PRIVATE
31
32   PUBLIC   usr_def_zgr        ! called by domzgr.F90
33
34  !! * Substitutions
35#  include "vectopt_loop_substitute.h90"
36   !!----------------------------------------------------------------------
37   !! NEMO/OPA 4.0 , NEMO Consortium (2016)
38   !! $Id: domzgr.F90 6624 2016-05-26 08:59:48Z gm $
39   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
40   !!----------------------------------------------------------------------
41CONTAINS             
42
43   SUBROUTINE usr_def_zgr( ld_zco  , ld_zps  , ld_sco  , ld_isfcav,    &   ! type of vertical coordinate
44      &                    pdept_1d, pdepw_1d, pe3t_1d , pe3w_1d  ,    &   ! 1D reference vertical coordinate
45      &                    pdept , pdepw ,                             &   ! 3D t & w-points depth
46      &                    pe3t  , pe3u  , pe3v , pe3f ,               &   ! vertical scale factors
47      &                    pe3w  , pe3uw , pe3vw,                      &   !     -      -      -
48      &                    k_top  , k_bot    )                             ! top & bottom ocean level
49      !!---------------------------------------------------------------------
50      !!              ***  ROUTINE usr_def_zgr  ***
51      !!
52      !! ** Purpose :   User defined the vertical coordinates
53      !!
54      !!----------------------------------------------------------------------
55      LOGICAL                   , INTENT(out) ::   ld_zco, ld_zps, ld_sco      ! vertical coordinate flags
56      LOGICAL                   , INTENT(out) ::   ld_isfcav                   ! under iceshelf cavity flag
57      REAL(wp), DIMENSION(:)    , INTENT(out) ::   pdept_1d, pdepw_1d          ! 1D grid-point depth     [m]
58      REAL(wp), DIMENSION(:)    , INTENT(out) ::   pe3t_1d , pe3w_1d           ! 1D grid-point depth     [m]
59      REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   pdept, pdepw                ! grid-point depth        [m]
60      REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   pe3t , pe3u , pe3v , pe3f   ! vertical scale factors  [m]
61      REAL(wp), DIMENSION(:,:,:), INTENT(out) ::   pe3w , pe3uw, pe3vw         ! i-scale factors
62      INTEGER , DIMENSION(:,:)  , INTENT(out) ::   k_top, k_bot                ! first & last ocean level
63      !
64      INTEGER  ::   ji, jj, jk        ! dummy indices
65      INTEGER  ::   ik                ! local integers
66      REAL(wp) ::   zfact, z1_jpkm1   ! local scalar
67      REAL(wp) ::   ze3min            ! local scalar
68      REAL(wp), DIMENSION(jpi,jpj) ::   zht, zhu, z2d   ! 2D workspace
69      !!----------------------------------------------------------------------
70      !
71      IF(lwp) WRITE(numout,*)
72      IF(lwp) WRITE(numout,*) 'usr_def_zgr : OVERFLOW configuration (z-coordinate closed box ocean without cavities)'
73      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
74      !
75      !
76      ! type of vertical coordinate
77      ! ---------------------------
78      ! already set in usrdef_nam.F90 by reading the namusr_def namelist except for ISF
79      ld_isfcav = .FALSE.
80      !
81      !
82      ! Build the vertical coordinate system
83      ! ------------------------------------
84      !
85      !                       !==  UNmasked meter bathymetry  ==!
86      !
87      ! western continental shelf (500m deep) and eastern deep ocean (2000m deep)
88      ! with a hyperbolic tangent transition centered at 40km
89      ! NB: here glamt is in kilometers
90      !
91      zht(:,:) = + (  500. + 0.5 * 1500. * ( 1.0 + tanh( (glamt(:,:) - 40.) / 7. ) )  )
92      !
93      ! at u-point: recompute from glamu (see usrdef_hgr.F90) to avoid the masking when using lbc_lnk
94      zfact = rn_dx * 1.e-3         ! conversion in km
95      DO ji = 1, jpi
96         zhu(ji,:) = + (  500. + 0.5 * 1500. * ( 1.0 + tanh( ( zfact * REAL( ji-1 + nimpp-1 , wp ) - 40.) / 7. ) )  )
97      END DO
98      !
99      CALL zgr_z1d( pdept_1d, pdepw_1d, pe3t_1d , pe3w_1d )   ! Reference z-coordinate system
100      !
101      !
102      !                       !==  top masked level bathymetry  ==!  (all coordinates)
103      !
104      ! no ocean cavities : top ocean level is ONE, except over land
105      ! the ocean basin surrounded by land (1 grid-point) set through lbc_lnk call as jperio=0
106      z2d(:,:) = 1._wp                    ! surface ocean is the 1st level
107      CALL lbc_lnk( z2d, 'T', 1. )        ! closed basin
108      k_top(:,:) = z2d(:,:)
109      !
110      !                             
111      !
112      IF ( ln_sco ) THEN      !==  s-coordinate  ==!   (terrain-following coordinate)
113         !
114         k_bot(:,:) = jpkm1 * k_top(:,:)  !* bottom ocean = jpk-1 (here use k_top as a land mask)
115         !
116         !                                !* terrain-following coordinate with e3.(k)=cst)
117         !                                !  OVERFLOW case : identical with j-index (T=V, U=F)
118         z1_jpkm1 = 1._wp / REAL( jpkm1 , wp)
119         DO jk = 1, jpk
120            pdept(:,:,jk) = zht(:,:) * z1_jpkm1 * ( REAL( jk   , wp ) - 0.5_wp )
121            pdepw(:,:,jk) = zht(:,:) * z1_jpkm1 * ( REAL( jk-1 , wp )          )
122            pe3t (:,:,jk) = zht(:,:) * z1_jpkm1
123            pe3u (:,:,jk) = zhu(:,:) * z1_jpkm1
124            pe3v (:,:,jk) = zht(:,:) * z1_jpkm1
125            pe3f (:,:,jk) = zhu(:,:) * z1_jpkm1
126            pe3w (:,:,jk) = zht(:,:) * z1_jpkm1
127            pe3uw(:,:,jk) = zhu(:,:) * z1_jpkm1
128            pe3vw(:,:,jk) = zht(:,:) * z1_jpkm1
129         END DO     
130      ENDIF
131      !
132      !
133      IF ( ln_zco ) THEN      !==  z-coordinate  ==!   (step-like topography)
134         !
135         !                                !* bottom ocean compute from the depth of grid-points
136         k_bot(:,:) = jpkm1 * k_top(:,:)     ! here use k_top as a land mask
137         DO jk = 1, jpkm1
138            WHERE( pdept_1d(jk) < zht(:,:) .AND. zht(:,:) <= pdept_1d(jk+1) )   k_bot(:,:) = jk * k_top(:,:)
139         END DO
140         !                                !* horizontally uniform coordinate (reference z-co everywhere)
141         DO jk = 1, jpk
142            pdept(:,:,jk) = pdept_1d(jk)
143            pdepw(:,:,jk) = pdepw_1d(jk)
144            pe3t (:,:,jk) = pe3t_1d (jk)
145            pe3u (:,:,jk) = pe3t_1d (jk)
146            pe3v (:,:,jk) = pe3t_1d (jk)
147            pe3f (:,:,jk) = pe3t_1d (jk)
148            pe3w (:,:,jk) = pe3w_1d (jk)
149            pe3uw(:,:,jk) = pe3w_1d (jk)
150            pe3vw(:,:,jk) = pe3w_1d (jk)
151         END DO
152      ENDIF
153      !
154      !
155      IF ( ln_zps ) THEN      !==  zps-coordinate  ==!   (partial bottom-steps)
156         !
157         ze3min = 0.1_wp * rn_dz
158         IF(lwp) WRITE(numout,*) '   minimum thickness of the partial cells = 10 % of e3 = ', ze3min
159         !
160         !
161         !                                !* bottom ocean compute from the depth of grid-points
162         k_bot(:,:) = jpkm1
163         DO jk = jpkm1, 1, -1
164            WHERE( zht(:,:) <= pdepw_1d(jk) + ze3min )   k_bot(:,:) = jk-1
165         END DO
166         !
167         !                                !* vertical coordinate system
168         DO jk = 1, jpk                      ! initialization to the reference z-coordinate
169            pdept(:,:,jk) = pdept_1d(jk)
170            pdepw(:,:,jk) = pdepw_1d(jk)
171            pe3t (:,:,jk) = pe3t_1d (jk)
172            pe3u (:,:,jk) = pe3t_1d (jk)
173            pe3v (:,:,jk) = pe3t_1d (jk)
174            pe3f (:,:,jk) = pe3t_1d (jk)
175            pe3w (:,:,jk) = pe3w_1d (jk)
176            pe3uw(:,:,jk) = pe3w_1d (jk)
177            pe3vw(:,:,jk) = pe3w_1d (jk)
178         END DO
179         DO jj = 1, jpj                      ! bottom scale factors and depth at T- and W-points
180            DO ji = 1, jpi
181               ik = k_bot(ji,jj)
182               IF( ik /= jpkm1 ) THEN                 ! last level ==> ref 1d z-coordinate
183                  pdepw(ji,jj,ik+1) = MIN( zht(ji,jj) , pdepw_1d(ik+1) )
184                  pe3t (ji,jj,ik  ) = pdepw(ji,jj,ik+1) - pdepw(ji,jj,ik)
185                  pe3t (ji,jj,ik+1) = pe3t (ji,jj,ik) 
186                  !
187                  pdept(ji,jj,ik  ) = pdept(ji,jj,ik-1) + pe3t(ji,jj,ik) * 0.5_wp
188                  pe3w (ji,jj,ik+1) = pdepw(ji,jj,ik+1) - pdepw(ji,jj,ik)
189               ENDIF
190            END DO
191         END DO         
192         !                                   ! bottom scale factors and depth at  U-, V-, UW and VW-points
193         !                                   ! usually Computed as the minimum of neighbooring scale factors
194         pe3u (:,:,:) = pe3t(:,:,:)          ! HERE OVERFLOW configuration :
195         pe3v (:,:,:) = pe3t(:,:,:)          !    e3 increases with i-index and identical with j-index
196         pe3f (:,:,:) = pe3t(:,:,:)          !    so e3 minimum of (i,i+1) points is (i) point (idem in j-direction)
197         pe3uw(:,:,:) = pe3w(:,:,:)          !   
198         pe3vw(:,:,:) = pe3w(:,:,:)          !    ==>>  no need of lbc_lnk calls
199         !     
200      ENDIF
201      !
202   END SUBROUTINE usr_def_zgr
203
204
205   SUBROUTINE zgr_z1d( pdept_1d, pdepw_1d, pe3t_1d , pe3w_1d )   ! 1D reference vertical coordinate
206      !!----------------------------------------------------------------------
207      !!                   ***  ROUTINE zgr_z1d  ***
208      !!
209      !! ** Purpose :   set the depth of model levels and the resulting
210      !!      vertical scale factors.
211      !!
212      !! ** Method  :   z-coordinate system (use in all type of coordinate)
213      !!      The depth of model levels is defined from an analytical
214      !!      function the derivative of which gives the scale factors.
215      !!      both depth and scale factors only depend on k (1d arrays).
216      !!              w-level: pdepw_1d  = pdep(k)
217      !!                       pe3w_1d(k) = dk(pdep)(k)     = e3(k)
218      !!              t-level: pdept_1d  = pdep(k+0.5)
219      !!                       pe3t_1d(k) = dk(pdep)(k+0.5) = e3(k+0.5)
220      !!
221      !!            ===    Here constant vertical resolution   ===
222      !!
223      !! ** Action  : - pdept_1d, pdepw_1d : depth of T- and W-point (m)
224      !!              - pe3t_1d , pe3w_1d  : scale factors at T- and W-levels (m)
225      !!----------------------------------------------------------------------
226      REAL(wp), DIMENSION(:), INTENT(out) ::   pdept_1d, pdepw_1d   ! 1D grid-point depth        [m]
227      REAL(wp), DIMENSION(:), INTENT(out) ::   pe3t_1d , pe3w_1d    ! 1D vertical scale factors  [m]
228      !
229      INTEGER  ::   jk       ! dummy loop indices
230      REAL(wp) ::   zt, zw   ! local scalar
231      !!----------------------------------------------------------------------
232      !
233      IF(lwp) THEN                         ! Parameter print
234         WRITE(numout,*)
235         WRITE(numout,*) '    zgr_z1d : Reference vertical z-coordinates: uniform dz = ', rn_dz
236         WRITE(numout,*) '    ~~~~~~~'
237      ENDIF
238      !
239      ! Reference z-coordinate (depth - scale factor at T- and W-points)   ! Madec & Imbard 1996 function
240      ! ----------------------
241      DO jk = 1, jpk
242         zw = REAL( jk , wp )
243         zt = REAL( jk , wp ) + 0.5_wp
244         pdepw_1d(jk) =    rn_dz *   REAL( jk-1 , wp )
245         pdept_1d(jk) =    rn_dz * ( REAL( jk-1 , wp ) + 0.5_wp )
246         pe3w_1d (jk) =    rn_dz
247         pe3t_1d (jk) =    rn_dz
248      END DO
249      !
250      IF(lwp) THEN                        ! control print
251         WRITE(numout,*)
252         WRITE(numout,*) '              Reference 1D z-coordinate depth and scale factors:'
253         WRITE(numout, "(9x,' level  gdept_1d  gdepw_1d  e3t_1d   e3w_1d  ')" )
254         WRITE(numout, "(10x, i4, 4f9.2)" ) ( jk, pdept_1d(jk), pdepw_1d(jk), pe3t_1d(jk), pe3w_1d(jk), jk = 1, jpk )
255      ENDIF
256      !
257   END SUBROUTINE zgr_z1d
258   
259   !!======================================================================
260END MODULE usrdef_zgr
Note: See TracBrowser for help on using the repository browser.