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 NEMO/trunk/src/OCE/DYN – NEMO

source: NEMO/trunk/src/OCE/DYN/wet_dry.F90 @ 10069

Last change on this file since 10069 was 10069, checked in by nicolasmartin, 6 years ago

Fix mistakes of previous commit on SVN keywords property

  • Property svn:keywords set to Id
File size: 19.2 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_init      : initialisation of wetting and drying
17   !!   wad_lmt       : horizontal flux limiter and limited velocity when wetting and drying happens
18   !!   wad_lmt_bt    : same as wad_lmt for the barotropic stepping (dynspg_ts)
19   !!----------------------------------------------------------------------
20   USE oce            ! ocean dynamics and tracers
21   USE dom_oce        ! ocean space and time domain
22   USE sbc_oce  , ONLY: ln_rnf   ! surface boundary condition: ocean
23   USE sbcrnf         ! river runoff
24   !
25   USE in_out_manager ! I/O manager
26   USE lbclnk         ! ocean lateral boundary conditions (or mpp link)
27   USE lib_mpp        ! MPP library
28   USE timing         ! timing of the main modules
29
30   IMPLICIT NONE
31   PRIVATE
32
33   !!----------------------------------------------------------------------
34   !! critical depths,filters, limiters,and masks for  Wetting and Drying
35   !! ---------------------------------------------------------------------
36
37   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) ::   wdmask   !: u- and v- limiter
38   !                                                           !  (can include negative depths)
39   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:) ::   wdramp, wdrampu, wdrampv !: for hpg limiting
40
41   LOGICAL,  PUBLIC  ::   ln_wd_il    !: Wetting/drying il activation switch (T:on,F:off)
42   LOGICAL,  PUBLIC  ::   ln_wd_dl    !: Wetting/drying dl activation switch (T:on,F:off)
43   REAL(wp), PUBLIC  ::   rn_wdmin0   !: depth at which wetting/drying starts
44   REAL(wp), PUBLIC  ::   rn_wdmin1   !: minimum water depth on dried cells
45   REAL(wp), PUBLIC  ::   r_rn_wdmin1 !: 1/minimum water depth on dried cells
46   REAL(wp), PUBLIC  ::   rn_wdmin2   !: tolerance of minimum water depth on dried cells
47   REAL(wp), PUBLIC  ::   rn_wdld     !: land elevation below which wetting/drying will be considered
48   INTEGER , PUBLIC  ::   nn_wdit     !: maximum number of iteration for W/D limiter
49   LOGICAL,  PUBLIC  ::   ln_wd_dl_bc !: DL scheme: True implies 3D velocities are set to the barotropic values at points
50                                      !: where the flow is from wet points on less than half the barotropic sub-steps 
51   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)     
52   REAL(wp), PUBLIC  ::   ssh_ref     !: height of z=0 with respect to the geoid;
53
54   LOGICAL,  PUBLIC  ::   ll_wd       !: Wetting/drying activation switch if either ln_wd_il or ln_wd_dl
55
56   PUBLIC   wad_init                  ! initialisation routine called by step.F90
57   PUBLIC   wad_lmt                   ! routine called by sshwzv.F90
58   PUBLIC   wad_lmt_bt                ! routine called by dynspg_ts.F90
59
60   !! * Substitutions
61#  include "vectopt_loop_substitute.h90"
62   !!----------------------------------------------------------------------
63CONTAINS
64
65   SUBROUTINE wad_init
66      !!----------------------------------------------------------------------
67      !!                     ***  ROUTINE wad_init  ***
68      !!                   
69      !! ** Purpose :   read wetting and drying namelist and print the variables.
70      !!
71      !! ** input   : - namwad namelist
72      !!----------------------------------------------------------------------
73      INTEGER  ::   ios, ierr   ! Local integer
74      !!
75      NAMELIST/namwad/ ln_wd_il, ln_wd_dl   , rn_wdmin0, rn_wdmin1, rn_wdmin2, rn_wdld,   &
76         &             nn_wdit , ln_wd_dl_bc, ln_wd_dl_rmp
77      !!----------------------------------------------------------------------
78      !
79      REWIND( numnam_ref )              ! Namelist namwad in reference namelist : Parameters for Wetting/Drying
80      READ  ( numnam_ref, namwad, IOSTAT = ios, ERR = 905)
81905   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namwad in reference namelist', .TRUE.) 
82      REWIND( numnam_cfg )              ! Namelist namwad in configuration namelist : 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      ENDIF
102      IF( .NOT. ln_read_cfg ) THEN
103         IF(lwp) WRITE(numout,*) '      No configuration file so seting ssh_ref to zero  '
104         ssh_ref=0._wp
105      ENDIF
106
107      r_rn_wdmin1 = 1 / rn_wdmin1
108      ll_wd = .FALSE.
109      IF( ln_wd_il .OR. ln_wd_dl ) THEN
110         ll_wd = .TRUE.
111         ALLOCATE( wdmask(jpi,jpj),   STAT=ierr )
112         ALLOCATE( wdramp(jpi,jpj), wdrampu(jpi,jpj), wdrampv(jpi,jpj), STAT=ierr ) 
113         IF( ierr /= 0 ) CALL ctl_stop('STOP', 'wad_init : Array allocation error')
114      ENDIF
115      !
116   END SUBROUTINE wad_init
117
118
119   SUBROUTINE wad_lmt( sshb1, sshemp, z2dt )
120      !!----------------------------------------------------------------------
121      !!                  ***  ROUTINE wad_lmt  ***
122      !!                   
123      !! ** Purpose :   generate flux limiters for wetting/drying
124      !!
125      !! ** Method  : - Prevent negative depth occurring (Not ready for Agrif)
126      !!
127      !! ** Action  : - calculate flux limiter and W/D flag
128      !!----------------------------------------------------------------------
129      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   sshb1        !!gm DOCTOR names: should start with p !
130      REAL(wp), DIMENSION(:,:), INTENT(in   ) ::   sshemp
131      REAL(wp)                , INTENT(in   ) ::   z2dt
132      !
133      INTEGER  ::   ji, jj, jk, jk1     ! dummy loop indices
134      INTEGER  ::   jflag               ! local scalar
135      REAL(wp) ::   zcoef, zdep1, zdep2 ! local scalars
136      REAL(wp) ::   zzflxp, zzflxn      ! local scalars
137      REAL(wp) ::   zdepwd              ! local scalar, always wet cell depth
138      REAL(wp) ::   ztmp                ! local scalars
139      REAL(wp),  DIMENSION(jpi,jpj) ::   zwdlmtu, zwdlmtv   ! W/D flux limiters
140      REAL(wp),  DIMENSION(jpi,jpj) ::   zflxp  ,  zflxn    ! local 2D workspace
141      REAL(wp),  DIMENSION(jpi,jpj) ::   zflxu  ,  zflxv    ! local 2D workspace
142      REAL(wp),  DIMENSION(jpi,jpj) ::   zflxu1 , zflxv1    ! local 2D workspace
143      !!----------------------------------------------------------------------
144      IF( ln_timing )   CALL timing_start('wad_lmt')      !
145      !
146      DO jk = 1, jpkm1
147         un(:,:,jk) = un(:,:,jk) * zwdlmtu(:,:) 
148         vn(:,:,jk) = vn(:,:,jk) * zwdlmtv(:,:) 
149      END DO
150      jflag  = 0
151      zdepwd = 50._wp      ! maximum depth on which that W/D could possibly happen
152      !
153      zflxp(:,:)   = 0._wp
154      zflxn(:,:)   = 0._wp
155      zflxu(:,:)   = 0._wp
156      zflxv(:,:)   = 0._wp
157      !
158      zwdlmtu(:,:) = 1._wp
159      zwdlmtv(:,:) = 1._wp
160      !
161      DO jk = 1, jpkm1     ! Horizontal Flux in u and v direction
162         zflxu(:,:) = zflxu(:,:) + e3u_n(:,:,jk) * un(:,:,jk) * umask(:,:,jk)
163         zflxv(:,:) = zflxv(:,:) + e3v_n(:,:,jk) * vn(:,:,jk) * vmask(:,:,jk)
164      END DO
165      zflxu(:,:) = zflxu(:,:) * e2u(:,:)
166      zflxv(:,:) = zflxv(:,:) * e1v(:,:)
167      !
168      wdmask(:,:) = 1._wp
169      DO jj = 2, jpj
170         DO ji = 2, jpi 
171            !
172            IF( tmask(ji,jj,1)        < 0.5_wp )   CYCLE    ! we don't care about land cells
173            IF( ht_0(ji,jj) - ssh_ref > zdepwd )   CYCLE    ! and cells which are unlikely to dry
174            !
175            zflxp(ji,jj) = MAX( zflxu(ji,jj) , 0._wp ) - MIN( zflxu(ji-1,jj  ) , 0._wp )   &
176               &         + MAX( zflxv(ji,jj) , 0._wp ) - MIN( zflxv(ji,  jj-1) , 0._wp ) 
177            zflxn(ji,jj) = MIN( zflxu(ji,jj) , 0._wp ) - MAX( zflxu(ji-1,jj  ) , 0._wp )   &
178               &         + MIN( zflxv(ji,jj) , 0._wp ) - MAX( zflxv(ji,  jj-1) , 0._wp ) 
179            !
180            zdep2 = ht_0(ji,jj) + sshb1(ji,jj) - rn_wdmin1
181            IF( zdep2 <= 0._wp ) THEN     ! add more safty, but not necessary
182               sshb1(ji,jj) = rn_wdmin1 - ht_0(ji,jj)
183               IF(zflxu(ji,  jj) > 0._wp) zwdlmtu(ji  ,jj) = 0._wp
184               IF(zflxu(ji-1,jj) < 0._wp) zwdlmtu(ji-1,jj) = 0._wp
185               IF(zflxv(ji,  jj) > 0._wp) zwdlmtv(ji  ,jj) = 0._wp
186               IF(zflxv(ji,jj-1) < 0._wp) zwdlmtv(ji,jj-1) = 0._wp 
187               wdmask(ji,jj) = 0._wp
188            END IF
189         END DO
190      END DO
191      !
192      !           ! HPG limiter from jholt
193      wdramp(:,:) = min((ht_0(:,:) + sshb1(:,:) - rn_wdmin1)/(rn_wdmin0 - rn_wdmin1),1.0_wp)
194      !jth assume don't need a lbc_lnk here
195      DO jj = 1, jpjm1
196         DO ji = 1, jpim1
197            wdrampu(ji,jj) = MIN( wdramp(ji,jj) , wdramp(ji+1,jj) )
198            wdrampv(ji,jj) = MIN( wdramp(ji,jj) , wdramp(ji,jj+1) )
199         END DO
200      END DO
201      !           ! end HPG limiter
202      !
203      !
204      DO jk1 = 1, nn_wdit + 1      !==  start limiter iterations  ==!
205         !
206         zflxu1(:,:) = zflxu(:,:) * zwdlmtu(:,:)
207         zflxv1(:,:) = zflxv(:,:) * zwdlmtv(:,:)
208         jflag = 0     ! flag indicating if any further iterations are needed
209         !
210         DO jj = 2, jpj
211            DO ji = 2, jpi 
212               IF( tmask(ji, jj, 1) < 0.5_wp )   CYCLE
213               IF( ht_0(ji,jj)      > zdepwd )   CYCLE
214               !
215               ztmp = e1e2t(ji,jj)
216               !
217               zzflxp = MAX( zflxu1(ji,jj) , 0._wp ) - MIN( zflxu1(ji-1,jj  ) , 0._wp)   &
218                  &   + MAX( zflxv1(ji,jj) , 0._wp ) - MIN( zflxv1(ji,  jj-1) , 0._wp) 
219               zzflxn = MIN( zflxu1(ji,jj) , 0._wp ) - MAX( zflxu1(ji-1,jj  ) , 0._wp)   &
220                  &   + MIN( zflxv1(ji,jj) , 0._wp ) - MAX( zflxv1(ji,  jj-1) , 0._wp) 
221               !
222               zdep1 = (zzflxp + zzflxn) * z2dt / ztmp
223               zdep2 = ht_0(ji,jj) + sshb1(ji,jj) - rn_wdmin1 - z2dt * sshemp(ji,jj)
224               !
225               IF( zdep1 > zdep2 ) THEN
226                  wdmask(ji, jj) = 0._wp
227                  zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zflxp(ji,jj) * z2dt )
228                  !zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zzflxp * z2dt )
229                  ! flag if the limiter has been used but stop flagging if the only
230                  ! changes have zeroed the coefficient since further iterations will
231                  ! not change anything
232                  IF( zcoef > 0._wp ) THEN   ;   jflag = 1 
233                  ELSE                       ;   zcoef = 0._wp
234                  ENDIF
235                  IF( jk1 > nn_wdit )   zcoef = 0._wp
236                  IF( zflxu1(ji  ,jj  ) > 0._wp )   zwdlmtu(ji  ,jj  ) = zcoef
237                  IF( zflxu1(ji-1,jj  ) < 0._wp )   zwdlmtu(ji-1,jj  ) = zcoef
238                  IF( zflxv1(ji  ,jj  ) > 0._wp )   zwdlmtv(ji  ,jj  ) = zcoef
239                  IF( zflxv1(ji  ,jj-1) < 0._wp )   zwdlmtv(ji  ,jj-1) = zcoef
240               ENDIF
241            END DO
242         END DO
243         CALL lbc_lnk_multi( zwdlmtu, 'U', 1., zwdlmtv, 'V', 1. )
244         !
245         IF( lk_mpp )   CALL mpp_max(jflag)   !max over the global domain
246         !
247         IF( jflag == 0 )   EXIT
248         !
249      END DO  ! jk1 loop
250      !
251      DO jk = 1, jpkm1
252         un(:,:,jk) = un(:,:,jk) * zwdlmtu(:,:) 
253         vn(:,:,jk) = vn(:,:,jk) * zwdlmtv(:,:) 
254      END DO
255      un_b(:,:) = un_b(:,:) * zwdlmtu(:, :)
256      vn_b(:,:) = vn_b(:,:) * zwdlmtv(:, :)
257      !
258!!gm TO BE SUPPRESSED ?  these lbc_lnk are useless since zwdlmtu and zwdlmtv are defined everywhere !
259      CALL lbc_lnk_multi( un  , 'U', -1., vn  , 'V', -1. )
260      CALL lbc_lnk_multi( un_b, 'U', -1., vn_b, 'V', -1. )
261!!gm
262      !
263      IF(jflag == 1 .AND. lwp)   WRITE(numout,*) 'Need more iterations in wad_lmt!!!'
264      !
265      !IF( ln_rnf      )   CALL sbc_rnf_div( hdivn )          ! runoffs (update hdivn field)
266      !
267      IF( ln_timing )   CALL timing_stop('wad_lmt')      !
268      !
269   END SUBROUTINE wad_lmt
270
271
272   SUBROUTINE wad_lmt_bt( zflxu, zflxv, sshn_e, zssh_frc, rdtbt )
273      !!----------------------------------------------------------------------
274      !!                  ***  ROUTINE wad_lmt  ***
275      !!                   
276      !! ** Purpose :   limiting flux in the barotropic stepping (dynspg_ts)
277      !!
278      !! ** Method  : - Prevent negative depth occurring (Not ready for Agrif)
279      !!
280      !! ** Action  : - calculate flux limiter and W/D flag
281      !!----------------------------------------------------------------------
282      REAL(wp)                , INTENT(in   ) ::   rdtbt    ! ocean time-step index
283      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   zflxu,  zflxv, sshn_e, zssh_frc 
284      !
285      INTEGER  ::   ji, jj, jk, jk1     ! dummy loop indices
286      INTEGER  ::   jflag               ! local integer
287      REAL(wp) ::   z2dt
288      REAL(wp) ::   zcoef, zdep1, zdep2 ! local scalars
289      REAL(wp) ::   zzflxp, zzflxn      ! local scalars
290      REAL(wp) ::   zdepwd              ! local scalar, always wet cell depth
291      REAL(wp) ::   ztmp                ! local scalars
292      REAL(wp), DIMENSION(jpi,jpj) ::   zwdlmtu, zwdlmtv         !: W/D flux limiters
293      REAL(wp), DIMENSION(jpi,jpj) ::   zflxp,  zflxn            ! local 2D workspace
294      REAL(wp), DIMENSION(jpi,jpj) ::   zflxu1, zflxv1           ! local 2D workspace
295      !!----------------------------------------------------------------------
296      IF( ln_timing )   CALL timing_start('wad_lmt_bt')      !
297      !
298      jflag  = 0
299      zdepwd = 50._wp   ! maximum depth that ocean cells can have W/D processes
300      !
301      z2dt = rdtbt   
302      !
303      zflxp(:,:)   = 0._wp
304      zflxn(:,:)   = 0._wp
305      zwdlmtu(:,:) = 1._wp
306      zwdlmtv(:,:) = 1._wp
307      !
308      DO jj = 2, jpj      ! Horizontal Flux in u and v direction
309         DO ji = 2, jpi 
310            !
311            IF( tmask(ji, jj, 1 ) < 0.5_wp) CYCLE   ! we don't care about land cells
312            IF( ht_0(ji,jj) > zdepwd )      CYCLE   ! and cells which are unlikely to dry
313            !
314            zflxp(ji,jj) = MAX( zflxu(ji,jj) , 0._wp ) - MIN( zflxu(ji-1,jj  ) , 0._wp )   &
315               &         + MAX( zflxv(ji,jj) , 0._wp ) - MIN( zflxv(ji,  jj-1) , 0._wp ) 
316            zflxn(ji,jj) = MIN( zflxu(ji,jj) , 0._wp ) - MAX( zflxu(ji-1,jj  ) , 0._wp )   &
317               &         + MIN( zflxv(ji,jj) , 0._wp ) - MAX( zflxv(ji,  jj-1) , 0._wp ) 
318            !
319            zdep2 = ht_0(ji,jj) + sshn_e(ji,jj) - rn_wdmin1
320            IF( zdep2 <= 0._wp ) THEN  !add more safety, but not necessary
321              sshn_e(ji,jj) = rn_wdmin1 - ht_0(ji,jj)
322              IF( zflxu(ji  ,jj  ) > 0._wp)   zwdlmtu(ji  ,jj  ) = 0._wp
323              IF( zflxu(ji-1,jj  ) < 0._wp)   zwdlmtu(ji-1,jj  ) = 0._wp
324              IF( zflxv(ji  ,jj  ) > 0._wp)   zwdlmtv(ji  ,jj  ) = 0._wp
325              IF( zflxv(ji  ,jj-1) < 0._wp)   zwdlmtv(ji  ,jj-1) = 0._wp 
326            ENDIF
327         END DO
328      END DO
329      !
330      DO jk1 = 1, nn_wdit + 1      !! start limiter iterations
331         !
332         zflxu1(:,:) = zflxu(:,:) * zwdlmtu(:,:)
333         zflxv1(:,:) = zflxv(:,:) * zwdlmtv(:,:)
334         jflag = 0     ! flag indicating if any further iterations are needed
335         !
336         DO jj = 2, jpj
337            DO ji = 2, jpi 
338               !
339               IF( tmask(ji, jj, 1 ) < 0.5_wp )   CYCLE
340               IF( ht_0(ji,jj)       > zdepwd )   CYCLE
341               !
342               ztmp = e1e2t(ji,jj)
343               !
344               zzflxp = max(zflxu1(ji,jj), 0._wp) - min(zflxu1(ji-1,jj),   0._wp)   &
345                  &   + max(zflxv1(ji,jj), 0._wp) - min(zflxv1(ji,  jj-1), 0._wp) 
346               zzflxn = min(zflxu1(ji,jj), 0._wp) - max(zflxu1(ji-1,jj),   0._wp)   &
347                  &   + min(zflxv1(ji,jj), 0._wp) - max(zflxv1(ji,  jj-1), 0._wp) 
348         
349               zdep1 = (zzflxp + zzflxn) * z2dt / ztmp
350               zdep2 = ht_0(ji,jj) + sshn_e(ji,jj) - rn_wdmin1 - z2dt * zssh_frc(ji,jj)
351         
352               IF(zdep1 > zdep2) THEN
353                 zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zflxp(ji,jj) * z2dt )
354                 !zcoef = ( ( zdep2 - rn_wdmin2 ) * ztmp - zzflxn * z2dt ) / ( zzflxp * z2dt )
355                 ! flag if the limiter has been used but stop flagging if the only
356                 ! changes have zeroed the coefficient since further iterations will
357                 ! not change anything
358                 IF( zcoef > 0._wp ) THEN
359                    jflag = 1 
360                 ELSE
361                    zcoef = 0._wp
362                 ENDIF
363                 IF(jk1 > nn_wdit) zcoef = 0._wp
364                 IF(zflxu1(ji,  jj) > 0._wp) zwdlmtu(ji  ,jj) = zcoef
365                 IF(zflxu1(ji-1,jj) < 0._wp) zwdlmtu(ji-1,jj) = zcoef
366                 IF(zflxv1(ji,  jj) > 0._wp) zwdlmtv(ji  ,jj) = zcoef
367                 IF(zflxv1(ji,jj-1) < 0._wp) zwdlmtv(ji,jj-1) = zcoef
368               END IF
369            END DO ! ji loop
370         END DO  ! jj loop
371         !
372         CALL lbc_lnk_multi( zwdlmtu, 'U', 1., zwdlmtv, 'V', 1. )
373         !
374         IF(lk_mpp) CALL mpp_max(jflag)   !max over the global domain
375         !
376         IF(jflag == 0)   EXIT
377         !
378      END DO  ! jk1 loop
379      !
380      zflxu(:,:) = zflxu(:,:) * zwdlmtu(:, :) 
381      zflxv(:,:) = zflxv(:,:) * zwdlmtv(:, :) 
382      !
383!!gm THIS lbc_lnk is useless since it is already done at the end of the jk1-loop
384      CALL lbc_lnk_multi( zflxu, 'U', -1., zflxv, 'V', -1. )
385!!gm end
386      !
387      IF( jflag == 1 .AND. lwp )   WRITE(numout,*) 'Need more iterations in wad_lmt_bt!!!'
388      !
389      !IF( ln_rnf      )   CALL sbc_rnf_div( hdivn )          ! runoffs (update hdivn field)
390      !
391      IF( ln_timing )   CALL timing_stop('wad_lmt_bt')      !
392      !
393   END SUBROUTINE wad_lmt_bt
394
395   !!==============================================================================
396END MODULE wet_dry
Note: See TracBrowser for help on using the repository browser.