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