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