Index: /branches/2013/dev_r4050_NOC_WaD/NEMOGCM/NEMO/OPA_SRC/DOM/dom_oce.F90
===================================================================
--- /branches/2013/dev_r4050_NOC_WaD/NEMOGCM/NEMO/OPA_SRC/DOM/dom_oce.F90	(revision 4354)
+++ /branches/2013/dev_r4050_NOC_WaD/NEMOGCM/NEMO/OPA_SRC/DOM/dom_oce.F90	(revision 4355)
@@ -54,4 +54,10 @@
    REAL(wp), PUBLIC ::   rdtmax          !: maximum time step on tracers
    REAL(wp), PUBLIC ::   rdth            !: depth variation of tracer step
+
+
+   LOGICAL,  PUBLIC ::   ln_wad     =  .FALSE.  !: turn on wetting/drying (T: on, F: off)
+   REAL(wp), PUBLIC ::   rn_wadmin  = 0.10_wp   !: minimum water depth on dried cells
+   REAL(wp), PUBLIC ::   rn_wadmine = 0.01_wp   !: tolerrance of minimum water depth on dried cells
+   REAL(wp), PUBLIC ::   rn_landele = 20.0_wp   !: land elevation below which wetting/drying will be considered
 
    !                                                  !!! associated variables
Index: /branches/2013/dev_r4050_NOC_WaD/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90
===================================================================
--- /branches/2013/dev_r4050_NOC_WaD/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90	(revision 4354)
+++ /branches/2013/dev_r4050_NOC_WaD/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90	(revision 4355)
@@ -81,4 +81,5 @@
       !
                              CALL dom_nam      ! read namelist ( namrun, namdom, namcla )
+                             CALL dom_wad      ! read namelist ( namwad ), wetting/drying
                              CALL dom_clo      ! Closed seas and lake
                              CALL dom_hgr      ! Horizontal mesh
@@ -269,4 +270,32 @@
       !
    END SUBROUTINE dom_nam
+
+   SUBROUTINE dom_wad
+      !!----------------------------------------------------------------------
+      !!                     ***  ROUTINE dom_wad  ***
+      !!                    
+      !! ** Purpose :   read wetting/drying namelists and print the variables.
+      !!
+      !! ** input   : - namwad namelist
+      !!----------------------------------------------------------------------
+      USE ioipsl
+      NAMELIST/namwad/ ln_wad , rn_wadmin, rn_wadmine, rn_landele
+      !!----------------------------------------------------------------------
+
+      REWIND( numnam )              ! Namelist namrun : parameters of the run
+      READ  ( numnam, namwad )
+      !
+      IF(lwp) THEN                  ! control print
+         WRITE(numout,*)
+         WRITE(numout,*) 'dom_wad  : set up wetting/drying '
+         WRITE(numout,*) '~~~~~~~ '
+         WRITE(numout,*) '   Namelist namwad'
+         WRITE(numout,*) '      turn on(off) wetting/drying     ln_wad        = ', ln_wad
+         WRITE(numout,*) '      minimal dry cell depth          rn_ wadmin    = ', rn_wadmin 
+         WRITE(numout,*) '      tollerance of minimual depth    rn_wadmine    = ', rn_wadmine
+         WRITE(numout,*) ' threshold elevation for land masking rn_landele    = ', rn_landele
+      ENDIF
+
+   END SUBROUTINE dom_wad
 
 
Index: /branches/2013/dev_r4050_NOC_WaD/NEMOGCM/NEMO/OPA_SRC/DYN/wadlmt.F90
===================================================================
--- /branches/2013/dev_r4050_NOC_WaD/NEMOGCM/NEMO/OPA_SRC/DYN/wadlmt.F90	(revision 4355)
+++ /branches/2013/dev_r4050_NOC_WaD/NEMOGCM/NEMO/OPA_SRC/DYN/wadlmt.F90	(revision 4355)
@@ -0,0 +1,165 @@
+
+MODULE wadlmt
+   !!==============================================================================
+   !!                       ***  MODULE  wadlmt  ***
+   !! compute the flux limiter for water depth positivity in !wetting/drying case
+   !!==============================================================================
+   !! History :   
+   !!  NEMO      3.5  ! 2014-01  ((H.Liu)  Original code
+   !!----------------------------------------------------------------------
+
+   !!----------------------------------------------------------------------
+   !!   wad_lmt    : Compute the horizontal flux limiter and the limited velocity
+   !!                when wetting and drying happens 
+   !!----------------------------------------------------------------------
+   USE oce             ! ocean dynamics and tracers
+   USE dom_oce         ! ocean space and time domain
+   USE sbc_oce, ONLY : ln_rnf   ! surface boundary condition: ocean
+   USE sbcrnf          ! river runoff 
+   USE cla             ! cross land advection             (cla_div routine)
+   USE in_out_manager  ! I/O manager
+   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
+   USE lib_mpp         ! MPP library
+   USE wrk_nemo        ! Memory Allocation
+   USE timing          ! Timing
+
+   IMPLICIT NONE
+   PRIVATE
+
+
+   PUBLIC   wad_lmt    ! routine called by step.F90
+
+   SUBROUTINE wad_lmt( kt )
+      !!----------------------------------------------------------------------
+      !!                  ***  ROUTINE wad_lmt  ***
+      !!                    
+      !! ** Purpose :   generate flux limiters for wetting/drying
+      !!
+      !! ** Method  : - Prevent negative depth occurring (Not ready for Agrif) 
+      !!
+      !! ** Action  : - update: uwdlmt(:,:), vwdlmt(:,:) 
+      !!----------------------------------------------------------------------
+      INTEGER, INTENT(in) ::   kt       ! ocean time-step index
+      !
+      INTEGER  ::   ji, jj, jk, jk1     ! dummy loop indices
+      INTEGER  ::   zflag, z2dt         ! local scalar
+      REAL(wp) ::   zflxp, zflxn        ! local scalars
+      REAL(wp) ::   zcoef, zdep1, zdep2 ! local scalars
+      REAL(wp) ::   ztmp                ! local scalars
+      REAL(wp), POINTER,  DIMENSION(:,:) ::   zflxu,  zflxv   ! specific 2D workspace
+      REAL(wp), POINTER,  DIMENSION(:,:) ::   zflxu1, zflxv1  ! specific 2D workspace
+      REAL(wp), POINTER,  DIMENSION(:,:) ::   uwdlmt, vwdlmt  !: W/D flux limiters
+      !!----------------------------------------------------------------------
+      !
+
+      IF( nn_timing == 1 )  CALL timing_start('wad_lmt')
+
+      IF(ln_wad) THEN
+
+        CALL wrk_alloc( jpi, jpj, zflxu, zflxv, zflxu1, zflxv1, uwdlmt, vwdlmt)
+        !
+       
+        IF(lwp) WRITE(numout,*)
+        IF(lwp) WRITE(numout,*) 'wad_lmt : wetting/drying limiters and velocity limiting'
+       
+        z2dt = 2. * rdt                                 ! Euler or leap-frog time step 
+        IF( neuler == 0 .AND. kt == nit000 )  z2dt = rdt
+       
+       
+        zflxu(:,:)  = 0._wp
+        zflxv(:,:)  = 0._wp
+        uwdlmt(:,:) = 1._wp
+        vwdlmt(:,:) = 1._wp
+       
+       
+        ! Horizontal Flux in u and v direction
+        DO jk = 1, jpkm1  
+           DO jj = 1, jpjm1
+              DO ji = 1, jpim1
+                 zflxu(ji,jj) = zflxu(ji,jj) + fse3u(ji,jj,jk) * un(ji,jj,jk) * umask(ji,jj,jk)
+                 zflxv(ji,jj) = zflxv(ji,jj) + fse3v(ji,jj,jk) * vn(ji,jj,jk) * vmask(ji,jj,jk)
+              END DO  
+           END DO  
+        END DO
+       
+        zflxu(:,:) = zflxu(:,:) * e2u(:,:)
+        zflxv(:,:) = zflxv(:,:) * e1v(:,:)
+       
+        DO jk1 = 1, 10
+       
+        zflag = 0     ! flag marking that if we need further iteration of W/D limiters?
+       
+        zflxu1(:,:) = zflxu(:,:) * uwdlmt(:,:)
+        zflxv1(:,:) = zflxv(:,:) * vwdlmt(:,:)
+       
+       
+        DO jj = 2, jpjm1
+           DO ji = fs_2, fs_jpim1   ! vector opt.
+              zflxp = max(zflxu1(ji,jj), 0._wp) - min(zflxu1(ji-1,jj), 0._wp) + &
+                    & max(zflxv1(ji,jj), 0._wp) - min(zflxv1(ji,jj-1), 0._wp) 
+              zflxn = min(zflxu1(ji,jj), 0._wp) - max(zflxu1(ji-1,jj), 0._wp) + &
+                    & min(zflxv1(ji,jj), 0._wp) - max(zflxv1(ji,jj-1), 0._wp) 
+       
+              IF(tmask(ji, jj, 1) < 0.5_wp) CYCLE 
+       
+              !IF( .NOT. AGRIF_Root() ) THEN    !Wetting and Drying but not for AGRIF
+                 IF (((nbondi ==  1).OR.(nbondi == 2)).AND.(ji == nlci-1)) CYCLE ! east
+                 IF (((nbondi == -1).OR.(nbondi == 2)).AND.(ji == 2     )) CYCLE ! west 
+                 IF (((nbondj ==  1).OR.(nbondj == 2)).AND.(jj == nlcj-1)) CYCLE ! north
+                 IF (((nbondj == -1).OR.(nbondj == 2)).AND.(jj == 2     )) CYCLE ! south
+              !END IF
+       
+              
+              ztmp = e1t(ji,jj) * e2t(ji,jj)
+       
+              zdep1 = (zflxp + zflxn) * z2dt / ztmp
+              zdep2 = bathy(ji,jj) + sshb(ji,jj) - rn_wadmin 
+       
+              IF(zdep2 < 0._wp) THEN
+                !WRITE(numout,*) 'depth less than minimum depth, cell(ji,jj):', ji,jj
+                zdep2 = 0._wp
+                sshb(ji,jj) = rn_wad_dep1 - bathy(ji,jj)
+              END IF
+       
+              IF(zdep1 > zdep2) THEN
+                zflag = 1
+                zcoef = ( ( zdep2 + rn_wadmine ) * ztmp + zflxn * z2dt ) / ( zflxp * z2dt )
+                IF(zflxu1(ji,  jj) >= 0._wp) uwdlmt(ji,jj  ) = zcoef
+                IF(zflxu1(ji-1,jj) <  0._wp) uwdlmt(ji-1,jj) = zcoef
+                IF(zflxv1(ji,  jj) >= 0._wp) vwdlmt(ji,jj  ) = zcoef
+                IF(zflxv1(ji,jj-1) <  0._wp) vwdlmt(ji,jj-1) = zcoef
+              END IF
+           END DO ! ji loop
+        END DO  ! jj loop
+       
+        CALL lbc_lnk( uwdlmt, 'U', 1. )   ;   CALL lbc_lnk( vwdlmt , 'V', 1. ) 
+       
+        IF(zflag == 0) EXIT
+       
+        END DO  ! jk1 loop
+       
+        DO jk = 1, jpkm1  
+           DO jj = 2, jpjm1
+              DO ji = fs_2, fs_jpim1   ! vector opt.
+               un(ji,jj,jk)   = 0.5_wp * ( sign(1._wp, un(ji,jj,jk)   + 1._wp ) * uwdlmt(ji,jj) 
+               vn(ji,jj,jk)   = 0.5_wp * ( sign(1._wp, vn(ji,jj,jk)   + 1._wp ) * vwdlmt(ji,jj) 
+               un(ji-1,jj,jk) = 0.5_wp * ( 1._wp - sign(1._wp, un(ji-1,jj,jk) ) * uwdlmt(ji-1,jj) 
+               vn(ji,jj-1,jk) = 0.5_wp * ( 1._wp - sign(1._wp, vn(ji,jj-1,jk) ) * vwdlmt(ji,jj-1) 
+              END DO
+           END DO
+        END DO
+       
+        IF(zflag == 1 .AND. lwp) WRITE(numout,*) 'Need more iterations in wad_lmt!!!'
+       
+        !IF( ln_rnf      )   CALL sbc_rnf_div( hdivn )          ! runoffs (update hdivn field)
+        !IF( nn_cla == 1 )   CALL cla_div    ( kt )             ! Cross Land Advection (update hdivn field)
+        !
+        !
+        CALL wrk_dealloc( jpi, jpj, zflxu, zflxv, zflxu1, zflxv1, uwdlmt, vwdlmt)
+      !
+      END IF
+
+      IF( nn_timing == 1 )  CALL timing_stop('wad_lmt')
+   END SUBROUTINE wad_lmt
+   !!======================================================================
+END MODULE wadlmt
Index: /branches/2013/dev_r4050_NOC_WaD/NEMOGCM/NEMO/OPA_SRC/oce.F90
===================================================================
--- /branches/2013/dev_r4050_NOC_WaD/NEMOGCM/NEMO/OPA_SRC/oce.F90	(revision 4354)
+++ /branches/2013/dev_r4050_NOC_WaD/NEMOGCM/NEMO/OPA_SRC/oce.F90	(revision 4355)
@@ -55,4 +55,5 @@
    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   snwice_fmass       !: time evolution of mass of snow+ice               [Kg/m2/s]
 
+
    !!----------------------------------------------------------------------
    !! NEMO/OPA 4.0 , NEMO Consortium (2011)
@@ -66,5 +67,5 @@
       !!                   ***  FUNCTION oce_alloc  ***
       !!----------------------------------------------------------------------
-      INTEGER :: ierr(3)
+      INTEGER :: ierr(4)
       !!----------------------------------------------------------------------
       !
Index: /branches/2013/dev_r4050_NOC_WaD/NEMOGCM/NEMO/OPA_SRC/step.F90
===================================================================
--- /branches/2013/dev_r4050_NOC_WaD/NEMOGCM/NEMO/OPA_SRC/step.F90	(revision 4354)
+++ /branches/2013/dev_r4050_NOC_WaD/NEMOGCM/NEMO/OPA_SRC/step.F90	(revision 4355)
@@ -87,4 +87,9 @@
 
       !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
+      !  Ocean dynamics : Wetting/Drying Limiter
+      !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
+                         CALL wad_lmt( kstp )         ! wetting/drying limiter
+
+      !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
       ! Update data, open boundaries, surface boundary condition (including sea-ice)
       !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
Index: /branches/2013/dev_r4050_NOC_WaD/NEMOGCM/NEMO/OPA_SRC/step_oce.F90
===================================================================
--- /branches/2013/dev_r4050_NOC_WaD/NEMOGCM/NEMO/OPA_SRC/step_oce.F90	(revision 4354)
+++ /branches/2013/dev_r4050_NOC_WaD/NEMOGCM/NEMO/OPA_SRC/step_oce.F90	(revision 4355)
@@ -116,4 +116,5 @@
    USE trcstp           ! passive tracer time-stepping      (trc_stp routine)
 #endif
+   USE wadlmt           ! wetting/drying limiter            (wad_lmt routine)
    !!----------------------------------------------------------------------
    !! NEMO/OPA 3.3 , NEMO Consortium (2010)
