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.
dommsk.F90 in NEMO/branches/2019/dev_r11613_ENHANCE-04_namelists_as_internalfiles/src/OCE/DOM – NEMO

source: NEMO/branches/2019/dev_r11613_ENHANCE-04_namelists_as_internalfiles/src/OCE/DOM/dommsk.F90 @ 11671

Last change on this file since 11671 was 11671, checked in by acc, 5 years ago

Branch 2019/dev_r11613_ENHANCE-04_namelists_as_internalfiles. Final, non-substantive changes to complete this branch. These changes remove all REWIND statements on the old namelist fortran units (now character variables for internal files). These changes have been left until last since they are easily repeated via a script and it may be preferable to use the previous revision for merge purposes and reapply these last changes separately. This branch has been fully SETTE tested.
  • Property svn:keywords set to Id
File size: 14.1 KB
Line 
1MODULE dommsk
2   !!======================================================================
3   !!                       ***  MODULE dommsk   ***
4   !! Ocean initialization : domain land/sea mask
5   !!======================================================================
6   !! History :  OPA  ! 1987-07  (G. Madec)  Original code
7   !!            6.0  ! 1993-03  (M. Guyon)  symetrical conditions (M. Guyon)
8   !!            7.0  ! 1996-01  (G. Madec)  suppression of common work arrays
9   !!             -   ! 1996-05  (G. Madec)  mask computed from tmask
10   !!            8.0  ! 1997-02  (G. Madec)  mesh information put in domhgr.F
11   !!            8.1  ! 1997-07  (G. Madec)  modification of kbat and fmask
12   !!             -   ! 1998-05  (G. Roullet)  free surface
13   !!            8.2  ! 2000-03  (G. Madec)  no slip accurate
14   !!             -   ! 2001-09  (J.-M. Molines)  Open boundaries
15   !!   NEMO     1.0  ! 2002-08  (G. Madec)  F90: Free form and module
16   !!             -   ! 2005-11  (V. Garnier) Surface pressure gradient organization
17   !!            3.2  ! 2009-07  (R. Benshila) Suppression of rigid-lid option
18   !!            3.6  ! 2015-05  (P. Mathiot) ISF: add wmask,wumask and wvmask
19   !!            4.0  ! 2016-06  (G. Madec, S. Flavoni)  domain configuration / user defined interface
20   !!----------------------------------------------------------------------
21
22   !!----------------------------------------------------------------------
23   !!   dom_msk       : compute land/ocean mask
24   !!----------------------------------------------------------------------
25   USE oce            ! ocean dynamics and tracers
26   USE dom_oce        ! ocean space and time domain
27   USE usrdef_fmask   ! user defined fmask
28   USE bdy_oce        ! open boundary
29   !
30   USE in_out_manager ! I/O manager
31   USE iom            ! IOM library
32   USE lbclnk         ! ocean lateral boundary conditions (or mpp link)
33   USE lib_mpp        ! Massively Parallel Processing library
34
35   IMPLICIT NONE
36   PRIVATE
37
38   PUBLIC   dom_msk    ! routine called by inidom.F90
39
40   !                            !!* Namelist namlbc : lateral boundary condition *
41   REAL(wp)        :: rn_shlat   ! type of lateral boundary condition on velocity
42   LOGICAL, PUBLIC :: ln_vorlat  !  consistency of vorticity boundary condition
43   !                                            with analytical eqs.
44
45   !! * Substitutions
46#  include "vectopt_loop_substitute.h90"
47   !!----------------------------------------------------------------------
48   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
49   !! $Id$
50   !! Software governed by the CeCILL license (see ./LICENSE)
51   !!----------------------------------------------------------------------
52CONTAINS
53
54   SUBROUTINE dom_msk( k_top, k_bot )
55      !!---------------------------------------------------------------------
56      !!                 ***  ROUTINE dom_msk  ***
57      !!
58      !! ** Purpose :   Compute land/ocean mask arrays at tracer points, hori-
59      !!      zontal velocity points (u & v), vorticity points (f) points.
60      !!
61      !! ** Method  :   The ocean/land mask  at t-point is deduced from ko_top
62      !!      and ko_bot, the indices of the fist and last ocean t-levels which
63      !!      are either defined in usrdef_zgr or read in zgr_read.
64      !!                The velocity masks (umask, vmask, wmask, wumask, wvmask)
65      !!      are deduced from a product of the two neighboring tmask.
66      !!                The vorticity mask (fmask) is deduced from tmask taking
67      !!      into account the choice of lateral boundary condition (rn_shlat) :
68      !!         rn_shlat = 0, free slip  (no shear along the coast)
69      !!         rn_shlat = 2, no slip  (specified zero velocity at the coast)
70      !!         0 < rn_shlat < 2, partial slip   | non-linear velocity profile
71      !!         2 < rn_shlat, strong slip        | in the lateral boundary layer
72      !!
73      !!      tmask_i : interior ocean mask at t-point, i.e. excluding duplicated
74      !!                rows/lines due to cyclic or North Fold boundaries as well
75      !!                as MPP halos.
76      !!      tmask_h : halo mask at t-point, i.e. excluding duplicated rows/lines
77      !!                due to cyclic or North Fold boundaries as well as MPP halos.
78      !!
79      !! ** Action :   tmask, umask, vmask, wmask, wumask, wvmask : land/ocean mask
80      !!                         at t-, u-, v- w, wu-, and wv-points (=0. or 1.)
81      !!               fmask   : land/ocean mask at f-point (=0., or =1., or
82      !!                         =rn_shlat along lateral boundaries)
83      !!               tmask_i : interior ocean mask
84      !!               tmask_h : halo mask
85      !!               ssmask , ssumask, ssvmask, ssfmask : 2D ocean mask
86      !!----------------------------------------------------------------------
87      INTEGER, DIMENSION(:,:), INTENT(in) ::   k_top, k_bot   ! first and last ocean level
88      !
89      INTEGER  ::   ji, jj, jk     ! dummy loop indices
90      INTEGER  ::   iif, iil       ! local integers
91      INTEGER  ::   ijf, ijl       !   -       -
92      INTEGER  ::   iktop, ikbot   !   -       -
93      INTEGER  ::   ios, inum
94      REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   zwf   ! 2D workspace
95      !!
96      NAMELIST/namlbc/ rn_shlat, ln_vorlat
97      NAMELIST/nambdy/ ln_bdy ,nb_bdy, ln_coords_file, cn_coords_file,         &
98         &             ln_mask_file, cn_mask_file, cn_dyn2d, nn_dyn2d_dta,     &
99         &             cn_dyn3d, nn_dyn3d_dta, cn_tra, nn_tra_dta,             &
100         &             ln_tra_dmp, ln_dyn3d_dmp, rn_time_dmp, rn_time_dmp_out, &
101         &             cn_ice, nn_ice_dta,                                     &
102         &             ln_vol, nn_volctl, nn_rimwidth
103      !!---------------------------------------------------------------------
104      !
105      READ  ( numnam_ref, namlbc, IOSTAT = ios, ERR = 901 )
106901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namlbc in reference namelist' )
107      READ  ( numnam_cfg, namlbc, IOSTAT = ios, ERR = 902 )
108902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namlbc in configuration namelist' )
109      IF(lwm) WRITE ( numond, namlbc )
110     
111      IF(lwp) THEN                  ! control print
112         WRITE(numout,*)
113         WRITE(numout,*) 'dommsk : ocean mask '
114         WRITE(numout,*) '~~~~~~'
115         WRITE(numout,*) '   Namelist namlbc'
116         WRITE(numout,*) '      lateral momentum boundary cond.    rn_shlat  = ',rn_shlat
117         WRITE(numout,*) '      consistency with analytical form   ln_vorlat = ',ln_vorlat 
118      ENDIF
119      !
120      IF(lwp) WRITE(numout,*)
121      IF     (      rn_shlat == 0.               ) THEN   ;   IF(lwp) WRITE(numout,*) '   ==>>>   ocean lateral  free-slip'
122      ELSEIF (      rn_shlat == 2.               ) THEN   ;   IF(lwp) WRITE(numout,*) '   ==>>>   ocean lateral  no-slip'
123      ELSEIF ( 0. < rn_shlat .AND. rn_shlat < 2. ) THEN   ;   IF(lwp) WRITE(numout,*) '   ==>>>   ocean lateral  partial-slip'
124      ELSEIF ( 2. < rn_shlat                     ) THEN   ;   IF(lwp) WRITE(numout,*) '   ==>>>   ocean lateral  strong-slip'
125      ELSE
126         CALL ctl_stop( 'dom_msk: wrong value for rn_shlat (i.e. a negalive value). We stop.' )
127      ENDIF
128
129      !  Ocean/land mask at t-point  (computed from ko_top and ko_bot)
130      ! ----------------------------
131      !
132      tmask(:,:,:) = 0._wp
133      DO jj = 1, jpj
134         DO ji = 1, jpi
135            iktop = k_top(ji,jj)
136            ikbot = k_bot(ji,jj)
137            IF( iktop /= 0 ) THEN       ! water in the column
138               tmask(ji,jj,iktop:ikbot  ) = 1._wp
139            ENDIF
140         END DO 
141      END DO
142      !
143      ! the following call is mandatory
144      ! it masks boundaries (bathy=0) where needed depending on the configuration (closed, periodic...) 
145      CALL lbc_lnk( 'dommsk', tmask  , 'T', 1._wp )      ! Lateral boundary conditions
146
147     ! Mask corrections for bdy (read in mppini2)
148      READ  ( numnam_ref, nambdy, IOSTAT = ios, ERR = 903)
149903   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nambdy in reference namelist' )
150      READ  ( numnam_cfg, nambdy, IOSTAT = ios, ERR = 904 )
151904   IF( ios >  0 )   CALL ctl_nam ( ios , 'nambdy in configuration namelist' )
152      ! ------------------------
153      IF ( ln_bdy .AND. ln_mask_file ) THEN
154         CALL iom_open( cn_mask_file, inum )
155         CALL iom_get ( inum, jpdom_data, 'bdy_msk', bdytmask(:,:) )
156         CALL iom_close( inum )
157         DO jk = 1, jpkm1
158            DO jj = 1, jpj
159               DO ji = 1, jpi
160                  tmask(ji,jj,jk) = tmask(ji,jj,jk) * bdytmask(ji,jj)
161               END DO
162            END DO
163         END DO
164      ENDIF
165         
166      !  Ocean/land mask at u-, v-, and f-points   (computed from tmask)
167      ! ----------------------------------------
168      ! NB: at this point, fmask is designed for free slip lateral boundary condition
169      DO jk = 1, jpk
170         DO jj = 1, jpjm1
171            DO ji = 1, fs_jpim1   ! vector loop
172               umask(ji,jj,jk) = tmask(ji,jj  ,jk) * tmask(ji+1,jj  ,jk)
173               vmask(ji,jj,jk) = tmask(ji,jj  ,jk) * tmask(ji  ,jj+1,jk)
174            END DO
175            DO ji = 1, jpim1      ! NO vector opt.
176               fmask(ji,jj,jk) = tmask(ji,jj  ,jk) * tmask(ji+1,jj  ,jk)   &
177                  &            * tmask(ji,jj+1,jk) * tmask(ji+1,jj+1,jk)
178            END DO
179         END DO
180      END DO
181      CALL lbc_lnk_multi( 'dommsk', umask, 'U', 1., vmask, 'V', 1., fmask, 'F', 1. )      ! Lateral boundary conditions
182 
183      ! Ocean/land mask at wu-, wv- and w points    (computed from tmask)
184      !-----------------------------------------
185      wmask (:,:,1) = tmask(:,:,1)     ! surface
186      wumask(:,:,1) = umask(:,:,1)
187      wvmask(:,:,1) = vmask(:,:,1)
188      DO jk = 2, jpk                   ! interior values
189         wmask (:,:,jk) = tmask(:,:,jk) * tmask(:,:,jk-1)
190         wumask(:,:,jk) = umask(:,:,jk) * umask(:,:,jk-1)   
191         wvmask(:,:,jk) = vmask(:,:,jk) * vmask(:,:,jk-1)
192      END DO
193
194
195      ! Ocean/land column mask at t-, u-, and v-points   (i.e. at least 1 wet cell in the vertical)
196      ! ----------------------------------------------
197      ssmask (:,:) = MAXVAL( tmask(:,:,:), DIM=3 )
198      ssumask(:,:) = MAXVAL( umask(:,:,:), DIM=3 )
199      ssvmask(:,:) = MAXVAL( vmask(:,:,:), DIM=3 )
200
201
202      ! Interior domain mask  (used for global sum)
203      ! --------------------
204      !
205      iif = nn_hls   ;   iil = nlci - nn_hls + 1
206      ijf = nn_hls   ;   ijl = nlcj - nn_hls + 1
207      !
208      !                          ! halo mask : 0 on the halo and 1 elsewhere
209      tmask_h(:,:) = 1._wp                 
210      tmask_h( 1 :iif,   :   ) = 0._wp      ! first columns
211      tmask_h(iil:jpi,   :   ) = 0._wp      ! last  columns (including mpp extra columns)
212      tmask_h(   :   , 1 :ijf) = 0._wp      ! first rows
213      tmask_h(   :   ,ijl:jpj) = 0._wp      ! last  rows (including mpp extra rows)
214      !
215      !                          ! north fold mask
216      tpol(1:jpiglo) = 1._wp 
217      fpol(1:jpiglo) = 1._wp
218      IF( jperio == 3 .OR. jperio == 4 ) THEN      ! T-point pivot
219         tpol(jpiglo/2+1:jpiglo) = 0._wp
220         fpol(     1    :jpiglo) = 0._wp
221         IF( mjg(nlej) == jpjglo ) THEN                  ! only half of the nlcj-1 row for tmask_h
222            DO ji = iif+1, iil-1
223               tmask_h(ji,nlej-1) = tmask_h(ji,nlej-1) * tpol(mig(ji))
224            END DO
225         ENDIF
226      ENDIF
227      !
228      IF( jperio == 5 .OR. jperio == 6 ) THEN      ! F-point pivot
229         tpol(     1    :jpiglo) = 0._wp
230         fpol(jpiglo/2+1:jpiglo) = 0._wp
231      ENDIF
232      !
233      !                          ! interior mask : 2D ocean mask x halo mask
234      tmask_i(:,:) = ssmask(:,:) * tmask_h(:,:)
235
236
237      ! Lateral boundary conditions on velocity (modify fmask)
238      ! --------------------------------------- 
239      IF( rn_shlat /= 0 ) THEN      ! Not free-slip lateral boundary condition
240         !
241         ALLOCATE( zwf(jpi,jpj) )
242         !
243         DO jk = 1, jpk
244            zwf(:,:) = fmask(:,:,jk)         
245            DO jj = 2, jpjm1
246               DO ji = fs_2, fs_jpim1   ! vector opt.
247                  IF( fmask(ji,jj,jk) == 0._wp ) THEN
248                     fmask(ji,jj,jk) = rn_shlat * MIN( 1._wp , MAX( zwf(ji+1,jj), zwf(ji,jj+1),   &
249                        &                                           zwf(ji-1,jj), zwf(ji,jj-1)  )  )
250                  ENDIF
251               END DO
252            END DO
253            DO jj = 2, jpjm1
254               IF( fmask(1,jj,jk) == 0._wp ) THEN
255                  fmask(1  ,jj,jk) = rn_shlat * MIN( 1._wp , MAX( zwf(2,jj), zwf(1,jj+1), zwf(1,jj-1) ) )
256               ENDIF
257               IF( fmask(jpi,jj,jk) == 0._wp ) THEN
258                  fmask(jpi,jj,jk) = rn_shlat * MIN( 1._wp , MAX( zwf(jpi,jj+1), zwf(jpim1,jj), zwf(jpi,jj-1) ) )
259               ENDIF
260            END DO         
261            DO ji = 2, jpim1
262               IF( fmask(ji,1,jk) == 0._wp ) THEN
263                  fmask(ji, 1 ,jk) = rn_shlat * MIN( 1._wp , MAX( zwf(ji+1,1), zwf(ji,2), zwf(ji-1,1) ) )
264               ENDIF
265               IF( fmask(ji,jpj,jk) == 0._wp ) THEN
266                  fmask(ji,jpj,jk) = rn_shlat * MIN( 1._wp , MAX( zwf(ji+1,jpj), zwf(ji-1,jpj), zwf(ji,jpjm1) ) )
267               ENDIF
268            END DO
269#if defined key_agrif 
270            IF( .NOT. AGRIF_Root() ) THEN
271               IF ((nbondi ==  1).OR.(nbondi == 2)) fmask(nlci-1 , :     ,jk) = 0.e0      ! east
272               IF ((nbondi == -1).OR.(nbondi == 2)) fmask(1      , :     ,jk) = 0.e0      ! west
273               IF ((nbondj ==  1).OR.(nbondj == 2)) fmask(:      ,nlcj-1 ,jk) = 0.e0      ! north
274               IF ((nbondj == -1).OR.(nbondj == 2)) fmask(:      ,1      ,jk) = 0.e0      ! south
275            ENDIF 
276#endif
277         END DO
278         !
279         DEALLOCATE( zwf )
280         !
281         CALL lbc_lnk( 'dommsk', fmask, 'F', 1._wp )      ! Lateral boundary conditions on fmask
282         !
283         ! CAUTION : The fmask may be further modified in dyn_vor_init ( dynvor.F90 ) depending on ln_vorlat
284         !
285      ENDIF
286     
287      ! User defined alteration of fmask (use to reduce ocean transport in specified straits)
288      ! --------------------------------
289      !
290      CALL usr_def_fmask( cn_cfg, nn_cfg, fmask )
291      !
292   END SUBROUTINE dom_msk
293   
294   !!======================================================================
295END MODULE dommsk
Note: See TracBrowser for help on using the repository browser.