ROMS
Loading...
Searching...
No Matches
tl_obc_volcons.F
Go to the documentation of this file.
1#include "cppdefs.h"
3#ifdef TANGENT
4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group Andrew M. Moore !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This routine computes tangent linear integral mass flux "obc_flux" !
13! across the open boundaries, which is needed to enforce global mass !
14! conservation constraint. !
15! !
16!=======================================================================
17!
18 implicit none
19!
20 PRIVATE
21 PUBLIC :: tl_obc_flux_tile
22 PUBLIC :: tl_set_duv_bc_tile
23!
24 CONTAINS
25!
26!***********************************************************************
27 SUBROUTINE tl_obc_flux (ng, tile, kinp)
28!***********************************************************************
29!
30 USE mod_param
31 USE mod_grid
32 USE mod_ocean
33!
34! Imported variable declarations.
35!
36 integer, intent(in) :: ng, tile, kinp
37!
38! Local variable declarations.
39!
40# include "tile.h"
41!
42 CALL tl_obc_flux_tile (ng, tile, &
43 & lbi, ubi, lbj, ubj, &
44 & imins, imaxs, jmins, jmaxs, &
45 & kinp, &
46# ifdef MASKING
47 & grid(ng) % umask, &
48 & grid(ng) % vmask, &
49# endif
50 & grid(ng) % h, &
51 & grid(ng) % tl_h, &
52 & grid(ng) % om_v, &
53 & grid(ng) % on_u, &
54 & ocean(ng) % ubar, &
55 & ocean(ng) % vbar, &
56 & ocean(ng) % zeta, &
57 & ocean(ng) % tl_ubar, &
58 & ocean(ng) % tl_vbar, &
59 & ocean(ng) % tl_zeta)
60
61 RETURN
62 END SUBROUTINE tl_obc_flux
63!
64!***********************************************************************
65 SUBROUTINE tl_obc_flux_tile (ng, tile, &
66 & LBi, UBi, LBj, UBj, &
67 & IminS, ImaxS, JminS, JmaxS, &
68 & kinp, &
69# ifdef MASKING
70 & umask, vmask, &
71# endif
72 & h, tl_h, &
73 & om_v, on_u, &
74 & ubar, vbar, zeta, &
75 & tl_ubar, tl_vbar, tl_zeta)
76!***********************************************************************
77!
78 USE mod_param
79 USE mod_parallel
80 USE mod_scalars
81
82# ifdef DISTRIBUTE
83!
84 USE distribute_mod, ONLY : mp_reduce
85# endif
86!
87! Imported variable declarations.
88!
89 integer, intent(in) :: ng, tile
90 integer, intent(in) :: lbi, ubi, lbj, ubj
91 integer, intent(in) :: imins, imaxs, jmins, jmaxs
92 integer, intent(in) :: kinp
93!
94# ifdef ASSUMED_SHAPE
95# ifdef MASKING
96 real(r8), intent(in) :: umask(lbi:,lbj:)
97 real(r8), intent(in) :: vmask(lbi:,lbj:)
98# endif
99 real(r8), intent(in) :: h(lbi:,lbj:)
100 real(r8), intent(in) :: tl_h(lbi:,lbj:)
101 real(r8), intent(in) :: om_v(lbi:,lbj:)
102 real(r8), intent(in) :: on_u(lbi:,lbj:)
103
104 real(r8), intent(in) :: ubar(lbi:,lbj:,:)
105 real(r8), intent(in) :: vbar(lbi:,lbj:,:)
106 real(r8), intent(in) :: zeta(lbi:,lbj:,:)
107 real(r8), intent(in) :: tl_ubar(lbi:,lbj:,:)
108 real(r8), intent(in) :: tl_vbar(lbi:,lbj:,:)
109 real(r8), intent(in) :: tl_zeta(lbi:,lbj:,:)
110# else
111# ifdef MASKING
112 real(r8), intent(in) :: umask(lbi:ubi,lbj:ubj)
113 real(r8), intent(in) :: vmask(lbi:ubi,lbj:ubj)
114# endif
115 real(r8), intent(in) :: h(lbi:ubi,lbj:ubj)
116 real(r8), intent(in) :: tl_h(lbi:ubi,lbj:ubj)
117 real(r8), intent(in) :: om_v(lbi:ubi,lbj:ubj)
118 real(r8), intent(in) :: on_u(lbi:ubi,lbj:ubj)
119
120 real(r8), intent(in) :: ubar(lbi:ubi,lbj:ubj,:)
121 real(r8), intent(in) :: vbar(lbi:ubi,lbj:ubj,:)
122 real(r8), intent(in) :: zeta(lbi:ubi,lbj:ubj,:)
123 real(r8), intent(in) :: tl_ubar(lbi:ubi,lbj:ubj,:)
124 real(r8), intent(in) :: tl_vbar(lbi:ubi,lbj:ubj,:)
125 real(r8), intent(in) :: tl_zeta(lbi:ubi,lbj:ubj,:)
126# endif
127!
128! Local variable declarations.
129!
130 integer :: nsub, i, j
131
132 real(r8) :: cff, my_area, my_flux
133 real(r8) :: tl_cff, tl_my_area, tl_my_flux
134
135# ifdef DISTRIBUTE
136 real(r8), dimension(2) :: rbuffer
137 character (len=3), dimension(2) :: op_handle
138# endif
139
140# include "set_bounds.h"
141!
142!-----------------------------------------------------------------------
143! Compute open segments cross-section area and mass flux.
144!-----------------------------------------------------------------------
145!
146!^ my_area=0.0_r8
147!^ my_flux=0.0_r8
148!^
149 tl_my_area=0.0_r8
150 tl_my_flux=0.0_r8
151!
152 IF (tl_volcons(iwest,ng)) THEN
153 IF (domain(ng)%Western_Edge(tile)) THEN
154 DO j=jstr,jend
155 cff=0.5_r8*on_u(istr,j)* &
156 & (zeta(istr-1,j,kinp)+h(istr-1,j)+ &
157 & zeta(istr ,j,kinp)+h(istr ,j))
158 tl_cff=0.5_r8*on_u(istr,j)* &
159 & (tl_zeta(istr-1,j,kinp)+tl_h(istr-1,j)+ &
160 & tl_zeta(istr ,j,kinp)+tl_h(istr-1,j))
161# ifdef MASKING
162 cff=cff*umask(istr,j)
163 tl_cff=tl_cff*umask(istr,j)
164# endif
165!^ my_area=my_area+cff
166!^
167 tl_my_area=tl_my_area+tl_cff
168!^ my_flux=my_flux+cff*ubar(Istr,j,kinp)
169!^
170 tl_my_flux=tl_my_flux+ &
171 & tl_cff*ubar(istr,j,kinp)+ &
172 & cff*tl_ubar(istr,j,kinp)
173 END DO
174 END IF
175 END IF
176
177 IF (tl_volcons(ieast,ng)) THEN
178 IF (domain(ng)%Eastern_Edge(tile)) THEN
179 DO j=jstr,jend
180 cff=0.5_r8*on_u(iend+1,j)* &
181 & (zeta(iend ,j,kinp)+h(iend ,j)+ &
182 & zeta(iend+1,j,kinp)+h(iend+1,j))
183 tl_cff=0.5_r8*on_u(iend+1,j)* &
184 & (tl_zeta(iend ,j,kinp)+tl_h(iend ,j)+ &
185 & tl_zeta(iend+1,j,kinp)+tl_h(iend+1,j))
186# ifdef MASKING
187 cff=cff*umask(iend+1,j)
188 tl_cff=tl_cff*umask(iend+1,j)
189# endif
190!^ my_area=my_area+cff
191!^
192 tl_my_area=tl_my_area+tl_cff
193!^ my_flux=my_flux-cff*ubar(Iend+1,j,kinp)
194!^
195 tl_my_flux=tl_my_flux- &
196 & tl_cff*ubar(iend+1,j,kinp)- &
197 & cff*tl_ubar(iend+1,j,kinp)
198 END DO
199 END IF
200 END IF
201
202 IF (tl_volcons(isouth,ng)) THEN
203 IF (domain(ng)%Southern_Edge(tile)) THEN
204 DO i=istr,iend
205 cff=0.5_r8*om_v(i,jstr)* &
206 & (zeta(i,jstr-1,kinp)+h(i,jstr-1)+ &
207 & zeta(i,jstr ,kinp)+h(i,jstr ))
208 tl_cff=0.5_r8*om_v(i,jstr)* &
209 & (tl_zeta(i,jstr-1,kinp)+tl_h(i,jstr-1)+ &
210 & tl_zeta(i,jstr ,kinp)+tl_h(i,jstr ))
211# ifdef MASKING
212 cff=cff*vmask(i,jstr)
213 tl_cff=tl_cff*vmask(i,jstr)
214# endif
215!^ my_area=my_area+cff
216!^
217 tl_my_area=tl_my_area+tl_cff
218!^ my_flux=my_flux+cff*vbar(i,JstrV-1,kinp)
219!^
220 tl_my_flux=tl_my_flux+ &
221 & tl_cff*vbar(i,jstrv-1,kinp)+ &
222 & cff*tl_vbar(i,jstrv-1,kinp)
223 END DO
224 END IF
225 END IF
226
227 IF (tl_volcons(inorth,ng)) THEN
228 IF (domain(ng)%Northern_Edge(tile)) THEN
229 DO i=istr,iend
230 cff=0.5_r8*om_v(i,jend+1)* &
231 & (zeta(i,jend ,kinp)+h(i,jend )+ &
232 & zeta(i,jend+1,kinp)+h(i,jend+1))
233 tl_cff=0.5_r8*om_v(i,jend+1)* &
234 & (tl_zeta(i,jend ,kinp)+tl_h(i,jend )+ &
235 & tl_zeta(i,jend+1,kinp)+tl_h(i,jend+1))
236# ifdef MASKING
237 cff=cff*vmask(i,jend+1)
238 tl_cff=tl_cff*vmask(i,jend+1)
239# endif
240!^ my_area=my_area+cff
241!^
242 tl_my_area=tl_my_area+tl_cff
243!^ my_flux=my_flux-cff*vbar(i,Jend+1,kinp)
244!^
245 tl_my_flux=tl_my_flux- &
246 & tl_cff*vbar(i,jend+1,kinp)- &
247 & cff*tl_vbar(i,jend+1,kinp)
248 END DO
249 END IF
250 END IF
251!
252!-----------------------------------------------------------------------
253! Perform global summation and compute correction velocity.
254!-----------------------------------------------------------------------
255!
256 IF (any(tl_volcons(:,ng))) THEN
257# ifdef DISTRIBUTE
258 nsub=1 ! distributed-memory
259# else
260 IF (domain(ng)%SouthWest_Corner(tile).and. &
261 & domain(ng)%NorthEast_Corner(tile)) THEN
262 nsub=1 ! non-tiled application
263 ELSE
264 nsub=ntilex(ng)*ntilee(ng) ! tiled application
265 END IF
266# endif
267!$OMP CRITICAL (TL_OBC_VOLUME)
268 IF (tile_count.eq.0) THEN
269!^ bc_flux=0.0_r8
270!^ bc_area=0.0_r8
271!^
272 tl_bc_flux=0.0_r8
273 tl_bc_area=0.0_r8
274 END IF
275!^ bc_area=bc_area+my_area
276!^ bc_flux=bc_flux+my_flux
277!^
278 tl_bc_area=tl_bc_area+tl_my_area
279 tl_bc_flux=tl_bc_flux+tl_my_flux
281 IF (tile_count.eq.nsub) THEN
282 tile_count=0
283# ifdef DISTRIBUTE
284!^ rbuffer(1)=bc_area
285!^ rbuffer(2)=bc_flux
286!^
287 rbuffer(1)=tl_bc_area
288 rbuffer(2)=tl_bc_flux
289 op_handle(1)='SUM'
290 op_handle(2)='SUM'
291 CALL mp_reduce (ng, itlm, 2, rbuffer, op_handle)
292!^ bc_area=rbuffer(1)
293!^ bc_flux=rbuffer(2)
294!^
295 tl_bc_area=rbuffer(1)
296 tl_bc_flux=rbuffer(2)
297# endif
298!^ ubar_xs=bc_flux/bc_area
299!^
302 END IF
303!$OMP END CRITICAL (TL_OBC_VOLUME)
304 END IF
305
306 RETURN
307 END SUBROUTINE tl_obc_flux_tile
308!
309!***********************************************************************
310 SUBROUTINE tl_set_duv_bc_tile (ng, tile, &
311 & LBi, UBi, LBj, UBj, &
312 & IminS, ImaxS, JminS, JmaxS, &
313 & kinp, &
314# ifdef MASKING
315 & umask, vmask, &
316# endif
317 & om_v, on_u, &
318 & ubar, vbar, &
319 & tl_ubar, tl_vbar, &
320 & Drhs, Duon, Dvom, &
321 & tl_Drhs, tl_Duon, tl_Dvom)
322!***********************************************************************
323!
324 USE mod_param
325 USE mod_scalars
326# ifdef DISTRIBUTE
327!
329# endif
330!
331! Imported variable declarations.
332!
333 integer, intent(in) :: ng, tile
334 integer, intent(in) :: lbi, ubi, lbj, ubj
335 integer, intent(in) :: imins, imaxs, jmins, jmaxs
336 integer, intent(in) :: kinp
337!
338# ifdef ASSUMED_SHAPE
339# ifdef MASKING
340 real(r8), intent(in) :: umask(lbi:,lbj:)
341 real(r8), intent(in) :: vmask(lbi:,lbj:)
342# endif
343 real(r8), intent(in) :: om_v(lbi:,lbj:)
344 real(r8), intent(in) :: on_u(lbi:,lbj:)
345 real(r8), intent(in) :: ubar(lbi:,lbj:,:)
346 real(r8), intent(in) :: vbar(lbi:,lbj:,:)
347 real(r8), intent(in) :: drhs(imins:,jmins:)
348 real(r8), intent(in) :: duon(imins:,jmins:)
349 real(r8), intent(in) :: dvom(imins:,jmins:)
350 real(r8), intent(in) :: tl_ubar(lbi:,lbj:,:)
351 real(r8), intent(in) :: tl_vbar(lbi:,lbj:,:)
352 real(r8), intent(in) :: tl_drhs(imins:,jmins:)
353
354 real(r8), intent(inout) :: tl_duon(imins:,jmins:)
355 real(r8), intent(inout) :: tl_dvom(imins:,jmins:)
356# else
357# ifdef MASKING
358 real(r8), intent(in) :: umask(lbi:ubi,lbj:ubj)
359 real(r8), intent(in) :: vmask(lbi:ubi,lbj:ubj)
360# endif
361 real(r8), intent(in) :: om_v(lbi:ubi,lbj:ubj)
362 real(r8), intent(in) :: on_u(lbi:ubi,lbj:ubj)
363 real(r8), intent(in) :: ubar(lbi:ubi,lbj:ubj,:)
364 real(r8), intent(in) :: vbar(lbi:ubi,lbj:ubj,:)
365 real(r8), intent(in) :: drhs(imins:imaxs,jmins:jmaxs)
366 real(r8), intent(in) :: duon(imins:imaxs,jmins:jmaxs)
367 real(r8), intent(in) :: dvom(imins:imaxs,jmins:jmaxs)
368 real(r8), intent(in) :: tl_ubar(lbi:ubi,lbj:ubj,:)
369 real(r8), intent(in) :: tl_vbar(lbi:ubi,lbj:ubj,:)
370 real(r8), intent(in) :: tl_drhs(imins:,jmins:)
371
372 real(r8), intent(inout) :: tl_duon(imins:imaxs,jmins:jmaxs)
373 real(r8), intent(inout) :: tl_dvom(imins:imaxs,jmins:jmaxs)
374# endif
375!
376! Local variable declarations.
377!
378 integer :: i, j
379
380# include "set_bounds.h"
381!
382!-----------------------------------------------------------------------
383! Set vertically integrated mass fluxes "Duon" and "Dvom" along
384! the open boundaries in such a way that the integral volume is
385! conserved. This is done by applying "ubar_xs" correction to
386! the velocities.
387!-----------------------------------------------------------------------
388!
389# ifdef DISTRIBUTE
390# define I_RANGE IstrU,MIN(Iend+1,Lm(ng))
391# define J_RANGE JstrV,MIN(Jend+1,Mm(ng))
392# else
393# define I_RANGE MAX(2,IstrU-1),MIN(Iend+1,Lm(ng))
394# define J_RANGE MAX(2,JstrV-1),MIN(Jend+1,Mm(ng))
395# endif
396
397 IF (tl_volcons(iwest,ng)) THEN
398 IF (domain(ng)%Western_Edge(tile)) THEN
399 DO j=-2+j_range+1
400!^ Duon(Istr,j)=0.5_r8*(Drhs(Istr,j)+Drhs(Istr-1,j))* &
401!^ & (ubar(Istr,j,kinp)-ubar_xs)* &
402!^ & on_u(Istr,j)
403!^
404 tl_duon(istr,j)=0.5_r8* &
405 & ((tl_drhs(istr,j)+tl_drhs(istr-1,j))* &
406 & (ubar(istr,j,kinp)-ubar_xs)+ &
407 & (drhs(istr,j)+drhs(istr-1,j))* &
408 & (tl_ubar(istr,j,kinp)-tl_ubar_xs))* &
409 & on_u(istr,j)
410# ifdef MASKING
411!^ Duon(Istr,j)=Duon(Istr,j)*umask(Istr,j)
412!^
413 tl_duon(istr,j)=tl_duon(istr,j)*umask(istr,j)
414# endif
415 END DO
416 END IF
417 END IF
418
419 IF (tl_volcons(ieast,ng)) THEN
420 IF (domain(ng)%Eastern_Edge(tile)) THEN
421 DO j=-2+j_range+1
422!^ Duon(Iend+1,j)=0.5_r8*(Drhs(Iend+1,j)+Drhs(Iend,j))* &
423!^ & (ubar(Iend+1,j,kinp)+ubar_xs)* &
424!^ & on_u(Iend+1,j)
425!^
426 tl_duon(iend+1,j)=0.5_r8* &
427 & ((tl_drhs(iend+1,j)+tl_drhs(iend,j))* &
428 & (ubar(iend+1,j,kinp)+ubar_xs)+ &
429 & (drhs(iend+1,j)+drhs(iend,j))* &
430 & (tl_ubar(iend+1,j,kinp)+tl_ubar_xs))* &
431 & on_u(iend+1,j)
432# ifdef MASKING
433!^ Duon(Iend+1,j)=Duon(Iend+1,j)*umask(Iend+1,j)
434!^
435 tl_duon(iend+1,j)=tl_duon(iend+1,j)*umask(iend+1,j)
436# endif
437 END DO
438 END IF
439 END IF
440
441 IF (tl_volcons(isouth,ng)) THEN
442 IF (domain(ng)%Southern_Edge(tile)) THEN
443 DO i=-2+i_range+1
444!^ Dvom(i,Jstr)=0.5_r8*(Drhs(i,Jstr)+Drhs(i,Jstr-1))* &
445!^ & (vbar(i,Jstr,kinp)-ubar_xs)* &
446!^ & om_v(i,Jstr)
447!^
448 tl_dvom(i,jstr)=0.5_r8* &
449 & ((tl_drhs(i,jstr)+tl_drhs(i,jstr-1))* &
450 & (vbar(i,jstr,kinp)-ubar_xs)+ &
451 & (drhs(i,jstr)+drhs(i,jstr-1))* &
452 & (tl_vbar(i,jstr,kinp)-tl_ubar_xs))* &
453 & om_v(i,jstr)
454# ifdef MASKING
455!^ Dvom(i,Jstr)=Dvom(i,Jstr)*vmask(i,Jstr)
456!^
457 tl_dvom(i,jstr)=tl_dvom(i,jstr)*vmask(i,jstr)
458# endif
459 END DO
460 END IF
461 END IF
462
463 IF (tl_volcons(inorth,ng)) THEN
464 IF (domain(ng)%Northern_Edge(tile)) THEN
465 DO i=-2+i_range+1
466!^ Dvom(i,Jend+1)=0.5_r8*(Drhs(i,Jend+1)+Drhs(i,Jend))* &
467!^ & (vbar(i,Jend+1,kinp)+ubar_xs)* &
468!^ & om_v(i,Jend+1)
469!^
470 tl_dvom(i,jend+1)=0.5_r8* &
471 & ((tl_drhs(i,jend+1)+tl_drhs(i,jend))* &
472 & (vbar(i,jend+1,kinp)+ubar_xs)+ &
473 & (drhs(i,jend+1)+drhs(i,jend))* &
474 & (tl_vbar(i,jend+1,kinp)+tl_ubar_xs))* &
475 & om_v(i,jend+1)
476# ifdef MASKING
477!^ Dvom(i,Jend+1)=Dvom(i,Jend+1)*vmask(i,Jend+1)
478!^
479 tl_dvom(i,jend+1)=tl_dvom(i,jend+1)*vmask(i,jend+1)
480# endif
481 END DO
482 END IF
483 END IF
484
485# ifdef DISTRIBUTE
486!
487! Do a special exchange to avoid having three ghost points for high
488! order numerical stencil.
489!
490 IF (tl_volcons(iwest,ng).or.tl_volcons(ieast,ng)) THEN
491!^ CALL mp_exchange2d (ng, tile, iNLM, 1, &
492!^ & IminS, ImaxS, JminS, JmaxS, &
493!^ & NghostPoints, &
494!^ & EWperiodic(ng), NSperiodic(ng), &
495!^ & Duon)
496!^
497 CALL mp_exchange2d (ng, tile, itlm, 1, &
498 & imins, imaxs, jmins, jmaxs, &
499 & nghostpoints, &
500 & ewperiodic(ng), nsperiodic(ng), &
501 & tl_duon)
502 END IF
503
504 IF (tl_volcons(isouth,ng).or.tl_volcons(inorth,ng)) THEN
505!^ CALL mp_exchange2d (ng, tile, iNLM, 1, &
506!^ & IminS, ImaxS, JminS, JmaxS, &
507!^ & NghostPoints, &
508!^ & EWperiodic(ng), NSperiodic(ng), &
509!^ & Dvom)
510!^
511 CALL mp_exchange2d (ng, tile, itlm, 1, &
512 & imins, imaxs, jmins, jmaxs, &
513 & nghostpoints, &
514 & ewperiodic(ng), nsperiodic(ng), &
515 & tl_dvom)
516 END IF
517# endif
518
519# undef I_RANGE
520# undef J_RANGE
521
522 RETURN
523 END SUBROUTINE tl_set_duv_bc_tile
524!
525!***********************************************************************
526 SUBROUTINE tl_conserve_mass_tile (ng, tile, &
527 & LBi, UBi, LBj, UBj, &
528 & IminS, ImaxS, JminS, JmaxS, &
529 & kinp, &
530# ifdef MASKING
531 & umask, vmask, &
532# endif
533 & tl_ubar, tl_vbar)
534!***********************************************************************
535!
536 USE mod_param
537 USE mod_scalars
538!
539! Imported variable declarations.
540!
541 integer, intent(in) :: ng, tile
542 integer, intent(in) :: LBi, UBi, LBj, UBj
543 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
544 integer, intent(in) :: kinp
545!
546# ifdef ASSUMED_SHAPE
547# ifdef MASKING
548 real(r8), intent(in) :: umask(LBi:,LBj:)
549 real(r8), intent(in) :: vmask(LBi:,LBj:)
550# endif
551 real(r8), intent(inout) :: tl_ubar(LBi:,LBj:,:)
552 real(r8), intent(inout) :: tl_vbar(LBi:,LBj:,:)
553# else
554# ifdef MASKING
555 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
556 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
557# endif
558 real(r8), intent(inout) :: tl_ubar(LBi:UBi,LBj:UBj,:)
559 real(r8), intent(inout) :: tl_vbar(LBi:UBi,LBj:UBj,:)
560# endif
561!
562! Local variable declarations.
563!
564 integer :: i, j
565
566# include "set_bounds.h"
567!
568!-----------------------------------------------------------------------
569! Corrects velocities across the open boundaries to enforce global
570! mass conservation constraint.
571!-----------------------------------------------------------------------
572!
573 IF (tl_volcons(iwest,ng)) THEN
574 IF (domain(ng)%Western_Edge(tile)) THEN
575 DO j=jstr,jend
576!^ ubar(Istr,j,kinp)=ubar(Istr,j,kinp)-ubar_xs
577!^
578 tl_ubar(istr,j,kinp)=tl_ubar(istr,j,kinp)-tl_ubar_xs
579# ifdef MASKING
580!^ ubar(Istr,j,kinp)=ubar(Istr,j,kinp)*umask(Istr,j)
581!^
582 tl_ubar(istr,j,kinp)=tl_ubar(istr,j,kinp)* &
583 & umask(istr,j)
584# endif
585 END DO
586 END IF
587 END IF
588
589 IF (tl_volcons(ieast,ng)) THEN
590 IF (domain(ng)%Eastern_Edge(tile)) THEN
591 DO j=jstr,jend
592!^ ubar(Iend+1,j,kinp)=ubar(Iend+1,j,kinp)+ubar_xs
593!^
594 tl_ubar(iend+1,j,kinp)=tl_ubar(iend+1,j,kinp)+tl_ubar_xs
595# ifdef MASKING
596!^ ubar(Iend+1,j,kinp)=ubar(Iend+1,j,kinp)*umask(Iend+1,j)
597!^
598 tl_ubar(iend+1,j,kinp)=tl_ubar(iend+1,j,kinp)* &
599 & umask(iend+1,j)
600# endif
601 END DO
602 END IF
603 END IF
604
605 IF (tl_volcons(isouth,ng)) THEN
606 IF (domain(ng)%Southern_Edge(tile)) THEN
607 DO i=istr,iend
608!^ vbar(i,Jstr,kinp)=(vbar(i,Jstr,kinp)-ubar_xs)
609!^
610 tl_vbar(i,jstr,kinp)=(tl_vbar(i,jstr,kinp)-tl_ubar_xs)
611# ifdef MASKING
612!^ vbar(i,Jstr,kinp)=vbar(i,Jstr,kinp)*vmask(i,Jstr)
613!^
614 tl_vbar(i,jstr,kinp)=tl_vbar(i,jstr,kinp)* &
615 & vmask(i,jstr)
616# endif
617 END DO
618 END IF
619 END IF
620
621 IF (tl_volcons(inorth,ng)) THEN
622 IF (domain(ng)%Northern_Edge(tile)) THEN
623 DO i=istr,iend
624!^ vbar(i,Jend+1,kinp)=(vbar(i,Jend+1,kinp)+ubar_xs)
625!^
626 tl_vbar(i,jend+1,kinp)=(tl_vbar(i,jend+1,kinp)+tl_ubar_xs)
627# ifdef MASKING
628!^ vbar(i,Jend+1,kinp)=vbar(i,Jend+1,kinp)*vmask(i,Jend+1)
629!^
630 tl_vbar(i,jend+1,kinp)=tl_vbar(i,jend+1,kinp)* &
631 & vmask(i,jend+1)
632# endif
633 END DO
634 END IF
635 END IF
636
637 RETURN
638 END SUBROUTINE tl_conserve_mass_tile
639#endif
640 END MODULE tl_obc_volcons_mod
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
integer tile_count
integer, dimension(:), allocatable ntilex
Definition mod_param.F:685
integer nghostpoints
Definition mod_param.F:710
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer, dimension(:), allocatable ntilee
Definition mod_param.F:686
integer, parameter itlm
Definition mod_param.F:663
real(dp) ubar_xs
real(dp) tl_bc_area
logical, dimension(:), allocatable ewperiodic
integer, parameter iwest
logical, dimension(:), allocatable nsperiodic
real(dp) bc_area
real(dp) tl_ubar_xs
integer, parameter isouth
logical, dimension(:,:), allocatable tl_volcons
real(dp) tl_bc_flux
integer, parameter ieast
integer, parameter inorth
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine, public tl_obc_flux_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, kinp, umask, vmask, h, tl_h, om_v, on_u, ubar, vbar, zeta, tl_ubar, tl_vbar, tl_zeta)
subroutine, public tl_set_duv_bc_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, kinp, umask, vmask, om_v, on_u, ubar, vbar, tl_ubar, tl_vbar, drhs, duon, dvom, tl_drhs, tl_duon, tl_dvom)
subroutine tl_conserve_mass_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, kinp, umask, vmask, tl_ubar, tl_vbar)
subroutine tl_obc_flux(ng, tile, kinp)