Ignore:
Timestamp:
2015-10-26T15:49:40+01:00 (5 years ago)
Author:
cetlod
Message:

merge the simplification branch onto the trunk, see ticket #1612

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/NEMOGCM/NEMO/OPA_SRC/DOM/dommsk.F90

    r5552 r5836  
    1717   !!             -   ! 2005-11  (V. Garnier) Surface pressure gradient organization 
    1818   !!            3.2  ! 2009-07  (R. Benshila) Suppression of rigid-lid option 
     19   !!            3.6  ! 2015-05  (P. Mathiot) ISF: add wmask,wumask and wvmask 
    1920   !!---------------------------------------------------------------------- 
    2021 
    2122   !!---------------------------------------------------------------------- 
    2223   !!   dom_msk        : compute land/ocean mask 
    23    !!   dom_msk_nsa    : update land/ocean mask when no-slip accurate option is used. 
    2424   !!---------------------------------------------------------------------- 
    2525   USE oce             ! ocean dynamics and tracers 
     
    3636 
    3737   PUBLIC   dom_msk         ! routine called by inidom.F90 
    38    PUBLIC   dom_msk_alloc   ! routine called by nemogcm.F90 
    3938 
    4039   !                            !!* Namelist namlbc : lateral boundary condition * 
     
    4241   LOGICAL, PUBLIC :: ln_vorlat  !  consistency of vorticity boundary condition  
    4342   !                                            with analytical eqs. 
    44  
    45  
    46    INTEGER, ALLOCATABLE, SAVE, DIMENSION(:,:) ::  icoord ! Workspace for dom_msk_nsa() 
    4743 
    4844   !! * Substitutions 
     
    5450   !!---------------------------------------------------------------------- 
    5551CONTAINS 
    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  
    6952 
    7053   SUBROUTINE dom_msk 
     
    129112      !!               tmask_i  : interior ocean mask 
    130113      !!---------------------------------------------------------------------- 
    131       ! 
    132       INTEGER  ::   ji, jj, jk      ! dummy loop indices 
     114      INTEGER  ::   ji, jj, jk               ! dummy loop indices 
    133115      INTEGER  ::   iif, iil, ii0, ii1, ii   ! local integers 
    134116      INTEGER  ::   ijf, ijl, ij0, ij1       !   -       - 
     
    199181      END DO   
    200182 
    201 !!gm  ???? 
    202 #if defined key_zdfkpp 
    203       IF( cp_cfg == 'orca' ) THEN 
    204          IF( jp_cfg == 2 )   THEN       ! land point on Bab el Mandeb zonal section 
    205             ij0 =  87   ;   ij1 =  88 
    206             ii0 = 160   ;   ii1 = 161 
    207             tmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) = 0._wp 
    208          ELSE 
    209             IF(lwp) WRITE(numout,*) 
    210             IF(lwp) WRITE(numout,cform_war) 
    211             IF(lwp) WRITE(numout,*) 
    212             IF(lwp) WRITE(numout,*)'          A mask must be applied on Bab el Mandeb strait' 
    213             IF(lwp) WRITE(numout,*)'          in case of ORCAs configurations' 
    214             IF(lwp) WRITE(numout,*)'          This is a problem which is not yet solved' 
    215             IF(lwp) WRITE(numout,*) 
    216          ENDIF 
    217       ENDIF 
    218 #endif 
    219 !!gm end 
    220  
    221183      ! Interior domain mask (used for global sum) 
    222184      ! -------------------- 
     
    284246      ! 3. Ocean/land mask at wu-, wv- and w points  
    285247      !---------------------------------------------- 
    286       wmask (:,:,1) = tmask(:,:,1) ! ???????? 
    287       wumask(:,:,1) = umask(:,:,1) ! ???????? 
    288       wvmask(:,:,1) = vmask(:,:,1) ! ???????? 
    289       DO jk=2,jpk 
    290          wmask (:,:,jk)=tmask(:,:,jk) * tmask(:,:,jk-1) 
    291          wumask(:,:,jk)=umask(:,:,jk) * umask(:,:,jk-1)    
    292          wvmask(:,:,jk)=vmask(:,:,jk) * vmask(:,:,jk-1) 
     248      wmask (:,:,1) = tmask(:,:,1)     ! surface 
     249      wumask(:,:,1) = umask(:,:,1) 
     250      wvmask(:,:,1) = vmask(:,:,1) 
     251      DO jk = 2, jpk                   ! interior values 
     252         wmask (:,:,jk) = tmask(:,:,jk) * tmask(:,:,jk-1) 
     253         wumask(:,:,jk) = umask(:,:,jk) * umask(:,:,jk-1)    
     254         wvmask(:,:,jk) = vmask(:,:,jk) * vmask(:,:,jk-1) 
    293255      END DO 
    294256 
     
    339301      ENDIF 
    340302 
    341  
    342       ! mask for second order calculation of vorticity 
    343       ! ---------------------------------------------- 
    344       CALL dom_msk_nsa 
    345  
    346        
    347303      ! Lateral boundary conditions on velocity (modify fmask) 
    348304      ! ---------------------------------------      
     
    377333      IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN   ! ORCA_R2 configuration 
    378334         !                                                 ! Increased lateral friction near of some straits 
    379          IF( nn_cla == 0 ) THEN 
    380             !                                ! Gibraltar strait  : partial slip (fmask=0.5) 
    381             ij0 = 101   ;   ij1 = 101 
    382             ii0 = 139   ;   ii1 = 140   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  0.5_wp 
    383             ij0 = 102   ;   ij1 = 102 
    384             ii0 = 139   ;   ii1 = 140   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  0.5_wp 
    385             ! 
    386             !                                ! Bab el Mandeb : partial slip (fmask=1) 
    387             ij0 =  87   ;   ij1 =  88 
    388             ii0 = 160   ;   ii1 = 160   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  1._wp 
    389             ij0 =  88   ;   ij1 =  88 
    390             ii0 = 159   ;   ii1 = 159   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  1._wp 
    391             ! 
    392          ENDIF 
     335         !                                ! Gibraltar strait  : partial slip (fmask=0.5) 
     336         ij0 = 101   ;   ij1 = 101 
     337         ii0 = 139   ;   ii1 = 140   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  0.5_wp 
     338         ij0 = 102   ;   ij1 = 102 
     339         ii0 = 139   ;   ii1 = 140   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  0.5_wp 
     340         ! 
     341         !                                ! Bab el Mandeb : partial slip (fmask=1) 
     342         ij0 =  87   ;   ij1 =  88 
     343         ii0 = 160   ;   ii1 = 160   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  1._wp 
     344         ij0 =  88   ;   ij1 =  88 
     345         ii0 = 159   ;   ii1 = 159   ;   fmask( mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) , 1:jpk ) =  1._wp 
     346         ! 
    393347         !                                ! Danish straits  : strong slip (fmask > 2) 
    394348! We keep this as an example but it is instable in this case  
     
    500454      ! 
    501455   END SUBROUTINE dom_msk 
    502  
    503 #if defined key_noslip_accurate 
    504    !!---------------------------------------------------------------------- 
    505    !!   'key_noslip_accurate' :         accurate no-slip boundary condition 
    506    !!---------------------------------------------------------------------- 
    507     
    508    SUBROUTINE dom_msk_nsa 
    509       !!--------------------------------------------------------------------- 
    510       !!                 ***  ROUTINE dom_msk_nsa  *** 
    511       !!  
    512       !! ** Purpose : 
    513       !! 
    514       !! ** Method  : 
    515       !! 
    516       !! ** Action : 
    517       !!---------------------------------------------------------------------- 
    518       INTEGER  ::   ji, jj, jk, jl      ! dummy loop indices 
    519       INTEGER  ::   ine, inw, ins, inn, itest, ierror, iind, ijnd 
    520       REAL(wp) ::   zaa 
    521       !!--------------------------------------------------------------------- 
    522       ! 
    523       IF( nn_timing == 1 )  CALL timing_start('dom_msk_nsa') 
    524       ! 
    525       IF(lwp) WRITE(numout,*) 
    526       IF(lwp) WRITE(numout,*) 'dom_msk_nsa : noslip accurate boundary condition' 
    527       IF(lwp) WRITE(numout,*) '~~~~~~~~~~~   using Schchepetkin and O Brian scheme' 
    528       IF( lk_mpp )   CALL ctl_stop( ' mpp version is not yet implemented' ) 
    529  
    530       ! mask for second order calculation of vorticity 
    531       ! ---------------------------------------------- 
    532       ! noslip boundary condition: fmask=1  at convex corner, store 
    533       ! index of straight coast meshes ( 'west', refering to a coast, 
    534       ! means west of the ocean, aso) 
    535        
    536       DO jk = 1, jpk 
    537          DO jl = 1, 4 
    538             npcoa(jl,jk) = 0 
    539             DO ji = 1, 2*(jpi+jpj) 
    540                nicoa(ji,jl,jk) = 0 
    541                njcoa(ji,jl,jk) = 0 
    542             END DO 
    543          END DO 
    544       END DO 
    545        
    546       IF( jperio == 2 ) THEN 
    547          WRITE(numout,*) ' ' 
    548          WRITE(numout,*) ' symetric boundary conditions need special' 
    549          WRITE(numout,*) ' treatment not implemented. we stop.' 
    550          STOP 
    551       ENDIF 
    552        
    553       ! convex corners 
    554        
    555       DO jk = 1, jpkm1 
    556          DO jj = 1, jpjm1 
    557             DO ji = 1, jpim1 
    558                zaa = tmask(ji  ,jj,jk) + tmask(ji  ,jj+1,jk)   & 
    559                   &+ tmask(ji+1,jj,jk) + tmask(ji+1,jj+1,jk) 
    560                IF( ABS(zaa-3._wp) <= 0.1_wp )   fmask(ji,jj,jk) = 1._wp 
    561             END DO 
    562          END DO 
    563       END DO 
    564  
    565       ! north-south straight coast 
    566  
    567       DO jk = 1, jpkm1 
    568          inw = 0 
    569          ine = 0 
    570          DO jj = 2, jpjm1 
    571             DO ji = 2, jpim1 
    572                zaa = tmask(ji+1,jj,jk) + tmask(ji+1,jj+1,jk) 
    573                IF( ABS(zaa-2._wp) <= 0.1_wp .AND. fmask(ji,jj,jk) == 0._wp ) THEN 
    574                   inw = inw + 1 
    575                   nicoa(inw,1,jk) = ji 
    576                   njcoa(inw,1,jk) = jj 
    577                   IF( nprint == 1 ) WRITE(numout,*) ' west  : ', jk, inw, ji, jj 
    578                ENDIF 
    579                zaa = tmask(ji,jj,jk) + tmask(ji,jj+1,jk) 
    580                IF( ABS(zaa-2._wp) <= 0.1_wp .AND. fmask(ji,jj,jk) == 0._wp ) THEN 
    581                   ine = ine + 1 
    582                   nicoa(ine,2,jk) = ji 
    583                   njcoa(ine,2,jk) = jj 
    584                   IF( nprint == 1 ) WRITE(numout,*) ' east  : ', jk, ine, ji, jj 
    585                ENDIF 
    586             END DO 
    587          END DO 
    588          npcoa(1,jk) = inw 
    589          npcoa(2,jk) = ine 
    590       END DO 
    591  
    592       ! west-east straight coast 
    593  
    594       DO jk = 1, jpkm1 
    595          ins = 0 
    596          inn = 0 
    597          DO jj = 2, jpjm1 
    598             DO ji =2, jpim1 
    599                zaa = tmask(ji,jj+1,jk) + tmask(ji+1,jj+1,jk) 
    600                IF( ABS(zaa-2._wp) <= 0.1_wp .AND. fmask(ji,jj,jk) == 0._wp ) THEN 
    601                   ins = ins + 1 
    602                   nicoa(ins,3,jk) = ji 
    603                   njcoa(ins,3,jk) = jj 
    604                   IF( nprint == 1 ) WRITE(numout,*) ' south : ', jk, ins, ji, jj 
    605                ENDIF 
    606                zaa = tmask(ji+1,jj,jk) + tmask(ji,jj,jk) 
    607                IF( ABS(zaa-2._wp) <= 0.1_wp .AND. fmask(ji,jj,jk) == 0._wp ) THEN 
    608                   inn = inn + 1 
    609                   nicoa(inn,4,jk) = ji 
    610                   njcoa(inn,4,jk) = jj 
    611                   IF( nprint == 1 ) WRITE(numout,*) ' north : ', jk, inn, ji, jj 
    612                ENDIF 
    613             END DO 
    614          END DO 
    615          npcoa(3,jk) = ins 
    616          npcoa(4,jk) = inn 
    617       END DO 
    618  
    619       itest = 2 * ( jpi + jpj ) 
    620       DO jk = 1, jpk 
    621          IF( npcoa(1,jk) > itest .OR. npcoa(2,jk) > itest .OR.   & 
    622              npcoa(3,jk) > itest .OR. npcoa(4,jk) > itest ) THEN 
    623              
    624             WRITE(ctmp1,*) ' level jk = ',jk 
    625             WRITE(ctmp2,*) ' straight coast index arraies are too small.:' 
    626             WRITE(ctmp3,*) ' npe, npw, nps, npn = ', npcoa(1,jk), npcoa(2,jk),   & 
    627                 &                                     npcoa(3,jk), npcoa(4,jk) 
    628             WRITE(ctmp4,*) ' 2*(jpi+jpj) = ',itest,'. we stop.' 
    629             CALL ctl_stop( ctmp1, ctmp2, ctmp3, ctmp4 ) 
    630         ENDIF 
    631       END DO 
    632  
    633       ierror = 0 
    634       iind = 0 
    635       ijnd = 0 
    636       IF( nperio == 1 .OR. nperio == 4 .OR. nperio == 6 )   iind = 2 
    637       IF( nperio == 3 .OR. nperio == 4 .OR. nperio == 5 .OR. nperio == 6 )   ijnd = 2 
    638       DO jk = 1, jpk 
    639          DO jl = 1, npcoa(1,jk) 
    640             IF( nicoa(jl,1,jk)+3 > jpi+iind ) THEN 
    641                ierror = ierror+1 
    642                icoord(ierror,1) = nicoa(jl,1,jk) 
    643                icoord(ierror,2) = njcoa(jl,1,jk) 
    644                icoord(ierror,3) = jk 
    645             ENDIF 
    646          END DO 
    647          DO jl = 1, npcoa(2,jk) 
    648             IF(nicoa(jl,2,jk)-2 < 1-iind ) THEN 
    649                ierror = ierror + 1 
    650                icoord(ierror,1) = nicoa(jl,2,jk) 
    651                icoord(ierror,2) = njcoa(jl,2,jk) 
    652                icoord(ierror,3) = jk 
    653             ENDIF 
    654          END DO 
    655          DO jl = 1, npcoa(3,jk) 
    656             IF( njcoa(jl,3,jk)+3 > jpj+ijnd ) THEN 
    657                ierror = ierror + 1 
    658                icoord(ierror,1) = nicoa(jl,3,jk) 
    659                icoord(ierror,2) = njcoa(jl,3,jk) 
    660                icoord(ierror,3) = jk 
    661             ENDIF 
    662          END DO 
    663          DO jl = 1, npcoa(4,jk) 
    664             IF( njcoa(jl,4,jk)-2 < 1) THEN 
    665                ierror=ierror + 1 
    666                icoord(ierror,1) = nicoa(jl,4,jk) 
    667                icoord(ierror,2) = njcoa(jl,4,jk) 
    668                icoord(ierror,3) = jk 
    669             ENDIF 
    670          END DO 
    671       END DO 
    672        
    673       IF( ierror > 0 ) THEN 
    674          IF(lwp) WRITE(numout,*) 
    675          IF(lwp) WRITE(numout,*) '              Problem on lateral conditions' 
    676          IF(lwp) WRITE(numout,*) '                 Bad marking off at points:' 
    677          DO jl = 1, ierror 
    678             IF(lwp) WRITE(numout,*) 'Level:',icoord(jl,3),   & 
    679                &                  '  Point(',icoord(jl,1),',',icoord(jl,2),')' 
    680          END DO 
    681          CALL ctl_stop( 'We stop...' ) 
    682       ENDIF 
    683       ! 
    684       IF( nn_timing == 1 )  CALL timing_stop('dom_msk_nsa') 
    685       ! 
    686    END SUBROUTINE dom_msk_nsa 
    687  
    688 #else 
    689    !!---------------------------------------------------------------------- 
    690    !!   Default option :                                      Empty routine 
    691    !!---------------------------------------------------------------------- 
    692    SUBROUTINE dom_msk_nsa        
    693    END SUBROUTINE dom_msk_nsa 
    694 #endif 
    695456    
    696457   !!====================================================================== 
Note: See TracChangeset for help on using the changeset viewer.