source: NEMO/trunk/tests/ISOMIP/MY_SRC/usrdef_zgr.F90 @ 12377

Last change on this file since 12377 was 12377, checked in by acc, 7 months ago

The big one. Merging all 2019 developments from the option 1 branch back onto the trunk.

This changeset reproduces 2019/dev_r11943_MERGE_2019 on the trunk using a 2-URL merge
onto a working copy of the trunk. I.e.:

svn merge —ignore-ancestry \

svn+ssh://acc@forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/NEMO/trunk \
svn+ssh://acc@forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/NEMO/branches/2019/dev_r11943_MERGE_2019 ./

The —ignore-ancestry flag avoids problems that may otherwise arise from the fact that
the merge history been trunk and branch may have been applied in a different order but
care has been taken before this step to ensure that all applicable fixes and updates
are present in the merge branch.

The trunk state just before this step has been branched to releases/release-4.0-HEAD
and that branch has been immediately tagged as releases/release-4.0.2. Any fixes
or additions in response to tickets on 4.0, 4.0.1 or 4.0.2 should be done on
releases/release-4.0-HEAD. From now on future 'point' releases (e.g. 4.0.2) will
remain unchanged with periodic releases as needs demand. Note release-4.0-HEAD is a
transitional naming convention. Future full releases, say 4.2, will have a release-4.2
branch which fulfills this role and the first point release (e.g. 4.2.0) will be made
immediately following the release branch creation.

2020 developments can be started from any trunk revision later than this one.

  • Property svn:keywords set to Id
File size: 12.4 KB
Line 
1MODULE usrdef_zgr
2   !!======================================================================
3   !!                     ***  MODULE usrdef_zgr  ***
4   !!
5   !!                       ===  ISOMIP case  ===
6   !!
7   !! user defined :  vertical coordinate system of a user configuration
8   !!======================================================================
9   !! History :  4.0  ! 2016-08  (G. Madec,   S. Flavoni)  Original code
10   !!                 ! 2017-02  (P. Mathiot, S. Flavoni)  Adapt code to ISOMIP case
11   !!----------------------------------------------------------------------
12
13   !!----------------------------------------------------------------------
14   !!   usr_def_zgr   : user defined vertical coordinate system (required)
15   !!       zgr_z1d   : reference 1D z-coordinate
16   !!---------------------------------------------------------------------
17   USE oce            ! ocean variables
18   USE dom_oce ,  ONLY: mj0   , mj1   , nimpp , njmpp   ! ocean space and time domain
19   USE dom_oce ,  ONLY: glamt , gphit                   ! ocean space and time domain
20   USE usrdef_nam     ! User defined : namelist variables
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   !!----------------------------------------------------------------------
33   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
34   !! $Id$
35   !! Software governed by the CeCILL license (see ./LICENSE)
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(in   ) ::   ld_zco, ld_zps, ld_sco      ! vertical coordinate flags ( read in namusr_def )
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  ::   ji , jj, jk       ! dummy indices
61      INTEGER  ::   ij0, ij1          ! dummy indices 
62      INTEGER  ::   ik                ! local integers
63      REAL(wp) ::   zfact, z1_jpkm1   ! local scalar
64      REAL(wp) ::   ze3min, zdepth    ! local scalar
65      REAL(wp), DIMENSION(jpi,jpj) ::   zht  , zhu         ! bottom depth
66      REAL(wp), DIMENSION(jpi,jpj) ::   zhisf, zhisfu      ! top depth
67      REAL(wp), DIMENSION(jpi,jpj) ::   zmsk 
68      REAL(wp), DIMENSION(jpi,jpj) ::   z2d                ! 2d workspace
69      !!----------------------------------------------------------------------
70      !
71      IF(lwp) WRITE(numout,*)
72      IF(lwp) WRITE(numout,*) 'usr_def_zgr : ISOMIP configuration (z(ps)- or s-coordinate closed box ocean without cavities)'
73      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~'
74      !
75      !
76      ! type of vertical coordinate
77      ! ---------------------------
78      ! set in usrdef_nam.F90 by reading the namusr_def namelist except for ISF
79      ld_isfcav = .TRUE.       ! ISF Ice Shelves Flag
80      !
81      !
82      ! Build the vertical coordinate system
83      ! ------------------------------------
84      !
85      !                       !==  isfdraft  ==!
86      !
87      ! the ocean basin surrounded by land (1 grid-point) set through lbc_lnk call as jperio=0
88      z2d(:,:) = 1._wp                    ! surface ocean is the 1st level
89      CALL lbc_lnk( 'usrdef_zgr', z2d, 'T', 1. )        ! closed basin since jperio = 0 (see userdef_nam.F90)
90      zmsk(:,:) = NINT( z2d(:,:) )
91      !
92      !
93      zht  (:,:) = rbathy 
94      zhisf(:,:) = 200._wp
95      ij0 = 1 ; ij1 = 40
96      DO jj = mj0(ij0), mj1(ij1)
97         zhisf(:,jj)=700.0_wp-(gphit(:,jj)+80.0_wp)*125.0_wp
98      END DO
99      zhisf(:,:) = zhisf(:,:) * zmsk(:,:)
100      !
101      CALL zgr_z1d( pdept_1d, pdepw_1d, pe3t_1d , pe3w_1d )   ! Reference z-coordinate system
102      !
103      !                       !==  top masked level bathymetry  ==!  (all coordinates)
104      !
105      IF ( ld_zps ) THEN      !==  zps-coordinate  ==!   (partial bottom-steps)
106         !
107         ze3min = 0.1_wp * rn_e3
108         IF(lwp) WRITE(numout,*) '   minimum thickness of the partial cells = 10 % of e3 = ', ze3min
109         !
110         !                                !* bottom ocean compute from the depth of grid-points
111         k_bot(:,:) = jpkm1
112         DO jk = jpkm1, 1, -1
113            WHERE( zht(:,:) < pdepw_1d(jk) + ze3min )   k_bot(:,:) = jk-1
114         END DO
115         !                                !* top ocean compute from the depth of grid-points
116         k_top(:,:) = 1                   !
117         DO jk = 2, jpkm1
118            zdepth = pdepw_1d(jk+1) - ze3min
119            WHERE( zhisf(:,:) > 0.0 .AND. zhisf(:,:) >= zdepth )   k_top(:,:) = (jk + 1) 
120         END DO
121         !
122         !                                   !* vertical coordinate system
123         DO jk = 1, jpk                      ! initialization to the reference z-coordinate
124            pdept(:,:,jk) = pdept_1d(jk)
125            pdepw(:,:,jk) = pdepw_1d(jk)
126            pe3t (:,:,jk) = pe3t_1d (jk)
127            pe3u (:,:,jk) = pe3t_1d (jk)
128            pe3v (:,:,jk) = pe3t_1d (jk)
129            pe3f (:,:,jk) = pe3t_1d (jk)
130            pe3w (:,:,jk) = pe3w_1d (jk)
131            pe3uw(:,:,jk) = pe3w_1d (jk)
132            pe3vw(:,:,jk) = pe3w_1d (jk)
133         END DO
134         DO jj = 1, jpj                      ! top scale factors and depth at T- and W-points
135            DO ji = 1, jpi
136               ik = k_top(ji,jj)
137               IF ( ik > 2 ) THEN
138                  ! pdeptw at the interface
139                  pdepw(ji,jj,ik  ) = MAX( zhisf(ji,jj) , pdepw(ji,jj,ik) )
140                  ! e3t in both side of the interface
141                  pe3t (ji,jj,ik  ) = pdepw(ji,jj,ik+1) - pdepw(ji,jj,ik)
142                  ! pdept in both side of the interface (from previous e3t)
143                  pdept(ji,jj,ik  ) = pdepw(ji,jj,ik  ) + pe3t (ji,jj,ik  ) * 0.5_wp
144                  pdept(ji,jj,ik-1) = pdepw(ji,jj,ik  ) - pe3t (ji,jj,ik  ) * 0.5_wp
145                  ! pe3w on both side of the interface
146                  pe3w (ji,jj,ik+1) = pdept(ji,jj,ik+1) - pdept(ji,jj,ik  )
147                  pe3w (ji,jj,ik  ) = pdept(ji,jj,ik  ) - pdept(ji,jj,ik-1)
148                  ! e3t into the ice shelf
149                  pe3t (ji,jj,ik-1) = pdepw(ji,jj,ik  ) - pdepw(ji,jj,ik-1)
150                  pe3w (ji,jj,ik-1) = pdept(ji,jj,ik-1) - pdept(ji,jj,ik-2)
151               END IF
152            END DO
153         END DO         
154         DO jj = 1, jpj                      ! bottom scale factors and depth at T- and W-points
155            DO ji = 1, jpi
156               ik = k_bot(ji,jj)
157               pdepw(ji,jj,ik+1) = MIN( zht(ji,jj) , pdepw_1d(ik+1) )
158               pe3t (ji,jj,ik  ) = pdepw(ji,jj,ik+1) - pdepw(ji,jj,ik)
159               pe3t (ji,jj,ik+1) = pe3t (ji,jj,ik  ) 
160               !
161               pdept(ji,jj,ik  ) = pdepw(ji,jj,ik  ) + pe3t (ji,jj,ik  ) * 0.5_wp
162               pdept(ji,jj,ik+1) = pdepw(ji,jj,ik+1) + pe3t (ji,jj,ik+1) * 0.5_wp
163               pe3w (ji,jj,ik+1) = pdept(ji,jj,ik+1) - pdept(ji,jj,ik)
164            END DO
165         END DO         
166         !                                   ! bottom scale factors and depth at  U-, V-, UW and VW-points
167         pe3u (:,:,:) = pe3t(:,:,:)
168         pe3uw(:,:,:) = pe3w(:,:,:)
169         DO jk = 1, jpk                      ! Computed as the minimum of neighbooring scale factors
170            DO jj = 1, jpjm1
171               DO ji = 1, jpi
172                  pe3v (ji,jj,jk) = MIN( pe3t(ji,jj,jk), pe3t(ji,jj+1,jk) )
173                  pe3vw(ji,jj,jk) = MIN( pe3w(ji,jj,jk), pe3w(ji,jj+1,jk) )
174                  pe3f (ji,jj,jk) = pe3v(ji,jj,jk)
175               END DO
176            END DO
177         END DO
178         CALL lbc_lnk( 'usrdef_zgr', pe3v , 'V', 1._wp )   ;   CALL lbc_lnk( 'usrdef_zgr', pe3vw, 'V', 1._wp )
179         CALL lbc_lnk( 'usrdef_zgr', pe3f , 'F', 1._wp )
180         DO jk = 1, jpk
181            ! set to z-scale factor if zero (i.e. along closed boundaries) because of lbclnk
182            WHERE( pe3u (:,:,jk) == 0._wp )   pe3u (:,:,jk) = pe3t_1d(jk)
183            WHERE( pe3v (:,:,jk) == 0._wp )   pe3v (:,:,jk) = pe3t_1d(jk)
184            WHERE( pe3f (:,:,jk) == 0._wp )   pe3f (:,:,jk) = pe3t_1d(jk)
185            WHERE( pe3uw(:,:,jk) == 0._wp )   pe3uw(:,:,jk) = pe3w_1d(jk)
186            WHERE( pe3vw(:,:,jk) == 0._wp )   pe3vw(:,:,jk) = pe3w_1d(jk)
187         END DO
188         !
189      ENDIF
190      !
191   END SUBROUTINE usr_def_zgr
192
193
194   SUBROUTINE zgr_z1d( pdept_1d, pdepw_1d, pe3t_1d , pe3w_1d )   ! 1D reference vertical coordinate
195      !!----------------------------------------------------------------------
196      !!                   ***  ROUTINE zgr_z1d  ***
197      !!
198      !! ** Purpose :   set the depth of model levels and the resulting
199      !!      vertical scale factors.
200      !!
201      !! ** Method  :   1D z-coordinate system (use in all type of coordinate)
202      !!       The depth of model levels is set from dep(k), an analytical function:
203      !!                   w-level: depw_1d  = dep(k)
204      !!                   t-level: dept_1d  = dep(k+0.5)
205      !!       The scale factors are the discrete derivative of the depth:
206      !!                   e3w_1d(jk) = dk[ dept_1d ]
207      !!                   e3t_1d(jk) = dk[ depw_1d ]
208      !!
209      !!            ===    Here constant vertical resolution   ===
210      !!
211      !! ** Action  : - pdept_1d, pdepw_1d : depth of T- and W-point (m)
212      !!              - pe3t_1d , pe3w_1d  : scale factors at T- and W-levels (m)
213      !!----------------------------------------------------------------------
214      REAL(wp), DIMENSION(:), INTENT(out) ::   pdept_1d, pdepw_1d   ! 1D grid-point depth        [m]
215      REAL(wp), DIMENSION(:), INTENT(out) ::   pe3t_1d , pe3w_1d    ! 1D vertical scale factors  [m]
216      !
217      INTEGER  ::   jk       ! dummy loop indices
218      REAL(wp) ::   zt, zw   ! local scalar
219      !!----------------------------------------------------------------------
220      !
221      IF(lwp) THEN                         ! Parameter print
222         WRITE(numout,*)
223         WRITE(numout,*) '    zgr_z1d : Reference vertical z-coordinates: uniform dz = ', rn_e3
224         WRITE(numout,*) '    ~~~~~~~'
225      ENDIF
226      !
227      ! Reference z-coordinate (depth - scale factor at T- and W-points)   ! Madec & Imbard 1996 function
228      ! ----------------------
229      DO jk = 1, jpk
230         zw = REAL( jk , wp )
231         zt = REAL( jk , wp ) + 0.5_wp
232         pdepw_1d(jk) =    rn_e3 *   REAL( jk-1 , wp )
233         pdept_1d(jk) =    rn_e3 * ( REAL( jk-1 , wp ) + 0.5_wp )
234         pe3w_1d (jk) =    rn_e3
235         pe3t_1d (jk) =    rn_e3
236      END DO
237      !
238      IF(lwp) THEN                        ! control print
239         WRITE(numout,*)
240         WRITE(numout,*) '              Reference 1D z-coordinate depth and scale factors:'
241         WRITE(numout, "(9x,' level  gdept_1d  gdepw_1d  e3t_1d   e3w_1d  ')" )
242         WRITE(numout, "(10x, i4, 4f9.2)" ) ( jk, pdept_1d(jk), pdepw_1d(jk), pe3t_1d(jk), pe3w_1d(jk), jk = 1, jpk )
243      ENDIF
244      !
245   END SUBROUTINE zgr_z1d
246   
247   !!======================================================================
248END MODULE usrdef_zgr
Note: See TracBrowser for help on using the repository browser.