ROMS
Loading...
Searching...
No Matches
ad_v3dbc_im.F
Go to the documentation of this file.
1#include "cppdefs.h"
3#if defined ADJOINT && defined SOLVE3D
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 adjoint lateral boundary conditions for total !
13! 3D V-velocity. It updates the specified "nout" time index. !
14! !
15! BASIC STATE variables needed: v !
16! !
17!=======================================================================
18!
19 implicit none
20
21 PRIVATE
22 PUBLIC :: ad_v3dbc, ad_v3dbc_tile
23
24 CONTAINS
25!
26!***********************************************************************
27 SUBROUTINE ad_v3dbc (ng, tile, nout)
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, nout
37!
38! Local variable declarations.
39!
40# include "tile.h"
41!
42 CALL ad_v3dbc_tile (ng, tile, &
43 & lbi, ubi, lbj, ubj, n(ng), &
44 & imins, imaxs, jmins, jmaxs, &
45 & nstp(ng), nout, &
46 & ocean(ng) % ad_v)
47 RETURN
48 END SUBROUTINE ad_v3dbc
49
50!
51!***********************************************************************
52 SUBROUTINE ad_v3dbc_tile (ng, tile, &
53 & LBi, UBi, LBj, UBj, UBk, &
54 & IminS, ImaxS, JminS, JmaxS, &
55 & nstp, nout, &
56 & ad_v)
57!***********************************************************************
58!
59 USE mod_param
60 USE mod_boundary
61 USE mod_clima
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, ubk
70 integer, intent(in) :: imins, imaxs, jmins, jmaxs
71 integer, intent(in) :: nstp, nout
72!
73# ifdef ASSUMED_SHAPE
74 real(r8), intent(inout) :: ad_v(lbi:,lbj:,:,:)
75# else
76 real(r8), intent(inout) :: ad_v(lbi:ubi,lbj:ubj,ubk,2)
77# endif
78!
79! Local variable declarations.
80!
81 integer :: jmin, jmax
82 integer :: i, j, k
83
84 real(r8) :: ce, cx, cff
85 real(r8) :: obc_in, obc_out, tau
86 real(r8) :: adfac
87
88 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_grad
89
90# include "set_bounds.h"
91!
92!-----------------------------------------------------------------------
93! Initialize adjoint private variables.
94!-----------------------------------------------------------------------
95!
96 ad_grad(lbi:ubi,lbj:ubj)=0.0_r8
97!
98!-----------------------------------------------------------------------
99! Boundary corners.
100!-----------------------------------------------------------------------
101!
102 IF (.not.(ewperiodic(ng).or.nsperiodic(ng))) THEN
103 IF (domain(ng)%NorthEast_Corner(tile)) THEN
104 IF (lbc_apply(ng)%north(iend+1).and. &
105 & lbc_apply(ng)%east (jend+1)) THEN
106 DO k=1,n(ng)
107!^ tl_v(Iend+1,Jend+1,k,nout)=0.5_r8* &
108!^ & (tl_v(Iend+1,Jend ,k,nout)+ &
109!^ & tl_v(Iend ,Jend+1,k,nout))
110!^
111 adfac=0.5_r8*ad_v(iend+1,jend+1,k,nout)
112 ad_v(iend+1,jend ,k,nout)=ad_v(iend+1,jend ,k,nout)+ &
113 & adfac
114 ad_v(iend ,jend+1,k,nout)=ad_v(iend ,jend+1,k,nout)+ &
115 & adfac
116 ad_v(iend+1,jend+1,k,nout)=0.0_r8
117 END DO
118 END IF
119 END IF
120 IF (domain(ng)%NorthWest_Corner(tile)) THEN
121 IF (lbc_apply(ng)%north(istr-1).and. &
122 & lbc_apply(ng)%west (jend+1)) THEN
123 DO k=1,n(ng)
124!^ tl_v(Istr-1,Jend+1,k,nout)=0.5_r8* &
125!^ & (tl_v(Istr-1,Jend ,k,nout)+ &
126!^ & tl_v(Istr ,Jend+1,k,nout))
127!^
128 adfac=0.5_r8*ad_v(istr-1,jend+1,k,nout)
129 ad_v(istr-1,jend ,k,nout)=ad_v(istr-1,jend ,k,nout)+ &
130 & adfac
131 ad_v(istr ,jend+1,k,nout)=ad_v(istr ,jend+1,k,nout)+ &
132 & adfac
133 ad_v(istr-1,jend+1,k,nout)=0.0_r8
134 END DO
135 END IF
136 END IF
137 IF (domain(ng)%SouthEast_Corner(tile)) THEN
138 IF (lbc_apply(ng)%south(iend+1).and. &
139 & lbc_apply(ng)%east (jstr )) THEN
140 DO k=1,n(ng)
141!^ tl_v(Iend+1,Jstr,k,nout)=0.5_r8* &
142!^ & (tl_v(Iend ,Jstr ,k,nout)+ &
143!^ & tl_v(Iend+1,Jstr+1,k,nout))
144!^
145 adfac=0.5_r8*ad_v(iend+1,jstr,k,nout)
146 ad_v(iend ,jstr ,k,nout)=ad_v(iend ,jstr ,k,nout)+ &
147 & adfac
148 ad_v(iend+1,jstr+1,k,nout)=ad_v(iend+1,jstr+1,k,nout)+ &
149 & adfac
150 ad_v(iend+1,jstr ,k,nout)=0.0_r8
151 END DO
152 END IF
153 END IF
154 IF (domain(ng)%SouthWest_Corner(tile)) THEN
155 IF (lbc_apply(ng)%south(istr-1).and. &
156 & lbc_apply(ng)%west (jstr )) THEN
157 DO k=1,n(ng)
158!^ tl_v(Istr-1,Jstr,k,nout)=0.5_r8* &
159!^ & (tl_v(Istr ,Jstr ,k,nout)+ &
160!^ & tl_v(Istr-1,Jstr+1,k,nout))
161!^
162 adfac=0.5_r8*ad_v(istr-1,jstr,k,nout)
163 ad_v(istr ,jstr ,k,nout)=ad_v(istr ,jstr ,k,nout)+ &
164 & adfac
165 ad_v(istr-1,jstr+1,k,nout)=ad_v(istr-1,jstr+1,k,nout)+ &
166 & adfac
167 ad_v(istr-1,jstr ,k,nout)=0.0_r8
168 END DO
169 END IF
170 END IF
171 END IF
172!
173!-----------------------------------------------------------------------
174! Lateral boundary conditions at the eastern edge.
175!-----------------------------------------------------------------------
176!
177 IF (domain(ng)%Eastern_Edge(tile)) THEN
178!
179! Eastern edge, implicit upstream radiation condition.
180!
181 IF (ad_lbc(ieast,isvvel,ng)%radiation) THEN
182 IF (iic(ng).ne.0) THEN
183 DO k=1,n(ng)
184 DO j=jstrv,jend
185 IF (lbc_apply(ng)%east(j)) THEN
186# if defined CELERITY_READ && defined FORWARD_READ
187 IF (ad_lbc(ieast,isvvel,ng)%nudging) THEN
188 IF (lnudgem3clm(ng)) THEN
189 obc_out=0.5_r8* &
190 & (clima(ng)%M3nudgcof(iend+1,j-1,k)+ &
191 & clima(ng)%M3nudgcof(iend+1,j ,k))
192 obc_in =obcfac(ng)*obc_out
193 ELSE
194 obc_out=m3obc_out(ng,ieast)
195 obc_in =m3obc_in(ng,ieast)
196 END IF
197 IF (boundary(ng)%v_east_Cx(j,k).lt.0.0_r8) THEN
198 tau=obc_in
199 ELSE
200 tau=obc_out
201 END IF
202 tau=tau*dt(ng)
203 END IF
204 cx=boundary(ng)%v_east_Cx(j,k)
205# ifdef RADIATION_2D
206 ce=boundary(ng)%v_east_Ce(j,k)
207# else
208 ce=0.0_r8
209# endif
210 cff=boundary(ng)%v_east_C2(j,k)
211# endif
212# ifdef MASKING
213!^ tl_v(Iend+1,j,k,nout)=tl_v(Iend+1,j,k,nout)* &
214!^ & GRID(ng)%vmask(Iend+1,j)
215!^
216 ad_v(iend+1,j,k,nout)=ad_v(iend+1,j,k,nout)* &
217 & grid(ng)%vmask(iend+1,j)
218# endif
219 IF (ad_lbc(ieast,isvvel,ng)%nudging) THEN
220!^ tl_v(Iend+1,j,k,nout)=tl_v(Iend+1,j,k,nout)- &
221!^ & tau*tl_v(Iend+1,j,k,nstp)
222!^
223 ad_v(iend+1,j,k,nstp)=ad_v(iend+1,j,k,nstp)- &
224 & tau*ad_v(iend+1,j,k,nout)
225 END IF
226!^ tl_v(Iend+1,j,k,nout)=(cff*tl_v(Iend+1,j,k,nstp)+ &
227!^ & Cx *tl_v(Iend ,j,k,nout)- &
228!^ & MAX(Ce,0.0_r8)* &
229!^ & tl_grad(Iend+1,j-1)- &
230!^ & MIN(Ce,0.0_r8)* &
231!^ & tl_grad(Iend+1,j ))/ &
232!^ & (cff+Cx)
233!^
234 adfac=ad_v(iend+1,j,k,nout)/(cff+cx)
235 ad_grad(iend+1,j-1)=ad_grad(iend+1,j-1)- &
236 & max(ce,0.0_r8)*adfac
237 ad_grad(iend+1,j )=ad_grad(iend+1,j )- &
238 & min(ce,0.0_r8)*adfac
239 ad_v(iend ,j,k,nout)=ad_v(iend ,j,k,nout)+cx *adfac
240 ad_v(iend+1,j,k,nstp)=ad_v(iend+1,j,k,nstp)+cff*adfac
241 ad_v(iend+1,j,k,nout)=0.0_r8
242 END IF
243 END DO
244 END DO
245 END IF
246!
247! Eastern edge, clamped boundary condition.
248!
249 ELSE IF (ad_lbc(ieast,isvvel,ng)%clamped) THEN
250 DO k=1,n(ng)
251 DO j=jstrv,jend
252 IF (lbc_apply(ng)%east(j)) THEN
253# ifdef MASKING
254!^ tl_v(Iend+1,j,k,nout)=tl_v(Iend+1,j,k,nout)* &
255!^ & GRID(ng)%vmask(Iend+1,j)
256!^
257 ad_v(iend+1,j,k,nout)=ad_v(iend+1,j,k,nout)* &
258 & grid(ng)%vmask(iend+1,j)
259# endif
260# ifdef ADJUST_BOUNDARY
261 IF (lobc(ieast,isvvel,ng)) THEN
262!^ tl_v(Iend+1,j,k,nout)=BOUNDARY(ng)%tl_v_east(j,k)
263!^
264 boundary(ng)%ad_v_east(j,k)= &
265 & boundary(ng)%ad_v_east(j,k)+ &
266 & ad_v(iend+1,j,k,nout)
267 ad_v(iend+1,j,k,nout)=0.0_r8
268 ELSE
269!^ tl_v(Iend+1,j,k,nout)=0.0_r8
270!^
271 ad_v(iend+1,j,k,nout)=0.0_r8
272 END IF
273# else
274!^ tl_v(Iend+1,j,k,nout)=0.0_r8
275!^
276 ad_v(iend+1,j,k,nout)=0.0_r8
277# endif
278 END IF
279 END DO
280 END DO
281!
282! Eastern edge, gradient boundary condition.
283!
284 ELSE IF (ad_lbc(ieast,isvvel,ng)%gradient) THEN
285 DO k=1,n(ng)
286 DO j=jstrv,jend
287 IF (lbc_apply(ng)%east(j)) THEN
288# ifdef MASKING
289!^ tl_v(Iend+1,j,k,nout)=tl_v(Iend+1,j,k,nout)* &
290!^ & GRID(ng)%vmask(Iend+1,j)
291!^
292 ad_v(iend+1,j,k,nout)=ad_v(iend+1,j,k,nout)* &
293 & grid(ng)%vmask(iend+1,j)
294# endif
295!^ tl_v(Iend+1,j,k,nout)=tl_v(Iend,j,k,nout)
296!^
297 ad_v(iend ,j,k,nout)=ad_v(iend ,j,k,nout)+ &
298 & ad_v(iend+1,j,k,nout)
299 ad_v(iend+1,j,k,nout)=0.0_r8
300 END IF
301 END DO
302 END DO
303!
304! Eastern edge, closed boundary condition: free slip (gamma2=1) or
305! no slip (gamma2=-1).
306!
307 ELSE IF (ad_lbc(ieast,isvvel,ng)%closed) THEN
308 IF (nsperiodic(ng)) THEN
309 jmin=jstrv
310 jmax=jend
311 ELSE
312 jmin=jstr
313 jmax=jendr
314 END IF
315 DO k=1,n(ng)
316 DO j=jmin,jmax
317 IF (lbc_apply(ng)%east(j)) THEN
318# ifdef MASKING
319!^ tl_v(Iend+1,j,k,nout)=tl_v(Iend+1,j,k,nout)* &
320!^ & GRID(ng)%vmask(Iend+1,j)
321!^
322 ad_v(iend+1,j,k,nout)=ad_v(iend+1,j,k,nout)* &
323 & grid(ng)%vmask(iend+1,j)
324# endif
325!^ tl_v(Iend+1,j,k,nout)=gamma2(ng)*tl_v(Iend,j,k,nout)
326!^
327 ad_v(iend ,j,k,nout)=ad_v(iend ,j,k,nout)+ &
328 & gamma2(ng)*ad_v(iend+1,j,k,nout)
329 ad_v(iend+1,j,k,nout)=0.0_r8
330 END IF
331 END DO
332 END DO
333 END IF
334 END IF
335!
336!-----------------------------------------------------------------------
337! Lateral boundary conditions at the western edge.
338!-----------------------------------------------------------------------
339!
340 IF (domain(ng)%Western_Edge(tile)) THEN
341!
342! Western edge, implicit upstream radiation condition.
343!
344 IF (ad_lbc(iwest,isvvel,ng)%radiation) THEN
345 IF (iic(ng).ne.0) THEN
346 DO k=1,n(ng)
347 DO j=jstrv,jend
348 IF (lbc_apply(ng)%west(j)) THEN
349# if defined CELERITY_READ && defined FORWARD_READ
350 IF (ad_lbc(iwest,isvvel,ng)%nudging) THEN
351 IF (lnudgem3clm(ng)) THEN
352 obc_out=0.5_r8* &
353 & (clima(ng)%M3nudgcof(istr-1,j-1,k)+ &
354 & clima(ng)%M3nudgcof(istr-1,j ,k))
355 obc_in =obcfac(ng)*obc_out
356 ELSE
357 obc_out=m3obc_out(ng,iwest)
358 obc_in =m3obc_in(ng,iwest)
359 END IF
360 IF (boundary(ng)%v_west_Cx(j,k).lt.0.0_r8) THEN
361 tau=obc_in
362 ELSE
363 tau=obc_out
364 END IF
365 tau=tau*dt(ng)
366 END IF
367 cx=boundary(ng)%v_west_Cx(j,k)
368# ifdef RADIATION_2D
369 ce=boundary(ng)%v_west_Ce(j,k)
370# else
371 ce=0.0_r8
372# endif
373 cff=boundary(ng)%v_west_C2(j,k)
374# endif
375# ifdef MASKING
376!^ tl_v(Istr-1,j,k,nout)=tl_v(Istr-1,j,k,nout)* &
377!^ & GRID(ng)%vmask(Istr-1,j)
378!^
379 ad_v(istr-1,j,k,nout)=ad_v(istr-1,j,k,nout)* &
380 & grid(ng)%vmask(istr-1,j)
381# endif
382 IF (ad_lbc(iwest,isvvel,ng)%nudging) THEN
383!^ tl_v(Istr-1,j,k,nout)=tl_v(Istr-1,j,k,nout)- &
384!^ & tau*tl_v(Istr-1,j,k,nstp)
385!^
386 ad_v(istr-1,j,k,nstp)=ad_v(istr-1,j,k,nstp)- &
387 & tau*ad_v(istr-1,j,k,nout)
388 END IF
389!^ tl_v(Istr-1,j,k,nout)=(cff*tl_v(Istr-1,j,k,nstp)+ &
390!^ & Cx *tl_v(Istr ,j,k,nout)- &
391!^ & MAX(Ce,0.0_r8)* &
392!^ & tl_grad(Istr-1,j-1)- &
393!^ & MIN(Ce,0.0_r8)* &
394!^ & tl_grad(Istr-1,j ))/ &
395!^ & (cff+Cx)
396!^
397 adfac=ad_v(istr-1,j,k,nout)/(cff+cx)
398 ad_grad(istr-1,j-1)=ad_grad(istr-1,j-1)- &
399 & max(ce,0.0_r8)*adfac
400 ad_grad(istr-1,j )=ad_grad(istr-1,j )- &
401 & min(ce,0.0_r8)*adfac
402 ad_v(istr-1,j,k,nstp)=ad_v(istr-1,j,k,nstp)+cff*adfac
403 ad_v(istr ,j,k,nout)=ad_v(istr ,j,k,nout)+cx *adfac
404 ad_v(istr-1,j,k,nout)=0.0_r8
405 END IF
406 END DO
407 END DO
408 END IF
409!
410! Western edge, clamped boundary condition.
411!
412 ELSE IF (ad_lbc(iwest,isvvel,ng)%clamped) THEN
413 DO k=1,n(ng)
414 DO j=jstrv,jend
415 IF (lbc_apply(ng)%west(j)) THEN
416# ifdef MASKING
417!^ tl_v(Istr-1,j,k,nout)=tl_v(Istr-1,j,k,nout)* &
418!^ & GRID(ng)%vmask(Istr-1,j)
419!^
420 ad_v(istr-1,j,k,nout)=ad_v(istr-1,j,k,nout)* &
421 & grid(ng)%vmask(istr-1,j)
422# endif
423# ifdef ADJUST_BOUNDARY
424 IF (lobc(iwest,isvvel,ng)) THEN
425!^ tl_v(Istr-1,j,k,nout)=BOUNDARY(ng)%tl_v_west(j,k)
426!^
427 boundary(ng)%ad_v_west(j,k)= &
428 & boundary(ng)%ad_v_west(j,k)+ &
429 & ad_v(istr-1,j,k,nout)
430 ad_v(istr-1,j,k,nout)=0.0_r8
431 ELSE
432!^ tl_v(Istr-1,j,k,nout)=0.0_r8
433!^
434 ad_v(istr-1,j,k,nout)=0.0_r8
435 END IF
436# else
437!^ tl_v(Istr-1,j,k,nout)=0.0_r8
438!^
439 ad_v(istr-1,j,k,nout)=0.0_r8
440# endif
441 END IF
442 END DO
443 END DO
444!
445! Western edge, gradient boundary condition.
446!
447 ELSE IF (ad_lbc(iwest,isvvel,ng)%gradient) THEN
448 DO k=1,n(ng)
449 DO j=jstrv,jend
450 IF (lbc_apply(ng)%west(j)) THEN
451# ifdef MASKING
452!^ tl_v(Istr-1,j,k,nout)=tl_v(Istr-1,j,k,nout)* &
453!^ & GRID(ng)%vmask(Istr-1,j)
454!^
455 ad_v(istr-1,j,k,nout)=ad_v(istr-1,j,k,nout)* &
456 & grid(ng)%vmask(istr-1,j)
457# endif
458!^ tl_v(Istr-1,j,k,nout)=tl_v(Istr,j,k,nout)
459!^
460 ad_v(istr ,j,k,nout)=ad_v(istr ,j,k,nout)+ &
461 & ad_v(istr-1,j,k,nout)
462 ad_v(istr-1,j,k,nout)=0.0_r8
463 END IF
464 END DO
465 END DO
466!
467! Western edge, closed boundary condition: free slip (gamma2=1) or
468! no slip (gamma2=-1).
469!
470 ELSE IF (ad_lbc(iwest,isvvel,ng)%closed) THEN
471 IF (nsperiodic(ng)) THEN
472 jmin=jstrv
473 jmax=jend
474 ELSE
475 jmin=jstr
476 jmax=jendr
477 END IF
478 DO k=1,n(ng)
479 DO j=jmin,jmax
480 IF (lbc_apply(ng)%west(j)) THEN
481# ifdef MASKING
482!^ tl_v(Istr-1,j,k,nout)=tl_v(Istr-1,j,k,nout)* &
483!^ & GRID(ng)%vmask(Istr-1,j)
484!^
485 ad_v(istr-1,j,k,nout)=ad_v(istr-1,j,k,nout)* &
486 & grid(ng)%vmask(istr-1,j)
487# endif
488!^ tl_v(Istr-1,j,k,nout)=gamma2(ng)*tl_v(Istr,j,k,nout)
489!^
490 ad_v(istr ,j,k,nout)=ad_v(istr ,j,k,nout)+ &
491 & gamma2(ng)*ad_v(istr-1,j,k,nout)
492 ad_v(istr-1,j,k,nout)=0.0
493 END IF
494 END DO
495 END DO
496 END IF
497 END IF
498!
499!-----------------------------------------------------------------------
500! Lateral boundary conditions at the northern edge.
501!-----------------------------------------------------------------------
502!
503 IF (domain(ng)%Northern_Edge(tile)) THEN
504!
505! Northern edge, implicit upstream radiation condition.
506!
507 IF (ad_lbc(inorth,isvvel,ng)%radiation) THEN
508 IF (iic(ng).ne.0) THEN
509 DO k=1,n(ng)
510 DO i=istr,iend
511 IF (lbc_apply(ng)%north(i)) THEN
512# if defined CELERITY_READ && defined FORWARD_READ
513 IF (ad_lbc(inorth,isvvel,ng)%nudging) THEN
514 IF (lnudgem3clm(ng)) THEN
515 obc_out=0.5_r8* &
516 & (clima(ng)%M3nudgcof(i,jend ,k)+ &
517 & clima(ng)%M3nudgcof(i,jend+1,k))
518 obc_in =obcfac(ng)*obc_out
519 ELSE
520 obc_out=m3obc_out(ng,inorth)
521 obc_in =m3obc_in(ng,inorth)
522 END IF
523 IF (boundary(ng)%v_south_Ce(i,k).lt.0.0_r8) THEN
524 tau=obc_in
525 ELSE
526 tau=obc_out
527 END IF
528 tau=tau*dt(ng)
529 END IF
530# ifdef RADIATION_2D
531 cx=boundary(ng)%v_south_Cx(i,k)
532# else
533 cx=0.0_r8
534# endif
535 ce=boundary(ng)%v_south_Ce(i,k)
536 cff=boundary(ng)%v_south_C2(i,k)
537# endif
538# ifdef MASKING
539!^ tl_v(i,Jend+1,k,nout)=tl_v(i,Jend+1,k,nout)* &
540!^ & GRID(ng)%vmask(i,Jend+1)
541!^
542 ad_v(i,jend+1,k,nout)=ad_v(i,jend+1,k,nout)* &
543 & grid(ng)%vmask(i,jend+1)
544# endif
545 IF (ad_lbc(inorth,isvvel,ng)%nudging) THEN
546!^ tl_v(i,Jend+1,k,nout)=tl_v(i,Jend+1,k,nout)- &
547!^ & tau*tl_v(i,Jend+1,k,nstp)
548!^
549 ad_v(i,jend+1,k,nstp)=ad_v(i,jend+1,k,nstp)- &
550 & tau*ad_v(i,jend+1,k,nout)
551 END IF
552!^ tl_v(i,Jend+1,k,nout)=(cff*tl_v(i,Jend+1,k,nstp)+ &
553!^ & Ce *tl_v(i,Jend ,k,nout)- &
554!^ & MAX(Cx,0.0_r8)* &
555!^ & tl_grad(i ,Jend+1)- &
556!^ & MIN(Cx,0.0_r8)* &
557!^ & tl_grad(i+1,Jend+1))/ &
558!^ & (cff+Ce)
559!^
560 adfac=ad_v(i,jend+1,k,nout)/(cff+ce)
561 ad_grad(i ,jend+1)=ad_grad(i ,jend+1)- &
562 & max(cx,0.0_r8)*adfac
563 ad_grad(i+1,jend+1)=ad_grad(i+1,jend+1)- &
564 & min(cx,0.0_r8)*adfac
565 ad_v(i,jend ,k,nstp)=ad_v(i,jend ,k,nstp)+ce *adfac
566 ad_v(i,jend+1,k,nstp)=ad_v(i,jend+1,k,nstp)+cff*adfac
567 ad_v(i,jend+1,k,nout)=0.0_r8
568 END IF
569 END DO
570 END DO
571 END IF
572!
573! Northern edge, clamped boundary condition.
574!
575 ELSE IF (ad_lbc(inorth,isvvel,ng)%clamped) THEN
576 DO k=1,n(ng)
577 DO i=istr,iend
578 IF (lbc_apply(ng)%north(i)) THEN
579# ifdef MASKING
580!^ tl_v(i,Jend+1,k,nout)=tl_v(i,Jend+1,k,nout)* &
581!^ & GRID(ng)%vmask(i,Jend+1)
582!^
583 ad_v(i,jend+1,k,nout)=ad_v(i,jend+1,k,nout)* &
584 & grid(ng)%vmask(i,jend+1)
585# endif
586# ifdef ADJUST_BOUNDARY
587 IF (lobc(inorth,isvvel,ng)) THEN
588!^ tl_v(i,Jend+1,k,nout)=BOUNDARY(ng)%tl_v_north(i,k)
589!^
590 boundary(ng)%ad_v_north(i,k)=boundary(ng)% &
591 & ad_v_north(i,k)+ &
592 & ad_v(i,jend+1,k,nout)
593 ad_v(i,jend+1,k,nout)=0.0_r8
594 ELSE
595!^ tl_v(i,Jend+1,k,nout)=0.0_r8
596!^
597 ad_v(i,jend+1,k,nout)=0.0_r8
598 END IF
599# else
600!^ tl_v(i,Jend+1,k,nout)=0.0_r8
601!^
602 ad_v(i,jend+1,k,nout)=0.0_r8
603# endif
604 END IF
605 END DO
606 END DO
607!
608! Northern edge, gradient boundary condition.
609!
610 ELSE IF (ad_lbc(inorth,isvvel,ng)%gradient) THEN
611 DO k=1,n(ng)
612 DO i=istr,iend
613 IF (lbc_apply(ng)%north(i)) THEN
614# ifdef MASKING
615!^ tl_v(i,Jend+1,k,nout)=tl_v(i,Jend+1,k,nout)* &
616!^ & GRID(ng)%vmask(i,Jend+1)
617!^
618 ad_v(i,jend+1,k,nout)=ad_v(i,jend+1,k,nout)* &
619 & grid(ng)%vmask(i,jend+1)
620# endif
621!^ tl_v(i,Jend+1,k,nout)=tl_v(i,Jend,k,nout)
622!^
623 ad_v(i,jend ,k,nout)=ad_v(i,jend ,k,nout)+ &
624 & ad_v(i,jend+1,k,nout)
625 ad_v(i,jend+1,k,nout)=0.0_r8
626 END IF
627 END DO
628 END DO
629!
630! Northern edge, closed boundary condition.
631!
632 ELSE IF (ad_lbc(inorth,isvvel,ng)%closed) THEN
633 DO k=1,n(ng)
634 DO i=istr,iend
635 IF (lbc_apply(ng)%north(i)) THEN
636!^ tl_v(i,Jend+1,k,nout)=0.0_r8
637!^
638 ad_v(i,jend+1,k,nout)=0.0_r8
639 END IF
640 END DO
641 END DO
642 END IF
643 END IF
644!
645!-----------------------------------------------------------------------
646! Lateral boundary conditions at the southern edge.
647!-----------------------------------------------------------------------
648!
649 IF (domain(ng)%Southern_Edge(tile)) THEN
650!
651! Southern edge, implicit upstream radiation condition.
652!
653 IF (ad_lbc(isouth,isvvel,ng)%radiation) THEN
654 IF (iic(ng).ne.0) THEN
655 DO k=1,n(ng)
656 DO i=istr,iend
657 IF (lbc_apply(ng)%south(i)) THEN
658# if defined CELERITY_READ && defined FORWARD_READ
659 IF (ad_lbc(isouth,isvvel,ng)%nudging) THEN
660 IF (lnudgem3clm(ng)) THEN
661 obc_out=0.5_r8* &
662 & (clima(ng)%M3nudgcof(i,jstr-1,k)+ &
663 & clima(ng)%M3nudgcof(i,jstr ,k))
664 obc_in =obcfac(ng)*obc_out
665 ELSE
666 obc_out=m3obc_out(ng,isouth)
667 obc_in =m3obc_in(ng,isouth)
668 END IF
669 IF (boundary(ng)%v_south_Ce(i,k).lt.0.0_r8) THEN
670 tau=obc_in
671 ELSE
672 tau=obc_out
673 END IF
674 tau=tau*dt(ng)
675 END IF
676# ifdef RADIATION_2D
677 cx=boundary(ng)%v_south_Cx(i,k)
678# else
679 cx=0.0_r8
680# endif
681 ce=boundary(ng)%v_south_Ce(i,k)
682 cff=boundary(ng)%v_south_C2(i,k)
683# endif
684# ifdef MASKING
685!^ tl_v(i,Jstr,k,nout)=tl_v(i,Jstr,k,nout)* &
686!^ & GRID(ng)%vmask(i,Jstr)
687!^
688 ad_v(i,jstr,k,nout)=ad_v(i,jstr,k,nout)* &
689 & grid(ng)%vmask(i,jstr)
690# endif
691 IF (ad_lbc(isouth,isvvel,ng)%nudging) THEN
692!^ tl_v(i,Jstr,k,nout)=tl_v(i,Jstr,k,nout)- &
693!^ & tau*tl_v(i,Jstr,k,nstp)
694!^
695 ad_v(i,jstr,k,nstp)=ad_v(i,jstr,k,nstp)- &
696 & tau*ad_v(i,jstr,k,nout)
697 END IF
698!^ tl_v(i,Jstr,k,nout)=(cff*tl_v(i,Jstr ,k,nstp)+ &
699!^ & Ce *tl_v(i,Jstr+1,k,nout)- &
700!^ & MAX(Cx,0.0_r8)* &
701!^ & tl_grad(i ,Jstr)- &
702!^ & MIN(Cx,0.0_r8)* &
703!^ & tl_grad(i+1,Jstr))/ &
704!^ & (cff+Ce)
705!^
706 adfac=ad_v(i,jstr,k,nout)/(cff+ce)
707 ad_grad(i ,jstr)=ad_grad(i ,jstr)- &
708 & max(cx,0.0_r8)*adfac
709 ad_grad(i+1,jstr)=ad_grad(i+1,jstr)- &
710 & min(cx,0.0_r8)*adfac
711 ad_v(i,jstr ,k,nstp)=ad_v(i,jstr ,k,nstp)+cff*adfac
712 ad_v(i,jstr+1,k,nout)=ad_v(i,jstr+1,k,nout)+ce *adfac
713 ad_v(i,jstr ,k,nout)=0.0_r8
714 END IF
715 END DO
716 END DO
717 END IF
718!
719! Southern edge, clamped boundary condition.
720!
721 ELSE IF (ad_lbc(isouth,isvvel,ng)%clamped) THEN
722 DO k=1,n(ng)
723 DO i=istr,iend
724 IF (lbc_apply(ng)%south(i)) THEN
725# ifdef MASKING
726!^ tl_v(i,Jstr,k,nout)=tl_v(i,Jstr,k,nout)* &
727!^ & GRID(ng)%vmask(i,Jstr)
728!^
729 ad_v(i,jstr,k,nout)=ad_v(i,jstr,k,nout)* &
730 & grid(ng)%vmask(i,jstr)
731# endif
732# ifdef ADJUST_BOUNDARY
733 IF (lobc(isouth,isvvel,ng)) THEN
734!^ tl_v(i,Jstr,k,nout)=BOUNDARY(ng)%tl_v_south(i,k)
735!^
736 boundary(ng)%ad_v_south(i,k)=boundary(ng)% &
737 & ad_v_south(i,k)+ &
738 & ad_v(i,jstr,k,nout)
739 ad_v(i,jstr,k,nout)=0.0_r8
740 ELSE
741!^ tl_v(i,Jstr,k,nout)=0.0_r8
742!^
743 ad_v(i,jstr,k,nout)=0.0_r8
744 END IF
745# else
746!^ tl_v(i,Jstr,k,nout)=0.0_r8
747!^
748 ad_v(i,jstr,k,nout)=0.0_r8
749# endif
750 END IF
751 END DO
752 END DO
753!
754! Southern edge, gradient boundary condition.
755!
756 ELSE IF (ad_lbc(isouth,isvvel,ng)%gradient) THEN
757 DO k=1,n(ng)
758 DO i=istr,iend
759 IF (lbc_apply(ng)%south(i)) THEN
760# ifdef MASKING
761!^ tl_v(i,Jstr,k,nout)=tl_v(i,Jstr,k,nout)* &
762!^ & GRID(ng)%vmask(i,Jstr)
763!^
764 ad_v(i,jstr ,k,nout)=ad_v(i,jstr ,k,nout)* &
765 & grid(ng)%vmask(i,jstr)
766# endif
767!^ tl_v(i,Jstr,k,nout)=tl_v(i,Jstr+1,k,nout)
768!^
769 ad_v(i,jstr+1,k,nout)=ad_v(i,jstr+1,k,nout)+ &
770 & ad_v(i,jstr ,k,nout)
771 ad_v(i,jstr ,k,nout)=0.0_r8
772 END IF
773 END DO
774 END DO
775!
776! Southern edge, closed boundary condition.
777!
778 ELSE IF (ad_lbc(isouth,isvvel,ng)%closed) THEN
779 DO k=1,n(ng)
780 DO i=istr,iend
781 IF (lbc_apply(ng)%south(i)) THEN
782!^ tl_v(i,Jstr,k,nout)=0.0_r8
783!^
784 ad_v(i,jstr,k,nout)=0.0_r8
785 END IF
786 END DO
787 END DO
788 END IF
789 END IF
790
791 RETURN
792 END SUBROUTINE ad_v3dbc_tile
793#endif
794 END MODULE ad_v3dbc_mod
subroutine, public ad_v3dbc_tile(ng, tile, lbi, ubi, lbj, ubj, ubk, imins, imaxs, jmins, jmaxs, nstp, nout, ad_v)
Definition ad_v3dbc_im.F:57
subroutine, public ad_v3dbc(ng, tile, nout)
Definition ad_v3dbc_im.F:28
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_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
integer isvvel
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
type(t_lbc), dimension(:,:,:), allocatable ad_lbc
Definition mod_param.F:378
integer, dimension(:), allocatable n
Definition mod_param.F:479
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
real(dp), dimension(:,:), allocatable m3obc_out
integer, dimension(:), allocatable iic
real(dp), dimension(:), allocatable dt
logical, dimension(:,:,:), allocatable lobc
logical, dimension(:), allocatable ewperiodic
integer, parameter iwest
logical, dimension(:), allocatable nsperiodic
real(dp), dimension(:), allocatable obcfac
real(r8), dimension(:), allocatable gamma2
logical, dimension(:), allocatable lnudgem3clm
integer, parameter isouth
integer, parameter ieast
integer, parameter inorth
real(dp), dimension(:,:), allocatable m3obc_in
integer, dimension(:), allocatable nstp