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