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.
Changeset 6225 for branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/OPA_SRC/DOM/dommsk.F90 – NEMO

Ignore:
Timestamp:
2016-01-08T10:35:19+01:00 (8 years ago)
Author:
jamesharle
Message:

Update MPP_BDY_UPDATE branch to be consistent with head of trunk

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2014/dev_r4704_NOC5_MPP_BDY_UPDATE/NEMOGCM/NEMO/OPA_SRC/DOM/dommsk.F90

    r4624 r6225  
    77   !!            6.0  ! 1993-03  (M. Guyon)  symetrical conditions (M. Guyon) 
    88   !!            7.0  ! 1996-01  (G. Madec)  suppression of common work arrays 
    9    !!             -   ! 1996-05  (G. Madec)  mask computed from tmask and sup- 
    10    !!                 !                      pression of the double computation of bmask 
     9   !!             -   ! 1996-05  (G. Madec)  mask computed from tmask 
    1110   !!            8.0  ! 1997-02  (G. Madec)  mesh information put in domhgr.F 
    1211   !!            8.1  ! 1997-07  (G. Madec)  modification of mbathy and fmask 
     
    1716   !!             -   ! 2005-11  (V. Garnier) Surface pressure gradient organization 
    1817   !!            3.2  ! 2009-07  (R. Benshila) Suppression of rigid-lid option 
     18   !!            3.6  ! 2015-05  (P. Mathiot) ISF: add wmask,wumask and wvmask 
    1919   !!---------------------------------------------------------------------- 
    2020 
    2121   !!---------------------------------------------------------------------- 
    2222   !!   dom_msk        : compute land/ocean mask 
    23    !!   dom_msk_nsa    : update land/ocean mask when no-slip accurate option is used. 
    2423   !!---------------------------------------------------------------------- 
    2524   USE oce             ! ocean dynamics and tracers 
    2625   USE dom_oce         ! ocean space and time domain 
     26   ! 
    2727   USE in_out_manager  ! I/O manager 
    2828   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
    29    USE lib_mpp 
    30    USE dynspg_oce      ! choice/control of key cpp for surface pressure gradient 
     29   USE lib_mpp         ! 
    3130   USE wrk_nemo        ! Memory allocation 
    3231   USE timing          ! Timing 
     
    3534   PRIVATE 
    3635 
    37    PUBLIC   dom_msk         ! routine called by inidom.F90 
    38    PUBLIC   dom_msk_alloc   ! routine called by nemogcm.F90 
     36   PUBLIC   dom_msk    ! routine called by inidom.F90 
    3937 
    4038   !                            !!* Namelist namlbc : lateral boundary condition * 
     
    4341   !                                            with analytical eqs. 
    4442 
    45  
    46    INTEGER, ALLOCATABLE, SAVE, DIMENSION(:,:) ::  icoord ! Workspace for dom_msk_nsa() 
    47  
    4843   !! * Substitutions 
    4944#  include "vectopt_loop_substitute.h90" 
     
    5449   !!---------------------------------------------------------------------- 
    5550CONTAINS 
    56     
    57    INTEGER FUNCTION dom_msk_alloc() 
    58       !!--------------------------------------------------------------------- 
    59       !!                 ***  FUNCTION dom_msk_alloc  *** 
    60       !!--------------------------------------------------------------------- 
    61       dom_msk_alloc = 0 
    62 #if defined key_noslip_accurate 
    63       ALLOCATE(icoord(jpi*jpj*jpk,3), STAT=dom_msk_alloc) 
    64 #endif 
    65       IF( dom_msk_alloc /= 0 )   CALL ctl_warn('dom_msk_alloc: failed to allocate icoord array') 
    66       ! 
    67    END FUNCTION dom_msk_alloc 
    68  
    6951 
    7052   SUBROUTINE dom_msk 
     
    7355      !! 
    7456      !! ** Purpose :   Compute land/ocean mask arrays at tracer points, hori- 
    75       !!      zontal velocity points (u & v), vorticity points (f) and baro- 
    76       !!      tropic stream function  points (b). 
     57      !!      zontal velocity points (u & v), vorticity points (f) points. 
    7758      !! 
    7859      !! ** Method  :   The ocean/land mask is computed from the basin bathy- 
     
    9273      !!                1. IF mbathy( ji ,jj) and mbathy( ji ,jj+1) 
    9374      !!                  and mbathy(ji+1,jj) and mbathy(ji+1,jj+1) >= jk. 
    94       !!      b-point : the same definition as for f-point of the first ocean 
    95       !!                level (surface level) but with 0 along coastlines. 
    9675      !!      tmask_i : interior ocean mask at t-point, i.e. excluding duplicated 
    9776      !!                rows/lines due to cyclic or North Fold boundaries as well 
     
    10786      !! 
    10887      !!      N.B. If nperio not equal to 0, the land/ocean mask arrays 
    109       !!      are defined with the proper value at lateral domain boundaries, 
    110       !!      but bmask. indeed, bmask defined the domain over which the 
    111       !!      barotropic stream function is computed. this domain cannot 
    112       !!      contain identical columns because the matrix associated with 
    113       !!      the barotropic stream function equation is then no more inverti- 
    114       !!      ble. therefore bmask is set to 0 along lateral domain boundaries 
    115       !!      even IF nperio is not zero. 
     88      !!      are defined with the proper value at lateral domain boundaries. 
    11689      !! 
    11790      !!      In case of open boundaries (lk_bdy=T): 
    11891      !!        - tmask is set to 1 on the points to be computed bay the open 
    11992      !!          boundaries routines. 
    120       !!        - bmask is  set to 0 on the open boundaries. 
    12193      !! 
    12294      !! ** Action :   tmask    : land/ocean mask at t-point (=0. or 1.) 
     
    12597      !!               fmask    : land/ocean mask at f-point (=0. or 1.) 
    12698      !!                          =rn_shlat along lateral boundaries 
    127       !!               bmask    : land/ocean mask at barotropic stream 
    128       !!                          function point (=0. or 1.) and set to 0 along lateral boundaries 
    12999      !!               tmask_i  : interior ocean mask 
    130100      !!---------------------------------------------------------------------- 
    131       ! 
    132       INTEGER  ::   ji, jj, jk      ! dummy loop indices 
     101      INTEGER  ::   ji, jj, jk               ! dummy loop indices 
    133102      INTEGER  ::   iif, iil, ii0, ii1, ii   ! local integers 
    134103      INTEGER  ::   ijf, ijl, ij0, ij1       !   -       - 
    135104      INTEGER  ::   ios 
     105      INTEGER  ::   isrow                    ! index for ORCA1 starting row 
    136106      INTEGER , POINTER, DIMENSION(:,:) ::  imsk 
    137107      REAL(wp), POINTER, DIMENSION(:,:) ::  zwf 
     
    184154         END DO   
    185155      END DO   
    186  
    187 !!gm  ???? 
    188 #if defined key_zdfkpp 
    189       IF( cp_cfg == 'orca' ) THEN 
    190          IF( jp_cfg == 2 )   THEN       ! land point on Bab el Mandeb zonal section 
    191             ij0 =  87   ;   ij1 =  88 
    192             ii0 = 160   ;   ii1 = 161 
    193             tmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) = 0._wp 
    194          ELSE 
    195             IF(lwp) WRITE(numout,*) 
    196             IF(lwp) WRITE(numout,cform_war) 
    197             IF(lwp) WRITE(numout,*) 
    198             IF(lwp) WRITE(numout,*)'          A mask must be applied on Bab el Mandeb strait' 
    199             IF(lwp) WRITE(numout,*)'          in case of ORCAs configurations' 
    200             IF(lwp) WRITE(numout,*)'          This is a problem which is not yet solved' 
    201             IF(lwp) WRITE(numout,*) 
    202          ENDIF 
    203       ENDIF 
    204 #endif 
    205 !!gm end 
     156       
     157      ! (ISF) define barotropic mask and mask the ice shelf point 
     158      ssmask(:,:)=tmask(:,:,1) ! at this stage ice shelf is not masked 
     159       
     160      DO jk = 1, jpk 
     161         DO jj = 1, jpj 
     162            DO ji = 1, jpi 
     163               IF( REAL( misfdep(ji,jj) - jk, wp ) - 0.1_wp >= 0._wp )   THEN 
     164                  tmask(ji,jj,jk) = 0._wp 
     165               END IF 
     166            END DO   
     167         END DO   
     168      END DO   
    206169 
    207170      ! Interior domain mask (used for global sum) 
    208171      ! -------------------- 
    209       tmask_i(:,:) = tmask(:,:,1) 
     172      tmask_i(:,:) = ssmask(:,:)            ! (ISH) tmask_i = 1 even on the ice shelf 
     173 
     174      tmask_h(:,:) = 1._wp                 ! 0 on the halo and 1 elsewhere 
    210175      iif = jpreci                         ! ??? 
    211176      iil = nlci - jpreci + 1 
     
    213178      ijl = nlcj - jprecj + 1 
    214179 
    215       tmask_i( 1 :iif,   :   ) = 0._wp      ! first columns 
    216       tmask_i(iil:jpi,   :   ) = 0._wp      ! last  columns (including mpp extra columns) 
    217       tmask_i(   :   , 1 :ijf) = 0._wp      ! first rows 
    218       tmask_i(   :   ,ijl:jpj) = 0._wp      ! last  rows (including mpp extra rows) 
     180      tmask_h( 1 :iif,   :   ) = 0._wp      ! first columns 
     181      tmask_h(iil:jpi,   :   ) = 0._wp      ! last  columns (including mpp extra columns) 
     182      tmask_h(   :   , 1 :ijf) = 0._wp      ! first rows 
     183      tmask_h(   :   ,ijl:jpj) = 0._wp      ! last  rows (including mpp extra rows) 
    219184 
    220185      ! north fold mask 
     
    227192         IF( mjg(nlej) == jpjglo ) THEN                  ! only half of the nlcj-1 row 
    228193            DO ji = iif+1, iil-1 
    229                tmask_i(ji,nlej-1) = tmask_i(ji,nlej-1) * tpol(mig(ji)) 
     194               tmask_h(ji,nlej-1) = tmask_h(ji,nlej-1) * tpol(mig(ji)) 
    230195            END DO 
    231196         ENDIF 
    232197      ENDIF 
     198      
     199      tmask_i(:,:) = tmask_i(:,:) * tmask_h(:,:) 
     200 
    233201      IF( jperio == 5 .OR. jperio == 6 ) THEN      ! F-point pivot 
    234202         tpol(     1    :jpiglo) = 0._wp 
     
    250218         END DO 
    251219      END DO 
    252       CALL lbc_lnk( umask, 'U', 1._wp )      ! Lateral boundary conditions 
    253       CALL lbc_lnk( vmask, 'V', 1._wp ) 
    254       CALL lbc_lnk( fmask, 'F', 1._wp ) 
    255  
    256  
    257       ! 4. ocean/land mask for the elliptic equation 
    258       ! -------------------------------------------- 
    259       bmask(:,:) = tmask(:,:,1)       ! elliptic equation is written at t-point 
    260       ! 
    261       !                               ! Boundary conditions 
    262       !                                    ! cyclic east-west : bmask must be set to 0. on rows 1 and jpi 
    263       IF( nperio == 1 .OR. nperio == 4 .OR. nperio == 6 ) THEN 
    264          bmask( 1 ,:) = 0._wp 
    265          bmask(jpi,:) = 0._wp 
    266       ENDIF 
    267       IF( nperio == 2 ) THEN               ! south symmetric :  bmask must be set to 0. on row 1 
    268          bmask(:, 1 ) = 0._wp 
    269       ENDIF 
    270       !                                    ! north fold :  
    271       IF( nperio == 3 .OR. nperio == 4 ) THEN   ! T-pt pivot : bmask set to 0. on row jpj and on half jpjglo-1 row 
    272          DO ji = 1, jpi                       
    273             ii = ji + nimpp - 1 
    274             bmask(ji,jpj-1) = bmask(ji,jpj-1) * tpol(ii) 
    275             bmask(ji,jpj  ) = 0._wp 
    276          END DO 
    277       ENDIF 
    278       IF( nperio == 5 .OR. nperio == 6 ) THEN   ! F-pt pivot and T-pt elliptic eq. : bmask set to 0. on row jpj 
    279          bmask(:,jpj) = 0._wp 
    280       ENDIF 
    281       ! 
    282       IF( lk_mpp ) THEN                    ! mpp specificities 
    283          !                                      ! bmask is set to zero on the overlap region 
    284          IF( nbondi /= -1 .AND. nbondi /= 2 )   bmask(  1 :jpreci,:) = 0._wp 
    285          IF( nbondi /=  1 .AND. nbondi /= 2 )   bmask(nlci:jpi   ,:) = 0._wp 
    286          IF( nbondj /= -1 .AND. nbondj /= 2 )   bmask(:,  1 :jprecj) = 0._wp 
    287          IF( nbondj /=  1 .AND. nbondj /= 2 )   bmask(:,nlcj:jpj   ) = 0._wp 
    288          ! 
    289          IF( npolj == 3 .OR. npolj == 4 ) THEN  ! north fold : bmask must be set to 0. on rows jpj-1 and jpj 
    290             DO ji = 1, nlci 
    291                ii = ji + nimpp - 1 
    292                bmask(ji,nlcj-1) = bmask(ji,nlcj-1) * tpol(ii) 
    293                bmask(ji,nlcj  ) = 0._wp 
    294             END DO 
    295          ENDIF 
    296          IF( npolj == 5 .OR. npolj == 6 ) THEN  ! F-pt pivot and T-pt elliptic eq. : bmask set to 0. on row jpj 
    297             DO ji = 1, nlci 
    298                bmask(ji,nlcj  ) = 0._wp 
    299             END DO 
    300          ENDIF 
    301       ENDIF 
    302  
    303  
    304       ! mask for second order calculation of vorticity 
    305       ! ---------------------------------------------- 
    306       CALL dom_msk_nsa 
    307  
    308        
     220      ! (ISF) MIN(1,SUM(umask)) is here to check if you have effectively at least 1 wet cell at u point 
     221      DO jj = 1, jpjm1 
     222         DO ji = 1, fs_jpim1   ! vector loop 
     223            ssumask(ji,jj)  = ssmask(ji,jj) * ssmask(ji+1,jj  )  * MIN(1._wp,SUM(umask(ji,jj,:))) 
     224            ssvmask(ji,jj)  = ssmask(ji,jj) * ssmask(ji  ,jj+1)  * MIN(1._wp,SUM(vmask(ji,jj,:))) 
     225         END DO 
     226         DO ji = 1, jpim1      ! NO vector opt. 
     227            ssfmask(ji,jj) =  ssmask(ji,jj  ) * ssmask(ji+1,jj  )   & 
     228               &            * ssmask(ji,jj+1) * ssmask(ji+1,jj+1) * MIN(1._wp,SUM(fmask(ji,jj,:))) 
     229         END DO 
     230      END DO 
     231      CALL lbc_lnk( umask  , 'U', 1._wp )      ! Lateral boundary conditions 
     232      CALL lbc_lnk( vmask  , 'V', 1._wp ) 
     233      CALL lbc_lnk( fmask  , 'F', 1._wp ) 
     234      CALL lbc_lnk( ssumask, 'U', 1._wp )      ! Lateral boundary conditions 
     235      CALL lbc_lnk( ssvmask, 'V', 1._wp ) 
     236      CALL lbc_lnk( ssfmask, 'F', 1._wp ) 
     237 
     238      ! 3. Ocean/land mask at wu-, wv- and w points  
     239      !---------------------------------------------- 
     240      wmask (:,:,1) = tmask(:,:,1)     ! surface 
     241      wumask(:,:,1) = umask(:,:,1) 
     242      wvmask(:,:,1) = vmask(:,:,1) 
     243      DO jk = 2, jpk                   ! interior values 
     244         wmask (:,:,jk) = tmask(:,:,jk) * tmask(:,:,jk-1) 
     245         wumask(:,:,jk) = umask(:,:,jk) * umask(:,:,jk-1)    
     246         wvmask(:,:,jk) = vmask(:,:,jk) * vmask(:,:,jk-1) 
     247      END DO 
     248 
    309249      ! Lateral boundary conditions on velocity (modify fmask) 
    310250      ! ---------------------------------------      
     
    339279      IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN   ! ORCA_R2 configuration 
    340280         !                                                 ! Increased lateral friction near of some straits 
    341          IF( nn_cla == 0 ) THEN 
    342             !                                ! Gibraltar strait  : partial slip (fmask=0.5) 
    343             ij0 = 101   ;   ij1 = 101 
    344             ii0 = 139   ;   ii1 = 140   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  0.5_wp 
    345             ij0 = 102   ;   ij1 = 102 
    346             ii0 = 139   ;   ii1 = 140   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  0.5_wp 
    347             ! 
    348             !                                ! Bab el Mandeb : partial slip (fmask=1) 
    349             ij0 =  87   ;   ij1 =  88 
    350             ii0 = 160   ;   ii1 = 160   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  1._wp 
    351             ij0 =  88   ;   ij1 =  88 
    352             ii0 = 159   ;   ii1 = 159   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  1._wp 
    353             ! 
    354          ENDIF 
     281         !                                ! Gibraltar strait  : partial slip (fmask=0.5) 
     282         ij0 = 101   ;   ij1 = 101 
     283         ii0 = 139   ;   ii1 = 140   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  0.5_wp 
     284         ij0 = 102   ;   ij1 = 102 
     285         ii0 = 139   ;   ii1 = 140   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  0.5_wp 
     286         ! 
     287         !                                ! Bab el Mandeb : partial slip (fmask=1) 
     288         ij0 =  87   ;   ij1 =  88 
     289         ii0 = 160   ;   ii1 = 160   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  1._wp 
     290         ij0 =  88   ;   ij1 =  88 
     291         ii0 = 159   ;   ii1 = 159   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  1._wp 
     292         ! 
    355293         !                                ! Danish straits  : strong slip (fmask > 2) 
    356294! We keep this as an example but it is instable in this case  
     
    364302      IF( cp_cfg == "orca" .AND. jp_cfg == 1 ) THEN   ! ORCA R1 configuration 
    365303         !                                                 ! Increased lateral friction near of some straits 
     304         ! This dirty section will be suppressed by simplification process: 
     305         ! all this will come back in input files 
     306         ! Currently these hard-wired indices relate to configuration with 
     307         ! extend grid (jpjglo=332) 
     308         ! 
     309         isrow = 332 - jpjglo 
     310         ! 
    366311         IF(lwp) WRITE(numout,*) 
    367312         IF(lwp) WRITE(numout,*) '   orca_r1: increase friction near the following straits : ' 
    368313         IF(lwp) WRITE(numout,*) '      Gibraltar ' 
    369          ii0 = 283   ;   ii1 = 284        ! Gibraltar Strait  
    370          ij0 = 200   ;   ij1 = 200   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1), 1:jpk ) = 2._wp   
     314         ii0 = 282           ;   ii1 = 283        ! Gibraltar Strait  
     315         ij0 = 241 - isrow   ;   ij1 = 241 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp   
    371316 
    372317         IF(lwp) WRITE(numout,*) '      Bhosporus ' 
    373          ii0 = 314   ;   ii1 = 315        ! Bhosporus Strait  
    374          ij0 = 208   ;   ij1 = 208   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1), 1:jpk ) = 2._wp   
     318         ii0 = 314           ;   ii1 = 315        ! Bhosporus Strait  
     319         ij0 = 248 - isrow   ;   ij1 = 248 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp   
    375320 
    376321         IF(lwp) WRITE(numout,*) '      Makassar (Top) ' 
    377          ii0 =  48   ;   ii1 =  48        ! Makassar Strait (Top)  
    378          ij0 = 149   ;   ij1 = 150   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1), 1:jpk ) = 3._wp   
     322         ii0 =  48           ;   ii1 =  48        ! Makassar Strait (Top)  
     323         ij0 = 189 - isrow   ;   ij1 = 190 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 3._wp   
    379324 
    380325         IF(lwp) WRITE(numout,*) '      Lombok ' 
    381          ii0 =  44   ;   ii1 =  44        ! Lombok Strait  
    382          ij0 = 124   ;   ij1 = 125   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1), 1:jpk ) = 2._wp   
     326         ii0 =  44           ;   ii1 =  44        ! Lombok Strait  
     327         ij0 = 164 - isrow   ;   ij1 = 165 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp   
    383328 
    384329         IF(lwp) WRITE(numout,*) '      Ombai ' 
    385          ii0 =  53   ;   ii1 =  53        ! Ombai Strait  
    386          ij0 = 124   ;   ij1 = 125   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1), 1:jpk ) = 2._wp   
     330         ii0 =  53           ;   ii1 =  53        ! Ombai Strait  
     331         ij0 = 164 - isrow   ;   ij1 = 165 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp   
    387332 
    388333         IF(lwp) WRITE(numout,*) '      Timor Passage ' 
    389          ii0 =  56   ;   ii1 =  56        ! Timor Passage  
    390          ij0 = 124   ;   ij1 = 125   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1), 1:jpk ) = 2._wp   
     334         ii0 =  56           ;   ii1 =  56        ! Timor Passage  
     335         ij0 = 164 - isrow   ;   ij1 = 165 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 2._wp   
    391336 
    392337         IF(lwp) WRITE(numout,*) '      West Halmahera ' 
    393          ii0 =  58   ;   ii1 =  58        ! West Halmahera Strait  
    394          ij0 = 141   ;   ij1 = 142   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1), 1:jpk ) = 3._wp   
     338         ii0 =  58           ;   ii1 =  58        ! West Halmahera Strait  
     339         ij0 = 181 - isrow   ;   ij1 = 182 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 3._wp   
    395340 
    396341         IF(lwp) WRITE(numout,*) '      East Halmahera ' 
    397          ii0 =  55   ;   ii1 =  55        ! East Halmahera Strait  
    398          ij0 = 141   ;   ij1 = 142   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1), 1:jpk ) = 3._wp   
     342         ii0 =  55           ;   ii1 =  55        ! East Halmahera Strait  
     343         ij0 = 181 - isrow   ;   ij1 = 182 - isrow   ;   fmask( mi0(ii0):mi1(ii1),mj0(ij0):mj1(ij1),1:jpk ) = 3._wp   
    399344         ! 
    400345      ENDIF 
    401346      ! 
    402347      CALL lbc_lnk( fmask, 'F', 1._wp )      ! Lateral boundary conditions on fmask 
    403  
     348      ! 
    404349      ! CAUTION : The fmask may be further modified in dyn_vor_init ( dynvor.F90 ) 
    405              
    406       IF( nprint == 1 .AND. lwp ) THEN      ! Control print 
    407          imsk(:,:) = INT( tmask_i(:,:) ) 
    408          WRITE(numout,*) ' tmask_i : ' 
    409          CALL prihin( imsk(:,:), jpi, jpj, 1, jpi, 1,   & 
    410                &                           1, jpj, 1, 1, numout) 
    411          WRITE (numout,*) 
    412          WRITE (numout,*) ' dommsk: tmask for each level' 
    413          WRITE (numout,*) ' ----------------------------' 
    414          DO jk = 1, jpk 
    415             imsk(:,:) = INT( tmask(:,:,jk) ) 
    416  
    417             WRITE(numout,*) 
    418             WRITE(numout,*) ' level = ',jk 
    419             CALL prihin( imsk(:,:), jpi, jpj, 1, jpi, 1,   & 
    420                &                              1, jpj, 1, 1, numout) 
    421          END DO 
    422          WRITE(numout,*) 
    423          WRITE(numout,*) ' dom_msk: vmask for each level' 
    424          WRITE(numout,*) ' -----------------------------' 
    425          DO jk = 1, jpk 
    426             imsk(:,:) = INT( vmask(:,:,jk) ) 
    427             WRITE(numout,*) 
    428             WRITE(numout,*) ' level = ',jk 
    429             CALL prihin( imsk(:,:), jpi, jpj, 1, jpi, 1,   & 
    430                &                              1, jpj, 1, 1, numout) 
    431          END DO 
    432          WRITE(numout,*) 
    433          WRITE(numout,*) ' dom_msk: fmask for each level' 
    434          WRITE(numout,*) ' -----------------------------' 
    435          DO jk = 1, jpk 
    436             imsk(:,:) = INT( fmask(:,:,jk) ) 
    437             WRITE(numout,*) 
    438             WRITE(numout,*) ' level = ',jk 
    439             CALL prihin( imsk(:,:), jpi, jpj, 1, jpi, 1,   & 
    440                &                              1, jpj, 1, 1, numout ) 
    441          END DO 
    442          WRITE(numout,*) 
    443          WRITE(numout,*) ' dom_msk: bmask ' 
    444          WRITE(numout,*) ' ---------------' 
    445          WRITE(numout,*) 
    446          imsk(:,:) = INT( bmask(:,:) ) 
    447          CALL prihin( imsk(:,:), jpi, jpj, 1, jpi, 1,   & 
    448             &                              1, jpj, 1, 1, numout ) 
    449       ENDIF 
    450350      ! 
    451351      CALL wrk_dealloc( jpi, jpj, imsk ) 
     
    455355      ! 
    456356   END SUBROUTINE dom_msk 
    457  
    458 #if defined key_noslip_accurate 
    459    !!---------------------------------------------------------------------- 
    460    !!   'key_noslip_accurate' :         accurate no-slip boundary condition 
    461    !!---------------------------------------------------------------------- 
    462     
    463    SUBROUTINE dom_msk_nsa 
    464       !!--------------------------------------------------------------------- 
    465       !!                 ***  ROUTINE dom_msk_nsa  *** 
    466       !!  
    467       !! ** Purpose : 
    468       !! 
    469       !! ** Method  : 
    470       !! 
    471       !! ** Action : 
    472       !!---------------------------------------------------------------------- 
    473       INTEGER  ::   ji, jj, jk, jl      ! dummy loop indices 
    474       INTEGER  ::   ine, inw, ins, inn, itest, ierror, iind, ijnd 
    475       REAL(wp) ::   zaa 
    476       !!--------------------------------------------------------------------- 
    477       ! 
    478       IF( nn_timing == 1 )  CALL timing_start('dom_msk_nsa') 
    479       ! 
    480       IF(lwp) WRITE(numout,*) 
    481       IF(lwp) WRITE(numout,*) 'dom_msk_nsa : noslip accurate boundary condition' 
    482       IF(lwp) WRITE(numout,*) '~~~~~~~~~~~   using Schchepetkin and O Brian scheme' 
    483       IF( lk_mpp )   CALL ctl_stop( ' mpp version is not yet implemented' ) 
    484  
    485       ! mask for second order calculation of vorticity 
    486       ! ---------------------------------------------- 
    487       ! noslip boundary condition: fmask=1  at convex corner, store 
    488       ! index of straight coast meshes ( 'west', refering to a coast, 
    489       ! means west of the ocean, aso) 
    490        
    491       DO jk = 1, jpk 
    492          DO jl = 1, 4 
    493             npcoa(jl,jk) = 0 
    494             DO ji = 1, 2*(jpi+jpj) 
    495                nicoa(ji,jl,jk) = 0 
    496                njcoa(ji,jl,jk) = 0 
    497             END DO 
    498          END DO 
    499       END DO 
    500        
    501       IF( jperio == 2 ) THEN 
    502          WRITE(numout,*) ' ' 
    503          WRITE(numout,*) ' symetric boundary conditions need special' 
    504          WRITE(numout,*) ' treatment not implemented. we stop.' 
    505          STOP 
    506       ENDIF 
    507        
    508       ! convex corners 
    509        
    510       DO jk = 1, jpkm1 
    511          DO jj = 1, jpjm1 
    512             DO ji = 1, jpim1 
    513                zaa = tmask(ji  ,jj,jk) + tmask(ji  ,jj+1,jk)   & 
    514                   &+ tmask(ji+1,jj,jk) + tmask(ji+1,jj+1,jk) 
    515                IF( ABS(zaa-3._wp) <= 0.1_wp )   fmask(ji,jj,jk) = 1._wp 
    516             END DO 
    517          END DO 
    518       END DO 
    519  
    520       ! north-south straight coast 
    521  
    522       DO jk = 1, jpkm1 
    523          inw = 0 
    524          ine = 0 
    525          DO jj = 2, jpjm1 
    526             DO ji = 2, jpim1 
    527                zaa = tmask(ji+1,jj,jk) + tmask(ji+1,jj+1,jk) 
    528                IF( ABS(zaa-2._wp) <= 0.1_wp .AND. fmask(ji,jj,jk) == 0._wp ) THEN 
    529                   inw = inw + 1 
    530                   nicoa(inw,1,jk) = ji 
    531                   njcoa(inw,1,jk) = jj 
    532                   IF( nprint == 1 ) WRITE(numout,*) ' west  : ', jk, inw, ji, jj 
    533                ENDIF 
    534                zaa = tmask(ji,jj,jk) + tmask(ji,jj+1,jk) 
    535                IF( ABS(zaa-2._wp) <= 0.1_wp .AND. fmask(ji,jj,jk) == 0._wp ) THEN 
    536                   ine = ine + 1 
    537                   nicoa(ine,2,jk) = ji 
    538                   njcoa(ine,2,jk) = jj 
    539                   IF( nprint == 1 ) WRITE(numout,*) ' east  : ', jk, ine, ji, jj 
    540                ENDIF 
    541             END DO 
    542          END DO 
    543          npcoa(1,jk) = inw 
    544          npcoa(2,jk) = ine 
    545       END DO 
    546  
    547       ! west-east straight coast 
    548  
    549       DO jk = 1, jpkm1 
    550          ins = 0 
    551          inn = 0 
    552          DO jj = 2, jpjm1 
    553             DO ji =2, jpim1 
    554                zaa = tmask(ji,jj+1,jk) + tmask(ji+1,jj+1,jk) 
    555                IF( ABS(zaa-2._wp) <= 0.1_wp .AND. fmask(ji,jj,jk) == 0._wp ) THEN 
    556                   ins = ins + 1 
    557                   nicoa(ins,3,jk) = ji 
    558                   njcoa(ins,3,jk) = jj 
    559                   IF( nprint == 1 ) WRITE(numout,*) ' south : ', jk, ins, ji, jj 
    560                ENDIF 
    561                zaa = tmask(ji+1,jj,jk) + tmask(ji,jj,jk) 
    562                IF( ABS(zaa-2._wp) <= 0.1_wp .AND. fmask(ji,jj,jk) == 0._wp ) THEN 
    563                   inn = inn + 1 
    564                   nicoa(inn,4,jk) = ji 
    565                   njcoa(inn,4,jk) = jj 
    566                   IF( nprint == 1 ) WRITE(numout,*) ' north : ', jk, inn, ji, jj 
    567                ENDIF 
    568             END DO 
    569          END DO 
    570          npcoa(3,jk) = ins 
    571          npcoa(4,jk) = inn 
    572       END DO 
    573  
    574       itest = 2 * ( jpi + jpj ) 
    575       DO jk = 1, jpk 
    576          IF( npcoa(1,jk) > itest .OR. npcoa(2,jk) > itest .OR.   & 
    577              npcoa(3,jk) > itest .OR. npcoa(4,jk) > itest ) THEN 
    578              
    579             WRITE(ctmp1,*) ' level jk = ',jk 
    580             WRITE(ctmp2,*) ' straight coast index arraies are too small.:' 
    581             WRITE(ctmp3,*) ' npe, npw, nps, npn = ', npcoa(1,jk), npcoa(2,jk),   & 
    582                 &                                     npcoa(3,jk), npcoa(4,jk) 
    583             WRITE(ctmp4,*) ' 2*(jpi+jpj) = ',itest,'. we stop.' 
    584             CALL ctl_stop( ctmp1, ctmp2, ctmp3, ctmp4 ) 
    585         ENDIF 
    586       END DO 
    587  
    588       ierror = 0 
    589       iind = 0 
    590       ijnd = 0 
    591       IF( nperio == 1 .OR. nperio == 4 .OR. nperio == 6 )   iind = 2 
    592       IF( nperio == 3 .OR. nperio == 4 .OR. nperio == 5 .OR. nperio == 6 )   ijnd = 2 
    593       DO jk = 1, jpk 
    594          DO jl = 1, npcoa(1,jk) 
    595             IF( nicoa(jl,1,jk)+3 > jpi+iind ) THEN 
    596                ierror = ierror+1 
    597                icoord(ierror,1) = nicoa(jl,1,jk) 
    598                icoord(ierror,2) = njcoa(jl,1,jk) 
    599                icoord(ierror,3) = jk 
    600             ENDIF 
    601          END DO 
    602          DO jl = 1, npcoa(2,jk) 
    603             IF(nicoa(jl,2,jk)-2 < 1-iind ) THEN 
    604                ierror = ierror + 1 
    605                icoord(ierror,1) = nicoa(jl,2,jk) 
    606                icoord(ierror,2) = njcoa(jl,2,jk) 
    607                icoord(ierror,3) = jk 
    608             ENDIF 
    609          END DO 
    610          DO jl = 1, npcoa(3,jk) 
    611             IF( njcoa(jl,3,jk)+3 > jpj+ijnd ) THEN 
    612                ierror = ierror + 1 
    613                icoord(ierror,1) = nicoa(jl,3,jk) 
    614                icoord(ierror,2) = njcoa(jl,3,jk) 
    615                icoord(ierror,3) = jk 
    616             ENDIF 
    617          END DO 
    618          DO jl = 1, npcoa(4,jk) 
    619             IF( njcoa(jl,4,jk)-2 < 1) THEN 
    620                ierror=ierror + 1 
    621                icoord(ierror,1) = nicoa(jl,4,jk) 
    622                icoord(ierror,2) = njcoa(jl,4,jk) 
    623                icoord(ierror,3) = jk 
    624             ENDIF 
    625          END DO 
    626       END DO 
    627        
    628       IF( ierror > 0 ) THEN 
    629          IF(lwp) WRITE(numout,*) 
    630          IF(lwp) WRITE(numout,*) '              Problem on lateral conditions' 
    631          IF(lwp) WRITE(numout,*) '                 Bad marking off at points:' 
    632          DO jl = 1, ierror 
    633             IF(lwp) WRITE(numout,*) 'Level:',icoord(jl,3),   & 
    634                &                  '  Point(',icoord(jl,1),',',icoord(jl,2),')' 
    635          END DO 
    636          CALL ctl_stop( 'We stop...' ) 
    637       ENDIF 
    638       ! 
    639       IF( nn_timing == 1 )  CALL timing_stop('dom_msk_nsa') 
    640       ! 
    641    END SUBROUTINE dom_msk_nsa 
    642  
    643 #else 
    644    !!---------------------------------------------------------------------- 
    645    !!   Default option :                                      Empty routine 
    646    !!---------------------------------------------------------------------- 
    647    SUBROUTINE dom_msk_nsa        
    648    END SUBROUTINE dom_msk_nsa 
    649 #endif 
    650357    
    651358   !!====================================================================== 
Note: See TracChangeset for help on using the changeset viewer.