ROMS
Loading...
Searching...
No Matches
tl_zetabc.F
Go to the documentation of this file.
1#include "cppdefs.h"
3#ifdef TANGENT
4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group Andrew M. Moore !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This subroutine sets tangent linear lateral boundary conditions for !
13! free-surface. It updates the specified "kout" time index. !
14! !
15! BASIC STATE variables fields needed: zeta !
16! !
17!=======================================================================
18!
19 implicit none
20!
21 PRIVATE
22 PUBLIC :: tl_zetabc, tl_zetabc_tile
23!
24 CONTAINS
25!
26!***********************************************************************
27 SUBROUTINE tl_zetabc (ng, tile, kout)
28!***********************************************************************
29!
30 USE mod_param
31 USE mod_ocean
32 USE mod_stepping
33!
34! Imported variable declarations.
35!
36 integer, intent(in) :: ng, tile, kout
37!
38! Local variable declarations.
39!
40# include "tile.h"
41!
42 CALL tl_zetabc_tile (ng, tile, &
43 & lbi, ubi, lbj, ubj, &
44 & imins, imaxs, jmins, jmaxs, &
45 & krhs(ng), kstp(ng), kout, &
46 & ocean(ng) % zeta, &
47 & ocean(ng) % tl_zeta)
48
49 RETURN
50 END SUBROUTINE tl_zetabc
51!
52!***********************************************************************
53 SUBROUTINE tl_zetabc_tile (ng, tile, &
54 & LBi, UBi, LBj, UBj, &
55 & IminS, ImaxS, JminS, JmaxS, &
56 & krhs, kstp, kout, &
57 & zeta, tl_zeta)
58!***********************************************************************
59!
60 USE mod_param
61 USE mod_boundary
62 USE mod_grid
63 USE mod_ncparam
64 USE mod_scalars
65!
66! Imported variable declarations.
67!
68 integer, intent(in) :: ng, tile
69 integer, intent(in) :: lbi, ubi, lbj, ubj
70 integer, intent(in) :: imins, imaxs, jmins, jmaxs
71 integer, intent(in) :: krhs, kstp, kout
72!
73# ifdef ASSUMED_SHAPE
74 real(r8), intent(in) :: zeta(lbi:,lbj:,:)
75
76 real(r8), intent(inout) :: tl_zeta(lbi:,lbj:,:)
77# else
78 real(r8), intent(in) :: zeta(lbi:ubi,lbj:ubj,3)
79
80 real(r8), intent(inout) :: tl_zeta(lbi:ubi,lbj:ubj,3)
81# endif
82!
83! Local variable declarations.
84!
85 integer :: i, j, know
86
87 real(r8) :: ce, cx
88 real(r8) :: cff, cff1, cff2, dt2d, tau
89
90 real(r8) :: tl_ce, tl_cx
91 real(r8) :: tl_cff1, tl_cff2
92
93 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_grad
94
95# include "set_bounds.h"
96!
97!-----------------------------------------------------------------------
98! Set time-indices
99!-----------------------------------------------------------------------
100!
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!
112!-----------------------------------------------------------------------
113! Lateral boundary conditions at the western edge.
114!-----------------------------------------------------------------------
115!
116 IF (domain(ng)%Western_Edge(tile)) THEN
117!
118! Western edge, implicit upstream radiation condition.
119!
120 IF (tl_lbc(iwest,isfsur,ng)%radiation) THEN
121 IF (iic(ng).ne.0) 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!^
126 tl_grad(istr-1,j)=0.0_r8
127 END DO
128 DO j=jstr,jend
129 IF (lbc_apply(ng)%west(j)) THEN
130# if defined CELERITY_READ && defined FORWARD_READ
131 IF (tl_lbc(iwest,isfsur,ng)%nudging) THEN
132 IF (boundary(ng)%zeta_west_Cx(j).eq.0.0_r8) THEN
133 tau=fsobc_in(ng,iwest)
134 ELSE
135 tau=fsobc_out(ng,iwest)
136 END IF
137 tau=tau*dt2d
138 END IF
139 cx=boundary(ng)%zeta_west_Cx(j)
140# ifdef RADIATION_2D
141 ce=boundary(ng)%zeta_west_Ce(j)
142# else
143 ce=0.0_r8
144# endif
145 cff=boundary(ng)%zeta_west_C2(j)
146# endif
147!^ zeta(Istr-1,j,kout)=(cff*zeta(Istr-1,j,know)+ &
148!^ & Cx *zeta(Istr ,j,kout)- &
149!^ & MAX(Ce,0.0_r8)*grad(Istr-1,j )- &
150!^ & MIN(Ce,0.0_r8)*grad(Istr-1,j+1))/ &
151!^ & (cff+Cx)
152!^
153 tl_zeta(istr-1,j,kout)=(cff*tl_zeta(istr-1,j,know)+ &
154 & cx *tl_zeta(istr ,j,kout)- &
155 & max(ce,0.0_r8)* &
156 & tl_grad(istr-1,j )- &
157 & min(ce,0.0_r8)* &
158 & tl_grad(istr-1,j+1))/ &
159 & (cff+cx)
160
161 IF (tl_lbc(iwest,isfsur,ng)%nudging) THEN
162!^ zeta(Istr-1,j,kout)=zeta(Istr-1,j,kout)+ &
163!^ & tau*(BOUNDARY(ng)%zeta_west(j)- &
164!^ & zeta(Istr-1,j,know))
165!^
166 tl_zeta(istr-1,j,kout)=tl_zeta(istr-1,j,kout)- &
167 & tau*tl_zeta(istr-1,j,know)
168 END IF
169# ifdef MASKING
170!^ zeta(Istr-1,j,kout)=zeta(Istr-1,j,kout)* &
171!^ & GRID(ng)%rmask(Istr-1,j)
172!^
173 tl_zeta(istr-1,j,kout)=tl_zeta(istr-1,j,kout)* &
174 & grid(ng)%rmask(istr-1,j)
175# endif
176 END IF
177 END DO
178 END IF
179!
180! Western edge, explicit Chapman boundary condition.
181!
182 ELSE IF (tl_lbc(iwest,isfsur,ng)%Chapman_explicit) THEN
183 DO j=jstr,jend
184 IF (lbc_apply(ng)%west(j)) THEN
185 cff=dt2d*grid(ng)%pm(istr,j)
186 cff1=sqrt(g*(grid(ng)%h(istr,j)+ &
187 & zeta(istr,j,know)))
188 tl_cff1=0.5_r8*g*(grid(ng)%tl_h(istr,j)+ &
189 & tl_zeta(istr,j,know))/cff1
190 cx=cff*cff1
191 tl_cx=cff*tl_cff1
192!^ zeta(Istr-1,j,kout)=(1.0_r8-Cx)*zeta(Istr-1,j,know)+ &
193!^ & Cx*zeta(Istr,j,know)
194!^
195 tl_zeta(istr-1,j,kout)=(1.0_r8-cx)*tl_zeta(istr-1,j,know)+&
196 & tl_cx*(zeta(istr-1,j,know)+ &
197 & zeta(istr ,j,know))+ &
198 & cx*tl_zeta(istr,j,know)
199# ifdef MASKING
200!^ zeta(Istr-1,j,kout)=zeta(Istr-1,j,kout)* &
201!^ & GRID(ng)%rmask(Istr-1,j)
202!^
203 tl_zeta(istr-1,j,kout)=tl_zeta(istr-1,j,kout)* &
204 & grid(ng)%rmask(istr-1,j)
205# endif
206 END IF
207 END DO
208!
209! Western edge, implicit Chapman boundary condition.
210!
211 ELSE IF (tl_lbc(iwest,isfsur,ng)%Chapman_implicit) THEN
212 DO j=jstr,jend
213 IF (lbc_apply(ng)%west(j)) THEN
214 cff=dt2d*grid(ng)%pm(istr,j)
215 cff1=sqrt(g*(grid(ng)%h(istr,j)+ &
216 & zeta(istr,j,know)))
217 tl_cff1=0.5_r8*g*(grid(ng)%tl_h(istr,j)+ &
218 & tl_zeta(istr,j,know))/cff1
219 cx=cff*cff1
220 tl_cx=cff*tl_cff1
221 cff2=1.0_r8/(1.0_r8+cx)
222 tl_cff2=-cff2*cff2*tl_cx
223!^ zeta(Istr-1,j,kout)=cff2*(zeta(Istr-1,j,know)+ &
224!^ & Cx*zeta(Istr,j,kout))
225!^
226 tl_zeta(istr-1,j,kout)=tl_cff2*(zeta(istr-1,j,know)+ &
227 & cx*zeta(istr,j,kout))+ &
228 & cff2*(tl_zeta(istr-1,j,know)+ &
229 & tl_cx*zeta(istr,j,kout)+ &
230 & cx*tl_zeta(istr,j,kout))
231# ifdef MASKING
232!^ zeta(Istr-1,j,kout)=zeta(Istr-1,j,kout)* &
233!^ & GRID(ng)%rmask(Istr-1,j)
234!^
235 tl_zeta(istr-1,j,kout)=tl_zeta(istr-1,j,kout)* &
236 & grid(ng)%rmask(istr-1,j)
237# endif
238 END IF
239 END DO
240!
241! Western edge, clamped boundary condition.
242!
243 ELSE IF (tl_lbc(iwest,isfsur,ng)%clamped) THEN
244 DO j=jstr,jend
245 IF (lbc_apply(ng)%west(j)) THEN
246!^ zeta(Istr-1,j,kout)=BOUNDARY(ng)%zeta_west(j)
247!^
248# ifdef ADJUST_BOUNDARY
249 IF (lobc(iwest,isfsur,ng)) THEN
250 tl_zeta(istr-1,j,kout)=boundary(ng)%tl_zeta_west(j)
251 ELSE
252 tl_zeta(istr-1,j,kout)=0.0_r8
253 END IF
254# else
255 tl_zeta(istr-1,j,kout)=0.0_r8
256# endif
257# ifdef MASKING
258!^ zeta(Istr-1,j,kout)=zeta(Istr-1,j,kout)* &
259!^ & GRID(ng)%rmask(Istr-1,j)
260!^
261 tl_zeta(istr-1,j,kout)=tl_zeta(istr-1,j,kout)* &
262 & grid(ng)%rmask(istr-1,j)
263# endif
264 END IF
265 END DO
266!
267! Western edge, gradient boundary condition.
268!
269 ELSE IF (tl_lbc(iwest,isfsur,ng)%gradient) THEN
270 DO j=jstr,jend
271 IF (lbc_apply(ng)%west(j)) THEN
272!^ zeta(Istr-1,j,kout)=zeta(Istr,j,kout)
273!^
274 tl_zeta(istr-1,j,kout)=tl_zeta(istr,j,kout)
275# ifdef MASKING
276!^ zeta(Istr-1,j,kout)=zeta(Istr-1,j,kout)* &
277!^ & GRID(ng)%rmask(Istr-1,j)
278!^
279 tl_zeta(istr-1,j,kout)=tl_zeta(istr-1,j,kout)* &
280 & grid(ng)%rmask(istr-1,j)
281# endif
282 END IF
283 END DO
284!
285! Western edge, closed boundary condition.
286!
287 ELSE IF (tl_lbc(iwest,isfsur,ng)%closed) THEN
288 DO j=jstr,jend
289 IF (lbc_apply(ng)%west(j)) THEN
290!^ zeta(Istr-1,j,kout)=zeta(Istr,j,kout)
291!^
292 tl_zeta(istr-1,j,kout)=tl_zeta(istr,j,kout)
293# ifdef MASKING
294!^ zeta(Istr-1,j,kout)=zeta(Istr-1,j,kout)* &
295!^ & GRID(ng)%rmask(Istr-1,j)
296!^
297 tl_zeta(istr-1,j,kout)=tl_zeta(istr-1,j,kout)* &
298 & grid(ng)%rmask(istr-1,j)
299# endif
300 END IF
301 END DO
302 END IF
303 END IF
304!
305!-----------------------------------------------------------------------
306! Lateral boundary conditions at the eastern edge.
307!-----------------------------------------------------------------------
308!
309 IF (domain(ng)%Eastern_Edge(tile)) THEN
310!
311! Eastern edge, implicit upstream radiation condition.
312!
313 IF (tl_lbc(ieast,isfsur,ng)%radiation) THEN
314 IF (iic(ng).ne.0) THEN
315 DO j=jstr,jend+1
316!^ grad(Iend+1,j)=zeta(Iend+1,j ,know)- &
317!^ & zeta(Iend+1,j-1,know)
318!^
319 tl_grad(iend+1,j)=0.0_r8
320 END DO
321 DO j=jstr,jend
322 IF (lbc_apply(ng)%east(j)) THEN
323# if defined CELERITY_READ && defined FORWARD_READ
324 IF (tl_lbc(ieast,isfsur,ng)%nudging) THEN
325 IF (boundary(ng)%zeta_east_Cx(j).eq.0.0_r8) THEN
326 tau=fsobc_in(ieast)
327 ELSE
328 tau=fsobc_out(ieast)
329 END IF
330 tau=tau*dt2d
331 END IF
332 cx=boundary(ng)%zeta_east_Cx(j)
333# ifdef RADIATION_2D
334 ce=boundary(ng)%zeta_east_Ce(j)
335# else
336 ce=0.0_r8
337# endif
338 cff=boundary(ng)%zeta_east_C2(j)
339# endif
340!^ zeta(Iend+1,j,kout)=(cff*zeta(Iend+1,j,know)+ &
341!^ & Cx *zeta(Iend ,j,kout)- &
342!^ & MAX(Ce,0.0_r8)*grad(Iend+1,j )- &
343!^ & MIN(Ce,0.0_r8)*grad(Iend+1,j+1))/ &
344!^ & (cff+Cx)
345!^
346 tl_zeta(iend+1,j,kout)=(cff*tl_zeta(iend+1,j,know)+ &
347 & cx *tl_zeta(iend ,j,kout)- &
348 & max(ce,0.0_r8)* &
349 & tl_grad(iend+1,j )- &
350 & min(ce,0.0_r8)* &
351 & tl_grad(iend+1,j+1))/ &
352 & (cff+cx)
353
354 IF (tl_lbc(ieast,isfsur,ng)%nudging) THEN
355!^ zeta(Iend+1,j,kout)=zeta(Iend+1,j,kout)+ &
356!^ & tau*(BOUNDARY(ng)%zeta_east(j)- &
357!^ & zeta(Iend+1,j,know))
358!^
359 tl_zeta(iend+1,j,kout)=tl_zeta(iend+1,j,kout)- &
360 & tau*tl_zeta(iend+1,j,know)
361 END IF
362# ifdef MASKING
363!^ zeta(Iend+1,j,kout)=zeta(Iend+1,j,kout)* &
364!^ & GRID(ng)%rmask(Iend+1,j)
365!^
366 tl_zeta(iend+1,j,kout)=tl_zeta(iend+1,j,kout)* &
367 & grid(ng)%rmask(iend+1,j)
368# endif
369 END IF
370 END DO
371 END IF
372!
373! Eastern edge, explicit Chapman boundary condition.
374!
375 ELSE IF (tl_lbc(ieast,isfsur,ng)%Chapman_explicit) THEN
376 DO j=jstr,jend
377 IF (lbc_apply(ng)%east(j)) THEN
378 cff=dt2d*grid(ng)%pm(iend,j)
379 cff1=sqrt(g*(grid(ng)%h(iend,j)+ &
380 & zeta(iend,j,know)))
381 tl_cff1=0.5_r8*g*(grid(ng)%tl_h(iend,j)+ &
382 & tl_zeta(iend,j,know))/cff1
383 cx=cff*cff1
384 tl_cx=cff*tl_cff1
385!^ zeta(Iend+1,j,kout)=(1.0_r8-Cx)*zeta(Iend+1,j,know)+ &
386!^ & Cx*zeta(Iend,j,know)
387!^
388 tl_zeta(iend+1,j,kout)=(1.0_r8-cx)*tl_zeta(iend+1,j,know)+&
389 & tl_cx*(zeta(iend+1,j,know)+ &
390 & zeta(iend ,j,know))+ &
391 & cx*tl_zeta(iend,j,know)
392# ifdef MASKING
393!^ zeta(Iend+1,j,kout)=zeta(Iend+1,j,kout)* &
394!^ & GRID(ng)%rmask(Iend+1,j)
395!^
396 tl_zeta(iend+1,j,kout)=tl_zeta(iend+1,j,kout)* &
397 & grid(ng)%rmask(iend+1,j)
398# endif
399 END IF
400 END DO
401!
402! Eastern edge, implicit Chapman boundary condition.
403!
404 ELSE IF (tl_lbc(ieast,isfsur,ng)%Chapman_implicit) THEN
405 DO j=jstr,jend
406 IF (lbc_apply(ng)%east(j)) THEN
407 cff=dt2d*grid(ng)%pm(iend,j)
408 cff1=sqrt(g*(grid(ng)%h(iend,j)+ &
409 & zeta(iend,j,know)))
410 tl_cff1=0.5_r8*g*(grid(ng)%tl_h(iend,j)+ &
411 & tl_zeta(iend,j,know))/cff1
412 cx=cff*cff1
413 tl_cx=cff*tl_cff1
414 cff2=1.0_r8/(1.0_r8+cx)
415 tl_cff2=-cff2*cff2*tl_cx
416!^ zeta(Iend+1,j,kout)=cff2*(zeta(Iend+1,j,know)+ &
417!^ & Cx*zeta(Iend,j,kout))
418!^
419 tl_zeta(iend+1,j,kout)=tl_cff2*(zeta(iend+1,j,know)+ &
420 & cx*zeta(iend,j,kout))+ &
421 & cff2*(tl_zeta(iend+1,j,know)+ &
422 & tl_cx*zeta(iend,j,kout)+ &
423 & cx*tl_zeta(iend,j,kout))
424# ifdef MASKING
425!^ zeta(Iend+1,j,kout)=zeta(Iend+1,j,kout)* &
426!^ & GRID(ng)%rmask(Iend+1,j)
427!^
428 tl_zeta(iend+1,j,kout)=tl_zeta(iend+1,j,kout)* &
429 & grid(ng)%rmask(iend+1,j)
430# endif
431 END IF
432 END DO
433!
434! Eastern edge, clamped boundary condition.
435!
436 ELSE IF (tl_lbc(ieast,isfsur,ng)%clamped) THEN
437 DO j=jstr,jend
438 IF (lbc_apply(ng)%east(j)) THEN
439!^ zeta(Iend+1,j,kout)=BOUNDARY(ng)%zeta_east(j)
440!^
441# ifdef ADJUST_BOUNDARY
442 IF (lobc(ieast,isfsur,ng)) THEN
443 tl_zeta(iend+1,j,kout)=boundary(ng)%tl_zeta_east(j)
444 ELSE
445 tl_zeta(iend+1,j,kout)=0.0_r8
446 END IF
447# else
448 tl_zeta(iend+1,j,kout)=0.0_r8
449# endif
450# ifdef MASKING
451!^ zeta(Iend+1,j,kout)=zeta(Iend+1,j,kout)* &
452!^ & GRID(ng)%rmask(Iend+1,j)
453!^
454 tl_zeta(iend+1,j,kout)=tl_zeta(iend+1,j,kout)* &
455 & grid(ng)%rmask(iend+1,j)
456# endif
457 END IF
458 END DO
459!
460! Eastern edge, gradient boundary condition.
461!
462 ELSE IF (tl_lbc(ieast,isfsur,ng)%gradient) THEN
463 DO j=jstr,jend
464 IF (lbc_apply(ng)%east(j)) THEN
465!^ zeta(Iend+1,j,kout)=zeta(Iend,j,kout)
466!^
467 tl_zeta(iend+1,j,kout)=tl_zeta(iend,j,kout)
468# ifdef MASKING
469!^ zeta(Iend+1,j,kout)=zeta(Iend+1,j,kout)* &
470!^ & GRID(ng)%rmask(Iend+1,j)
471!^
472 tl_zeta(iend+1,j,kout)=tl_zeta(iend+1,j,kout)* &
473 & grid(ng)%rmask(iend+1,j)
474# endif
475 END IF
476 END DO
477!
478! Eastern edge, closed boundary condition.
479!
480 ELSE IF (tl_lbc(ieast,isfsur,ng)%closed) THEN
481 DO j=jstr,jend
482 IF (lbc_apply(ng)%east(j)) THEN
483!^ zeta(Iend+1,j,kout)=zeta(Iend,j,kout)
484!^
485 tl_zeta(iend+1,j,kout)=tl_zeta(iend,j,kout)
486# ifdef MASKING
487!^ zeta(Iend+1,j,kout)=zeta(Iend+1,j,kout)* &
488!^ & GRID(ng)%rmask(Iend+1,j)
489!^
490 tl_zeta(iend+1,j,kout)=tl_zeta(iend+1,j,kout)* &
491 & grid(ng)%rmask(iend+1,j)
492# endif
493 END IF
494 END DO
495 END IF
496 END IF
497!
498!-----------------------------------------------------------------------
499! Lateral boundary conditions at the southern edge.
500!-----------------------------------------------------------------------
501!
502 IF (domain(ng)%Southern_Edge(tile)) THEN
503!
504! Southern edge, implicit upstream radiation condition.
505!
506 IF (tl_lbc(isouth,isfsur,ng)%radiation) THEN
507 IF (iic(ng).ne.0) THEN
508 DO i=istr,iend+1
509!^ grad(i,Jstr)=zeta(i ,Jstr,know)- &
510!^ & zeta(i-1,Jstr,know)
511!^
512 tl_grad(i,jstr)=0.0_r8
513 END DO
514 DO i=istr,iend
515 IF (lbc_apply(ng)%south(i)) THEN
516# if defined CELERITY_READ && defined FORWARD_READ
517 IF (tl_lbc(isouth,isfsur,ng)%nudging) THEN
518 IF (boundary(ng)%zeta_south_Ce(i).eq.0.0_r8) THEN
519 tau=fsobc_in(ng,isouth)
520 ELSE
521 tau=fsobc_out(ng,isouth)
522 END IF
523 tau=tau*dt2d
524 END IF
525# ifdef RADIATION_2D
526 cx=boundary(ng)%zeta_south_Cx(i)
527# else
528 cx=0.0_r8
529# endif
530 ce=boundary(ng)%zeta_south_Ce(i)
531 cff=boundary(ng)%zeta_south_C2(i)
532# endif
533!^ zeta(i,Jstr-1,kout)=(cff*zeta(i,Jstr-1,know)+ &
534!^ & Ce *zeta(i,Jstr ,kout)- &
535!^ & MAX(Cx,0.0_r8)*grad(i ,Jstr)- &
536!^ & MIN(Cx,0.0_r8)*grad(i+1,Jstr))/ &
537!^ & (cff+Ce)
538!^
539 tl_zeta(i,jstr-1,kout)=(cff*tl_zeta(i,jstr-1,know)+ &
540 & ce *tl_zeta(i,jstr ,kout)- &
541 & max(cx,0.0_r8)* &
542 & tl_grad(i ,jstr-1)- &
543 & min(cx,0.0_r8)* &
544 & tl_grad(i+1,jstr-1))/ &
545 & (cff+ce)
546
547 IF (tl_lbc(isouth,isfsur,ng)%nudging) THEN
548!^ zeta(i,Jstr-1,kout)=zeta(i,Jstr-1,kout)+ &
549!^ & tau*(BOUNDARY(ng)%zeta_south(i)- &
550!^ & zeta(i,Jstr-1,know))
551!^
552 tl_zeta(i,jstr-1,kout)=tl_zeta(i,jstr-1,kout)- &
553 & tau*tl_zeta(i,jstr-1,know)
554 END IF
555# ifdef MASKING
556!^ zeta(i,Jstr-1,kout)=zeta(i,Jstr-1,kout)* &
557!^ & GRID(ng)%rmask(i,Jstr-1)
558!^
559 tl_zeta(i,jstr-1,kout)=tl_zeta(i,jstr-1,kout)* &
560 & grid(ng)%rmask(i,jstr-1)
561# endif
562 END IF
563 END DO
564 END IF
565!
566! Southern edge, explicit Chapman boundary condition.
567!
568 ELSE IF (tl_lbc(isouth,isfsur,ng)%Chapman_explicit) THEN
569 DO i=istr,iend
570 IF (lbc_apply(ng)%south(i)) THEN
571 cff=dt2d*grid(ng)%pn(i,jstr)
572 cff1=sqrt(g*(grid(ng)%h(i,jstr)+ &
573 & zeta(i,jstr,know)))
574 tl_cff1=0.5_r8*g*(grid(ng)%tl_h(i,jstr)+ &
575 & tl_zeta(i,jstr,know))/cff1
576 ce=cff*cff1
577 tl_ce=cff*tl_cff1
578!^ zeta(i,Jstr-1,kout)=(1.0_r8-Ce)*zeta(i,Jstr-1,know)+ &
579!^ & Ce*zeta(i,Jstr,know)
580!^
581 tl_zeta(i,jstr-1,kout)=(1.0_r8-ce)*tl_zeta(i,jstr-1,know)+&
582 & tl_ce*(zeta(i,jstr-1,know)+ &
583 & zeta(i,jstr ,know))+ &
584 & ce*tl_zeta(i,jstr,know)
585# ifdef MASKING
586!^ zeta(i,Jstr-1,kout)=zeta(i,Jstr-1,kout)* &
587!^ & GRID(ng)%rmask(i,Jstr-1)
588!^
589 tl_zeta(i,jstr-1,kout)=tl_zeta(i,jstr-1,kout)* &
590 & grid(ng)%rmask(i,jstr-1)
591# endif
592 END IF
593 END DO
594!
595! Southern edge, implicit Chapman boundary condition.
596!
597 ELSE IF (tl_lbc(isouth,isfsur,ng)%Chapman_implicit) THEN
598 DO i=istr,iend
599 IF (lbc_apply(ng)%south(i)) THEN
600 cff=dt2d*grid(ng)%pn(i,jstr)
601 cff1=sqrt(g*(grid(ng)%h(i,jstr)+ &
602 & zeta(i,jstr,know)))
603 tl_cff1=0.5_r8*g*(grid(ng)%tl_h(i,jstr)+ &
604 & tl_zeta(i,jstr,know))/cff1
605 ce=cff*cff1
606 tl_ce=cff*tl_cff1
607 cff2=1.0_r8/(1.0_r8+ce)
608 tl_cff2=-cff2*cff2*tl_ce
609!^ zeta(i,Jstr-1,kout)=cff2*(zeta(i,Jstr-1,know)+ &
610!^ & Ce*zeta(i,Jstr,kout))
611!^
612 tl_zeta(i,jstr-1,kout)=tl_cff2*(zeta(i,jstr-1,know)+ &
613 & ce*zeta(i,jstr,kout))+ &
614 & cff2*(tl_zeta(i,jstr-1,know)+ &
615 & tl_ce*zeta(i,jstr,kout)+ &
616 & ce*tl_zeta(i,jstr,kout))
617# ifdef MASKING
618!^ zeta(i,Jstr-1,kout)=zeta(i,Jstr-1,kout)* &
619!^ & GRID(ng)%rmask(i,Jstr-1)
620!^
621 tl_zeta(i,jstr-1,kout)=tl_zeta(i,jstr-1,kout)* &
622 & grid(ng)%rmask(i,jstr-1)
623# endif
624 END IF
625 END DO
626!
627! Southern edge, clamped boundary condition.
628!
629 ELSE IF (tl_lbc(isouth,isfsur,ng)%clamped) THEN
630 DO i=istr,iend
631 IF (lbc_apply(ng)%south(i)) THEN
632!^ zeta(i,Jstr-1,kout)=BOUNDARY(ng)%zeta_south(i)
633!^
634# ifdef ADJUST_BOUNDARY
635 IF (lobc(isouth,isfsur,ng)) THEN
636 tl_zeta(i,jstr-1,kout)=boundary(ng)%tl_zeta_south(i)
637 ELSE
638 tl_zeta(i,jstr-1,kout)=0.0_r8
639 END IF
640# else
641 tl_zeta(i,jstr-1,kout)=0.0_r8
642# endif
643# ifdef MASKING
644!^ zeta(i,Jstr-1,kout)=zeta(i,Jstr-1,kout)* &
645!^ & GRID(ng)%rmask(i,Jstr-1)
646!^
647 tl_zeta(i,jstr-1,kout)=tl_zeta(i,jstr-1,kout)* &
648 & grid(ng)%rmask(i,jstr-1)
649# endif
650 END IF
651 END DO
652!
653! Southern edge, gradient boundary condition.
654!
655 ELSE IF (tl_lbc(isouth,isfsur,ng)%gradient) THEN
656 DO i=istr,iend
657 IF (lbc_apply(ng)%south(i)) THEN
658!^ zeta(i,Jstr-1,kout)=zeta(i,Jstr,kout)
659!^
660 tl_zeta(i,jstr-1,kout)=tl_zeta(i,jstr,kout)
661# ifdef MASKING
662!^ zeta(i,Jstr-1,kout)=zeta(i,Jstr-1,kout)* &
663!^ & GRID(ng)%rmask(i,Jstr-1)
664!^
665 tl_zeta(i,jstr-1,kout)=tl_zeta(i,jstr-1,kout)* &
666 & grid(ng)%rmask(i,jstr-1)
667# endif
668 END IF
669 END DO
670!
671! Southern edge, closed boundary condition.
672!
673 ELSE IF (tl_lbc(isouth,isfsur,ng)%closed) THEN
674 DO i=istr,iend
675 IF (lbc_apply(ng)%south(i)) THEN
676!^ zeta(i,Jstr-1,kout)=zeta(i,Jstr,kout)
677!^
678 tl_zeta(i,jstr-1,kout)=tl_zeta(i,jstr,kout)
679# ifdef MASKING
680!^ zeta(i,Jstr-1,kout)=zeta(i,Jstr-1,kout)* &
681!^ & GRID(ng)%rmask(i,Jstr-1)
682!^
683 tl_zeta(i,jstr-1,kout)=tl_zeta(i,jstr-1,kout)* &
684 & grid(ng)%rmask(i,jstr-1)
685# endif
686 END IF
687 END DO
688 END IF
689 END IF
690!
691!-----------------------------------------------------------------------
692! Lateral boundary conditions at the northern edge.
693!-----------------------------------------------------------------------
694!
695 IF (domain(ng)%Northern_Edge(tile)) THEN
696!
697! Northern edge, implicit upstream radiation condition.
698!
699 IF (tl_lbc(inorth,isfsur,ng)%radiation) THEN
700 IF (iic(ng).ne.0) THEN
701 DO i=istr,iend+1
702!^ grad(i,Jend+1)=zeta(i ,Jend+1,know)- &
703!^ & zeta(i-1,Jend+1,know)
704!^
705 tl_grad(i,jend+1)=0.0_r8
706 END DO
707 DO i=istr,iend
708 IF (lbc_apply(ng)%north(i)) THEN
709# if defined CELERITY_READ && defined FORWARD_READ
710 IF (tl_lbc(inorth,isfsur,ng)%nudging) THEN
711 IF (boundary(ng)%zeta_north_Ce(i).eq.0.0_r8) THEN
712 tau=fsobc_in(ng,inorth)
713 ELSE
714 tau=fsobc_out(ng,inorth)
715 END IF
716 tau=tau*dt2d
717 END IF
718# ifdef RADIATION_2D
719 cx=boundary(ng)%zeta_north_Cx(i)
720# else
721 cx=0.0_r8
722# endif
723 ce=boundary(ng)%zeta_north_Ce(i)
724 cff=boundary(ng)%zeta_north_C2(i)
725# endif
726!^ zeta(i,Jend+1,kout)=(cff*zeta(i,Jend+1,know)+ &
727!^ & Ce *zeta(i,Jend ,kout)- &
728!^ & MAX(Cx,0.0_r8)*grad(i ,Jend+1)- &
729!^ & MIN(Cx,0.0_r8)*grad(i+1,Jend+1))/ &
730!^ & (cff+Ce)
731!^
732 tl_zeta(i,jend+1,kout)=(cff*tl_zeta(i,jend+1,know)+ &
733 & ce *tl_zeta(i,jend ,kout)- &
734 & max(cx,0.0_r8)* &
735 & tl_grad(i ,jend+1)- &
736 & min(cx,0.0_r8)* &
737 & tl_grad(i+1,jend+1))/ &
738 & (cff+ce)
739
740 IF (tl_lbc(inorth,isfsur,ng)%nudging) THEN
741!^ zeta(i,Jend+1,kout)=zeta(i,Jend+1,kout)+ &
742!^ & tau*(BOUNDARY(ng)%zeta_north(i)- &
743!^ & zeta(i,Jend+1,know))
744!^
745 tl_zeta(i,jend+1,kout)=tl_zeta(i,jend+1,kout)- &
746 & tau*tl_zeta(i,jend+1,know)
747 END IF
748# ifdef MASKING
749!^ zeta(i,Jend+1,kout)=zeta(i,Jend+1,kout)* &
750!^ & GRID(ng)%rmask(i,Jend+1)
751!^
752 tl_zeta(i,jend+1,kout)=tl_zeta(i,jend+1,kout)* &
753 & grid(ng)%rmask(i,jend+1)
754# endif
755 END IF
756 END DO
757 END IF
758!
759! Northern edge, explicit Chapman boundary condition.
760!
761 ELSE IF (tl_lbc(inorth,isfsur,ng)%Chapman_explicit) THEN
762 DO i=istr,iend
763 IF (lbc_apply(ng)%north(i)) THEN
764 cff=dt2d*grid(ng)%pn(i,jend)
765 cff1=sqrt(g*(grid(ng)%h(i,jend)+ &
766 & zeta(i,jend,know)))
767 tl_cff1=0.5_r8*g*(grid(ng)%tl_h(i,jend)+ &
768 & tl_zeta(i,jend,know))/cff1
769 ce=cff*cff1
770 tl_ce=cff*tl_cff1
771!^ zeta(i,Jend+1,kout)=(1.0_r8-Ce)*zeta(i,Jend+1,know)+ &
772!^ & Ce*zeta(i,Jend,know)
773!^
774 tl_zeta(i,jend+1,kout)=(1.0_r8-ce)*tl_zeta(i,jend+1,know)+&
775 & tl_ce*(zeta(i,jend+1,know)+ &
776 & zeta(i,jend ,know))+ &
777 & ce*tl_zeta(i,jend,know)
778# ifdef MASKING
779!^ zeta(i,Jend+1,kout)=zeta(i,Jend+1,kout)* &
780!^ & GRID(ng)%rmask(i,Jend+1)
781!^
782 tl_zeta(i,jend+1,kout)=tl_zeta(i,jend+1,kout)* &
783 & grid(ng)%rmask(i,jend+1)
784# endif
785 END IF
786 END DO
787!
788! Northern edge, implicit Chapman boundary condition.
789!
790 ELSE IF (tl_lbc(inorth,isfsur,ng)%Chapman_implicit) THEN
791 DO i=istr,iend
792 IF (lbc_apply(ng)%north(i)) THEN
793 cff=dt2d*grid(ng)%pn(i,jend)
794 cff1=sqrt(g*(grid(ng)%h(i,jend)+ &
795 & zeta(i,jend,know)))
796 tl_cff1=0.5_r8*g*(grid(ng)%tl_h(i,jend)+ &
797 & tl_zeta(i,jend,know))/cff1
798 ce=cff*cff1
799 tl_ce=cff*tl_cff1
800 cff2=1.0_r8/(1.0_r8+ce)
801 tl_cff2=-cff2*cff2*tl_ce
802!^ zeta(i,Jend+1,kout)=cff2*(zeta(i,Jend+1,know)+ &
803!^ & Ce*zeta(i,Jend,kout))
804!^
805 tl_zeta(i,jend+1,kout)=tl_cff2*(zeta(i,jend+1,know)+ &
806 & ce*zeta(i,jend,kout))+ &
807 & cff2*(tl_zeta(i,jend+1,know)+ &
808 & tl_ce*zeta(i,jend,kout)+ &
809 & ce*tl_zeta(i,jend,kout))
810# ifdef MASKING
811!^ zeta(i,Jend+1,kout)=zeta(i,Jend+1,kout)* &
812!^ & GRID(ng)%rmask(i,Jend+1)
813!^
814 tl_zeta(i,jend+1,kout)=tl_zeta(i,jend+1,kout)* &
815 & grid(ng)%rmask(i,jend+1)
816# endif
817 END IF
818 END DO
819!
820! Northern edge, clamped boundary condition.
821!
822 ELSE IF (tl_lbc(inorth,isfsur,ng)%clamped) THEN
823 DO i=istr,iend
824 IF (lbc_apply(ng)%north(i)) THEN
825!^ zeta(i,Jend+1,kout)=BOUNDARY(ng)%zeta_north(i)
826!^
827# ifdef ADJUST_BOUNDARY
828 IF (lobc(inorth,isfsur,ng)) THEN
829 tl_zeta(i,jend+1,kout)=boundary(ng)%tl_zeta_north(i)
830 ELSE
831 tl_zeta(i,jend+1,kout)=0.0_r8
832 END IF
833# else
834 tl_zeta(i,jend+1,kout)=0.0_r8
835# endif
836# ifdef MASKING
837!^ zeta(i,Jend+1,kout)=zeta(i,Jend+1,kout)* &
838!^ & GRID(ng)%rmask(i,Jend+1)
839!^
840 tl_zeta(i,jend+1,kout)=tl_zeta(i,jend+1,kout)* &
841 & grid(ng)%rmask(i,jend+1)
842# endif
843 END IF
844 END DO
845!
846! Northern edge, gradient boundary condition.
847!
848 ELSE IF (tl_lbc(inorth,isfsur,ng)%gradient) THEN
849 DO i=istr,iend
850 IF (lbc_apply(ng)%north(i)) THEN
851!^ zeta(i,Jend+1,kout)=zeta(i,Jend,kout)
852!^
853 tl_zeta(i,jend+1,kout)=tl_zeta(i,jend,kout)
854# ifdef MASKING
855!^ zeta(i,Jend+1,kout)=zeta(i,Jend+1,kout)* &
856!^ & GRID(ng)%rmask(i,Jend+1)
857!^
858 tl_zeta(i,jend+1,kout)=tl_zeta(i,jend+1,kout)* &
859 & grid(ng)%rmask(i,jend+1)
860# endif
861 END IF
862 END DO
863!
864! Northern edge, closed boundary condition.
865!
866 ELSE IF (tl_lbc(inorth,isfsur,ng)%closed) THEN
867 DO i=istr,iend
868 IF (lbc_apply(ng)%north(i)) THEN
869!^ zeta(i,Jend+1,kout)=zeta(i,Jend,kout)
870!^
871 tl_zeta(i,jend+1,kout)=tl_zeta(i,jend,kout)
872# ifdef MASKING
873!^ zeta(i,Jend+1,kout)=zeta(i,Jend+1,kout)* &
874!^ & GRID(ng)%rmask(i,Jend+1)
875!^
876 tl_zeta(i,jend+1,kout)=tl_zeta(i,jend+1,kout)* &
877 & grid(ng)%rmask(i,jend+1)
878# endif
879 END IF
880 END DO
881 END IF
882 END IF
883!
884!-----------------------------------------------------------------------
885! Boundary corners.
886!-----------------------------------------------------------------------
887!
888 IF (.not.(ewperiodic(ng).or.nsperiodic(ng))) THEN
889 IF (domain(ng)%SouthWest_Corner(tile)) THEN
890 IF (lbc_apply(ng)%south(istr-1).and. &
891 & lbc_apply(ng)%west (jstr-1)) THEN
892!^ zeta(Istr-1,Jstr-1,kout)=0.5_r8*(zeta(Istr ,Jstr-1,kout)+ &
893!^ & zeta(Istr-1,Jstr ,kout))
894!^
895 tl_zeta(istr-1,jstr-1,kout)=0.5_r8* &
896 & (tl_zeta(istr ,jstr-1,kout)+ &
897 & tl_zeta(istr-1,jstr ,kout))
898 END IF
899 END IF
900 IF (domain(ng)%SouthEast_Corner(tile)) THEN
901 IF (lbc_apply(ng)%south(iend+1).and. &
902 & lbc_apply(ng)%east (jstr-1)) THEN
903!^ zeta(Iend+1,Jstr-1,kout)=0.5_r8*(zeta(Iend ,Jstr-1,kout)+ &
904!^ & zeta(Iend+1,Jstr ,kout))
905!^
906 tl_zeta(iend+1,jstr-1,kout)=0.5_r8* &
907 & (tl_zeta(iend ,jstr-1,kout)+ &
908 & tl_zeta(iend+1,jstr ,kout))
909 END IF
910 END IF
911 IF (domain(ng)%NorthWest_Corner(tile)) THEN
912 IF (lbc_apply(ng)%north(istr-1).and. &
913 & lbc_apply(ng)%west (jend+1)) THEN
914!^ zeta(Istr-1,Jend+1,kout)=0.5_r8*(zeta(Istr-1,Jend ,kout)+ &
915!^ & zeta(Istr ,Jend+1,kout))
916!^
917 tl_zeta(istr-1,jend+1,kout)=0.5_r8* &
918 & (tl_zeta(istr-1,jend ,kout)+ &
919 & tl_zeta(istr ,jend+1,kout))
920 END IF
921 END IF
922 IF (domain(ng)%NorthEast_Corner(tile)) THEN
923 IF (lbc_apply(ng)%north(iend+1).and. &
924 & lbc_apply(ng)%east (jend+1)) THEN
925!^ zeta(Iend+1,Jend+1,kout)=0.5_r8*(zeta(Iend+1,Jend ,kout)+ &
926!^ & zeta(Iend ,Jend+1,kout))
927!^
928 tl_zeta(iend+1,jend+1,kout)=0.5_r8* &
929 & (tl_zeta(iend+1,jend ,kout)+ &
930 & tl_zeta(iend ,jend+1,kout))
931 END IF
932 END IF
933 END IF
934
935# if defined WET_DRY
936!
937!-----------------------------------------------------------------------
938! Ensure that water level on boundary cells is above bed elevation.
939!-----------------------------------------------------------------------
940!
941 cff=dcrit(ng)-eps
942 IF (.not.ewperiodic(ng)) THEN
943 IF (domain(ng)%Western_Edge(tile)) THEN
944 DO j=jstr,jend
945 IF (lbc_apply(ng)%west(j)) THEN
946 IF (zeta(istr-1,j,kout).le. &
947 & (dcrit(ng)-grid(ng)%h(istr-1,j))) THEN
948!^ zeta(Istr-1,j,kout)=cff-GRID(ng)%h(Istr-1,j)
949!^
950 tl_zeta(istr-1,j,kout)=-grid(ng)%tl_h(istr-1,j)
951 END IF
952 END IF
953 END DO
954 END IF
955 IF (domain(ng)%Eastern_Edge(tile)) THEN
956 DO j=jstr,jend
957 IF (lbc_apply(ng)%east(j)) THEN
958 IF (zeta(iend+1,j,kout).le. &
959 & (dcrit(ng)-grid(ng)%h(iend+1,j))) THEN
960!^ zeta(Iend+1,j,kout)=cff-GRID(ng)%h(Iend+1,j)
961!^
962 tl_zeta(iend+1,j,kout)=-grid(ng)%tl_h(iend+1,j)
963 END IF
964 END IF
965 END DO
966 END IF
967 END IF
968
969 IF (.not.nsperiodic(ng)) THEN
970 IF (domain(ng)%Southern_Edge(tile)) THEN
971 DO i=istr,iend
972 IF (lbc_apply(ng)%south(i)) THEN
973 IF (zeta(i,jstr-1,kout).le. &
974 & (dcrit(ng)-grid(ng)%h(i,jstr-1))) THEN
975!^ zeta(i,Jstr-1,kout)=cff-GRID(ng)%h(i,Jstr-1)
976!^
977 tl_zeta(i,jstr-1,kout)=-grid(ng)%tl_h(i,jstr-1)
978 END IF
979 END IF
980 END DO
981 END IF
982 IF (domain(ng)%Northern_Edge(tile)) THEN
983 DO i=istr,iend
984 IF (lbc_apply(ng)%north(i)) THEN
985 IF (zeta(i,jend+1,kout).le. &
986 & (dcrit(ng)-grid(ng)%h(i,jend+1))) THEN
987!^ zeta(i,Jend+1,kout)=cff-GRID(ng)%h(i,Jend+1)
988!^
989 tl_zeta(i,jend+1,kout)=-grid(ng)%tl_h(i,jend+1)
990 END IF
991 END IF
992 END DO
993 END IF
994 END IF
995
996 IF (.not.(ewperiodic(ng).or.nsperiodic(ng))) THEN
997 IF (domain(ng)%SouthWest_Corner(tile)) THEN
998 IF (lbc_apply(ng)%south(istr-1).and. &
999 & lbc_apply(ng)%west (jstr-1)) THEN
1000 IF (zeta(istr-1,jstr-1,kout).le. &
1001 & (dcrit(ng)-grid(ng)%h(istr-1,jstr-1))) THEN
1002!^ zeta(Istr-1,Jstr-1,kout)=cff-GRID(ng)%h(Istr-1,Jstr-1)
1003!^
1004 tl_zeta(istr-1,jstr-1,kout)=-grid(ng)%tl_h(istr-1,jstr-1)
1005 END IF
1006 END IF
1007 END IF
1008 IF (domain(ng)%SouthEast_Corner(tile)) THEN
1009 IF (lbc_apply(ng)%south(iend+1).and. &
1010 & lbc_apply(ng)%east (jstr-1)) THEN
1011 IF (zeta(iend+1,jstr-1,kout).le. &
1012 & (dcrit(ng)-grid(ng)%h(iend+1,jstr-1))) THEN
1013!^ zeta(Iend+1,Jstr-1,kout)=cff-GRID(ng)%h(Iend+1,Jstr-1)
1014!^
1015 tl_zeta(iend+1,jstr-1,kout)=-grid(ng)%tl_h(iend+1,jstr-1)
1016 END IF
1017 END IF
1018 END IF
1019 IF (domain(ng)%NorthWest_Corner(tile)) THEN
1020 IF (lbc_apply(ng)%north(istr-1).and. &
1021 & lbc_apply(ng)%west (jend+1)) THEN
1022 IF (zeta(istr-1,jend+1,kout).le. &
1023 & (dcrit(ng)-grid(ng)%h(istr-1,jend+1))) THEN
1024!^ zeta(Istr-1,Jend+1,kout)=cff-GRID(ng)%h(Istr-1,Jend+1)
1025!^
1026 tl_zeta(istr-1,jend+1,kout)=-grid(ng)%tl_h(istr-1,jend+1)
1027 END IF
1028 END IF
1029 END IF
1030 IF (domain(ng)%NorthEast_Corner(tile)) THEN
1031 IF (lbc_apply(ng)%north(iend+1).and. &
1032 & lbc_apply(ng)%east (jend+1)) THEN
1033 IF (zeta(iend+1,jend+1,kout).le. &
1034 & (dcrit(ng)-grid(ng)%h(iend+1,jend+1))) THEN
1035!^ zeta(Iend+1,Jend+1,kout)=cff-GRID(ng)%h(Iend+1,Jend+1)
1036!^
1037 tl_zeta(iend+1,jend+1,kout)=-grid(ng)%tl_h(iend+1,jend+1)
1038 END IF
1039 END IF
1040 END IF
1041 END IF
1042# endif
1043
1044 RETURN
1045 END SUBROUTINE tl_zetabc_tile
1046#endif
1047 END MODULE tl_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 tl_lbc
Definition mod_param.F:379
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer, dimension(:), allocatable iic
real(r8), dimension(:), allocatable dcrit
logical, dimension(:,:,:), allocatable lobc
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, public tl_zetabc(ng, tile, kout)
Definition tl_zetabc.F:28
subroutine, public tl_zetabc_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, krhs, kstp, kout, zeta, tl_zeta)
Definition tl_zetabc.F:58