65
66
67
68
69 integer, intent(in) :: ng, tile
70 integer, intent(in) :: LBi, UBi, LBj, UBj
71 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
72 integer, intent(in) :: krhs, kstp, kout
73
74#ifdef ASSUMED_SHAPE
75 real(r8), intent(inout) :: zeta(LBi:,LBj:,:)
76#else
77 real(r8), intent(inout) :: zeta(LBi:UBi,LBj:UBj,:)
78#endif
79
80
81
82 integer :: i, j, know
83
84 real(r8), parameter :: eps =1.0e-20_r8
85
86 real(r8) :: Ce, Cx
87 real(r8) :: cff, cff1, cff2, dt2d, dZde, dZdt, dZdx, tau
88
89 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: grad
90
91#include "set_bounds.h"
92
93
94
95
96
97#if defined STEP2D_FB_AB3_AM4 || defined STEP2D_FB_LF_AM3
98 know=kstp
99 dt2d=dtfast(ng)
100#else
101 IF (first_2d_step) THEN
102 know=krhs
103 dt2d=dtfast(ng)
104 ELSE IF (predictor_2d_step(ng)) THEN
105 know=krhs
106 dt2d=2.0_r8*dtfast(ng)
107 ELSE
108 know=kstp
109 dt2d=dtfast(ng)
110 END IF
111#endif
112
113
114
115
116
117 IF (domain(ng)%Western_Edge(tile)) THEN
118
119
120
121 IF (lbc(iwest,isfsur,ng)%radiation) THEN
122 DO j=jstr,jend+1
123 grad(istr-1,j)=zeta(istr-1,j ,know)- &
124 & zeta(istr-1,j-1,know)
125#ifdef MASKING
126 grad(istr-1,j)=grad(istr-1,j)*grid(ng)%vmask(istr-1,j)
127#endif
128 grad(istr,j)=zeta(istr,j ,know)- &
129 & zeta(istr,j-1,know)
130#ifdef MASKING
131 grad(istr,j)=grad(istr,j)*grid(ng)%vmask(istr,j)
132#endif
133 END DO
134 DO j=jstr,jend
135 IF (lbc_apply(ng)%west(j)) THEN
136 dzdt=zeta(istr,j,know)-zeta(istr ,j,kout)
137 dzdx=zeta(istr,j,kout)-zeta(istr+1,j,kout)
138
139 IF (lbc(iwest,isfsur,ng)%nudging) THEN
140 IF ((dzdt*dzdx).lt.0.0_r8) THEN
141 tau=fsobc_in(ng,iwest)
142 ELSE
143 tau=fsobc_out(ng,iwest)
144 END IF
145 tau=tau*dt2d
146 END IF
147
148 IF ((dzdt*dzdx).lt.0.0_r8) dzdt=0.0_r8
149 IF ((dzdt*(grad(istr,j)+grad(istr,j+1))).gt.0.0_r8) THEN
150 dzde=grad(istr,j )
151 ELSE
152 dzde=grad(istr,j+1)
153 END IF
154 cff=max(dzdx*dzdx+dzde*dzde,eps)
155 cx=dzdt*dzdx
156#ifdef RADIATION_2D
157 ce=min(cff,max(dzdt*dzde,-cff))
158#else
159 ce=0.0_r8
160#endif
161#if defined CELERITY_WRITE && defined FORWARD_WRITE
162 boundary(ng)%zeta_west_Cx(j)=cx
163 boundary(ng)%zeta_west_Ce(j)=ce
164 boundary(ng)%zeta_west_C2(j)=cff
165#endif
166 zeta(istr-1,j,kout)=(cff*zeta(istr-1,j,know)+ &
167 & cx *zeta(istr ,j,kout)- &
168 & max(ce,0.0_r8)*grad(istr-1,j )- &
169 & min(ce,0.0_r8)*grad(istr-1,j+1))/ &
170 & (cff+cx)
171
172 IF (lbc(iwest,isfsur,ng)%nudging) THEN
173 zeta(istr-1,j,kout)=zeta(istr-1,j,kout)+ &
174 & tau*(boundary(ng)%zeta_west(j)- &
175 & zeta(istr-1,j,know))
176 END IF
177#ifdef MASKING
178 zeta(istr-1,j,kout)=zeta(istr-1,j,kout)* &
179 & grid(ng)%rmask(istr-1,j)
180#endif
181 END IF
182 END DO
183
184
185
186 ELSE IF (lbc(iwest,isfsur,ng)%Chapman_explicit) THEN
187 DO j=jstr,jend
188 IF (lbc_apply(ng)%west(j)) THEN
189 cff=dt2d*grid(ng)%pm(istr,j)
190#ifdef WET_DRY
191 cff1=sqrt(g*(max(grid(ng)%h(istr,j)+ &
192 & zeta(istr,j,know),dcrit(ng))))
193#else
194 cff1=sqrt(g*(grid(ng)%h(istr,j)+ &
195 & zeta(istr,j,know)))
196#endif
197 cx=cff*cff1
198 zeta(istr-1,j,kout)=(1.0_r8-cx)*zeta(istr-1,j,know)+ &
199 & cx*zeta(istr,j,know)
200#ifdef MASKING
201 zeta(istr-1,j,kout)=zeta(istr-1,j,kout)* &
202 & grid(ng)%rmask(istr-1,j)
203#endif
204 END IF
205 END DO
206
207
208
209 ELSE IF (lbc(iwest,isfsur,ng)%Chapman_implicit) THEN
210 DO j=jstr,jend
211 IF (lbc_apply(ng)%west(j)) THEN
212 cff=dt2d*grid(ng)%pm(istr,j)
213#ifdef WET_DRY
214 cff1=sqrt(g*(max(grid(ng)%h(istr,j)+ &
215 & zeta(istr,j,know),dcrit(ng))))
216#else
217 cff1=sqrt(g*(grid(ng)%h(istr,j)+ &
218 & zeta(istr,j,know)))
219#endif
220 cx=cff*cff1
221 cff2=1.0_r8/(1.0_r8+cx)
222 zeta(istr-1,j,kout)=cff2*(zeta(istr-1,j,know)+ &
223 & cx*zeta(istr,j,kout))
224#ifdef MASKING
225 zeta(istr-1,j,kout)=zeta(istr-1,j,kout)* &
226 & grid(ng)%rmask(istr-1,j)
227#endif
228 END IF
229 END DO
230
231
232
233 ELSE IF (lbc(iwest,isfsur,ng)%clamped) THEN
234 DO j=jstr,jend
235 IF (lbc_apply(ng)%west(j)) THEN
236 zeta(istr-1,j,kout)=boundary(ng)%zeta_west(j)
237#ifdef MASKING
238 zeta(istr-1,j,kout)=zeta(istr-1,j,kout)* &
239 & grid(ng)%rmask(istr-1,j)
240#endif
241 END IF
242 END DO
243
244
245
246 ELSE IF (lbc(iwest,isfsur,ng)%gradient) THEN
247 DO j=jstr,jend
248 IF (lbc_apply(ng)%west(j)) THEN
249 zeta(istr-1,j,kout)=zeta(istr,j,kout)
250#ifdef MASKING
251 zeta(istr-1,j,kout)=zeta(istr-1,j,kout)* &
252 & grid(ng)%rmask(istr-1,j)
253#endif
254 END IF
255 END DO
256
257
258
259 ELSE IF (lbc(iwest,isfsur,ng)%closed) THEN
260 DO j=jstr,jend
261 IF (lbc_apply(ng)%west(j)) THEN
262 zeta(istr-1,j,kout)=zeta(istr,j,kout)
263#ifdef MASKING
264 zeta(istr-1,j,kout)=zeta(istr-1,j,kout)* &
265 & grid(ng)%rmask(istr-1,j)
266#endif
267 END IF
268 END DO
269 END IF
270 END IF
271
272
273
274
275
276 IF (domain(ng)%Eastern_Edge(tile)) THEN
277
278
279
280 IF (lbc(ieast,isfsur,ng)%radiation) THEN
281 DO j=jstr,jend+1
282 grad(iend ,j)=zeta(iend ,j ,know)- &
283 & zeta(iend ,j-1,know)
284#ifdef MASKING
285 grad(iend ,j)=grad(iend ,j)*grid(ng)%vmask(iend ,j)
286#endif
287 grad(iend+1,j)=zeta(iend+1,j ,know)- &
288 & zeta(iend+1,j-1,know)
289#ifdef MASKING
290 grad(iend+1,j)=grad(iend+1,j)*grid(ng)%vmask(iend+1,j)
291#endif
292 END DO
293 DO j=jstr,jend
294 IF (lbc_apply(ng)%east(j)) THEN
295 dzdt=zeta(iend,j,know)-zeta(iend ,j,kout)
296 dzdx=zeta(iend,j,kout)-zeta(iend-1,j,kout)
297
298 IF (lbc(ieast,isfsur,ng)%nudging) THEN
299 IF ((dzdt*dzdx).lt.0.0_r8) THEN
300 tau=fsobc_in(ng,ieast)
301 ELSE
302 tau=fsobc_out(ng,ieast)
303 END IF
304 tau=tau*dt2d
305 END IF
306
307 IF ((dzdt*dzdx).lt.0.0_r8) dzdt=0.0_r8
308 IF ((dzdt*(grad(iend,j)+grad(iend,j+1))).gt.0.0_r8) THEN
309 dzde=grad(iend,j )
310 ELSE
311 dzde=grad(iend,j+1)
312 END IF
313 cff=max(dzdx*dzdx+dzde*dzde,eps)
314 cx=dzdt*dzdx
315#ifdef RADIATION_2D
316 ce=min(cff,max(dzdt*dzde,-cff))
317#else
318 ce=0.0_r8
319#endif
320#if defined CELERITY_WRITE && defined FORWARD_WRITE
321 boundary(ng)%zeta_east_Cx(j)=cx
322 boundary(ng)%zeta_east_Ce(j)=ce
323 boundary(ng)%zeta_east_C2(j)=cff
324#endif
325 zeta(iend+1,j,kout)=(cff*zeta(iend+1,j,know)+ &
326 & cx *zeta(iend ,j,kout)- &
327 & max(ce,0.0_r8)*grad(iend+1,j )- &
328 & min(ce,0.0_r8)*grad(iend+1,j+1))/ &
329 & (cff+cx)
330
331 IF (lbc(ieast,isfsur,ng)%nudging) THEN
332 zeta(iend+1,j,kout)=zeta(iend+1,j,kout)+ &
333 & tau*(boundary(ng)%zeta_east(j)- &
334 & zeta(iend+1,j,know))
335 END IF
336#ifdef MASKING
337 zeta(iend+1,j,kout)=zeta(iend+1,j,kout)* &
338 & grid(ng)%rmask(iend+1,j)
339#endif
340 END IF
341 END DO
342
343
344
345 ELSE IF (lbc(ieast,isfsur,ng)%Chapman_explicit) THEN
346 DO j=jstr,jend
347 IF (lbc_apply(ng)%east(j)) THEN
348 cff=dt2d*grid(ng)%pm(iend,j)
349#ifdef WET_DRY
350 cff1=sqrt(g*(max(grid(ng)%h(iend,j)+ &
351 & zeta(iend,j,know),dcrit(ng))))
352#else
353 cff1=sqrt(g*(grid(ng)%h(iend,j)+ &
354 & zeta(iend,j,know)))
355#endif
356 cx=cff*cff1
357 zeta(iend+1,j,kout)=(1.0_r8-cx)*zeta(iend+1,j,know)+ &
358 & cx*zeta(iend,j,know)
359#ifdef MASKING
360 zeta(iend+1,j,kout)=zeta(iend+1,j,kout)* &
361 & grid(ng)%rmask(iend+1,j)
362#endif
363 END IF
364 END DO
365
366
367
368 ELSE IF (lbc(ieast,isfsur,ng)%Chapman_implicit) THEN
369 DO j=jstr,jend
370 IF (lbc_apply(ng)%east(j)) THEN
371 cff=dt2d*grid(ng)%pm(iend,j)
372#ifdef WET_DRY
373 cff1=sqrt(g*(max(grid(ng)%h(iend,j)+ &
374 & zeta(iend,j,know),dcrit(ng))))
375#else
376 cff1=sqrt(g*(grid(ng)%h(iend,j)+ &
377 & zeta(iend,j,know)))
378#endif
379 cx=cff*cff1
380 cff2=1.0_r8/(1.0_r8+cx)
381 zeta(iend+1,j,kout)=cff2*(zeta(iend+1,j,know)+ &
382 & cx*zeta(iend,j,kout))
383#ifdef MASKING
384 zeta(iend+1,j,kout)=zeta(iend+1,j,kout)* &
385 & grid(ng)%rmask(iend+1,j)
386#endif
387 END IF
388 END DO
389
390
391
392 ELSE IF (lbc(ieast,isfsur,ng)%clamped) THEN
393 DO j=jstr,jend
394 IF (lbc_apply(ng)%east(j)) THEN
395 zeta(iend+1,j,kout)=boundary(ng)%zeta_east(j)
396#ifdef MASKING
397 zeta(iend+1,j,kout)=zeta(iend+1,j,kout)* &
398 & grid(ng)%rmask(iend+1,j)
399#endif
400 END IF
401 END DO
402
403
404
405 ELSE IF (lbc(ieast,isfsur,ng)%gradient) THEN
406 DO j=jstr,jend
407 IF (lbc_apply(ng)%east(j)) THEN
408 zeta(iend+1,j,kout)=zeta(iend,j,kout)
409#ifdef MASKING
410 zeta(iend+1,j,kout)=zeta(iend+1,j,kout)* &
411 & grid(ng)%rmask(iend+1,j)
412#endif
413 END IF
414 END DO
415
416
417
418 ELSE IF (lbc(ieast,isfsur,ng)%closed) THEN
419 DO j=jstr,jend
420 IF (lbc_apply(ng)%east(j)) THEN
421 zeta(iend+1,j,kout)=zeta(iend,j,kout)
422#ifdef MASKING
423 zeta(iend+1,j,kout)=zeta(iend+1,j,kout)* &
424 & grid(ng)%rmask(iend+1,j)
425#endif
426 END IF
427 END DO
428 END IF
429 END IF
430
431
432
433
434
435 IF (domain(ng)%Southern_Edge(tile)) THEN
436
437
438
439 IF (lbc(isouth,isfsur,ng)%radiation) THEN
440 DO i=istr,iend+1
441 grad(i,jstr )=zeta(i ,jstr,know)- &
442 & zeta(i-1,jstr,know)
443#ifdef MASKING
444 grad(i,jstr )=grad(i,jstr )*grid(ng)%umask(i,jstr )
445#endif
446 grad(i,jstr-1)=zeta(i ,jstr-1,know)- &
447 & zeta(i-1,jstr-1,know)
448#ifdef MASKING
449 grad(i,jstr-1)=grad(i,jstr-1)*grid(ng)%umask(i,jstr-1)
450#endif
451 END DO
452 DO i=istr,iend
453 IF (lbc_apply(ng)%south(i)) THEN
454 dzdt=zeta(i,jstr,know)-zeta(i,jstr ,kout)
455 dzde=zeta(i,jstr,kout)-zeta(i,jstr-1,kout)
456
457 IF (lbc(isouth,isfsur,ng)%nudging) THEN
458 IF ((dzdt*dzde).lt.0.0_r8) THEN
459 tau=fsobc_in(ng,isouth)
460 ELSE
461 tau=fsobc_out(ng,isouth)
462 END IF
463 tau=tau*dt2d
464 END IF
465
466 IF ((dzdt*dzde).lt.0.0_r8) dzdt=0.0_r8
467 IF ((dzdt*(grad(i,jstr)+grad(i+1,jstr))).gt.0.0_r8) THEN
468 dzdx=grad(i ,jstr)
469 ELSE
470 dzdx=grad(i+1,jstr)
471 END IF
472 cff=max(dzdx*dzdx+dzde*dzde,eps)
473#ifdef RADIATION_2D
474 cx=min(cff,max(dzdt*dzdx,-cff))
475#else
476 cx=0.0_r8
477#endif
478 ce=dzdt*dzde
479#if defined CELERITY_WRITE && defined FORWARD_WRITE
480 boundary(ng)%zeta_south_Cx(i)=cx
481 boundary(ng)%zeta_south_Ce(i)=ce
482 boundary(ng)%zeta_south_C2(i)=cff
483#endif
484 zeta(i,jstr-1,kout)=(cff*zeta(i,jstr-1,know)+ &
485 & ce *zeta(i,jstr ,kout)- &
486 & max(cx,0.0_r8)*grad(i ,jstr)- &
487 & min(cx,0.0_r8)*grad(i+1,jstr))/ &
488 & (cff+ce)
489
490 IF (lbc(isouth,isfsur,ng)%nudging) THEN
491 zeta(i,jstr-1,kout)=zeta(i,jstr-1,kout)+ &
492 & tau*(boundary(ng)%zeta_south(i)- &
493 & zeta(i,jstr-1,know))
494 END IF
495#ifdef MASKING
496 zeta(i,jstr-1,kout)=zeta(i,jstr-1,kout)* &
497 & grid(ng)%rmask(i,jstr-1)
498#endif
499 END IF
500 END DO
501
502
503
504 ELSE IF (lbc(isouth,isfsur,ng)%Chapman_explicit) THEN
505 DO i=istr,iend
506 IF (lbc_apply(ng)%south(i)) THEN
507 cff=dt2d*grid(ng)%pn(i,jstr)
508#ifdef WET_DRY
509 cff1=sqrt(g*(max(grid(ng)%h(i,jstr)+ &
510 & zeta(i,jstr,know),dcrit(ng))))
511#else
512 cff1=sqrt(g*(grid(ng)%h(i,jstr)+ &
513 & zeta(i,jstr,know)))
514#endif
515 ce=cff*cff1
516 zeta(i,jstr-1,kout)=(1.0_r8-ce)*zeta(i,jstr-1,know)+ &
517 & ce*zeta(i,jstr,know)
518#ifdef MASKING
519 zeta(i,jstr-1,kout)=zeta(i,jstr-1,kout)* &
520 & grid(ng)%rmask(i,jstr-1)
521#endif
522 END IF
523 END DO
524
525
526
527 ELSE IF (lbc(isouth,isfsur,ng)%Chapman_implicit) THEN
528 DO i=istr,iend
529 IF (lbc_apply(ng)%south(i)) THEN
530 cff=dt2d*grid(ng)%pn(i,jstr)
531#ifdef WET_DRY
532 cff1=sqrt(g*(max(grid(ng)%h(i,jstr)+ &
533 & zeta(i,jstr,know),dcrit(ng))))
534#else
535 cff1=sqrt(g*(grid(ng)%h(i,jstr)+ &
536 & zeta(i,jstr,know)))
537#endif
538 ce=cff*cff1
539 cff2=1.0_r8/(1.0_r8+ce)
540 zeta(i,jstr-1,kout)=cff2*(zeta(i,jstr-1,know)+ &
541 & ce*zeta(i,jstr,kout))
542#ifdef MASKING
543 zeta(i,jstr-1,kout)=zeta(i,jstr-1,kout)* &
544 & grid(ng)%rmask(i,jstr-1)
545#endif
546 END IF
547 END DO
548
549
550
551 ELSE IF (lbc(isouth,isfsur,ng)%clamped) THEN
552 DO i=istr,iend
553 IF (lbc_apply(ng)%south(i)) THEN
554 zeta(i,jstr-1,kout)=boundary(ng)%zeta_south(i)
555#ifdef MASKING
556 zeta(i,jstr-1,kout)=zeta(i,jstr-1,kout)* &
557 & grid(ng)%rmask(i,jstr-1)
558#endif
559 END IF
560 END DO
561
562
563
564 ELSE IF (lbc(isouth,isfsur,ng)%gradient) THEN
565 DO i=istr,iend
566 IF (lbc_apply(ng)%south(i)) THEN
567 zeta(i,jstr-1,kout)=zeta(i,jstr,kout)
568#ifdef MASKING
569 zeta(i,jstr-1,kout)=zeta(i,jstr-1,kout)* &
570 & grid(ng)%rmask(i,jstr-1)
571#endif
572 END IF
573 END DO
574
575
576
577 ELSE IF (lbc(isouth,isfsur,ng)%closed) THEN
578 DO i=istr,iend
579 IF (lbc_apply(ng)%south(i)) THEN
580 zeta(i,jstr-1,kout)=zeta(i,jstr,kout)
581#ifdef MASKING
582 zeta(i,jstr-1,kout)=zeta(i,jstr-1,kout)* &
583 & grid(ng)%rmask(i,jstr-1)
584#endif
585 END IF
586 END DO
587 END IF
588 END IF
589
590
591
592
593
594 IF (domain(ng)%Northern_Edge(tile)) THEN
595
596
597
598 IF (lbc(inorth,isfsur,ng)%radiation) THEN
599 DO i=istr,iend+1
600 grad(i,jend )=zeta(i ,jend ,know)- &
601 & zeta(i-1,jend ,know)
602#ifdef MASKING
603 grad(i,jend )=grad(i,jend )*grid(ng)%umask(i,jend )
604#endif
605 grad(i,jend+1)=zeta(i ,jend+1,know)- &
606 & zeta(i-1,jend+1,know)
607#ifdef MASKING
608 grad(i,jend+1)=grad(i,jend+1)*grid(ng)%umask(i,jend+1)
609#endif
610 END DO
611 DO i=istr,iend
612 IF (lbc_apply(ng)%north(i)) THEN
613 dzdt=zeta(i,jend,know)-zeta(i,jend ,kout)
614 dzde=zeta(i,jend,kout)-zeta(i,jend-1,kout)
615
616 IF (lbc(inorth,isfsur,ng)%nudging) THEN
617 IF ((dzdt*dzde).lt.0.0_r8) THEN
618 tau=fsobc_in(ng,inorth)
619 ELSE
620 tau=fsobc_out(ng,inorth)
621 END IF
622 tau=tau*dt2d
623 END IF
624
625 IF ((dzdt*dzde).lt.0.0_r8) dzdt=0.0_r8
626 IF ((dzdt*(grad(i,jend)+grad(i+1,jend))).gt.0.0_r8) THEN
627 dzdx=grad(i ,jend)
628 ELSE
629 dzdx=grad(i+1,jend)
630 END IF
631 cff=max(dzdx*dzdx+dzde*dzde,eps)
632#ifdef RADIATION_2D
633 cx=min(cff,max(dzdt*dzdx,-cff))
634#else
635 cx=0.0_r8
636#endif
637 ce=dzdt*dzde
638#if defined CELERITY_WRITE && defined FORWARD_WRITE
639 boundary(ng)%zeta_north_Cx(i)=cx
640 boundary(ng)%zeta_north_Ce(i)=ce
641 boundary(ng)%zeta_north_C2(i)=cff
642#endif
643 zeta(i,jend+1,kout)=(cff*zeta(i,jend+1,know)+ &
644 & ce *zeta(i,jend ,kout)- &
645 & max(cx,0.0_r8)*grad(i ,jend+1)- &
646 & min(cx,0.0_r8)*grad(i+1,jend+1))/ &
647 & (cff+ce)
648
649 IF (lbc(inorth,isfsur,ng)%nudging) THEN
650 zeta(i,jend+1,kout)=zeta(i,jend+1,kout)+ &
651 & tau*(boundary(ng)%zeta_north(i)- &
652 & zeta(i,jend+1,know))
653 END IF
654#ifdef MASKING
655 zeta(i,jend+1,kout)=zeta(i,jend+1,kout)* &
656 & grid(ng)%rmask(i,jend+1)
657#endif
658 END IF
659 END DO
660
661
662
663 ELSE IF (lbc(inorth,isfsur,ng)%Chapman_explicit) THEN
664 DO i=istr,iend
665 IF (lbc_apply(ng)%north(i)) THEN
666 cff=dt2d*grid(ng)%pn(i,jend)
667#ifdef WET_DRY
668 cff1=sqrt(g*(max(grid(ng)%h(i,jend)+ &
669 & zeta(i,jend,know),dcrit(ng))))
670#else
671 cff1=sqrt(g*(grid(ng)%h(i,jend)+ &
672 & zeta(i,jend,know)))
673#endif
674 ce=cff*cff1
675 zeta(i,jend+1,kout)=(1.0_r8-ce)*zeta(i,jend+1,know)+ &
676 & ce*zeta(i,jend,know)
677#ifdef MASKING
678 zeta(i,jend+1,kout)=zeta(i,jend+1,kout)* &
679 & grid(ng)%rmask(i,jend+1)
680#endif
681 END IF
682 END DO
683
684
685
686 ELSE IF (lbc(inorth,isfsur,ng)%Chapman_implicit) THEN
687 DO i=istr,iend
688 IF (lbc_apply(ng)%north(i)) THEN
689 cff=dt2d*grid(ng)%pn(i,jend)
690#ifdef WET_DRY
691 cff1=sqrt(g*(max(grid(ng)%h(i,jend)+ &
692 & zeta(i,jend,know),dcrit(ng))))
693#else
694 cff1=sqrt(g*(grid(ng)%h(i,jend)+ &
695 & zeta(i,jend,know)))
696#endif
697 ce=cff*cff1
698 cff2=1.0_r8/(1.0_r8+ce)
699 zeta(i,jend+1,kout)=cff2*(zeta(i,jend+1,know)+ &
700 & ce*zeta(i,jend,kout))
701#ifdef MASKING
702 zeta(i,jend+1,kout)=zeta(i,jend+1,kout)* &
703 & grid(ng)%rmask(i,jend+1)
704#endif
705 END IF
706 END DO
707
708
709
710 ELSE IF (lbc(inorth,isfsur,ng)%clamped) THEN
711 DO i=istr,iend
712 IF (lbc_apply(ng)%north(i)) THEN
713 zeta(i,jend+1,kout)=boundary(ng)%zeta_north(i)
714#ifdef MASKING
715 zeta(i,jend+1,kout)=zeta(i,jend+1,kout)* &
716 & grid(ng)%rmask(i,jend+1)
717#endif
718 END IF
719 END DO
720
721
722
723 ELSE IF (lbc(inorth,isfsur,ng)%gradient) THEN
724 DO i=istr,iend
725 IF (lbc_apply(ng)%north(i)) THEN
726 zeta(i,jend+1,kout)=zeta(i,jend,kout)
727#ifdef MASKING
728 zeta(i,jend+1,kout)=zeta(i,jend+1,kout)* &
729 & grid(ng)%rmask(i,jend+1)
730#endif
731 END IF
732 END DO
733
734
735
736 ELSE IF (lbc(inorth,isfsur,ng)%closed) THEN
737 DO i=istr,iend
738 IF (lbc_apply(ng)%north(i)) THEN
739 zeta(i,jend+1,kout)=zeta(i,jend,kout)
740#ifdef MASKING
741 zeta(i,jend+1,kout)=zeta(i,jend+1,kout)* &
742 & grid(ng)%rmask(i,jend+1)
743#endif
744 END IF
745 END DO
746 END IF
747 END IF
748
749
750
751
752
753 IF (.not.(ewperiodic(ng).or.nsperiodic(ng))) THEN
754 IF (domain(ng)%SouthWest_Corner(tile)) THEN
755 IF (lbc_apply(ng)%south(istr-1).and. &
756 & lbc_apply(ng)%west (jstr-1)) THEN
757 zeta(istr-1,jstr-1,kout)=0.5_r8*(zeta(istr ,jstr-1,kout)+ &
758 & zeta(istr-1,jstr ,kout))
759 END IF
760 END IF
761 IF (domain(ng)%SouthEast_Corner(tile)) THEN
762 IF (lbc_apply(ng)%south(iend+1).and. &
763 & lbc_apply(ng)%east (jstr-1)) THEN
764 zeta(iend+1,jstr-1,kout)=0.5_r8*(zeta(iend ,jstr-1,kout)+ &
765 & zeta(iend+1,jstr ,kout))
766 END IF
767 END IF
768 IF (domain(ng)%NorthWest_Corner(tile)) THEN
769 IF (lbc_apply(ng)%north(istr-1).and. &
770 & lbc_apply(ng)%west (jend+1)) THEN
771 zeta(istr-1,jend+1,kout)=0.5_r8*(zeta(istr-1,jend ,kout)+ &
772 & zeta(istr ,jend+1,kout))
773 END IF
774 END IF
775 IF (domain(ng)%NorthEast_Corner(tile)) THEN
776 IF (lbc_apply(ng)%north(iend+1).and. &
777 & lbc_apply(ng)%east (jend+1)) THEN
778 zeta(iend+1,jend+1,kout)=0.5_r8*(zeta(iend+1,jend ,kout)+ &
779 & zeta(iend ,jend+1,kout))
780 END IF
781 END IF
782 END IF
783
784#if defined WET_DRY
785
786
787
788
789
790 cff=dcrit(ng)-eps
791 IF (.not.ewperiodic(ng)) THEN
792 IF (domain(ng)%Western_Edge(tile)) THEN
793 DO j=jstr,jend
794 IF (lbc_apply(ng)%west(j)) THEN
795 IF (zeta(istr-1,j,kout).le. &
796 & (dcrit(ng)-grid(ng)%h(istr-1,j))) THEN
797 zeta(istr-1,j,kout)=cff-grid(ng)%h(istr-1,j)
798 END IF
799 END IF
800 END DO
801 END IF
802 IF (domain(ng)%Eastern_Edge(tile)) THEN
803 DO j=jstr,jend
804 IF (lbc_apply(ng)%east(j)) THEN
805 IF (zeta(iend+1,j,kout).le. &
806 & (dcrit(ng)-grid(ng)%h(iend+1,j))) THEN
807 zeta(iend+1,j,kout)=cff-grid(ng)%h(iend+1,j)
808 END IF
809 END IF
810 END DO
811 END IF
812 END IF
813
814 IF (.not.nsperiodic(ng)) THEN
815 IF (domain(ng)%Southern_Edge(tile)) THEN
816 DO i=istr,iend
817 IF (lbc_apply(ng)%south(i)) THEN
818 IF (zeta(i,jstr-1,kout).le. &
819 & (dcrit(ng)-grid(ng)%h(i,jstr-1))) THEN
820 zeta(i,jstr-1,kout)=cff-grid(ng)%h(i,jstr-1)
821 END IF
822 END IF
823 END DO
824 END IF
825 IF (domain(ng)%Northern_Edge(tile)) THEN
826 DO i=istr,iend
827 IF (lbc_apply(ng)%north(i)) THEN
828 IF (zeta(i,jend+1,kout).le. &
829 & (dcrit(ng)-grid(ng)%h(i,jend+1))) THEN
830 zeta(i,jend+1,kout)=cff-grid(ng)%h(i,jend+1)
831 END IF
832 END IF
833 END DO
834 END IF
835 END IF
836
837 IF (.not.(ewperiodic(ng).or.nsperiodic(ng))) THEN
838 IF (domain(ng)%SouthWest_Corner(tile)) THEN
839 IF (lbc_apply(ng)%south(istr-1).and. &
840 & lbc_apply(ng)%west (jstr-1)) THEN
841 IF (zeta(istr-1,jstr-1,kout).le. &
842 & (dcrit(ng)-grid(ng)%h(istr-1,jstr-1))) THEN
843 zeta(istr-1,jstr-1,kout)=cff-grid(ng)%h(istr-1,jstr-1)
844 END IF
845 END IF
846 END IF
847 IF (domain(ng)%SouthEast_Corner(tile)) THEN
848 IF (lbc_apply(ng)%south(iend+1).and. &
849 & lbc_apply(ng)%east (jstr-1)) THEN
850 IF (zeta(iend+1,jstr-1,kout).le. &
851 & (dcrit(ng)-grid(ng)%h(iend+1,jstr-1))) THEN
852 zeta(iend+1,jstr-1,kout)=cff-grid(ng)%h(iend+1,jstr-1)
853 END IF
854 END IF
855 END IF
856 IF (domain(ng)%NorthWest_Corner(tile)) THEN
857 IF (lbc_apply(ng)%north(istr-1).and. &
858 & lbc_apply(ng)%west (jend+1)) THEN
859 IF (zeta(istr-1,jend+1,kout).le. &
860 & (dcrit(ng)-grid(ng)%h(istr-1,jend+1))) THEN
861 zeta(istr-1,jend+1,kout)=cff-grid(ng)%h(istr-1,jend+1)
862 END IF
863 END IF
864 END IF
865 IF (domain(ng)%NorthEast_Corner(tile)) THEN
866 IF (lbc_apply(ng)%north(iend+1).and. &
867 & lbc_apply(ng)%east (jend+1)) THEN
868 IF (zeta(iend+1,jend+1,kout).le. &
869 & (dcrit(ng)-grid(ng)%h(iend+1,jend+1))) THEN
870 zeta(iend+1,jend+1,kout)=cff-grid(ng)%h(iend+1,jend+1)
871 END IF
872 END IF
873 END IF
874 END IF
875#endif
876
877 RETURN