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/trunk/src/OCE/DOM – NEMO

source: NEMO/trunk/src/OCE/DOM/dommsk.F90 @ 12377

Last change on this file since 12377 was 12377, checked in by acc, 4 years 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: 13.9 KB
RevLine 
[3]1MODULE dommsk
[1566]2   !!======================================================================
[3]3   !!                       ***  MODULE dommsk   ***
4   !! Ocean initialization : domain land/sea mask
[1566]5   !!======================================================================
6   !! History :  OPA  ! 1987-07  (G. Madec)  Original code
[2528]7   !!            6.0  ! 1993-03  (M. Guyon)  symetrical conditions (M. Guyon)
8   !!            7.0  ! 1996-01  (G. Madec)  suppression of common work arrays
[6140]9   !!             -   ! 1996-05  (G. Madec)  mask computed from tmask
[2528]10   !!            8.0  ! 1997-02  (G. Madec)  mesh information put in domhgr.F
[7646]11   !!            8.1  ! 1997-07  (G. Madec)  modification of kbat and fmask
[1566]12   !!             -   ! 1998-05  (G. Roullet)  free surface
[2528]13   !!            8.2  ! 2000-03  (G. Madec)  no slip accurate
[1566]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
[5836]18   !!            3.6  ! 2015-05  (P. Mathiot) ISF: add wmask,wumask and wvmask
[7646]19   !!            4.0  ! 2016-06  (G. Madec, S. Flavoni)  domain configuration / user defined interface
[1566]20   !!----------------------------------------------------------------------
[3]21
22   !!----------------------------------------------------------------------
[7646]23   !!   dom_msk       : compute land/ocean mask
[3]24   !!----------------------------------------------------------------------
[7646]25   USE oce            ! ocean dynamics and tracers
26   USE dom_oce        ! ocean space and time domain
27   USE usrdef_fmask   ! user defined fmask
[9124]28   USE bdy_oce        ! open boundary
29   !
[7646]30   USE in_out_manager ! I/O manager
[9600]31   USE iom            ! IOM library
[7646]32   USE lbclnk         ! ocean lateral boundary conditions (or mpp link)
33   USE lib_mpp        ! Massively Parallel Processing library
[3]34
35   IMPLICIT NONE
36   PRIVATE
37
[6140]38   PUBLIC   dom_msk    ! routine called by inidom.F90
[3]39
[1601]40   !                            !!* Namelist namlbc : lateral boundary condition *
[4147]41   REAL(wp)        :: rn_shlat   ! type of lateral boundary condition on velocity
42   LOGICAL, PUBLIC :: ln_vorlat  !  consistency of vorticity boundary condition
[3294]43   !                                            with analytical eqs.
[2715]44
[3]45   !! * Substitutions
[12377]46#  include "do_loop_substitute.h90"
[1566]47   !!----------------------------------------------------------------------
[9598]48   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
[7753]49   !! $Id$
[10068]50   !! Software governed by the CeCILL license (see ./LICENSE)
[1566]51   !!----------------------------------------------------------------------
[3]52CONTAINS
[2715]53
[7646]54   SUBROUTINE dom_msk( k_top, k_bot )
[3]55      !!---------------------------------------------------------------------
56      !!                 ***  ROUTINE dom_msk  ***
57      !!
58      !! ** Purpose :   Compute land/ocean mask arrays at tracer points, hori-
[6125]59      !!      zontal velocity points (u & v), vorticity points (f) points.
[3]60      !!
[7646]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) :
[1601]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
[3]72      !!
[7646]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.
[3]78      !!
[7646]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
[3]86      !!----------------------------------------------------------------------
[7646]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
[9019]94      REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   zwf   ! 2D workspace
[1601]95      !!
[3294]96      NAMELIST/namlbc/ rn_shlat, ln_vorlat
[7646]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, &
[9657]101         &             cn_ice, nn_ice_dta,                                     &
[11536]102         &             ln_vol, nn_volctl, nn_rimwidth
[3]103      !!---------------------------------------------------------------------
[3294]104      !
[4147]105      READ  ( numnam_ref, namlbc, IOSTAT = ios, ERR = 901 )
[11536]106901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namlbc in reference namelist' )
[4147]107      READ  ( numnam_cfg, namlbc, IOSTAT = ios, ERR = 902 )
[11536]108902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namlbc in configuration namelist' )
[4624]109      IF(lwm) WRITE ( numond, namlbc )
[1566]110     
111      IF(lwp) THEN                  ! control print
[3]112         WRITE(numout,*)
113         WRITE(numout,*) 'dommsk : ocean mask '
114         WRITE(numout,*) '~~~~~~'
[1566]115         WRITE(numout,*) '   Namelist namlbc'
[3294]116         WRITE(numout,*) '      lateral momentum boundary cond.    rn_shlat  = ',rn_shlat
117         WRITE(numout,*) '      consistency with analytical form   ln_vorlat = ',ln_vorlat 
[3]118      ENDIF
[9169]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'
[1601]125      ELSE
[9527]126         CALL ctl_stop( 'dom_msk: wrong value for rn_shlat (i.e. a negalive value). We stop.' )
[3]127      ENDIF
128
[7646]129      !  Ocean/land mask at t-point  (computed from ko_top and ko_bot)
130      ! ----------------------------
[1566]131      !
[7753]132      tmask(:,:,:) = 0._wp
[12377]133      DO_2D_11_11
134         iktop = k_top(ji,jj)
135         ikbot = k_bot(ji,jj)
136         IF( iktop /= 0 ) THEN       ! water in the column
137            tmask(ji,jj,iktop:ikbot  ) = 1._wp
138         ENDIF
139      END_2D
[11233]140      !
141      ! the following call is mandatory
142      ! it masks boundaries (bathy=0) where needed depending on the configuration (closed, periodic...) 
[10425]143      CALL lbc_lnk( 'dommsk', tmask  , 'T', 1._wp )      ! Lateral boundary conditions
[3]144
[7646]145     ! Mask corrections for bdy (read in mppini2)
146      READ  ( numnam_ref, nambdy, IOSTAT = ios, ERR = 903)
[11536]147903   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nambdy in reference namelist' )
[7646]148      READ  ( numnam_cfg, nambdy, IOSTAT = ios, ERR = 904 )
[11536]149904   IF( ios >  0 )   CALL ctl_nam ( ios , 'nambdy in configuration namelist' )
[7646]150      ! ------------------------
151      IF ( ln_bdy .AND. ln_mask_file ) THEN
[9600]152         CALL iom_open( cn_mask_file, inum )
153         CALL iom_get ( inum, jpdom_data, 'bdy_msk', bdytmask(:,:) )
154         CALL iom_close( inum )
[12377]155         DO_3D_11_11( 1, jpkm1 )
156            tmask(ji,jj,jk) = tmask(ji,jj,jk) * bdytmask(ji,jj)
157         END_3D
[3]158      ENDIF
[7646]159         
160      !  Ocean/land mask at u-, v-, and f-points   (computed from tmask)
161      ! ----------------------------------------
162      ! NB: at this point, fmask is designed for free slip lateral boundary condition
[3]163      DO jk = 1, jpk
164         DO jj = 1, jpjm1
[12377]165            DO ji = 1, jpim1   ! vector loop
[3]166               umask(ji,jj,jk) = tmask(ji,jj  ,jk) * tmask(ji+1,jj  ,jk)
167               vmask(ji,jj,jk) = tmask(ji,jj  ,jk) * tmask(ji  ,jj+1,jk)
[1271]168            END DO
[1694]169            DO ji = 1, jpim1      ! NO vector opt.
[3]170               fmask(ji,jj,jk) = tmask(ji,jj  ,jk) * tmask(ji+1,jj  ,jk)   &
[62]171                  &            * tmask(ji,jj+1,jk) * tmask(ji+1,jj+1,jk)
[3]172            END DO
173         END DO
174      END DO
[10425]175      CALL lbc_lnk_multi( 'dommsk', umask, 'U', 1., vmask, 'V', 1., fmask, 'F', 1. )      ! Lateral boundary conditions
[7646]176 
177      ! Ocean/land mask at wu-, wv- and w points    (computed from tmask)
178      !-----------------------------------------
[7753]179      wmask (:,:,1) = tmask(:,:,1)     ! surface
180      wumask(:,:,1) = umask(:,:,1)
181      wvmask(:,:,1) = vmask(:,:,1)
[5836]182      DO jk = 2, jpk                   ! interior values
[7753]183         wmask (:,:,jk) = tmask(:,:,jk) * tmask(:,:,jk-1)
184         wumask(:,:,jk) = umask(:,:,jk) * umask(:,:,jk-1)   
185         wvmask(:,:,jk) = vmask(:,:,jk) * vmask(:,:,jk-1)
[5120]186      END DO
[3]187
[7646]188
189      ! Ocean/land column mask at t-, u-, and v-points   (i.e. at least 1 wet cell in the vertical)
190      ! ----------------------------------------------
191      ssmask (:,:) = MAXVAL( tmask(:,:,:), DIM=3 )
192      ssumask(:,:) = MAXVAL( umask(:,:,:), DIM=3 )
193      ssvmask(:,:) = MAXVAL( vmask(:,:,:), DIM=3 )
194
195
196      ! Interior domain mask  (used for global sum)
197      ! --------------------
198      !
[9019]199      iif = nn_hls   ;   iil = nlci - nn_hls + 1
200      ijf = nn_hls   ;   ijl = nlcj - nn_hls + 1
[7646]201      !
202      !                          ! halo mask : 0 on the halo and 1 elsewhere
[7753]203      tmask_h(:,:) = 1._wp                 
[7646]204      tmask_h( 1 :iif,   :   ) = 0._wp      ! first columns
205      tmask_h(iil:jpi,   :   ) = 0._wp      ! last  columns (including mpp extra columns)
206      tmask_h(   :   , 1 :ijf) = 0._wp      ! first rows
207      tmask_h(   :   ,ijl:jpj) = 0._wp      ! last  rows (including mpp extra rows)
208      !
209      !                          ! north fold mask
210      tpol(1:jpiglo) = 1._wp 
211      fpol(1:jpiglo) = 1._wp
212      IF( jperio == 3 .OR. jperio == 4 ) THEN      ! T-point pivot
213         tpol(jpiglo/2+1:jpiglo) = 0._wp
214         fpol(     1    :jpiglo) = 0._wp
215         IF( mjg(nlej) == jpjglo ) THEN                  ! only half of the nlcj-1 row for tmask_h
216            DO ji = iif+1, iil-1
217               tmask_h(ji,nlej-1) = tmask_h(ji,nlej-1) * tpol(mig(ji))
218            END DO
219         ENDIF
220      ENDIF
221      !
222      IF( jperio == 5 .OR. jperio == 6 ) THEN      ! F-point pivot
223         tpol(     1    :jpiglo) = 0._wp
224         fpol(jpiglo/2+1:jpiglo) = 0._wp
225      ENDIF
226      !
227      !                          ! interior mask : 2D ocean mask x halo mask
[7753]228      tmask_i(:,:) = ssmask(:,:) * tmask_h(:,:)
[7646]229
230
[3]231      ! Lateral boundary conditions on velocity (modify fmask)
[7646]232      ! --------------------------------------- 
233      IF( rn_shlat /= 0 ) THEN      ! Not free-slip lateral boundary condition
234         !
[9019]235         ALLOCATE( zwf(jpi,jpj) )
[7646]236         !
237         DO jk = 1, jpk
[7753]238            zwf(:,:) = fmask(:,:,jk)         
[12377]239            DO_2D_00_00
240               IF( fmask(ji,jj,jk) == 0._wp ) THEN
241                  fmask(ji,jj,jk) = rn_shlat * MIN( 1._wp , MAX( zwf(ji+1,jj), zwf(ji,jj+1),   &
242                     &                                           zwf(ji-1,jj), zwf(ji,jj-1)  )  )
243               ENDIF
244            END_2D
[7646]245            DO jj = 2, jpjm1
246               IF( fmask(1,jj,jk) == 0._wp ) THEN
247                  fmask(1  ,jj,jk) = rn_shlat * MIN( 1._wp , MAX( zwf(2,jj), zwf(1,jj+1), zwf(1,jj-1) ) )
[3]248               ENDIF
[7646]249               IF( fmask(jpi,jj,jk) == 0._wp ) THEN
250                  fmask(jpi,jj,jk) = rn_shlat * MIN( 1._wp , MAX( zwf(jpi,jj+1), zwf(jpim1,jj), zwf(jpi,jj-1) ) )
251               ENDIF
252            END DO         
253            DO ji = 2, jpim1
254               IF( fmask(ji,1,jk) == 0._wp ) THEN
255                  fmask(ji, 1 ,jk) = rn_shlat * MIN( 1._wp , MAX( zwf(ji+1,1), zwf(ji,2), zwf(ji-1,1) ) )
256               ENDIF
257               IF( fmask(ji,jpj,jk) == 0._wp ) THEN
258                  fmask(ji,jpj,jk) = rn_shlat * MIN( 1._wp , MAX( zwf(ji+1,jpj), zwf(ji-1,jpj), zwf(ji,jpjm1) ) )
259               ENDIF
[3]260            END DO
[9023]261#if defined key_agrif 
[9090]262            IF( .NOT. AGRIF_Root() ) THEN
263               IF ((nbondi ==  1).OR.(nbondi == 2)) fmask(nlci-1 , :     ,jk) = 0.e0      ! east
264               IF ((nbondi == -1).OR.(nbondi == 2)) fmask(1      , :     ,jk) = 0.e0      ! west
265               IF ((nbondj ==  1).OR.(nbondj == 2)) fmask(:      ,nlcj-1 ,jk) = 0.e0      ! north
266               IF ((nbondj == -1).OR.(nbondj == 2)) fmask(:      ,1      ,jk) = 0.e0      ! south
267            ENDIF 
[9023]268#endif
[3]269         END DO
[5836]270         !
[9019]271         DEALLOCATE( zwf )
[5836]272         !
[10425]273         CALL lbc_lnk( 'dommsk', fmask, 'F', 1._wp )      ! Lateral boundary conditions on fmask
[3]274         !
[7646]275         ! CAUTION : The fmask may be further modified in dyn_vor_init ( dynvor.F90 ) depending on ln_vorlat
[5506]276         !
[2528]277      ENDIF
[7646]278     
279      ! User defined alteration of fmask (use to reduce ocean transport in specified straits)
280      ! --------------------------------
[2528]281      !
[7646]282      CALL usr_def_fmask( cn_cfg, nn_cfg, fmask )
[6140]283      !
[3]284   END SUBROUTINE dom_msk
285   
286   !!======================================================================
287END MODULE dommsk
Note: See TracBrowser for help on using the repository browser.