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 7514 for branches – NEMO

Changeset 7514 for branches


Ignore:
Timestamp:
2016-12-20T15:40:04+01:00 (7 years ago)
Author:
acc
Message:

Branch dev_merge_2016. Tidying up wetting and drying test cases. Only usrdef routines are now required in the MY_SRC directory and all test cases will run. Still some robustness issues to resolve with test case 6

Location:
branches/2016/dev_merge_2016/NEMOGCM
Files:
8 deleted
11 edited

Legend:

Unmodified
Added
Removed
  • branches/2016/dev_merge_2016/NEMOGCM/CONFIG/WAD_TEST_CASES/EXP00/README.py

    r7412 r7514  
    1 ncks -d y,17,17,1 WAD_1ts_00010101_00010101_grid_T.nc sshtime.nc 
    2 ncks -4 -A -v gdepw_0,bathy -d y,17,17,1 -d z,9,9,1 -C mesh_mask_wad1.nc sshtime.nc 
     1ncks -O -d y,17,17,1 WAD_1ts_00010101_00010101_grid_T.nc sshtime.nc 
     2ncks -4 -A -v gdepw_0,ht_wd -d y,17,17,1 -d z,10,10,1 -C mesh_mask.nc sshtime.nc 
     3 
     4python2.7 matpoly2.py sshtime.nc 
     5animate wadfr*.png 
  • branches/2016/dev_merge_2016/NEMOGCM/CONFIG/WAD_TEST_CASES/EXP00/matpoly2.py

    r7412 r7514  
    2121fw = Dataset(pathin) 
    2222ssh = fw.variables['sossheig'][:,:,:] 
    23 #bat = fw.variables['gdepw_0'][0,0,:,:] 
    24 bat = fw.variables['bathy'][0,:,:] 
     23bat = fw.variables['gdepw_0'][0,0,:,:] 
     24#bat = fw.variables['ht_wd'][0,:,:] 
    2525#vot = fw.variables['votemper'][:,0,:,:] 
    2626vot = fw.variables['sosaline'][:,:,:] 
     
    4949 dy   = np.int(t/24.0) 
    5050 
    51  tfac = 12.0/3600.0 
     51 tfac = 18.0/3600.0 
    5252 t24  = np.int(np.mod(t*tfac,24)) 
    5353 dy   = np.int(t*tfac/24.0) 
     
    112112 p.set_array(np.array(colors)) 
    113113 
     114 ax.set_ylim([-10., 4.0]) 
     115 ax.set_xlim([0., 50.0]) 
    114116 ax.add_collection(p) 
    115117 ax.plot(xs,ys, '--', color='black', ms=10) 
  • branches/2016/dev_merge_2016/NEMOGCM/CONFIG/WAD_TEST_CASES/MY_SRC/usrdef_istate.F90

    r7467 r7514  
    1414   !!---------------------------------------------------------------------- 
    1515   USE par_oce        ! ocean space and time domain 
    16    USE dom_oce , ONLY : mi0, mig, mjg, glamt, ht_0 
     16   USE dom_oce , ONLY : mi0, mig, mjg, glamt, gphit, ht_0 
    1717   USE phycst         ! physical constants 
    1818   USE wet_dry        ! Wetting and drying 
     
    7474      !!---------------------------------------------------------------------- 
    7575      ! 
    76       ! Uniform T & S in all test cases 
     76      ! Uniform T & S in most test cases 
    7777      pts(:,:,:,jp_tem) = 10._wp 
    78       !tsb(:,:,:,jp_tem) = 10._wp 
    7978      pts(:,:,:,jp_sal) = 35._wp 
    80       !tsb(:,:,:,jp_sal) = 35._wp 
    8179      SELECT CASE ( nn_cfg )  
    8280         !                                        ! ==================== 
     
    9088            do ji = 1,jpi 
    9189             pssh(ji,:) = ( -5.5_wp + 5.5_wp*glamt(ji,1)/50._wp)*ptmask(ji,:,1) 
    92              !write(*,*) 'SSH ',jpk,mig(ji),pssh(ji,1),ht_0(ji,1),-5.5_wp+5.5_wp*glamt(ji,1)/50._wp,ptmask(ji,1,1) 
    9390            end do 
    9491            !                                     ! ==================== 
     
    10198            ! 
    10299            do ji = 1,jpi 
    103              pssh(ji,:) = ( -5.5_wp + 3.9_wp*FLOAT(jpiglo - mig(ji))/FLOAT(jpiglo-1))*ptmask(ji,:,1) 
     100             pssh(ji,:) = ( -5.5_wp + 6.5_wp*glamt(ji,1)/50._wp)*ptmask(ji,:,1) 
    104101            end do 
    105102            !                                     ! ==================== 
     
    112109            ! 
    113110            do ji = 1,jpi 
    114              pssh(ji,:) = ( -7.5_wp + 6.9_wp*FLOAT(jpiglo - mig(ji))/FLOAT(jpiglo-1))*ptmask(ji,:,1) 
     111             pssh(ji,:) = ( -5.5_wp + 5.5_wp*glamt(ji,1)/50._wp)*ptmask(ji,:,1) 
    115112            end do 
    116113 
     
    125122            ! 
    126123            DO ji = 1, jpi 
    127                zi = MAX(1.0-FLOAT((mig(ji)-25)**2)/400.0, 0.0 ) 
     124               zi = MAX(1.0-((glamt(ji,1)-25._wp)**2)/400.0, 0.0 ) 
    128125               DO jj = 1, jpj 
    129                   zj = MAX(1.0-FLOAT((mjg(jj)-17)**2)/144.0, 0.0 ) 
     126                  zj = MAX(1.0-((gphit(1,jj)-17._wp)**2)/144.0, 0.0 ) 
    130127                  pssh(ji,jj) = -8.5_wp + 8.5_wp*zi*zj 
    131128               END DO 
     
    141138            IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 
    142139            ! 
    143             ! Needed rn_wdmin2 increased to 0.01 for this case? 
    144140            do ji = 1,jpi 
    145              pssh(ji,:) = ( -5.5_wp + 9.0_wp*FLOAT(mig(ji))/FLOAT(jpiglo-1))*ptmask(ji,:,1) 
     141             pssh(ji,:) = ( -5.5_wp + 5.5_wp*glamt(ji,1)/50._wp)*ptmask(ji,:,1) 
    146142            end do 
    147143 
     
    166162             !pssh(ji,:) = ( -4.5_wp + 8.0_wp*FLOAT(jpiglo - mig(ji))/FLOAT(jpiglo-1))*ptmask(ji,:,1) 
    167163             !6e 
    168              pssh(ji,:) = ( -3.5_wp + 7.0_wp*FLOAT(jpiglo - mig(ji))/FLOAT(jpiglo-1))*ptmask(ji,:,1) 
     164             pssh(ji,:) = ( -3.5_wp + 5.5_wp*(50._wp-glamt(ji,1))/50._wp)*ptmask(ji,:,1) 
    169165             !6f 
    170166             !pssh(ji,:) = ( 0.5_wp + 3.75_wp*FLOAT(jpiglo - mig(ji))/FLOAT(jpiglo-1))*ptmask(ji,:,1) 
     
    173169            do ji = mi0(jpiglo/2), mi0(jpiglo) 
    174170             pts(ji,:,:,jp_sal) = 30._wp 
    175              !tsb(ji,:,:,jp_sal) = 30._wp 
    176171            end do 
    177172            ! 
     
    191186            IF( ht_wd(ji,jj) + pssh(ji,jj) < rn_wdmin1 ) THEN 
    192187               pssh(ji,jj) = ptmask(ji,jj,1)*( rn_wdmin1 - ht_wd(ji,jj) ) 
    193                !IF (jj.eq.4) write(*,*) 'SSH2 ',mig(ji),pssh(ji,4),ht_0(ji,4),-5.5_wp+5.5_wp*glamt(ji,4)/50._wp,ptmask(ji,4,1) 
    194188            ENDIF 
    195189         end do 
    196190      end do 
    197       !sshb = sshn 
    198       !ssha = sshn 
    199191      ! 
    200       !     CALL wad_istate                      ! WAD test configuration : start from pre-defined T-S fields and initial ssh slope 
    201       !    
    202192   END SUBROUTINE usr_def_istate 
    203193 
  • branches/2016/dev_merge_2016/NEMOGCM/CONFIG/WAD_TEST_CASES/MY_SRC/usrdef_zgr.F90

    r7467 r7514  
    106106                 zi = MIN((glamt(ji,1) - 10.0)/40.0, 1.0 ) 
    107107                 zht(ji,:) = MAX(zbathy*zi, 1.5)  
    108                  IF(lwp)write(numout,*) 'ZDTA ',ji,zi,zht(ji,1) 
    109108               END DO 
    110109               zht(mi0(1):mi1(1),:) = -4._wp 
     
    122121               DO ji = 1, jpi 
    123122                 zi = MAX(1.0-((glamt(ji,1)-25.0)**2)/484.0, 0.0 ) 
    124                  zi = MAX(1.0-((glamt(ji,1)-25.0)**2)/484.0, -2.0 ) 
    125                  zht(ji,:) = MAX(zbathy*zi, -20.0) 
    126                  IF(lwp)write(numout,*) 'ZDTA ',ji,zi,zht(ji,1) 
    127                END DO 
    128                zht(mi0(1):mi1(2),:) = -4._wp 
    129                zht(mi0(jpiglo-1):mi1(jpiglo),:) = -4._wp 
     123                 zht(ji,:) = MAX(zbathy*zi, 1.5) 
     124               END DO 
     125               zht(mi0(1):mi1(1),:) = -4._wp 
     126               zht(mi0(jpiglo):mi1(jpiglo),:) = -4._wp 
    130127               zht(:,mj0(1):mj1(1)) = -4._wp 
    131128               zht(:,mj0(jpjglo):mj1(jpjglo)) = -4._wp 
     
    142139               DO jj = 1, jpj 
    143140                 zj = MAX(1.0-((gphit(1,jj)-17.0)**2)/196.0, -2.0 ) 
    144                  zht(ji,jj) = MAX(zbathy*zi*zj, -2.0) 
    145                END DO 
    146                  IF(lwp)write(numout,*) 'ZDTA ',ji,zi,zht(ji,1) 
     141                 zht(ji,jj) = MAX(zbathy*zi*zj, 1.5) 
     142               END DO 
    147143               END DO 
    148144               zht(mi0(1):mi1(1),:) = -4._wp 
    149                zht(mi0(jpiglo-1):mi1(jpiglo),:) = -4._wp 
     145               zht(mi0(jpiglo):mi1(jpiglo),:) = -4._wp 
    150146               zht(:,mj0(1):mj1(1)) = -4._wp 
    151147               zht(:,mj0(jpjglo):mj1(jpjglo)) = -4._wp 
     
    160156               DO ji = 1, jpi 
    161157                 zi = MIN(glamt(ji,1)/45.0, 1.0 ) 
    162                  zht(ji,:) = MAX(zbathy*zi, 0.5) 
     158                 zht(ji,:) = MAX(zbathy*zi, 1.5) 
    163159                 IF( glamt(ji,1) >= 46.0 ) THEN 
    164160                   zht(ji,:) = 10.0 
     
    170166                 ELSE IF( glamt(ji,1) >= 4.0 .AND. glamt(ji,1) < 15.0 ) THEN 
    171167                   zi = 2.0/11.0 
    172                    zht(ji,:) = MAX(2.5 - zi*(16.0-glamt(ji,1)), 0.5) 
     168                   zht(ji,:) = MAX(2.5 - zi*(16.0-glamt(ji,1)), 1.5) 
    173169                 ELSE IF( glamt(ji,1) >= 1.0 .AND. glamt(ji,1) < 4.0 ) THEN 
    174                    zht(ji,:) = 0.5 
     170                   zht(ji,:) = 1.5 
    175171                 ENDIF 
    176                  IF(lwp)write(numout,*) 'ZDTA ',ji,zi,zht(ji,1) 
    177172               END DO 
    178173               !                                    ! =========================== 
     
    192187                 zi = MAX(1.0-((glamt(ji,1)-25.0)**2)/484.0, -2.0 ) 
    193188                 zj = 0.95*MAX(EXP(-1.0*((glamt(ji,1)-25.0)**2)/32.0) , 0.0 ) 
    194                  zht(ji,:) = MAX(zbathy*(zi-zj), -2.0) 
    195                  IF(lwp)write(numout,*) 'ZDTA ',ji,zi,zht(ji,1) 
     189                 zht(ji,:) = MAX(zbathy*(zi-zj), 1.0) 
    196190               END DO 
    197191               zht(mi0(1):mi1(1),:) = -4._wp 
     
    214208         zhu(ji,:) = 0.5_wp * ( zht(ji,:) + zht(ji+1,:) ) 
    215209      END DO 
    216       IF(lwp) write(numout,*) 'E3a zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:)) 
    217210      CALL lbc_lnk( zhu, 'U', 1. )     ! boundary condition: this mask the surrounding grid-points 
    218211      !                                ! ==>>>  set by hand non-zero value on first/last columns & rows  
    219       IF(lwp) write(numout,*) 'E3b zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:)) 
    220212      DO ji = mi0(1), mi1(1)              ! first row of global domain only 
    221213         zhu(ji,:) = zht(1,:) 
     
    224216         zhu(ji,:) = zht(jpi,:) 
    225217      END DO 
    226       IF(lwp) write(numout,*) 'E3c zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:)) 
    227218      ! at v-point: averaging zht 
    228219      zhv = 0._wp 
     
    231222      END DO 
    232223      CALL lbc_lnk( zhv, 'V', 1. )     ! boundary condition: this mask the surrounding grid-points 
    233       IF(lwp) write(numout,*) 'E3d zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:)) 
    234224      !avoid the zero depth on T- (U-,V-,F-) points 
    235225        DO jj = 1, jpj 
     
    243233          END DO 
    244234        END DO 
    245       IF(lwp) write(numout,*) 'E3e zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:)) 
    246235      !      
    247236      CALL zgr_z1d( pdept_1d, pdepw_1d, pe3t_1d , pe3w_1d )   ! Reference z-coordinate system 
     
    274263            DO ji = 1, jpi 
    275264              IF( zht(ji,jj) <= -(rn_wdld - rn_wdmin2)) THEN 
    276                 !WRITE(numout,*) 'KBOT ',ji,jj,zht(ji,jj),rn_wdld,rn_wdmin2 
    277265                k_bot(ji,jj) = 0 
    278266                zht(ji,jj) = rn_wdmin1 * REAL( jpkm1 , wp) * 0.5_wp 
     
    280268                k_bot(ji,jj) = 2 
    281269                zht(ji,jj) = rn_wdmin1 * REAL( jpkm1 , wp) * 0.5_wp 
    282                 !WRITE(numout,*) 'KBOT1 ',ji,jj,zht(ji,jj),rn_wdmin1 
    283270              ENDIF 
    284271           END DO 
     
    287274         zhu(ji,:) = 0.5_wp * ( zht(ji,:) + zht(ji+1,:) ) 
    288275      END DO 
    289       IF(lwp) write(numout,*) 'E3f zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:)) 
    290276      !                                ! ==>>>  set by hand non-zero value on 
    291277      !                                first/last columns & rows  
    292       IF(lwp) write(numout,*) 'E3f Zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:)) 
    293278      CALL lbc_lnk( zhu, 'U', 1. )     ! boundary condition: this mask the surrounding grid-points 
    294       IF(lwp) write(numout,*) 'E3F zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:)) 
    295279      DO ji = mi0(1), mi1(1)              ! first row of global domain only 
    296280         zhu(ji,:) = zht(ji,:) 
     
    300284      END DO 
    301285      ! at v-point: averaging zht 
    302       IF(lwp) write(numout,*) 'E3g zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:)) 
    303286      zhv = 0._wp 
    304287      DO jj = 1, jpjm1 
    305288         zhv(:,jj) = 0.5_wp * ( zht(:,jj) + zht(:,jj+1) ) 
    306289      END DO 
    307       IF(lwp) write(numout,*) 'E3h zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:)) 
    308290      CALL lbc_lnk( zhv, 'V', 1. ) 
    309       IF(lwp) write(numout,*) 'E3H zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:)) 
    310291      DO jj = mj0(1), mj1(1)   ! first  row of global domain only 
    311292         zhv(:,jj) = zht(:,jj) 
     
    359340         CALL lbc_lnk( pe3vw, 'V', 1. ) 
    360341         CALL lbc_lnk( pe3f , 'F', 1. ) 
    361       WHERE( pe3t (:,:,:) == 0._wp )   pe3t (:,:,:) = 1._wp 
    362       WHERE( pe3u (:,:,:) == 0._wp )   pe3u (:,:,:) = 1._wp 
    363       WHERE( pe3v (:,:,:) == 0._wp )   pe3v (:,:,:) = 1._wp 
    364       WHERE( pe3f (:,:,:) == 0._wp )   pe3f (:,:,:) = 1._wp 
    365       WHERE( pe3w (:,:,:) == 0._wp )   pe3w (:,:,:) = 1._wp 
    366       WHERE( pe3uw(:,:,:) == 0._wp )   pe3uw(:,:,:) = 1._wp 
    367       WHERE( pe3vw(:,:,:) == 0._wp )   pe3vw(:,:,:) = 1._wp 
    368       IF(lwp) write(numout,*) 'E3 TU ',SUM(0._wp/pe3t(:,:,:)), MINVAL(pe3t(:,:,:)), rn_wdmin1 ; CALL FLUSH(numout) 
    369       IF(lwp) write(numout,*) 'E3 tu ',SUM(0._wp/pe3u(:,:,:)), MINVAL(pe3u(:,:,:)) 
    370       IF(lwp) write(numout,*) 'E3 tv ',SUM(0._wp/pe3v(:,:,:)), MINVAL(pe3v(:,:,:)) 
    371       IF(lwp) write(numout,*) 'E3 tf ',SUM(0._wp/pe3f(:,:,:)), MINVAL(pe3f(:,:,:)) 
    372       IF(lwp) write(numout,*) 'E3 zu ', MINVAL(zht(:,:)), MINVAL(zhu(:,:)), MINVAL(zhv(:,:)) 
     342         WHERE( pe3t (:,:,:) == 0._wp )   pe3t (:,:,:) = 1._wp 
     343         WHERE( pe3u (:,:,:) == 0._wp )   pe3u (:,:,:) = 1._wp 
     344         WHERE( pe3v (:,:,:) == 0._wp )   pe3v (:,:,:) = 1._wp 
     345         WHERE( pe3f (:,:,:) == 0._wp )   pe3f (:,:,:) = 1._wp 
     346         WHERE( pe3w (:,:,:) == 0._wp )   pe3w (:,:,:) = 1._wp 
     347         WHERE( pe3uw(:,:,:) == 0._wp )   pe3uw(:,:,:) = 1._wp 
     348         WHERE( pe3vw(:,:,:) == 0._wp )   pe3vw(:,:,:) = 1._wp 
    373349      ENDIF 
    374350      ! 
  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90

    r7421 r7514  
    4141   USE zdfddm          ! vertical  physics: double diffusion 
    4242   USE diahth          ! thermocline diagnostics 
     43   USE wet_dry         ! wetting and drying 
    4344   ! 
    4445   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
     
    153154 
    154155      CALL iom_put( "ssh" , sshn )                 ! sea surface height 
     156      IF( iom_use("wetdep") )   &                  ! wet depth 
     157         CALL iom_put( "wetdep" , ht_wd(:,:) + sshn(:,:) ) 
    155158       
    156159      CALL iom_put( "toce", tsn(:,:,:,jp_tem) )    ! 3D temperature 
  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90

    r7490 r7514  
    3939   USE domc1d         ! 1D configuration: column location 
    4040   USE dyncor_c1d     ! 1D configuration: Coriolis term    (cor_c1d routine) 
    41    ! 
     41   USE wet_dry        ! wetting and drying 
    4242   ! 
    4343   USE in_out_manager ! I/O manager 
     
    664664      ENDIF 
    665665      ! 
     666      IF( ln_wd ) THEN              ! wetting and drying domain 
     667         CALL iom_rstput( 0, 0, inum, 'ht_0'   , ht_0   , ktype = jp_r8 ) 
     668         CALL iom_rstput( 0, 0, inum, 'ht_wd'  , ht_wd  , ktype = jp_r8 ) 
     669      ENDIF 
     670      ! 
    666671      !                                ! ============================ 
    667672      !                                !        close the files  
  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DOM/domvvl.F90

    r7430 r7514  
    2424   USE sbc_oce         ! ocean surface boundary condition 
    2525   USE wet_dry         ! wetting and drying 
     26   USE usrdef_istate   ! user defined initial state (wad only) 
    2627   USE restart         ! ocean restart 
    2728   ! 
     
    876877            ! 
    877878            IF( ln_wd .AND. ( cn_cfg == 'wad' ) ) THEN 
    878                ! 
    879                CALL wad_istate                  ! WAD test configuration : start from  
     879              ! Wetting and drying test case 
     880              CALL usr_def_istate( gdept_b, tmask, tsb, ub, vb, sshb  ) 
     881                       tsn  (:,:,:,:) = tsb (:,:,:,:)       ! set now values from to before ones 
     882                       sshn (:,:)     = sshb(:,:) 
     883                       un   (:,:,:)   = ub  (:,:,:) 
     884                       vn   (:,:,:)   = vb  (:,:,:) 
    880885                                                ! uniform T-S fields and initial ssh slope 
    881886               ! needs to be called here and in istate which is called later. 
     
    893898                DO ji = 1, jpi 
    894899                  IF( e3t_0(ji,jj,1) <= 0.5_wp * rn_wdmin1 ) THEN 
     900                     ! potential bug 
     901                     ! Warning this assumes 2 layers only over wetting locations. needs investigating 
    895902                     e3t_b(ji,jj,:) = 0.5_wp * rn_wdmin1 
    896903                     e3t_n(ji,jj,:) = 0.5_wp * rn_wdmin1 
    897904                     e3t_a(ji,jj,:) = 0.5_wp * rn_wdmin1 
    898                      sshb(ji,jj) = rn_wdmin1 - ht_0(ji,jj)           !!gm I don't understand that ! 
    899                      sshn(ji,jj) = rn_wdmin1 - ht_0(ji,jj) 
    900                      ssha(ji,jj) = rn_wdmin1 - ht_0(ji,jj) 
     905                     sshb(ji,jj) = rn_wdmin1 - ht_wd(ji,jj)           !!gm I don't understand that ! 
     906                     sshn(ji,jj) = rn_wdmin1 - ht_wd(ji,jj) 
     907                     ssha(ji,jj) = rn_wdmin1 - ht_wd(ji,jj) 
    901908                  ENDIF 
    902909                ENDDO 
  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DOM/domwri.F90

    r7421 r7514  
    1818   USE dom_oce         ! ocean space and time domain 
    1919   USE phycst ,   ONLY :   rsmall 
     20   USE wet_dry,   ONLY :   ln_wd, ht_wd 
    2021   ! 
    2122   USE in_out_manager  ! I/O manager 
     
    194195      IF( ln_sco ) THEN                                         ! s-coordinate stiffness 
    195196         CALL dom_stiff( zprt ) 
    196          CALL iom_rstput( 0, 0, inum, 'stiffness', zprt )      !    ! Max. grid stiffness ratio 
     197         CALL iom_rstput( 0, 0, inum, 'stiffness', zprt )       ! Max. grid stiffness ratio 
     198      ENDIF 
     199      ! 
     200      IF( ln_wd ) THEN                                          ! wetting and drying domain 
     201         CALL iom_rstput( 0, 0, inum, 'ht_0'   , ht_0   , ktype = jp_r8 ) 
     202         CALL iom_rstput( 0, 0, inum, 'ht_wd'  , ht_wd  , ktype = jp_r8 ) 
    197203      ENDIF 
    198204      !                                     ! ============================ 
  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DYN/dynhpg.F90

    r7437 r7514  
    432432      INTEGER  ::   ji, jj, jk, jii, jjj                 ! dummy loop indices 
    433433      REAL(wp) ::   zcoef0, zuap, zvap, znad, ztmp       ! temporary scalars 
    434       LOGICAL  ::   ll_tmp1, ll_tmp2, ll_tmp3            ! local logical variables 
     434      LOGICAL  ::   ll_tmp1, ll_tmp2                     ! local logical variables 
    435435      REAL(wp), POINTER, DIMENSION(:,:,:) ::  zhpi, zhpj 
    436436      REAL(wp), POINTER, DIMENSION(:,:)   ::  zcpx, zcpy !W/D pressure filter 
     
    438438      ! 
    439439      CALL wrk_alloc( jpi,jpj,jpk, zhpi, zhpj ) 
    440       IF(ln_wd) CALL wrk_alloc( jpi,jpj, zcpx, zcpy ) 
     440      IF( ln_wd ) CALL wrk_alloc( jpi,jpj, zcpx, zcpy ) 
    441441      ! 
    442442      IF( kt == nit000 ) THEN 
     
    451451      ENDIF 
    452452      ! 
    453       IF(ln_wd) THEN 
     453      IF( ln_wd ) THEN 
    454454        DO jj = 2, jpjm1 
    455455           DO ji = 2, jpim1  
    456              ll_tmp1 = MIN(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-ht_0(ji,jj), -ht_0(ji+1,jj))  
    457              ll_tmp2 = MAX(sshn(ji,jj) + ht_0(ji,jj), sshn(ji+1,jj) + ht_0(ji+1,jj)) > rn_wdmin1 + rn_wdmin2 
    458              ll_tmp3 = MAX(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-ht_0(ji,jj), -ht_0(ji+1,jj)) + & 
    459                                                        & rn_wdmin1 + rn_wdmin2 
    460  
    461              IF(ll_tmp1.AND.ll_tmp2) THEN 
     456             ll_tmp1 = MIN(   sshn(ji,jj)               ,   sshn(ji+1,jj) ) >                & 
     457                  &    MAX( -ht_wd(ji,jj)               , -ht_wd(ji+1,jj) ) .AND.            & 
     458                  &    MAX(   sshn(ji,jj) + ht_wd(ji,jj),   sshn(ji+1,jj) + ht_wd(ji+1,jj) ) & 
     459                  &                                                         > rn_wdmin1 + rn_wdmin2 
     460             ll_tmp2 = MAX(   sshn(ji,jj)               ,   sshn(ji+1,jj) ) >                & 
     461                  &    MAX( -ht_wd(ji,jj)               , -ht_wd(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 
     462 
     463             IF(ll_tmp1) THEN 
    462464               zcpx(ji,jj) = 1.0_wp 
    463              ELSE IF(ll_tmp3) THEN 
    464                ! no worries about sshn(ji+1,jj)-sshn(ji,jj) = 0, it won't happen ! here 
    465                zcpx(ji,jj) = ABS((sshn(ji+1,jj) + ht_0(ji+1,jj) - sshn(ji,jj) - ht_0(ji,jj)) / & 
    466                            &     (sshn(ji+1,jj) - sshn(ji,jj))) 
     465             ELSE IF(ll_tmp2) THEN 
     466               ! no worries about  sshn(ji+1,jj) -  sshn(ji  ,jj) = 0, it won't happen ! here 
     467               zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + ht_wd(ji+1,jj) - sshn(ji,jj) - ht_wd(ji,jj)) & 
     468                           &    / (sshn(ji+1,jj) -  sshn(ji  ,jj)) ) 
    467469             ELSE 
    468470               zcpx(ji,jj) = 0._wp 
    469471             END IF 
    470472       
    471              ll_tmp1 = MIN(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-ht_0(ji,jj), -ht_0(ji,jj+1))  
    472              ll_tmp2 = MAX(sshn(ji,jj) + ht_0(ji,jj), sshn(ji,jj+1) + ht_0(ji,jj+1)) > rn_wdmin1 + rn_wdmin2 
    473              ll_tmp3 = MAX(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-ht_0(ji,jj), -ht_0(ji,jj+1)) + & 
    474                                                        & rn_wdmin1 + rn_wdmin2 
    475  
    476              IF(ll_tmp1.AND.ll_tmp2) THEN 
     473             ll_tmp1 = MIN(   sshn(ji,jj)               ,   sshn(ji,jj+1) ) >                & 
     474                  &    MAX( -ht_wd(ji,jj)               , -ht_wd(ji,jj+1) ) .AND.            & 
     475                  &    MAX(   sshn(ji,jj) + ht_wd(ji,jj),   sshn(ji,jj+1) + ht_wd(ji,jj+1) ) & 
     476                  &                                                         > rn_wdmin1 + rn_wdmin2 
     477             ll_tmp2 = MAX(   sshn(ji,jj)               ,   sshn(ji,jj+1) ) >                & 
     478                  &    MAX( -ht_wd(ji,jj)               , -ht_wd(ji,jj+1) ) + rn_wdmin1 + rn_wdmin2 
     479 
     480             IF(ll_tmp1) THEN 
    477481               zcpy(ji,jj) = 1.0_wp 
    478              ELSE IF(ll_tmp3) THEN 
    479                ! no worries about sshn(ji,jj+1)-sshn(ji,jj) = 0, it won't happen ! here 
    480                zcpy(ji,jj) = ABS((sshn(ji,jj+1) + ht_0(ji,jj+1) - sshn(ji,jj) - ht_0(ji,jj)) / & 
    481                            &     (sshn(ji,jj+1) - sshn(ji,jj))) 
     482             ELSE IF(ll_tmp2) THEN 
     483               ! no worries about  sshn(ji,jj+1) -  sshn(ji,jj  ) = 0, it won't happen ! here 
     484               zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + ht_wd(ji,jj+1) - sshn(ji,jj) - ht_wd(ji,jj)) & 
     485                           &    / (sshn(ji,jj+1) -  sshn(ji,jj  )) ) 
    482486             ELSE 
    483487               zcpy(ji,jj) = 0._wp 
     
    486490        END DO 
    487491        CALL lbc_lnk( zcpx, 'U', 1._wp )    ;   CALL lbc_lnk( zcpy, 'V', 1._wp ) 
    488       ENDIF 
    489  
     492      END IF 
    490493 
    491494      ! Surface value 
     
    504507 
    505508 
    506             IF(ln_wd) THEN 
     509            IF( ln_wd ) THEN 
    507510 
    508511              zhpi(ji,jj,1) = zhpi(ji,jj,1) * zcpx(ji,jj) 
     
    535538                  &           * ( gde3w_n(ji  ,jj+1,jk) - gde3w_n(ji,jj,jk) ) * r1_e2v(ji,jj) 
    536539 
    537                IF(ln_wd) THEN 
     540               IF( ln_wd ) THEN 
    538541                 zhpi(ji,jj,jk) = zhpi(ji,jj,jk) * zcpx(ji,jj) 
    539542                 zhpj(ji,jj,jk) = zhpj(ji,jj,jk) * zcpy(ji,jj)  
     
    550553      ! 
    551554      CALL wrk_dealloc( jpi,jpj,jpk, zhpi, zhpj ) 
    552       IF(ln_wd) CALL wrk_dealloc( jpi,jpj, zcpx, zcpy ) 
     555      IF( ln_wd ) CALL wrk_dealloc( jpi,jpj, zcpx, zcpy ) 
    553556      ! 
    554557   END SUBROUTINE hpg_sco 
     
    695698      CALL wrk_alloc( jpi, jpj, jpk, drhox, drhoy, drhoz, drhou, drhov, drhow ) 
    696699      CALL wrk_alloc( jpi, jpj, jpk, rho_i, rho_j, rho_k,  zhpi,  zhpj        ) 
    697       IF(ln_wd) CALL wrk_alloc( jpi,jpj, zcpx, zcpy ) 
    698       ! 
    699       ! 
    700       IF(ln_wd) THEN 
     700      IF( ln_wd ) CALL wrk_alloc( jpi,jpj, zcpx, zcpy ) 
     701      ! 
     702      ! 
     703      IF( ln_wd ) THEN 
    701704        DO jj = 2, jpjm1 
    702705           DO ji = 2, jpim1  
    703              ll_tmp1 = MIN(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-ht_0(ji,jj), -ht_0(ji+1,jj)) & 
    704                      & .and. MAX(sshn(ji,jj) + ht_0(ji,jj), sshn(ji+1,jj) + ht_0(ji+1,jj)) & 
    705                      &  > rn_wdmin1 + rn_wdmin2 
    706              ll_tmp2 = MAX(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-ht_0(ji,jj), -ht_0(ji+1,jj)) +& 
    707                                                        & rn_wdmin1 + rn_wdmin2 
     706             ll_tmp1 = MIN(   sshn(ji,jj)               ,   sshn(ji+1,jj) ) >                & 
     707                  &    MAX( -ht_wd(ji,jj)               , -ht_wd(ji+1,jj) ) .AND.            & 
     708                  &    MAX(   sshn(ji,jj) + ht_wd(ji,jj),   sshn(ji+1,jj) + ht_wd(ji+1,jj) ) & 
     709                  &                                                         > rn_wdmin1 + rn_wdmin2 
     710             ll_tmp2 = MAX(   sshn(ji,jj)               ,   sshn(ji+1,jj) ) >                & 
     711                  &    MAX( -ht_wd(ji,jj)               , -ht_wd(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 
    708712 
    709713             IF(ll_tmp1) THEN 
    710714               zcpx(ji,jj) = 1.0_wp 
    711715             ELSE IF(ll_tmp2) THEN 
    712                ! no worries about sshn(ji+1,jj)-sshn(ji,jj) = 0, it won't happen ! here 
    713                zcpx(ji,jj) = ABS((sshn(ji+1,jj) + ht_0(ji+1,jj) - sshn(ji,jj) - ht_0(ji,jj)) /& 
    714                            &     (sshn(ji+1,jj) - sshn(ji,jj))) 
     716               ! no worries about  sshn(ji+1,jj) -  sshn(ji  ,jj) = 0, it won't happen ! here 
     717               zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + ht_wd(ji+1,jj) - sshn(ji,jj) - ht_wd(ji,jj)) & 
     718                           &    / (sshn(ji+1,jj) -  sshn(ji  ,jj)) ) 
    715719             ELSE 
    716720               zcpx(ji,jj) = 0._wp 
    717721             END IF 
    718722       
    719              ll_tmp1 = MIN(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-ht_0(ji,jj), -ht_0(ji,jj+1)) & 
    720                      & .and. MAX(sshn(ji,jj) + ht_0(ji,jj), sshn(ji,jj+1) + ht_0(ji,jj+1)) & 
    721                      &  > rn_wdmin1 + rn_wdmin2 
    722              ll_tmp2 = MAX(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-ht_0(ji,jj), -ht_0(ji,jj+1)) +& 
    723                                                        & rn_wdmin1 + rn_wdmin2 
     723             ll_tmp1 = MIN(   sshn(ji,jj)               ,   sshn(ji,jj+1) ) >                & 
     724                  &    MAX( -ht_wd(ji,jj)               , -ht_wd(ji,jj+1) ) .AND.            & 
     725                  &    MAX(   sshn(ji,jj) + ht_wd(ji,jj),   sshn(ji,jj+1) + ht_wd(ji,jj+1) ) & 
     726                  &                                                         > rn_wdmin1 + rn_wdmin2 
     727             ll_tmp2 = MAX(   sshn(ji,jj)               ,   sshn(ji,jj+1) ) >                & 
     728                  &    MAX( -ht_wd(ji,jj)               , -ht_wd(ji,jj+1) ) + rn_wdmin1 + rn_wdmin2 
    724729 
    725730             IF(ll_tmp1) THEN 
    726731               zcpy(ji,jj) = 1.0_wp 
    727732             ELSE IF(ll_tmp2) THEN 
    728                ! no worries about sshn(ji,jj+1)-sshn(ji,jj) = 0, it won't happen ! here 
    729                zcpy(ji,jj) = ABS((sshn(ji,jj+1) + ht_0(ji,jj+1) - sshn(ji,jj) - ht_0(ji,jj)) /& 
    730                            &     (sshn(ji,jj+1) - sshn(ji,jj))) 
     733               ! no worries about  sshn(ji,jj+1) -  sshn(ji,jj  ) = 0, it won't happen ! here 
     734               zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + ht_wd(ji,jj+1) - sshn(ji,jj) - ht_wd(ji,jj)) & 
     735                           &    / (sshn(ji,jj+1) -  sshn(ji,jj  )) ) 
    731736             ELSE 
    732737               zcpy(ji,jj) = 0._wp 
     
    735740        END DO 
    736741        CALL lbc_lnk( zcpx, 'U', 1._wp )    ;   CALL lbc_lnk( zcpy, 'V', 1._wp ) 
    737       ENDIF 
    738  
     742      END IF 
    739743 
    740744      IF( kt == nit000 ) THEN 
     
    907911            zhpi(ji,jj,1) = ( rho_k(ji+1,jj  ,1) - rho_k(ji,jj,1) - rho_i(ji,jj,1) ) * r1_e1u(ji,jj) 
    908912            zhpj(ji,jj,1) = ( rho_k(ji  ,jj+1,1) - rho_k(ji,jj,1) - rho_j(ji,jj,1) ) * r1_e2v(ji,jj) 
    909             IF(ln_wd) THEN 
     913            IF( ln_wd ) THEN 
    910914              zhpi(ji,jj,1) = zhpi(ji,jj,1) * zcpx(ji,jj) 
    911915              zhpj(ji,jj,1) = zhpj(ji,jj,1) * zcpy(ji,jj)  
     
    930934                  &           + (  ( rho_k(ji,jj+1,jk) - rho_k(ji,jj,jk  ) )    & 
    931935                  &               -( rho_j(ji,jj  ,jk) - rho_j(ji,jj,jk-1) )  ) * r1_e2v(ji,jj) 
    932                IF(ln_wd) THEN 
     936               IF( ln_wd ) THEN 
    933937                 zhpi(ji,jj,jk) = zhpi(ji,jj,jk) * zcpx(ji,jj) 
    934938                 zhpj(ji,jj,jk) = zhpj(ji,jj,jk) * zcpy(ji,jj)  
     
    944948      CALL wrk_dealloc( jpi, jpj, jpk, drhox, drhoy, drhoz, drhou, drhov, drhow ) 
    945949      CALL wrk_dealloc( jpi, jpj, jpk, rho_i, rho_j, rho_k,  zhpi,  zhpj        ) 
    946       IF(ln_wd) CALL wrk_dealloc( jpi,jpj, zcpx, zcpy ) 
     950      IF( ln_wd ) CALL wrk_dealloc( jpi,jpj, zcpx, zcpy ) 
    947951      ! 
    948952   END SUBROUTINE hpg_djc 
     
    981985      CALL wrk_alloc( jpi,jpj,jpk,   zdept, zrhh ) 
    982986      CALL wrk_alloc( jpi,jpj,       zsshu_n, zsshv_n ) 
    983       IF(ln_wd) CALL wrk_alloc( jpi,jpj, zcpx, zcpy ) 
     987      IF( ln_wd ) CALL wrk_alloc( jpi,jpj, zcpx, zcpy ) 
    984988      ! 
    985989      IF( kt == nit000 ) THEN 
     
    994998      IF( ln_linssh )   znad = 0._wp 
    995999 
    996       IF(ln_wd) THEN 
     1000      IF( ln_wd ) THEN 
    9971001        DO jj = 2, jpjm1 
    9981002           DO ji = 2, jpim1  
    999              ll_tmp1 = MIN(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-ht_0(ji,jj), -ht_0(ji+1,jj)) & 
    1000                      & .and. MAX(sshn(ji,jj) + ht_0(ji,jj), sshn(ji+1,jj) + ht_0(ji+1,jj)) & 
    1001                      &  > rn_wdmin1 + rn_wdmin2 
    1002              ll_tmp2 = MAX(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-ht_0(ji,jj), -ht_0(ji+1,jj)) +& 
    1003                                                        & rn_wdmin1 + rn_wdmin2 
     1003             ll_tmp1 = MIN(   sshn(ji,jj)               ,   sshn(ji+1,jj) ) >                & 
     1004                  &    MAX( -ht_wd(ji,jj)               , -ht_wd(ji+1,jj) ) .AND.            & 
     1005                  &    MAX(   sshn(ji,jj) + ht_wd(ji,jj),   sshn(ji+1,jj) + ht_wd(ji+1,jj) ) & 
     1006                  &                                                         > rn_wdmin1 + rn_wdmin2 
     1007             ll_tmp2 = MAX(   sshn(ji,jj)               ,   sshn(ji+1,jj) ) >                & 
     1008                  &    MAX( -ht_wd(ji,jj)               , -ht_wd(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 
    10041009 
    10051010             IF(ll_tmp1) THEN 
    10061011               zcpx(ji,jj) = 1.0_wp 
    10071012             ELSE IF(ll_tmp2) THEN 
    1008                ! no worries about sshn(ji+1,jj)-sshn(ji,jj) = 0, it won't happen ! here 
    1009                zcpx(ji,jj) = ABS((sshn(ji+1,jj) + ht_0(ji+1,jj) - sshn(ji,jj) - ht_0(ji,jj)) /& 
    1010                            &     (sshn(ji+1,jj) - sshn(ji,jj))) 
     1013               ! no worries about  sshn(ji+1,jj) -  sshn(ji  ,jj) = 0, it won't happen ! here 
     1014               zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + ht_wd(ji+1,jj) - sshn(ji,jj) - ht_wd(ji,jj)) & 
     1015                           &    / (sshn(ji+1,jj) -  sshn(ji  ,jj)) ) 
    10111016             ELSE 
    10121017               zcpx(ji,jj) = 0._wp 
    10131018             END IF 
    10141019       
    1015              ll_tmp1 = MIN(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-ht_0(ji,jj), -ht_0(ji,jj+1)) & 
    1016                      & .and. MAX(sshn(ji,jj) + ht_0(ji,jj), sshn(ji,jj+1) + ht_0(ji,jj+1)) & 
    1017                      &  > rn_wdmin1 + rn_wdmin2 
    1018              ll_tmp2 = MAX(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-ht_0(ji,jj), -ht_0(ji,jj+1)) +& 
    1019                                                        & rn_wdmin1 + rn_wdmin2 
    1020  
    1021              IF(ll_tmp1.OR.ll_tmp2) THEN 
     1020             ll_tmp1 = MIN(   sshn(ji,jj)               ,   sshn(ji,jj+1) ) >                & 
     1021                  &    MAX( -ht_wd(ji,jj)               , -ht_wd(ji,jj+1) ) .AND.            & 
     1022                  &    MAX(   sshn(ji,jj) + ht_wd(ji,jj),   sshn(ji,jj+1) + ht_wd(ji,jj+1) ) & 
     1023                  &                                                         > rn_wdmin1 + rn_wdmin2 
     1024             ll_tmp2 = MAX(   sshn(ji,jj)               ,   sshn(ji,jj+1) ) >                & 
     1025                  &    MAX( -ht_wd(ji,jj)               , -ht_wd(ji,jj+1) ) + rn_wdmin1 + rn_wdmin2 
     1026 
     1027             IF(ll_tmp1) THEN 
    10221028               zcpy(ji,jj) = 1.0_wp 
    10231029             ELSE IF(ll_tmp2) THEN 
    1024                ! no worries about sshn(ji,jj+1)-sshn(ji,jj) = 0, it won't happen ! here 
    1025                zcpy(ji,jj) = ABS((sshn(ji,jj+1) + ht_0(ji,jj+1) - sshn(ji,jj) - ht_0(ji,jj)) /& 
    1026                            &     (sshn(ji,jj+1) - sshn(ji,jj))) 
     1030               ! no worries about  sshn(ji,jj+1) -  sshn(ji,jj  ) = 0, it won't happen ! here 
     1031               zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + ht_wd(ji,jj+1) - sshn(ji,jj) - ht_wd(ji,jj)) & 
     1032                           &    / (sshn(ji,jj+1) -  sshn(ji,jj  )) ) 
    10271033             ELSE 
    10281034               zcpy(ji,jj) = 0._wp 
     
    10311037        END DO 
    10321038        CALL lbc_lnk( zcpx, 'U', 1._wp )    ;   CALL lbc_lnk( zcpy, 'V', 1._wp ) 
    1033       ENDIF 
     1039      END IF 
    10341040 
    10351041      ! Clean 3-D work arrays 
     
    12151221                 zdpdx2 = zcoef0 * r1_e1u(ji,jj) * REAL(jis-jid, wp) * (zpwes + zpwed) 
    12161222               ENDIF 
    1217                IF(ln_wd) THEN 
     1223               IF( ln_wd ) THEN 
    12181224                  zdpdx1 = zdpdx1 * zcpx(ji,jj) 
    12191225                  zdpdx2 = zdpdx2 * zcpx(ji,jj) 
     
    12741280                  zdpdy2 = zcoef0 * r1_e2v(ji,jj) * REAL(jjs-jjd, wp) * (zpnss + zpnsd ) 
    12751281               ENDIF 
    1276                IF(ln_wd) THEN 
     1282               IF( ln_wd ) THEN 
    12771283                  zdpdy1 = zdpdy1 * zcpy(ji,jj) 
    12781284                  zdpdy2 = zdpdy2 * zcpy(ji,jj) 
     
    12891295      CALL wrk_dealloc( jpi,jpj,jpk,   zdept, zrhh ) 
    12901296      CALL wrk_dealloc( jpi,jpj,       zsshu_n, zsshv_n ) 
    1291       IF(ln_wd) CALL wrk_dealloc( jpi,jpj, zcpx, zcpy ) 
     1297      IF( ln_wd ) CALL wrk_dealloc( jpi,jpj, zcpx, zcpy ) 
    12921298      ! 
    12931299   END SUBROUTINE hpg_prj 
  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_ts.F90

    r7482 r7514  
    261261!!gm  
    262262!!             
    263 !!            IF ( .not. ln_sco ) THEN 
    264 !! 
    265 !! !!gm  agree the JC comment  : this should be done in a much clear way 
    266 !! 
    267 !! ! JC: It not clear yet what should be the depth at f-points over land in z-coordinate case 
    268 !! !     Set it to zero for the time being  
    269 !! !              IF( rn_hmin < 0._wp ) THEN    ;   jk = - INT( rn_hmin )                                      ! from a nb of level 
    270 !! !              ELSE                          ;   jk = MINLOC( gdepw_0, mask = gdepw_0 > rn_hmin, dim = 1 )  ! from a depth 
    271 !! !              ENDIF 
    272 !! !              zhf(:,:) = gdepw_0(:,:,jk+1) 
    273 !!             ELSE 
    274 !!               zhf(:,:) = hbatf(:,:) 
    275 !!            END IF 
    276 !! 
    277 !!            DO jj = 1, jpjm1 
    278 !!               zhf(:,jj) = zhf(:,jj) * (1._wp- umask(:,jj,1) * umask(:,jj+1,1)) 
    279 !!            END DO 
     263              IF ( .not. ln_sco ) THEN 
     264   
     265   !!gm  agree the JC comment  : this should be done in a much clear way 
     266   
     267   ! JC: It not clear yet what should be the depth at f-points over land in z-coordinate case 
     268   !     Set it to zero for the time being  
     269   !              IF( rn_hmin < 0._wp ) THEN    ;   jk = - INT( rn_hmin )                                      ! from a nb of level 
     270   !              ELSE                          ;   jk = MINLOC( gdepw_0, mask = gdepw_0 > rn_hmin, dim = 1 )  ! from a depth 
     271   !              ENDIF 
     272   !              zhf(:,:) = gdepw_0(:,:,jk+1) 
     273               ELSE 
     274                 !zhf(:,:) = hbatf(:,:) 
     275                 DO jj = 1, jpjm1 
     276                   DO ji = 1, jpim1 
     277                     zhf(ji,jj) = MAX( 0._wp,                                & 
     278                                & ( ht_0(ji  ,jj  )*tmask(ji  ,jj  ,1) +     & 
     279                                &   ht_0(ji+1,jj  )*tmask(ji+1,jj  ,1) +     & 
     280                                &   ht_0(ji  ,jj+1)*tmask(ji  ,jj+1,1) +     & 
     281                                &   ht_0(ji+1,jj+1)*tmask(ji+1,jj+1,1) ) /   & 
     282                                & ( tmask(ji  ,jj  ,1) + tmask(ji+1,jj  ,1) +& 
     283                                &   tmask(ji  ,jj+1,1) + tmask(ji+1,jj+1,1) +& 
     284                                &   rsmall  ) ) 
     285                   END DO 
     286                 END DO 
     287              END IF 
     288   
     289              DO jj = 1, jpjm1 
     290                 zhf(:,jj) = zhf(:,jj) * (1._wp- umask(:,jj,1) * umask(:,jj+1,1)) 
     291              END DO 
    280292!!gm end 
    281293 
     
    381393      IF( .NOT.ln_linssh ) THEN                 ! Variable volume : remove surface pressure gradient 
    382394        IF( ln_wd ) THEN                        ! Calculating and applying W/D gravity filters 
    383           DO jj = 2, jpjm1 
    384              DO ji = 2, jpim1 
    385                 ll_tmp1 = MIN(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-ht_0(ji,jj), -ht_0(ji+1,jj))   & 
    386                         & .and. MAX(sshn(ji,jj) + ht_0(ji,jj), sshn(ji+1,jj) + ht_0(ji+1,jj))   & 
    387                         &  > rn_wdmin1 + rn_wdmin2 
    388                 ll_tmp2 = MAX(sshn(ji,jj), sshn(ji+1,jj)) > MAX(-ht_0(ji,jj), -ht_0(ji+1,jj))   & 
    389                         &                                   + rn_wdmin1 + rn_wdmin2 
     395           DO jj = 2, jpjm1 
     396              DO ji = 2, jpim1  
     397                ll_tmp1 = MIN(   sshn(ji,jj)               ,   sshn(ji+1,jj) ) >                & 
     398                     &    MAX( -ht_wd(ji,jj)               , -ht_wd(ji+1,jj) ) .AND.            & 
     399                     &    MAX(   sshn(ji,jj) + ht_wd(ji,jj),   sshn(ji+1,jj) + ht_wd(ji+1,jj) ) & 
     400                     &                                                         > rn_wdmin1 + rn_wdmin2 
     401                ll_tmp2 = MAX(   sshn(ji,jj)               ,   sshn(ji+1,jj) ) >                & 
     402                     &    MAX( -ht_wd(ji,jj)               , -ht_wd(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 
     403    
    390404                IF(ll_tmp1) THEN 
    391                   zcpx(ji,jj)    = 1.0_wp 
    392                 ELSEIF(ll_tmp2) THEN 
    393                    ! no worries about sshn(ji+1,jj)-sshn(ji,jj) = 0, it won't happen here 
    394                   zcpx(ji,jj) = ABS((sshn(ji+1,jj) + ht_0(ji+1,jj) - sshn(ji,jj) - ht_0(ji,jj)) & 
    395                         &          /(sshn(ji+1,jj) - sshn(ji,jj))) 
     405                  zcpx(ji,jj) = 1.0_wp 
     406                ELSE IF(ll_tmp2) THEN 
     407                  ! no worries about  sshn(ji+1,jj) -  sshn(ji  ,jj) = 0, it won't happen ! here 
     408                  zcpx(ji,jj) = ABS( (sshn(ji+1,jj) + ht_wd(ji+1,jj) - sshn(ji,jj) - ht_wd(ji,jj)) & 
     409                              &    / (sshn(ji+1,jj) -  sshn(ji  ,jj)) ) 
    396410                ELSE 
    397                   zcpx(ji,jj)    = 0._wp 
     411                  zcpx(ji,jj) = 0._wp 
    398412                END IF 
    399  
    400                 ll_tmp1 = MIN(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-ht_0(ji,jj), -ht_0(ji,jj+1))   & 
    401                         & .and. MAX(sshn(ji,jj) + ht_0(ji,jj), sshn(ji,jj+1) + ht_0(ji,jj+1))   & 
    402                         &  > rn_wdmin1 + rn_wdmin2 
    403                 ll_tmp2 = MAX(sshn(ji,jj), sshn(ji,jj+1)) > MAX(-ht_0(ji,jj), -ht_0(ji,jj+1))   & 
    404                         &                                   + rn_wdmin1 + rn_wdmin2 
     413          
     414                ll_tmp1 = MIN(   sshn(ji,jj)               ,   sshn(ji,jj+1) ) >                & 
     415                     &    MAX( -ht_wd(ji,jj)               , -ht_wd(ji,jj+1) ) .AND.            & 
     416                     &    MAX(   sshn(ji,jj) + ht_wd(ji,jj),   sshn(ji,jj+1) + ht_wd(ji,jj+1) ) & 
     417                     &                                                         > rn_wdmin1 + rn_wdmin2 
     418                ll_tmp2 = MAX(   sshn(ji,jj)               ,   sshn(ji,jj+1) ) >                & 
     419                     &    MAX( -ht_wd(ji,jj)               , -ht_wd(ji,jj+1) ) + rn_wdmin1 + rn_wdmin2 
     420    
    405421                IF(ll_tmp1) THEN 
    406                    zcpy(ji,jj)    = 1.0_wp 
    407                 ELSEIF(ll_tmp2) THEN 
    408                    ! no worries about sshn(ji,jj+1)-sshn(ji,jj) = 0, it won't happen here 
    409                   zcpy(ji,jj) = ABS((sshn(ji,jj+1) + ht_0(ji,jj+1) - sshn(ji,jj) - ht_0(ji,jj)) & 
    410                         &          /(sshn(ji,jj+1) - sshn(ji,jj))) 
     422                  zcpy(ji,jj) = 1.0_wp 
     423                ELSE IF(ll_tmp2) THEN 
     424                  ! no worries about  sshn(ji,jj+1) -  sshn(ji,jj  ) = 0, it won't happen ! here 
     425                  zcpy(ji,jj) = ABS( (sshn(ji,jj+1) + ht_wd(ji,jj+1) - sshn(ji,jj) - ht_wd(ji,jj)) & 
     426                              &    / (sshn(ji,jj+1) -  sshn(ji,jj  )) ) 
    411427                ELSE 
    412                   zcpy(ji,jj)    = 0._wp 
    413                 ENDIF 
    414  
    415              END DO 
     428                  zcpy(ji,jj) = 0._wp 
     429                END IF 
     430              END DO 
    416431           END DO 
    417  
    418            CALL lbc_lnk( zcpx, 'U', 1._wp )    ;   CALL lbc_lnk( zcpy, 'V', 1._wp ) 
    419  
     432  
    420433           DO jj = 2, jpjm1 
    421434              DO ji = 2, jpim1 
    422                  zu_trd(ji,jj) = ( zu_trd(ji,jj) - grav * ( sshn(ji+1,jj  ) - sshn(ji  ,jj ) )   & 
    423                         &                        * r1_e1u(ji,jj) ) * zcpx(ji,jj)  
    424                  zv_trd(ji,jj) = ( zv_trd(ji,jj) - grav * ( sshn(ji  ,jj+1) - sshn(ji  ,jj ) )   & 
    425                         &                        * r1_e2v(ji,jj) ) * zcpy(ji,jj)  
     435                 zu_trd(ji,jj) = zu_trd(ji,jj) - grav * ( sshn(ji+1,jj  ) - sshn(ji  ,jj ) )   & 
     436                        &                        * r1_e1u(ji,jj) * zcpx(ji,jj) 
     437                 zv_trd(ji,jj) = zv_trd(ji,jj) - grav * ( sshn(ji  ,jj+1) - sshn(ji  ,jj ) )   & 
     438                        &                        * r1_e2v(ji,jj) * zcpy(ji,jj) 
    426439              END DO 
    427440           END DO 
     
    570583      ENDIF 
    571584 
    572       IF( ln_wd ) THEN      !preserve the positivity of water depth 
    573                           !ssh[b,n,a] should have already been processed for this 
    574          sshbb_e(:,:) = MAX(sshbb_e(:,:), rn_wdmin1 - ht_0(:,:)) 
    575          sshb_e(:,:)  = MAX(sshb_e(:,:) , rn_wdmin1 - ht_0(:,:)) 
    576       ENDIF 
    577585      ! 
    578586      IF (ln_bt_fw) THEN                  ! FORWARD integration: start from NOW fields                     
     
    649657            zhup2_e (:,:) = hu_0(:,:) + zwx(:,:)                ! Ocean depth at U- and V-points 
    650658            zhvp2_e (:,:) = hv_0(:,:) + zwy(:,:) 
    651             IF( ln_wd ) THEN 
    652               zhup2_e(:,:) = MAX(zhup2_e (:,:), rn_wdmin1) 
    653               zhvp2_e(:,:) = MAX(zhvp2_e (:,:), rn_wdmin1) 
    654             END IF 
    655659         ELSE 
    656660            zhup2_e (:,:) = hu_n(:,:) 
     
    705709         END DO 
    706710         ssha_e(:,:) = (  sshn_e(:,:) - rdtbt * ( zssh_frc(:,:) + zhdiv(:,:) )  ) * ssmask(:,:) 
    707          IF( ln_wd ) ssha_e(:,:) = MAX(ssha_e(:,:), rn_wdmin1 - ht_0(:,:))  
     711          
    708712         CALL lbc_lnk( ssha_e, 'T',  1._wp ) 
    709713 
     
    752756         IF( ln_wd ) THEN                   ! Calculating and applying W/D gravity filters 
    753757           DO jj = 2, jpjm1 
    754               DO ji = 2, jpim1 
    755                  ll_tmp1 = MIN( zsshp2_e(ji,jj), zsshp2_e(ji+1,jj) ) > MAX( -ht_0(ji,jj), -ht_0(ji+1,jj) ) & 
    756                         & .AND. MAX( zsshp2_e(ji,jj) + ht_0(ji,jj), zsshp2_e(ji+1,jj) + ht_0(ji+1,jj) )    & 
    757                         &                                  > rn_wdmin1 + rn_wdmin2 
    758                  ll_tmp2 = MAX( zsshp2_e(ji,jj), zsshp2_e(ji+1,jj) ) > MAX( -ht_0(ji,jj), -ht_0(ji+1,jj) ) & 
    759                         &                                  + rn_wdmin1 + rn_wdmin2 
    760                  IF(ll_tmp1) THEN 
    761                     zcpx(ji,jj) = 1._wp 
    762                  ELSE IF(ll_tmp2) THEN 
    763                     ! no worries about zsshp2_e(ji+1,jj)-zsshp2_e(ji,jj) = 0, it won't happen here 
    764                     zcpx(ji,jj) = ABS( (zsshp2_e(ji+1,jj) + ht_0(ji+1,jj) - zsshp2_e(ji,jj) - ht_0(ji,jj)) & 
    765                         &             / (zsshp2_e(ji+1,jj) - zsshp2_e(ji,jj)) ) 
    766                  ELSE 
    767                     zcpx(ji,jj)    = 0._wp 
    768                  END IF 
    769  
    770                  ll_tmp1 = MIN( zsshp2_e(ji,jj), zsshp2_e(ji,jj+1) ) > MAX( -ht_0(ji,jj), -ht_0(ji,jj+1) ) & 
    771                         & .AND. MAX( zsshp2_e(ji,jj) + ht_0(ji,jj), zsshp2_e(ji,jj+1) + ht_0(ji,jj+1) )    & 
    772                         &                                  > rn_wdmin1 + rn_wdmin2 
    773                  ll_tmp2 = MAX( zsshp2_e(ji,jj), zsshp2_e(ji,jj+1) ) > MAX( -ht_0(ji,jj), -ht_0(ji,jj+1) ) & 
    774                         &                                  + rn_wdmin1 + rn_wdmin2 
    775                  IF(ll_tmp1) THEN 
    776                     zcpy(ji,jj) = 1._wp 
    777                  ELSE IF(ll_tmp2) THEN 
    778                     ! no worries about zsshp2_e(ji,jj+1)-zsshp2_e(ji,jj) = 0, it won't happen here 
    779                     zcpy(ji,jj) = ABS( (zsshp2_e(ji,jj+1) + ht_0(ji,jj+1) - zsshp2_e(ji,jj) - ht_0(ji,jj)) & 
    780                         &             / (zsshp2_e(ji,jj+1) - zsshp2_e(ji,jj)) ) 
    781                  ELSE 
    782                     zcpy(ji,jj)    = 0._wp 
    783                  END IF 
     758              DO ji = 2, jpim1  
     759                ll_tmp1 = MIN( zsshp2_e(ji,jj)               , zsshp2_e(ji+1,jj) ) >                & 
     760                     &    MAX(   -ht_wd(ji,jj)               ,   -ht_wd(ji+1,jj) ) .AND.            & 
     761                     &    MAX( zsshp2_e(ji,jj) + ht_wd(ji,jj), zsshp2_e(ji+1,jj) + ht_wd(ji+1,jj) ) & 
     762                     &                                                             > rn_wdmin1 + rn_wdmin2 
     763                ll_tmp2 = MAX( zsshp2_e(ji,jj)               , zsshp2_e(ji+1,jj) ) >                & 
     764                     &    MAX(   -ht_wd(ji,jj)               ,   -ht_wd(ji+1,jj) ) + rn_wdmin1 + rn_wdmin2 
     765    
     766                IF(ll_tmp1) THEN 
     767                  zcpx(ji,jj) = 1.0_wp 
     768                ELSE IF(ll_tmp2) THEN 
     769                  ! no worries about  zsshp2_e(ji+1,jj) - zsshp2_e(ji  ,jj) = 0, it won't happen ! here 
     770                  zcpx(ji,jj) = ABS( (zsshp2_e(ji+1,jj) +    ht_wd(ji+1,jj) - zsshp2_e(ji,jj) - ht_wd(ji,jj)) & 
     771                              &    / (zsshp2_e(ji+1,jj) - zsshp2_e(ji  ,jj)) ) 
     772                ELSE 
     773                  zcpx(ji,jj) = 0._wp 
     774                END IF 
     775          
     776                ll_tmp1 = MIN( zsshp2_e(ji,jj)               , zsshp2_e(ji,jj+1) ) >                & 
     777                     &    MAX(   -ht_wd(ji,jj)               ,   -ht_wd(ji,jj+1) ) .AND.            & 
     778                     &    MAX( zsshp2_e(ji,jj) + ht_wd(ji,jj), zsshp2_e(ji,jj+1) + ht_wd(ji,jj+1) ) & 
     779                     &                                                             > rn_wdmin1 + rn_wdmin2 
     780                ll_tmp2 = MAX( zsshp2_e(ji,jj)               , zsshp2_e(ji,jj+1) ) >                & 
     781                     &    MAX(   -ht_wd(ji,jj)               ,   -ht_wd(ji,jj+1) ) + rn_wdmin1 + rn_wdmin2 
     782    
     783                IF(ll_tmp1) THEN 
     784                  zcpy(ji,jj) = 1.0_wp 
     785                ELSE IF(ll_tmp2) THEN 
     786                  ! no worries about  zsshp2_e(ji,jj+1) - zsshp2_e(ji,jj  ) = 0, it won't happen ! here 
     787                  zcpy(ji,jj) = ABS( (zsshp2_e(ji,jj+1) +    ht_wd(ji,jj+1) - zsshp2_e(ji,jj) - ht_wd(ji,jj)) & 
     788                              &    / (zsshp2_e(ji,jj+1) - zsshp2_e(ji,jj  )) ) 
     789                ELSE 
     790                  zcpy(ji,jj) = 0._wp 
     791                END IF 
    784792              END DO 
    785             END DO 
    786             CALL lbc_lnk( zcpx, 'U', 1._wp )    ;   CALL lbc_lnk( zcpy, 'V', 1._wp ) 
    787          ENDIF 
     793           END DO 
     794         END IF 
    788795         ! 
    789796         ! Compute associated depths at U and V points: 
     
    802809               END DO 
    803810            END DO 
    804  
    805             IF( ln_wd ) THEN 
    806               zhust_e(:,:) = MAX(zhust_e (:,:), rn_wdmin1 ) 
    807               zhvst_e(:,:) = MAX(zhvst_e (:,:), rn_wdmin1 ) 
    808             END IF 
    809811 
    810812         ENDIF 
     
    12431245            zxr2 = r1_e1t(ji,jj) * r1_e1t(ji,jj) 
    12441246            zyr2 = r1_e2t(ji,jj) * r1_e2t(ji,jj) 
    1245             zcu(ji,jj) = SQRT( grav * ht_0(ji,jj) * (zxr2 + zyr2) ) 
     1247            zcu(ji,jj) = SQRT( grav * MAX(ht_0(ji,jj),0._wp) * (zxr2 + zyr2) ) 
    12461248         END DO 
    12471249      END DO 
  • branches/2016/dev_merge_2016/NEMOGCM/NEMO/OPA_SRC/DYN/wet_dry.F90

    r7436 r7514  
    3232 
    3333   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) 
    3436 
    3537   LOGICAL,  PUBLIC  ::   ln_wd       !: Wetting/drying activation switch (T:on,F:off) 
     
    8587      ! 
    8688      IF(ln_wd) THEN 
    87          ALLOCATE( wdmask(jpi,jpj), STAT=ierr ) 
     89         ALLOCATE( wdmask(jpi,jpj), ht_wd(jpi,jpj),  STAT=ierr ) 
    8890         IF( ierr /= 0 ) CALL ctl_stop('STOP', 'wad_init : Array allocation error') 
    8991      ENDIF 
     
    123125      IF(ln_wd) THEN 
    124126 
    125         CALL wrk_alloc( jpi,jpj,  zflxp, zflxn, zflxu, zflxv, zflxu1, zflxv1 ) 
    126         CALL wrk_alloc( jpi,jpj,  zwdlmtu, zwdlmtv) 
     127        CALL wrk_alloc( jpi, jpj, zflxp, zflxn, zflxu, zflxv, zflxu1, zflxv1 ) 
     128        CALL wrk_alloc( jpi, jpj, zwdlmtu, zwdlmtv) 
    127129        ! 
    128130        
     
    159161           DO ji = 2, jpi  
    160162 
    161              IF( tmask(ji,jj,1) == 0._wp  )  CYCLE   ! we don't care about land cells 
    162              IF( ht_0 (ji,jj)   >  zdepwd )   CYCLE   ! and cells which will unlikely go dried out 
     163             IF( tmask(ji, jj, 1) < 0.5_wp ) CYCLE   ! we don't care about land cells 
     164             IF( ht_wd(ji,jj) > zdepwd )      CYCLE   ! and cells which are unlikely to dry 
    163165 
    164166              zflxp(ji,jj) = max(zflxu(ji,jj), 0._wp) - min(zflxu(ji-1,jj),   0._wp) + & 
     
    167169                           & min(zflxv(ji,jj), 0._wp) - max(zflxv(ji,  jj-1), 0._wp)  
    168170 
    169               zdep2 = ht_0(ji,jj) + sshb1(ji,jj) - rn_wdmin1 
    170               IF(zdep2 < 0._wp) THEN  !add more safty, but not necessary 
     171              zdep2 = ht_wd(ji,jj) + sshb1(ji,jj) - rn_wdmin1 
     172              IF(zdep2 .le. 0._wp) THEN  !add more safty, but not necessary 
    171173                !zdep2 = 0._wp 
    172                 sshb1(ji,jj) = rn_wdmin1 - ht_0(ji,jj) 
     174                sshb1(ji,jj) = rn_wdmin1 - ht_wd(ji,jj) 
     175                wdmask(ji,jj) = 0._wp 
    173176              END IF 
    174177           ENDDO 
     
    186189              DO ji = 2, jpi  
    187190         
    188                  wdmask(ji,jj) = 0 
    189                  IF( tmask(ji,jj,1) < 0.5_wp) CYCLE  
    190                  IF( ht_0(ji,jj) > zdepwd) CYCLE 
     191                 IF( tmask(ji, jj, 1) < 0.5_wp ) CYCLE  
     192                 IF( ht_wd(ji,jj) > zdepwd )      CYCLE 
    191193         
    192194                 ztmp = e1e2t(ji,jj) 
     
    198200           
    199201                 zdep1 = (zzflxp + zzflxn) * z2dt / ztmp 
    200                  zdep2 = ht_0(ji,jj) + sshb1(ji,jj) - rn_wdmin1 - z2dt * sshemp(ji,jj)  ! this one can be moved out of the loop 
     202                 zdep2 = ht_wd(ji,jj) + sshb1(ji,jj) - rn_wdmin1 - z2dt * sshemp(ji,jj) 
    201203           
    202204                 IF(zdep1 > zdep2) THEN 
    203205                   zflag = 1 
    204206                   wdmask(ji, jj) = 0 
    205                    zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zflxp(ji,jj) * z2dt ) 
     207                   !zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zflxp(ji,jj) * z2dt ) 
     208                   zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zzflxp * z2dt ) 
    206209                   zcoef = max(zcoef, 0._wp) 
    207210                   IF(jk1 > nn_wdit) zcoef = 0._wp 
     
    307310           DO ji = 2, jpi  
    308311 
    309              IF(tmask(ji,jj,1) < 0.5_wp) CYCLE   ! we don't care about land cells 
    310              IF(ht_0 (ji,jj)   > zdepwd) CYCLE       ! and cells which will unlikely go dried out 
     312             IF( tmask(ji, jj, 1 ) < 0.5_wp) CYCLE   ! we don't care about land cells 
     313             IF( ht_wd(ji,jj) > zdepwd )      CYCLE   ! and cells which are unlikely to dry 
    311314 
    312315              zflxp(ji,jj) = max(zflxu(ji,jj), 0._wp) - min(zflxu(ji-1,jj),   0._wp) + & 
     
    315318                           & min(zflxv(ji,jj), 0._wp) - max(zflxv(ji,  jj-1), 0._wp)  
    316319 
    317               zdep2 = ht_0(ji,jj) + sshn_e(ji,jj) - rn_wdmin1 
     320              zdep2 = ht_wd(ji,jj) + sshn_e(ji,jj) - rn_wdmin1 
    318321           ENDDO 
    319322        END DO 
     
    330333              DO ji = 2, jpi  
    331334         
    332                  IF(tmask(ji,jj,1) < 0.5_wp) CYCLE  
    333                  IF(ht_0 (ji,jj)   > zdepwd) CYCLE 
     335                 IF( tmask(ji, jj, 1 ) < 0.5_wp) CYCLE  
     336                 IF( ht_wd(ji,jj) > zdepwd )      CYCLE 
    334337         
    335338                 ztmp = e1e2t(ji,jj) 
     
    341344           
    342345                 zdep1 = (zzflxp + zzflxn) * z2dt / ztmp 
    343                  zdep2 = ht_0(ji,jj) + sshn_e(ji,jj) - rn_wdmin1   ! this one can be moved out of the loop 
     346                 zdep2 = ht_wd(ji,jj) + sshn_e(ji,jj) - rn_wdmin1 
    344347                 zdep2 = zdep2 - z2dt * zssh_frc(ji,jj) 
    345348           
    346349                 IF(zdep1 > zdep2) THEN 
    347350                   zflag = 1 
    348                    !wdmask(ji, jj) = 1 
    349                    zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zflxp(ji,jj) * z2dt ) 
     351                   !zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zflxp(ji,jj) * z2dt ) 
     352                   zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zzflxp * z2dt ) 
    350353                   zcoef = max(zcoef, 0._wp) 
    351354                   IF(jk1 > nn_wdit) zcoef = 0._wp 
     
    487490            do ji = 1,jpi 
    488491             !6a 
    489              sshn(ji,:) = ( -5.5_wp + 9.0_wp*FLOAT(jpiglo - mig(ji))/FLOAT(jpiglo-1))*tmask(ji,:,1) 
     492             !sshn(ji,:) = ( -5.5_wp + 9.0_wp*FLOAT(jpiglo - mig(ji))/FLOAT(jpiglo-1))*tmask(ji,:,1) 
    490493             !Some variations in initial slope that have been tested 
    491494             !6b 
     
    495498             !6d 
    496499             !sshn(ji,:) = ( -4.5_wp + 8.0_wp*FLOAT(jpiglo - mig(ji))/FLOAT(jpiglo-1))*tmask(ji,:,1) 
     500             !6e 
     501             sshn(ji,:) = ( -3.5_wp + 7.0_wp*FLOAT(jpiglo - mig(ji))/FLOAT(jpiglo-1))*tmask(ji,:,1) 
     502             !6f 
     503             !sshn(ji,:) = ( 0.5_wp + 3.75_wp*FLOAT(jpiglo - mig(ji))/FLOAT(jpiglo-1))*tmask(ji,:,1) 
    497504            end do 
    498  
     505            ! 
     506            do ji = mi0(jpiglo/2), mi0(jpiglo) 
     507             tsn(ji,:,:,jp_sal) = 30._wp 
     508             tsb(ji,:,:,jp_sal) = 30._wp 
     509            end do 
    499510            ! 
    500511            !                                    ! =========================== 
     
    511522      do jj = 1,jpj 
    512523         do ji = 1,jpi 
    513             IF( ht_0(ji,jj) + sshn(ji,jj) < rn_wdmin1 ) THEN 
    514                sshn(ji,jj) = tmask(ji,jj,1)*( rn_wdmin1 - ht_0(ji,jj) ) 
     524            IF( ht_wd(ji,jj) + sshn(ji,jj) < rn_wdmin1 ) THEN 
     525               sshn(ji,jj) = tmask(ji,jj,1)*( rn_wdmin1 - ht_wd(ji,jj) ) 
    515526            ENDIF 
    516527         end do 
Note: See TracChangeset for help on using the changeset viewer.