64
65
73
74
75
76 integer, intent(in) :: ng, tile
77 integer, intent(in) :: LBi, UBi, LBj, UBj
78 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
79 integer, intent(in) :: krhs, kstp, kout
80
81# ifdef ASSUMED_SHAPE
82 real(r8), intent(in) :: ubar(LBi:,LBj:,:)
83 real(r8), intent(in) :: vbar(LBi:,LBj:,:)
84 real(r8), intent(in) :: zeta(LBi:,LBj:,:)
85 real(r8), intent(in) :: tl_vbar(LBi:,LBj:,:)
86 real(r8), intent(in) :: tl_zeta(LBi:,LBj:,:)
87
88 real(r8), intent(inout) :: tl_ubar(LBi:,LBj:,:)
89# else
90 real(r8), intent(in) :: ubar(LBi:UBi,LBj:UBj,:)
91 real(r8), intent(in) :: vbar(LBi:UBi,LBj:UBj,:)
92 real(r8), intent(in) :: zeta(LBi:UBi,LBj:UBj,:)
93 real(r8), intent(in) :: tl_vbar(LBi:UBi,LBj:UBj,:)
94 real(r8), intent(in) :: tl_zeta(LBi:UBi,LBj:UBj,:)
95
96 real(r8), intent(inout) :: tl_ubar(LBi:UBi,LBj:UBj,:)
97# endif
98
99
100
101 integer :: Imin, Imax
102 integer :: i, j, know
103
104 real(r8) :: Ce, Cx, Zx
105 real(r8) :: bry_pgr, bry_cor, bry_str
106 real(r8) :: cff, cff1, cff2, cff3, dt2d
107 real(r8) :: obc_in, obc_out, tau
108# if defined ATM_PRESS && defined PRESS_COMPENSATE
109 real(r8) :: OneAtm, fac
110# endif
111
112 real(r8) :: tl_Ce, tl_Cx, tl_Zx
113 real(r8) :: tl_bry_pgr, tl_bry_cor, tl_bry_str, tl_bry_val
114 real(r8) :: tl_cff, tl_cff1, tl_cff2, tl_cff3
115
116 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_grad
117
118# include "set_bounds.h"
119
120
121
122
123
124 IF (first_2d_step) THEN
125 know=krhs
128 know=krhs
130 ELSE
131 know=kstp
133 END IF
134# if defined ATM_PRESS && defined PRESS_COMPENSATE
135 oneatm=1013.25_r8
136 fac=100.0_r8/(
g*
rho0)
137# endif
138
139
140
141
142
143 IF (
domain(ng)%Western_Edge(tile))
THEN
144
145
146
148 IF (
iic(ng).ne.0)
THEN
149 DO j=jstr,jend+1
150
151
152
153 tl_grad(istr,j)=0.0_r8
154 END DO
155 DO j=jstr,jend
157# if defined CELERITY_READ && defined FORWARD_READ
160 obc_out=0.5_r8* &
161 & (
clima(ng)%M2nudgcof(istr-1,j)+ &
162 &
clima(ng)%M2nudgcof(istr ,j))
163 obc_in =
obcfac(ng)*obc_out
164 ELSE
167 END IF
168 IF (
boundary(ng)%ubar_west_Cx(j).lt.0.0_r8)
THEN
169 tau=obc_in
170 ELSE
171 tau=obc_out
172 END IF
173 tau=tau*dt2d
174 END IF
176# ifdef RADIATION_2D
178# else
179 ce=0.0_r8
180# endif
182# endif
183
184
185
186
187
188
189 tl_ubar(istr,j,kout)=(cff*tl_ubar(istr ,j,know)+ &
190 & cx *tl_ubar(istr+1,j,kout)- &
191 & max(ce,0.0_r8)* &
192 & tl_grad(istr,j )- &
193 & min(ce,0.0_r8)* &
194 & tl_grad(istr,j+1))/ &
195 & (cff+cx)
196
198
199
200
201
202 tl_ubar(istr,j,kout)=tl_ubar(istr,j,kout)- &
203 & tau*tl_ubar(istr,j,know)
204 END IF
205# ifdef MASKING
206
207
208
209 tl_ubar(istr,j,kout)=tl_ubar(istr,j,kout)* &
210 &
grid(ng)%umask(istr,j)
211# endif
212 END IF
213 END DO
214 END IF
215
216
217
219 DO j=jstr,jend
221# if defined SSH_TIDES_NOT_YET && !defined UV_TIDES_NOT_YET
223 bry_pgr=-
g*(zeta(istr,j,know)- &
225 & 0.5_r8*
grid(ng)%pm(istr,j)
226 tl_bry_pgr=-
g*tl_zeta(istr,j,know)* &
227 & 0.5_r8*
grid(ng)%pm(istr,j)
228# ifdef ADJUST_BOUNDARY
230 tl_bry_pgr=tl_bry_pgr+ &
232 & 0.5_r8*
grid(ng)%pm(istr,j)
233 END IF
234# endif
235 ELSE
236 bry_pgr=-
g*(zeta(istr ,j,know)- &
237 & zeta(istr-1,j,know))* &
238 & 0.5_r8*(
grid(ng)%pm(istr-1,j)+ &
239 &
grid(ng)%pm(istr ,j))
240 tl_bry_pgr=-
g*(tl_zeta(istr ,j,know)- &
241 & tl_zeta(istr-1,j,know))* &
242 & 0.5_r8*(
grid(ng)%pm(istr-1,j)+ &
243 &
grid(ng)%pm(istr ,j))
244 END IF
245# ifdef UV_COR
246 bry_cor=0.125_r8*(vbar(istr-1,j ,know)+ &
247 & vbar(istr-1,j+1,know)+ &
248 & vbar(istr ,j ,know)+ &
249 & vbar(istr ,j+1,know))* &
250 & (
grid(ng)%f(istr-1,j)+ &
251 &
grid(ng)%f(istr ,j))
252 tl_bry_cor=0.125_r8*(tl_vbar(istr-1,j ,know)+ &
253 & tl_vbar(istr-1,j+1,know)+ &
254 & tl_vbar(istr ,j ,know)+ &
255 & tl_vbar(istr ,j+1,know))* &
256 & (
grid(ng)%f(istr-1,j)+ &
257 &
grid(ng)%f(istr ,j))
258# else
259 bry_cor=0.0_r8
260 tl_bry_cor=0.0_r8
261# endif
262 cff1=1.0_r8/(0.5_r8*(
grid(ng)%h(istr-1,j)+ &
263 & zeta(istr-1,j,know)+ &
264 &
grid(ng)%h(istr ,j)+ &
265 & zeta(istr ,j,know)))
266 tl_cff1=-cff1*cff1*0.5_r8*(
grid(ng)%tl_h(istr-1,j)+ &
267 & tl_zeta(istr-1,j,know)+ &
268 &
grid(ng)%tl_h(istr ,j)+ &
269 & tl_zeta(istr ,j,know))
270 bry_str=cff1*(
forces(ng)%sustr(istr,j)- &
271 &
forces(ng)%bustr(istr,j))
272 tl_bry_str=tl_cff1*(
forces(ng)%sustr(istr,j)- &
273 &
forces(ng)%bustr(istr,j))+ &
274 & cff1*(
forces(ng)%tl_sustr(istr,j)- &
275 &
forces(ng)%tl_bustr(istr,j))
276 cx=1.0_r8/sqrt(
g*0.5_r8*(
grid(ng)%h(istr-1,j)+ &
277 & zeta(istr-1,j,know)+ &
278 &
grid(ng)%h(istr ,j)+ &
279 & zeta(istr ,j,know)))
280 tl_cx=-cx*cx*cx*0.25_r8*
g*(
grid(ng)%tl_h(istr-1,j)+ &
281 & tl_zeta(istr-1,j,know)+ &
282 &
grid(ng)%tl_h(istr ,j)+ &
283 & tl_zeta(istr ,j,know))
284 cff2=
grid(ng)%om_u(istr,j)*cx
285 tl_cff2=
grid(ng)%om_u(istr,j)*tl_cx
286
287
288
289
290
291 tl_bry_val=tl_ubar(istr+1,j,know)+ &
292 & tl_cff2*(bry_pgr+ &
293 & bry_cor+ &
294 & bry_str)+ &
295 & cff2*(tl_bry_pgr+ &
296 & tl_bry_cor+ &
297 & tl_bry_str)
298# else
299
300
301# ifdef ADJUST_BOUNDARY
303 tl_bry_val=
boundary(ng)%tl_ubar_west(j)
304 ELSE
305 tl_bry_val=0.0_r8
306 END IF
307# else
308 tl_bry_val=0.0_r8
309# endif
310# endif
311 cff=1.0_r8/(0.5_r8*(
grid(ng)%h(istr-1,j)+ &
312 & zeta(istr-1,j,know)+ &
313 &
grid(ng)%h(istr ,j)+ &
314 & zeta(istr ,j,know)))
315 tl_cff=-cff*cff*(0.5_r8*(
grid(ng)%tl_h(istr-1,j)+ &
316 & tl_zeta(istr-1,j,know)+ &
317 &
grid(ng)%tl_h(istr ,j)+ &
318 & tl_zeta(istr ,j,know)))
320 tl_cx=0.5_r8*
g*tl_cff/cx
321# if defined ATM_PRESS && defined PRESS_COMPENSATE
322
323
324
325
326
327
328
329
330
331 tl_ubar(istr,j,kout)=tl_bry_val- &
332 & tl_cx* &
333 & (0.5_r8* &
334 & (zeta(istr-1,j,know)+ &
335 & zeta(istr ,j,know)+ &
336 & fac*(
forces(ng)%Pair(istr-1,j)+ &
337 &
forces(ng)%Pair(istr ,j)- &
338 & 2.0_r8*oneatm))- &
340 & cx* &
341 & (0.5_r8* &
342 & (tl_zeta(istr-1,j,know)+ &
343 & tl_zeta(istr ,j,know)))
344# else
345
346
347
348
349
350 tl_ubar(istr,j,kout)=tl_bry_val- &
351 & tl_cx* &
352 & (0.5_r8*(zeta(istr-1,j,know)+ &
353 & zeta(istr ,j,know))- &
355 & cx* &
356 & (0.5_r8*(tl_zeta(istr-1,j,know)+ &
357 & tl_zeta(istr ,j,know)))
358# endif
359# ifdef ADJUST_BOUNDARY
361 tl_ubar(istr,j,kout)=tl_ubar(istr,j,kout)+ &
363 END IF
364# endif
365# ifdef MASKING
366
367
368
369 tl_ubar(istr,j,kout)=tl_ubar(istr,j,kout)* &
370 &
grid(ng)%umask(istr,j)
371# endif
372 END IF
373 END DO
374
375
376
378 DO j=jstr,jend
380# if defined SSH_TIDES_NOT_YET && !defined UV_TIDES_NOT_YET
382 bry_pgr=-
g*(zeta(istr,j,know)- &
384 & 0.5_r8*
grid(ng)%pm(istr,j)
385 tl_bry_pgr=-
g*tl_zeta(istr,j,know)* &
386 & 0.5_r8*
grid(ng)%pm(istr,j)
387# ifdef ADJUST_BOUNDARY
389 tl_bry_pgr=tl_bry_pgr+ &
391 & 0.5_r8*
grid(ng)%pm(istr,j)
392 END IF
393# endif
394 ELSE
395 bry_pgr=-
g*(zeta(istr ,j,know)- &
396 & zeta(istr-1,j,know))* &
397 & 0.5_r8*(
grid(ng)%pm(istr-1,j)+ &
398 &
grid(ng)%pm(istr ,j))
399 tl_bry_pgr=-
g*(tl_zeta(istr ,j,know)- &
400 & tl_zeta(istr-1,j,know))* &
401 & 0.5_r8*(
grid(ng)%pm(istr-1,j)+ &
402 &
grid(ng)%pm(istr ,j))
403 END IF
404# ifdef UV_COR
405 bry_cor=0.125_r8*(vbar(istr-1,j ,know)+ &
406 & vbar(istr-1,j+1,know)+ &
407 & vbar(istr ,j ,know)+ &
408 & vbar(istr ,j+1,know))* &
409 & (
grid(ng)%f(istr-1,j)+ &
410 &
grid(ng)%f(istr ,j))
411 tl_bry_cor=0.125_r8*(tl_vbar(istr-1,j ,know)+ &
412 & tl_vbar(istr-1,j+1,know)+ &
413 & tl_vbar(istr ,j ,know)+ &
414 & tl_vbar(istr ,j+1,know))* &
415 & (
grid(ng)%f(istr-1,j)+ &
416 &
grid(ng)%f(istr ,j))
417# else
418 bry_cor=0.0_r8
419 tl_bry_cor=0.0_r8
420# endif
421 cff1=1.0_r8/(0.5_r8*(
grid(ng)%h(istr-1,j)+ &
422 & zeta(istr-1,j,know)+ &
423 &
grid(ng)%h(istr ,j)+ &
424 & zeta(istr ,j,know)))
425 tl_cff1=-cff1*cff1*0.5_r8*(
grid(ng)%tl_h(istr-1,j)+ &
426 & tl_zeta(istr-1,j,know)+ &
427 &
grid(ng)%tl_h(istr ,j)+ &
428 & tl_zeta(istr ,j,know))
429 bry_str=cff1*(
forces(ng)%sustr(istr,j)- &
430 &
forces(ng)%bustr(istr,j))
431 tl_bry_str=tl_cff1*(
forces(ng)%sustr(istr,j)- &
432 &
forces(ng)%bustr(istr,j))+ &
433 & cff1*(
forces(ng)%tl_sustr(istr,j)- &
434 &
forces(ng)%tl_bustr(istr,j))
435 cx=1.0_r8/sqrt(
g*0.5_r8*(
grid(ng)%h(istr-1,j)+ &
436 & zeta(istr-1,j,know)+ &
437 &
grid(ng)%h(istr ,j)+ &
438 & zeta(istr ,j,know)))
439 tl_cx=-cx*cx*cx*0.25_r8*
g*(
grid(ng)%tl_h(istr-1,j)+ &
440 & tl_zeta(istr-1,j,know)+ &
441 &
grid(ng)%tl_h(istr ,j)+ &
442 & tl_zeta(istr ,j,know))
443 cff2=
grid(ng)%om_u(istr,j)*cx
444 tl_cff2=
grid(ng)%om_u(istr,j)*tl_cx
445
446
447
448
449
450 tl_bry_val=tl_ubar(istr+1,j,know)+ &
451 & tl_cff2*(bry_pgr+ &
452 & bry_cor+ &
453 & bry_str)+ &
454 & cff2*(tl_bry_pgr+ &
455 & tl_bry_cor+ &
456 & tl_bry_str)
457# else
458
459
460# ifdef ADJUST_BOUNDARY
462 tl_bry_val=
boundary(ng)%tl_ubar_west(j)
463 ELSE
464 tl_bry_val=0.0_r8
465 END IF
466# else
467 tl_bry_val=0.0_r8
468# endif
469# endif
470# ifdef WET_DRY_NOT_YET
471 cff=0.5_r8*(
grid(ng)%h(istr-1,j)+ &
472 & zeta(istr-1,j,know)+ &
473 &
grid(ng)%h(istr ,j)+ &
474 & zeta(istr ,j,know))
475 tl_cff=0.5_r8*(
grid(ng)%tl_h(istr-1,j)+ &
476 & tl_zeta(istr-1,j,know)+ &
477 &
grid(ng)%tl_h(istr ,j)+ &
478 & tl_zeta(istr ,j,know))
479# else
480 cff=0.5_r8*(
grid(ng)%h(istr-1,j)+ &
481 &
grid(ng)%h(istr ,j))
482 tl_cff=0.5_r8*(
grid(ng)%tl_h(istr-1,j)+ &
483 &
grid(ng)%tl_h(istr ,j))
484# endif
486 tl_cff1=-0.5_r8*cff1*tl_cff/cff
487 cx=dt2d*cff1*cff*0.5_r8*(
grid(ng)%pm(istr-1,j)+ &
488 &
grid(ng)%pm(istr ,j))
489 tl_cx=dt2d*0.5_r8*(
grid(ng)%pm(istr-1,j)+ &
490 &
grid(ng)%pm(istr ,j))* &
491 & (cff1*tl_cff+ &
492 & tl_cff1*cff)
493 zx=(0.5_r8+cx)*zeta(istr ,j,know)+ &
494 & (0.5_r8-cx)*zeta(istr-1,j,know)
495 tl_zx=(0.5_r8+cx)*tl_zeta(istr ,j,know)+ &
496 & (0.5_r8-cx)*tl_zeta(istr-1,j,know)+ &
497 & tl_cx*(zeta(istr ,j,know)- &
498 & zeta(istr-1,j,know))
500 cff2=(1.0_r8-
co/cx)**2
501 tl_cff2=2.0_r8*cff2*
co*tl_cx/(cx*cx)
502 cff3=zeta(istr,j,kout)+ &
503 & cx*zeta(istr-1,j,know)- &
504 & (1.0_r8+cx)*zeta(istr,j,know)
505 tl_cff3=tl_zeta(istr,j,kout)+ &
506 & cx*tl_zeta(istr-1,j,know)+ &
507 & tl_cx*(zeta(istr-1,j,know)+ &
508 & zeta(istr ,j,know))- &
509 & (1.0_r8+cx)*tl_zeta(istr,j,know)
510 zx=zx+cff2*cff3
511 tl_zx=tl_zx+cff2*tl_cff3+ &
512 & tl_cff2*cff3
513 END IF
514
515
516
517
518
519
520 tl_ubar(istr,j,kout)=0.5_r8* &
521 & ((1.0_r8-cx)* &
522 & tl_ubar(istr,j,know)- &
523 & tl_cx*(ubar(istr ,j,know)- &
524 & ubar(istr+1,j,know))+ &
525 & cx*tl_ubar(istr+1,j,know)+ &
526 & tl_bry_val- &
527 & tl_cff1* &
529 & cff1*tl_zx)
530# ifdef ADJUST_BOUNDARY
532 tl_ubar(istr,j,kout)=tl_ubar(istr,j,kout)+ &
533 & 0.5_r8*cff1* &
535 END IF
536# endif
537# ifdef MASKING
538
539
540
541 tl_ubar(istr,j,kout)=tl_ubar(istr,j,kout)* &
542 &
grid(ng)%umask(istr,j)
543# endif
544 END IF
545 END DO
546
547
548
550 DO j=jstr,jend
552
553
554# ifdef ADJUST_BOUNDARY
556 tl_ubar(istr,j,kout)=
boundary(ng)%tl_ubar_west(j)
557 ELSE
558 tl_ubar(istr,j,kout)=0.0_r8
559 END IF
560# else
561 tl_ubar(istr,j,kout)=0.0_r8
562# endif
563# ifdef MASKING
564
565
566
567 tl_ubar(istr,j,kout)=tl_ubar(istr,j,kout)* &
568 &
grid(ng)%umask(istr,j)
569# endif
570 END IF
571 END DO
572
573
574
576 DO j=jstr,jend
578
579
580 tl_ubar(istr,j,kout)=tl_ubar(istr+1,j,kout)
581# ifdef MASKING
582
583
584
585 tl_ubar(istr,j,kout)=tl_ubar(istr,j,kout)* &
586 &
grid(ng)%umask(istr,j)
587# endif
588 END IF
589 END DO
590
591
592
594 DO j=jstr,jend
597
598
599
600
601 tl_bry_pgr=-
g*tl_zeta(istr,j,know)* &
602 & 0.5_r8*
grid(ng)%pm(istr,j)
603# ifdef ADJUST_BOUNDARY
605 tl_bry_pgr=tl_bry_pgr+ &
607 & 0.5_r8*
grid(ng)%pm(istr,j)
608 END IF
609# endif
610 ELSE
611
612
613
614
615
616 tl_bry_pgr=-
g*(tl_zeta(istr ,j,know)- &
617 & tl_zeta(istr-1,j,know))* &
618 & 0.5_r8*(
grid(ng)%pm(istr-1,j)+ &
619 &
grid(ng)%pm(istr ,j))
620 END IF
621# ifdef UV_COR
622
623
624
625
626
627
628
629 tl_bry_cor=0.125_r8*(tl_vbar(istr-1,j ,know)+ &
630 & tl_vbar(istr-1,j+1,know)+ &
631 & tl_vbar(istr ,j ,know)+ &
632 & tl_vbar(istr ,j+1,know))* &
633 & (
grid(ng)%f(istr-1,j)+ &
634 &
grid(ng)%f(istr ,j))
635# else
636
637
638 tl_bry_cor=0.0_r8
639# endif
640 cff=1.0_r8/(0.5_r8*(
grid(ng)%h(istr-1,j)+ &
641 & zeta(istr-1,j,know)+ &
642 &
grid(ng)%h(istr ,j)+ &
643 & zeta(istr ,j,know)))
644 tl_cff=-cff*cff*0.5_r8*(
grid(ng)%tl_h(istr-1,j)+ &
645 & tl_zeta(istr-1,j,know)+ &
646 &
grid(ng)%tl_h(istr ,j)+ &
647 & tl_zeta(istr ,j,know))
648
649
650
651 tl_bry_str=tl_cff*(
forces(ng)%sustr(istr,j)- &
652 &
forces(ng)%bustr(istr,j))+ &
653 & cff*(
forces(ng)%tl_sustr(istr,j)- &
654 &
forces(ng)%tl_bustr(istr,j))
655
656
657
658
659
660 tl_ubar(istr,j,kout)=tl_ubar(istr,j,know)+ &
661 & dt2d*(tl_bry_pgr+ &
662 & tl_bry_cor+ &
663 & tl_bry_str)
664# ifdef MASKING
665
666
667
668 tl_ubar(istr,j,kout)=tl_ubar(istr,j,kout)* &
669 &
grid(ng)%umask(istr,j)
670# endif
671 END IF
672 END DO
673
674
675
677 DO j=jstr,jend
679
680
681 tl_ubar(istr,j,kout)=0.0_r8
682 END IF
683 END DO
684 END IF
685 END IF
686
687
688
689
690
691 IF (
domain(ng)%Eastern_Edge(tile))
THEN
692
693
694
696 IF (
iic(ng).ne.0)
THEN
697 DO j=jstr,jend+1
698
699
700
701 tl_grad(iend+1,j)=0.0_r8
702 END DO
703 DO j=jstr,jend
705# if defined CELERITY_READ && defined FORWARD_READ
708 obc_out=0.5_r8* &
709 & (
clima(ng)%M2nudgcof(iend ,j)+ &
710 &
clima(ng)%M2nudgcof(iend+1,j))
711 obc_in =
obcfac(ng)*obc_out
712 ELSE
715 END IF
716 IF (
boundary(ng)%ubar_east_Cx(j).lt.0.0_r8)
THEN
717 tau=obc_in
718 ELSE
719 tau=obc_out
720 END IF
721 tau=tau*dt2d
722 END IF
724# ifdef RADIATION_2D
726# else
727 ce=0.0_r8
728# endif
730# endif
731
732
733
734
735
736
737 tl_ubar(iend+1,j,kout)=(cff*tl_ubar(iend+1,j,know)+ &
738 & cx *tl_ubar(iend ,j,kout)- &
739 & max(ce,0.0_r8)* &
740 & tl_grad(iend+1,j )- &
741 & min(ce,0.0_r8)* &
742 & tl_grad(iend+1,j+1))/ &
743 & (cff+cx)
744
746
747
748
749
750 tl_ubar(iend+1,j,kout)=tl_ubar(iend+1,j,kout)- &
751 & tau*tl_ubar(iend+1,j,know)
752 END IF
753# ifdef MASKING
754
755
756
757 tl_ubar(iend+1,j,kout)=tl_ubar(iend+1,j,kout)* &
758 &
grid(ng)%umask(iend+1,j)
759# endif
760 END IF
761 END DO
762 END IF
763
764
765
767 DO j=jstr,jend
769# if defined SSH_TIDES_NOT_YET && !defined UV_TIDES_NOT_YET
772 & zeta(iend,j,know))* &
773 & 0.5_r8*
grid(ng)%pm(iend,j)
774 tl_bry_pgr=
g*tl_zeta(iend,j,know)* &
775 & 0.5_r8*
grid(ng)%pm(iend,j)
776# ifdef ADJUST_BOUNDARY
778 tl_bry_pgr=tl_bry_pgr- &
780 & 0.5_r8*
grid(ng)%pm(iend,j)
781 END IF
782# endif
783 ELSE
784 bry_pgr=-
g*(zeta(iend+1,j,know)- &
785 & zeta(iend ,j,know))* &
786 & 0.5_r8*(
grid(ng)%pm(iend ,j)+ &
787 &
grid(ng)%pm(iend+1,j))
788 tl_bry_pgr=-
g*(tl_zeta(iend+1,j,know)- &
789 & tl_zeta(iend ,j,know))* &
790 & 0.5_r8*(
grid(ng)%pm(iend ,j)+ &
791 &
grid(ng)%pm(iend+1,j))
792 END IF
793# ifdef UV_COR
794 bry_cor=0.125_r8*(vbar(iend ,j ,know)+ &
795 & vbar(iend ,j+1,know)+ &
796 & vbar(iend+1,j ,know)+ &
797 & vbar(iend+1,j+1,know))* &
798 & (
grid(ng)%f(iend ,j)+ &
799 &
grid(ng)%f(iend+1,j))
800 tl_bry_cor=0.125_r8*(tl_vbar(iend ,j ,know)+ &
801 & tl_vbar(iend ,j+1,know)+ &
802 & tl_vbar(iend+1,j ,know)+ &
803 & tl_vbar(iend+1,j+1,know))* &
804 & (
grid(ng)%f(iend ,j)+ &
805 &
grid(ng)%f(iend+1,j))
806# else
807 bry_cor=0.0_r8
808 tl_bry_cor=0.0_r8
809# endif
810 cff1=1.0_r8/(0.5_r8*(
grid(ng)%h(iend ,j)+ &
811 & zeta(iend ,j,know)+ &
812 &
grid(ng)%h(iend+1,j)+ &
813 & zeta(iend+1,j,know)))
814 tl_cff1=-cff1*cff1*0.5_r8*(
grid(ng)%tl_h(iend ,j)+ &
815 & tl_zeta(iend ,j,know)+ &
816 &
grid(ng)%tl_h(iend+1,j)+ &
817 & tl_zeta(iend+1,j,know))
818 bry_str=cff1*(
forces(ng)%sustr(iend+1,j)- &
819 &
forces(ng)%bustr(iend+1,j))
820 tl_bry_str=tl_cff1*(
forces(ng)%sustr(iend+1,j)- &
821 &
forces(ng)%bustr(iend+1,j))+ &
822 & cff1*(
forces(ng)%tl_sustr(iend+1,j)- &
823 &
forces(ng)%tl_bustr(iend+1,j))
824 cx=1.0_r8/sqrt(
g*0.5_r8*(
grid(ng)%h(iend+1,j)+ &
825 & zeta(iend+1,j,know)+ &
826 &
grid(ng)%h(iend ,j)+ &
827 & zeta(iend ,j,know)))
828 tl_cx=-cx*cx*cx*0.25_r8*
g*(
grid(ng)%tl_h(iend+1,j)+ &
829 & tl_zeta(iend+1,j,know)+ &
830 &
grid(ng)%tl_h(iend ,j)+ &
831 & tl_zeta(iend ,j,know))
832 cff2=
grid(ng)%om_u(iend+1,j)*cx
833 tl_cff2=
grid(ng)%om_u(iend+1,j)*tl_cx
834
835
836
837
838
839 tl_bry_val=tl_ubar(iend,j,know)+ &
840 & tl_cff2*(bry_pgr+ &
841 & bry_cor+ &
842 & bry_str)+ &
843 & cff2*(tl_bry_pgr+ &
844 & tl_bry_cor+ &
845 & tl_bry_str)
846# else
847
848
849# ifdef ADJUST_BOUNDARY
851 tl_bry_val=
boundary(ng)%tl_ubar_east(j)
852 ELSE
853 tl_bry_val=0.0_r8
854 END IF
855# else
856 tl_bry_val=0.0_r8
857# endif
858# endif
859 cff=1.0_r8/(0.5_r8*(
grid(ng)%h(iend ,j)+ &
860 & zeta(iend ,j,know)+ &
861 &
grid(ng)%h(iend+1,j)+ &
862 & zeta(iend+1,j,know)))
863 tl_cff=-cff*cff*(0.5_r8*(
grid(ng)%tl_h(iend ,j)+ &
864 & tl_zeta(iend ,j,know)+ &
865 &
grid(ng)%tl_h(iend+1,j)+ &
866 & tl_zeta(iend+1,j,know)))
868 tl_cx=0.5_r8*
g*tl_cff/cx
869# if defined ATM_PRESS && defined PRESS_COMPENSATE
870
871
872
873
874
875
876
877
878
879 tl_ubar(iend+1,j,kout)=tl_bry_val+ &
880 & tl_cx* &
881 & (0.5_r8* &
882 & (zeta(iend ,j,know)+ &
883 & zeta(iend+1,j,know)+ &
884 & fac*(
forces(ng)%Pair(iend ,j)+ &
885 &
forces(ng)%Pair(iend+1,j)- &
886 & 2.0_r8*oneatm))- &
888 & cx* &
889 & (0.5_r8*(tl_zeta(iend ,j,know)+ &
890 & tl_zeta(iend+1,j,know)))
891# else
892
893
894
895
896
897 tl_ubar(iend+1,j,kout)=tl_bry_val+ &
898 & tl_cx* &
899 & (0.5_r8*(zeta(iend ,j,know)+ &
900 & zeta(iend+1,j,know))- &
902 & cx* &
903 & (0.5_r8*(tl_zeta(iend ,j,know)+ &
904 & tl_zeta(iend+1,j,know)))
905# endif
906# ifdef ADJUST_BOUNDARY
908 tl_ubar(iend+1,j,kout)=tl_ubar(iend+1,j,kout)- &
910 END IF
911# endif
912# ifdef MASKING
913
914
915
916 tl_ubar(iend+1,j,kout)=tl_ubar(iend+1,j,kout)* &
917 &
grid(ng)%umask(iend+1,j)
918# endif
919 END IF
920 END DO
921
922
923
925 DO j=jstr,jend
927# if defined SSH_TIDES_NOT_YET && !defined UV_TIDES_NOT_YET
930 & zeta(iend,j,know))* &
931 & 0.5_r8*
grid(ng)%pm(iend,j)
932 tl_bry_pgr=
g*tl_zeta(iend,j,know)* &
933 & 0.5_r8*
grid(ng)%pm(iend,j)
934# ifdef ADJUST_BOUNDARY
936 tl_bry_pgr=tl_bry_pgr- &
938 & 0.5_r8*
grid(ng)%pm(iend,j)
939 END IF
940# endif
941 ELSE
942 bry_pgr=-
g*(zeta(iend+1,j,know)- &
943 & zeta(iend ,j,know))* &
944 & 0.5_r8*(
grid(ng)%pm(iend ,j)+ &
945 &
grid(ng)%pm(iend+1,j))
946 tl_bry_pgr=-
g*(tl_zeta(iend+1,j,know)- &
947 & tl_zeta(iend ,j,know))* &
948 & 0.5_r8*(
grid(ng)%pm(iend ,j)+ &
949 &
grid(ng)%pm(iend+1,j))
950 END IF
951# ifdef UV_COR
952 bry_cor=0.125_r8*(vbar(iend ,j ,know)+ &
953 & vbar(iend ,j+1,know)+ &
954 & vbar(iend+1,j ,know)+ &
955 & vbar(iend+1,j+1,know))* &
956 & (
grid(ng)%f(iend ,j)+ &
957 &
grid(ng)%f(iend+1,j))
958 tl_bry_cor=0.125_r8*(tl_vbar(iend ,j ,know)+ &
959 & tl_vbar(iend ,j+1,know)+ &
960 & tl_vbar(iend+1,j ,know)+ &
961 & tl_vbar(iend+1,j+1,know))* &
962 & (
grid(ng)%f(iend ,j)+ &
963 &
grid(ng)%f(iend+1,j))
964# else
965 bry_cor=0.0_r8
966 tl_bry_cor=0.0_r8
967# endif
968 cff1=1.0_r8/(0.5_r8*(
grid(ng)%h(iend ,j)+ &
969 & zeta(iend ,j,know)+ &
970 &
grid(ng)%h(iend+1,j)+ &
971 & zeta(iend+1,j,know)))
972 tl_cff1=-cff1*cff1*0.5_r8*(
grid(ng)%tl_h(iend ,j)+ &
973 & tl_zeta(iend ,j,know)+ &
974 &
grid(ng)%tl_h(iend+1,j)+ &
975 & tl_zeta(iend+1,j,know))
976 bry_str=cff1*(
forces(ng)%sustr(iend+1,j)- &
977 &
forces(ng)%bustr(iend+1,j))
978 tl_bry_str=tl_cff1*(
forces(ng)%sustr(iend+1,j)- &
979 &
forces(ng)%bustr(iend+1,j))+ &
980 & cff1*(
forces(ng)%tl_sustr(iend+1,j)- &
981 &
forces(ng)%tl_bustr(iend+1,j))
982 cx=1.0_r8/sqrt(
g*0.5_r8*(
grid(ng)%h(iend+1,j)+ &
983 & zeta(iend+1,j,know)+ &
984 &
grid(ng)%h(iend ,j)+ &
985 & zeta(iend ,j,know)))
986 tl_cx=-cx*cx*cx*0.25_r8*
g*(
grid(ng)%tl_h(iend+1,j)+ &
987 & tl_zeta(iend+1,j,know)+ &
988 &
grid(ng)%tl_h(iend ,j)+ &
989 & tl_zeta(iend ,j,know))
990 cff2=
grid(ng)%om_u(iend+1,j)*cx
991 tl_cff2=
grid(ng)%om_u(iend+1,j)*tl_cx
992
993
994
995
996
997 tl_bry_val=tl_ubar(iend,j,know)+ &
998 & tl_cff2*(bry_pgr+ &
999 & bry_cor+ &
1000 & bry_str)+ &
1001 & cff2*(tl_bry_pgr+ &
1002 & tl_bry_cor+ &
1003 & tl_bry_str)
1004# else
1005
1006
1007# ifdef ADJUST_BOUNDARY
1009 tl_bry_val=
boundary(ng)%tl_ubar_east(j)
1010 ELSE
1011 tl_bry_val=0.0_r8
1012 END IF
1013# else
1014 tl_bry_val=0.0_r8
1015# endif
1016# endif
1017# ifdef WET_DRY_NOT_YET
1018 cff=0.5_r8*(
grid(ng)%h(iend ,j)+ &
1019 & zeta(iend ,j,know)+ &
1020 &
grid(ng)%h(iend+1,j)+ &
1021 & zeta(iend+1,j,know))
1022 tl_cff=0.5_r8*(
grid(ng)%tl_h(iend ,j)+ &
1023 & tl_zeta(iend ,j,know)+ &
1024 &
grid(ng)%tl_h(iend+1,j)+ &
1025 & tl_zeta(iend+1,j,know))
1026# else
1027 cff=0.5_r8*(
grid(ng)%h(iend ,j)+ &
1028 &
grid(ng)%h(iend+1,j))
1029 tl_cff=0.5_r8*(
grid(ng)%tl_h(iend ,j)+ &
1030 &
grid(ng)%tl_h(iend+1,j))
1031# endif
1033 tl_cff1=-0.5_r8*cff1*tl_cff/cff
1034 cx=dt2d*cff1*cff*0.5_r8*(
grid(ng)%pm(iend ,j)+ &
1035 &
grid(ng)%pm(iend+1,j))
1036 tl_cx=dt2d*0.5_r8*(
grid(ng)%pm(iend ,j)+ &
1037 &
grid(ng)%pm(iend+1,j))* &
1038 & (cff1*tl_cff+ &
1039 & tl_cff1*cff)
1040 zx=(0.5_r8+cx)*zeta(iend ,j,know)+ &
1041 & (0.5_r8-cx)*zeta(iend+1,j,know)
1042 tl_zx=(0.5_r8+cx)*tl_zeta(iend ,j,know)+ &
1043 & (0.5_r8-cx)*tl_zeta(iend+1,j,know)+ &
1044 & tl_cx*(zeta(iend ,j,know)- &
1045 & zeta(iend+1,j,know))
1047 cff2=(1.0_r8-
co/cx)**2
1048 tl_cff2=2.0_r8*cff2*
co*tl_cx/(cx*cx)
1049 cff3=zeta(iend,j,kout)+ &
1050 & cx*zeta(iend+1,j,know)- &
1051 & (1.0_r8+cx)*zeta(iend,j,know)
1052 tl_cff3=tl_zeta(iend,j,kout)+ &
1053 & cx*tl_zeta(iend+1,j,know)+ &
1054 & tl_cx*(zeta(iend ,j,know)+ &
1055 & zeta(iend+1,j,know))- &
1056 & (1.0_r8+cx)*tl_zeta(iend,j,know)
1057 zx=zx+cff2*cff3
1058 tl_zx=tl_zx+cff2*tl_cff3+ &
1059 & tl_cff2*cff3
1060 END IF
1061
1062
1063
1064
1065
1066
1067 tl_ubar(iend+1,j,kout)=0.5_r8* &
1068 & ((1.0_r8-cx)* &
1069 & tl_ubar(iend+1,j,know)+ &
1070 & tl_cx*(ubar(iend ,j,know)- &
1071 & ubar(iend+1,j,know))+ &
1072 & cx*tl_ubar(iend,j,know)+ &
1073 & tl_bry_val+ &
1074 & tl_cff1* &
1075 & (zx-
boundary(ng)%zeta_east(j))- &
1076 & cff1*tl_zx)
1077# ifdef ADJUST_BOUNDARY
1079 tl_ubar(iend+1,j,kout)=tl_ubar(iend+1,j,kout)- &
1080 & 0.5_r8*cff1* &
1082 END IF
1083# endif
1084# ifdef MASKING
1085
1086
1087
1088 tl_ubar(iend+1,j,kout)=tl_ubar(iend+1,j,kout)* &
1089 &
grid(ng)%umask(iend+1,j)
1090# endif
1091 END IF
1092 END DO
1093
1094
1095
1097 DO j=jstr,jend
1099
1100
1101# ifdef ADJUST_BOUNDARY
1103 tl_ubar(iend+1,j,kout)=
boundary(ng)%tl_ubar_east(j)
1104 ELSE
1105 tl_ubar(iend+1,j,kout)=0.0_r8
1106 END IF
1107# else
1108 tl_ubar(iend+1,j,kout)=0.0_r8
1109# endif
1110# ifdef MASKING
1111
1112
1113
1114 tl_ubar(iend+1,j,kout)=tl_ubar(iend+1,j,kout)* &
1115 &
grid(ng)%umask(iend+1,j)
1116# endif
1117 END IF
1118 END DO
1119
1120
1121
1123 DO j=jstr,jend
1125
1126
1127 tl_ubar(iend+1,j,kout)=tl_ubar(iend,j,kout)
1128# ifdef MASKING
1129
1130
1131
1132 tl_ubar(iend+1,j,kout)=tl_ubar(iend+1,j,kout)* &
1133 &
grid(ng)%umask(iend+1,j)
1134# endif
1135 END IF
1136 END DO
1137
1138
1139
1141 DO j=jstr,jend
1144
1145
1146
1147
1148 tl_bry_pgr=
g*tl_zeta(iend,j,know)* &
1149 & 0.5_r8*
grid(ng)%pm(iend,j)
1150# ifdef ADJUST_BOUNDARY
1152 tl_bry_pgr=tl_bry_pgr- &
1154 & 0.5_r8*
grid(ng)%pm(iend,j)
1155 END IF
1156# endif
1157 ELSE
1158
1159
1160
1161
1162
1163 tl_bry_pgr=-
g*(tl_zeta(iend+1,j,know)- &
1164 & tl_zeta(iend ,j,know))* &
1165 & 0.5_r8*(
grid(ng)%pm(iend ,j)+ &
1166 &
grid(ng)%pm(iend+1,j))
1167 END IF
1168# ifdef UV_COR
1169
1170
1171
1172
1173
1174
1175
1176 tl_bry_cor=0.125_r8*(tl_vbar(iend, j ,know)+ &
1177 & tl_vbar(iend ,j+1,know)+ &
1178 & tl_vbar(iend+1,j ,know)+ &
1179 & tl_vbar(iend+1,j+1,know))* &
1180 & (
grid(ng)%f(iend ,j)+ &
1181 &
grid(ng)%f(iend+1,j))
1182# else
1183
1184
1185 tl_bry_cor=0.0_r8
1186# endif
1187 cff=1.0_r8/(0.5_r8*(
grid(ng)%h(iend ,j)+ &
1188 & zeta(iend ,j,know)+ &
1189 &
grid(ng)%h(iend+1,j)+ &
1190 & zeta(iend+1,j,know)))
1191 tl_cff=-cff*cff*0.5_r8*(
grid(ng)%tl_h(iend ,j)+ &
1192 & tl_zeta(iend ,j,know)+ &
1193 &
grid(ng)%tl_h(iend+1,j)+ &
1194 & tl_zeta(iend+1,j,know))
1195
1196
1197
1198 tl_bry_str=tl_cff*(
forces(ng)%sustr(iend+1,j)- &
1199 &
forces(ng)%bustr(iend+1,j))+ &
1200 & cff*(
forces(ng)%tl_sustr(iend+1,j)- &
1201 &
forces(ng)%tl_bustr(iend+1,j))
1202
1203
1204
1205
1206
1207 tl_ubar(iend+1,j,kout)=tl_ubar(iend+1,j,know)+ &
1208 & dt2d*(tl_bry_pgr+ &
1209 & tl_bry_cor+ &
1210 & tl_bry_str)
1211# ifdef MASKING
1212
1213
1214
1215 tl_ubar(iend+1,j,kout)=tl_ubar(iend+1,j,kout)* &
1216 &
grid(ng)%umask(iend+1,j)
1217# endif
1218 END IF
1219 END DO
1220
1221
1222
1224 DO j=jstr,jend
1226
1227
1228 tl_ubar(iend+1,j,kout)=0.0_r8
1229 END IF
1230 END DO
1231 END IF
1232 END IF
1233
1234
1235
1236
1237
1238 IF (
domain(ng)%Southern_Edge(tile))
THEN
1239
1240
1241
1243 IF (
iic(ng).ne.0)
THEN
1244 DO i=istru-1,iend
1245
1246
1247
1248 tl_grad(i,jstr-1)=0.0_r8
1249 END DO
1250 DO i=istru,iend
1252# if defined CELERITY_READ && defined FORWARD_READ
1255 obc_out=0.5_r8* &
1256 & (
clima(ng)%M2nudgcof(i-1,jstr-1)+ &
1257 &
clima(ng)%M2nudgcof(i ,jstr-1))
1258 obc_in =
obcfac(ng)*obc_out
1259 ELSE
1262 END IF
1263 IF (
boundary(ng)%ubar_south_Ce(i).lt.0.0_r8)
THEN
1264 tau=obc_in
1265 ELSE
1266 tau=obc_out
1267 END IF
1268 tau=tau*dt2d
1269 END IF
1270# ifdef RADIATION_2D
1272# else
1273 cx=0.0_r8
1274# endif
1277# endif
1278
1279
1280
1281
1282
1283
1284 tl_ubar(i,jstr-1,kout)=(cff*tl_ubar(i,jstr-1,know)+ &
1285 & ce *tl_ubar(i,jstr ,kout)- &
1286 & max(cx,0.0_r8)* &
1287 & tl_grad(i-1,jstr-1)- &
1288 & min(cx,0.0_r8)* &
1289 & tl_grad(i ,jstr-1))/ &
1290 & (cff+ce)
1291
1293
1294
1295
1296
1297 tl_ubar(i,jstr-1,kout)=tl_ubar(i,jstr-1,kout)- &
1298 & tau*tl_ubar(i,jstr-1,know)
1299 END IF
1300# ifdef MASKING
1301
1302
1303
1304 tl_ubar(i,jstr-1,kout)=tl_ubar(i,jstr-1,kout)* &
1305 &
grid(ng)%umask(i,jstr-1)
1306# endif
1307 END IF
1308 END DO
1309 END IF
1310
1311
1312
1316 DO i=istru,iend
1318 cff=dt2d*0.5_r8*(
grid(ng)%pn(i-1,jstr)+ &
1319 &
grid(ng)%pn(i ,jstr))
1320 cff1=sqrt(
g*0.5_r8*(
grid(ng)%h(i-1,jstr)+ &
1321 & zeta(i-1,jstr,know)+ &
1322 &
grid(ng)%h(i ,jstr)+ &
1323 & zeta(i ,jstr,know)))
1324 tl_cff1=0.25_r8*
g*(
grid(ng)%tl_h(i-1,jstr)+ &
1325 & tl_zeta(i-1,jstr,know)+ &
1326 &
grid(ng)%tl_h(i ,jstr)+ &
1327 & tl_zeta(i ,jstr,know))/cff1
1328 ce=cff*cff1
1329 tl_ce=cff*tl_cff1
1330 cff2=1.0_r8/(1.0_r8+ce)
1331 tl_cff2=-cff2*cff2*tl_ce
1332
1333
1334
1335 tl_ubar(i,jstr-1,kout)=tl_cff2*(ubar(i,jstr-1,know)+ &
1336 & ce*ubar(i,jstr,kout))+ &
1337 & cff2*(tl_ubar(i,jstr-1,know)+ &
1338 & tl_ce*ubar(i,jstr,kout)+ &
1339 & ce*tl_ubar(i,jstr,kout))
1340# ifdef MASKING
1341
1342
1343
1344 tl_ubar(i,jstr-1,kout)=tl_ubar(i,jstr-1,kout)* &
1345 &
grid(ng)%umask(i,jstr-1)
1346# endif
1347 END IF
1348 END DO
1349
1350
1351
1353 DO i=istru,iend
1355
1356
1357# ifdef ADJUST_BOUNDARY
1359 tl_ubar(i,jstr-1,kout)=
boundary(ng)%tl_ubar_south(i)
1360 ELSE
1361 tl_ubar(i,jstr-1,kout)=0.0_r8
1362 END IF
1363# else
1364 tl_ubar(i,jstr-1,kout)=0.0_r8
1365# endif
1366# ifdef MASKING
1367
1368
1369
1370 tl_ubar(i,jstr-1,kout)=tl_ubar(i,jstr-1,kout)* &
1371 &
grid(ng)%umask(i,jstr-1)
1372# endif
1373 END IF
1374 END DO
1375
1376
1377
1379 DO i=istru,iend
1381
1382
1383 tl_ubar(i,jstr-1,kout)=tl_ubar(i,jstr,kout)
1384# ifdef MASKING
1385
1386
1387
1388 tl_ubar(i,jstr-1,kout)=tl_ubar(i,jstr-1,kout)* &
1389 &
grid(ng)%umask(i,jstr-1)
1390# endif
1391 END IF
1392 END DO
1393
1394
1395
1396
1399 imin=istru
1400 imax=iend
1401 ELSE
1402 imin=istr
1403 imax=iendr
1404 END IF
1405 DO i=imin,imax
1407
1408
1409 tl_ubar(i,jstr-1,kout)=
gamma2(ng)*tl_ubar(i,jstr,kout)
1410# ifdef MASKING
1411
1412
1413
1414 tl_ubar(i,jstr-1,kout)=tl_ubar(i,jstr-1,kout)* &
1415 &
grid(ng)%umask(i,jstr-1)
1416# endif
1417 END IF
1418 END DO
1419 END IF
1420 END IF
1421
1422
1423
1424
1425
1426 IF (
domain(ng)%Northern_Edge(tile))
THEN
1427
1428
1429
1431 IF (
iic(ng).ne.0)
THEN
1432 DO i=istru-1,iend
1433
1434
1435
1436 tl_grad(i,jend+1)=0.0_r8
1437 END DO
1438 DO i=istru,iend
1440# if defined CELERITY_READ && defined FORWARD_READ
1443 obc_out=0.5_r8* &
1444 & (
clima(ng)%M2nudgcof(i-1,jend+1)+ &
1445 &
clima(ng)%M2nudgcof(i ,jend+1))
1446 obc_in =
obcfac(ng)*obc_out
1447 ELSE
1450 END IF
1451 IF (
boundary(ng)%ubar_north_Ce(i).lt.0.0_r8)
THEN
1452 tau=obc_in
1453 ELSE
1454 tau=obc_out
1455 END IF
1456 tau=tau*dt2d
1457 END IF
1458# ifdef RADIATION_2D
1460# else
1461 cx=0.0_r8
1462# endif
1465# endif
1466
1467
1468
1469
1470
1471
1472 tl_ubar(i,jend+1,kout)=(cff*tl_ubar(i,jend+1,know)+ &
1473 & ce *tl_ubar(i,jend ,kout)- &
1474 & max(cx,0.0_r8)* &
1475 & tl_grad(i-1,jend+1)- &
1476 & min(cx,0.0_r8)* &
1477 & tl_grad(i ,jend+1))/ &
1478 & (cff+ce)
1479
1481
1482
1483
1484
1485 tl_ubar(i,jend+1,kout)=tl_ubar(i,jend+1,kout)- &
1486 & tau*tl_ubar(i,jend+1,know)
1487 END IF
1488# ifdef MASKING
1489
1490
1491
1492 tl_ubar(i,jend+1,kout)=tl_ubar(i,jend+1,kout)* &
1493 &
grid(ng)%umask(i,jend+1)
1494# endif
1495 END IF
1496 END DO
1497 END IF
1498
1499
1500
1504 DO i=istru,iend
1506 cff=dt2d*0.5_r8*(
grid(ng)%pn(i-1,jend)+ &
1507 &
grid(ng)%pn(i ,jend))
1508 cff1=sqrt(
g*0.5_r8*(
grid(ng)%h(i-1,jend)+ &
1509 & zeta(i-1,jend,know)+ &
1510 &
grid(ng)%h(i ,jend)+ &
1511 & zeta(i ,jend,know)))
1512 tl_cff1=0.25_r8*
g*(
grid(ng)%tl_h(i-1,jend)+ &
1513 & tl_zeta(i-1,jend,know)+ &
1514 &
grid(ng)%tl_h(i ,jend)+ &
1515 & tl_zeta(i ,jend,know))/cff1
1516 ce=cff*cff1
1517 tl_ce=cff*tl_cff1
1518 cff2=1.0_r8/(1.0_r8+ce)
1519 tl_cff2=-cff2*cff2*tl_ce
1520
1521
1522
1523 tl_ubar(i,jend+1,kout)=tl_cff2*(ubar(i,jend+1,know)+ &
1524 & ce*ubar(i,jend,kout))+ &
1525 & cff2*(tl_ubar(i,jend+1,know)+ &
1526 & tl_ce*ubar(i,jend,kout)+ &
1527 & ce*tl_ubar(i,jend,kout))
1528# ifdef MASKING
1529
1530
1531
1532 tl_ubar(i,jend+1,kout)=tl_ubar(i,jend+1,kout)* &
1533 &
grid(ng)%umask(i,jend+1)
1534# endif
1535 END IF
1536 END DO
1537
1538
1539
1541 DO i=istru,iend
1543
1544
1545# ifdef ADJUST_BOUNDARY
1547 tl_ubar(i,jend+1,kout)=
boundary(ng)%tl_ubar_north(i)
1548 ELSE
1549 tl_ubar(i,jend+1,kout)=0.0_r8
1550 END IF
1551# else
1552 tl_ubar(i,jend+1,kout)=0.0_r8
1553# endif
1554# ifdef MASKING
1555
1556
1557
1558 tl_ubar(i,jend+1,kout)=tl_ubar(i,jend+1,kout)* &
1559 &
grid(ng)%umask(i,jend+1)
1560# endif
1561 END IF
1562 END DO
1563
1564
1565
1567 DO i=istru,iend
1569
1570
1571 tl_ubar(i,jend+1,kout)=tl_ubar(i,jend,kout)
1572# ifdef MASKING
1573
1574
1575
1576 tl_ubar(i,jend+1,kout)=tl_ubar(i,jend+1,kout)* &
1577 &
grid(ng)%umask(i,jend+1)
1578# endif
1579 END IF
1580 END DO
1581
1582
1583
1584
1587 imin=istru
1588 imax=iend
1589 ELSE
1590 imin=istr
1591 imax=iendr
1592 END IF
1593 DO i=imin,imax
1595
1596
1597 tl_ubar(i,jend+1,kout)=
gamma2(ng)*tl_ubar(i,jend,kout)
1598
1599# ifdef MASKING
1600
1601
1602
1603 tl_ubar(i,jend+1,kout)=tl_ubar(i,jend+1,kout)* &
1604 &
grid(ng)%umask(i,jend+1)
1605# endif
1606 END IF
1607 END DO
1608 END IF
1609 END IF
1610
1611
1612
1613
1614
1616 IF (
domain(ng)%SouthWest_Corner(tile))
THEN
1619
1620
1621
1622 tl_ubar(istr,jstr-1,kout)=0.5_r8* &
1623 & (tl_ubar(istr+1,jstr-1,kout)+ &
1624 & tl_ubar(istr ,jstr ,kout))
1625 END IF
1626 END IF
1627 IF (
domain(ng)%SouthEast_Corner(tile))
THEN
1630
1631
1632
1633 tl_ubar(iend+1,jstr-1,kout)=0.5_r8* &
1634 & (tl_ubar(iend ,jstr-1,kout)+ &
1635 & tl_ubar(iend+1,jstr ,kout))
1636 END IF
1637 END IF
1638 IF (
domain(ng)%NorthWest_Corner(tile))
THEN
1641
1642
1643
1644 tl_ubar(istr,jend+1,kout)=0.5_r8* &
1645 & (tl_ubar(istr ,jend ,kout)+ &
1646 & tl_ubar(istr+1,jend+1,kout))
1647 END IF
1648 END IF
1649 IF (
domain(ng)%NorthEast_Corner(tile))
THEN
1652
1653
1654
1655 tl_ubar(iend+1,jend+1,kout)=0.5_r8* &
1656 & (tl_ubar(iend+1,jend ,kout)+ &
1657 & tl_ubar(iend ,jend+1,kout))
1658 END IF
1659 END IF
1660 END IF
1661
1662# if defined WET_DRY_NOT_YET
1663
1664
1665
1666
1667
1668
1669
1671 IF (
domain(ng)%Western_Edge(tile))
THEN
1672 DO j=jstr,jend
1675
1676
1677
1678
1679
1680
1681 END IF
1682 END DO
1683 END IF
1684 IF (
domain(ng)%Eastern_Edge(tile))
THEN
1685 DO j=jstr,jend
1688
1689
1690
1691
1692
1693
1694 END IF
1695 END DO
1696 END IF
1697 END IF
1698
1700 IF (
domain(ng)%Southern_Edge(tile))
THEN
1701 DO i=istru,iend
1704
1705
1706
1707
1708
1709
1710 END IF
1711 END DO
1712 END IF
1713 IF (
domain(ng)%Northern_Edge(tile))
THEN
1714 DO i=istr,iend
1717
1718
1719
1720
1721
1722
1723 END IF
1724 END DO
1725 END IF
1726 END IF
1727
1729 IF (
domain(ng)%SouthWest_Corner(tile))
THEN
1734
1735
1736
1737
1738
1739
1740 END IF
1741 END IF
1742 IF (
domain(ng)%SouthEast_Corner(tile))
THEN
1747
1748
1749
1750
1751
1752
1753 END IF
1754 END IF
1755 IF (
domain(ng)%NorthWest_Corner(tile))
THEN
1760
1761
1762
1763
1764
1765
1766 END IF
1767 END IF
1768 IF (
domain(ng)%NorthEast_Corner(tile))
THEN
1773
1774
1775
1776
1777
1778
1779 END IF
1780 END IF
1781 END IF
1782# endif
1783
1784 RETURN
type(t_boundary), dimension(:), allocatable boundary
type(t_apply), dimension(:), allocatable lbc_apply
type(t_clima), dimension(:), allocatable clima
type(t_forces), dimension(:), allocatable forces
type(t_grid), dimension(:), allocatable grid
type(t_lbc), dimension(:,:,:), allocatable tl_lbc
type(t_lbc), dimension(:,:,:), allocatable lbc
type(t_domain), dimension(:), allocatable domain
logical, dimension(:), allocatable lnudgem2clm
integer, dimension(:), allocatable iic
logical, dimension(:,:,:), allocatable lobc
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
logical, dimension(:), allocatable predictor_2d_step
real(dp), dimension(:), allocatable obcfac
real(r8), dimension(:), allocatable gamma2
integer, parameter isouth
real(dp), dimension(:), allocatable dtfast
real(dp), dimension(:,:), allocatable m2obc_out
integer, parameter inorth
real(dp), dimension(:,:), allocatable m2obc_in