Context Navigation

Back to Ticket #805

Ticket #805: v2dbc_im.F

File v2dbc_im.F, 51.6 KB (added by m.hadfield, 7 years ago)

Line
1	#include "cppdefs.h"
2
3	MODULE v2dbc_mod
4	!
5	!svn $Id$
6	!=======================================================================
7	! Copyright (c) 2002-2019 The ROMS/TOMS Group !
8	! Licensed under a MIT/X style license !
9	! See License_ROMS.txt Hernan G. Arango !
10	!========================================== Alexander F. Shchepetkin ===
11	! !
12	! This subroutine sets lateral boundary conditions for vertically !
13	! integrated V-velocity. !
14	! !
15	!=======================================================================
16	!
17	implicit none
18	!
19	PRIVATE
20	PUBLIC :: v2dbc, v2dbc_tile
21	!
22	CONTAINS
23	!
24	!***********************************************************************
25	SUBROUTINE v2dbc (ng, tile, kout)
26	!***********************************************************************
27	!
28	USE mod_param
29	USE mod_ocean
30	USE mod_stepping
31	!
32	! Imported variable declarations.
33	!
34	integer, intent(in) :: ng, tile, kout
35	!
36	! Local variable declarations.
37	!
38	#include "tile.h"
39	!
40	CALL v2dbc_tile (ng, tile, &
41	& LBi, UBi, LBj, UBj, &
42	& IminS, ImaxS, JminS, JmaxS, &
43	& krhs(ng), kstp(ng), kout, &
44	& OCEAN(ng) % ubar, &
45	& OCEAN(ng) % vbar, &
46	& OCEAN(ng) % zeta)
47
48	RETURN
49	END SUBROUTINE v2dbc
50	!
51	!***********************************************************************
52	SUBROUTINE v2dbc_tile (ng, tile, &
53	& LBi, UBi, LBj, UBj, &
54	& IminS, ImaxS, JminS, JmaxS, &
55	& krhs, kstp, kout, &
56	& ubar, vbar, zeta)
57	!***********************************************************************
58	!
59	USE mod_param
60	USE mod_boundary
61	USE mod_clima
62	USE mod_forces
63	USE mod_grid
64	USE mod_ncparam
65	#ifdef NESTING
66	USE mod_nesting
67	#endif
68	USE mod_scalars
69	!
70	! Imported variable declarations.
71	!
72	integer, intent(in) :: ng, tile
73	integer, intent(in) :: LBi, UBi, LBj, UBj
74	integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
75	integer, intent(in) :: krhs, kstp, kout
76	!
77	#ifdef ASSUMED_SHAPE
78	real(r8), intent(in) :: ubar(LBi:,LBj:,:)
79	real(r8), intent(in) :: zeta(LBi:,LBj:,:)
80
81	real(r8), intent(inout) :: vbar(LBi:,LBj:,:)
82	#else
83	real(r8), intent(in) :: ubar(LBi:UBi,LBj:UBj,3)
84	real(r8), intent(in) :: zeta(LBi:UBi,LBj:UBj,3)
85
86	real(r8), intent(inout) :: vbar(LBi:UBi,LBj:UBj,3)
87	#endif
88	!
89	! Local variable declarations.
90	!
91	integer :: Jmin, Jmax
92	integer :: i, j, know
93	#ifdef NESTING
94	integer :: Idg, Jdg, cr, dg, m, rg, tnew, told
95	#endif
96
97	real(r8), parameter :: eps = 1.0E-20_r8
98
99	real(r8) :: Ce, Cx, Ze
100	real(r8) :: bry_pgr, bry_cor, bry_str, bry_val
101	real(r8) :: cff, cff1, cff2, cff3, dt2d, dVde, dVdt, dVdx
102	real(r8) :: obc_in, obc_out, phi, tau
103
104	# if defined ATM_PRESS && defined PRESS_COMPENSATE
105	real(r8) :: OneAtm, fac
106	# endif
107
108	real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: grad
109
110	#include "set_bounds.h"
111	!
112	!-----------------------------------------------------------------------
113	! Set time-indices
114	!-----------------------------------------------------------------------
115	!
116	IF (FIRST_2D_STEP) THEN
117	know=krhs
118	dt2d=dtfast(ng)
119	ELSE IF (PREDICTOR_2D_STEP(ng)) THEN
120	know=krhs
121	dt2d=2.0_r8*dtfast(ng)
122	ELSE
123	know=kstp
124	dt2d=dtfast(ng)
125	END IF
126	# if defined ATM_PRESS && defined PRESS_COMPENSATE
127	OneAtm=1013.25_r8 ! 1 atm = 1013.25 mb
128	fac=100.0_r8/(g*rho0)
129	# endif
130	!
131	!-----------------------------------------------------------------------
132	! Lateral boundary conditions at the southern edge.
133	!-----------------------------------------------------------------------
134	!
135	IF (DOMAIN(ng)%Southern_Edge(tile)) THEN
136	!
137	! Southern edge, implicit upstream radiation condition.
138	!
139	IF (LBC(isouth,isVbar,ng)%radiation) THEN
140	DO i=Istr,Iend+1
141	grad(i,Jstr )=vbar(i ,Jstr ,know)- &
142	& vbar(i-1,Jstr ,know)
143	grad(i,Jstr+1)=vbar(i ,Jstr+1,know)- &
144	& vbar(i-1,Jstr+1,know)
145	END DO
146	DO i=Istr,Iend
147	IF (LBC_apply(ng)%south(i)) THEN
148	dVdt=vbar(i,Jstr+1,know)-vbar(i,Jstr+1,kout)
149	dVde=vbar(i,Jstr+1,kout)-vbar(i,Jstr+2,kout)
150
151	IF (LBC(isouth,isVbar,ng)%nudging) THEN
152	IF (LnudgeM2CLM(ng)) THEN
153	obc_out=0.5_r8* &
154	& (CLIMA(ng)%M2nudgcof(i,Jstr-1)+ &
155	& CLIMA(ng)%M2nudgcof(i,Jstr ))
156	obc_in =obcfac(ng)*obc_out
157	ELSE
158	obc_out=M2obc_out(ng,isouth)
159	obc_in =M2obc_in (ng,isouth)
160	END IF
161	IF ((dVdt*dVde).lt.0.0_r8) THEN
162	tau=obc_in
163	ELSE
164	tau=obc_out
165	END IF
166	#ifdef IMPLICIT_NUDGING
167	IF (tau.gt.0.0_r8) tau=1.0_r8/tau
168	#else
169	tau=tau*dt2d
170	#endif
171	END IF
172
173	IF ((dVdt*dVde).lt.0.0_r8) dVdt=0.0_r8
174	IF ((dVdt*(grad(i ,Jstr+1)+ &
175	& grad(i+1,Jstr+1))).gt.0.0_r8) THEN
176	dVdx=grad(i ,Jstr+1)
177	ELSE
178	dVdx=grad(i+1,Jstr+1)
179	END IF
180	cff=MAX(dVdxdVdx+dVdedVde,eps)
181	#ifdef RADIATION_2D
182	Cx=MIN(cff,MAX(dVdt*dVdx,-cff))
183	#else
184	Cx=0.0_r8
185	#endif
186	Ce=dVdt*dVde
187	#if defined CELERITY_WRITE && defined FORWARD_WRITE
188	BOUNDARY(ng)%vbar_south_Cx(i)=Cx
189	BOUNDARY(ng)%vbar_south_Ce(i)=Ce
190	BOUNDARY(ng)%vbar_south_C2(i)=cff
191	#endif
192	vbar(i,Jstr,kout)=(cff*vbar(i,Jstr ,know)+ &
193	& Ce *vbar(i,Jstr+1,kout)- &
194	& MAX(Cx,0.0_r8)*grad(i ,Jstr)- &
195	& MIN(Cx,0.0_r8)*grad(i+1,Jstr))/ &
196	& (cff+Ce)
197
198	IF (LBC(isouth,isVbar,ng)%nudging) THEN
199	#ifdef IMPLICIT_NUDGING
200	phi=dt(ng)/(tau+dt(ng))
201	vbar(i,Jstr,kout)=(1.0_r8-phi)*vbar(i,Jstr,kout)+ &
202	& phi*BOUNDARY(ng)%vbar_south(i)
203
204	#else
205	vbar(i,Jstr,kout)=vbar(i,Jstr,kout)+ &
206	& tau*(BOUNDARY(ng)%vbar_south(i)- &
207	& vbar(i,Jstr,know))
208	#endif
209	END IF
210	#ifdef MASKING
211	vbar(i,Jstr,kout)=vbar(i,Jstr,kout)* &
212	& GRID(ng)%vmask(i,Jstr)
213	#endif
214	END IF
215	END DO
216	!
217	! Southern edge, Flather boundary condition.
218	!
219	ELSE IF (LBC(isouth,isVbar,ng)%Flather) THEN
220	DO i=Istr,Iend
221	IF (LBC_apply(ng)%south(i)) THEN
222	#if defined SSH_TIDES && !defined UV_TIDES
223	IF (LBC(isouth,isFsur,ng)%acquire) THEN
224	bry_pgr=-g*(zeta(i,Jstr,know)- &
225	& BOUNDARY(ng)%zeta_south(i))* &
226	& 0.5_r8*GRID(ng)%pn(i,Jstr)
227	ELSE
228	bry_pgr=-g*(zeta(i,Jstr ,know)- &
229	& zeta(i,Jstr-1,know))* &
230	& 0.5_r8*(GRID(ng)%pn(i,Jstr-1)+ &
231	& GRID(ng)%pn(i,Jstr ))
232	END IF
233	# ifdef UV_COR
234	bry_cor=-0.125_r8*(ubar(i ,Jstr-1,know)+ &
235	& ubar(i+1,Jstr-1,know)+ &
236	& ubar(i ,Jstr ,know)+ &
237	& ubar(i+1,Jstr ,know))* &
238	& (GRID(ng)%f(i,Jstr-1)+ &
239	& GRID(ng)%f(i,Jstr ))
240	# else
241	bry_cor=0.0_r8
242	# endif
243	cff1=1.0_r8/(0.5_r8*(GRID(ng)%h(i,Jstr-1)+ &
244	& zeta(i,Jstr-1,know)+ &
245	& GRID(ng)%h(i,Jstr )+ &
246	& zeta(i,Jstr ,know)))
247	bry_str=cff1*(FORCES(ng)%svstr(i,Jstr)- &
248	& FORCES(ng)%bvstr(i,Jstr))
249	Ce=1.0_r8/SQRT(g0.5_r8(GRID(ng)%h(i,Jstr-1)+ &
250	& zeta(i,Jstr-1,know)+ &
251	& GRID(ng)%h(i,Jstr )+ &
252	& zeta(i,Jstr ,know)))
253	cff2=GRID(ng)%on_v(i,Jstr)*Ce
254	!! cff2=dt2d
255	bry_val=vbar(i,Jstr+1,know)+ &
256	& cff2*(bry_pgr+ &
257	& bry_cor+ &
258	& bry_str)
259	#else
260	bry_val=BOUNDARY(ng)%vbar_south(i)
261	#endif
262	cff=1.0_r8/(0.5_r8*(GRID(ng)%h(i,Jstr-1)+ &
263	& zeta(i,Jstr-1,know)+ &
264	& GRID(ng)%h(i,Jstr )+ &
265	& zeta(i,Jstr ,know)))
266	Ce=SQRT(g*cff)
267	# if defined ATM_PRESS && defined PRESS_COMPENSATE
268	vbar(i,Jstr,kout)= &
269	& bry_val-Ce(0.5_r8(zeta(i,Jstr-1,know)+ &
270	& zeta(i,Jstr ,know)+ &
271	& fac*(FORCES(ng)%Pair(i,Jstr-1)+ &
272	& FORCES(ng)%Pair(i,Jstr )- &
273	& 2.0_r8*OneAtm))- &
274	& BOUNDARY(ng)%zeta_south(i))
275	# else
276	vbar(i,Jstr,kout)=bry_val- &
277	& Ce(0.5_r8(zeta(i,Jstr-1,know)+ &
278	& zeta(i,Jstr ,know))- &
279	& BOUNDARY(ng)%zeta_south(i))
280	# endif
281	#ifdef MASKING
282	vbar(i,Jstr,kout)=vbar(i,Jstr,kout)* &
283	& GRID(ng)%vmask(i,Jstr)
284	#endif
285	END IF
286	END DO
287	!
288	! Southern edge, Shchepetkin boundary condition (Maison et al., 2010).
289	!
290	ELSE IF (LBC(isouth,isVbar,ng)%Shchepetkin) THEN
291	DO i=Istr,Iend
292	IF (LBC_apply(ng)%south(i)) THEN
293	#if defined SSH_TIDES && !defined UV_TIDES
294	IF (LBC(isouth,isFsur,ng)%acquire) THEN
295	bry_pgr=-g*(zeta(i,Jstr,know)- &
296	& BOUNDARY(ng)%zeta_south(i))* &
297	& 0.5_r8*GRID(ng)%pn(i,Jstr)
298	ELSE
299	bry_pgr=-g*(zeta(i,Jstr ,know)- &
300	& zeta(i,Jstr-1,know))* &
301	& 0.5_r8*(GRID(ng)%pn(i,Jstr-1)+ &
302	& GRID(ng)%pn(i,Jstr ))
303	END IF
304	# ifdef UV_COR
305	bry_cor=-0.125_r8*(ubar(i ,Jstr-1,know)+ &
306	& ubar(i+1,Jstr-1,know)+ &
307	& ubar(i ,Jstr ,know)+ &
308	& ubar(i+1,Jstr ,know))* &
309	& (GRID(ng)%f(i,Jstr-1)+ &
310	& GRID(ng)%f(i,Jstr ))
311	# else
312	bry_cor=0.0_r8
313	# endif
314	cff1=1.0_r8/(0.5_r8*(GRID(ng)%h(i,Jstr-1)+ &
315	& zeta(i,Jstr-1,know)+ &
316	& GRID(ng)%h(i,Jstr )+ &
317	& zeta(i,Jstr ,know)))
318	bry_str=cff1*(FORCES(ng)%svstr(i,Jstr)- &
319	& FORCES(ng)%bvstr(i,Jstr))
320	Ce=1.0_r8/SQRT(g0.5_r8(GRID(ng)%h(i,Jstr-1)+ &
321	& zeta(i,Jstr-1,know)+ &
322	& GRID(ng)%h(i,Jstr )+ &
323	& zeta(i,Jstr ,know)))
324	cff2=GRID(ng)%on_v(i,Jstr)*Ce
325	!! cff2=dt2d
326	bry_val=vbar(i,Jstr+1,know)+ &
327	& cff2*(bry_pgr+ &
328	& bry_cor+ &
329	& bry_str)
330	#else
331	bry_val=BOUNDARY(ng)%vbar_south(i)
332	#endif
333	cff=0.5_r8*(GRID(ng)%h(i,Jstr-1)+ &
334	& GRID(ng)%h(i,Jstr ))
335	cff1=SQRT(g/cff)
336	Ce=dt2dcff1cff0.5_r8(GRID(ng)%pn(i,Jstr-1)+ &
337	& GRID(ng)%pn(i,Jstr ))
338	Ze=(0.5_r8+Ce)*zeta(i,Jstr ,know)+ &
339	& (0.5_r8-Ce)*zeta(i,Jstr-1,know)
340	IF (Ce.gt.Co) THEN
341	cff2=(1.0_r8-Co/Ce)**2
342	cff3=zeta(i,Jstr,kout)+ &
343	& Ce*zeta(i,Jstr-1,know)- &
344	& (1.0_r8+Ce)*zeta(i,Jstr,know)
345	Ze=Ze+cff2*cff3
346	END IF
347	vbar(i,Jstr,kout)=0.5_r8* &
348	& ((1.0_r8-Ce)*vbar(i,Jstr,know)+ &
349	& Ce*vbar(i,Jstr+1,know)+ &
350	& bry_val- &
351	& cff1*(Ze-BOUNDARY(ng)%zeta_south(i)))
352	#ifdef MASKING
353	vbar(i,Jstr,kout)=vbar(i,Jstr,kout)* &
354	& GRID(ng)%vmask(i,Jstr)
355	#endif
356	END IF
357	END DO
358	!
359	! Southern edge, clamped boundary condition.
360	!
361	ELSE IF (LBC(isouth,isVbar,ng)%clamped) THEN
362	DO i=Istr,Iend
363	IF (LBC_apply(ng)%south(i)) THEN
364	vbar(i,Jstr,kout)=BOUNDARY(ng)%vbar_south(i)
365	#ifdef MASKING
366	vbar(i,Jstr,kout)=vbar(i,Jstr,kout)* &
367	& GRID(ng)%vmask(i,Jstr)
368	#endif
369	END IF
370	END DO
371	!
372	! Southern edge, gradient boundary condition.
373	!
374	ELSE IF (LBC(isouth,isVbar,ng)%gradient) THEN
375	DO i=Istr,Iend
376	IF (LBC_apply(ng)%south(i)) THEN
377	vbar(i,Jstr,kout)=vbar(i,Jstr+1,kout)
378	#ifdef MASKING
379	vbar(i,Jstr,kout)=vbar(i,Jstr,kout)* &
380	& GRID(ng)%vmask(i,Jstr)
381	#endif
382	END IF
383	END DO
384	!
385	! Southern edge, reduced-physics boundary condition.
386	!
387	ELSE IF (LBC(isouth,isVbar,ng)%reduced) THEN
388	DO i=Istr,Iend
389	IF (LBC_apply(ng)%south(i)) THEN
390	IF (LBC(isouth,isFsur,ng)%acquire) THEN
391	bry_pgr=-g*(zeta(i,Jstr,know)- &
392	& BOUNDARY(ng)%zeta_south(i))* &
393	& 0.5_r8*GRID(ng)%pn(i,Jstr)
394	ELSE
395	bry_pgr=-g*(zeta(i,Jstr ,know)- &
396	& zeta(i,Jstr-1,know))* &
397	& 0.5_r8*(GRID(ng)%pn(i,Jstr-1)+ &
398	& GRID(ng)%pn(i,Jstr ))
399	END IF
400	#ifdef UV_COR
401	bry_cor=-0.125_r8*(ubar(i ,Jstr-1,know)+ &
402	& ubar(i+1,Jstr-1,know)+ &
403	& ubar(i ,Jstr ,know)+ &
404	& ubar(i+1,Jstr ,know))* &
405	& (GRID(ng)%f(i,Jstr-1)+ &
406	& GRID(ng)%f(i,Jstr ))
407	#else
408	bry_cor=0.0_r8
409	#endif
410	cff=1.0_r8/(0.5_r8*(GRID(ng)%h(i,Jstr-1)+ &
411	& zeta(i,Jstr-1,know)+ &
412	& GRID(ng)%h(i,Jstr )+ &
413	& zeta(i,Jstr ,know)))
414	bry_str=cff*(FORCES(ng)%svstr(i,Jstr)- &
415	& FORCES(ng)%bvstr(i,Jstr))
416	vbar(i,Jstr,kout)=vbar(i,Jstr,know)+ &
417	& dt2d*(bry_pgr+ &
418	& bry_cor+ &
419	& bry_str)
420	#ifdef MASKING
421	vbar(i,Jstr,kout)=vbar(i,Jstr,kout)* &
422	& GRID(ng)%vmask(i,Jstr)
423	#endif
424	END IF
425	END DO
426	!
427	! Southern edge, closed boundary condition.
428	!
429	ELSE IF (LBC(isouth,isVbar,ng)%closed) THEN
430	DO i=Istr,Iend
431	IF (LBC_apply(ng)%south(i)) THEN
432	vbar(i,Jstr,kout)=0.0_r8
433	END IF
434	END DO
435
436	#ifdef NESTING
437	!
438	! If refinement grid and southern edge, impose mass flux from donor
439	! coaser grid for volume and mass conservation.
440	!
441	ELSE IF (LBC(isouth,isVbar,ng)%nested) THEN
442	DO cr=1,Ncontact
443	dg=Vcontact(cr)%donor_grid
444	rg=Vcontact(cr)%receiver_grid
445	IF (RefinedGrid(ng).and. &
446	& (rg.eq.ng).and.(DXmax(dg).gt.DXmax(rg))) THEN
447	told=3-RollingIndex(cr)
448	tnew=RollingIndex(cr)
449	DO i=Istr,Iend
450	m=BRY_CONTACT(isouth,cr)%C2Bindex(i)
451	Idg=Vcontact(cr)%Idg(m) ! for debugging
452	Jdg=Vcontact(cr)%Jdg(m) ! purposes
453	cff=0.5_r8GRID(ng)%om_v(i,Jstr) &
454	& (GRID(ng)%h(i,Jstr-1)+zeta(i,Jstr-1,kout)+ &
455	& GRID(ng)%h(i,Jstr )+zeta(i,Jstr ,kout))
456	cff1=GRID(ng)%om_v(i,Jstr)/REFINED(cr)%om_v(m)
457	bry_val=cff1*REFINED(cr)%DV_avg2(1,m,tnew)/cff
458	# ifdef MASKING
459	bry_val=bry_val*GRID(ng)%vmask(i,Jstr)
460	# endif
461	# ifdef NESTING_DEBUG
462	BRY_CONTACT(isouth,cr)%Mflux(i)=cff*bry_val
463	# endif
464	vbar(i,Jstr,kout)=bry_val
465	END DO
466	END IF
467	END DO
468	#endif
469	END IF
470	END IF
471	!
472	!-----------------------------------------------------------------------
473	! Lateral boundary conditions at the northern edge.
474	!-----------------------------------------------------------------------
475	!
476	IF (DOMAIN(ng)%Northern_Edge(tile)) THEN
477	!
478	! Northern edge, implicit upstream radiation condition.
479	!
480	IF (LBC(inorth,isVbar,ng)%radiation) THEN
481	DO i=Istr,Iend+1
482	grad(i,Jend )=vbar(i ,Jend ,know)- &
483	& vbar(i-1,Jend ,know)
484	grad(i,Jend+1)=vbar(i ,Jend+1,know)- &
485	& vbar(i-1,Jend+1,know)
486	END DO
487	DO i=Istr,Iend
488	IF (LBC_apply(ng)%north(i)) THEN
489	dVdt=vbar(i,Jend,know)-vbar(i,Jend ,kout)
490	dVde=vbar(i,Jend,kout)-vbar(i,Jend-1,kout)
491
492	IF (LBC(inorth,isVbar,ng)%nudging) THEN
493	IF (LnudgeM2CLM(ng)) THEN
494	obc_out=0.5_r8* &
495	& (CLIMA(ng)%M2nudgcof(i,Jend )+ &
496	& CLIMA(ng)%M2nudgcof(i,Jend+1))
497	obc_in =obcfac(ng)*obc_out
498	ELSE
499	obc_out=M2obc_out(ng,inorth)
500	obc_in =M2obc_in (ng,inorth)
501	END IF
502	IF ((dVdt*dVde).lt.0.0_r8) THEN
503	tau=obc_in
504	ELSE
505	tau=obc_out
506	END IF
507	#ifdef IMPLICIT_NUDGING
508	IF (tau.gt.0.0_r8) tau=1.0_r8/tau
509	#else
510	tau=tau*dt2d
511	#endif
512	END IF
513
514	IF ((dVdt*dVde).lt.0.0_r8) dVdt=0.0_r8
515	IF ((dVdt*(grad(i ,Jend)+ &
516	& grad(i+1,Jend))).gt.0.0_r8) THEN
517	dVdx=grad(i ,Jend)
518	ELSE
519	dVdx=grad(i+1,Jend)
520	END IF
521	cff=MAX(dVdxdVdx+dVdedVde,eps)
522	#ifdef RADIATION_2D
523	Cx=MIN(cff,MAX(dVdt*dVdx,-cff))
524	#else
525	Cx=0.0_r8
526	#endif
527	Ce=dVdt*dVde
528	#if defined CELERITY_WRITE && defined FORWARD_WRITE
529	BOUNDARY(ng)%vbar_north_Cx(i)=Cx
530	BOUNDARY(ng)%vbar_north_Ce(i)=Ce
531	BOUNDARY(ng)%vbar_north_C2(i)=cff
532	#endif
533	vbar(i,Jend+1,kout)=(cff*vbar(i,Jend+1,know)+ &
534	& Ce *vbar(i,Jend ,kout)- &
535	& MAX(Cx,0.0_r8)*grad(i ,Jend+1)- &
536	& MIN(Cx,0.0_r8)*grad(i+1,Jend+1))/ &
537	& (cff+Ce)
538
539	IF (LBC(inorth,isVbar,ng)%nudging) THEN
540	#ifdef IMPLICIT_NUDGING
541	phi=dt(ng)/(tau+dt(ng))
542	vbar(i,Jend+1,kout)=(1.0_r8-phi)*vbar(i,Jend+1,kout)+ &
543	& phi*BOUNDARY(ng)%vbar_north(i)
544
545	#else
546	vbar(i,Jend+1,kout)=vbar(i,Jend+1,kout)+ &
547	& tau*(BOUNDARY(ng)%vbar_north(i)- &
548	& vbar(i,Jend+1,know))
549	#endif
550	END IF
551	#ifdef MASKING
552	vbar(i,Jend+1,kout)=vbar(i,Jend+1,kout)* &
553	& GRID(ng)%vmask(i,Jend+1)
554	#endif
555	END IF
556	END DO
557	!
558	! Northern edge, Flather boundary condition.
559	!
560	ELSE IF (LBC(inorth,isVbar,ng)%Flather) THEN
561	DO i=Istr,Iend
562	IF (LBC_apply(ng)%north(i)) THEN
563	#if defined SSH_TIDES && !defined UV_TIDES
564	IF (LBC(inorth,isFsur,ng)%acquire) THEN
565	bry_pgr=-g*(BOUNDARY(ng)%zeta_north(i)- &
566	& zeta(i,Jend,know))* &
567	& 0.5_r8*GRID(ng)%pn(i,Jend)
568	ELSE
569	bry_pgr=-g*(zeta(i,Jend+1,know)- &
570	& zeta(i,Jend ,know))* &
571	& 0.5_r8*(GRID(ng)%pn(i,Jend )+ &
572	& GRID(ng)%pn(i,Jend+1))
573	END IF
574	# ifdef UV_COR
575	bry_cor=-0.125_r8*(ubar(i ,Jend ,know)+ &
576	& ubar(i+1,Jend ,know)+ &
577	& ubar(i ,Jend+1,know)+ &
578	& ubar(i+1,Jend+1,know))* &
579	& (GRID(ng)%f(i,Jend )+ &
580	& GRID(ng)%f(i,Jend+1))
581	# else
582	bry_cor=0.0_r8
583	# endif
584	cff1=1.0_r8/(0.5_r8*(GRID(ng)%h(i,Jend )+ &
585	& zeta(i,Jend ,know)+ &
586	& GRID(ng)%h(i,Jend+1)+ &
587	& zeta(i,Jend+1,know)))
588	bry_str=cff1*(FORCES(ng)%svstr(i,Jend+1)- &
589	& FORCES(ng)%bvstr(i,Jend+1))
590	Ce=1.0_r8/SQRT(g0.5_r8(GRID(ng)%h(i,Jend+1)+ &
591	& zeta(i,Jend+1,know)+ &
592	& GRID(ng)%h(i,Jend )+ &
593	& zeta(i,Jend ,know)))
594	cff2=GRID(ng)%on_v(i,Jend+1)*Ce
595	!! cff2=dt2d
596	bry_val=vbar(i,Jend,know)+ &
597	& cff2*(bry_pgr+ &
598	& bry_cor+ &
599	& bry_str)
600	#else
601	bry_val=BOUNDARY(ng)%vbar_north(i)
602	#endif
603	cff=1.0_r8/(0.5_r8*(GRID(ng)%h(i,Jend )+ &
604	& zeta(i,Jend ,know)+ &
605	& GRID(ng)%h(i,Jend+1)+ &
606	& zeta(i,Jend+1,know)))
607	Ce=SQRT(g*cff)
608	# if defined ATM_PRESS && defined PRESS_COMPENSATE
609	vbar(i,Jend+1,kout)= &
610	& bry_val+Ce(0.5_r8(zeta(i,Jend ,know)+ &
611	& zeta(i,Jend+1,know)+ &
612	& fac*(FORCES(ng)%Pair(i,Jend )+ &
613	& FORCES(ng)%Pair(i,Jend+1)- &
614	& 2.0_r8*OneAtm))- &
615	& BOUNDARY(ng)%zeta_north(i))
616	# else
617	vbar(i,Jend+1,kout)=bry_val+ &
618	& Ce(0.5_r8(zeta(i,Jend ,know)+ &
619	& zeta(i,Jend+1,know))- &
620	& BOUNDARY(ng)%zeta_north(i))
621	# endif
622	#ifdef MASKING
623	vbar(i,Jend+1,kout)=vbar(i,Jend+1,kout)* &
624	& GRID(ng)%vmask(i,Jend+1)
625	#endif
626	END IF
627	END DO
628	!
629	! Northern edge, Shchepetkin boundary condition (Maison et al., 2010).
630	!
631	ELSE IF (LBC(inorth,isVbar,ng)%Shchepetkin) THEN
632	DO i=Istr,Iend
633	IF (LBC_apply(ng)%north(i)) THEN
634	#if defined SSH_TIDES && !defined UV_TIDES
635	IF (LBC(inorth,isFsur,ng)%acquire) THEN
636	bry_pgr=-g*(BOUNDARY(ng)%zeta_north(i)- &
637	& zeta(i,Jend,know))* &
638	& 0.5_r8*GRID(ng)%pn(i,Jend)
639	ELSE
640	bry_pgr=-g*(zeta(i,Jend+1,know)- &
641	& zeta(i,Jend ,know))* &
642	& 0.5_r8*(GRID(ng)%pn(i,Jend )+ &
643	& GRID(ng)%pn(i,Jend+1))
644	END IF
645	# ifdef UV_COR
646	bry_cor=-0.125_r8*(ubar(i ,Jend ,know)+ &
647	& ubar(i+1,Jend ,know)+ &
648	& ubar(i ,Jend+1,know)+ &
649	& ubar(i+1,Jend+1,know))* &
650	& (GRID(ng)%f(i,Jend )+ &
651	& GRID(ng)%f(i,Jend+1))
652	# else
653	bry_cor=0.0_r8
654	# endif
655	cff1=1.0_r8/(0.5_r8*(GRID(ng)%h(i,Jend )+ &
656	& zeta(i,Jend ,know)+ &
657	& GRID(ng)%h(i,Jend+1)+ &
658	& zeta(i,Jend+1,know)))
659	bry_str=cff1*(FORCES(ng)%svstr(i,Jend+1)- &
660	& FORCES(ng)%bvstr(i,Jend+1))
661	Ce=1.0_r8/SQRT(g0.5_r8(GRID(ng)%h(i,Jend+1)+ &
662	& zeta(i,Jend+1,know)+ &
663	& GRID(ng)%h(i,Jend )+ &
664	& zeta(i,Jend ,know)))
665	cff2=GRID(ng)%on_v(i,Jend+1)*Ce
666	!! cff2=dt2d
667	bry_val=vbar(i,Jend,know)+ &
668	& cff2*(bry_pgr+ &
669	& bry_cor+ &
670	& bry_str)
671	#else
672	bry_val=BOUNDARY(ng)%vbar_north(i)
673	#endif
674	cff=0.5_r8*(GRID(ng)%h(i,Jend )+ &
675	& GRID(ng)%h(i,Jend+1))
676	cff1=SQRT(g/cff)
677	Ce=dt2dcff1cff0.5_r8(GRID(ng)%pn(i,Jend )+ &
678	& GRID(ng)%pn(i,Jend+1))
679	Ze=(0.5_r8+Ce)*zeta(i,Jend ,know)+ &
680	& (0.5_r8-Ce)*zeta(i,Jend+1,know)
681	IF (Ce.gt.Co) THEN
682	cff2=(1.0_r8-Co/Ce)**2
683	cff3=zeta(i,Jend,kout)+ &
684	& Ce*zeta(i,Jend+1,know)- &
685	& (1.0_r8+Ce)*zeta(i,Jend,know)
686	Ze=Ze+cff2*cff3
687	END IF
688	vbar(i,Jend+1,kout)=0.5_r8* &
689	& ((1.0_r8-Ce)*vbar(i,Jend+1,know)+ &
690	& Ce*vbar(i,Jend,know)+ &
691	& bry_val+ &
692	& cff1*(Ze-BOUNDARY(ng)%zeta_north(i)))
693	#ifdef MASKING
694	vbar(i,Jend+1,kout)=vbar(i,Jend+1,kout)* &
695	& GRID(ng)%vmask(i,Jend+1)
696	#endif
697	END IF
698	END DO
699	!
700	! Northern edge, clamped boundary condition.
701	!
702	ELSE IF (LBC(inorth,isVbar,ng)%clamped) THEN
703	DO i=Istr,Iend
704	IF (LBC_apply(ng)%north(i)) THEN
705	vbar(i,Jend+1,kout)=BOUNDARY(ng)%vbar_north(i)
706	#ifdef MASKING
707	vbar(i,Jend+1,kout)=vbar(i,Jend+1,kout)* &
708	& GRID(ng)%vmask(i,Jend+1)
709	#endif
710	END IF
711	END DO
712	!
713	! Northern edge, gradient boundary condition.
714	!
715	ELSE IF (LBC(inorth,isVbar,ng)%gradient) THEN
716	DO i=Istr,Iend
717	IF (LBC_apply(ng)%north(i)) THEN
718	vbar(i,Jend+1,kout)=vbar(i,Jend,kout)
719	#ifdef MASKING
720	vbar(i,Jend+1,kout)=vbar(i,Jend+1,kout)* &
721	& GRID(ng)%vmask(i,Jend+1)
722	#endif
723	END IF
724	END DO
725	!
726	! Northern edge, reduced-physics boundary condition.
727	!
728	ELSE IF (LBC(inorth,isVbar,ng)%reduced) THEN
729	DO i=Istr,Iend
730	IF (LBC_apply(ng)%north(i)) THEN
731	IF (LBC(inorth,isFsur,ng)%acquire) THEN
732	bry_pgr=-g*(BOUNDARY(ng)%zeta_north(i)- &
733	& zeta(i,Jend,know))* &
734	& 0.5_r8*GRID(ng)%pn(i,Jend)
735	ELSE
736	bry_pgr=-g*(zeta(i,Jend+1,know)- &
737	& zeta(i,Jend ,know))* &
738	& 0.5_r8*(GRID(ng)%pn(i,Jend )+ &
739	& GRID(ng)%pn(i,Jend+1))
740	END IF
741	#ifdef UV_COR
742	bry_cor=-0.125_r8*(ubar(i ,Jend ,know)+ &
743	& ubar(i+1,Jend ,know)+ &
744	& ubar(i ,Jend+1,know)+ &
745	& ubar(i+1,Jend+1,know))* &
746	& (GRID(ng)%f(i,Jend )+ &
747	& GRID(ng)%f(i,Jend+1))
748	#else
749	bry_cor=0.0_r8
750	#endif
751	cff=1.0_r8/(0.5_r8*(GRID(ng)%h(i,Jend )+ &
752	& zeta(i,Jend ,know)+ &
753	& GRID(ng)%h(i,Jend+1)+ &
754	& zeta(i,Jend+1,know)))
755	bry_str=cff*(FORCES(ng)%svstr(i,Jend+1)- &
756	& FORCES(ng)%bvstr(i,Jend+1))
757	vbar(i,Jend+1,kout)=vbar(i,Jend+1,know)+ &
758	& dt2d*(bry_pgr+ &
759	& bry_cor+ &
760	& bry_str)
761	#ifdef MASKING
762	vbar(i,Jend+1,kout)=vbar(i,Jend+1,kout)* &
763	& GRID(ng)%vmask(i,Jend+1)
764	#endif
765	END IF
766	END DO
767	!
768	! Northern edge, closed boundary condition.
769	!
770	ELSE IF (LBC(inorth,isVbar,ng)%closed) THEN
771	DO i=Istr,Iend
772	IF (LBC_apply(ng)%north(i)) THEN
773	vbar(i,Jend+1,kout)=0.0_r8
774	END IF
775	END DO
776
777	#ifdef NESTING
778	!
779	! If refinement grid and northern edge, impose mass flux from donor
780	! coaser grid for volume and mass conservation.
781	!
782	ELSE IF (LBC(inorth,isVbar,ng)%nested) THEN
783	DO cr=1,Ncontact
784	dg=Vcontact(cr)%donor_grid
785	rg=Vcontact(cr)%receiver_grid
786	IF (RefinedGrid(ng).and. &
787	& (rg.eq.ng).and.(DXmax(dg).gt.DXmax(rg))) THEN
788	told=3-RollingIndex(cr)
789	tnew=RollingIndex(cr)
790	DO i=Istr,Iend
791	m=BRY_CONTACT(inorth,cr)%C2Bindex(i)
792	Idg=Vcontact(cr)%Idg(m) ! for debugging
793	Jdg=Vcontact(cr)%Jdg(m) ! purposes
794	cff=0.5_r8GRID(ng)%om_v(i,Jend+1) &
795	& (GRID(ng)%h(i,Jend+1)+zeta(i,Jend+1,kout)+ &
796	& GRID(ng)%h(i,Jend )+zeta(i,Jend ,kout))
797	cff1=GRID(ng)%om_v(i,Jend+1)/REFINED(cr)%om_v(m)
798	bry_val=cff1*REFINED(cr)%DV_avg2(1,m,tnew)/cff
799	# ifdef MASKING
800	bry_val=bry_val*GRID(ng)%vmask(i,Jend+1)
801	# endif
802	# ifdef NESTING_DEBUG
803	BRY_CONTACT(inorth,cr)%Mflux(i)=cff*bry_val
804	# endif
805	vbar(i,Jend+1,kout)=bry_val
806	END DO
807	END IF
808	END DO
809	#endif
810	END IF
811	END IF
812	!
813	!-----------------------------------------------------------------------
814	! Lateral boundary conditions at the western edge.
815	!-----------------------------------------------------------------------
816	!
817	IF (DOMAIN(ng)%Western_Edge(tile)) THEN
818	!
819	! Western edge, implicit upstream radiation condition.
820	!
821	IF (LBC(iwest,isVbar,ng)%radiation) THEN
822	DO j=JstrV-1,Jend
823	grad(Istr-1,j)=vbar(Istr-1,j+1,know)- &
824	& vbar(Istr-1,j ,know)
825	grad(Istr ,j)=vbar(Istr ,j+1,know)- &
826	& vbar(Istr ,j ,know)
827	END DO
828	DO j=JstrV,Jend
829	IF (LBC_apply(ng)%west(j)) THEN
830	dVdt=vbar(Istr,j,know)-vbar(Istr ,j,kout)
831	dVdx=vbar(Istr,j,kout)-vbar(Istr+1,j,kout)
832
833	IF (LBC(iwest,isVbar,ng)%nudging) THEN
834	IF (LnudgeM2CLM(ng)) THEN
835	obc_out=0.5_r8* &
836	& (CLIMA(ng)%M2nudgcof(Istr-1,j-1)+ &
837	& CLIMA(ng)%M2nudgcof(Istr-1,j ))
838	obc_in =obcfac(ng)*obc_out
839	ELSE
840	obc_out=M2obc_out(ng,iwest)
841	obc_in =M2obc_in (ng,iwest)
842	END IF
843	IF ((dVdt*dVdx).lt.0.0_r8) THEN
844	tau=obc_in
845	ELSE
846	tau=obc_out
847	END IF
848	#ifdef IMPLICIT_NUDGING
849	IF (tau.gt.0.0_r8) tau=1.0_r8/tau
850	#else
851	tau=tau*dt2d
852	#endif
853	END IF
854
855	IF ((dVdt*dVdx).lt.0.0_r8) dVdt=0.0_r8
856	IF ((dVdt*(grad(Istr,j-1)+ &
857	& grad(Istr,j ))).gt.0.0_r8) THEN
858	dVde=grad(Istr,j-1)
859	ELSE
860	dVde=grad(Istr,j )
861	END IF
862	cff=MAX(dVdxdVdx+dVdedVde,eps)
863	Cx=dVdt*dVdx
864	#ifdef RADIATION_2D
865	Ce=MIN(cff,MAX(dVdt*dVde,-cff))
866	#else
867	Ce=0.0_r8
868	#endif
869	#if defined CELERITY_WRITE && defined FORWARD_WRITE
870	BOUNDARY(ng)%vbar_west_Cx(j)=Cx
871	BOUNDARY(ng)%vbar_west_Ce(j)=Ce
872	BOUNDARY(ng)%vbar_west_C2(j)=cff
873	#endif
874	vbar(Istr-1,j,kout)=(cff*vbar(Istr-1,j,know)+ &
875	& Cx *vbar(Istr ,j,kout)- &
876	& MAX(Ce,0.0_r8)*grad(Istr-1,j-1)- &
877	& MIN(Ce,0.0_r8)*grad(Istr-1,j ))/ &
878	& (cff+Cx)
879
880	IF (LBC(iwest,isVbar,ng)%nudging) THEN
881	#ifdef IMPLICIT_NUDGING
882	phi=dt(ng)/(tau+dt(ng))
883	vbar(Istr-1,j,kout)=(1.0_r8-phi)*vbar(Istr-1,j,kout)+ &
884	& phi*BOUNDARY(ng)%vbar_west(j)
885	#else
886	vbar(Istr-1,j,kout)=vbar(Istr-1,j,kout)+ &
887	& tau*(BOUNDARY(ng)%vbar_west(j)- &
888	& vbar(Istr-1,j,know))
889	#endif
890	END IF
891	#ifdef MASKING
892	vbar(Istr-1,j,kout)=vbar(Istr-1,j,kout)* &
893	& GRID(ng)%vmask(Istr-1,j)
894	#endif
895	END IF
896	END DO
897	!
898	! Western edge, Chapman boundary condition.
899	!
900	ELSE IF (LBC(iwest,isVbar,ng)%Flather.or. &
901	& LBC(iwest,isVbar,ng)%reduced.or. &
902	& LBC(iwest,isVbar,ng)%Shchepetkin) THEN
903	DO j=JstrV,Jend
904	IF (LBC_apply(ng)%west(j)) THEN
905	cff=dt2d0.5_r8(GRID(ng)%pm(Istr,j-1)+ &
906	& GRID(ng)%pm(Istr,j ))
907	cff1=SQRT(g0.5_r8(GRID(ng)%h(Istr,j-1)+ &
908	& zeta(Istr,j-1,know)+ &
909	& GRID(ng)%h(Istr,j )+ &
910	& zeta(Istr,j ,know)))
911	Cx=cff*cff1
912	cff2=1.0_r8/(1.0_r8+Cx)
913	vbar(Istr-1,j,kout)=cff2*(vbar(Istr-1,j,know)+ &
914	& Cx*vbar(Istr,j,kout))
915	#ifdef MASKING
916	vbar(Istr-1,j,kout)=vbar(Istr-1,j,kout)* &
917	& GRID(ng)%vmask(Istr-1,j)
918	#endif
919	END IF
920	END DO
921	!
922	! Western edge, clamped boundary condition.
923	!
924	ELSE IF (LBC(iwest,isVbar,ng)%clamped) THEN
925	DO j=JstrV,Jend
926	IF (LBC_apply(ng)%west(j)) THEN
927	vbar(Istr-1,j,kout)=BOUNDARY(ng)%vbar_west(j)
928	#ifdef MASKING
929	vbar(Istr-1,j,kout)=vbar(Istr-1,j,kout)* &
930	& GRID(ng)%vmask(Istr-1,j)
931	#endif
932	END IF
933	END DO
934	!
935	! Western edge, gradient boundary condition.
936	!
937	ELSE IF (LBC(iwest,isVbar,ng)%gradient) THEN
938	DO j=JstrV,Jend
939	IF (LBC_apply(ng)%west(j)) THEN
940	vbar(Istr-1,j,kout)=vbar(Istr,j,kout)
941	#ifdef MASKING
942	vbar(Istr-1,j,kout)=vbar(Istr-1,j,kout)* &
943	& GRID(ng)%vmask(Istr-1,j)
944	#endif
945	END IF
946	END DO
947	!
948	! Western edge, closed boundary condition: free slip (gamma2=1) or
949	! no slip (gamma2=-1).
950	!
951	ELSE IF (LBC(iwest,isVbar,ng)%closed) THEN
952	IF (NSperiodic(ng)) THEN
953	Jmin=JstrV
954	Jmax=Jend
955	ELSE
956	Jmin=Jstr
957	Jmax=JendR
958	END IF
959	DO j=Jmin,Jmax
960	IF (LBC_apply(ng)%west(j)) THEN
961	vbar(Istr-1,j,kout)=gamma2(ng)*vbar(Istr,j,kout)
962	#ifdef MASKING
963	vbar(Istr-1,j,kout)=vbar(Istr-1,j,kout)* &
964	& GRID(ng)%vmask(Istr-1,j)
965	#endif
966	END IF
967	END DO
968	END IF
969	END IF
970	!
971	!-----------------------------------------------------------------------
972	! Lateral boundary conditions at the eastern edge.
973	!-----------------------------------------------------------------------
974	!
975	IF (DOMAIN(ng)%Eastern_Edge(tile)) THEN
976	!
977	! Eastern edge, implicit upstream radiation condition.
978	!
979	IF (LBC(ieast,isVbar,ng)%radiation) THEN
980	DO j=JstrV-1,Jend
981	grad(Iend ,j)=vbar(Iend ,j+1,know)- &
982	& vbar(Iend ,j ,know)
983	grad(Iend+1,j)=vbar(Iend+1,j+1,know)- &
984	& vbar(Iend+1,j ,know)
985	END DO
986	DO j=JstrV,Jend
987	IF (LBC_apply(ng)%east(j)) THEN
988	dVdt=vbar(Iend,j,know)-vbar(Iend ,j,kout)
989	dVdx=vbar(Iend,j,kout)-vbar(Iend-1,j,kout)
990
991	IF (LBC(ieast,isVbar,ng)%nudging) THEN
992	IF (LnudgeM2CLM(ng)) THEN
993	obc_out=0.5_r8* &
994	& (CLIMA(ng)%M2nudgcof(Iend+1,j-1)+ &
995	& CLIMA(ng)%M2nudgcof(Iend+1,j ))
996	obc_in =obcfac(ng)*obc_out
997	ELSE
998	obc_out=M2obc_out(ng,ieast)
999	obc_in =M2obc_in (ng,ieast)
1000	END IF
1001	IF ((dVdt*dVdx).lt.0.0_r8) THEN
1002	tau=obc_in
1003	ELSE
1004	tau=obc_out
1005	END IF
1006	tau=tau*dt2d
1007	END IF
1008
1009	IF ((dVdt*dVdx).lt.0.0_r8) dVdt=0.0_r8
1010	IF ((dVdt*(grad(Iend,j-1)+ &
1011	& grad(Iend,j ))).gt.0.0_r8) THEN
1012	dVde=grad(Iend,j-1)
1013	ELSE
1014	dVde=grad(Iend,j )
1015	END IF
1016	cff=MAX(dVdxdVdx+dVdedVde,eps)
1017	Cx=dVdt*dVdx
1018	#ifdef RADIATION_2D
1019	Ce=MIN(cff,MAX(dVdt*dVde,-cff))
1020	#else
1021	Ce=0.0_r8
1022	#endif
1023	#if defined CELERITY_WRITE && defined FORWARD_WRITE
1024	BOUNDARY(ng)%vbar_east_Cx(j)=Cx
1025	BOUNDARY(ng)%vbar_east_Ce(j)=Ce
1026	BOUNDARY(ng)%vbar_east_C2(j)=cff
1027	#endif
1028	vbar(Iend+1,j,kout)=(cff*vbar(Iend+1,j,know)+ &
1029	& Cx *vbar(Iend ,j,kout)- &
1030	& MAX(Ce,0.0_r8)*grad(Iend+1,j-1)- &
1031	& MIN(Ce,0.0_r8)*grad(Iend+1,j ))/ &
1032	& (cff+Cx)
1033
1034	IF (LBC(ieast,isVbar,ng)%nudging) THEN
1035	vbar(Iend+1,j,kout)=vbar(Iend+1,j,kout)+ &
1036	& tau*(BOUNDARY(ng)%vbar_east(j)- &
1037	& vbar(Iend+1,j,know))
1038	END IF
1039	#ifdef MASKING
1040	vbar(Iend+1,j,kout)=vbar(Iend+1,j,kout)* &
1041	& GRID(ng)%vmask(Iend+1,j)
1042	#endif
1043	END IF
1044	END DO
1045	!
1046	! Eastern edge, Chapman boundary condition.
1047	!
1048	ELSE IF (LBC(ieast,isVbar,ng)%Flather.or. &
1049	& LBC(ieast,isVbar,ng)%reduced.or. &
1050	& LBC(ieast,isVbar,ng)%Shchepetkin) THEN
1051	DO j=JstrV,Jend
1052	IF (LBC_apply(ng)%east(j)) THEN
1053	cff=dt2d0.5_r8(GRID(ng)%pm(Iend,j-1)+ &
1054	& GRID(ng)%pm(Iend,j ))
1055	cff1=SQRT(g0.5_r8(GRID(ng)%h(Iend,j-1)+ &
1056	& zeta(Iend,j-1,know)+ &
1057	& GRID(ng)%h(Iend,j )+ &
1058	& zeta(Iend,j ,know)))
1059	Cx=cff*cff1
1060	cff2=1.0_r8/(1.0_r8+Cx)
1061	vbar(Iend+1,j,kout)=cff2*(vbar(Iend+1,j,know)+ &
1062	& Cx*vbar(Iend,j,kout))
1063	#ifdef MASKING
1064	vbar(Iend+1,j,kout)=vbar(Iend+1,j,kout)* &
1065	& GRID(ng)%vmask(Iend+1,j)
1066	#endif
1067	END IF
1068	END DO
1069	!
1070	! Eastern edge, clamped boundary condition.
1071	!
1072	ELSE IF (LBC(ieast,isVbar,ng)%clamped) THEN
1073	DO j=JstrV,Jend
1074	IF (LBC_apply(ng)%east(j)) THEN
1075	vbar(Iend+1,j,kout)=BOUNDARY(ng)%vbar_east(j)
1076	#ifdef MASKING
1077	vbar(Iend+1,j,kout)=vbar(Iend+1,j,kout)* &
1078	& GRID(ng)%vmask(Iend+1,j)
1079	#endif
1080	END IF
1081	END DO
1082	!
1083	! Eastern edge, gradient boundary condition.
1084	!
1085	ELSE IF (LBC(ieast,isVbar,ng)%gradient) THEN
1086	DO j=JstrV,Jend
1087	IF (LBC_apply(ng)%east(j)) THEN
1088	vbar(Iend+1,j,kout)=vbar(Iend,j,kout)
1089	#ifdef MASKING
1090	vbar(Iend+1,j,kout)=vbar(Iend+1,j,kout)* &
1091	& GRID(ng)%vmask(Iend+1,j)
1092	#endif
1093	END IF
1094	END DO
1095	!
1096	! Eastern edge, closed boundary condition: free slip (gamma2=1) or
1097	! no slip (gamma2=-1).
1098	!
1099	ELSE IF (LBC(ieast,isVbar,ng)%closed) THEN
1100	IF (NSperiodic(ng)) THEN
1101	Jmin=JstrV
1102	Jmax=Jend
1103	ELSE
1104	Jmin=Jstr
1105	Jmax=JendR
1106	END IF
1107	DO j=Jmin,Jmax
1108	IF (LBC_apply(ng)%east(j)) THEN
1109	vbar(Iend+1,j,kout)=gamma2(ng)*vbar(Iend,j,kout)
1110	#ifdef MASKING
1111	vbar(Iend+1,j,kout)=vbar(Iend+1,j,kout)* &
1112	& GRID(ng)%vmask(Iend+1,j)
1113	#endif
1114	END IF
1115	END DO
1116	END IF
1117	END IF
1118	!
1119	!-----------------------------------------------------------------------
1120	! Boundary corners.
1121	!-----------------------------------------------------------------------
1122	!
1123	IF (.not.(EWperiodic(ng).or.NSperiodic(ng))) THEN
1124	IF (DOMAIN(ng)%SouthWest_Corner(tile)) THEN
1125	IF (LBC_apply(ng)%south(Istr-1).and. &
1126	& LBC_apply(ng)%west (Jstr )) THEN
1127	vbar(Istr-1,Jstr,kout)=0.5_r8*(vbar(Istr ,Jstr ,kout)+ &
1128	& vbar(Istr-1,Jstr+1,kout))
1129	END IF
1130	END IF
1131	IF (DOMAIN(ng)%SouthEast_Corner(tile)) THEN
1132	IF (LBC_apply(ng)%south(Iend+1).and. &
1133	& LBC_apply(ng)%east (Jstr )) THEN
1134	vbar(Iend+1,Jstr,kout)=0.5_r8*(vbar(Iend ,Jstr ,kout)+ &
1135	& vbar(Iend+1,Jstr+1,kout))
1136	END IF
1137	END IF
1138	IF (DOMAIN(ng)%NorthWest_Corner(tile)) THEN
1139	IF (LBC_apply(ng)%north(Istr-1).and. &
1140	& LBC_apply(ng)%west (Jend+1)) THEN
1141	vbar(Istr-1,Jend+1,kout)=0.5_r8*(vbar(Istr-1,Jend ,kout)+ &
1142	& vbar(Istr ,Jend+1,kout))
1143	END IF
1144	END IF
1145	IF (DOMAIN(ng)%NorthEast_Corner(tile)) THEN
1146	IF (LBC_apply(ng)%north(Iend+1).and. &
1147	& LBC_apply(ng)%east (Jend+1)) THEN
1148	vbar(Iend+1,Jend+1,kout)=0.5_r8*(vbar(Iend+1,Jend ,kout)+ &
1149	& vbar(Iend ,Jend+1,kout))
1150	END IF
1151	END IF
1152	END IF
1153
1154	#if defined WET_DRY
1155	!
1156	!-----------------------------------------------------------------------
1157	! Impose wetting and drying conditions.
1158	!-----------------------------------------------------------------------
1159	!
1160	IF (.not.EWperiodic(ng)) THEN
1161	IF (DOMAIN(ng)%Western_Edge(tile)) THEN
1162	DO j=JstrV,Jend
1163	IF (LBC_apply(ng)%west(j)) THEN
1164	cff1=ABS(ABS(GRID(ng)%vmask_wet(Istr-1,j))-1.0_r8)
1165	cff2=0.5_r8+DSIGN(0.5_r8,vbar(Istr-1,j,kout))* &
1166	& GRID(ng)%vmask_wet(Istr-1,j)
1167	cff=0.5_r8GRID(ng)%vmask_wet(Istr-1,j)cff1+ &
1168	& cff2*(1.0_r8-cff1)
1169	vbar(Istr,j,kout)=vbar(Istr,j,kout)*cff
1170	END IF
1171	END DO
1172	END IF
1173	IF (DOMAIN(ng)%Eastern_Edge(tile)) THEN
1174	DO j=JstrV,Jend
1175	IF (LBC_apply(ng)%east(j)) THEN
1176	cff1=ABS(ABS(GRID(ng)%vmask_wet(Iend+1,j))-1.0_r8)
1177	cff2=0.5_r8+DSIGN(0.5_r8,vbar(Iend+1,j,kout))* &
1178	& GRID(ng)%vmask_wet(Iend+1,j)
1179	cff=0.5_r8GRID(ng)%vmask_wet(Iend+1,j)cff1+ &
1180	& cff2*(1.0_r8-cff1)
1181	vbar(Iend+1,j,kout)=vbar(Iend+1,j,kout)*cff
1182	END IF
1183	END DO
1184	END IF
1185	END IF
1186
1187	IF (.not.NSperiodic(ng)) THEN
1188	IF (DOMAIN(ng)%Southern_Edge(tile)) THEN
1189	DO i=Istr,Iend
1190	IF (LBC_apply(ng)%south(i)) THEN
1191	cff1=ABS(ABS(GRID(ng)%vmask_wet(i,Jstr))-1.0_r8)
1192	cff2=0.5_r8+DSIGN(0.5_r8,vbar(i,Jstr,kout))* &
1193	& GRID(ng)%vmask_wet(i,Jstr)
1194	cff=0.5_r8GRID(ng)%vmask_wet(i,Jstr)cff1+ &
1195	& cff2*(1.0_r8-cff1)
1196	vbar(i,Jstr,kout)=vbar(i,Jstr,kout)*cff
1197	END IF
1198	END DO
1199	END IF
1200	IF (DOMAIN(ng)%Northern_Edge(tile)) THEN
1201	DO i=Istr,Iend
1202	IF (LBC_apply(ng)%north(i)) THEN
1203	cff1=ABS(ABS(GRID(ng)%vmask_wet(i,Jend+1))-1.0_r8)
1204	cff2=0.5_r8+DSIGN(0.5_r8,vbar(i,Jend+1,kout))* &
1205	& GRID(ng)%vmask_wet(i,Jend+1)
1206	cff=0.5_r8GRID(ng)%vmask_wet(i,Jend+1)cff1+ &
1207	& cff2*(1.0_r8-cff1)
1208	vbar(i,Jend+1,kout)=vbar(i,Jend+1,kout)*cff
1209	END IF
1210	END DO
1211	END IF
1212	END IF
1213
1214	IF (.not.(EWperiodic(ng).or.NSperiodic(ng))) THEN
1215	IF (DOMAIN(ng)%SouthWest_Corner(tile)) THEN
1216	IF (LBC_apply(ng)%south(Istr-1).and. &
1217	& LBC_apply(ng)%west (Jstr )) THEN
1218	cff1=ABS(ABS(GRID(ng)%vmask_wet(Istr-1,Jstr))-1.0_r8)
1219	cff2=0.5_r8+DSIGN(0.5_r8,vbar(Istr-1,Jstr,kout))* &
1220	& GRID(ng)%vmask_wet(Istr-1,Jstr)
1221	cff=0.5_r8GRID(ng)%vmask_wet(Istr-1,Jstr)cff1+ &
1222	& cff2*(1.0_r8-cff1)
1223	vbar(Istr-1,Jstr,kout)=vbar(Istr-1,Jstr,kout)*cff
1224	END IF
1225	END IF
1226	IF (DOMAIN(ng)%SouthEast_Corner(tile)) THEN
1227	IF (LBC_apply(ng)%south(Iend+1).and. &
1228	& LBC_apply(ng)%east (Jstr )) THEN
1229	cff1=ABS(ABS(GRID(ng)%vmask_wet(Iend+1,Jstr))-1.0_r8)
1230	cff2=0.5_r8+DSIGN(0.5_r8,vbar(Iend+1,Jstr,kout))* &
1231	& GRID(ng)%vmask_wet(Iend+1,Jstr)
1232	cff=0.5_r8GRID(ng)%vmask_wet(Iend+1,Jstr)cff1+ &
1233	& cff2*(1.0_r8-cff1)
1234	vbar(Iend+1,Jstr,kout)=vbar(Iend+1,Jstr,kout)*cff
1235	END IF
1236	END IF
1237	IF (DOMAIN(ng)%NorthWest_Corner(tile)) THEN
1238	IF (LBC_apply(ng)%north(Istr-1).and. &
1239	& LBC_apply(ng)%west (Jend+1)) THEN
1240	cff1=ABS(ABS(GRID(ng)%vmask_wet(Istr-1,Jend+1))-1.0_r8)
1241	cff2=0.5_r8+DSIGN(0.5_r8,vbar(Istr-1,Jend+1,kout))* &
1242	& GRID(ng)%vmask_wet(Istr-1,Jend+1)
1243	cff=0.5_r8GRID(ng)%vmask_wet(Istr-1,Jend+1)cff1+ &
1244	& cff2*(1.0_r8-cff1)
1245	vbar(Istr-1,Jend+1,kout)=vbar(Istr-1,Jend+1,kout)*cff
1246	END IF
1247	END IF
1248	IF (DOMAIN(ng)%NorthEast_Corner(tile)) THEN
1249	IF (LBC_apply(ng)%north(Iend+1).and. &
1250	& LBC_apply(ng)%east (Jend+1)) THEN
1251	cff1=ABS(ABS(GRID(ng)%vmask_wet(Iend+1,Jend+1))-1.0_r8)
1252	cff2=0.5_r8+DSIGN(0.5_r8,vbar(Iend+1,Jend+1,kout))* &
1253	& GRID(ng)%vmask_wet(Iend+1,Jend+1)
1254	cff=0.5_r8GRID(ng)%vmask_wet(Iend+1,Jend+1)cff1+ &
1255	& cff2*(1.0_r8-cff1)
1256	vbar(Iend+1,Jend+1,kout)=vbar(Iend+1,Jend+1,kout)*cff
1257	END IF
1258	END IF
1259	END IF
1260	#endif
1261
1262	RETURN
1263
1264	END SUBROUTINE v2dbc_tile
1265	END MODULE v2dbc_mod

Download in other formats:

Original Format