ROMS
Loading...
Searching...
No Matches
zetabc.F
Go to the documentation of this file.
1#include "cppdefs.h"
2 MODULE zetabc_mod
3!
4!git $Id$
5!================================================== Hernan G. Arango ===
6! Copyright (c) 2002-2025 The ROMS Group !
7! Licensed under a MIT/X style license !
8! See License_ROMS.md !
9!=======================================================================
10! !
11! This routine sets lateral boundary conditions for free-surface. !
12#if defined STEP2D_FB_AB3_AM4 || defined STEP2D_FB_LF_AM3
13! Notice that "zetabc_local" is used for either the Forward-Backward !
14! AB3-AM4 or Forward-Backward LF-AM3 barotropic kernels where the !
15! boundary conditions are loaded into private array "zeta_new". !
16#endif
17! !
18!=======================================================================
19!
20 USE mod_param
21 USE mod_boundary
22 USE mod_grid
23 USE mod_ncparam
24 USE mod_ocean
25 USE mod_scalars
26 USE mod_stepping
27!
28 implicit none
29!
30 PRIVATE
31 PUBLIC :: zetabc_tile
32#if defined STEP2D_FB_AB3_AM4 || defined STEP2D_FB_LF_AM3
33 PUBLIC :: zetabc_local
34#endif
35!
36 CONTAINS
37!
38!***********************************************************************
39 SUBROUTINE zetabc (ng, tile, kout)
40!***********************************************************************
41!
42! Imported variable declarations.
43!
44 integer, intent(in) :: ng, tile, kout
45!
46! Local variable declarations.
47!
48# include "tile.h"
49!
50 CALL zetabc_tile (ng, tile, &
51 & lbi, ubi, lbj, ubj, &
52 & imins, imaxs, jmins, jmaxs, &
53 & krhs(ng), kstp(ng), kout, &
54 & ocean(ng) % zeta)
55
56 RETURN
57 END SUBROUTINE zetabc
58!
59!***********************************************************************
60 SUBROUTINE zetabc_tile (ng, tile, &
61 & LBi, UBi, LBj, UBj, &
62 & IminS, ImaxS, JminS, JmaxS, &
63 & krhs, kstp, kout, &
64 & zeta)
65!***********************************************************************
66!
67! Imported variable declarations.
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! Local variable declarations.
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! Set time-indices
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! Lateral boundary conditions at the western edge.
115!-----------------------------------------------------------------------
116!
117 IF (domain(ng)%Western_Edge(tile)) THEN
118!
119! Western edge, implicit upstream radiation condition.
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! Western edge, explicit Chapman boundary condition.
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! Western edge, implicit Chapman boundary condition.
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! Western edge, clamped boundary condition.
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! Western edge, gradient boundary condition.
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! Western edge, closed boundary condition.
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! Lateral boundary conditions at the eastern edge.
274!-----------------------------------------------------------------------
275!
276 IF (domain(ng)%Eastern_Edge(tile)) THEN
277!
278! Eastern edge, implicit upstream radiation condition.
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! Eastern edge, explicit Chapman boundary condition.
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! Eastern edge, implicit Chapman boundary condition.
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! Eastern edge, clamped boundary condition.
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! Eastern edge, gradient boundary condition.
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! Eastern edge, closed boundary condition.
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! Lateral boundary conditions at the southern edge.
433!-----------------------------------------------------------------------
434!
435 IF (domain(ng)%Southern_Edge(tile)) THEN
436!
437! Southern edge, implicit upstream radiation condition.
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! Southern edge, explicit Chapman boundary condition.
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! Southern edge, implicit Chapman boundary condition.
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! Southern edge, clamped boundary condition.
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! Southern edge, gradient boundary condition.
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! Southern edge, closed boundary condition.
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! Lateral boundary conditions at the northern edge.
592!-----------------------------------------------------------------------
593!
594 IF (domain(ng)%Northern_Edge(tile)) THEN
595!
596! Northern edge, implicit upstream radiation condition.
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! Northern edge, explicit Chapman boundary condition.
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! Northern edge, implicit Chapman boundary condition.
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! Northern edge, clamped boundary condition.
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! Northern edge, gradient boundary condition.
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! Northern edge, closed boundary condition.
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! Boundary corners.
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! Ensure that water level on boundary cells is above bed elevation.
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
878 END SUBROUTINE zetabc_tile
879
880#if defined STEP2D_FB_AB3_AM4 || defined STEP2D_FB_LF_AM3
881!
882!***********************************************************************
883 SUBROUTINE zetabc_local (ng, tile, &
884 & LBi, UBi, LBj, UBj, &
885 & IminS, ImaxS, JminS, JmaxS, &
886 & kstp, &
887 & zeta, &
888 & zeta_new)
889!***********************************************************************
890!
891! Imported variable declarations.
892!
893 integer, intent(in) :: ng, tile
894 integer, intent(in) :: lbi, ubi, lbj, ubj
895 integer, intent(in) :: imins, imaxs, jmins, jmaxs
896 integer, intent(in) :: kstp
897!
898# ifdef ASSUMED_SHAPE
899 real(r8), intent(in ) :: zeta(lbi:,lbj:,:)
900 real(r8), intent(inout) :: zeta_new(imins:,jmins:)
901# else
902 real(r8), intent(in ) :: zeta(lbi:ubi,lbj:ubj,:)
903 real(r8), intent(inout) :: zeta_new(imins:imaxs,jmins:jmaxs)
904# endif
905!
906! Local variable declarations.
907!
908 integer :: i, j
909!
910 real(r8), parameter :: eps =1.0e-20_r8
911
912 real(r8) :: ce, cx
913 real(r8) :: cff, cff1, cff2, dzde, dzdt, dzdx, tau
914
915 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: grad
916
917# include "set_bounds.h"
918!
919!-----------------------------------------------------------------------
920! Lateral boundary conditions at the western edge.
921!-----------------------------------------------------------------------
922!
923 IF (domain(ng)%Western_Edge(tile)) THEN
924!
925! Western edge, implicit upstream radiation condition.
926!
927 IF (lbc(iwest,isfsur,ng)%radiation) THEN
928 DO j=jstrv-1,jend+1
929 grad(istr-1,j)=zeta(istr-1,j ,kstp)- &
930 & zeta(istr-1,j-1,kstp)
931# ifdef MASKING
932 grad(istr-1,j)=grad(istr-1,j)*grid(ng)%vmask(istr-1,j)
933# endif
934 grad(istr,j)=zeta(istr,j ,kstp)- &
935 & zeta(istr,j-1,kstp)
936# ifdef MASKING
937 grad(istr,j)=grad(istr,j)*grid(ng)%vmask(istr,j)
938# endif
939 END DO
940 DO j=jstrv-1,jend
941 IF (lbc_apply(ng)%west(j)) THEN
942 dzdt=zeta(istr,j,kstp)-zeta_new(istr ,j)
943 dzdx=zeta_new(istr,j) -zeta_new(istr+1,j)
944
945 IF (lbc(iwest,isfsur,ng)%nudging) THEN
946 IF ((dzdt*dzdx).lt.0.0_r8) THEN
947 tau=fsobc_in(ng,iwest)
948 ELSE
949 tau=fsobc_out(ng,iwest)
950 END IF
951 tau=tau*dtfast(ng)
952 END IF
953
954 IF ((dzdt*dzdx).lt.0.0_r8) dzdt=0.0_r8
955 IF ((dzdt*(grad(istr,j)+grad(istr,j+1))).gt.0.0_r8) THEN
956 dzde=grad(istr,j )
957 ELSE
958 dzde=grad(istr,j+1)
959 END IF
960 cff=max(dzdx*dzdx+dzde*dzde,eps)
961 cx=dzdt*dzdx
962# ifdef RADIATION_2D
963 ce=min(cff,max(dzdt*dzde,-cff))
964# else
965 ce=0.0_r8
966# endif
967# if defined CELERITY_WRITE && defined FORWARD_WRITE
968 boundary(ng)%zeta_west_Cx(j)=cx
969 boundary(ng)%zeta_west_Ce(j)=ce
970 boundary(ng)%zeta_west_C2(j)=cff
971# endif
972 zeta_new(istr-1,j)=(cff*zeta(istr-1,j,kstp)+ &
973 & cx *zeta_new(istr,j)- &
974 & max(ce,0.0_r8)*grad(istr-1,j )- &
975 & min(ce,0.0_r8)*grad(istr-1,j+1))/ &
976 & (cff+cx)
977
978 IF (lbc(iwest,isfsur,ng)%nudging) THEN
979 zeta_new(istr-1,j)=zeta_new(istr-1,j)+ &
980 & tau*(boundary(ng)%zeta_west(j)- &
981 & zeta(istr-1,j,kstp))
982 END IF
983# ifdef MASKING
984 zeta_new(istr-1,j)=zeta_new(istr-1,j)* &
985 & grid(ng)%rmask(istr-1,j)
986# endif
987 END IF
988 END DO
989!
990! Western edge, explicit Chapman boundary condition.
991!
992 ELSE IF (lbc(iwest,isfsur,ng)%Chapman_explicit) THEN
993 DO j=jstrv-1,jend
994 IF (lbc_apply(ng)%west(j)) THEN
995 cff=dtfast(ng)*grid(ng)%pm(istr,j)
996# ifdef WET_DRY
997 cff1=sqrt(g*(max(grid(ng)%h(istr,j)+ &
998 & zeta(istr,j,kstp),dcrit(ng))))
999# else
1000 cff1=sqrt(g*(grid(ng)%h(istr,j)+ &
1001 & zeta(istr,j,kstp)))
1002# endif
1003 cx=cff*cff1
1004 zeta_new(istr-1,j)=(1.0_r8-cx)*zeta(istr-1,j,kstp)+ &
1005 & cx*zeta(istr,j,kstp)
1006# ifdef MASKING
1007 zeta_new(istr-1,j)=zeta_new(istr-1,j)* &
1008 & grid(ng)%rmask(istr-1,j)
1009# endif
1010 END IF
1011 END DO
1012!
1013! Western edge, implicit Chapman boundary condition.
1014!
1015 ELSE IF (lbc(iwest,isfsur,ng)%Chapman_implicit) THEN
1016 DO j=jstrv-1,jend
1017 IF (lbc_apply(ng)%west(j)) THEN
1018 cff=dtfast(ng)*grid(ng)%pm(istr,j)
1019# ifdef WET_DRY
1020 cff1=sqrt(g*(max(grid(ng)%h(istr,j)+ &
1021 & zeta(istr,j,kstp),dcrit(ng))))
1022# else
1023 cff1=sqrt(g*(grid(ng)%h(istr,j)+ &
1024 & zeta(istr,j,kstp)))
1025# endif
1026 cx=cff*cff1
1027 cff2=1.0_r8/(1.0_r8+cx)
1028 zeta_new(istr-1,j)=cff2*(zeta(istr-1,j,kstp)+ &
1029 & cx*zeta_new(istr,j))
1030# ifdef MASKING
1031 zeta_new(istr-1,j)=zeta_new(istr-1,j)* &
1032 & grid(ng)%rmask(istr-1,j)
1033# endif
1034 END IF
1035 END DO
1036!
1037! Western edge, clamped boundary condition.
1038!
1039 ELSE IF (lbc(iwest,isfsur,ng)%clamped) THEN
1040 DO j=jstrv-1,jend
1041 IF (lbc_apply(ng)%west(j)) THEN
1042 zeta_new(istr-1,j)=boundary(ng)%zeta_west(j)
1043# ifdef MASKING
1044 zeta_new(istr-1,j)=zeta_new(istr-1,j)* &
1045 & grid(ng)%rmask(istr-1,j)
1046# endif
1047 END IF
1048 END DO
1049!
1050! Western edge, gradient boundary condition.
1051!
1052 ELSE IF (lbc(iwest,isfsur,ng)%gradient) THEN
1053 DO j=jstrv-1,jend
1054 IF (lbc_apply(ng)%west(j)) THEN
1055 zeta_new(istr-1,j)=zeta_new(istr,j)
1056# ifdef MASKING
1057 zeta_new(istr-1,j)=zeta_new(istr-1,j)* &
1058 & grid(ng)%rmask(istr-1,j)
1059# endif
1060 END IF
1061 END DO
1062!
1063! Western edge, closed boundary condition.
1064!
1065 ELSE IF (lbc(iwest,isfsur,ng)%closed) THEN
1066 DO j=jstrv-1,jend
1067 IF (lbc_apply(ng)%west(j)) THEN
1068 zeta_new(istr-1,j)=zeta_new(istr,j)
1069# ifdef MASKING
1070 zeta_new(istr-1,j)=zeta_new(istr-1,j)* &
1071 & grid(ng)%rmask(istr-1,j)
1072# endif
1073 END IF
1074 END DO
1075 END IF
1076 END IF
1077!
1078!-----------------------------------------------------------------------
1079! Lateral boundary conditions at the eastern edge.
1080!-----------------------------------------------------------------------
1081!
1082 IF (domain(ng)%Eastern_Edge(tile)) THEN
1083!
1084! Eastern edge, implicit upstream radiation condition.
1085!
1086 IF (lbc(ieast,isfsur,ng)%radiation) THEN
1087 DO j=jstrv-1,jend+1
1088 grad(iend ,j)=zeta(iend ,j ,kstp)- &
1089 & zeta(iend ,j-1,kstp)
1090# ifdef MASKING
1091 grad(iend ,j)=grad(iend ,j)*grid(ng)%vmask(iend ,j)
1092# endif
1093 grad(iend+1,j)=zeta(iend+1,j ,kstp)- &
1094 & zeta(iend+1,j-1,kstp)
1095# ifdef MASKING
1096 grad(iend+1,j)=grad(iend+1,j)*grid(ng)%vmask(iend+1,j)
1097# endif
1098 END DO
1099 DO j=jstrv-1,jend
1100 IF (lbc_apply(ng)%east(j)) THEN
1101 dzdt=zeta(iend,j,kstp)-zeta_new(iend ,j)
1102 dzdx=zeta_new(iend,j) -zeta_new(iend-1,j)
1103
1104 IF (lbc(ieast,isfsur,ng)%nudging) THEN
1105 IF ((dzdt*dzdx).lt.0.0_r8) THEN
1106 tau=fsobc_in(ng,ieast)
1107 ELSE
1108 tau=fsobc_out(ng,ieast)
1109 END IF
1110 tau=tau*dtfast(ng)
1111 END IF
1112
1113 IF ((dzdt*dzdx).lt.0.0_r8) dzdt=0.0_r8
1114 IF ((dzdt*(grad(iend,j)+grad(iend,j+1))).gt.0.0_r8) THEN
1115 dzde=grad(iend,j )
1116 ELSE
1117 dzde=grad(iend,j+1)
1118 END IF
1119 cff=max(dzdx*dzdx+dzde*dzde,eps)
1120 cx=dzdt*dzdx
1121# ifdef RADIATION_2D
1122 ce=min(cff,max(dzdt*dzde,-cff))
1123# else
1124 ce=0.0_r8
1125# endif
1126# if defined CELERITY_WRITE && defined FORWARD_WRITE
1127 boundary(ng)%zeta_east_Cx(j)=cx
1128 boundary(ng)%zeta_east_Ce(j)=ce
1129 boundary(ng)%zeta_east_C2(j)=cff
1130# endif
1131 zeta_new(iend+1,j)=(cff*zeta(iend+1,j,kstp)+ &
1132 & cx *zeta_new(iend,j)- &
1133 & max(ce,0.0_r8)*grad(iend+1,j )- &
1134 & min(ce,0.0_r8)*grad(iend+1,j+1))/ &
1135 & (cff+cx)
1136
1137 IF (lbc(ieast,isfsur,ng)%nudging) THEN
1138 zeta_new(iend+1,j)=zeta_new(iend+1,j)+ &
1139 & tau*(boundary(ng)%zeta_east(j)- &
1140 & zeta(iend+1,j,kstp))
1141 END IF
1142# ifdef MASKING
1143 zeta_new(iend+1,j)=zeta_new(iend+1,j)* &
1144 & grid(ng)%rmask(iend+1,j)
1145# endif
1146 END IF
1147 END DO
1148!
1149! Eastern edge, explicit Chapman boundary condition.
1150!
1151 ELSE IF (lbc(ieast,isfsur,ng)%Chapman_explicit) THEN
1152 DO j=jstrv-1,jend
1153 IF (lbc_apply(ng)%east(j)) THEN
1154 cff=dtfast(ng)*grid(ng)%pm(iend,j)
1155# ifdef WET_DRY
1156 cff1=sqrt(g*(max(grid(ng)%h(iend,j)+ &
1157 & zeta(iend,j,kstp),dcrit(ng))))
1158# else
1159 cff1=sqrt(g*(grid(ng)%h(iend,j)+ &
1160 & zeta(iend,j,kstp)))
1161# endif
1162 cx=cff*cff1
1163 zeta_new(iend+1,j)=(1.0_r8-cx)*zeta(iend+1,j,kstp)+ &
1164 & cx*zeta(iend,j,kstp)
1165# ifdef MASKING
1166 zeta_new(iend+1,j)=zeta_new(iend+1,j)* &
1167 & grid(ng)%rmask(iend+1,j)
1168# endif
1169 END IF
1170 END DO
1171!
1172! Eastern edge, implicit Chapman boundary condition.
1173!
1174 ELSE IF (lbc(ieast,isfsur,ng)%Chapman_implicit) THEN
1175 DO j=jstrv-1,jend
1176 IF (lbc_apply(ng)%east(j)) THEN
1177 cff=dtfast(ng)*grid(ng)%pm(iend,j)
1178# ifdef WET_DRY
1179 cff1=sqrt(g*(max(grid(ng)%h(iend,j)+ &
1180 & zeta(iend,j,kstp),dcrit(ng))))
1181# else
1182 cff1=sqrt(g*(grid(ng)%h(iend,j)+ &
1183 & zeta(iend,j,kstp)))
1184# endif
1185 cx=cff*cff1
1186 cff2=1.0_r8/(1.0_r8+cx)
1187 zeta_new(iend+1,j)=cff2*(zeta(iend+1,j,kstp)+ &
1188 & cx*zeta_new(iend,j))
1189# ifdef MASKING
1190 zeta_new(iend+1,j)=zeta_new(iend+1,j)* &
1191 & grid(ng)%rmask(iend+1,j)
1192# endif
1193 END IF
1194 END DO
1195!
1196! Eastern edge, clamped boundary condition.
1197!
1198 ELSE IF (lbc(ieast,isfsur,ng)%clamped) THEN
1199 DO j=jstrv-1,jend
1200 IF (lbc_apply(ng)%east(j)) THEN
1201 zeta_new(iend+1,j)=boundary(ng)%zeta_east(j)
1202# ifdef MASKING
1203 zeta_new(iend+1,j)=zeta_new(iend+1,j)* &
1204 & grid(ng)%rmask(iend+1,j)
1205# endif
1206 END IF
1207 END DO
1208!
1209! Eastern edge, gradient boundary condition.
1210!
1211 ELSE IF (lbc(ieast,isfsur,ng)%gradient) THEN
1212 DO j=jstrv-1,jend
1213 IF (lbc_apply(ng)%east(j)) THEN
1214 zeta_new(iend+1,j)=zeta_new(iend,j)
1215# ifdef MASKING
1216 zeta_new(iend+1,j)=zeta_new(iend+1,j)* &
1217 & grid(ng)%rmask(iend+1,j)
1218# endif
1219 END IF
1220 END DO
1221!
1222! Eastern edge, closed boundary condition.
1223!
1224 ELSE IF (lbc(ieast,isfsur,ng)%closed) THEN
1225 DO j=jstrv-1,jend
1226 IF (lbc_apply(ng)%east(j)) THEN
1227 zeta_new(iend+1,j)=zeta_new(iend,j)
1228# ifdef MASKING
1229 zeta_new(iend+1,j)=zeta_new(iend+1,j)* &
1230 & grid(ng)%rmask(iend+1,j)
1231# endif
1232 END IF
1233 END DO
1234 END IF
1235 END IF
1236!
1237!-----------------------------------------------------------------------
1238! Lateral boundary conditions at the southern edge.
1239!-----------------------------------------------------------------------
1240!
1241 IF (domain(ng)%Southern_Edge(tile)) THEN
1242!
1243! Southern edge, implicit upstream radiation condition.
1244!
1245 IF (lbc(isouth,isfsur,ng)%radiation) THEN
1246 DO i=istru-1,iend+1
1247 grad(i,jstr )=zeta(i ,jstr,kstp)- &
1248 & zeta(i-1,jstr,kstp)
1249# ifdef MASKING
1250 grad(i,jstr )=grad(i,jstr )*grid(ng)%umask(i,jstr )
1251# endif
1252 grad(i,jstr-1)=zeta(i ,jstr-1,kstp)- &
1253 & zeta(i-1,jstr-1,kstp)
1254# ifdef MASKING
1255 grad(i,jstr-1)=grad(i,jstr-1)*grid(ng)%umask(i,jstr-1)
1256# endif
1257 END DO
1258 DO i=istru-1,iend
1259 IF (lbc_apply(ng)%south(i)) THEN
1260 dzdt=zeta(i,jstr,kstp)-zeta_new(i,jstr )
1261 dzde=zeta_new(i,jstr) -zeta_new(i,jstr-1)
1262
1263 IF (lbc(isouth,isfsur,ng)%nudging) THEN
1264 IF ((dzdt*dzde).lt.0.0_r8) THEN
1265 tau=fsobc_in(ng,isouth)
1266 ELSE
1267 tau=fsobc_out(ng,isouth)
1268 END IF
1269 tau=tau*dtfast(ng)
1270 END IF
1271
1272 IF ((dzdt*dzde).lt.0.0_r8) dzdt=0.0_r8
1273 IF ((dzdt*(grad(i,jstr)+grad(i+1,jstr))).gt.0.0_r8) THEN
1274 dzdx=grad(i ,jstr)
1275 ELSE
1276 dzdx=grad(i+1,jstr)
1277 END IF
1278 cff=max(dzdx*dzdx+dzde*dzde,eps)
1279# ifdef RADIATION_2D
1280 cx=min(cff,max(dzdt*dzdx,-cff))
1281# else
1282 cx=0.0_r8
1283# endif
1284 ce=dzdt*dzde
1285# if defined CELERITY_WRITE && defined FORWARD_WRITE
1286 boundary(ng)%zeta_south_Cx(i)=cx
1287 boundary(ng)%zeta_south_Ce(i)=ce
1288 boundary(ng)%zeta_south_C2(i)=cff
1289# endif
1290 zeta_new(i,jstr-1)=(cff*zeta(i,jstr-1,kstp)+ &
1291 & ce *zeta_new(i,jstr)- &
1292 & max(cx,0.0_r8)*grad(i ,jstr)- &
1293 & min(cx,0.0_r8)*grad(i+1,jstr))/ &
1294 & (cff+ce)
1295
1296 IF (lbc(isouth,isfsur,ng)%nudging) THEN
1297 zeta_new(i,jstr-1)=zeta_new(i,jstr-1)+ &
1298 & tau*(boundary(ng)%zeta_south(i)- &
1299 & zeta(i,jstr-1,kstp))
1300 END IF
1301# ifdef MASKING
1302 zeta_new(i,jstr-1)=zeta_new(i,jstr-1)* &
1303 & grid(ng)%rmask(i,jstr-1)
1304# endif
1305 END IF
1306 END DO
1307!
1308! Southern edge, explicit Chapman boundary condition.
1309!
1310 ELSE IF (lbc(isouth,isfsur,ng)%Chapman_explicit) THEN
1311 DO i=istru-1,iend
1312 IF (lbc_apply(ng)%south(i)) THEN
1313 cff=dtfast(ng)*grid(ng)%pn(i,jstr)
1314# ifdef WET_DRY
1315 cff1=sqrt(g*(max(grid(ng)%h(i,jstr)+ &
1316 & zeta(i,jstr,kstp),dcrit(ng))))
1317# else
1318 cff1=sqrt(g*(grid(ng)%h(i,jstr)+ &
1319 & zeta(i,jstr,kstp)))
1320# endif
1321 ce=cff*cff1
1322 zeta_new(i,jstr-1)=(1.0_r8-ce)*zeta(i,jstr-1,kstp)+ &
1323 & ce*zeta(i,jstr,kstp)
1324# ifdef MASKING
1325 zeta_new(i,jstr-1)=zeta_new(i,jstr-1)* &
1326 & grid(ng)%rmask(i,jstr-1)
1327# endif
1328 END IF
1329 END DO
1330!
1331! Southern edge, implicit Chapman boundary condition.
1332!
1333 ELSE IF (lbc(isouth,isfsur,ng)%Chapman_implicit) THEN
1334 DO i=istru-1,iend
1335 IF (lbc_apply(ng)%south(i)) THEN
1336 cff=dtfast(ng)*grid(ng)%pn(i,jstr)
1337# ifdef WET_DRY
1338 cff1=sqrt(g*(max(grid(ng)%h(i,jstr)+ &
1339 & zeta(i,jstr,kstp),dcrit(ng))))
1340# else
1341 cff1=sqrt(g*(grid(ng)%h(i,jstr)+ &
1342 & zeta(i,jstr,kstp)))
1343# endif
1344 ce=cff*cff1
1345 cff2=1.0_r8/(1.0_r8+ce)
1346 zeta_new(i,jstr-1)=cff2*(zeta(i,jstr-1,kstp)+ &
1347 & ce*zeta_new(i,jstr))
1348# ifdef MASKING
1349 zeta_new(i,jstr-1)=zeta_new(i,jstr-1)* &
1350 & grid(ng)%rmask(i,jstr-1)
1351# endif
1352 END IF
1353 END DO
1354!
1355! Southern edge, clamped boundary condition.
1356!
1357 ELSE IF (lbc(isouth,isfsur,ng)%clamped) THEN
1358 DO i=istru-1,iend
1359 IF (lbc_apply(ng)%south(i)) THEN
1360 zeta_new(i,jstr-1)=boundary(ng)%zeta_south(i)
1361# ifdef MASKING
1362 zeta_new(i,jstr-1)=zeta_new(i,jstr-1)* &
1363 & grid(ng)%rmask(i,jstr-1)
1364# endif
1365 END IF
1366 END DO
1367!
1368! Southern edge, gradient boundary condition.
1369!
1370 ELSE IF (lbc(isouth,isfsur,ng)%gradient) THEN
1371 DO i=istru-1,iend
1372 IF (lbc_apply(ng)%south(i)) THEN
1373 zeta_new(i,jstr-1)=zeta_new(i,jstr)
1374# ifdef MASKING
1375 zeta_new(i,jstr-1)=zeta_new(i,jstr-1)* &
1376 & grid(ng)%rmask(i,jstr-1)
1377# endif
1378 END IF
1379 END DO
1380!
1381! Southern edge, closed boundary condition.
1382!
1383 ELSE IF (lbc(isouth,isfsur,ng)%closed) THEN
1384 DO i=istru-1,iend
1385 IF (lbc_apply(ng)%south(i)) THEN
1386 zeta_new(i,jstr-1)=zeta_new(i,jstr)
1387# ifdef MASKING
1388 zeta_new(i,jstr-1)=zeta_new(i,jstr-1)* &
1389 & grid(ng)%rmask(i,jstr-1)
1390# endif
1391 END IF
1392 END DO
1393 END IF
1394 END IF
1395!
1396!-----------------------------------------------------------------------
1397! Lateral boundary conditions at the northern edge.
1398!-----------------------------------------------------------------------
1399!
1400 IF (domain(ng)%Northern_Edge(tile)) THEN
1401!
1402! Northern edge, implicit upstream radiation condition.
1403!
1404 IF (lbc(inorth,isfsur,ng)%radiation) THEN
1405 DO i=istru-1,iend+1
1406 grad(i,jend )=zeta(i ,jend ,kstp)- &
1407 & zeta(i-1,jend ,kstp)
1408# ifdef MASKING
1409 grad(i,jend )=grad(i,jend )*grid(ng)%umask(i,jend )
1410# endif
1411 grad(i,jend+1)=zeta(i ,jend+1,kstp)- &
1412 & zeta(i-1,jend+1,kstp)
1413# ifdef MASKING
1414 grad(i,jend+1)=grad(i,jend+1)*grid(ng)%umask(i,jend+1)
1415# endif
1416 END DO
1417 DO i=istru-1,iend
1418 IF (lbc_apply(ng)%north(i)) THEN
1419 dzdt=zeta(i,jend,kstp)-zeta_new(i,jend )
1420 dzde=zeta_new(i,jend) -zeta_new(i,jend-1)
1421
1422 IF (lbc(inorth,isfsur,ng)%nudging) THEN
1423 IF ((dzdt*dzde).lt.0.0_r8) THEN
1424 tau=fsobc_in(ng,inorth)
1425 ELSE
1426 tau=fsobc_out(ng,inorth)
1427 END IF
1428 tau=tau*dtfast(ng)
1429 END IF
1430
1431 IF ((dzdt*dzde).lt.0.0_r8) dzdt=0.0_r8
1432 IF ((dzdt*(grad(i,jend)+grad(i+1,jend))).gt.0.0_r8) THEN
1433 dzdx=grad(i ,jend)
1434 ELSE
1435 dzdx=grad(i+1,jend)
1436 END IF
1437 cff=max(dzdx*dzdx+dzde*dzde,eps)
1438# ifdef RADIATION_2D
1439 cx=min(cff,max(dzdt*dzdx,-cff))
1440# else
1441 cx=0.0_r8
1442# endif
1443 ce=dzdt*dzde
1444# if defined CELERITY_WRITE && defined FORWARD_WRITE
1445 boundary(ng)%zeta_north_Cx(i)=cx
1446 boundary(ng)%zeta_north_Ce(i)=ce
1447 boundary(ng)%zeta_north_C2(i)=cff
1448# endif
1449 zeta_new(i,jend+1)=(cff*zeta(i,jend+1,kstp)+ &
1450 & ce *zeta_new(i,jend)- &
1451 & max(cx,0.0_r8)*grad(i ,jend+1)- &
1452 & min(cx,0.0_r8)*grad(i+1,jend+1))/ &
1453 & (cff+ce)
1454
1455 IF (lbc(inorth,isfsur,ng)%nudging) THEN
1456 zeta_new(i,jend+1)=zeta_new(i,jend+1)+ &
1457 & tau*(boundary(ng)%zeta_north(i)- &
1458 & zeta(i,jend+1,kstp))
1459 END IF
1460# ifdef MASKING
1461 zeta_new(i,jend+1)=zeta_new(i,jend+1)* &
1462 & grid(ng)%rmask(i,jend+1)
1463# endif
1464 END IF
1465 END DO
1466!
1467! Northern edge, explicit Chapman boundary condition.
1468!
1469 ELSE IF (lbc(inorth,isfsur,ng)%Chapman_explicit) THEN
1470 DO i=istru-1,iend
1471 IF (lbc_apply(ng)%north(i)) THEN
1472 cff=dtfast(ng)*grid(ng)%pn(i,jend)
1473 cff1=sqrt(g*(grid(ng)%h(i,jend)+ &
1474 & zeta(i,jend,kstp)))
1475 ce=cff*cff1
1476 zeta_new(i,jend+1)=(1.0_r8-ce)*zeta(i,jend+1,kstp)+ &
1477 & ce*zeta(i,jend,kstp)
1478# ifdef MASKING
1479 zeta_new(i,jend+1)=zeta_new(i,jend+1)* &
1480 & grid(ng)%rmask(i,jend+1)
1481# endif
1482 END IF
1483 END DO
1484!
1485! Northern edge, implicit Chapman boundary condition.
1486!
1487 ELSE IF (lbc(inorth,isfsur,ng)%Chapman_implicit) THEN
1488 DO i=istru-1,iend
1489 IF (lbc_apply(ng)%north(i)) THEN
1490 cff=dtfast(ng)*grid(ng)%pn(i,jend)
1491# ifdef WET_DRY
1492 cff1=sqrt(g*(max(grid(ng)%h(i,jend)+ &
1493 & zeta(i,jend,kstp),dcrit(ng))))
1494# else
1495 cff1=sqrt(g*(grid(ng)%h(i,jend)+ &
1496 & zeta(i,jend,kstp)))
1497# endif
1498 ce=cff*cff1
1499 cff2=1.0_r8/(1.0_r8+ce)
1500 zeta_new(i,jend+1)=cff2*(zeta(i,jend+1,kstp)+ &
1501 & ce*zeta_new(i,jend))
1502# ifdef MASKING
1503 zeta_new(i,jend+1)=zeta_new(i,jend+1)* &
1504 & grid(ng)%rmask(i,jend+1)
1505# endif
1506 END IF
1507 END DO
1508!
1509! Northern edge, clamped boundary condition.
1510!
1511 ELSE IF (lbc(inorth,isfsur,ng)%clamped) THEN
1512 DO i=istru-1,iend
1513 IF (lbc_apply(ng)%north(i)) THEN
1514 zeta_new(i,jend+1)=boundary(ng)%zeta_north(i)
1515# ifdef MASKING
1516 zeta_new(i,jend+1)=zeta_new(i,jend+1)* &
1517 & grid(ng)%rmask(i,jend+1)
1518# endif
1519 END IF
1520 END DO
1521!
1522! Northern edge, gradient boundary condition.
1523!
1524 ELSE IF (lbc(inorth,isfsur,ng)%gradient) THEN
1525 DO i=istru-1,iend
1526 IF (lbc_apply(ng)%north(i)) THEN
1527 zeta_new(i,jend+1)=zeta_new(i,jend)
1528# ifdef MASKING
1529 zeta_new(i,jend+1)=zeta_new(i,jend+1)* &
1530 & grid(ng)%rmask(i,jend+1)
1531# endif
1532 END IF
1533 END DO
1534!
1535! Northern edge, closed boundary condition.
1536!
1537 ELSE IF (lbc(inorth,isfsur,ng)%closed) THEN
1538 DO i=istru-1,iend
1539 IF (lbc_apply(ng)%north(i)) THEN
1540 zeta_new(i,jend+1)=zeta_new(i,jend)
1541# ifdef MASKING
1542 zeta_new(i,jend+1)=zeta_new(i,jend+1)* &
1543 & grid(ng)%rmask(i,jend+1)
1544# endif
1545 END IF
1546 END DO
1547 END IF
1548 END IF
1549!
1550!-----------------------------------------------------------------------
1551! Boundary corners.
1552!-----------------------------------------------------------------------
1553!
1554 IF (.not.(ewperiodic(ng).or.nsperiodic(ng))) THEN
1555 IF (domain(ng)%SouthWest_Corner(tile)) THEN
1556 IF (lbc_apply(ng)%south(istr-1).and. &
1557 & lbc_apply(ng)%west (jstr-1)) THEN
1558 zeta_new(istr-1,jstr-1)=0.5_r8*(zeta_new(istr ,jstr-1)+ &
1559 & zeta_new(istr-1,jstr ))
1560 END IF
1561 END IF
1562 IF (domain(ng)%SouthEast_Corner(tile)) THEN
1563 IF (lbc_apply(ng)%south(iend+1).and. &
1564 & lbc_apply(ng)%east (jstr-1)) THEN
1565 zeta_new(iend+1,jstr-1)=0.5_r8*(zeta_new(iend ,jstr-1)+ &
1566 & zeta_new(iend+1,jstr ))
1567 END IF
1568 END IF
1569 IF (domain(ng)%NorthWest_Corner(tile)) THEN
1570 IF (lbc_apply(ng)%north(istr-1).and. &
1571 & lbc_apply(ng)%west (jend+1)) THEN
1572 zeta_new(istr-1,jend+1)=0.5_r8*(zeta_new(istr-1,jend )+ &
1573 & zeta_new(istr ,jend+1))
1574 END IF
1575 END IF
1576 IF (domain(ng)%NorthEast_Corner(tile)) THEN
1577 IF (lbc_apply(ng)%north(iend+1).and. &
1578 & lbc_apply(ng)%east (jend+1)) THEN
1579 zeta_new(iend+1,jend+1)=0.5_r8*(zeta_new(iend+1,jend )+ &
1580 & zeta_new(iend ,jend+1))
1581 END IF
1582 END IF
1583 END IF
1584
1585# if defined WET_DRY
1586!
1587!-----------------------------------------------------------------------
1588! Ensure that water level on boundary cells is above bed elevation.
1589!-----------------------------------------------------------------------
1590!
1591 cff=dcrit(ng)-eps
1592 IF (.not.ewperiodic(ng)) THEN
1593 IF (domain(ng)%Western_Edge(tile)) THEN
1594 DO j=jstrv-1,jend
1595 IF (lbc_apply(ng)%west(j)) THEN
1596 IF (zeta_new(istr-1,j).le. &
1597 & (dcrit(ng)-grid(ng)%h(istr-1,j))) THEN
1598 zeta_new(istr-1,j)=cff-grid(ng)%h(istr-1,j)
1599 END IF
1600 END IF
1601 END DO
1602 END IF
1603 IF (domain(ng)%Eastern_Edge(tile)) THEN
1604 DO j=jstrv-1,jend
1605 IF (lbc_apply(ng)%east(j)) THEN
1606 IF (zeta_new(iend+1,j).le. &
1607 & (dcrit(ng)-grid(ng)%h(iend+1,j))) THEN
1608 zeta_new(iend+1,j)=cff-grid(ng)%h(iend+1,j)
1609 END IF
1610 END IF
1611 END DO
1612 END IF
1613 END IF
1614!
1615 IF (.not.nsperiodic(ng)) THEN
1616 IF (domain(ng)%Southern_Edge(tile)) THEN
1617 DO i=istru-1,iend
1618 IF (lbc_apply(ng)%south(i)) THEN
1619 IF (zeta_new(i,jstr-1).le. &
1620 & (dcrit(ng)-grid(ng)%h(i,jstr-1))) THEN
1621 zeta_new(i,jstr-1)=cff-grid(ng)%h(i,jstr-1)
1622 END IF
1623 END IF
1624 END DO
1625 END IF
1626 IF (domain(ng)%Northern_Edge(tile)) THEN
1627 DO i=istru-1,iend
1628 IF (lbc_apply(ng)%north(i)) THEN
1629 IF (zeta_new(i,jend+1).le. &
1630 & (dcrit(ng)-grid(ng)%h(i,jend+1))) THEN
1631 zeta_new(i,jend+1)=cff-grid(ng)%h(i,jend+1)
1632 END IF
1633 END IF
1634 END DO
1635 END IF
1636 END IF
1637!
1638 IF (.not.(ewperiodic(ng).or.nsperiodic(ng))) THEN
1639 IF (domain(ng)%SouthWest_Corner(tile)) THEN
1640 IF (lbc_apply(ng)%south(istr-1).and. &
1641 & lbc_apply(ng)%west (jstr-1)) THEN
1642 IF (zeta_new(istr-1,jstr-1).le. &
1643 & (dcrit(ng)-grid(ng)%h(istr-1,jstr-1))) THEN
1644 zeta_new(istr-1,jstr-1)=cff-grid(ng)%h(istr-1,jstr-1)
1645 END IF
1646 END IF
1647 END IF
1648 IF (domain(ng)%SouthEast_Corner(tile)) THEN
1649 IF (lbc_apply(ng)%south(iend+1).and. &
1650 & lbc_apply(ng)%east (jstr-1)) THEN
1651 IF (zeta_new(iend+1,jstr-1).le. &
1652 & (dcrit(ng)-grid(ng)%h(iend+1,jstr-1))) THEN
1653 zeta_new(iend+1,jstr-1)=cff-grid(ng)%h(iend+1,jstr-1)
1654 END IF
1655 END IF
1656 END IF
1657 IF (domain(ng)%NorthWest_Corner(tile)) THEN
1658 IF (lbc_apply(ng)%north(istr-1).and. &
1659 & lbc_apply(ng)%west (jend+1)) THEN
1660 IF (zeta_new(istr-1,jend+1).le. &
1661 & (dcrit(ng)-grid(ng)%h(istr-1,jend+1))) THEN
1662 zeta_new(istr-1,jend+1)=cff-grid(ng)%h(istr-1,jend+1)
1663 END IF
1664 END IF
1665 END IF
1666 IF (domain(ng)%NorthEast_Corner(tile)) THEN
1667 IF (lbc_apply(ng)%north(iend+1).and. &
1668 & lbc_apply(ng)%east (jend+1)) THEN
1669 IF (zeta_new(iend+1,jend+1).le. &
1670 & (dcrit(ng)-grid(ng)%h(iend+1,jend+1))) THEN
1671 zeta_new(iend+1,jend+1)=cff-grid(ng)%h(iend+1,jend+1)
1672 END IF
1673 END IF
1674 END IF
1675 END IF
1676# endif
1677!
1678 RETURN
1679!
1680 END SUBROUTINE zetabc_local
1681#endif
1682!
1683 END MODULE zetabc_mod
type(t_boundary), dimension(:), allocatable boundary
type(t_apply), dimension(:), allocatable lbc_apply
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
integer isfsur
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
type(t_lbc), dimension(:,:,:), allocatable lbc
Definition mod_param.F:375
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
real(r8), dimension(:), allocatable dcrit
logical, dimension(:), allocatable ewperiodic
integer, parameter iwest
logical, dimension(:), allocatable nsperiodic
logical, dimension(:), allocatable predictor_2d_step
real(dp), dimension(:,:), allocatable fsobc_out
integer, parameter isouth
real(dp), dimension(:), allocatable dtfast
integer, parameter ieast
real(dp) g
integer, parameter inorth
real(dp), dimension(:,:), allocatable fsobc_in
integer, dimension(:), allocatable kstp
integer, dimension(:), allocatable krhs
subroutine zetabc(ng, tile, kout)
Definition zetabc.F:40
subroutine, public zetabc_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, krhs, kstp, kout, zeta)
Definition zetabc.F:65