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 7646 for trunk/NEMOGCM/NEMO/OPA_SRC/DYN/wet_dry.F90 – NEMO

Ignore:
Timestamp:
2017-02-06T10:25:03+01:00 (7 years ago)
Author:
timgraham
Message:

Merge of dev_merge_2016 into trunk. UPDATE TO ARCHFILES NEEDED for XIOS2.
LIM_SRC_s/limrhg.F90 to follow in next commit due to change of kind (I'm unable to do it in this commit).
Merged using the following steps:

1) svn merge --reintegrate svn+ssh://forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/trunk .
2) Resolve minor conflicts in sette.sh and namelist_cfg for ORCA2LIM3 (due to a change in trunk after branch was created)
3) svn commit
4) svn switch svn+ssh://forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/trunk
5) svn merge svn+ssh://forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/branches/2016/dev_merge_2016 .
6) At this stage I checked out a clean copy of the branch to compare against what is about to be committed to the trunk.
6) svn commit #Commit code to the trunk

In this commit I have also reverted a change to Fcheck_archfile.sh which was causing problems on the Paris machine.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/NEMOGCM/NEMO/OPA_SRC/DYN/wet_dry.F90

    r6152 r7646  
    1  
    21MODULE wet_dry 
    32   !!============================================================================== 
     
    76   !! only effects if wetting/drying is on (ln_wd == .true.) 
    87   !!============================================================================== 
    9    !! History :    
    10    !!  NEMO      3.6  ! 2014-09  ((H.Liu)  Original code 
     8   !! History :  3.6  ! 2014-09  ((H.Liu)  Original code 
    119   !!                 ! will add the runoff and periodic BC case later 
    1210   !!---------------------------------------------------------------------- 
     
    3331   !! --------------------------------------------------------------------- 
    3432 
    35    REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) ::   wduflt, wdvflt !: u- and v- filter 
    3633   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) ::   wdmask         !: u- and v- limiter  
     34   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) ::   ht_wd          !: wetting and drying t-pt depths 
     35                                                                     !  (can include negative depths) 
    3736 
    3837   LOGICAL,  PUBLIC  ::   ln_wd       !: Wetting/drying activation switch (T:on,F:off) 
     
    7776         WRITE(numout,*) 
    7877         WRITE(numout,*) 'wad_init : Wetting and drying initialization through namelist read' 
    79          WRITE(numout,*) '~~~~~~~ ' 
     78         WRITE(numout,*) '~~~~~~~~' 
    8079         WRITE(numout,*) '   Namelist namwad' 
    8180         WRITE(numout,*) '      Logical activation                 ln_wd      = ', ln_wd 
     
    8483         WRITE(numout,*) '      land elevation threshold         rn_wdld      = ', rn_wdld 
    8584         WRITE(numout,*) '      Max iteration for W/D limiter    nn_wdit      = ', nn_wdit 
    86        ENDIF 
    87  
     85      ENDIF 
     86      ! 
    8887      IF(ln_wd) THEN 
    89          ALLOCATE( wduflt(jpi,jpj), wdvflt(jpi,jpj), wdmask(jpi,jpj), STAT=ierr ) 
     88         ALLOCATE( wdmask(jpi,jpj), ht_wd(jpi,jpj), STAT=ierr ) 
    9089         IF( ierr /= 0 ) CALL ctl_stop('STOP', 'wad_init : Array allocation error') 
    9190      ENDIF 
     91      ! 
    9292   END SUBROUTINE wad_init 
     93 
    9394 
    9495   SUBROUTINE wad_lmt( sshb1, sshemp, z2dt ) 
     
    107108      ! 
    108109      INTEGER  ::   ji, jj, jk, jk1     ! dummy loop indices 
    109       INTEGER  ::   zflag               ! local scalar 
     110      INTEGER  ::   jflag               ! local scalar 
    110111      REAL(wp) ::   zcoef, zdep1, zdep2 ! local scalars 
    111112      REAL(wp) ::   zzflxp, zzflxn      ! local scalars 
     
    116117      REAL(wp), POINTER,  DIMENSION(:,:) ::   zflxu,  zflxv            ! local 2D workspace 
    117118      REAL(wp), POINTER,  DIMENSION(:,:) ::   zflxu1, zflxv1           ! local 2D workspace 
    118  
    119119      !!---------------------------------------------------------------------- 
    120120      ! 
     
    131131        !IF(lwp) WRITE(numout,*) 'wad_lmt : wetting/drying limiters and velocity limiting' 
    132132        
    133         zflag  = 0 
     133        jflag  = 0 
    134134        zdepwd = 50._wp   !maximum depth on which that W/D could possibly happen 
    135135 
     
    145145        ! Horizontal Flux in u and v direction 
    146146        DO jk = 1, jpkm1   
    147            DO jj = 1, jpjm1 
    148               DO ji = 1, jpim1 
     147           DO jj = 1, jpj 
     148              DO ji = 1, jpi 
    149149                 zflxu(ji,jj) = zflxu(ji,jj) + e3u_n(ji,jj,jk) * un(ji,jj,jk) * umask(ji,jj,jk) 
    150150                 zflxv(ji,jj) = zflxv(ji,jj) + e3v_n(ji,jj,jk) * vn(ji,jj,jk) * vmask(ji,jj,jk) 
     
    156156        zflxv(:,:) = zflxv(:,:) * e1v(:,:) 
    157157        
    158         DO jj = 2, jpjm1 
    159            DO ji = 2, jpim1  
    160  
    161              IF(tmask(ji, jj, 1) < 0.5_wp) CYCLE   ! we don't care about land cells 
    162              IF(bathy(ji,jj) > zdepwd) CYCLE       ! and cells which will unlikely go dried out 
     158        wdmask(:,:) = 1 
     159        DO jj = 2, jpj 
     160           DO ji = 2, jpi  
     161 
     162             IF( tmask(ji, jj, 1) < 0.5_wp ) CYCLE   ! we don't care about land cells 
     163             IF( ht_wd(ji,jj) > zdepwd )      CYCLE   ! and cells which are unlikely to dry 
    163164 
    164165              zflxp(ji,jj) = max(zflxu(ji,jj), 0._wp) - min(zflxu(ji-1,jj),   0._wp) + & 
     
    167168                           & min(zflxv(ji,jj), 0._wp) - max(zflxv(ji,  jj-1), 0._wp)  
    168169 
    169               zdep2 = bathy(ji,jj) + sshb1(ji,jj) - rn_wdmin1 
    170               IF(zdep2 < 0._wp) THEN  !add more safty, but not necessary 
    171                 !zdep2 = 0._wp 
    172                 sshb1(ji,jj) = rn_wdmin1 - bathy(ji,jj) 
     170              zdep2 = ht_wd(ji,jj) + sshb1(ji,jj) - rn_wdmin1 
     171              IF(zdep2 .le. 0._wp) THEN  !add more safty, but not necessary 
     172                sshb1(ji,jj) = rn_wdmin1 - ht_wd(ji,jj) 
     173                IF(zflxu(ji,  jj) > 0._wp) zwdlmtu(ji  ,jj) = 0._wp 
     174                IF(zflxu(ji-1,jj) < 0._wp) zwdlmtu(ji-1,jj) = 0._wp 
     175                IF(zflxv(ji,  jj) > 0._wp) zwdlmtv(ji  ,jj) = 0._wp 
     176                IF(zflxv(ji,jj-1) < 0._wp) zwdlmtv(ji,jj-1) = 0._wp  
     177                wdmask(ji,jj) = 0._wp 
    173178              END IF 
    174179           ENDDO 
     
    182187           zflxu1(:,:) = zflxu(:,:) * zwdlmtu(:,:) 
    183188           zflxv1(:,:) = zflxv(:,:) * zwdlmtv(:,:) 
    184            
    185            DO jj = 2, jpjm1 
    186               DO ji = 2, jpim1  
     189           jflag = 0     ! flag indicating if any further iterations are needed 
     190           
     191           DO jj = 2, jpj 
     192              DO ji = 2, jpi  
    187193         
    188                  wdmask(ji,jj) = 0 
    189                  IF(tmask(ji, jj, 1) < 0.5_wp) CYCLE  
    190                  IF(bathy(ji,jj) > zdepwd) CYCLE 
     194                 IF( tmask(ji, jj, 1) < 0.5_wp ) CYCLE  
     195                 IF( ht_wd(ji,jj) > zdepwd )      CYCLE 
    191196         
    192197                 ztmp = e1e2t(ji,jj) 
     
    198203           
    199204                 zdep1 = (zzflxp + zzflxn) * z2dt / ztmp 
    200                  zdep2 = bathy(ji,jj) + sshb1(ji,jj) - rn_wdmin1 - z2dt * sshemp(ji,jj)  ! this one can be moved out of the loop 
    201            
    202                  IF(zdep1 > zdep2) THEN 
    203                    zflag = 1 
    204                    wdmask(ji, jj) = 1 
     205                 zdep2 = ht_wd(ji,jj) + sshb1(ji,jj) - rn_wdmin1 - z2dt * sshemp(ji,jj) 
     206           
     207                 IF( zdep1 > zdep2 ) THEN 
     208                   wdmask(ji, jj) = 0 
    205209                   zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zflxp(ji,jj) * z2dt ) 
    206                    zcoef = max(zcoef, 0._wp) 
     210                   !zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zzflxp * z2dt ) 
     211                   ! flag if the limiter has been used but stop flagging if the only 
     212                   ! changes have zeroed the coefficient since further iterations will 
     213                   ! not change anything 
     214                   IF( zcoef > 0._wp ) THEN 
     215                      jflag = 1  
     216                   ELSE 
     217                      zcoef = 0._wp 
     218                   ENDIF 
    207219                   IF(jk1 > nn_wdit) zcoef = 0._wp 
    208220                   IF(zflxu1(ji,  jj) > 0._wp) zwdlmtu(ji  ,jj) = zcoef 
    209221                   IF(zflxu1(ji-1,jj) < 0._wp) zwdlmtu(ji-1,jj) = zcoef 
    210222                   IF(zflxv1(ji,  jj) > 0._wp) zwdlmtv(ji  ,jj) = zcoef 
    211                    IF(zflxv1(ji,jj-1) < 0._wp) zwdlmtv(ji-1,jj) = zcoef 
     223                   IF(zflxv1(ji,jj-1) < 0._wp) zwdlmtv(ji,jj-1) = zcoef 
    212224                 END IF 
    213225              END DO ! ji loop 
     
    217229           CALL lbc_lnk( zwdlmtv, 'V', 1. ) 
    218230 
    219            IF(lk_mpp) CALL mpp_max(zflag)   !max over the global domain 
    220  
    221            IF(zflag == 0) EXIT 
    222            
    223            zflag = 0     ! flag indicating if any further iteration is needed? 
     231           IF(lk_mpp) CALL mpp_max(jflag)   !max over the global domain 
     232 
     233           IF(jflag == 0) EXIT 
     234           
    224235        END DO  ! jk1 loop 
    225236        
     
    231242        CALL lbc_lnk( un, 'U', -1. ) 
    232243        CALL lbc_lnk( vn, 'V', -1. ) 
    233         
    234         IF(zflag == 1 .AND. lwp) WRITE(numout,*) 'Need more iterations in wad_lmt!!!' 
     244      ! 
     245        un_b(:,:) = un_b(:,:) * zwdlmtu(:, :) 
     246        vn_b(:,:) = vn_b(:,:) * zwdlmtv(:, :) 
     247        CALL lbc_lnk( un_b, 'U', -1. ) 
     248        CALL lbc_lnk( vn_b, 'V', -1. ) 
     249        
     250        IF(jflag == 1 .AND. lwp) WRITE(numout,*) 'Need more iterations in wad_lmt!!!' 
    235251        
    236252        !IF( ln_rnf      )   CALL sbc_rnf_div( hdivn )          ! runoffs (update hdivn field) 
     
    240256        CALL wrk_dealloc( jpi, jpj, zflxp, zflxn, zflxu, zflxv, zflxu1, zflxv1 ) 
    241257        CALL wrk_dealloc( jpi, jpj, zwdlmtu, zwdlmtv) 
    242       ! 
    243       END IF 
    244  
     258        ! 
     259      ENDIF 
     260      ! 
    245261      IF( nn_timing == 1 )  CALL timing_stop('wad_lmt') 
     262      ! 
    246263   END SUBROUTINE wad_lmt 
     264 
    247265 
    248266   SUBROUTINE wad_lmt_bt( zflxu, zflxv, sshn_e, zssh_frc, rdtbt ) 
     
    260278      ! 
    261279      INTEGER  ::   ji, jj, jk, jk1     ! dummy loop indices 
    262       INTEGER  ::   zflag         ! local scalar 
     280      INTEGER  ::   jflag               ! local scalar 
    263281      REAL(wp) ::   z2dt 
    264282      REAL(wp) ::   zcoef, zdep1, zdep2 ! local scalars 
     
    269287      REAL(wp), POINTER,  DIMENSION(:,:) ::   zflxp,  zflxn            ! local 2D workspace 
    270288      REAL(wp), POINTER,  DIMENSION(:,:) ::   zflxu1, zflxv1           ! local 2D workspace 
    271  
    272       !!---------------------------------------------------------------------- 
    273       ! 
    274  
     289      !!---------------------------------------------------------------------- 
     290      ! 
    275291      IF( nn_timing == 1 )  CALL timing_start('wad_lmt_bt') 
    276292 
     
    284300        !IF(lwp) WRITE(numout,*) 'wad_lmt_bt : wetting/drying limiters and velocity limiting' 
    285301        
    286         zflag  = 0 
     302        jflag  = 0 
    287303        zdepwd = 50._wp   !maximum depth that ocean cells can have W/D processes 
    288304 
     
    291307        zflxp(:,:)   = 0._wp 
    292308        zflxn(:,:)   = 0._wp 
    293         !zflxu(:,:)   = 0._wp 
    294         !zflxv(:,:)   = 0._wp 
    295309 
    296310        zwdlmtu(:,:)  = 1._wp 
     
    299313        ! Horizontal Flux in u and v direction 
    300314        
    301         !zflxu(:,:) = zflxu(:,:) * e2u(:,:) 
    302         !zflxv(:,:) = zflxv(:,:) * e1v(:,:) 
    303         
    304         DO jj = 2, jpjm1 
    305            DO ji = 2, jpim1  
    306  
    307              IF(tmask(ji, jj, 1) < 0.5_wp) CYCLE   ! we don't care about land cells 
    308              IF(bathy(ji,jj) > zdepwd) CYCLE       ! and cells which will unlikely go dried out 
     315        DO jj = 2, jpj 
     316           DO ji = 2, jpi  
     317 
     318             IF( tmask(ji, jj, 1 ) < 0.5_wp) CYCLE   ! we don't care about land cells 
     319             IF( ht_wd(ji,jj) > zdepwd )      CYCLE   ! and cells which are unlikely to dry 
    309320 
    310321              zflxp(ji,jj) = max(zflxu(ji,jj), 0._wp) - min(zflxu(ji-1,jj),   0._wp) + & 
     
    313324                           & min(zflxv(ji,jj), 0._wp) - max(zflxv(ji,  jj-1), 0._wp)  
    314325 
    315               zdep2 = bathy(ji,jj) + sshn_e(ji,jj) - rn_wdmin1 
    316               IF(zdep2 < 0._wp) THEN  !add more safty, but not necessary 
    317                 !zdep2 = 0._wp 
    318                sshn_e(ji,jj) = rn_wdmin1 - bathy(ji,jj) 
     326              zdep2 = ht_wd(ji,jj) + sshn_e(ji,jj) - rn_wdmin1 
     327              IF(zdep2 .le. 0._wp) THEN  !add more safety, but not necessary 
     328                sshn_e(ji,jj) = rn_wdmin1 - ht_wd(ji,jj) 
     329                IF(zflxu(ji,  jj) > 0._wp) zwdlmtu(ji  ,jj) = 0._wp 
     330                IF(zflxu(ji-1,jj) < 0._wp) zwdlmtu(ji-1,jj) = 0._wp 
     331                IF(zflxv(ji,  jj) > 0._wp) zwdlmtv(ji  ,jj) = 0._wp 
     332                IF(zflxv(ji,jj-1) < 0._wp) zwdlmtv(ji,jj-1) = 0._wp  
    319333              END IF 
    320334           ENDDO 
     
    328342           zflxu1(:,:) = zflxu(:,:) * zwdlmtu(:,:) 
    329343           zflxv1(:,:) = zflxv(:,:) * zwdlmtv(:,:) 
    330            
    331            DO jj = 2, jpjm1 
    332               DO ji = 2, jpim1  
     344           jflag = 0     ! flag indicating if any further iterations are needed 
     345           
     346           DO jj = 2, jpj 
     347              DO ji = 2, jpi  
    333348         
    334                  wdmask(ji,jj) = 0 
    335                  IF(tmask(ji, jj, 1) < 0.5_wp) CYCLE  
    336                  IF(bathy(ji,jj) > zdepwd) CYCLE 
     349                 IF( tmask(ji, jj, 1 ) < 0.5_wp) CYCLE  
     350                 IF( ht_wd(ji,jj) > zdepwd )      CYCLE 
    337351         
    338352                 ztmp = e1e2t(ji,jj) 
     
    344358           
    345359                 zdep1 = (zzflxp + zzflxn) * z2dt / ztmp 
    346                  zdep2 = bathy(ji,jj) + sshn_e(ji,jj) - rn_wdmin1   ! this one can be moved out of the loop 
    347                  zdep2 = zdep2 - z2dt * zssh_frc(ji,jj) 
     360                 zdep2 = ht_wd(ji,jj) + sshn_e(ji,jj) - rn_wdmin1 - z2dt * zssh_frc(ji,jj) 
    348361           
    349362                 IF(zdep1 > zdep2) THEN 
    350                    zflag = 1 
    351                    !wdmask(ji, jj) = 1 
    352363                   zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zflxp(ji,jj) * z2dt ) 
    353                    zcoef = max(zcoef, 0._wp) 
     364                   !zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zzflxp * z2dt ) 
     365                   ! flag if the limiter has been used but stop flagging if the only 
     366                   ! changes have zeroed the coefficient since further iterations will 
     367                   ! not change anything 
     368                   IF( zcoef > 0._wp ) THEN 
     369                      jflag = 1  
     370                   ELSE 
     371                      zcoef = 0._wp 
     372                   ENDIF 
    354373                   IF(jk1 > nn_wdit) zcoef = 0._wp 
    355374                   IF(zflxu1(ji,  jj) > 0._wp) zwdlmtu(ji  ,jj) = zcoef 
    356375                   IF(zflxu1(ji-1,jj) < 0._wp) zwdlmtu(ji-1,jj) = zcoef 
    357376                   IF(zflxv1(ji,  jj) > 0._wp) zwdlmtv(ji  ,jj) = zcoef 
    358                    IF(zflxv1(ji,jj-1) < 0._wp) zwdlmtv(ji-1,jj) = zcoef 
     377                   IF(zflxv1(ji,jj-1) < 0._wp) zwdlmtv(ji,jj-1) = zcoef 
    359378                 END IF 
    360379              END DO ! ji loop 
     
    364383           CALL lbc_lnk( zwdlmtv, 'V', 1. ) 
    365384 
    366            IF(lk_mpp) CALL mpp_max(zflag)   !max over the global domain 
    367  
    368            IF(zflag == 0) EXIT 
    369            
    370            zflag = 0     ! flag indicating if any further iteration is needed? 
     385           IF(lk_mpp) CALL mpp_max(jflag)   !max over the global domain 
     386 
     387           IF(jflag == 0) EXIT 
     388           
    371389        END DO  ! jk1 loop 
    372390        
     
    377395        CALL lbc_lnk( zflxv, 'V', -1. ) 
    378396        
    379         IF(zflag == 1 .AND. lwp) WRITE(numout,*) 'Need more iterations in wad_lmt_bt!!!' 
     397        IF(jflag == 1 .AND. lwp) WRITE(numout,*) 'Need more iterations in wad_lmt_bt!!!' 
    380398        
    381399        !IF( ln_rnf      )   CALL sbc_rnf_div( hdivn )          ! runoffs (update hdivn field) 
     
    385403        CALL wrk_dealloc( jpi, jpj, zflxp, zflxn, zflxu1, zflxv1 ) 
    386404        CALL wrk_dealloc( jpi, jpj, zwdlmtu, zwdlmtv) 
    387       ! 
     405        ! 
    388406      END IF 
    389407 
    390408      IF( nn_timing == 1 )  CALL timing_stop('wad_lmt') 
    391409   END SUBROUTINE wad_lmt_bt 
     410 
     411   !!============================================================================== 
    392412END MODULE wet_dry 
Note: See TracChangeset for help on using the changeset viewer.