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.
wet_dry.F90 in branches/UKMO/ROMS_WAD_7832/NEMOGCM/NEMO/OPA_SRC/DYN – NEMO

source: branches/UKMO/ROMS_WAD_7832/NEMOGCM/NEMO/OPA_SRC/DYN/wet_dry.F90 @ 8981

Last change on this file since 8981 was 8981, checked in by deazer, 7 years ago

Deal with rn_ssh_ref, this is not a namelist variable but read in from cfg file
change to ssh_ref, set to zero if no cfg file such as test cases
ht_wd_0 looks is obsolete.

File size: 19.5 KB
Line 
1MODULE wet_dry
2
3   !! includes updates to namelist namwad for diagnostic outputs of ROMS wetting and drying
4
5   !!==============================================================================
6   !!                       ***  MODULE  wet_dry  ***
7   !! Wetting and drying includes initialisation routine and routines to
8   !! compute and apply flux limiters and preserve water depth positivity
9   !! only effects if wetting/drying is on (ln_wd_il == .true. or ln_wd_dl==.true. )
10   !!==============================================================================
11   !! History :  3.6  ! 2014-09  ((H.Liu)  Original code
12   !!                 ! will add the runoff and periodic BC case later
13   !!----------------------------------------------------------------------
14
15   !!----------------------------------------------------------------------
16   !!   wad_lmt    : Compute the horizontal flux limiter and the limited velocity
17   !!                when wetting and drying happens
18   !!----------------------------------------------------------------------
19   USE oce             ! ocean dynamics and tracers
20   USE dom_oce         ! ocean space and time domain
21   USE sbc_oce, ONLY : ln_rnf   ! surface boundary condition: ocean
22   USE sbcrnf          ! river runoff
23   USE in_out_manager  ! I/O manager
24   USE lbclnk          ! ocean lateral boundary conditions (or mpp link)
25   USE lib_mpp         ! MPP library
26   USE wrk_nemo        ! Memory Allocation
27   USE timing          ! Timing
28
29   IMPLICIT NONE
30   PRIVATE
31
32   !!----------------------------------------------------------------------
33   !! critical depths,filters, limiters,and masks for  Wetting and Drying
34   !! ---------------------------------------------------------------------
35
36   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) ::   wdmask         !: u- and v- limiter  (can include negative depths)
37   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) ::   wdramp, wdrampu, wdrampv !: for hpg limiting
38
39   LOGICAL,  PUBLIC  ::   ln_wd_il    !: Wetting/drying il activation switch (T:on,F:off)
40   LOGICAL,  PUBLIC  ::   ln_wd_dl    !: Wetting/drying dl activation switch (T:on,F:off)
41   REAL(wp), PUBLIC  ::   rn_wdmin0   !: depth at which wetting/drying starts
42   REAL(wp), PUBLIC  ::   rn_wdmin1   !: minimum water depth on dried cells
43   REAL(wp), PUBLIC  ::   r_rn_wdmin1 !: 1/minimum water depth on dried cells
44   REAL(wp), PUBLIC  ::   rn_wdmin2   !: tolerance of minimum water depth on dried cells
45   REAL(wp), PUBLIC  ::   rn_wdld     !: land elevation below which wetting/drying will be considered
46   INTEGER , PUBLIC  ::   nn_wdit     !: maximum number of iteration for W/D limiter
47   LOGICAL,  PUBLIC  ::   ln_wd_dl_bc !: DL scheme: True implies 3D velocities are set to the barotropic values at points
48                                      !: where the flow is from wet points on less than half the barotropic sub-steps 
49   LOGICAL,  PUBLIC  ::  ln_wd_dl_rmp !: use a ramp for the dl flux limiter between 2 rn_wdmin1 and rn_wdmin1 (rather than a cut-off at rn_wdmin1)     
50   REAL(wp), PUBLIC  ::   ssh_ref     !: height of z=0 with respect to the geoid;
51
52   LOGICAL,  PUBLIC  ::   ll_wd       !: Wetting/drying activation switch if either ln_wd_il or ln_wd_dl
53
54   PUBLIC   wad_init                  ! initialisation routine called by step.F90
55   PUBLIC   wad_lmt                   ! routine called by sshwzv.F90
56   PUBLIC   wad_lmt_bt                ! routine called by dynspg_ts.F90
57
58   !! * Substitutions
59#  include "vectopt_loop_substitute.h90"
60CONTAINS
61
62   SUBROUTINE wad_init
63      !!----------------------------------------------------------------------
64      !!                     ***  ROUTINE wad_init  ***
65      !!                   
66      !! ** Purpose :   read wetting and drying namelist and print the variables.
67      !!
68      !! ** input   : - namwad namelist
69      !!----------------------------------------------------------------------
70      NAMELIST/namwad/ ln_wd_il, ln_wd_dl, rn_wdmin0, rn_wdmin1, rn_wdmin2, rn_wdld, &
71                     & nn_wdit, ln_wd_dl_bc, ln_wd_dl_rmp
72
73      INTEGER  ::   ios                 ! Local integer output status for namelist read
74      INTEGER  ::   ierr                ! Local integer status array allocation
75      !!----------------------------------------------------------------------
76
77      REWIND( numnam_ref )              ! Namelist namwad in reference namelist
78                                        ! : Parameters for Wetting/Drying
79      READ  ( numnam_ref, namwad, IOSTAT = ios, ERR = 905)
80905   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namwad in reference namelist', .TRUE.) 
81      REWIND( numnam_cfg )              ! Namelist namwad in configuration namelist
82                                        ! : Parameters for Wetting/Drying
83      READ  ( numnam_cfg, namwad, IOSTAT = ios, ERR = 906)
84906   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namwad in configuration namelist', .TRUE. )
85      IF(lwm) WRITE ( numond, namwad )
86
87      IF(lwp) THEN                  ! control print
88         WRITE(numout,*)
89         WRITE(numout,*) 'wad_init : Wetting and drying initialization through namelist read'
90         WRITE(numout,*) '~~~~~~~~'
91         WRITE(numout,*) '   Namelist namwad'
92         WRITE(numout,*) '      Logical for Iter Lim wd option   ln_wd_il     = ', ln_wd_il
93         WRITE(numout,*) '      Logical for Dir. Lim wd option   ln_wd_dl     = ', ln_wd_dl
94         WRITE(numout,*) '      Depth at which wet/drying starts rn_wdmin0    = ', rn_wdmin0
95         WRITE(numout,*) '      Minimum wet depth on dried cells rn_wdmin1    = ', rn_wdmin1
96         WRITE(numout,*) '      Tolerance of min wet depth       rn_wdmin2    = ', rn_wdmin2
97         WRITE(numout,*) '      land elevation threshold         rn_wdld      = ', rn_wdld
98         WRITE(numout,*) '      Max iteration for W/D limiter    nn_wdit      = ', nn_wdit
99         WRITE(numout,*) '      T => baroclinic u,v=0 at dry pts: ln_wd_dl_bc = ', ln_wd_dl_bc     
100         WRITE(numout,*) '      use a ramp for rwd limiter:  ln_wd_dl_rwd_rmp = ', ln_wd_dl_rmp
101 
102      ENDIF
103      IF( .NOT. ln_read_cfg ) THEN
104         WRITE(numout,*) '      No configuration file so seting ssh_ref to zero  '
105          ssh_ref=0.0
106      ENDIF
107
108      r_rn_wdmin1=1/rn_wdmin1
109      ll_wd = .FALSE.
110      IF(ln_wd_il .OR. ln_wd_dl) THEN
111         ll_wd = .TRUE.
112         ALLOCATE( wdmask(jpi,jpj),   STAT=ierr )
113         ALLOCATE( wdramp(jpi,jpj), wdrampu(jpi,jpj), wdrampv(jpi,jpj), STAT=ierr ) 
114         IF( ierr /= 0 ) CALL ctl_stop('STOP', 'wad_init : Array allocation error')
115      ENDIF
116      !
117   END SUBROUTINE wad_init
118
119
120   SUBROUTINE wad_lmt( sshb1, sshemp, z2dt )
121      !!----------------------------------------------------------------------
122      !!                  ***  ROUTINE wad_lmt  ***
123      !!                   
124      !! ** Purpose :   generate flux limiters for wetting/drying
125      !!
126      !! ** Method  : - Prevent negative depth occurring (Not ready for Agrif)
127      !!
128      !! ** Action  : - calculate flux limiter and W/D flag
129      !!----------------------------------------------------------------------
130      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   sshb1
131      REAL(wp), DIMENSION(:,:), INTENT(in)    ::   sshemp
132      REAL(wp), INTENT(in) :: z2dt
133      !
134      INTEGER  ::   ji, jj, jk, jk1     ! dummy loop indices
135      INTEGER  ::   jflag               ! local scalar
136      REAL(wp) ::   zcoef, zdep1, zdep2 ! local scalars
137      REAL(wp) ::   zzflxp, zzflxn      ! local scalars
138      REAL(wp) ::   zdepwd              ! local scalar, always wet cell depth
139      REAL(wp) ::   ztmp                ! local scalars
140      REAL(wp), POINTER,  DIMENSION(:,:) ::   zwdlmtu, zwdlmtv         !: W/D flux limiters
141      REAL(wp), POINTER,  DIMENSION(:,:) ::   zflxp,  zflxn            ! local 2D workspace
142      REAL(wp), POINTER,  DIMENSION(:,:) ::   zflxu,  zflxv            ! local 2D workspace
143      REAL(wp), POINTER,  DIMENSION(:,:) ::   zflxu1, zflxv1           ! local 2D workspace
144      !!----------------------------------------------------------------------
145      !
146
147      IF( nn_timing == 1 )  CALL timing_start('wad_lmt')
148
149
150      CALL wrk_alloc( jpi, jpj, zflxp, zflxn, zflxu, zflxv, zflxu1, zflxv1 )
151      CALL wrk_alloc( jpi, jpj, zwdlmtu, zwdlmtv)
152      !
153       
154      jflag  = 0
155      zdepwd = 50._wp   !maximum depth on which that W/D could possibly happen
156
157       
158      zflxp(:,:)   = 0._wp
159      zflxn(:,:)   = 0._wp
160      zflxu(:,:)   = 0._wp
161      zflxv(:,:)   = 0._wp
162
163      zwdlmtu(:,:)  = 1._wp
164      zwdlmtv(:,:)  = 1._wp
165       
166      ! Horizontal Flux in u and v direction
167      DO jk = 1, jpkm1 
168         DO jj = 1, jpj
169            DO ji = 1, jpi
170               zflxu(ji,jj) = zflxu(ji,jj) + e3u_n(ji,jj,jk) * un(ji,jj,jk) * umask(ji,jj,jk)
171               zflxv(ji,jj) = zflxv(ji,jj) + e3v_n(ji,jj,jk) * vn(ji,jj,jk) * vmask(ji,jj,jk)
172            END DO 
173         END DO 
174      END DO
175       
176      zflxu(:,:) = zflxu(:,:) * e2u(:,:)
177      zflxv(:,:) = zflxv(:,:) * e1v(:,:)
178       
179      wdmask(:,:) = 1
180      DO jj = 2, jpj
181         DO ji = 2, jpi 
182
183            IF( tmask(ji, jj, 1) < 0.5_wp ) CYCLE              ! we don't care about land cells
184            IF( ht_0(ji,jj) - ssh_ref > zdepwd ) CYCLE   ! and cells which are unlikely to dry
185
186            zflxp(ji,jj) = max(zflxu(ji,jj), 0._wp) - min(zflxu(ji-1,jj),   0._wp) + &
187                         & max(zflxv(ji,jj), 0._wp) - min(zflxv(ji,  jj-1), 0._wp) 
188            zflxn(ji,jj) = min(zflxu(ji,jj), 0._wp) - max(zflxu(ji-1,jj),   0._wp) + &
189                         & min(zflxv(ji,jj), 0._wp) - max(zflxv(ji,  jj-1), 0._wp) 
190
191            zdep2 = ht_0(ji,jj) + sshb1(ji,jj) - rn_wdmin1
192            IF(zdep2 .le. 0._wp) THEN  !add more safty, but not necessary
193               sshb1(ji,jj) = rn_wdmin1 - ht_0(ji,jj)
194               IF(zflxu(ji,  jj) > 0._wp) zwdlmtu(ji  ,jj) = 0._wp
195               IF(zflxu(ji-1,jj) < 0._wp) zwdlmtu(ji-1,jj) = 0._wp
196               IF(zflxv(ji,  jj) > 0._wp) zwdlmtv(ji  ,jj) = 0._wp
197               IF(zflxv(ji,jj-1) < 0._wp) zwdlmtv(ji,jj-1) = 0._wp 
198               wdmask(ji,jj) = 0._wp
199            END IF
200         ENDDO
201      END DO
202
203
204! HPG limiter from jholt
205      wdramp(:,:) = min((ht_0(:,:) + sshb1(:,:) - rn_wdmin1)/(rn_wdmin0 - rn_wdmin1),1.0_wp)
206!jth assume don't need a lbc_lnk here
207      DO jj = 1, jpjm1
208         DO ji = 1, jpim1
209            wdrampu(ji,jj) = min(wdramp(ji,jj),wdramp(ji+1,jj))
210            wdrampv(ji,jj) = min(wdramp(ji,jj),wdramp(ji,jj+1))
211         END DO
212      END DO
213! end HPG limiter
214
215
216     
217        !! start limiter iterations
218      DO jk1 = 1, nn_wdit + 1
219       
220         
221         zflxu1(:,:) = zflxu(:,:) * zwdlmtu(:,:)
222         zflxv1(:,:) = zflxv(:,:) * zwdlmtv(:,:)
223         jflag = 0     ! flag indicating if any further iterations are needed
224         
225         DO jj = 2, jpj
226            DO ji = 2, jpi 
227       
228               IF( tmask(ji, jj, 1) < 0.5_wp ) CYCLE
229               IF( ht_0(ji,jj) > zdepwd )      CYCLE
230       
231               ztmp = e1e2t(ji,jj)
232
233               zzflxp = max(zflxu1(ji,jj), 0._wp) - min(zflxu1(ji-1,jj),   0._wp) + &
234                      & max(zflxv1(ji,jj), 0._wp) - min(zflxv1(ji,  jj-1), 0._wp) 
235               zzflxn = min(zflxu1(ji,jj), 0._wp) - max(zflxu1(ji-1,jj),   0._wp) + &
236                      & min(zflxv1(ji,jj), 0._wp) - max(zflxv1(ji,  jj-1), 0._wp) 
237         
238               zdep1 = (zzflxp + zzflxn) * z2dt / ztmp
239               zdep2 = ht_0(ji,jj) + sshb1(ji,jj) - rn_wdmin1 - z2dt * sshemp(ji,jj)
240         
241               IF( zdep1 > zdep2 ) THEN
242                  wdmask(ji, jj) = 0
243                  zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zflxp(ji,jj) * z2dt )
244                  !zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zzflxp * z2dt )
245                  ! flag if the limiter has been used but stop flagging if the only
246                  ! changes have zeroed the coefficient since further iterations will
247                  ! not change anything
248                  IF( zcoef > 0._wp ) THEN
249                     jflag = 1 
250                  ELSE
251                     zcoef = 0._wp
252                  ENDIF
253                  IF(jk1 > nn_wdit) zcoef = 0._wp
254                  IF(zflxu1(ji,  jj) > 0._wp) zwdlmtu(ji  ,jj) = zcoef
255                  IF(zflxu1(ji-1,jj) < 0._wp) zwdlmtu(ji-1,jj) = zcoef
256                  IF(zflxv1(ji,  jj) > 0._wp) zwdlmtv(ji  ,jj) = zcoef
257                  IF(zflxv1(ji,jj-1) < 0._wp) zwdlmtv(ji,jj-1) = zcoef
258               END IF
259            END DO ! ji loop
260         END DO  ! jj loop
261
262         CALL lbc_lnk( zwdlmtu, 'U', 1. )
263         CALL lbc_lnk( zwdlmtv, 'V', 1. )
264
265         IF(lk_mpp) CALL mpp_max(jflag)   !max over the global domain
266
267         IF(jflag == 0) EXIT
268         
269      END DO  ! jk1 loop
270       
271      DO jk = 1, jpkm1
272         un(:,:,jk) = un(:,:,jk) * zwdlmtu(:, :) 
273         vn(:,:,jk) = vn(:,:,jk) * zwdlmtv(:, :) 
274      END DO
275
276      CALL lbc_lnk( un, 'U', -1. )
277      CALL lbc_lnk( vn, 'V', -1. )
278        !
279      un_b(:,:) = un_b(:,:) * zwdlmtu(:, :)
280      vn_b(:,:) = vn_b(:,:) * zwdlmtv(:, :)
281      CALL lbc_lnk( un_b, 'U', -1. )
282      CALL lbc_lnk( vn_b, 'V', -1. )
283       
284      IF(jflag == 1 .AND. lwp) WRITE(numout,*) 'Need more iterations in wad_lmt!!!'
285       
286      !IF( ln_rnf      )   CALL sbc_rnf_div( hdivn )          ! runoffs (update hdivn field)
287      !IF( nn_cla == 1 )   CALL cla_div    ( kt )             ! Cross Land Advection (update hdivn field)
288      !
289      !
290      CALL wrk_dealloc( jpi, jpj, zflxp, zflxn, zflxu, zflxv, zflxu1, zflxv1 )
291      CALL wrk_dealloc( jpi, jpj, zwdlmtu, zwdlmtv)
292      !
293      !
294      IF( nn_timing == 1 )  CALL timing_stop('wad_lmt')
295      !
296   END SUBROUTINE wad_lmt
297
298
299   SUBROUTINE wad_lmt_bt( zflxu, zflxv, sshn_e, zssh_frc, rdtbt )
300      !!----------------------------------------------------------------------
301      !!                  ***  ROUTINE wad_lmt  ***
302      !!                   
303      !! ** Purpose :   limiting flux in the barotropic stepping (dynspg_ts)
304      !!
305      !! ** Method  : - Prevent negative depth occurring (Not ready for Agrif)
306      !!
307      !! ** Action  : - calculate flux limiter and W/D flag
308      !!----------------------------------------------------------------------
309      REAL(wp), INTENT(in)    ::   rdtbt    ! ocean time-step index
310      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   zflxu,  zflxv, sshn_e, zssh_frc 
311      !
312      INTEGER  ::   ji, jj, jk, jk1     ! dummy loop indices
313      INTEGER  ::   jflag               ! local scalar
314      REAL(wp) ::   z2dt
315      REAL(wp) ::   zcoef, zdep1, zdep2 ! local scalars
316      REAL(wp) ::   zzflxp, zzflxn      ! local scalars
317      REAL(wp) ::   zdepwd              ! local scalar, always wet cell depth
318      REAL(wp) ::   ztmp                ! local scalars
319      REAL(wp), POINTER,  DIMENSION(:,:) ::   zwdlmtu, zwdlmtv         !: W/D flux limiters
320      REAL(wp), POINTER,  DIMENSION(:,:) ::   zflxp,  zflxn            ! local 2D workspace
321      REAL(wp), POINTER,  DIMENSION(:,:) ::   zflxu1, zflxv1           ! local 2D workspace
322      !!----------------------------------------------------------------------
323      !
324      IF( nn_timing == 1 )  CALL timing_start('wad_lmt_bt')
325
326
327      CALL wrk_alloc( jpi, jpj, zflxp, zflxn, zflxu1, zflxv1 )
328      CALL wrk_alloc( jpi, jpj, zwdlmtu, zwdlmtv)
329      !
330       
331      !IF(lwp) WRITE(numout,*)
332      !IF(lwp) WRITE(numout,*) 'wad_lmt_bt : wetting/drying limiters and velocity limiting'
333       
334      jflag  = 0
335      zdepwd = 50._wp   !maximum depth that ocean cells can have W/D processes
336
337      z2dt = rdtbt   
338       
339      zflxp(:,:)   = 0._wp
340      zflxn(:,:)   = 0._wp
341
342      zwdlmtu(:,:)  = 1._wp
343      zwdlmtv(:,:)  = 1._wp
344       
345        ! Horizontal Flux in u and v direction
346       
347      DO jj = 2, jpj
348         DO ji = 2, jpi 
349
350           IF( tmask(ji, jj, 1 ) < 0.5_wp) CYCLE   ! we don't care about land cells
351           IF( ht_0(ji,jj) > zdepwd )      CYCLE   ! and cells which are unlikely to dry
352
353            zflxp(ji,jj) = max(zflxu(ji,jj), 0._wp) - min(zflxu(ji-1,jj),   0._wp) + &
354                         & max(zflxv(ji,jj), 0._wp) - min(zflxv(ji,  jj-1), 0._wp) 
355            zflxn(ji,jj) = min(zflxu(ji,jj), 0._wp) - max(zflxu(ji-1,jj),   0._wp) + &
356                         & min(zflxv(ji,jj), 0._wp) - max(zflxv(ji,  jj-1), 0._wp) 
357
358            zdep2 = ht_0(ji,jj) + sshn_e(ji,jj) - rn_wdmin1
359            IF(zdep2 .le. 0._wp) THEN  !add more safety, but not necessary
360              sshn_e(ji,jj) = rn_wdmin1 - ht_0(ji,jj)
361              IF(zflxu(ji,  jj) > 0._wp) zwdlmtu(ji  ,jj) = 0._wp
362              IF(zflxu(ji-1,jj) < 0._wp) zwdlmtu(ji-1,jj) = 0._wp
363              IF(zflxv(ji,  jj) > 0._wp) zwdlmtv(ji  ,jj) = 0._wp
364              IF(zflxv(ji,jj-1) < 0._wp) zwdlmtv(ji,jj-1) = 0._wp 
365            END IF
366         ENDDO
367      END DO
368
369     
370      !! start limiter iterations
371      DO jk1 = 1, nn_wdit + 1
372       
373         
374         zflxu1(:,:) = zflxu(:,:) * zwdlmtu(:,:)
375         zflxv1(:,:) = zflxv(:,:) * zwdlmtv(:,:)
376         jflag = 0     ! flag indicating if any further iterations are needed
377         
378         DO jj = 2, jpj
379            DO ji = 2, jpi 
380       
381               IF( tmask(ji, jj, 1 ) < 0.5_wp) CYCLE
382               IF( ht_0(ji,jj) > zdepwd )      CYCLE
383       
384               ztmp = e1e2t(ji,jj)
385
386               zzflxp = max(zflxu1(ji,jj), 0._wp) - min(zflxu1(ji-1,jj),   0._wp) + &
387                      & max(zflxv1(ji,jj), 0._wp) - min(zflxv1(ji,  jj-1), 0._wp) 
388               zzflxn = min(zflxu1(ji,jj), 0._wp) - max(zflxu1(ji-1,jj),   0._wp) + &
389                      & min(zflxv1(ji,jj), 0._wp) - max(zflxv1(ji,  jj-1), 0._wp) 
390         
391               zdep1 = (zzflxp + zzflxn) * z2dt / ztmp
392               zdep2 = ht_0(ji,jj) + sshn_e(ji,jj) - rn_wdmin1 - z2dt * zssh_frc(ji,jj)
393         
394               IF(zdep1 > zdep2) THEN
395                 zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zflxp(ji,jj) * z2dt )
396                 !zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zzflxp * z2dt )
397                 ! flag if the limiter has been used but stop flagging if the only
398                 ! changes have zeroed the coefficient since further iterations will
399                 ! not change anything
400                 IF( zcoef > 0._wp ) THEN
401                    jflag = 1 
402                 ELSE
403                    zcoef = 0._wp
404                 ENDIF
405                 IF(jk1 > nn_wdit) zcoef = 0._wp
406                 IF(zflxu1(ji,  jj) > 0._wp) zwdlmtu(ji  ,jj) = zcoef
407                 IF(zflxu1(ji-1,jj) < 0._wp) zwdlmtu(ji-1,jj) = zcoef
408                 IF(zflxv1(ji,  jj) > 0._wp) zwdlmtv(ji  ,jj) = zcoef
409                 IF(zflxv1(ji,jj-1) < 0._wp) zwdlmtv(ji,jj-1) = zcoef
410               END IF
411            END DO ! ji loop
412         END DO  ! jj loop
413
414         CALL lbc_lnk( zwdlmtu, 'U', 1. )
415         CALL lbc_lnk( zwdlmtv, 'V', 1. )
416
417         IF(lk_mpp) CALL mpp_max(jflag)   !max over the global domain
418
419         IF(jflag == 0) EXIT
420         
421      END DO  ! jk1 loop
422       
423      zflxu(:,:) = zflxu(:,:) * zwdlmtu(:, :) 
424      zflxv(:,:) = zflxv(:,:) * zwdlmtv(:, :) 
425
426      CALL lbc_lnk( zflxu, 'U', -1. )
427      CALL lbc_lnk( zflxv, 'V', -1. )
428       
429      IF(jflag == 1 .AND. lwp) WRITE(numout,*) 'Need more iterations in wad_lmt_bt!!!'
430       
431      !IF( ln_rnf      )   CALL sbc_rnf_div( hdivn )          ! runoffs (update hdivn field)
432      !IF( nn_cla == 1 )   CALL cla_div    ( kt )             ! Cross Land Advection (update hdivn field)
433      !
434      !
435      CALL wrk_dealloc( jpi, jpj, zflxp, zflxn, zflxu1, zflxv1 )
436      CALL wrk_dealloc( jpi, jpj, zwdlmtu, zwdlmtv)
437      !
438
439      IF( nn_timing == 1 )  CALL timing_stop('wad_lmt')
440   END SUBROUTINE wad_lmt_bt
441
442   !!==============================================================================
443END MODULE wet_dry
Note: See TracBrowser for help on using the repository browser.