57
58
64
65
66
67 integer, intent(in) :: ng, tile
68 integer, intent(in) :: LBi, UBi, LBj, UBj
69 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
70 integer, intent(in) :: krhs, kstp, kout
71
72# ifdef ASSUMED_SHAPE
73 real(r8), intent(in) :: zeta(LBi:,LBj:,:)
74
75 real(r8), intent(inout) :: ad_zeta(LBi:,LBj:,:)
76# else
77 real(r8), intent(in) :: zeta(LBi:UBi,LBj:UBj,3)
78
79 real(r8), intent(inout) :: ad_zeta(LBi:UBi,LBj:UBj,3)
80# endif
81
82
83
84 integer :: i, j, know
85
86 real(r8) :: Ce, Cx
87 real(r8) :: cff, cff1, cff2, dt2d, tau
88
89 real(r8) :: ad_Ce, ad_Cx
90 real(r8) :: ad_cff1, ad_cff2
91 real(r8) :: adfac
92
93 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_grad
94
95# include "set_bounds.h"
96
97
98
99
100
101 ad_ce=0.0_r8
102 ad_cx=0.0_r8
103 ad_cff1=0.0_r8
104 ad_cff2=0.0_r8
105
106 ad_grad(lbi:ubi,lbj:ubj)=0.0_r8
107
108
109
110
111
112 IF (first_2d_step) THEN
113 know=krhs
116 know=krhs
118 ELSE
119 know=kstp
121 END IF
122
123# if defined WET_DRY
124
125
126
127
128
130
132 IF (
domain(ng)%NorthEast_Corner(tile))
THEN
135 IF (zeta(iend+1,jend+1,kout).le. &
136 & (
dcrit(ng)-
grid(ng)%h(iend+1,jend+1)))
THEN
137
138
139 grid(ng)%ad_h(iend+1,jend+1)=
grid(ng)%ad_h(iend+1,jend+1)-&
140 & ad_zeta(iend+1,jend+1,kout)
141 ad_zeta(iend+1,jend+1,kout)=0.0_r8
142 END IF
143 END IF
144 END IF
145 IF (
domain(ng)%NorthWest_Corner(tile))
THEN
148 IF (zeta(istr-1,jend+1,kout).le. &
149 & (
dcrit(ng)-
grid(ng)%h(istr-1,jend+1)))
THEN
150
151
152 grid(ng)%ad_h(istr-1,jend+1)=
grid(ng)%ad_h(istr-1,jend+1)-&
153 & ad_zeta(istr-1,jend+1,kout)
154 ad_zeta(istr-1,jend+1,kout)=0.0_r8
155 END IF
156 END IF
157 END IF
158 IF (
domain(ng)%SouthEast_Corner(tile))
THEN
161 IF (zeta(iend+1,jstr-1,kout).le. &
162 & (
dcrit(ng)-
grid(ng)%h(iend+1,jstr-1)))
THEN
163
164
165 grid(ng)%ad_h(iend+1,jstr-1)=
grid(ng)%ad_h(iend+1,jstr-1)-&
166 & ad_zeta(iend+1,jstr-1,kout)
167 tl_zeta(iend+1,jstr-1,kout)=0.0_r8
168 END IF
169 END IF
170 END IF
171 IF (
domain(ng)%SouthWest_Corner(tile))
THEN
174 IF (zeta(istr-1,jstr-1,kout).le. &
175 & (
dcrit(ng)-
grid(ng)%h(istr-1,jstr-1)))
THEN
176
177
178 grid(ng)%ad_h(istr-1,jstr-1)=
grid(ng)%ad_h(istr-1,jstr-1)-&
179 & ad_zeta(istr-1,jstr-1,kout)
180 ad_zeta(istr-1,jstr-1,kout)=0.0_r8
181 END IF
182 END IF
183 END IF
184 END IF
185
187 IF (
domain(ng)%Northern_Edge(tile))
THEN
188 DO i=istr,iend
190 IF (zeta(i,jend+1,kout).le. &
192
193
194 grid(ng)%ad_h(i,jend+1)=
grid(ng)%ad_h(i,jend+1)- &
195 & ad_zeta(i,jend+1,kout)
196 ad_zeta(i,jend+1,kout)=0.0_r8
197 END IF
198 END IF
199 END DO
200 END IF
201 IF (
domain(ng)%Southern_Edge(tile))
THEN
202 DO i=istr,iend
204 IF (zeta(i,jstr-1,kout).le. &
206
207
208 grid(ng)%ad_h(i,jstr-1)=
grid(ng)%ad_h(i,jstr-1)- &
209 & ad_zeta(i,jstr-1,kout)
210 ad_zeta(i,jstr-1,kout)=0.0_r8
211 END IF
212 END IF
213 END DO
214 END IF
215 END IF
216
218 IF (
domain(ng)%Eastern_Edge(tile))
THEN
219 DO j=jstr,jend
221 IF (zeta(iend+1,j,kout).le. &
223
224
225 grid(ng)%ad_h(iend+1,j)=
grid(ng)%ad_h(iend+1,j)- &
226 & ad_zeta(iend+1,j,kout)
227 ad_zeta(iend+1,j,kout)=0.0_r8
228 END IF
229 END IF
230 END DO
231 END IF
232 IF (
domain(ng)%Western_Edge(tile))
THEN
233 DO j=jstr,jend
235 IF (zeta(istr-1,j,kout).le. &
237
238
239 grid(ng)%ad_h(istr-1,j)=
grid(ng)%ad_h(istr-1,j)- &
240 & ad_zeta(istr-1,j,kout)
241 ad_zeta(istr-1,j,kout)=0.0_r8
242 END IF
243 END IF
244 END DO
245 END IF
246 END IF
247# endif
248
249
250
251
252
254 IF (
domain(ng)%NorthEast_Corner(tile))
THEN
257
258
259
260
261 adfac=0.5_r8*ad_zeta(iend+1,jend+1,kout)
262 ad_zeta(iend ,jend+1,kout)=ad_zeta(iend ,jend+1,kout)+ &
263 & adfac
264 ad_zeta(iend+1,jend ,kout)=ad_zeta(iend+1,jend ,kout)+ &
265 & adfac
266 ad_zeta(iend+1,jend+1,kout)=0.0_r8
267 END IF
268 END IF
269 IF (
domain(ng)%NorthWest_Corner(tile))
THEN
272
273
274
275
276 adfac=0.5_r8*ad_zeta(istr-1,jend+1,kout)
277 ad_zeta(istr-1,jend ,kout)=ad_zeta(istr-1,jend ,kout)+ &
278 & adfac
279 ad_zeta(istr ,jend+1,kout)=ad_zeta(istr ,jend+1,kout)+ &
280 & adfac
281 ad_zeta(istr-1,jend+1,kout)=0.0_r8
282 END IF
283 END IF
284 IF (
domain(ng)%SouthEast_Corner(tile))
THEN
287
288
289
290
291 adfac=0.5_r8*ad_zeta(iend+1,jstr-1,kout)
292 ad_zeta(iend ,jstr-1,kout)=ad_zeta(iend ,jstr-1,kout)+ &
293 & adfac
294 ad_zeta(iend+1,jstr ,kout)=ad_zeta(iend+1,jstr ,kout)+ &
295 & adfac
296 ad_zeta(iend+1,jstr-1,kout)=0.0_r8
297 END IF
298 END IF
299 IF (
domain(ng)%SouthWest_Corner(tile))
THEN
302
303
304
305
306 adfac=0.5_r8*ad_zeta(istr-1,jstr-1,kout)
307 ad_zeta(istr ,jstr-1,kout)=ad_zeta(istr ,jstr-1,kout)+ &
308 & adfac
309 ad_zeta(istr-1,jstr ,kout)=ad_zeta(istr-1,jstr ,kout)+ &
310 & adfac
311 ad_zeta(istr-1,jstr-1,kout)=0.0_r8
312 END IF
313 END IF
314 END IF
315
316
317
318
319
320 IF (
domain(ng)%Northern_Edge(tile))
THEN
321
322
323
325 IF (
iic(ng).ne.0)
THEN
326 DO i=istr,iend
328# if defined CELERITY_READ && defined FORWARD_READ
330 IF (
boundary(ng)%zeta_north_Ce(i).eq.0.0_r8)
THEN
332 ELSE
334 END IF
335 tau=tau*dt2d
336 END IF
337# ifdef RADIATION_2D
339# else
340 cx=0.0_r8
341# endif
344# endif
345# ifdef MASKING
346
347
348
349 ad_zeta(i,jend+1,kout)=ad_zeta(i,jend+1,kout)* &
350 &
grid(ng)%rmask(i,jend+1)
351# endif
353
354
355
356 ad_zeta(i,jend+1,know)=ad_zeta(i,jend+1,know)- &
357 & tau*ad_zeta(i,jend+1,kout)
358 END IF
359
360
361
362
363
364
365
366
367 adfac=ad_zeta(i,jend+1,kout)/(cff+ce)
368 ad_grad(i ,jend+1)=ad_grad(i ,jend+1)- &
369 & max(cx,0.0_r8)*adfac
370 ad_grad(i+1,jend+1)=ad_grad(i+1,jend+1)- &
371 & min(cx,0.0_r8)*adfac
372 ad_zeta(i,jend ,kout)=ad_zeta(i,jend ,kout)+ce *adfac
373 ad_zeta(i,jend+1,know)=ad_zeta(i,jend+1,know)+cff*adfac
374 ad_zeta(i,jend+1,kout)=0.0_r8
375 END IF
376 END DO
377 END IF
378
379
380
382 DO i=istr,iend
384 cff=dt2d*
grid(ng)%pn(i,jend)
385 cff1=sqrt(
g*(
grid(ng)%h(i,jend)+ &
386 & zeta(i,jend,know)))
387 ce=cff*cff1
388# ifdef MASKING
389
390
391
392 ad_zeta(i,jend+1,kout)=ad_zeta(i,jend+1,kout)* &
393 &
grid(ng)%rmask(i,jend+1)
394# endif
395
396
397
398
399
400 ad_zeta(i,jend+1,know)=ad_zeta(i,jend+1,know)+ &
401 & (1.0_r8-ce)*ad_zeta(i,jend+1,kout)
402 ad_ce=ad_ce+(zeta(i,jend+1,know)+ &
403 & zeta(i,jend ,know))*ad_zeta(i,jend+1,kout)
404 ad_zeta(i,jend,know)=ad_zeta(i,jend,know)+ &
405 & ce*ad_zeta(i,jend+1,kout)
406 ad_zeta(i,jend+1,kout)=0.0_r8
407
408
409 ad_cff1=ad_cff1+cff*ad_ce
410 ad_ce=0.0_r8
411
412
413
414 adfac=0.5_r8*
g*ad_cff1/cff1
415 grid(ng)%ad_h(i,jend)=
grid(ng)%ad_h(i,jend)+adfac
416 ad_zeta(i,jend,know)=ad_zeta(i,jend,know)+adfac
417 ad_cff1=0.0_r8
418 END IF
419 END DO
420
421
422
424 DO i=istr,iend
426 cff=dt2d*
grid(ng)%pn(i,jend)
427 cff1=sqrt(
g*(
grid(ng)%h(i,jend)+ &
428 & zeta(i,jend,know)))
429 ce=cff*cff1
430 cff2=1.0_r8/(1.0_r8+ce)
431# ifdef MASKING
432
433
434
435 ad_zeta(i,jend+1,kout)=ad_zeta(i,jend+1,kout)* &
436 &
grid(ng)%rmask(i,jend+1)
437# endif
438
439
440
441
442
443
444 adfac=cff2*ad_zeta(i,jend+1,kout)
445 ad_zeta(i,jend ,kout)=ad_zeta(i,jend ,kout)+ce*adfac
446 ad_zeta(i,jend+1,know)=ad_zeta(i,jend+1,know)+adfac
447 ad_ce=ad_ce+zeta(i,jend,kout)*adfac
448 ad_cff2=ad_cff2+ &
449 & (zeta(i,jend+1,know)+ &
450 & ce*zeta(i,jend,kout))*ad_zeta(i,jend+1,kout)
451 ad_zeta(i,jend+1,kout)=0.0_r8
452
453
454 ad_ce=ad_ce-cff2*cff2*ad_cff2
455 ad_cff2=0.0_r8
456
457
458 ad_cff1=ad_cff1+cff*ad_ce
459 ad_ce=0.0_r8
460
461
462
463 adfac=0.5_r8*
g*ad_cff1/cff1
464 grid(ng)%ad_h(i,jend)=
grid(ng)%ad_h(i,jend)+adfac
465 ad_zeta(i,jend,know)=ad_zeta(i,jend,know)+adfac
466 ad_cff1=0.0_r8
467 END IF
468 END DO
469
470
471
473 DO i=istr,iend
475# ifdef MASKING
476
477
478
479 ad_zeta(i,jend+1,kout)=ad_zeta(i,jend+1,kout)* &
480 &
grid(ng)%rmask(i,jend+1)
481# endif
482# ifdef ADJUST_BOUNDARY
484
485
487 & ad_zeta_north(i)+ &
488 & ad_zeta(i,jend+1,kout)
489 ad_zeta(i,jend+1,kout)=0.0_r8
490 ELSE
491
492
493 ad_zeta(i,jend+1,kout)=0.0_r8
494 END IF
495# else
496
497
498 ad_zeta(i,jend+1,kout)=0.0_r8
499# endif
500 END IF
501 END DO
502
503
504
506 DO i=istr,iend
508# ifdef MASKING
509
510
511
512 ad_zeta(i,jend+1,kout)=ad_zeta(i,jend+1,kout)* &
513 &
grid(ng)%rmask(i,jend+1)
514# endif
515
516
517 ad_zeta(i,jend ,kout)=ad_zeta(i,jend ,kout)+ &
518 & ad_zeta(i,jend+1,kout)
519 ad_zeta(i,jend+1,kout)=0.0_r8
520 END IF
521 END DO
522
523
524
526 DO i=istr,iend
528# ifdef MASKING
529
530
531
532 ad_zeta(i,jend+1,kout)=ad_zeta(i,jend+1,kout)* &
533 &
grid(ng)%rmask(i,jend+1)
534# endif
535
536
537 ad_zeta(i,jend ,kout)=ad_zeta(i,jend ,kout)+ &
538 & ad_zeta(i,jend+1,kout)
539 ad_zeta(i,jend+1,kout)=0.0_r8
540 END IF
541 END DO
542 END IF
543 END IF
544
545
546
547
548
549 IF (
domain(ng)%Southern_Edge(tile))
THEN
550
551
552
554 IF (
iic(ng).ne.0)
THEN
555 DO i=istr,iend
557# if defined CELERITY_READ && defined FORWARD_READ
559 IF (
boundary(ng)%zeta_south_Ce(i).eq.0.0_r8)
THEN
561 ELSE
563 END IF
564 tau=tau*dt2d
565 END IF
566# ifdef RADIATION_2D
568# else
569 cx=0.0_r8
570# endif
573# endif
574# ifdef MASKING
575
576
577
578 ad_zeta(i,jstr-1,kout)=ad_zeta(i,jstr-1,kout)* &
579 &
grid(ng)%rmask(i,jstr-1)
580# endif
582
583
584
585 ad_zeta(i,jstr-1,know)=ad_zeta(i,jstr-1,know)- &
586 & tau*ad_zeta(i,jstr-1,kout)
587 END IF
588
589
590
591
592
593
594
595
596 adfac=ad_zeta(i,jstr-1,kout)/(cff+ce)
597 ad_grad(i ,jstr-1)=ad_grad(i ,jstr-1)- &
598 & max(cx,0.0_r8)*adfac
599 ad_grad(i+1,jstr-1)=ad_grad(i+1,jstr-1)- &
600 & min(cx,0.0_r8)*adfac
601 ad_zeta(i,jstr-1,know)=ad_zeta(i,jstr-1,know)+cff*adfac
602 ad_zeta(i,jstr ,kout)=ad_zeta(i,jstr ,kout)+ce *adfac
603 ad_zeta(i,jstr-1,kout)=0.0_r8
604 END IF
605 END DO
606 END IF
607
608
609
611 DO i=istr,iend
613 cff=dt2d*
grid(ng)%pn(i,jstr)
614 cff1=sqrt(
g*(
grid(ng)%h(i,jstr)+ &
615 & zeta(i,jstr,know)))
616 ce=cff*cff1
617# ifdef MASKING
618
619
620
621 ad_zeta(i,jstr-1,kout)=ad_zeta(i,jstr-1,kout)* &
622 &
grid(ng)%rmask(i,jstr-1)
623# endif
624
625
626
627
628
629 ad_zeta(i,jstr-1,know)=ad_zeta(i,jstr-1,know)+ &
630 & (1.0_r8-ce)*ad_zeta(i,jstr-1,kout)
631 ad_ce=ad_ce+(zeta(i,jstr-1,know)+ &
632 & zeta(i,jstr ,know))*ad_zeta(i,jstr-1,kout)
633 ad_zeta(i,jstr,know)=ad_zeta(i,jstr,know)+ &
634 & ce*ad_zeta(i,jstr-1,kout)
635 ad_zeta(i,jstr-1,kout)=0.0_r8
636
637
638 ad_cff1=ad_cff1+cff*ad_ce
639 ad_ce=0.0_r8
640
641
642
643 adfac=0.5_r8*
g*ad_cff1/cff1
644 grid(ng)%ad_h(i,jstr)=
grid(ng)%ad_h(i,jstr)+adfac
645 ad_zeta(i,jstr,know)=ad_zeta(i,jstr,know)+adfac
646 ad_cff1=0.0_r8
647 END IF
648 END DO
649
650
651
653 DO i=istr,iend
655 cff=dt2d*
grid(ng)%pn(i,jstr)
656 cff1=sqrt(
g*(
grid(ng)%h(i,jstr)+ &
657 & zeta(i,jstr,know)))
658 ce=cff*cff1
659 cff2=1.0_r8/(1.0_r8+ce)
660# ifdef MASKING
661
662
663
664 ad_zeta(i,jstr-1,kout)=ad_zeta(i,jstr-1,kout)* &
665 &
grid(ng)%rmask(i,jstr-1)
666# endif
667
668
669
670
671
672
673 adfac=cff2*ad_zeta(i,jstr-1,kout)
674 ad_zeta(i,jstr-1,know)=ad_zeta(i,jstr-1,know)+adfac
675 ad_zeta(i,jstr ,kout)=ad_zeta(i,jstr ,kout)+ce*adfac
676 ad_ce=ad_ce+zeta(i,jstr,kout)*adfac
677 ad_cff2=ad_cff2+ &
678 & (zeta(i,jstr-1,know)+ &
679 & ce*zeta(i,jstr,kout))*ad_zeta(i,jstr-1,kout)
680 ad_zeta(i,jstr-1,kout)=0.0_r8
681
682
683 ad_ce=ad_ce-cff2*cff2*ad_cff2
684 ad_cff2=0.0_r8
685
686
687 ad_cff1=ad_cff1+cff*ad_ce
688 ad_ce=0.0_r8
689
690
691
692 adfac=0.5_r8*
g*ad_cff1/cff1
693 grid(ng)%ad_h(i,jstr)=
grid(ng)%ad_h(i,jstr)+adfac
694 ad_zeta(i,jstr,know)=ad_zeta(i,jstr,know)+adfac
695 ad_cff1=0.0_r8
696 END IF
697 END DO
698
699
700
702 DO i=istr,iend
704# ifdef MASKING
705
706
707
708 ad_zeta(i,jstr-1,kout)=ad_zeta(i,jstr-1,kout)* &
709 &
grid(ng)%rmask(i,jstr-1)
710# endif
711# ifdef ADJUST_BOUNDARY
713
714
716 & ad_zeta_south(i)+ &
717 & ad_zeta(i,jstr-1,kout)
718 ad_zeta(i,jstr-1,kout)=0.0_r8
719 ELSE
720
721
722 ad_zeta(i,jstr-1,kout)=0.0_r8
723 END IF
724# else
725
726
727 ad_zeta(i,jstr-1,kout)=0.0_r8
728# endif
729 END IF
730 END DO
731
732
733
735 DO i=istr,iend
737# ifdef MASKING
738
739
740
741 ad_zeta(i,jstr-1,kout)=ad_zeta(i,jstr-1,kout)* &
742 &
grid(ng)%rmask(i,jstr-1)
743# endif
744
745
746 ad_zeta(i,jstr ,kout)=ad_zeta(i,jstr ,kout)+ &
747 & ad_zeta(i,jstr-1,kout)
748 ad_zeta(i,jstr-1,kout)=0.0_r8
749 END IF
750 END DO
751
752
753
755 DO i=istr,iend
757# ifdef MASKING
758
759
760
761 ad_zeta(i,jstr-1,kout)=ad_zeta(i,jstr-1,kout)* &
762 &
grid(ng)%rmask(i,jstr-1)
763# endif
764
765
766 ad_zeta(i,jstr ,kout)=ad_zeta(i,jstr ,kout)+ &
767 & ad_zeta(i,jstr-1,kout)
768 ad_zeta(i,jstr-1,kout)=0.0_r8
769 END IF
770 END DO
771 END IF
772 END IF
773
774
775
776
777
778 IF (
domain(ng)%Eastern_Edge(tile))
THEN
779
780
781
783 IF (
iic(ng).ne.0)
THEN
784 DO j=jstr,jend
786# if defined CELERITY_READ && defined FORWARD_READ
788 IF (
boundary(ng)%zeta_east_Cx(j).eq.0.0_r8)
THEN
790 ELSE
792 END IF
793 tau=tau*dt2d
794 END IF
796# ifdef RADIATION_2D
798# else
799 ce=0.0_r8
800# endif
802# endif
803# ifdef MASKING
804
805
806
807 ad_zeta(iend+1,j,kout)=ad_zeta(iend+1,j,kout)* &
808 &
grid(ng)%rmask(iend+1,j)
809# endif
811
812
813
814 ad_zeta(iend+1,j,know)=ad_zeta(iend+1,j,know)- &
815 & tau*ad_zeta(iend+1,j,kout)
816 END IF
817
818
819
820
821
822
823
824
825 adfac=ad_zeta(iend+1,j,kout)/(cff+cx)
826 ad_grad(iend+1,j )=ad_grad(iend+1,j )- &
827 & max(ce,0.0_r8)*adfac
828 ad_grad(iend+1,j+1)=ad_grad(iend+1,j+1)- &
829 & min(ce,0.0_r8)*adfac
830 ad_zeta(iend ,j,kout)=ad_zeta(iend ,j,kout)+cx *adfac
831 ad_zeta(iend+1,j,know)=ad_zeta(iend+1,j,know)+cff*adfac
832 ad_zeta(iend+1,j,kout)=0.0_r8
833 END IF
834 END DO
835 END IF
836
837
838
840 DO j=jstr,jend
842 cff=dt2d*
grid(ng)%pm(iend,j)
843 cff1=sqrt(
g*(
grid(ng)%h(iend,j)+ &
844 & zeta(iend,j,know)))
845 cx=cff*cff1
846# ifdef MASKING
847
848
849
850 ad_zeta(iend+1,j,kout)=ad_zeta(iend+1,j,kout)* &
851 &
grid(ng)%rmask(iend+1,j)
852# endif
853
854
855
856
857
858 ad_zeta(iend+1,j,know)=ad_zeta(iend+1,j,know)+ &
859 & (1.0_r8-cx)*ad_zeta(iend+1,j,kout)
860 ad_cx=ad_cx+(zeta(iend+1,j,know)+ &
861 & zeta(iend ,j,know))*ad_zeta(iend+1,j,kout)
862 ad_zeta(iend,j,know)=ad_zeta(iend,j,know)+ &
863 & cx*ad_zeta(iend+1,j,kout)
864 ad_zeta(iend+1,j,kout)=0.0_r8
865
866
867 ad_cff1=ad_cff1+cff*ad_cx
868 ad_cx=0.0_r8
869
870
871
872 adfac=0.5_r8*
g*ad_cff1/cff1
873 grid(ng)%ad_h(iend,j)=
grid(ng)%ad_h(iend,j)+adfac
874 ad_zeta(iend,j,know)=ad_zeta(iend,j,know)+adfac
875 ad_cff1=0.0_r8
876 END IF
877 END DO
878
879
880
882 DO j=jstr,jend
884 cff=dt2d*
grid(ng)%pm(iend,j)
885 cff1=sqrt(
g*(
grid(ng)%h(iend,j)+ &
886 & zeta(iend,j,know)))
887 cx=cff*cff1
888 cff2=1.0_r8/(1.0_r8+cx)
889# ifdef MASKING
890
891
892
893 ad_zeta(iend+1,j,kout)=ad_zeta(iend+1,j,kout)* &
894 &
grid(ng)%rmask(iend+1,j)
895# endif
896
897
898
899
900
901
902 adfac=cff2*ad_zeta(iend+1,j,kout)
903 ad_zeta(iend ,j,kout)=ad_zeta(iend ,j,kout)+cx*adfac
904 ad_zeta(iend+1,j,know)=ad_zeta(iend+1,j,know)+adfac
905 ad_cx=ad_cx+zeta(iend,j,kout)*adfac
906 ad_cff2=ad_cff2+ &
907 & (zeta(iend+1,j,know)+ &
908 & cx*zeta(iend,j,kout))*ad_zeta(iend+1,j,kout)
909 ad_zeta(iend+1,j,kout)=0.0_r8
910
911
912 ad_cx=ad_cx-cff2*cff2*ad_cff2
913 ad_cff2=0.0_r8
914
915
916 ad_cff1=ad_cff1+cff*ad_cx
917 ad_cx=0.0_r8
918
919
920
921 adfac=0.5_r8*
g*ad_cff1/cff1
922 grid(ng)%ad_h(iend,j)=
grid(ng)%ad_h(iend,j)+adfac
923 ad_zeta(iend,j,know)=ad_zeta(iend,j,know)+adfac
924 ad_cff1=0.0_r8
925 END IF
926 END DO
927
928
929
931 DO j=jstr,jend
933# ifdef MASKING
934
935
936
937 ad_zeta(iend+1,j,kout)=ad_zeta(iend+1,j,kout)* &
938 &
grid(ng)%rmask(iend+1,j)
939# endif
940# ifdef ADJUST_BOUNDARY
942
943
945 & ad_zeta_east(j)+ &
946 & ad_zeta(iend+1,j,kout)
947 ad_zeta(iend+1,j,kout)=0.0_r8
948 ELSE
949
950
951 ad_zeta(iend+1,j,kout)=0.0_r8
952 END IF
953# else
954
955
956 ad_zeta(iend+1,j,kout)=0.0_r8
957# endif
958 END IF
959 END DO
960
961
962
964 DO j=jstr,jend
966# ifdef MASKING
967
968
969
970 ad_zeta(iend+1,j,kout)=ad_zeta(iend+1,j,kout)* &
971 &
grid(ng)%rmask(iend+1,j)
972# endif
973
974
975 ad_zeta(iend ,j,kout)=ad_zeta(iend ,j,kout)+ &
976 & ad_zeta(iend+1,j,kout)
977 ad_zeta(iend+1,j,kout)=0.0_r8
978 END IF
979 END DO
980
981
982
984 DO j=jstr,jend
986# ifdef MASKING
987
988
989
990 ad_zeta(iend+1,j,kout)=ad_zeta(iend+1,j,kout)* &
991 &
grid(ng)%rmask(iend+1,j)
992# endif
993
994
995 ad_zeta(iend ,j,kout)=ad_zeta(iend ,j,kout)+ &
996 & ad_zeta(iend+1,j,kout)
997 ad_zeta(iend+1,j,kout)=0.0_r8
998 END IF
999 END DO
1000 END IF
1001 END IF
1002
1003
1004
1005
1006
1007 IF (
domain(ng)%Western_Edge(tile))
THEN
1008
1009
1010
1012 IF (
iic(ng).ne.0)
THEN
1013 DO j=jstr,jend
1015# if defined CELERITY_READ && defined FORWARD_READ
1017 IF (
boundary(ng)%zeta_west_Cx(j).eq.0.0_r8)
THEN
1019 ELSE
1021 END IF
1022 tau=tau*dt2d
1023 END IF
1025# ifdef RADIATION_2D
1027# else
1028 ce=0.0_r8
1029# endif
1031# endif
1032# ifdef MASKING
1033
1034
1035
1036 ad_zeta(istr-1,j,kout)=ad_zeta(istr-1,j,kout)* &
1037 &
grid(ng)%rmask(istr-1,j)
1038# endif
1040
1041
1042
1043 ad_zeta(istr-1,j,know)=ad_zeta(istr-1,j,know)- &
1044 & tau*ad_zeta(istr-1,j,kout)
1045
1046 END IF
1047
1048
1049
1050
1051
1052
1053
1054
1055 adfac=ad_zeta(istr-1,j,kout)/(cff+cx)
1056 ad_grad(istr-1,j )=ad_grad(istr-1,j )- &
1057 & max(ce,0.0_r8)*adfac
1058 ad_grad(istr-1,j+1)=ad_grad(istr-1,j+1)- &
1059 & min(ce,0.0_r8)*adfac
1060 ad_zeta(istr-1,j,know)=ad_zeta(istr-1,j,know)+cff*adfac
1061 ad_zeta(istr ,j,kout)=ad_zeta(istr ,j,kout)+cx *adfac
1062 ad_zeta(istr-1,j,kout)=0.0_r8
1063 END IF
1064 END DO
1065 END IF
1066
1067
1068
1070 DO j=jstr,jend
1072 cff=dt2d*
grid(ng)%pm(istr,j)
1073 cff1=sqrt(
g*(
grid(ng)%h(istr,j)+ &
1074 & zeta(istr,j,know)))
1075 cx=cff*cff1
1076# ifdef MASKING
1077
1078
1079
1080 ad_zeta(istr-1,j,kout)=ad_zeta(istr-1,j,kout)* &
1081 &
grid(ng)%rmask(istr-1,j)
1082# endif
1083
1084
1085
1086
1087
1088 ad_zeta(istr-1,j,know)=ad_zeta(istr-1,j,know)+ &
1089 & (1.0_r8-cx)*ad_zeta(istr-1,j,kout)
1090 ad_cx=ad_cx+(zeta(istr-1,j,know)+ &
1091 & zeta(istr ,j,know))*ad_zeta(istr-1,j,kout)
1092 ad_zeta(istr,j,know)=ad_zeta(istr,j,know)+ &
1093 & cx*ad_zeta(istr-1,j,kout)
1094 ad_zeta(istr-1,j,kout)=0.0_r8
1095
1096
1097 ad_cff1=ad_cff1+cff*ad_cx
1098 ad_cx=0.0_r8
1099
1100
1101
1102 adfac=0.5_r8*
g*ad_cff1/cff1
1103 grid(ng)%ad_h(istr,j)=
grid(ng)%ad_h(istr,j)+adfac
1104 ad_zeta(istr,j,know)=ad_zeta(istr,j,know)+adfac
1105 ad_cff1=0.0_r8
1106 END IF
1107 END DO
1108
1109
1110
1112 DO j=jstr,jend
1114 cff=dt2d*
grid(ng)%pm(istr,j)
1115 cff1=sqrt(
g*(
grid(ng)%h(istr,j)+ &
1116 & zeta(istr,j,know)))
1117 cx=cff*cff1
1118 cff2=1.0_r8/(1.0_r8+cx)
1119# ifdef MASKING
1120
1121
1122
1123 ad_zeta(istr-1,j,kout)=ad_zeta(istr-1,j,kout)* &
1124 &
grid(ng)%rmask(istr-1,j)
1125# endif
1126
1127
1128
1129
1130
1131
1132 adfac=cff2*ad_zeta(istr-1,j,kout)
1133 ad_zeta(istr-1,j,know)=ad_zeta(istr-1,j,know)+adfac
1134 ad_zeta(istr ,j,kout)=ad_zeta(istr ,j,kout)+cx*adfac
1135 ad_cx=ad_cx+zeta(istr,j,kout)*adfac
1136 ad_cff2=ad_cff2+ &
1137 & (zeta(istr-1,j,know)+ &
1138 & cx*zeta(istr,j,kout))*ad_zeta(istr-1,j,kout)
1139 ad_zeta(istr-1,j,kout)=0.0_r8
1140
1141
1142 ad_cx=ad_cx-cff2*cff2*ad_cff2
1143 ad_cff2=0.0_r8
1144
1145
1146 ad_cff1=ad_cff1+cff*ad_cx
1147 ad_cx=0.0_r8
1148
1149
1150
1151 adfac=0.5_r8*
g*ad_cff1/cff1
1152 grid(ng)%ad_h(istr,j)=
grid(ng)%ad_h(istr,j)+adfac
1153 ad_zeta(istr,j,know)=ad_zeta(istr,j,know)+adfac
1154 ad_cff1=0.0_r8
1155 END IF
1156 END DO
1157
1158
1159
1161 DO j=jstr,jend
1163# ifdef MASKING
1164
1165
1166
1167 ad_zeta(istr-1,j,kout)=ad_zeta(istr-1,j,kout)* &
1168 &
grid(ng)%rmask(istr-1,j)
1169# endif
1170# ifdef ADJUST_BOUNDARY
1172
1173
1175 & ad_zeta_west(j)+ &
1176 & ad_zeta(istr-1,j,kout)
1177 ad_zeta(istr-1,j,kout)=0.0_r8
1178 ELSE
1179
1180
1181 ad_zeta(istr-1,j,kout)=0.0_r8
1182 END IF
1183# else
1184
1185
1186 ad_zeta(istr-1,j,kout)=0.0_r8
1187# endif
1188 END IF
1189 END DO
1190
1191
1192
1194 DO j=jstr,jend
1196# ifdef MASKING
1197
1198
1199
1200 ad_zeta(istr-1,j,kout)=ad_zeta(istr-1,j,kout)* &
1201 &
grid(ng)%rmask(istr-1,j)
1202# endif
1203
1204
1205 ad_zeta(istr ,j,kout)=ad_zeta(istr ,j,kout)+ &
1206 & ad_zeta(istr-1,j,kout)
1207 ad_zeta(istr-1,j,kout)=0.0_r8
1208 END IF
1209 END DO
1210
1211
1212
1214 DO j=jstr,jend
1216# ifdef MASKING
1217
1218
1219
1220 ad_zeta(istr-1,j,kout)=ad_zeta(istr-1,j,kout)* &
1221 &
grid(ng)%rmask(istr-1,j)
1222# endif
1223
1224
1225 ad_zeta(istr ,j,kout)=ad_zeta(istr ,j,kout)+ &
1226 & ad_zeta(istr-1,j,kout)
1227 ad_zeta(istr-1,j,kout)=0.0_r8
1228 END IF
1229 END DO
1230 END IF
1231 END IF
1232
1233 RETURN
type(t_boundary), dimension(:), allocatable boundary
type(t_apply), dimension(:), allocatable lbc_apply
type(t_grid), dimension(:), allocatable grid
type(t_lbc), dimension(:,:,:), allocatable ad_lbc
type(t_domain), dimension(:), allocatable domain
integer, dimension(:), allocatable iic
real(r8), dimension(:), allocatable dcrit
logical, dimension(:,:,:), allocatable lobc
logical, dimension(:), allocatable ewperiodic
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 inorth
real(dp), dimension(:,:), allocatable fsobc_in