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 5066 for branches/UKMO/2014_Surge_Modelling/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90 – NEMO

Ignore:
Timestamp:
2015-02-06T17:02:20+01:00 (9 years ago)
Author:
rfurner
Message:

added current state of wetting and drying code to test...note it does not work

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/UKMO/2014_Surge_Modelling/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90

    r4754 r5066  
    1717   !!            3.4  ! 2012-08  (J. Siddorn) added Siddorn and Furner stretching function 
    1818   !!            3.4  ! 2012-12  (R. Bourdalle-Badie and G. Reffray)  modify C1D case   
     19   !!            3.6  ! 2014-09  (H. Liu) Modifications for Wetting/Drying 
    1920   !!---------------------------------------------------------------------- 
    2021 
     
    11211122      REAL(wp) ::   zrmax, ztaper   ! temporary scalars 
    11221123      REAL(wp) ::   zrfact 
     1124      REAL(wp) ::   zsmth 
    11231125      ! 
    11241126      REAL(wp), POINTER, DIMENSION(:,:  ) :: ztmpi1, ztmpi2, ztmpj1, ztmpj2 
     
    11761178      bathy(:,:) = MIN( rn_sbot_max, bathy(:,:) ) 
    11771179 
    1178       DO jj = 1, jpj 
    1179          DO ji = 1, jpi 
    1180            IF( bathy(ji,jj) > 0._wp )   bathy(ji,jj) = MAX( rn_sbot_min, bathy(ji,jj) ) 
    1181          END DO 
    1182       END DO 
     1180      IF (.NOT. ln_wd) THEN 
     1181        DO jj = 1, jpj 
     1182           DO ji = 1, jpi 
     1183             IF( bathy(ji,jj) > 0._wp )   bathy(ji,jj) = MAX( rn_sbot_min, bathy(ji,jj) ) 
     1184           END DO 
     1185        END DO 
     1186      ENDIF 
    11831187      !                                        ! ============================= 
    11841188      !                                        ! Define the envelop bathymetry   (hbatt) 
     
    11871191      zenv(:,:) = bathy(:,:) 
    11881192      ! 
    1189       ! set first land point adjacent to a wet cell to sbot_min as this needs to be included in smoothing 
    1190       DO jj = 1, jpj 
    1191          DO ji = 1, jpi 
    1192            IF( bathy(ji,jj) == 0._wp ) THEN 
    1193              iip1 = MIN( ji+1, jpi ) 
    1194              ijp1 = MIN( jj+1, jpj ) 
    1195              iim1 = MAX( ji-1, 1 ) 
    1196              ijm1 = MAX( jj-1, 1 ) 
    1197              IF( (bathy(iip1,jj) + bathy(iim1,jj) + bathy(ji,ijp1) + bathy(ji,ijm1) +              & 
    1198         &         bathy(iip1,ijp1) + bathy(iim1,ijm1) + bathy(iip1,ijp1) + bathy(iim1,ijm1)) > 0._wp ) THEN 
    1199                zenv(ji,jj) = rn_sbot_min 
     1193      IF (.NOT. ln_wd) THEN 
     1194        ! set first land point adjacent to a wet cell to sbot_min as this needs to be included in smoothing 
     1195        DO jj = 1, jpj 
     1196           DO ji = 1, jpi 
     1197             IF( bathy(ji,jj) == 0._wp ) THEN 
     1198               iip1 = MIN( ji+1, jpi ) 
     1199               ijp1 = MIN( jj+1, jpj ) 
     1200               iim1 = MAX( ji-1, 1 ) 
     1201               ijm1 = MAX( jj-1, 1 ) 
     1202               IF( (bathy(iip1,jj) + bathy(iim1,jj) + bathy(ji,ijp1) + bathy(ji,ijm1) +              & 
     1203          &         bathy(iip1,ijp1) + bathy(iim1,ijm1) + bathy(iip1,ijp1) + bathy(iim1,ijm1)) > 0._wp ) THEN 
     1204                 zenv(ji,jj) = rn_sbot_min 
     1205               ENDIF 
    12001206             ENDIF 
    1201            ENDIF 
    1202          END DO 
    1203       END DO 
     1207           END DO 
     1208        END DO 
     1209      END IF 
    12041210      ! apply lateral boundary condition   CAUTION: keep the value when the lbc field is zero 
    12051211      CALL lbc_lnk( zenv, 'T', 1._wp, 'no0' ) 
     
    12411247         zri(:,:) = 0._wp 
    12421248         zrj(:,:) = 0._wp 
     1249         zsmth = 0._wp 
    12431250         DO jj = 1, nlcj 
    12441251            DO ji = 1, nlci 
    12451252               iip1 = MIN( ji+1, nlci )      ! force zri = 0 on last line (ji=ncli+1 to jpi) 
    12461253               ijp1 = MIN( jj+1, nlcj )      ! force zrj = 0 on last raw  (jj=nclj+1 to jpj) 
    1247                IF( (zenv(ji,jj) > 0._wp) .AND. (zenv(iip1,jj) > 0._wp)) THEN 
     1254               IF( (zenv(ji,jj) > zsmth) .AND. (zenv(iip1,jj) > zsmth)) THEN 
    12481255                  zri(ji,jj) = ( zenv(iip1,jj  ) - zenv(ji,jj) ) / ( zenv(iip1,jj  ) + zenv(ji,jj) ) 
    12491256               END IF 
    1250                IF( (zenv(ji,jj) > 0._wp) .AND. (zenv(ji,ijp1) > 0._wp)) THEN 
     1257               IF( (zenv(ji,jj) > zsmth) .AND. (zenv(ji,ijp1) > zsmth)) THEN 
    12511258                  zrj(ji,jj) = ( zenv(ji  ,ijp1) - zenv(ji,jj) ) / ( zenv(ji  ,ijp1) + zenv(ji,jj) ) 
    12521259               END IF 
     
    12711278      END DO                                                !     End loop     ! 
    12721279      !                                                     ! ================ ! 
    1273       DO jj = 1, jpj 
    1274          DO ji = 1, jpi 
    1275             zenv(ji,jj) = MAX( zenv(ji,jj), rn_sbot_min ) ! set all points to avoid undefined scale value warnings 
    1276          END DO 
    1277       END DO 
     1280      IF(ln_wd) THEN 
     1281        DO jj = 1, jpj 
     1282           DO ji = 1, jpi 
     1283              zenv(ji,jj) = MAX( zenv(ji,jj), -rn_wdld )    ! fill out land bathy data 
     1284           END DO 
     1285        END DO 
     1286      ELSE 
     1287        DO jj = 1, jpj 
     1288           DO ji = 1, jpi 
     1289              zenv(ji,jj) = MAX( zenv(ji,jj), rn_sbot_min ) ! set all points to avoid undefined scale value warnings 
     1290           END DO 
     1291        END DO 
     1292      END IF 
    12781293      ! 
    12791294      ! Envelope bathymetry saved in hbatt 
     
    13051320      IF(lwp) THEN 
    13061321         WRITE(numout,*) 
    1307          WRITE(numout,*) ' zgr_sco: minimum depth of the envelop topography set to : ', rn_sbot_min 
     1322         IF (.NOT. ln_wd ) THEN 
     1323           WRITE(numout,*) ' zgr_sco: minimum depth of the envelop topography set to : ', rn_sbot_min 
     1324         ELSE 
     1325           WRITE(numout,*) ' zgr_sco: minimum positive depth of the envelop topography set to : ', rn_sbot_min 
     1326           WRITE(numout,*) ' zgr_sco: minimum negative depth of the envelop topography set to : ', -rn_wdld 
     1327         ENDIF 
    13081328      ENDIF 
    13091329      hbatu(:,:) = rn_sbot_min 
     
    13181338        END DO 
    13191339      END DO 
     1340 
     1341      IF(ln_wd) THEN               !avoid the zero depth on T- (U-,V-,F-) points 
     1342        DO jj = 1, jpj 
     1343          DO ji = 1, jpi 
     1344            IF(ABS(hbatt(ji,jj)) < rn_wdmin1) & 
     1345              & hbatt(ji,jj) = SIGN(1._wp, hbatt(ji,jj)) * rn_wdmin1 
     1346            IF(ABS(hbatu(ji,jj)) < rn_wdmin1) & 
     1347              & hbatu(ji,jj) = SIGN(1._wp, hbatu(ji,jj)) * rn_wdmin1 
     1348            IF(ABS(hbatv(ji,jj)) < rn_wdmin1) & 
     1349              & hbatv(ji,jj) = SIGN(1._wp, hbatv(ji,jj)) * rn_wdmin1 
     1350            IF(ABS(hbatf(ji,jj)) < rn_wdmin1) & 
     1351              & hbatf(ji,jj) = SIGN(1._wp, hbatf(ji,jj)) * rn_wdmin1 
     1352          END DO 
     1353        END DO 
     1354      END IF 
    13201355      !  
    13211356      ! Apply lateral boundary condition 
     
    13251360         DO ji = 1, jpi 
    13261361            IF( hbatu(ji,jj) == 0._wp ) THEN 
     1362               !No worries about the following line when ln_wd == .true. 
    13271363               IF( zhbat(ji,jj) == 0._wp )   hbatu(ji,jj) = rn_sbot_min 
    13281364               IF( zhbat(ji,jj) /= 0._wp )   hbatu(ji,jj) = zhbat(ji,jj) 
     
    13501386 
    13511387!!bug:  key_helsinki a verifer 
    1352       hift(:,:) = MIN( hift(:,:), hbatt(:,:) ) 
    1353       hifu(:,:) = MIN( hifu(:,:), hbatu(:,:) ) 
    1354       hifv(:,:) = MIN( hifv(:,:), hbatv(:,:) ) 
    1355       hiff(:,:) = MIN( hiff(:,:), hbatf(:,:) ) 
     1388 
     1389      IF (.NOT. ln_wd) THEN 
     1390        hift(:,:) = MIN( hift(:,:), hbatt(:,:) ) 
     1391        hifu(:,:) = MIN( hifu(:,:), hbatu(:,:) ) 
     1392        hifv(:,:) = MIN( hifv(:,:), hbatv(:,:) ) 
     1393        hiff(:,:) = MIN( hiff(:,:), hbatf(:,:) ) 
     1394      END IF 
    13561395 
    13571396      IF( nprint == 1 .AND. lwp )   THEN 
     
    13951434      CALL lbc_lnk( e3vw_0, 'V', 1._wp ) 
    13961435 
    1397       fsdepw(:,:,:) = gdepw_0 (:,:,:) 
    1398       fsde3w(:,:,:) = gdep3w_0(:,:,:) 
    1399       ! 
    1400       where (e3t_0   (:,:,:).eq.0.0)  e3t_0(:,:,:) = 1.0 
    1401       where (e3u_0   (:,:,:).eq.0.0)  e3u_0(:,:,:) = 1.0 
    1402       where (e3v_0   (:,:,:).eq.0.0)  e3v_0(:,:,:) = 1.0 
    1403       where (e3f_0   (:,:,:).eq.0.0)  e3f_0(:,:,:) = 1.0 
    1404       where (e3w_0   (:,:,:).eq.0.0)  e3w_0(:,:,:) = 1.0 
    1405       where (e3uw_0  (:,:,:).eq.0.0)  e3uw_0(:,:,:) = 1.0 
    1406       where (e3vw_0  (:,:,:).eq.0.0)  e3vw_0(:,:,:) = 1.0 
     1436      !fsdepw(:,:,:) = gdepw_0 (:,:,:) 
     1437      !fsde3w(:,:,:) = gdep3w_0(:,:,:) 
     1438      ! 
     1439      IF (.NOT. ln_wd) THEN 
     1440        where (e3t_0   (:,:,:).eq.0.0)  e3t_0(:,:,:) = 1.0 
     1441        where (e3u_0   (:,:,:).eq.0.0)  e3u_0(:,:,:) = 1.0 
     1442        where (e3v_0   (:,:,:).eq.0.0)  e3v_0(:,:,:) = 1.0 
     1443        where (e3f_0   (:,:,:).eq.0.0)  e3f_0(:,:,:) = 1.0 
     1444        where (e3w_0   (:,:,:).eq.0.0)  e3w_0(:,:,:) = 1.0 
     1445        where (e3uw_0  (:,:,:).eq.0.0)  e3uw_0(:,:,:) = 1.0 
     1446        where (e3vw_0  (:,:,:).eq.0.0)  e3vw_0(:,:,:) = 1.0 
     1447      END IF 
    14071448 
    14081449#if defined key_agrif 
     
    14491490               IF( scobot(ji,jj) >= fsdept(ji,jj,jk) )   mbathy(ji,jj) = MAX( 2, jk ) 
    14501491#endif 
    1451                IF( scobot(ji,jj) == 0._wp            )   mbathy(ji,jj) = 0 
    1452             END DO 
     1492            END DO 
     1493            IF(ln_wd) THEN 
     1494              IF( scobot(ji,jj) <= -(rn_wdld - rn_wdmin2)) THEN 
     1495                mbathy(ji,jj) = 0 
     1496              ELSEIF(scobot(ji,jj) <= rn_wdmin1) THEN 
     1497#if key_surge 
     1498                mbathy(ji,jj) = 1 
     1499#else 
     1500                mbathy(ji,jj) = 2 
     1501#endif 
     1502              ENDIF 
     1503            ELSE 
     1504              IF( scobot(ji,jj) == 0._wp            )   mbathy(ji,jj) = 0 
     1505            END IF 
    14531506         END DO 
    14541507      END DO 
     
    15081561         DO jj = 1, jpj 
    15091562 
    1510             IF( hbatt(ji,jj) > 0._wp) THEN 
     1563            IF( bathy(ji,jj) > 0._wp) THEN 
    15111564               DO jk = 1, mbathy(ji,jj) 
    15121565                 ! check coordinate is monotonically increasing 
    1513                  IF (fse3w(ji,jj,jk) <= 0._wp .OR. fse3t(ji,jj,jk) <= 0._wp ) THEN 
    1514                     WRITE(ctmp1,*) 'ERROR zgr_sco :   e3w   or e3t   =< 0  at point (i,j,k)= ', ji, jj, jk 
    1515                     WRITE(numout,*) 'ERROR zgr_sco :   e3w   or e3t   =< 0  at point (i,j,k)= ', ji, jj, jk 
     1566                 !RF test... 
     1567                 !IF (fse3w(ji,jj,jk) <= 0._wp .OR. fse3t(ji,jj,jk) <= 0._wp ) THEN 
     1568                 !   WRITE(ctmp1,*) 'ERROR zgr_sco :   e3w   or e3t   =< 0  at point (i,j,k)= ', ji, jj, jk 
     1569                 !   WRITE(numout,*) 'ERROR zgr_sco :   e3w   or e3t   =< 0  at point (i,j,k)= ', ji, jj, jk 
     1570                 IF (fse3w(ji,jj,jk) < 0._wp .OR. fse3t(ji,jj,jk) < 0._wp ) THEN 
     1571                    WRITE(ctmp1,*) 'ERROR zgr_sco :   e3w   or e3t   < 0  at point (i,j,k)= ', ji, jj, jk 
     1572                    WRITE(numout,*) 'ERROR zgr_sco :   e3w   or e3t   < 0  at point (i,j,k)= ', ji, jj, jk 
    15161573                    WRITE(numout,*) 'e3w',fse3w(ji,jj,:) 
    15171574                    WRITE(numout,*) 'e3t',fse3t(ji,jj,:) 
     
    15721629      INTEGER  ::   ji, jj, jk           ! dummy loop argument 
    15731630      REAL(wp) ::   zcoeft, zcoefw   ! temporary scalars 
     1631      REAL(wp) ::   ztmpu,  ztmpv,  ztmpf 
     1632      REAL(wp) ::   ztmpu1, ztmpv1, ztmpf1 
    15741633      ! 
    15751634      REAL(wp), POINTER, DIMENSION(:,:,:) :: z_gsigw3, z_gsigt3, z_gsi3w3 
     
    16251684      DO ji = 1, jpim1 
    16261685         DO jj = 1, jpjm1 
     1686            ! extended for Wetting/Drying case 
     1687            ztmpu = hbatt(ji,jj)+hbatt(ji+1,jj) 
     1688            ztmpv = hbatt(ji,jj)+hbatt(ji,jj+1) 
     1689            ztmpf = hbatt(ji,jj)+hbatt(ji+1,jj)+hbatt(ji,jj+1)+hbatt(ji+1,jj+1) 
     1690            ztmpu1 = hbatt(ji,jj)*hbatt(ji+1,jj) 
     1691            ztmpv1 = hbatt(ji,jj)*hbatt(ji,jj+1) 
     1692            ztmpf1 = MIN(hbatt(ji,jj), hbatt(ji+1,jj), hbatt(ji,jj+1), hbatt(ji+1,jj+1)) * & 
     1693                   & MAX(hbatt(ji,jj), hbatt(ji+1,jj), hbatt(ji,jj+1), hbatt(ji+1,jj+1)) 
    16271694            DO jk = 1, jpk 
    1628                z_esigtu3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigt3(ji+1,jj,jk) )   & 
    1629                   &              / ( hbatt(ji,jj)+hbatt(ji+1,jj) ) 
    1630                z_esigtv3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji,jj+1)*z_esigt3(ji,jj+1,jk) )   & 
    1631                   &              / ( hbatt(ji,jj)+hbatt(ji,jj+1) ) 
    1632                z_esigtf3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigt3(ji+1,jj,jk)     & 
    1633                   &                + hbatt(ji,jj+1)*z_esigt3(ji,jj+1,jk)+hbatt(ji+1,jj+1)*z_esigt3(ji+1,jj+1,jk) )   & 
    1634                   &              / ( hbatt(ji,jj)+hbatt(ji+1,jj)+hbatt(ji,jj+1)+hbatt(ji+1,jj+1) ) 
    1635                z_esigwu3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigw3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigw3(ji+1,jj,jk) )   & 
    1636                   &              / ( hbatt(ji,jj)+hbatt(ji+1,jj) ) 
    1637                z_esigwv3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigw3(ji,jj,jk)+hbatt(ji,jj+1)*z_esigw3(ji,jj+1,jk) )   & 
    1638                   &              / ( hbatt(ji,jj)+hbatt(ji,jj+1) ) 
     1695               IF((ztmpu1 < 0._wp.OR.ABS(ztmpu) < rn_wdmin1).AND.ln_wd) THEN 
     1696                 z_esigtu3(ji,jj,jk) = 0.5_wp * ( z_esigt3(ji,jj,jk) + z_esigt3(ji+1,jj,jk) ) 
     1697                 z_esigwu3(ji,jj,jk) = 0.5_wp * ( z_esigw3(ji,jj,jk) + z_esigw3(ji+1,jj,jk) ) 
     1698               ELSE 
     1699                 z_esigtu3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigt3(ji+1,jj,jk) ) & 
     1700                        &              / ztmpu 
     1701                 z_esigwu3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigw3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigw3(ji+1,jj,jk) ) & 
     1702                        &              / ztmpu 
     1703               END IF 
     1704 
     1705               IF((ztmpv1 < 0._wp.OR.ABS(ztmpv) < rn_wdmin1).AND.ln_wd) THEN 
     1706                 z_esigtv3(ji,jj,jk) = 0.5_wp * ( z_esigt3(ji,jj,jk) + z_esigt3(ji,jj+1,jk) ) 
     1707                 z_esigwv3(ji,jj,jk) = 0.5_wp * ( z_esigw3(ji,jj,jk) + z_esigw3(ji,jj+1,jk) ) 
     1708               ELSE 
     1709                 z_esigtv3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji,jj+1)*z_esigt3(ji,jj+1,jk) ) & 
     1710                        &              / ztmpv 
     1711                 z_esigwv3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigw3(ji,jj,jk)+hbatt(ji,jj+1)*z_esigw3(ji,jj+1,jk) ) & 
     1712                        &              / ztmpv 
     1713               END IF 
     1714 
     1715               IF((ztmpf1 < 0._wp.OR.ABS(ztmpf) < rn_wdmin1).AND.ln_wd) THEN 
     1716                 z_esigtf3(ji,jj,jk) = 0.25_wp * ( z_esigt3(ji,jj,jk)   + z_esigt3(ji+1,jj,jk)  & 
     1717                    &                            + z_esigt3(ji,jj+1,jk) + z_esigt3(ji+1,jj+1,jk) ) 
     1718               ELSE 
     1719                 z_esigtf3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigt3(ji+1,jj,jk)  & 
     1720                    &                + hbatt(ji,jj+1)*z_esigt3(ji,jj+1,jk)+hbatt(ji+1,jj+1)*z_esigt3(ji+1,jj+1,jk) ) & 
     1721                    &              / ztmpf 
     1722               END IF 
     1723 
    16391724               ! 
    16401725               e3t_0(ji,jj,jk) = ( (hbatt(ji,jj)-rn_hc)*z_esigt3 (ji,jj,jk) + rn_hc/REAL(jpkm1,wp) ) 
     
    16761761      REAL(wp) ::   zsmth               ! smoothing around critical depth 
    16771762      REAL(wp) ::   zzs, zzb           ! Surface and bottom cell thickness in sigma space 
     1763      REAL(wp) ::   ztmpu, ztmpv, ztmpf 
     1764      REAL(wp) ::   ztmpu1, ztmpv1, ztmpf1 
    16781765      ! 
    16791766      REAL(wp), POINTER, DIMENSION(:,:,:) :: z_gsigw3, z_gsigt3, z_gsi3w3 
     
    17541841        DO jj=1,jpj-1 
    17551842 
    1756           DO jk = 1, jpk 
    1757                 z_esigtu3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigt3(ji+1,jj,jk) ) / & 
    1758                                     ( hbatt(ji,jj)+hbatt(ji+1,jj) ) 
    1759                 z_esigtv3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji,jj+1)*z_esigt3(ji,jj+1,jk) ) / & 
    1760                                     ( hbatt(ji,jj)+hbatt(ji,jj+1) ) 
    1761                 z_esigtf3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigt3(ji+1,jj,jk) +  & 
    1762                                       hbatt(ji,jj+1)*z_esigt3(ji,jj+1,jk)+hbatt(ji+1,jj+1)*z_esigt3(ji+1,jj+1,jk) ) / & 
    1763                                     ( hbatt(ji,jj)+hbatt(ji+1,jj)+hbatt(ji,jj+1)+hbatt(ji+1,jj+1) ) 
    1764                 z_esigwu3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigw3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigw3(ji+1,jj,jk) ) / & 
    1765                                     ( hbatt(ji,jj)+hbatt(ji+1,jj) ) 
    1766                 z_esigwv3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigw3(ji,jj,jk)+hbatt(ji,jj+1)*z_esigw3(ji,jj+1,jk) ) / & 
    1767                                     ( hbatt(ji,jj)+hbatt(ji,jj+1) ) 
     1843           ! extend to suit for Wetting/Drying case 
     1844           ztmpu = hbatt(ji,jj)+hbatt(ji+1,jj) 
     1845           ztmpv = hbatt(ji,jj)+hbatt(ji,jj+1) 
     1846           ztmpf = hbatt(ji,jj)+hbatt(ji+1,jj)+hbatt(ji,jj+1)+hbatt(ji+1,jj+1) 
     1847           ztmpu1 = hbatt(ji,jj)*hbatt(ji+1,jj) 
     1848           ztmpv1 = hbatt(ji,jj)*hbatt(ji,jj+1) 
     1849           ztmpf1 = MIN(hbatt(ji,jj), hbatt(ji+1,jj), hbatt(ji,jj+1), hbatt(ji+1,jj+1)) * & 
     1850                  & MAX(hbatt(ji,jj), hbatt(ji+1,jj), hbatt(ji,jj+1), hbatt(ji+1,jj+1)) 
     1851           DO jk = 1, jpk 
     1852              IF((ztmpu1 < 0._wp.OR.ABS(ztmpu) < rn_wdmin1).AND.ln_wd) THEN 
     1853                z_esigtu3(ji,jj,jk) = 0.5_wp * ( z_esigt3(ji,jj,jk) + z_esigt3(ji+1,jj,jk) ) 
     1854                z_esigwu3(ji,jj,jk) = 0.5_wp * ( z_esigw3(ji,jj,jk) + z_esigw3(ji+1,jj,jk) ) 
     1855              ELSE 
     1856                z_esigtu3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigt3(ji+1,jj,jk) ) & 
     1857                       &              / ztmpu 
     1858                z_esigwu3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigw3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigw3(ji+1,jj,jk) ) & 
     1859                       &              / ztmpu 
     1860              END IF 
     1861 
     1862              IF((ztmpv1 < 0._wp.OR.ABS(ztmpv) < rn_wdmin1).AND.ln_wd) THEN 
     1863                z_esigtv3(ji,jj,jk) = 0.5_wp * ( z_esigt3(ji,jj,jk) + z_esigt3(ji,jj+1,jk) ) 
     1864                z_esigwv3(ji,jj,jk) = 0.5_wp * ( z_esigw3(ji,jj,jk) + z_esigw3(ji,jj+1,jk) ) 
     1865              ELSE 
     1866                z_esigtv3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji,jj+1)*z_esigt3(ji,jj+1,jk) ) & 
     1867                       &              / ztmpv 
     1868                z_esigwv3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigw3(ji,jj,jk)+hbatt(ji,jj+1)*z_esigw3(ji,jj+1,jk) ) & 
     1869                       &              / ztmpv 
     1870              END IF 
     1871 
     1872              IF((ztmpf1 < 0._wp.OR.ABS(ztmpf) < rn_wdmin1).AND.ln_wd) THEN 
     1873                z_esigtf3(ji,jj,jk) = 0.25_wp * ( z_esigt3(ji,jj,jk)   + z_esigt3(ji+1,jj,jk)  & 
     1874                   &                            + z_esigt3(ji,jj+1,jk) + z_esigt3(ji+1,jj+1,jk) ) 
     1875              ELSE 
     1876                z_esigtf3(ji,jj,jk) = ( hbatt(ji,jj)*z_esigt3(ji,jj,jk)+hbatt(ji+1,jj)*z_esigt3(ji+1,jj,jk)  & 
     1877                   &                + hbatt(ji,jj+1)*z_esigt3(ji,jj+1,jk)+hbatt(ji+1,jj+1)*z_esigt3(ji+1,jj+1,jk) ) & 
     1878                   &              / ztmpf 
     1879              END IF 
    17681880 
    17691881             e3t_0(ji,jj,jk)=(scosrf(ji,jj)+hbatt(ji,jj))*z_esigt3(ji,jj,jk) 
Note: See TracChangeset for help on using the changeset viewer.