ROMS
Loading...
Searching...
No Matches
rp_v2dbc_im.F
Go to the documentation of this file.
1#include "cppdefs.h"
3#ifdef TL_IOMS
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 representers tangent linear lateral boundary !
13! conditions for vertically integrated V-velocity. It updates the !
14! specified "kout" index. !
15! !
16! BASIC STATE variables needed: vbar !
17! !
18!=======================================================================
19!
20 implicit none
21!
22 PRIVATE
23 PUBLIC :: rp_v2dbc, rp_v2dbc_tile
24!
25 CONTAINS
26!
27!***********************************************************************
28 SUBROUTINE rp_v2dbc (ng, tile, kout)
29!***********************************************************************
30!
31 USE mod_param
32 USE mod_ocean
33 USE mod_stepping
34!
35! Imported variable declarations.
36!
37 integer, intent(in) :: ng, tile, kout
38!
39! Local variable declarations.
40!
41# include "tile.h"
42!
43 CALL rp_v2dbc_tile (ng, tile, &
44 & lbi, ubi, lbj, ubj, &
45 & imins, imaxs, jmins, jmaxs, &
46 & krhs(ng), kstp(ng), kout, &
47 & ocean(ng) % ubar, &
48 & ocean(ng) % vbar, &
49 & ocean(ng) % zeta, &
50 & ocean(ng) % tl_ubar, &
51 & ocean(ng) % tl_vbar, &
52 & ocean(ng) % tl_zeta)
53
54 RETURN
55 END SUBROUTINE rp_v2dbc
56!
57!***********************************************************************
58 SUBROUTINE rp_v2dbc_tile (ng, tile, &
59 & LBi, UBi, LBj, UBj, &
60 & IminS, ImaxS, JminS, JmaxS, &
61 & krhs, kstp, kout, &
62 & ubar, vbar, zeta, &
63 & tl_ubar, tl_vbar, tl_zeta)
64!***********************************************************************
65!
66 USE mod_param
67 USE mod_boundary
68 USE mod_clima
69 USE mod_forces
70 USE mod_grid
71 USE mod_ncparam
72 USE mod_scalars
73!
74! Imported variable declarations.
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_ubar(lbi:,lbj:,:)
86 real(r8), intent(in) :: tl_zeta(lbi:,lbj:,:)
87
88 real(r8), intent(inout) :: tl_vbar(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_ubar(lbi:ubi,lbj:ubj,:)
94 real(r8), intent(in) :: tl_zeta(lbi:ubi,lbj:ubj,:)
95
96 real(r8), intent(inout) :: tl_vbar(lbi:ubi,lbj:ubj,:)
97# endif
98!
99! Local variable declarations.
100!
101 integer :: jmin, jmax
102 integer :: i, j, know
103
104 real(r8) :: ce, cx, ze
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_ze
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! Set time-indices
122!-----------------------------------------------------------------------
123!
124 IF (first_2d_step) THEN
125 know=krhs
126 dt2d=dtfast(ng)
127 ELSE IF (predictor_2d_step(ng)) THEN
128 know=krhs
129 dt2d=2.0_r8*dtfast(ng)
130 ELSE
131 know=kstp
132 dt2d=dtfast(ng)
133 END IF
134# if defined ATM_PRESS && defined PRESS_COMPENSATE
135 oneatm=1013.25_r8 ! 1 atm = 1013.25 mb
136 fac=100.0_r8/(g*rho0)
137# endif
138!
139!-----------------------------------------------------------------------
140! Lateral boundary conditions at the southern edge.
141!-----------------------------------------------------------------------
142!
143 IF (domain(ng)%Southern_Edge(tile)) THEN
144!
145! Southern edge, implicit upstream radiation condition.
146!
147 IF (tl_lbc(isouth,isvbar,ng)%radiation) THEN
148 IF (iic(ng).ne.0) THEN
149 DO i=istr,iend+1
150!^ grad(i,Jstr)=vbar(i ,Jstr,know)- &
151!^ & vbar(i-1,Jstr,know)
152!^
153 tl_grad(i,jstr)=0.0_r8
154 END DO
155 DO i=istr,iend
156 IF (lbc_apply(ng)%south(i)) THEN
157# if defined CELERITY_READ && defined FORWARD_READ
158 IF (tl_lbc(isouth,isvbar,ng)%nudging) THEN
159 IF (lnudgem2clm(ng)) THEN
160 obc_out=0.5_r8* &
161 & (clima(ng)%M2nudgcof(i,jstr-1)+ &
162 & clima(ng)%M2nudgcof(i,jstr ))
163 obc_in =obcfac(ng)*obc_out
164 ELSE
165 obc_out=m2obc_out(ng,isouth)
166 obc_in =m2obc_in(ng,isouth)
167 END IF
168 IF (boundary(ng)%vbar_south_Ce(i).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
175# ifdef RADIATION_2D
176 cx=boundary(ng)%vbar_south_Cx(i)
177# else
178 cx=0.0_r8
179# endif
180 ce=boundary(ng)%vbar_south_Ce(i)
181 cff=boundary(ng)%vbar_south_C2(i)
182# endif
183!^ vbar(i,Jstr,kout)=(cff*vbar(i,Jstr ,know)+ &
184!^ & Ce *vbar(i,Jstr+1,kout)- &
185!^ & MAX(Cx,0.0_r8)*grad(i ,Jstr)- &
186!^ & MIN(Cx,0.0_r8)*grad(i+1,Jstr))/ &
187!^ & (cff+Ce)
188!^
189 tl_vbar(i,jstr,kout)=(cff*tl_vbar(i,jstr ,know)+ &
190 & ce *tl_vbar(i,jstr+1,kout)- &
191 & max(cx,0.0_r8)* &
192 & tl_grad(i ,jstr)- &
193 & min(cx,0.0_r8)* &
194 & tl_grad(i+1,jstr))/ &
195 & (cff+ce)
196
197 IF (tl_lbc(isouth,isvbar,ng)%nudging) THEN
198!^ vbar(i,Jstr,kout)=vbar(i,Jstr,kout)+ &
199!^ & tau*(BOUNDARY(ng)%vbar_south(i)- &
200!^ & vbar(i,Jstr,know))
201!^
202 tl_vbar(i,jstr,kout)=tl_vbar(i,jstr,kout)- &
203 & tau*tl_vbar(i,jstr,know)
204 END IF
205# ifdef MASKING
206!^ vbar(i,Jstr,kout)=vbar(i,Jstr,kout)* &
207!^ & GRID(ng)%vmask(i,Jstr)
208!^
209 tl_vbar(i,jstr,kout)=tl_vbar(i,jstr,kout)* &
210 & grid(ng)%vmask(i,jstr)
211# endif
212 END IF
213 END DO
214 END IF
215!
216! Southern edge, Flather boundary condition.
217!
218 ELSE IF (tl_lbc(isouth,isvbar,ng)%Flather) THEN
219 DO i=istr,iend
220 IF (lbc_apply(ng)%south(i)) THEN
221# if defined SSH_TIDES_NOT_YET && !defined UV_TIDES_NOT_YET
222 IF (tl_lbc(isouth,isfsur,ng)%acquire) THEN
223 bry_pgr=-g*(zeta(i,jstr,know)- &
224 & boundary(ng)%zeta_south(i))* &
225 & 0.5_r8*grid(ng)%pn(i,jstr)
226 tl_bry_pgr=-g*(tl_zeta(i,jstr,know)- &
227 & boundary(ng)%tl_zeta_south(i))* &
228 & 0.5_r8*grid(ng)%pn(i,jstr)
229 ELSE
230 bry_pgr=-g*(zeta(i,jstr ,know)- &
231 & zeta(i,jstr-1,know))* &
232 & 0.5_r8*(grid(ng)%pn(i,jstr-1)+ &
233 & grid(ng)%pn(i,jstr ))
234 tl_bry_pgr=-g*(tl_zeta(i,jstr ,know)- &
235 & tl_zeta(i,jstr-1,know))* &
236 & 0.5_r8*(grid(ng)%pn(i,jstr-1)+ &
237 & grid(ng)%pn(i,jstr ))
238 END IF
239# ifdef UV_COR
240 bry_cor=-0.125_r8*(ubar(i ,jstr-1,know)+ &
241 & ubar(i+1,jstr-1,know)+ &
242 & ubar(i ,jstr ,know)+ &
243 & ubar(i+1,jstr ,know))* &
244 & (grid(ng)%f(i,jstr-1)+ &
245 & grid(ng)%f(i,jstr ))
246 tl_bry_cor=-0.125_r8*(tl_ubar(i ,jstr-1,know)+ &
247 & tl_ubar(i+1,jstr-1,know)+ &
248 & tl_ubar(i ,jstr ,know)+ &
249 & tl_ubar(i+1,jstr ,know))* &
250 & (grid(ng)%f(i,jstr-1)+ &
251 & grid(ng)%f(i,jstr ))
252# else
253 bry_cor=0.0_r8
254 tl_bry_cor=0.0_r8
255# endif
256 cff1=1.0_r8/(0.5_r8*(grid(ng)%h(i,jstr-1)+ &
257 & zeta(i,jstr-1,know)+ &
258 & grid(ng)%h(i,jstr )+ &
259 & zeta(i,jstr ,know)))
260 tl_cff1=-cff1*cff1*(0.5_r8*(grid(ng)%tl_h(i,jstr-1)+ &
261 & tl_zeta(i,jstr-1,know)+ &
262 & grid(ng)%tl_h(i,jstr )+ &
263 & tl_zeta(i,jstr ,know)))+ &
264# ifdef TL_IOMS
265 & 2.0_r8*cff1
266# endif
267 bry_str=cff1*(forces(ng)%svstr(i,jstr)- &
268 & forces(ng)%bvstr(i,jstr))
269 tl_bry_str=tl_cff1*(forces(ng)%svstr(i,jstr)- &
270 & forces(ng)%bvstr(i,jstr))+ &
271 & cff1*(forces(ng)%tl_svstr(i,jstr)- &
272 & forces(ng)%tl_bvstr(i,jstr))- &
273# ifdef TL_IOMS
274 & bry_str
275# endif
276 ce=1.0_r8/sqrt(g*0.5_r8*(grid(ng)%h(i,jstr-1)+ &
277 & zeta(i,jstr-1,know)+ &
278 & grid(ng)%h(i,jstr )+ &
279 & zeta(i,jstr ,know)))
280 tl_ce=-ce*ce*ce*0.25_r8*g*(grid(ng)%tl_h(i,jstr-1)+ &
281 & tl_zeta(i,jstr-1,know)+ &
282 & grid(ng)%tl_h(i,jstr )+ &
283 & tl_zeta(i,jstr ,know))+ &
284# ifdef TL_IOMS
285 & g*0.5_r8*ce*ce*ce*(grid(ng)%h(i,jstr-1)+ &
286 & zeta(i,jstr-1,know)+ &
287 & grid(ng)%h(i,jstr )+ &
288 & zeta(i,jstr ,know))
289# endif
290 cff2=grid(ng)%on_v(i,jstr)*ce
291 tl_cff2=grid(ng)%on_v(i,jstr)*tl_ce
292!^ bry_val=vbar(i,Jstr+1,know)+ &
293!^ & cff2*(bry_pgr+ &
294!^ & bry_cor+ &
295!^ & bry_str)
296!^
297 tl_bry_val=tl_vbar(i,jstr+1,know)+ &
298 & tl_cff2*(bry_pgr+ &
299 & bry_cor+ &
300 & bry_str)+ &
301 & cff2*(tl_bry_pgr+ &
302 & tl_bry_cor+ &
303 & tl_bry_str)- &
304# ifdef TL_IOMS
305 & cff2*(bry_pgr+bry_cor+bry_str)
306# endif
307# else
308!^ bry_val=BOUNDARY(ng)%vbar_south(i)
309!^
310 tl_bry_val=boundary(ng)%tl_vbar_south(i)
311# endif
312 cff=1.0_r8/(0.5_r8*(grid(ng)%h(i,jstr-1)+ &
313 & zeta(i,jstr-1,know)+ &
314 & grid(ng)%h(i,jstr )+ &
315 & zeta(i,jstr ,know)))
316 tl_cff=-cff*cff*(0.5_r8*(grid(ng)%tl_h(i,jstr-1)+ &
317 & tl_zeta(i,jstr-1,know)+ &
318 & grid(ng)%tl_h(i,jstr )+ &
319 & tl_zeta(i,jstr ,know)))+ &
320# ifdef TL_IOMS
321 & 2.0_r8*cff
322# endif
323 ce=sqrt(g*cff)
324 tl_ce=0.5_r8*g*tl_cff/ce+ &
325# ifdef TL_IOMS
326 & 0.5_r8*ce
327# endif
328# if defined ATM_PRESS && defined PRESS_COMPENSATE
329!^ vbar(i,Jstr,kout)=bry_val- &
330!^ & Ce*(0.5_r8* &
331!^ & (zeta(i,Jstr-1,know)+ &
332!^ & zeta(i,Jstr ,know)+ &
333!^ & fac*(FORCES(ng)%Pair(i,Jstr-1)+ &
334!^ & FORCES(ng)%Pair(i,Jstr )- &
335!^ & 2.0_r8*OneAtm))- &
336!^ & BOUNDARY(ng)%zeta_south(i))
337!^
338 tl_vbar(i,jstr,kout)=tl_bry_val- &
339 & tl_ce* &
340 & (0.5_r8* &
341 & (zeta(i,jstr-1,know)+ &
342 & zeta(i,jstr ,know)+ &
343 & fac*(forces(ng)%Pair(i,jstr-1)+ &
344 & forces(ng)%Pair(i,jstr )- &
345 & 2.0_r8*oneatm))- &
346 & boundary(ng)%zeta_south(i))- &
347 & ce* &
348 & (0.5_r8*(tl_zeta(i,jstr-1,know)+ &
349 & tl_zeta(i,jstr ,know)))+ &
350# ifdef TL_IOMS
351 & ce* &
352 & (0.5_r8* &
353 & (zeta(i,jstr-1,know)+ &
354 & zeta(i,jstr ,know)+ &
355 & fac*(forces(ng)%Pair(i,jstr-1)+ &
356 & forces(ng)%Pair(i,jstr )- &
357 & 2.0_r8*oneatm))- &
358 & boundary(ng)%zeta_south(i))
359# endif
360# else
361!^ vbar(i,Jstr,kout)=bry_val- &
362!^ & Ce*(0.5_r8*(zeta(i,Jstr-1,know)+ &
363!^ & zeta(i,Jstr ,know))- &
364!^ & BOUNDARY(ng)%zeta_south(i))
365!^
366 tl_vbar(i,jstr,kout)=tl_bry_val- &
367 & tl_ce* &
368 & (0.5_r8*(zeta(i,jstr-1,know)+ &
369 & zeta(i,jstr ,know))- &
370 & boundary(ng)%zeta_south(i))- &
371 & ce* &
372 & (0.5_r8*(tl_zeta(i,jstr-1,know)+ &
373 & tl_zeta(i,jstr ,know))- &
374 & boundary(ng)%tl_zeta_south(i))+ &
375# ifdef TL_IOMS
376 & ce* &
377 & (0.5_r8*(zeta(i,jstr-1,know)+ &
378 & zeta(i,jstr ,know))- &
379 & boundary(ng)%zeta_south(i))
380# endif
381# endif
382# ifdef MASKING
383!^ vbar(i,Jstr,kout)=vbar(i,Jstr,kout)* &
384!^ & GRID(ng)%vmask(i,Jstr)
385!^
386 tl_vbar(i,jstr,kout)=tl_vbar(i,jstr,kout)* &
387 & grid(ng)%vmask(i,jstr)
388# endif
389 END IF
390 END DO
391!
392! Southern edge, Shchepetkin boundary condition (Maison et al., 2010).
393!
394 ELSE IF (tl_lbc(isouth,isvbar,ng)%Shchepetkin) THEN
395 DO i=istr,iend
396 IF (lbc_apply(ng)%south(i)) THEN
397# if defined SSH_TIDES_NOT_YET && !defined UV_TIDES_NOT_YET
398 IF (tl_lbc(isouth,isfsur,ng)%acquire) THEN
399 bry_pgr=-g*(zeta(i,jstr,know)- &
400 & boundary(ng)%zeta_south(i))* &
401 & 0.5_r8*grid(ng)%pn(i,jstr)
402 tl_bry_pgr=-g*(tl_zeta(i,jstr,know)- &
403 & boundary(ng)%tl_zeta_south(i))* &
404 & 0.5_r8*grid(ng)%pn(i,jstr)
405 ELSE
406 bry_pgr=-g*(zeta(i,jstr ,know)- &
407 & zeta(i,jstr-1,know))* &
408 & 0.5_r8*(grid(ng)%pn(i,jstr-1)+ &
409 & grid(ng)%pn(i,jstr ))
410 tl_bry_pgr=-g*(tl_zeta(i,jstr ,know)- &
411 & tl_zeta(i,jstr-1,know))* &
412 & 0.5_r8*(grid(ng)%pn(i,jstr-1)+ &
413 & grid(ng)%pn(i,jstr ))
414 END IF
415# ifdef UV_COR
416 bry_cor=-0.125_r8*(ubar(i ,jstr-1,know)+ &
417 & ubar(i+1,jstr-1,know)+ &
418 & ubar(i ,jstr ,know)+ &
419 & ubar(i+1,jstr ,know))* &
420 & (grid(ng)%f(i,jstr-1)+ &
421 & grid(ng)%f(i,jstr ))
422 tl_bry_cor=-0.125_r8*(tl_ubar(i ,jstr-1,know)+ &
423 & tl_ubar(i+1,jstr-1,know)+ &
424 & tl_ubar(i ,jstr ,know)+ &
425 & tl_ubar(i+1,jstr ,know))* &
426 & (grid(ng)%f(i,jstr-1)+ &
427 & grid(ng)%f(i,jstr ))
428# else
429 bry_cor=0.0_r8
430 tl_bry_cor=0.0_r8
431# endif
432 cff1=1.0_r8/(0.5_r8*(grid(ng)%h(i,jstr-1)+ &
433 & zeta(i,jstr-1,know)+ &
434 & grid(ng)%h(i,jstr )+ &
435 & zeta(i,jstr ,know)))
436 tl_cff1=-cff1*cff1*(0.5_r8*(grid(ng)%tl_h(i,jstr-1)+ &
437 & tl_zeta(i,jstr-1,know)+ &
438 & grid(ng)%tl_h(i,jstr )+ &
439 & tl_zeta(i,jstr ,know)))+ &
440# ifdef TL_IOMS
441 & 2.0_r8*cff1
442# endif
443 bry_str=cff1*(forces(ng)%svstr(i,jstr)- &
444 & forces(ng)%bvstr(i,jstr))
445 tl_bry_str=tl_cff1*(forces(ng)%svstr(i,jstr)- &
446 & forces(ng)%bvstr(i,jstr))+ &
447 & cff1*(forces(ng)%tl_svstr(i,jstr)- &
448 & forces(ng)%tl_bvstr(i,jstr))- &
449# ifdef TL_IOMS
450 & bry_str
451# endif
452 ce=1.0_r8/sqrt(g*0.5_r8*(grid(ng)%h(i,jstr-1)+ &
453 & zeta(i,jstr-1,know)+ &
454 & grid(ng)%h(i,jstr )+ &
455 & zeta(i,jstr ,know)))
456 tl_ce=-ce*ce*ce*0.25_r8*g*(grid(ng)%tl_h(i,jstr-1)+ &
457 & tl_zeta(i,jstr-1,know)+ &
458 & grid(ng)%tl_h(i,jstr )+ &
459 & tl_zeta(i,jstr ,know))+ &
460# ifdef TL_IOMS
461 & g*0.5_r8*ce*ce*ce*(grid(ng)%h(i,jstr-1)+ &
462 & zeta(i,jstr-1,know)+ &
463 & grid(ng)%h(i,jstr )+ &
464 & zeta(i,jstr ,know))
465# endif
466 cff2=grid(ng)%on_v(i,jstr)*ce
467 tl_cff2=grid(ng)%on_v(i,jstr)*tl_ce
468!^ bry_val=vbar(i,Jstr+1,know)+ &
469!^ & cff2*(bry_pgr+ &
470!^ & bry_cor+ &
471!^ & bry_str)
472!^
473 tl_bry_val=tl_vbar(i,jstr+1,know)+ &
474 & tl_cff2*(bry_pgr+ &
475 & bry_cor+ &
476 & bry_str)+ &
477 & cff2*(tl_bry_pgr+ &
478 & tl_bry_cor+ &
479 & tl_bry_str)- &
480# ifdef TL_IOMS
481 & cff2*(bry_pgr+bry_cor+bry_str)
482# endif
483# else
484!^ bry_val=BOUNDARY(ng)%vbar_south(i)
485!^
486 tl_bry_val=boundary(ng)%tl_vbar_south(i)
487# endif
488# ifdef WET_DRY_NOT_YET
489 cff=0.5_r8*(grid(ng)%h(i,jstr-1)+ &
490 & zeta(i,jstr-1,know)+ &
491 & grid(ng)%h(i,jstr )+ &
492 & zeta(i,jstr ,know))
493 tl_cff=0.5_r8*(grid(ng)%tl_h(i,jstr-1)+ &
494 & tl_zeta(i,jstr-1,know)+ &
495 & grid(ng)%tl_h(i,jstr )+ &
496 & tl_zeta(i,jstr ,know))
497# else
498 cff=0.5_r8*(grid(ng)%h(i,jstr-1)+ &
499 & grid(ng)%h(i,jstr ))
500 tl_cff=0.5_r8*(grid(ng)%tl_h(i,jstr-1)+ &
501 & grid(ng)%tl_h(i,jstr ))
502# endif
503 cff1=sqrt(g/cff)
504 tl_cff1=-0.5_r8*cff1*tl_cff/cff+ &
505# ifdef TL_IOMS
506 & 0.5_r8*cff1
507# endif
508 ce=dt2d*cff1*cff*0.5_r8*(grid(ng)%pn(i,jstr-1)+ &
509 & grid(ng)%pn(i,jstr ))
510 tl_ce=dt2d*0.5_r8*(grid(ng)%pn(i,jstr-1)+ &
511 & grid(ng)%pn(i,jstr ))* &
512 & (cff1*tl_cff+ &
513 & tl_cff1*cff)- &
514# ifdef TL_IOMS
515 & ce
516# endif
517 ze=(0.5_r8+ce)*zeta(i,jstr ,know)+ &
518 & (0.5_r8-ce)*zeta(i,jstr-1,know)
519 tl_ze=(0.5_r8+ce)*tl_zeta(i,jstr ,know)+ &
520 & (0.5_r8-ce)*tl_zeta(i,jstr-1,know)+ &
521 & tl_ce*(zeta(i,jstr ,know)- &
522 & zeta(i,jstr-1,know))- &
523# ifdef TL_IOMS
524 & ze ! HGA check
525# endif
526 IF (ce.gt.co) THEN
527 cff2=(1.0_r8-co/ce)**2
528 tl_cff2=2.0_r8*cff2*co*tl_ce/(ce*ce)- &
529# ifdef TL_IOMS
530 & cff2 ! HGA check
531# endif
532 cff3=zeta(i,jstr,kout)+ &
533 & ce*zeta(i,jstr-1,know)- &
534 & (1.0_r8+ce)*zeta(i,jstr,know)
535 tl_cff3=tl_zeta(i,jstr,kout)+ &
536 & ce*tl_zeta(i,jstr-1,know)+ &
537 & tl_ce*(zeta(i,jstr-1,know)+ &
538 & zeta(i,jstr ,know))- &
539 & (1.0_r8+ce)*tl_zeta(i,jstr,know)- &
540# ifdef TL_IOMS
541 & ce*zeta(i,jstr-1,know)+ &
542 & (1.0_r8+ce)*zeta(i,jstr,know) ! HGA check
543# endif
544 ze=ze+cff2*cff3
545 tl_ze=tl_ze+cff2*tl_cff3+ &
546 & tl_cff2*cff3- &
547# ifdef TL_IOMS
548 & cff2*cff3 ! HGA check
549# endif
550 END IF
551!^ vbar(i,Jstr,kout)=0.5_r8* &
552!^ & ((1.0_r8-Ce)*vbar(i,Jstr,know)+ &
553!^ & Ce*vbar(i,Jstr+1,know)+ &
554!^ & bry_val- &
555!^ & cff1*(Ze-BOUNDARY(ng)%zeta_south(i)))
556!^
557 tl_vbar(i,jstr,kout)=0.5_r8* &
558 & ((1.0_r8-ce)* &
559 & tl_vbar(i,jstr,know)- &
560 & tl_ce*(vbar(i,jstr ,know)- &
561 & vbar(i,jstr+1,know))+ &
562 & ce*tl_vbar(i,jstr+1,know)+ &
563 & tl_bry_val- &
564 & tl_cff1* &
565 & (ze-boundary(ng)%zeta_south(i))- &
566 & cff1*tl_ze)- &
567# ifdef TL_IOMS
568 & 0.5_r8* &
569 & ((1.0_r8-ce)*vbar(i,jstr,know)+ &
570 & ce*vbar(i,jstr+1,know)- &
571 & cff1* &
572 & (ze-boundary(ng)%zeta_south(i)))
573!! ! HGA check
574# endif
575# ifdef ADJUST_BOUNDARY
576 IF (lobc(isouth,isvbar,ng)) THEN
577 tl_vbar(i,jstr,kout)=tl_vbar(i,jstr,kout)+ &
578 & 0.5_r8*cff1* &
579 & boundary(ng)%tl_zeta_south(i)
580 END IF
581# endif
582# ifdef MASKING
583!^ vbar(i,Jstr,kout)=vbar(i,Jstr,kout)* &
584!^ & GRID(ng)%vmask(i,Jstr)
585!^
586 tl_vbar(i,jstr,kout)=tl_vbar(i,jstr,kout)* &
587 & grid(ng)%vmask(i,jstr)
588# endif
589 END IF
590 END DO
591!
592! Southern edge, clamped boundary condition.
593!
594 ELSE IF (tl_lbc(isouth,isvbar,ng)%clamped) THEN
595 DO i=istr,iend
596 IF (lbc_apply(ng)%south(i)) THEN
597!^ vbar(i,Jstr,kout)=BOUNDARY(ng)%vbar_south(i)
598!^
599 tl_vbar(i,jstr,kout)=boundary(ng)%tl_vbar_south(i)
600# ifdef MASKING
601!^ vbar(i,Jstr,kout)=vbar(i,Jstr,kout)* &
602!^ & GRID(ng)%vmask(i,Jstr)
603!^
604 tl_vbar(i,jstr,kout)=tl_vbar(i,jstr,kout)* &
605 & grid(ng)%vmask(i,jstr)
606# endif
607 END IF
608 END DO
609!
610! Southern edge, gradient boundary condition.
611!
612 ELSE IF (tl_lbc(isouth,isvbar,ng)%gradient) THEN
613 DO i=istr,iend
614 IF (lbc_apply(ng)%south(i)) THEN
615!^ vbar(i,Jstr,kout)=vbar(i,Jstr+1,kout)
616!^
617 tl_vbar(i,jstr,kout)=tl_vbar(i,jstr+1,kout)
618# ifdef MASKING
619!^ vbar(i,Jstr,kout)=vbar(i,Jstr,kout)* &
620!^ & GRID(ng)%vmask(i,Jstr)
621!^
622 tl_vbar(i,jstr,kout)=tl_vbar(i,jstr,kout)* &
623 & grid(ng)%vmask(i,jstr)
624# endif
625 END IF
626 END DO
627!
628! Southern edge, reduced-physics boundary condition.
629!
630 ELSE IF (tl_lbc(isouth,isvbar,ng)%reduced) THEN
631 DO i=istr,iend
632 IF (lbc_apply(ng)%south(i)) THEN
633 IF (tl_lbc(isouth,isfsur,ng)%acquire) THEN
634!^ bry_pgr=-g*(zeta(i,Jstr,know)- &
635!^ & BOUNDARY(ng)%zeta_south(i))* &
636!^ & 0.5_r8*GRID(ng)%pn(i,Jstr)
637!^
638 tl_bry_pgr=-g*(tl_zeta(i,jstr,know)- &
639 & boundary(ng)%tl_zeta_south(i))* &
640 & 0.5_r8*grid(ng)%pn(i,jstr)
641 ELSE
642!^ bry_pgr=-g*(zeta(i,Jstr,know)- &
643!^ & zeta(i,Jstr-1,know))* &
644!^ & 0.5_r8*(GRID(ng)%pn(i,Jstr-1)+ &
645!^ & GRID(ng)%pn(i,Jstr ))
646!^
647 tl_bry_pgr=-g*(tl_zeta(i,jstr ,know)- &
648 & tl_zeta(i,jstr-1,know))* &
649 & 0.5_r8*(grid(ng)%pn(i,jstr-1)+ &
650 & grid(ng)%pn(i,jstr ))
651 END IF
652# ifdef UV_COR
653!^ bry_cor=-0.125_r8*(ubar(i ,Jstr-1,know)+ &
654!^ & ubar(i+1,Jstr-1,know)+ &
655!^ & ubar(i ,Jstr ,know)+ &
656!^ & ubar(i+1,Jstr ,know))* &
657!^ & (GRID(ng)%f(i,Jstr-1)+ &
658!^ & GRID(ng)%f(i,Jstr ))
659!^
660 tl_bry_cor=-0.125_r8*(tl_ubar(i ,jstr-1,know)+ &
661 & tl_ubar(i+1,jstr-1,know)+ &
662 & tl_ubar(i ,jstr ,know)+ &
663 & tl_ubar(i+1,jstr ,know))* &
664 & (grid(ng)%f(i,jstr-1)+ &
665 & grid(ng)%f(i,jstr ))
666# else
667!^ bry_cor=0.0_r8
668!^
669 tl_bry_cor=0.0_r8
670# endif
671 cff=1.0_r8/(0.5_r8*(grid(ng)%h(i,jstr-1)+ &
672 & zeta(i,jstr-1,know)+ &
673 & grid(ng)%h(i,jstr )+ &
674 & zeta(i,jstr ,know)))
675 tl_cff=-cff*cff*0.5_r8*(grid(ng)%tl_h(i,jstr-1)+ &
676 & tl_zeta(i,jstr-1,know)+ &
677 & grid(ng)%tl_h(i,jstr )+ &
678 & tl_zeta(i,jstr ,know))+ &
679# ifdef TL_IOMS
680 & 2.0_r8*cff
681# endif
682!^ bry_str=cff*(FORCES(ng)%svstr(i,Jstr)- &
683!^ & FORCES(ng)%bvstr(i,Jstr))
684!^
685 tl_bry_str=tl_cff*(forces(ng)%svstr(i,jstr)- &
686 & forces(ng)%bvstr(i,jstr))+ &
687 & cff*(forces(ng)%tl_svstr(i,jstr)- &
688 & forces(ng)%tl_bvstr(i,jstr))- &
689# ifdef TL_IOMS
690 & cff*(forces(ng)%svstr(i,jstr)- &
691 & forces(ng)%bvstr(i,jstr))
692# endif
693!^ vbar(i,Jstr,kout)=vbar(i,Jstr,know)+ &
694!^ & dt2d*(bry_pgr+ &
695!^ & bry_cor+ &
696!^ & bry_str)
697!^
698 tl_vbar(i,jstr,kout)=tl_vbar(i,jstr,know)+ &
699 & dt2d*(tl_bry_pgr+ &
700 & tl_bry_cor+ &
701 & tl_bry_str)
702# ifdef MASKING
703!^ vbar(i,Jstr,kout)=vbar(i,Jstr,kout)* &
704!^ & GRID(ng)%vmask(i,Jstr)
705!^
706 tl_vbar(i,jstr,kout)=tl_vbar(i,jstr,kout)* &
707 & grid(ng)%vmask(i,jstr)
708# endif
709 END IF
710 END DO
711!
712! Southern edge, closed boundary condition.
713!
714 ELSE IF (tl_lbc(isouth,isvbar,ng)%closed) THEN
715 DO i=istr,iend
716 IF (lbc_apply(ng)%south(i)) THEN
717!^ vbar(i,Jstr,kout)=0.0_r8
718!^
719 tl_vbar(i,jstr,kout)=0.0_r8
720 END IF
721 END DO
722 END IF
723 END IF
724!
725!-----------------------------------------------------------------------
726! Lateral boundary conditions at the northern edge.
727!-----------------------------------------------------------------------
728!
729 IF (domain(ng)%Northern_Edge(tile)) THEN
730!
731! Northern edge, implicit upstream radiation condition.
732!
733 IF (tl_lbc(inorth,isvbar,ng)%radiation) THEN
734 IF (iic(ng).ne.0) THEN
735 DO i=istr,iend+1
736!^ grad(i,Jend+1)=vbar(i ,Jend+1,know)- &
737!^ & vbar(i-1,Jend+1,know)
738!^
739 tl_grad(i,jend+1)=0.0_r8
740 END DO
741 DO i=istr,iend
742 IF (lbc_apply(ng)%north(i)) THEN
743# if defined CELERITY_READ && defined FORWARD_READ
744 IF (tl_lbc(inorth,isvbar,ng)%nudging) THEN
745 IF (lnudgem2clm(ng)) THEN
746 obc_out=0.5_r8* &
747 & (clima(ng)%M2nudgcof(i,jend )+ &
748 & clima(ng)%M2nudgcof(i,jend+1))
749 obc_in =obcfac(ng)*obc_out
750 ELSE
751 obc_out=m2obc_out(ng,inorth)
752 obc_in =m2obc_in(ng,inorth)
753 END IF
754 IF (boundary(ng)%vbar_north_Ce(i).lt.0.0_r8) THEN
755 tau=obc_in
756 ELSE
757 tau=obc_out
758 END IF
759 tau=tau*dt2d
760 END IF
761# ifdef RADIATION_2D
762 cx=boundary(ng)%vbar_north_Cx(i)
763# else
764 cx=0.0_r8
765# endif
766 ce=boundary(ng)%vbar_north_Ce(i)
767 cff=boundary(ng)%vbar_north_C2(i)
768# endif
769!^ vbar(i,Jend+1,kout)=(cff*vbar(i,Jend+1,know)+ &
770!^ & Ce *vbar(i,Jend ,kout)- &
771!^ & MAX(Cx,0.0_r8)*grad(i ,Jend+1)- &
772!^ & MIN(Cx,0.0_r8)*grad(i+1,Jend+1))/ &
773!^ & (cff+Ce)
774!^
775 tl_vbar(i,jend+1,kout)=(cff*tl_vbar(i,jend+1,know)+ &
776 & ce *tl_vbar(i,jend ,kout)- &
777 & max(cx,0.0_r8)* &
778 & tl_grad(i ,jend+1)- &
779 & min(cx,0.0_r8)* &
780 & tl_grad(i+1,jend+1))/ &
781 & (cff+ce)
782
783 IF (tl_lbc(inorth,isvbar,ng)%nudging) THEN
784!^ vbar(i,Jend+1,kout)=vbar(i,Jend+1,kout)+ &
785!^ & tau*(BOUNDARY(ng)%vbar_north(i)- &
786!^ & vbar(i,Jend+1,know))
787!^
788 tl_vbar(i,jend+1,kout)=tl_vbar(i,jend+1,kout)- &
789 & tau*tl_vbar(i,jend+1,know)
790 END IF
791# ifdef MASKING
792!^ vbar(i,Jend+1,kout)=vbar(i,Jend+1,kout)* &
793!^ & GRID(ng)%vmask(i,Jend+1)
794!^
795 tl_vbar(i,jend+1,kout)=tl_vbar(i,jend+1,kout)* &
796 & grid(ng)%vmask(i,jend+1)
797# endif
798 END IF
799 END DO
800 END IF
801!
802! Northern edge, Flather boundary condition.
803!
804 ELSE IF (tl_lbc(inorth,isvbar,ng)%Flather) THEN
805 DO i=istr,iend
806 IF (lbc_apply(ng)%north(i)) THEN
807# if defined SSH_TIDES_NOT_YET && !defined UV_TIDES_NOT_YET
808 IF (tl_lbc(inorth,isfsur,ng)%acquire) THEN
809 bry_pgr=-g*(boundary(ng)%zeta_north(i)- &
810 & zeta(i,jend,know))* &
811 & 0.5_r8*grid(ng)%pn(i,jend)
812 tl_bry_pgr=-g*(boundary(ng)%tl_zeta_north(i)- &
813 & tl_zeta(i,jend,know))* &
814 & 0.5_r8*grid(ng)%pn(i,jend)
815 ELSE
816 bry_pgr=-g*(zeta(i,jend+1,know)- &
817 & zeta(i,jend ,know))* &
818 & 0.5_r8*(grid(ng)%pn(i,jend )+ &
819 & grid(ng)%pn(i,jend+1))
820 tl_bry_pgr=-g*(tl_zeta(i,jend+1,know)- &
821 & tl_zeta(i,jend ,know))* &
822 & 0.5_r8*(grid(ng)%pn(i,jend )+ &
823 & grid(ng)%pn(i,jend+1))
824 END IF
825# ifdef UV_COR
826 bry_cor=-0.125_r8*(ubar(i ,jend ,know)+ &
827 & ubar(i+1,jend ,know)+ &
828 & ubar(i ,jend+1,know)+ &
829 & ubar(i+1,jend+1,know))* &
830 & (grid(ng)%f(i,jend )+ &
831 & grid(ng)%f(i,jend+1))
832 tl_bry_cor=-0.125_r8*(tl_ubar(i ,jend ,know)+ &
833 & tl_ubar(i+1,jend ,know)+ &
834 & tl_ubar(i ,jend+1,know)+ &
835 & tl_ubar(i+1,jend+1,know))* &
836 & (grid(ng)%f(i,jend )+ &
837 & grid(ng)%f(i,jend+1))
838# else
839 bry_cor=0.0_r8
840 tl_bry_cor=0.0_r8
841# endif
842 cff1=1.0_r8/(0.5_r8*(grid(ng)%h(i,jend )+ &
843 & zeta(i,jend ,know)+ &
844 & grid(ng)%h(i,jend+1)+ &
845 & zeta(i,jend+1,know)))
846 tl_cff1=-cff1*cff1*0.5_r8*(grid(ng)%tl_h(i,jend )+ &
847 & tl_zeta(i,jend ,know)+ &
848 & grid(ng)%tl_h(i,jend+1)+ &
849 & tl_zeta(i,jend+1,know))+ &
850# ifdef TL_IOMS
851 & 2.0_r8*cff1
852# endif
853 bry_str=cff1*(forces(ng)%svstr(i,jend+1)- &
854 & forces(ng)%bvstr(i,jend+1))
855 tl_bry_str=tl_cff1*(forces(ng)%svstr(i,jend+1)- &
856 & forces(ng)%bvstr(i,jend+1))+ &
857 & cff1*(forces(ng)%tl_svstr(i,jend+1)- &
858 & forces(ng)%tl_bvstr(i,jend+1))- &
859# ifdef TL_IOMS
860 & bry_str
861# endif
862 ce=1.0_r8/sqrt(g*0.5_r8*(grid(ng)%h(i,jend+1)+ &
863 & zeta(i,jend+1,know)+ &
864 & grid(ng)%h(i,jend )+ &
865 & zeta(i,jend ,know)))
866 tl_ce=-ce*ce*ce*0.25_r8*g*(grid(ng)%tl_h(i,jend+1)+ &
867 & tl_zeta(i,jend+1,know)+ &
868 & grid(ng)%tl_h(i,jend )+ &
869 & tl_zeta(i,jend ,know))+ &
870# ifdef TL_IOMS
871 & g*0.5_r8*ce*ce*ce*(grid(ng)%h(i,jend+1)+ &
872 & zeta(i,jend+1,know)+ &
873 & grid(ng)%h(i,jend )+ &
874 & zeta(i,jend ,know))
875# endif
876 cff2=grid(ng)%on_v(i,jend+1)*ce
877 tl_cff2=grid(ng)%on_v(i,jend+1)*tl_ce
878!^ bry_val=vbar(i,Jend,know)+ &
879!^ & cff2*(bry_pgr+ &
880!^ & bry_cor+ &
881!^ & bry_str)
882!^
883 tl_bry_val=tl_vbar(i,jend,know)+ &
884 & tl_cff2*(bry_pgr+ &
885 & bry_cor+ &
886 & bry_str)+ &
887 & cff2*(tl_bry_pgr+ &
888 & tl_bry_cor+ &
889 & tl_bry_str)- &
890# ifdef TL_IOMS
891 & cff2*(bry_pgr+bry_cor+bry_str)
892# endif
893# else
894!^ bry_val=BOUNDARY(ng)%vbar_north(i)
895!^
896 tl_bry_val=boundary(ng)%tl_vbar_north(i)
897# endif
898 cff=1.0_r8/(0.5_r8*(grid(ng)%h(i,jend )+ &
899 & zeta(i,jend ,know)+ &
900 & grid(ng)%h(i,jend+1)+ &
901 & zeta(i,jend+1,know)))
902 tl_cff=-cff*cff*(0.5_r8*(grid(ng)%tl_h(i,jend )+ &
903 & tl_zeta(i,jend ,know)+ &
904 & grid(ng)%tl_h(i,jend+1)+ &
905 & tl_zeta(i,jend+1,know)))+ &
906# ifdef TL_IOMS
907 & 2.0_r8*cff
908# endif
909 ce=sqrt(g*cff)
910 tl_ce=0.5_r8*g*tl_cff/ce+ &
911# ifdef TL_IOMS
912 & 0.5_r8*ce
913# endif
914# if defined ATM_PRESS && defined PRESS_COMPENSATE
915!^ vbar(i,Jend+1,kout)=bry_val+ &
916!^ & Ce*(0.5_r8* &
917!^ & (zeta(i,Jend ,know)+ &
918!^ & zeta(i,Jend+1,know)+ &
919!^ & fac*(FORCES(ng)%Pair(i,Jend )+ &
920!^ & FORCES(ng)%Pair(i,Jend+1)- &
921!^ & 2.0_r8*OneAtm))- &
922!^ & BOUNDARY(ng)%zeta_north(i))
923!^
924 tl_vbar(i,jend+1,kout)=tl_bry_val+ &
925 & tl_ce* &
926 & (0.5_r8* &
927 & (zeta(i,jend ,know)+ &
928 & zeta(i,jend+1,know)+ &
929 & fac*(forces(ng)%Pair(i,jend )+ &
930 & forces(ng)%Pair(i,jend+1)- &
931 & 2.0_r8*oneatm))- &
932 & boundary(ng)%zeta_north(i))+ &
933 & ce* &
934 & (0.5_r8*(tl_zeta(i,jend ,know)+ &
935 & tl_zeta(i,jend+1,know)))- &
936# ifdef TL_IOMS
937 & ce* &
938 & (0.5_r8* &
939 & (zeta(i,jend ,know)+ &
940 & zeta(i,jend+1,know)+ &
941 & fac*(forces(ng)%Pair(i,jend )+ &
942 & forces(ng)%Pair(i,jend+1)- &
943 & 2.0_r8*oneatm))- &
944 & boundary(ng)%zeta_north(i))
945# endif
946# else
947!^ vbar(i,Jend+1,kout)=bry_val+ &
948!^ & Ce*(0.5_r8*(zeta(i,Jend ,know)+ &
949!^ & zeta(i,Jend+1,know))- &
950!^ & BOUNDARY(ng)%zeta_north(i))
951!^
952 tl_vbar(i,jend+1,kout)=tl_bry_val+ &
953 & tl_ce* &
954 & (0.5_r8*(zeta(i,jend ,know)+ &
955 & zeta(i,jend+1,know))- &
956 & boundary(ng)%zeta_north(i))+ &
957 & ce* &
958 & (0.5_r8*(tl_zeta(i,jend ,know)+ &
959 & tl_zeta(i,jend+1,know))- &
960 & boundary(ng)%tl_zeta_north(i))- &
961# ifdef TL_IOMS
962 & ce* &
963 & (0.5_r8*(zeta(i,jend ,know)+ &
964 & zeta(i,jend+1,know))- &
965 & boundary(ng)%zeta_north(i))
966# endif
967# endif
968# ifdef MASKING
969!^ vbar(i,Jend+1,kout)=vbar(i,Jend+1,kout)* &
970!^ & GRID(ng)%vmask(i,Jend+1)
971!^
972 tl_vbar(i,jend+1,kout)=tl_vbar(i,jend+1,kout)* &
973 & grid(ng)%vmask(i,jend+1)
974# endif
975 END IF
976 END DO
977!
978! Northern edge, Shchepetkin boundary condition (Maison et al., 2010).
979!
980 ELSE IF (lbc(inorth,isvbar,ng)%Shchepetkin) THEN
981 DO i=istr,iend
982 IF (lbc_apply(ng)%north(i)) THEN
983# if defined SSH_TIDES_NOT_YET && !defined UV_TIDES_NOT_YET
984 IF (tl_lbc(inorth,isfsur,ng)%acquire) THEN
985 bry_pgr=-g*(boundary(ng)%zeta_north(i)- &
986 & zeta(i,jend,know))* &
987 & 0.5_r8*grid(ng)%pn(i,jend)
988 tl_bry_pgr=-g*(boundary(ng)%tl_zeta_north(i)- &
989 & tl_zeta(i,jend,know))* &
990 & 0.5_r8*grid(ng)%pn(i,jend)
991 ELSE
992 bry_pgr=-g*(zeta(i,jend+1,know)- &
993 & zeta(i,jend ,know))* &
994 & 0.5_r8*(grid(ng)%pn(i,jend )+ &
995 & grid(ng)%pn(i,jend+1))
996 tl_bry_pgr=-g*(tl_zeta(i,jend+1,know)- &
997 & tl_zeta(i,jend ,know))* &
998 & 0.5_r8*(grid(ng)%pn(i,jend )+ &
999 & grid(ng)%pn(i,jend+1))
1000 END IF
1001# ifdef UV_COR
1002 bry_cor=-0.125_r8*(ubar(i ,jend ,know)+ &
1003 & ubar(i+1,jend ,know)+ &
1004 & ubar(i ,jend+1,know)+ &
1005 & ubar(i+1,jend+1,know))* &
1006 & (grid(ng)%f(i,jend )+ &
1007 & grid(ng)%f(i,jend+1))
1008 tl_bry_cor=-0.125_r8*(tl_ubar(i ,jend ,know)+ &
1009 & tl_ubar(i+1,jend ,know)+ &
1010 & tl_ubar(i ,jend+1,know)+ &
1011 & tl_ubar(i+1,jend+1,know))* &
1012 & (grid(ng)%f(i,jend )+ &
1013 & grid(ng)%f(i,jend+1))
1014# else
1015 bry_cor=0.0_r8
1016 tl_bry_cor=0.0_r8
1017# endif
1018 cff1=1.0_r8/(0.5_r8*(grid(ng)%h(i,jend )+ &
1019 & zeta(i,jend ,know)+ &
1020 & grid(ng)%h(i,jend+1)+ &
1021 & zeta(i,jend+1,know)))
1022 tl_cff1=-cff1*cff1*0.5_r8*(grid(ng)%tl_h(i,jend )+ &
1023 & tl_zeta(i,jend ,know)+ &
1024 & grid(ng)%tl_h(i,jend+1)+ &
1025 & tl_zeta(i,jend+1,know))+ &
1026# ifdef TL_IOMS
1027 & 2.0_r8*cff1
1028# endif
1029 bry_str=cff1*(forces(ng)%svstr(i,jend+1)- &
1030 & forces(ng)%bvstr(i,jend+1))
1031 tl_bry_str=tl_cff1*(forces(ng)%svstr(i,jend+1)- &
1032 & forces(ng)%bvstr(i,jend+1))+ &
1033 & cff1*(forces(ng)%tl_svstr(i,jend+1)- &
1034 & forces(ng)%tl_bvstr(i,jend+1))- &
1035# ifdef TL_IOMS
1036 & bry_str
1037# endif
1038 ce=1.0_r8/sqrt(g*0.5_r8*(grid(ng)%h(i,jend+1)+ &
1039 & zeta(i,jend+1,know)+ &
1040 & grid(ng)%h(i,jend )+ &
1041 & zeta(i,jend ,know)))
1042 tl_ce=-ce*ce*ce*0.25_r8*g*(grid(ng)%tl_h(i,jend+1)+ &
1043 & tl_zeta(i,jend+1,know)+ &
1044 & grid(ng)%tl_h(i,jend )+ &
1045 & tl_zeta(i,jend ,know))+ &
1046# ifdef TL_IOMS
1047 & g*0.5_r8*ce*ce*ce*(grid(ng)%h(i,jend+1)+ &
1048 & zeta(i,jend+1,know)+ &
1049 & grid(ng)%h(i,jend )+ &
1050 & zeta(i,jend ,know))
1051# endif
1052 cff2=grid(ng)%on_v(i,jend+1)*ce
1053 tl_cff2=grid(ng)%on_v(i,jend+1)*tl_ce
1054!^ bry_val=vbar(i,Jend,know)+ &
1055!^ & cff2*(bry_pgr+ &
1056!^ & bry_cor+ &
1057!^ & bry_str)
1058!^
1059 tl_bry_val=tl_vbar(i,jend,know)+ &
1060 & tl_cff2*(bry_pgr+ &
1061 & bry_cor+ &
1062 & bry_str)+ &
1063 & cff2*(tl_bry_pgr+ &
1064 & tl_bry_cor+ &
1065 & tl_bry_str)- &
1066# ifdef TL_IOMS
1067 & cff2*(bry_pgr+bry_cor+bry_str)
1068# endif
1069# else
1070!^ bry_val=BOUNDARY(ng)%vbar_north(i)
1071!^
1072 tl_bry_val=boundary(ng)%tl_vbar_north(i)
1073# endif
1074# ifdef WET_DRY_NOT_YET
1075 cff=0.5_r8*(grid(ng)%h(i,jend )+ &
1076 & zeta(i,jend ,know)+ &
1077 & grid(ng)%h(i,jend+1)+ &
1078 & zeta(i,jend+1,know))
1079 tl_cff=0.5_r8*(grid(ng)%tl_h(i,jend )+ &
1080 & tl_zeta(i,jend ,know)+ &
1081 & grid(ng)%tl_h(i,jend+1)+ &
1082 & tl_zeta(i,jend+1,know))
1083# else
1084
1085 cff=0.5_r8*(grid(ng)%h(i,jend )+ &
1086 & grid(ng)%h(i,jend+1))
1087 tl_cff=0.5_r8*(grid(ng)%tl_h(i,jend )+ &
1088 & grid(ng)%tl_h(i,jend+1))
1089# endif
1090 cff1=sqrt(g/cff)
1091 tl_cff1=-0.5_r8*cff1*tl_cff/cff+ &
1092# ifdef TL_IOMS
1093 & 0.5_r8*cff1
1094# endif
1095 ce=dt2d*cff1*cff*0.5_r8*(grid(ng)%pn(i,jend )+ &
1096 & grid(ng)%pn(i,jend+1))
1097 tl_ce=dt2d*0.5_r8*(grid(ng)%pn(i,jend )+ &
1098 & grid(ng)%pn(i,jend+1))* &
1099 & (cff1*tl_cff+ &
1100 & tl_cff1*cff)- &
1101# ifdef TL_IOMS
1102 & ce
1103# endif
1104 ze=(0.5_r8+ce)*zeta(i,jend ,know)+ &
1105 & (0.5_r8-ce)*zeta(i,jend+1,know)
1106 tl_ze=(0.5_r8+ce)*tl_zeta(i,jend ,know)+ &
1107 & (0.5_r8-ce)*tl_zeta(i,jend+1,know)+ &
1108 & tl_ce*(zeta(i,jend ,know)- &
1109 & zeta(i,jend+1,know))- &
1110# ifdef TL_IOMS
1111 & ze ! HGA check
1112# endif
1113 IF (ce.gt.co) THEN
1114 cff2=(1.0_r8-co/ce)**2
1115 tl_cff2=2.0_r8*cff2*co*tl_ce/(ce*ce)- &
1116# ifdef TL_IOMS
1117 & cff2 ! HGA check
1118# endif
1119 cff3=zeta(i,jend,kout)+ &
1120 & ce*zeta(i,jend+1,know)- &
1121 & (1.0_r8+ce)*zeta(i,jend,know)
1122 tl_cff3=tl_zeta(i,jend,kout)+ &
1123 & ce*tl_zeta(i,jend+1,know)+ &
1124 & tl_ce*(zeta(i,jend ,know)+ &
1125 & zeta(i,jend+1,know))- &
1126 & (1.0_r8+ce)*tl_zeta(i,jend,know)- &
1127# ifdef TL_IOMS
1128 & ce*zeta(i,jend+1,know)+ &
1129 & (1.0_r8+ce)*zeta(i,jend,know) ! HGA check
1130# endif
1131 ze=ze+cff2*cff3
1132 tl_ze=tl_ze+cff2*tl_cff3+ &
1133 & tl_cff2*cff3- &
1134# ifdef TL_IOMS
1135 & cff2*cff3 ! HGA check
1136# endif
1137 END IF
1138!^ vbar(i,Jend+1,kout)=0.5_r8* &
1139!^ & ((1.0_r8-Ce)*vbar(i,Jend+1,know)+ &
1140!^ & Ce*vbar(i,Jend,know)+ &
1141!^ & bry_val+ &
1142!^ & cff1*(Ze-BOUNDARY(ng)%zeta_north(i)))
1143!^
1144 tl_vbar(i,jend+1,kout)=0.5_r8* &
1145 & ((1.0_r8-ce)* &
1146 & tl_vbar(i,jend+1,know)+ &
1147 & tl_ce*(vbar(i,jend ,know)- &
1148 & vbar(i,jend+1,know))+ &
1149 & ce*tl_vbar(i,jend,know)+ &
1150 & tl_bry_val+ &
1151 & tl_cff1* &
1152 & (ze-boundary(ng)%zeta_north(i))- &
1153 & cff1*tl_ze)- &
1154# ifdef TL_IOMS
1155 & 0.5_r8* &
1156 & ((1.0_r8-ce)*vbar(i,jend+1,know)+ &
1157 & ce*vbar(i,jend,know)+ &
1158 & cff1* &
1159 & (ze-boundary(ng)%zeta_north(i)))
1160!! ! HGA check
1161# endif
1162# ifdef ADJUST_BOUNDARY
1163 IF (lobc(inorth,isvbar,ng)) THEN
1164 tl_vbar(i,jend+1,kout)=tl_vbar(i,jend+1,kout)- &
1165 & 0.5_r8*cff1* &
1166 & ce*boundary(ng)%tl_zeta_north(i)
1167 END IF
1168# endif
1169# ifdef MASKING
1170!^ vbar(i,Jend+1,kout)=vbar(i,Jend+1,kout)* &
1171!^ & GRID(ng)%vmask(i,Jend+1)
1172!^
1173 tl_vbar(i,jend+1,kout)=tl_vbar(i,jend+1,kout)* &
1174 & grid(ng)%vmask(i,jend+1)
1175# endif
1176 END IF
1177 END DO
1178!
1179! Northern edge, clamped boundary condition.
1180!
1181 ELSE IF (tl_lbc(inorth,isvbar,ng)%clamped) THEN
1182 DO i=istr,iend
1183 IF (lbc_apply(ng)%north(i)) THEN
1184!^ vbar(i,Jend+1,kout)=BOUNDARY(ng)%vbar_north(i)
1185# ifdef MASKING
1186!^ vbar(i,Jend+1,kout)=vbar(i,Jend+1,kout)* &
1187!^ & GRID(ng)%vmask(i,Jend+1)
1188# endif
1189!^
1190 tl_vbar(i,jend+1,kout)=boundary(ng)%tl_vbar_north(i)
1191# ifdef MASKING
1192 tl_vbar(i,jend+1,kout)=tl_vbar(i,jend+1,kout)* &
1193 & grid(ng)%vmask(i,jend+1)
1194# endif
1195 END IF
1196 END DO
1197!
1198! Northern edge, gradient boundary condition.
1199!
1200 ELSE IF (tl_lbc(inorth,isvbar,ng)%gradient) THEN
1201 DO i=istr,iend
1202 IF (lbc_apply(ng)%north(i)) THEN
1203!^ vbar(i,Jend+1,kout)=vbar(i,Jend,kout)
1204!^
1205 tl_vbar(i,jend+1,kout)=tl_vbar(i,jend,kout)
1206# ifdef MASKING
1207!^ vbar(i,Jend+1,kout)=vbar(i,Jend+1,kout)* &
1208!^ & GRID(ng)%vmask(i,Jend+1)
1209!^
1210 tl_vbar(i,jend+1,kout)=tl_vbar(i,jend+1,kout)* &
1211 & grid(ng)%vmask(i,jend+1)
1212# endif
1213 END IF
1214 END DO
1215!
1216! Northern edge, reduced-physics boundary condition.
1217!
1218 ELSE IF (tl_lbc(inorth,isvbar,ng)%reduced) THEN
1219 DO i=istr,iend
1220 IF (lbc_apply(ng)%north(i)) THEN
1221 IF (tl_lbc(inorth,isfsur,ng)%acquire) THEN
1222!^ bry_pgr=-g*(BOUNDARY(ng)%zeta_north(i)- &
1223!^ & zeta(i,Jend,know))* &
1224!^ & 0.5_r8*GRID(ng)%pn(i,Jend)
1225!^
1226 tl_bry_pgr=-g*(boundary(ng)%tl_zeta_north(i)- &
1227 & tl_zeta(i,jend,know))* &
1228 & 0.5_r8*grid(ng)%pn(i,jend)
1229 ELSE
1230!^ bry_pgr=-g*(zeta(i,Jend+1,know)- &
1231!^ & zeta(i,Jend ,know))* &
1232!^ & 0.5_r8*(GRID(ng)%pn(i,Jend )+ &
1233!^ & GRID(ng)%pn(i,Jend+1))
1234!^
1235 tl_bry_pgr=-g*(tl_zeta(i,jend+1,know)- &
1236 & tl_zeta(i,jend ,know))* &
1237 & 0.5_r8*(grid(ng)%pn(i,jend )+ &
1238 & grid(ng)%pn(i,jend+1))
1239 END IF
1240# ifdef UV_COR
1241!^ bry_cor=-0.125_r8*(ubar(i ,Jend ,know)+ &
1242!^ & ubar(i+1,Jend ,know)+ &
1243!^ & ubar(i ,Jend+1,know)+ &
1244!^ & ubar(i+1,Jend+1,know))* &
1245!^ & (GRID(ng)%f(i,Jend )+ &
1246!^ & GRID(ng)%f(i,Jend+1))
1247!^
1248 tl_bry_cor=-0.125_r8*(tl_ubar(i ,jend ,know)+ &
1249 & tl_ubar(i+1,jend ,know)+ &
1250 & tl_ubar(i ,jend+1,know)+ &
1251 & tl_ubar(i+1,jend+1,know))* &
1252 & (grid(ng)%f(i,jend )+ &
1253 & grid(ng)%f(i,jend+1))
1254# else
1255!^ bry_cor=0.0_r8
1256!^
1257 tl_bry_cor=0.0_r8
1258# endif
1259 cff=1.0_r8/(0.5_r8*(grid(ng)%h(i,jend )+ &
1260 & zeta(i,jend ,know)+ &
1261 & grid(ng)%h(i,jend+1)+ &
1262 & zeta(i,jend+1,know)))
1263 tl_cff=-cff*cff*0.5_r8*(grid(ng)%tl_h(i,jend )+ &
1264 & tl_zeta(i,jend ,know)+ &
1265 & grid(ng)%tl_h(i,jend+1)+ &
1266 & tl_zeta(i,jend+1,know))+ &
1267# ifdef TL_IOMS
1268 & 2.0_r8*cff
1269# endif
1270!^ bry_str=cff*(FORCES(ng)%svstr(i,Jend+1)- &
1271!^ & FORCES(ng)%bvstr(i,Jend+1))
1272!^
1273 tl_bry_str=tl_cff*(forces(ng)%svstr(i,jend+1)- &
1274 & forces(ng)%bvstr(i,jend+1))+ &
1275 & cff*(forces(ng)%tl_svstr(i,jend+1)- &
1276 & forces(ng)%tl_bvstr(i,jend+1))- &
1277# ifdef TL_IOMS
1278 & cff*(forces(ng)%svstr(i,jend+1)- &
1279 & forces(ng)%bvstr(i,jend+1))
1280# endif
1281!^ vbar(i,Jend+1,kout)=vbar(i,Jend+1,know)+ &
1282!^ & dt2d*(bry_pgr+ &
1283!^ & bry_cor+ &
1284!^ & bry_str)
1285!^
1286 tl_vbar(i,jend+1,kout)=tl_vbar(i,jend+1,know)+ &
1287 & dt2d*(tl_bry_pgr+ &
1288 & tl_bry_cor+ &
1289 & tl_bry_str)
1290# ifdef MASKING
1291!^ vbar(i,Jend+1,kout)=vbar(i,Jend+1,kout)* &
1292!^ & GRID(ng)%vmask(i,Jend+1)
1293!^
1294 tl_vbar(i,jend+1,kout)=tl_vbar(i,jend+1,kout)* &
1295 & grid(ng)%vmask(i,jend+1)
1296# endif
1297 END IF
1298 END DO
1299!
1300! Northern edge, closed boundary condition.
1301!
1302 ELSE IF (tl_lbc(inorth,isvbar,ng)%closed) THEN
1303 DO i=istr,iend
1304 IF (lbc_apply(ng)%north(i)) THEN
1305!^ vbar(i,Jend+1,kout)=0.0_r8
1306!^
1307 tl_vbar(i,jend+1,kout)=0.0_r8
1308 END IF
1309 END DO
1310 END IF
1311 END IF
1312!
1313!-----------------------------------------------------------------------
1314! Lateral boundary conditions at the western edge.
1315!-----------------------------------------------------------------------
1316!
1317 IF (domain(ng)%Western_Edge(tile)) THEN
1318!
1319! Western edge, implicit upstream radiation condition.
1320!
1321 IF (tl_lbc(iwest,isvbar,ng)%radiation) THEN
1322 IF (iic(ng).ne.0) THEN
1323 DO j=jstrv-1,jend
1324!^ grad(Istr-1,j)=vbar(Istr-1,j+1,know)- &
1325!^ & vbar(Istr-1,j ,know)
1326!^
1327 tl_grad(istr-1,j)=0.0_r8
1328 END DO
1329 DO j=jstrv,jend
1330 IF (lbc_apply(ng)%west(j)) THEN
1331# if defined CELERITY_READ && defined FORWARD_READ
1332 IF (tl_lbc(iwest,isvbar,ng)%nudging) THEN
1333 IF (lnudgem2clm(ng)) THEN
1334 obc_out=0.5_r8* &
1335 & (clima(ng)%M2nudgcof(istr-1,j-1)+ &
1336 & clima(ng)%M2nudgcof(istr-1,j ))
1337 obc_in =obcfac(ng)*obc_out
1338 ELSE
1339 obc_out=m2obc_out(ng,iwest)
1340 obc_in =m2obc_in(ng,iwest)
1341 END IF
1342 IF (boundary(ng)%vbar_west_Cx(j).lt.0.0_r8) THEN
1343 tau=obc_in
1344 ELSE
1345 tau=obc_out
1346 END IF
1347 tau=tau*dt2d
1348 END IF
1349 cx=boundary(ng)%vbar_west_Cx(j)
1350# ifdef RADIATION_2D
1351 ce=boundary(ng)%vbar_west_Ce(j)
1352# else
1353 ce=0.0_r8
1354# endif
1355 cff=boundary(ng)%vbar_west_C2(j)
1356# endif
1357!^ vbar(Istr-1,j,kout)=(cff*vbar(Istr-1,j,know)+ &
1358!^ & Cx *vbar(Istr ,j,kout)- &
1359!^ & MAX(Ce,0.0_r8)*grad(Istr-1,j-1)- &
1360!^ & MIN(Ce,0.0_r8)*grad(Istr-1,j ))/ &
1361!^ & (cff+Cx)
1362!^
1363 tl_vbar(istr-1,j,kout)=(cff*tl_vbar(istr-1,j,know)+ &
1364 & cx *tl_vbar(1,j,kout)- &
1365 & max(ce,0.0_r8)* &
1366 & tl_grad(istr-1,j-1)- &
1367 & min(ce,0.0_r8)* &
1368 & tl_grad(istr-1,j ))/ &
1369 & (cff+cx)
1370
1371 IF (tl_lbc(iwest,isvbar,ng)%nudging) THEN
1372!^ vbar(Istr-1,j,kout)=vbar(Istr-1,j,kout)+ &
1373!^ & tau*(BOUNDARY(ng)%vbar_west(j)- &
1374!^ & vbar(Istr-1,j,know))
1375!^
1376 tl_vbar(istr-1,j,kout)=tl_vbar(istr-1,j,kout)- &
1377 & tau*tl_vbar(1,j,know)
1378 END IF
1379# ifdef MASKING
1380!^ vbar(Istr-1,j,kout)=vbar(Istr-1,j,kout)* &
1381!^ & GRID(ng)%vmask(Istr-1,j)
1382!^
1383 tl_vbar(istr-1,j,kout)=tl_vbar(istr-1,j,kout)* &
1384 & grid(ng)%vmask(istr-1,j)
1385# endif
1386 END IF
1387 END DO
1388 END IF
1389!
1390! Western edge, Chapman boundary condition.
1391!
1392 ELSE IF (tl_lbc(iwest,isvbar,ng)%Flather.or. &
1393 & tl_lbc(iwest,isvbar,ng)%reduced.or. &
1394 & tl_lbc(iwest,isvbar,ng)%Shchepetkin) THEN
1395 DO j=jstrv,jend
1396 IF (lbc_apply(ng)%west(j)) THEN
1397 cff=dt2d*0.5_r8*(grid(ng)%pm(istr,j-1)+ &
1398 & grid(ng)%pm(istr,j ))
1399 cff1=sqrt(g*0.5_r8*(grid(ng)%h(istr,j-1)+ &
1400 & zeta(istr,j-1,know)+ &
1401 & grid(ng)%h(istr,j )+ &
1402 & zeta(istr,j ,know)))
1403 tl_cff1=0.25_r8*g*(grid(ng)%tl_h(istr,j-1)+ &
1404 & tl_zeta(istr,j-1,know)+ &
1405 & grid(ng)%tl_h(istr,j )+ &
1406 & tl_zeta(istr,j ,know))/cff1+ &
1407# ifdef TL_IOMS
1408 & 0.5_r8*cff1
1409# endif
1410 cx=cff*cff1
1411 tl_cx=cff*tl_cff1
1412 cff2=1.0_r8/(1.0_r8+cx)
1413 tl_cff2=-cff2*cff2*tl_cx+ &
1414# ifdef TL_IOMS
1415 & cff2*cff2*(1.0_r8+2.0_r8*cx)
1416# endif
1417!^ vbar(Istr-1,j,kout)=cff2*(vbar(Istr-1,j,know)+ &
1418!^ & Cx*vbar(Istr,j,kout))
1419!^
1420 tl_vbar(istr-1,j,kout)=tl_cff2*(vbar(istr-1,j,know)+ &
1421 & cx*vbar(istr,j,kout))+ &
1422 & cff2*(tl_vbar(istr-1,j,know)+ &
1423 & tl_cx*vbar(istr,j,kout)+ &
1424 & cx*tl_vbar(istr,j,kout))- &
1425# ifdef TL_IOMS
1426 & cff2*(vbar(istr-1,j,know)+ &
1427 & 2.0_r8*cx*vbar(istr,j,kout))
1428# endif
1429# ifdef MASKING
1430!^ vbar(Istr-1,j,kout)=vbar(Istr-1,j,kout)* &
1431!^ & GRID(ng)%vmask(Istr-1,j)
1432!^
1433 tl_vbar(istr-1,j,kout)=tl_vbar(istr-1,j,kout)* &
1434 & grid(ng)%vmask(istr-1,j)
1435# endif
1436 END IF
1437 END DO
1438!
1439! Western edge, clamped boundary condition.
1440!
1441 ELSE IF (tl_lbc(iwest,isvbar,ng)%clamped) THEN
1442 DO j=jstrv,jend
1443 IF (lbc_apply(ng)%west(j)) THEN
1444!^ vbar(Istr-1,j,kout)=BOUNDARY(ng)%vbar_west(j)
1445!^
1446 tl_vbar(istr-1,j,kout)=boundary(ng)%tl_vbar_west(j)
1447# ifdef MASKING
1448!^ vbar(Istr-1,j,kout)=vbar(Istr-1,j,kout)* &
1449!^ & GRID(ng)%vmask(Istr-1,j)
1450!^
1451 tl_vbar(istr-1,j,kout)=tl_vbar(istr-1,j,kout)* &
1452 & grid(ng)%vmask(istr-1,j)
1453# endif
1454 END IF
1455 END DO
1456!
1457! Western edge, gradient boundary condition.
1458!
1459 ELSE IF (tl_lbc(iwest,isvbar,ng)%gradient) THEN
1460 DO j=jstrv,jend
1461 IF (lbc_apply(ng)%west(j)) THEN
1462!^ vbar(Istr-1,j,kout)=vbar(Istr,j,kout)
1463!^
1464 tl_vbar(istr-1,j,kout)=tl_vbar(istr,j,kout)
1465# ifdef MASKING
1466!^ vbar(Istr-1,j,kout)=vbar(Istr-1,j,kout)* &
1467!^ & GRID(ng)%vmask(Istr-1,j)
1468!^
1469 tl_vbar(istr-1,j,kout)=tl_vbar(istr-1,j,kout)* &
1470 & grid(ng)%vmask(istr-1,j)
1471# endif
1472 END IF
1473 END DO
1474!
1475! Western edge, closed boundary condition: free slip (gamma2=1) or
1476! no slip (gamma2=-1).
1477!
1478 ELSE IF (tl_lbc(iwest,isvbar,ng)%closed) THEN
1479 IF (nsperiodic(ng)) THEN
1480 jmin=jstrv
1481 jmax=jend
1482 ELSE
1483 jmin=jstr
1484 jmax=jendr
1485 END IF
1486 DO j=jmin,jmax
1487 IF (lbc_apply(ng)%west(j)) THEN
1488!^ vbar(Istr-1,j,kout)=gamma2(ng)*vbar(Istr,j,kout)
1489!^
1490 tl_vbar(istr-1,j,kout)=gamma2(ng)*tl_vbar(istr,j,kout)
1491# ifdef MASKING
1492!^ vbar(Istr-1,j,kout)=vbar(Istr-1,j,kout)* &
1493!^ & GRID(ng)%vmask(Istr-1,j)
1494!^
1495 tl_vbar(istr-1,j,kout)=tl_vbar(istr-1,j,kout)* &
1496 & grid(ng)%vmask(istr-1,j)
1497# endif
1498 END IF
1499 END DO
1500 END IF
1501 END IF
1502!
1503!-----------------------------------------------------------------------
1504! Lateral boundary conditions at the eastern edge.
1505!-----------------------------------------------------------------------
1506!
1507 IF (domain(ng)%Eastern_Edge(tile)) THEN
1508!
1509! Eastern edge, implicit upstream radiation condition.
1510!
1511 IF (tl_lbc(ieast,isvbar,ng)%radiation) THEN
1512 IF (iic(ng).ne.0) THEN
1513 DO j=jstrv-1,jend
1514!^ grad(Iend+1,j)=vbar(Iend+1,j+1,know)- &
1515!^ & vbar(Iend+1,j ,know)
1516!^
1517 tl_grad(iend+1,j)=0.0_r8
1518 END DO
1519 DO j=jstrv,jend
1520 IF (lbc_apply(ng)%east(j)) THEN
1521# if defined CELERITY_READ && defined FORWARD_READ
1522 IF (tl_lbc(ieast,isvbar,ng)%nudging) THEN
1523 IF (lnudgem2clm(ng)) THEN
1524 obc_out=0.5_r8* &
1525 & (clima(ng)%M2nudgcof(iend+1,j-1)+ &
1526 & clima(ng)%M2nudgcof(iend+1,j ))
1527 obc_in =obcfac(ng)*obc_out
1528 ELSE
1529 obc_out=m2obc_out(ng,ieast)
1530 obc_in =m2obc_in(ng,ieast)
1531 END IF
1532 IF (boundary(ng)%vbar_east_Cx(j).lt.0.0_r8) THEN
1533 tau=obc_in
1534 ELSE
1535 tau=obc_out
1536 END IF
1537 tau=tau*dt2d
1538 END IF
1539 cx=boundary(ng)%vbar_east_Cx(j)
1540# ifdef RADIATION_2D
1541 ce=boundary(ng)%vbar_east_Ce(j)
1542# else
1543 ce=0.0_r8
1544# endif
1545 cff=boundary(ng)%vbar_east_C2(j)
1546# endif
1547!^ vbar(Iend+1,j,kout)=(cff*vbar(Iend+1,j,know)+ &
1548!^ & Cx *vbar(Iend ,j,kout)- &
1549!^ & MAX(Ce,0.0_r8)*grad(Iend+1,j-1)- &
1550!^ & MIN(Ce,0.0_r8)*grad(Iend+1,j ))/ &
1551!^ & (cff+Cx)
1552!^
1553 tl_vbar(iend+1,j,kout)=(cff*tl_vbar(iend+1,j,know)+ &
1554 & cx *tl_vbar(iend ,j,kout)- &
1555 & max(ce,0.0_r8)* &
1556 & tl_grad(iend+1,j-1)- &
1557 & min(ce,0.0_r8)* &
1558 & tl_grad(iend+1,j ))/ &
1559 & (cff+cx)
1560
1561 IF (tl_lbc(ieast,isvbar,ng)%nudging) THEN
1562!^ vbar(Iend+1,j,kout)=vbar(Iend+1,j,kout)+ &
1563!^ & tau*(BOUNDARY(ng)%vbar_east(j)- &
1564!^ & vbar(Iend+1,j,know))
1565!^
1566 tl_vbar(iend+1,j,kout)=tl_vbar(iend+1,j,kout)- &
1567 & tau*tl_vbar(iend+1,j,know)
1568 END IF
1569# ifdef MASKING
1570!^ vbar(Iend+1,j,kout)=vbar(Iend+1,j,kout)* &
1571!^ & GRID(ng)%vmask(Iend+1,j)
1572!^
1573 tl_vbar(iend+1,j,kout)=tl_vbar(iend+1,j,kout)* &
1574 & grid(ng)%vmask(iend+1,j)
1575# endif
1576 END IF
1577 END DO
1578 END IF
1579!
1580! Eastern edge, Chapman boundary condition.
1581!
1582 ELSE IF (tl_lbc(ieast,isvbar,ng)%Flather.or. &
1583 & tl_lbc(ieast,isvbar,ng)%reduced.or. &
1584 & tl_lbc(ieast,isvbar,ng)%Shchepetkin) THEN
1585 DO j=jstrv,jend
1586 IF (lbc_apply(ng)%east(j)) THEN
1587 cff=dt2d*0.5_r8*(grid(ng)%pm(iend,j-1)+ &
1588 & grid(ng)%pm(iend,j ))
1589 cff1=sqrt(g*0.5_r8*(grid(ng)%h(iend,j-1)+ &
1590 & zeta(iend,j-1,know)+ &
1591 & grid(ng)%h(iend,j )+ &
1592 & zeta(iend,j ,know)))
1593 tl_cff1=0.25_r8*g*(grid(ng)%tl_h(iend,j-1)+ &
1594 & tl_zeta(iend,j-1,know)+ &
1595 & grid(ng)%tl_h(iend,j )+ &
1596 & tl_zeta(iend,j ,know))/cff1+ &
1597# ifdef TL_IOMS
1598 & 0.5_r8*cff1
1599# endif
1600 cx=cff*cff1
1601 tl_cx=cff*tl_cff1
1602 cff2=1.0_r8/(1.0_r8+cx)
1603 tl_cff2=-cff2*cff2*tl_cx+ &
1604# ifdef TL_IOMS
1605 & cff2*cff2*(1.0_r8+2.0_r8*cx)
1606# endif
1607!^ vbar(Iend+1,j,kout)=cff2*(vbar(Iend+1,j,know)+ &
1608!^ & Cx*vbar(Iend,j,kout))
1609!^
1610 tl_vbar(iend+1,j,kout)=tl_cff2*(vbar(iend+1,j,know)+ &
1611 & cx*vbar(iend,j,kout))+ &
1612 & cff2*(tl_vbar(iend+1,j,know)+ &
1613 & tl_cx*vbar(iend,j,kout)+ &
1614 & cx*tl_vbar(iend,j,kout))- &
1615# ifdef TL_IOMS
1616 & cff2*(vbar(iend+1,j,know)+ &
1617 & 2.0_r8*cx*vbar(iend,j,kout))
1618# endif
1619# ifdef MASKING
1620!^ vbar(Iend+1,j,kout)=vbar(Iend+1,j,kout)* &
1621!^ & GRID(ng)%vmask(Iend+1,j)
1622!^
1623 tl_vbar(iend+1,j,kout)=tl_vbar(iend+1,j,kout)* &
1624 & grid(ng)%vmask(iend+1,j)
1625# endif
1626 END IF
1627 END DO
1628!
1629! Eastern edge, clamped boundary condition.
1630!
1631 ELSE IF (tl_lbc(ieast,isvbar,ng)%clamped) THEN
1632 DO j=jstrv,jend
1633 IF (lbc_apply(ng)%east(j)) THEN
1634!^ vbar(Iend+1,j,kout)=BOUNDARY(ng)%vbar_east(j)
1635!^
1636 tl_vbar(iend+1,j,kout)=boundary(ng)%tl_vbar_east(j)
1637# ifdef MASKING
1638!^ vbar(Iend+1,j,kout)=vbar(Iend+1,j,kout)* &
1639!^ & GRID(ng)%vmask(Iend+1,j)
1640!^
1641 tl_vbar(iend+1,j,kout)=tl_vbar(iend+1,j,kout)* &
1642 & grid(ng)%vmask(iend+1,j)
1643# endif
1644 END IF
1645 END DO
1646!
1647! Eastern edge, gradient boundary condition.
1648!
1649 ELSE IF (tl_lbc(ieast,isvbar,ng)%gradient) THEN
1650 DO j=jstrv,jend
1651 IF (lbc_apply(ng)%east(j)) THEN
1652!^ vbar(Iend+1,j,kout)=vbar(Iend,j,kout)
1653!^
1654 tl_vbar(iend+1,j,kout)=tl_vbar(iend,j,kout)
1655# ifdef MASKING
1656!^ vbar(Iend+1,j,kout)=vbar(Iend+1,j,kout)* &
1657!^ & GRID(ng)%vmask(Iend+1,j)
1658!^
1659 tl_vbar(iend+1,j,kout)=tl_vbar(iend+1,j,kout)* &
1660 & grid(ng)%vmask(iend+1,j)
1661# endif
1662 END IF
1663 END DO
1664!
1665! Eastern edge, closed boundary condition: free slip (gamma2=1) or
1666! no slip (gamma2=-1).
1667!
1668 ELSE IF (tl_lbc(ieast,isvbar,ng)%closed) THEN
1669 IF (nsperiodic(ng)) THEN
1670 jmin=jstrv
1671 jmax=jend
1672 ELSE
1673 jmin=jstr
1674 jmax=jendr
1675 END IF
1676 DO j=jmin,jmax
1677 IF (lbc_apply(ng)%east(j)) THEN
1678!^ vbar(Iend+1,j,kout)=gamma2(ng)*vbar(Iend,j,kout)
1679!^
1680 tl_vbar(iend+1,j,kout)=gamma2(ng)*tl_vbar(iend,j,kout)
1681# ifdef MASKING
1682!^ vbar(Iend+1,j,kout)=vbar(Iend+1,j,kout)* &
1683!^ & GRID(ng)%vmask(Iend+1,j)
1684!^
1685 tl_vbar(iend+1,j,kout)=tl_vbar(iend+1,j,kout)* &
1686 & grid(ng)%vmask(iend+1,j)
1687# endif
1688 END IF
1689 END DO
1690 END IF
1691 END IF
1692!
1693!-----------------------------------------------------------------------
1694! Boundary corners.
1695!-----------------------------------------------------------------------
1696!
1697 IF (.not.(ewperiodic(ng).or.nsperiodic(ng))) THEN
1698 IF (domain(ng)%SouthWest_Corner(tile)) THEN
1699 IF (lbc_apply(ng)%south(istr-1).and. &
1700 & lbc_apply(ng)%west (jstr )) THEN
1701!^ vbar(Istr-1,Jstr,kout)=0.5_r8*(vbar(Istr ,Jstr ,kout)+ &
1702!^ & vbar(Istr-1,Jstr+1,kout))
1703!^
1704 tl_vbar(istr-1,jstr,kout)=0.5_r8* &
1705 & (tl_vbar(istr ,jstr ,kout)+ &
1706 & tl_vbar(istr-1,jstr+1,kout))
1707 END IF
1708 END IF
1709 IF (domain(ng)%SouthEast_Corner(tile)) THEN
1710 IF (lbc_apply(ng)%south(iend+1).and. &
1711 & lbc_apply(ng)%east (jstr )) THEN
1712!^ vbar(Iend+1,Jstr,kout)=0.5_r8*(vbar(Iend ,Jstr ,kout)+ &
1713!^ & vbar(Iend+1,Jstr+1,kout))
1714!^
1715 tl_vbar(iend+1,jstr,kout)=0.5_r8* &
1716 & (tl_vbar(iend ,jstr ,kout)+ &
1717 & tl_vbar(iend+1,jstr+1,kout))
1718 END IF
1719 END IF
1720 IF (domain(ng)%NorthWest_Corner(tile)) THEN
1721 IF (lbc_apply(ng)%north(istr-1).and. &
1722 & lbc_apply(ng)%west (jend+1)) THEN
1723!^ vbar(Istr-1,Jend+1,kout)=0.5_r8*(vbar(Istr-1,Jend ,kout)+ &
1724!^ & vbar(Istr ,Jend+1,kout))
1725!^
1726 tl_vbar(istr-1,jend+1,kout)=0.5_r8* &
1727 & (tl_vbar(istr-1,jend ,kout)+ &
1728 & tl_vbar(istr ,jend+1,kout))
1729 END IF
1730 END IF
1731 IF (domain(ng)%NorthEast_Corner(tile)) THEN
1732 IF (lbc_apply(ng)%north(iend+1).and. &
1733 & lbc_apply(ng)%east (jend+1)) THEN
1734!^ vbar(Iend+1,Jend+1,kout)=0.5_r8*(vbar(Iend+1,Jend ,kout)+ &
1735!^ & vbar(Iend ,Jend+1,kout))
1736!^
1737 tl_vbar(iend+1,jend+1,kout)=0.5_r8* &
1738 & (tl_vbar(iend+1,jend ,kout)+ &
1739 & tl_vbar(iend ,jend+1,kout))
1740 END IF
1741 END IF
1742 END IF
1743
1744# if defined WET_DRY_NOT_YET
1745!
1746!-----------------------------------------------------------------------
1747! Impose wetting and drying conditions.
1748!
1749! HGA: need RPM code here.
1750!-----------------------------------------------------------------------
1751!
1752 IF (.not.ewperiodic(ng)) THEN
1753 IF (domain(ng)%Western_Edge(tile)) THEN
1754 DO j=jstrv,jend
1755 IF (lbc_apply(ng)%west(j).or. &
1756 & lbc(iwest,isvbar,ng)%nested) THEN
1757!^ cff1=ABS(ABS(GRID(ng)%vmask_wet(Istr-1,j))-1.0_r8)
1758!^ cff2=0.5_r8+DSIGN(0.5_r8,vbar(Istr-1,j,kout))* &
1759!^ & GRID(ng)%vmask_wet(Istr-1,j)
1760!^ cff=0.5_r8*GRID(ng)%vmask_wet(Istr-1,j)*cff1+ &
1761!^ & cff2*(1.0_r8-cff1)
1762!^ vbar(Istr,j,kout)=vbar(Istr,j,kout)*cff
1763 END IF
1764 END DO
1765 END IF
1766 IF (domain(ng)%Eastern_Edge(tile)) THEN
1767 DO j=jstrv,jend
1768 IF (lbc_apply(ng)%east(j).or. &
1769 & lbc(ieast,isvbar,ng)%nested) THEN
1770!^ cff1=ABS(ABS(GRID(ng)%vmask_wet(Iend+1,j))-1.0_r8)
1771!^ cff2=0.5_r8+DSIGN(0.5_r8,vbar(Iend+1,j,kout))* &
1772!^ & GRID(ng)%vmask_wet(Iend+1,j)
1773!^ cff=0.5_r8*GRID(ng)%vmask_wet(Iend+1,j)*cff1+ &
1774!^ & cff2*(1.0_r8-cff1)
1775!^ vbar(Iend+1,j,kout)=vbar(Iend+1,j,kout)*cff
1776 END IF
1777 END DO
1778 END IF
1779 END IF
1780!
1781 IF (.not.nsperiodic(ng)) THEN
1782 IF (domain(ng)%Southern_Edge(tile)) THEN
1783 DO i=istr,iend
1784 IF (lbc_apply(ng)%south(i).or. &
1785 & lbc(isouth,isvbar,ng)%nested) THEN
1786!^ cff1=ABS(ABS(GRID(ng)%vmask_wet(i,Jstr))-1.0_r8)
1787!^ cff2=0.5_r8+DSIGN(0.5_r8,vbar(i,Jstr,kout))* &
1788!^ & GRID(ng)%vmask_wet(i,Jstr)
1789!^ cff=0.5_r8*GRID(ng)%vmask_wet(i,Jstr)*cff1+ &
1790!^ & cff2*(1.0_r8-cff1)
1791!^ vbar(i,Jstr,kout)=vbar(i,Jstr,kout)*cff
1792 END IF
1793 END DO
1794 END IF
1795 IF (domain(ng)%Northern_Edge(tile)) THEN
1796 DO i=istr,iend
1797 IF (lbc_apply(ng)%north(i).or. &
1798 & lbc(inorth,isvbar,ng)%nested) THEN
1799!^ cff1=ABS(ABS(GRID(ng)%vmask_wet(i,Jend+1))-1.0_r8)
1800!^ cff2=0.5_r8+DSIGN(0.5_r8,vbar(i,Jend+1,kout))* &
1801!^ & GRID(ng)%vmask_wet(i,Jend+1)
1802!^ cff=0.5_r8*GRID(ng)%vmask_wet(i,Jend+1)*cff1+ &
1803!^ & cff2*(1.0_r8-cff1)
1804!^ vbar(i,Jend+1,kout)=vbar(i,Jend+1,kout)*cff
1805 END IF
1806 END DO
1807 END IF
1808 END IF
1809!
1810 IF (.not.(ewperiodic(ng).or.nsperiodic(ng))) THEN
1811 IF (domain(ng)%SouthWest_Corner(tile)) THEN
1812 IF ((lbc_apply(ng)%south(istr-1).and. &
1813 & lbc_apply(ng)%west (jstr )).or. &
1814 & (lbc(iwest,isvbar,ng)%nested.and. &
1815 & lbc(isouth,isvbar,ng)%nested)) THEN
1816!^ cff1=ABS(ABS(GRID(ng)%vmask_wet(Istr-1,Jstr))-1.0_r8)
1817!^ cff2=0.5_r8+DSIGN(0.5_r8,vbar(Istr-1,Jstr,kout))* &
1818!^ & GRID(ng)%vmask_wet(Istr-1,Jstr)
1819!^ cff=0.5_r8*GRID(ng)%vmask_wet(Istr-1,Jstr)*cff1+ &
1820!^ & cff2*(1.0_r8-cff1)
1821!^ vbar(Istr-1,Jstr,kout)=vbar(Istr-1,Jstr,kout)*cff
1822 END IF
1823 END IF
1824 IF (domain(ng)%SouthEast_Corner(tile)) THEN
1825 IF ((lbc_apply(ng)%south(iend+1).and. &
1826 & lbc_apply(ng)%east (jstr )).or. &
1827 & (lbc(ieast,isvbar,ng)%nested.and. &
1828 & lbc(isouth,isvbar,ng)%nested)) THEN
1829!^ cff1=ABS(ABS(GRID(ng)%vmask_wet(Iend+1,Jstr))-1.0_r8)
1830!^ cff2=0.5_r8+DSIGN(0.5_r8,vbar(Iend+1,Jstr,kout))* &
1831!^ & GRID(ng)%vmask_wet(Iend+1,Jstr)
1832!^ cff=0.5_r8*GRID(ng)%vmask_wet(Iend+1,Jstr)*cff1+ &
1833!^ & cff2*(1.0_r8-cff1)
1834!^ vbar(Iend+1,Jstr,kout)=vbar(Iend+1,Jstr,kout)*cff
1835 END IF
1836 END IF
1837 IF (domain(ng)%NorthWest_Corner(tile)) THEN
1838 IF ((lbc_apply(ng)%north(istr-1).and. &
1839 & lbc_apply(ng)%west (jend+1)).or. &
1840 & (lbc(iwest,isvbar,ng)%nested.and. &
1841 & lbc(inorth,isvbar,ng)%nested)) THEN
1842!^ cff1=ABS(ABS(GRID(ng)%vmask_wet(Istr-1,Jend+1))-1.0_r8)
1843!^ cff2=0.5_r8+DSIGN(0.5_r8,vbar(Istr-1,Jend+1,kout))* &
1844!^ & GRID(ng)%vmask_wet(Istr-1,Jend+1)
1845!^ cff=0.5_r8*GRID(ng)%vmask_wet(Istr-1,Jend+1)*cff1+ &
1846!^ & cff2*(1.0_r8-cff1)
1847!^ vbar(Istr-1,Jend+1,kout)=vbar(Istr-1,Jend+1,kout)*cff
1848 END IF
1849 END IF
1850 IF (domain(ng)%NorthEast_Corner(tile)) THEN
1851 IF ((lbc_apply(ng)%north(iend+1).and. &
1852 & lbc_apply(ng)%east (jend+1)).or. &
1853 & (lbc(ieast,isvbar,ng)%nested.and. &
1854 & lbc(inorth,isvbar,ng)%nested)) THEN
1855!^ cff1=ABS(ABS(GRID(ng)%vmask_wet(Iend+1,Jend+1))-1.0_r8)
1856!^ cff2=0.5_r8+DSIGN(0.5_r8,vbar(Iend+1,Jend+1,kout))* &
1857!^ & GRID(ng)%vmask_wet(Iend+1,Jend+1)
1858!^ cff=0.5_r8*GRID(ng)%vmask_wet(Iend+1,Jend+1)*cff1+ &
1859!^ & cff2*(1.0_r8-cff1)
1860!^ vbar(Iend+1,Jend+1,kout)=vbar(Iend+1,Jend+1,kout)*cff
1861 END IF
1862 END IF
1863 END IF
1864# endif
1865!
1866 RETURN
1867 END SUBROUTINE rp_v2dbc_tile
1868#endif
1869 END MODULE rp_v2dbc_mod
type(t_boundary), dimension(:), allocatable boundary
type(t_apply), dimension(:), allocatable lbc_apply
type(t_clima), dimension(:), allocatable clima
Definition mod_clima.F:153
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
integer isvbar
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_lbc), dimension(:,:,:), allocatable lbc
Definition mod_param.F:375
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
real(r8) co
logical, dimension(:), allocatable lnudgem2clm
integer, dimension(:), allocatable iic
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 obcfac
real(r8), dimension(:), allocatable gamma2
integer, parameter isouth
real(dp), dimension(:), allocatable dtfast
real(dp), dimension(:,:), allocatable m2obc_out
integer, parameter ieast
real(dp) g
real(dp) rho0
integer, parameter inorth
real(dp), dimension(:,:), allocatable m2obc_in
integer, dimension(:), allocatable kstp
integer, dimension(:), allocatable krhs
subroutine, public rp_v2dbc(ng, tile, kout)
Definition rp_v2dbc_im.F:29
subroutine, public rp_v2dbc_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, krhs, kstp, kout, ubar, vbar, zeta, tl_ubar, tl_vbar, tl_zeta)
Definition rp_v2dbc_im.F:64