Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Blame
Revision Log

source: trunk/v3dbc.F @ 2

Last change on this file since 2 was 2, checked in by pinsard, 17 years ago
initial import from official romsagrif distribution without any svn and CVS directories
File size: 16.7 KB

Line
1	!
2	! $Id: v3dbc.F,v 1.5 2005/11/07 10:38:47 pmarches Exp $
3	!
4	#ifndef CHILD
5	!
6	# include "cppdefs.h"
7	# ifdef SOLVE3D
8	subroutine v3dbc_tile (Istr,Iend,Jstr,Jend,grad)
9	# ifdef AGRIF
10	use AGRIF_Util
11	integer Istr,Iend,Jstr,Jend
12	real grad(PRIVATE_2D_SCRATCH_ARRAY)
13	if (AGRIF_Root()) then
14	call v3dbc_parent_tile (Istr,Iend,Jstr,Jend,grad)
15	else
16	call v3dbc_child_tile (Istr,Iend,Jstr,Jend,grad)
17	c call v3dbc_interp_tile(Istr,Iend,Jstr,Jend)
18	endif
19	return
20	end
21	!
22	! PARENT
23	!
24	subroutine v3dbc_parent_tile (Istr,Iend,Jstr,Jend,grad)
25	# endif
26	!
27	! Set lateral boundary conditions for ETA-component velocity
28	! v(:,:,:,nnew) for the parent grid.
29	!
30	# endif /* SOLVE3D */
31	#else
32	# ifdef SOLVE3D
33	!
34	! CHILD
35	!
36	subroutine v3dbc_child_tile (Istr,Iend,Jstr,Jend,grad)
37	!
38	! Set lateral boundary conditions for ETA-component velocity
39	! v(:,:,:,nnew) for the child grid.
40	!
41	# endif /* SOLVE3D */
42	#endif /* CHILD */
43	#ifdef SOLVE3D
44	!
45	! Common Code
46	!
47	# include "set_obc_definitions.h"
48	!
49	implicit none
50	# include "param.h"
51	# include "grid.h"
52	# include "ocean3d.h"
53	# include "climat.h"
54	# include "scalars.h"
55	# include "boundary.h"
56	integer Istr,Iend,Jstr,Jend, i,j,k
57	real grad(PRIVATE_2D_SCRATCH_ARRAY), cff,eps,
58	& cx,cy, dft,dfx,dfy, tau,tau_in,tau_out
59	parameter (eps=1.E-20)
60	!
61	# include "compute_auxiliary_bounds.h"
62	!
63	! Interpolations of the parent values to get vbry_east or vclm
64	!
65	# ifdef CHILD
66	call v3dbc_interp_tile(Istr,Iend,Jstr,Jend)
67	# endif
68	!
69	# if defined M3_FRC_BRY \|\| defined M3NUDGING
70	tau_in=dt*tauM_in
71	tau_out=dt*tauM_out
72	# endif
73	!
74	# ifndef NS_COM_PERIODIC
75	!
76	!====================================================================
77	! SOUTHERN BC
78	!====================================================================
79	if (SOUTHERN_EDGE) then
80	# ifdef OBC_COM_SOUTH
81	# ifdef OBC_COM_M3ORLANSKI
82	do k=1,N ! Southern edge radiation
83	do i=Istr,Iend+1 ! ======== ==== =========
84	grad(i,Jstr )=(v(i,Jstr ,k,nstp)-v(i-1,Jstr ,k,nstp))
85	# ifdef MASKING
86	& *pmask(i,Jstr)
87	# endif
88	grad(i,Jstr+1)=(v(i,Jstr+1,k,nstp)-v(i-1,Jstr+1,k,nstp))
89	# ifdef MASKING
90	& *pmask(i,Jstr+1)
91	# endif
92	enddo
93	do i=Istr,Iend
94	dft=v(i,Jstr+1,k,nstp)-v(i,Jstr+1,k,nnew)
95	dfx=v(i,Jstr+1,k,nnew)-v(i,Jstr+2,k,nnew)
96
97	if (dfx*dft .lt. 0.) then
98	dft=0. ! <-- cancel cx, if inflow
99	# if defined M3_FRC_BRY \|\| defined M3NUDGING
100	tau=tau_in
101	else
102	tau=tau_out
103	# endif
104	endif
105
106	if (dft*(grad(i,Jstr+1)+grad(i+1,Jstr+1)) .gt. 0.) then
107	dfy=grad(i,Jstr+1)
108	else
109	dfy=grad(i+1,Jstr+1)
110	endif
111
112	# ifdef OBC_COM_RAD_NORMAL
113	dfy=0.
114	# endif
115	cff=max(dfxdfx+dfydfy, eps)
116	cx=dft*dfx
117	# ifdef OBC_COM_RAD_NPO
118	cy=0.
119	# else
120	cy=min(cff,max(dft*dfy,-cff))
121	# endif
122
123	v(i,Jstr,k,nnew)=( cff*v(i,Jstr,k,nstp)
124	& +cx*v(i,Jstr+1,k,nnew)
125	& -max(cy,0.)*grad(i ,Jstr)
126	& -min(cy,0.)*grad(i+1,Jstr)
127	& )/(cff+cx)
128	# if defined M3_FRC_BRY \|\| defined M3NUDGING
129	v(i,Jstr,k,nnew)=(1.-tau)v(i,Jstr,k,nnew)+tau
130	# ifdef M3_FRC_BRY
131	& vbry_south(i,k)
132	# else
133	& vclm(i,Jstr,k)
134	# endif
135	# endif
136	# ifdef MASKING
137	v(i,Jstr,k,nnew)=v(i,Jstr,k,nnew)*vmask(i,Jstr)
138	# endif
139	enddo
140	enddo
141	# elif defined OBC_COM_M3SPECIFIED
142	! Southern edge Specified BC
143	! ======== ==== ========= ==
144	do k=1,N
145	do i=Istr,Iend
146	# ifdef M3_FRC_BRY
147	v(i,Jstr,k,nnew)=vbry_south(i,k) ! specified
148	# else
149	v(i,Jstr,k,nnew)=vclm(i,Jstr,k)
150	# endif
151	# ifdef MASKING
152	& *vmask(i,Jstr)
153	# endif
154	enddo
155	enddo
156	# else
157	! Southern edge gradient BC
158	! ======== ==== ======== ==
159	do k=1,N
160	do i=Istr,Iend
161	v(i,Jstr,k,nnew)=v(i,Jstr+1,k,nnew) ! gradient (default)
162	# ifdef MASKING
163	& *vmask(i,Jstr)
164	# endif
165	enddo
166	enddo
167	# endif
168	# else
169	do k=1,N ! Southern edge closed
170	do i=Istr,Iend ! ======== ==== ======
171	v(i,Jstr,k,nnew)=0. ! (no-flux: default)
172	enddo
173	enddo
174	# endif /* OBC_COM_SOUTH */
175	endif !<-- SOUTHERN_EDGE
176	!
177	!====================================================================
178	! NORTHERN BC
179	!====================================================================
180	if (NORTHERN_EDGE) then
181	# ifdef OBC_COM_NORTH
182	# ifdef OBC_COM_M3ORLANSKI
183	do k=1,N ! Northern edge radiation
184	do i=Istr,Iend+1 ! ======== ==== =========
185	grad(i,Jend )=(v(i,Jend ,k,nstp)-v(i-1,Jend ,k,nstp))
186	# ifdef MASKING
187	& *pmask(i,Jend)
188	# endif
189	grad(i,Jend+1)=(v(i,Jend+1,k,nstp)-v(i-1,Jend+1,k,nstp))
190	# ifdef MASKING
191	& *pmask(i,Jend+1)
192	# endif
193	enddo
194	do i=Istr,Iend
195	dft=v(i,Jend,k,nstp)-v(i,Jend ,k,nnew)
196	dfx=v(i,Jend,k,nnew)-v(i,Jend-1,k,nnew)
197
198	if (dfx*dft .lt. 0.) then
199	dft=0. ! <-- cancel cx, if inflow
200	# if defined M3_FRC_BRY \|\| defined M3NUDGING
201	tau=tau_in
202	else
203	tau=tau_out
204	# endif
205	endif
206
207	if (dft*(grad(i,Jend)+grad(i+1,Jend)) .gt. 0.) then
208	dfy=grad(i,Jend)
209	else
210	dfy=grad(i+1,Jend)
211	endif
212
213	# ifdef OBC_COM_RAD_NORMAL
214	dfy=0.
215	# endif
216	cff=max(dfxdfx+dfydfy, eps)
217	cx=dft*dfx
218	# ifdef OBC_COM_RAD_NPO
219	cy=0.
220	# else
221	cy=min(cff,max(dft*dfy,-cff))
222	# endif
223
224	v(i,Jend+1,k,nnew)=( cff*v(i,Jend+1,k,nstp)
225	& +cx*v(i,Jend,k,nnew)
226	& -max(cy,0.)*grad(i ,Jend+1)
227	& -min(cy,0.)*grad(i+1,Jend+1)
228	& )/(cff+cx)
229	# if defined M3_FRC_BRY \|\| defined M3NUDGING
230	v(i,Jend+1,k,nnew)=(1.-tau)v(i,Jend+1,k,nnew)+tau
231	# ifdef M3_FRC_BRY
232	& vbry_north(i,k)
233	# else
234	& vclm(i,Jend+1,k)
235	# endif
236	# endif
237	# ifdef MASKING
238	v(i,Jend+1,k,nnew)=v(i,Jend+1,k,nnew)*vmask(i,Jend+1)
239	# endif
240	enddo
241	enddo
242	!
243	# elif defined OBC_COM_M3SPECIFIED
244	! Northern edge Specified BC
245	! ======== ==== ========= ==
246	do k=1,N
247	do i=Istr,Iend
248	# ifdef M3_FRC_BRY
249	v(i,Jend+1,k,nnew)=vbry_north(i,k) ! specified
250	# else
251	v(i,Jend+1,k,nnew)=vclm(i,Jend+1,k)
252	# endif
253	# ifdef MASKING
254	& *vmask(i,Jend+1)
255	# endif
256	enddo
257	enddo
258	# else
259	do k=1,N
260	do i=Istr,Iend
261	! Northern edge gradient BC
262	! ======== ==== ======== ==
263	v(i,Jend+1,k,nnew)=v(i,Jend,k,nnew) ! gradient (default)
264	# ifdef MASKING
265	& *vmask(i,Jend+1)
266	# endif
267	enddo
268	enddo
269	# endif
270	# else
271	do k=1,N ! Northern edge closed
272	do i=Istr,Iend ! ======== ==== ======
273	v(i,Jend+1,k,nnew)=0. ! (no-flux: default)
274	enddo
275	enddo
276	# endif
277	endif !<-- NORTHERN_EDGE
278	# endif /* !NS_COM_PERIODIC */
279
280	# ifndef EW_COM_PERIODIC
281	!
282	!====================================================================
283	! WESTERN BC
284	!====================================================================
285	if (WESTERN_EDGE) then
286	# ifdef OBC_COM_WEST
287	# ifdef OBC_COM_M3ORLANSKI
288	do k=1,N ! Western edge radiation
289	do j=JstrV-1,Jend ! ======= ==== =========
290	grad(Istr-1,j)=v(Istr-1,j+1,k,nstp)-v(Istr-1,j,k,nstp)
291	grad(Istr ,j)=v(Istr ,j+1,k,nstp)-v(Istr ,j,k,nstp)
292	enddo
293	do j=JstrV,Jend
294	dft=v(Istr,j,k,nstp)-v(Istr ,j,k,nnew)
295	dfx=v(Istr,j,k,nnew)-v(Istr+1,j,k,nnew)
296
297	if (dfx*dft .lt. 0.) then
298	dft=0. ! <-- cancel cx, if inflow
299	# if defined M3_FRC_BRY \|\| defined M3NUDGING
300	tau=tau_in
301	else
302	tau=tau_out
303	# endif
304	endif
305
306	if (dft*(grad(Istr,j-1)+grad(Istr,j)) .gt. 0.) then
307	dfy=grad(Istr,j-1)
308	else
309	dfy=grad(Istr,j )
310	endif
311
312	# ifdef OBC_COM_RAD_NORMAL
313	dfy=0.
314	# endif
315	cff=max(dfxdfx+dfydfy, eps)
316	cx=dft*dfx
317	# ifdef OBC_COM_RAD_NPO
318	cy=0.
319	# else
320	cy=min(cff,max(dft*dfy,-cff))
321	# endif
322
323	v(Istr-1,j,k,nnew)=( cff*v(Istr-1,j,k,nstp)
324	& +cx*v(Istr,j,k,nnew)
325	& -max(cy,0.)*grad(Istr-1,j-1)
326	& -min(cy,0.)*grad(Istr-1,j )
327	& )/(cff+cx)
328	# if defined M3_FRC_BRY \|\| defined M3NUDGING
329	v(Istr-1,j,k,nnew)=(1.-tau)*v(Istr-1,j,k,nnew)
330	# ifdef M3_FRC_BRY
331	& +tau*vbry_west(j,k)
332	# else
333	& +tau*vclm(Istr-1,j,k)
334	# endif
335	# endif
336	# ifdef MASKING
337	v(Istr-1,j,k,nnew)=v(Istr-1,j,k,nnew)*vmask(Istr-1,j)
338	# endif
339	enddo
340	enddo
341	!
342	# elif defined OBC_COM_M3SPECIFIED
343	! Western edge Specified BC
344	! ======= ==== ========= ==
345	do k=1,N
346	do j=JstrV,Jend
347	# ifdef M3_FRC_BRY
348	v(Istr-1,j,k,nnew)=vbry_west(j,k) ! specified
349	# else
350	v(Istr-1,j,k,nnew)=vclm(Istr-1,j,k)
351	# endif
352	# ifdef MASKING
353	& *vmask(Istr-1,j)
354	# endif
355	enddo
356	enddo
357	# else
358	! Western edge gradient BC
359	! ======= ==== ======== ==
360	do k=1,N
361	do j=JstrV,Jend
362	v(Istr-1,j,k,nnew)=v(Istr,j,k,nnew) ! gradient (default)
363	# ifdef MASKING
364	& *vmask(Istr-1,j)
365	# endif
366	enddo
367	enddo
368	# endif
369	# else
370	# ifdef NS_COM_PERIODIC
371	# define J_RANGE JstrV,Jend
372	# else
373	# define J_RANGE Jstr,JendR
374	# endif
375	do k=1,N ! Wall: free-slip (gamma2=+1)
376	do j=J_RANGE ! ===== no-slip (gamma2=-1)
377	v(Istr-1,j,k,nnew)=gamma2*v(Istr,j,k,nnew)
378	# ifdef MASKING
379	& *vmask(Istr-1,j)
380	# endif
381	enddo
382	enddo
383	# undef J_RANGE
384	# endif
385	endif !<-- WESTERN_EDGE
386	!
387	!====================================================================
388	! EASTERN BC
389	!====================================================================
390	if (EASTERN_EDGE) then
391	# ifdef OBC_COM_EAST
392	# ifdef OBC_COM_M3ORLANSKI
393	do k=1,N ! Eastern edge radiation
394	do j=JstrV-1,Jend ! ======= ==== =========
395	grad(Iend ,j)=v(Iend ,j+1,k,nstp)-v(Iend ,j,k,nstp)
396	grad(Iend+1,j)=v(Iend+1,j+1,k,nstp)-v(Iend+1,j,k,nstp)
397	enddo
398	do j=JstrV,Jend
399	dft=v(Iend,j,k,nstp)-v(Iend ,j,k,nnew)
400	dfx=v(Iend,j,k,nnew)-v(Iend-1,j,k,nnew)
401
402	if (dfx*dft .lt. 0.) then
403	dft=0. ! <-- cancel cx, if inflow
404	# if defined M3_FRC_BRY \|\| defined M3NUDGING
405	tau=tau_in
406	else
407	tau=tau_out
408	# endif
409	endif
410
411	if (dft*(grad(Iend,j-1)+grad(Iend,j)) .gt. 0.) then
412	dfy=grad(Iend,j-1)
413	else
414	dfy=grad(Iend,j )
415	endif
416
417	# ifdef OBC_COM_RAD_NORMAL
418	dfy=0.
419	# endif
420	cff=max(dfxdfx+dfydfy, eps)
421	cx=dft*dfx
422	# ifdef OBC_COM_RAD_NPO
423	cy=0.
424	# else
425	cy=min(cff,max(dft*dfy,-cff))
426	# endif
427
428	v(Iend+1,j,k,nnew)=( cff*v(Iend+1,j,k,nstp)
429	& +cx*v(Iend,j,k,nnew)
430	& -max(cy,0.)*grad(Iend+1,j-1)
431	& -min(cy,0.)*grad(Iend+1,j )
432	& )/(cff+cx)
433	# if defined M3_FRC_BRY \|\| defined M3NUDGING
434	v(Iend+1,j,k,nnew)=(1.-tau)*v(Iend+1,j,k,nnew)
435	# ifdef M3_FRC_BRY
436	& +tau*vbry_east(j,k)
437	# else
438	& +tau*vclm(Iend+1,j,k)
439	# endif
440	# endif
441	# ifdef MASKING
442	v(Iend+1,j,k,nnew)=v(Iend+1,j,k,nnew)*vmask(Iend+1,j)
443	# endif
444	enddo
445	enddo
446	!
447	# elif defined OBC_COM_M3SPECIFIED
448	! Eastern edge Specified BC
449	! ======= ==== ========= ==
450	do k=1,N
451	do j=Jstr,Jend
452	# ifdef M3_FRC_BRY
453	v(Iend+1,j,k,nnew)=vbry_east(j,k) ! specified
454	# else
455	v(Iend+1,j,k,nnew)=vclm(Iend+1,j,k)
456	# endif
457	# ifdef MASKING
458	& *vmask(Iend+1,j)
459	# endif
460	enddo
461	enddo
462	# else
463	! Eastern edge gradient BC
464	! ======= ==== ======== ==
465	do k=1,N
466	do j=Jstr,Jend
467	v(Iend+1,j,k,nnew)=v(Iend,j,k,nnew) ! gradient (default)
468	# ifdef MASKING
469	& *vmask(Iend+1,j)
470	# endif
471	enddo
472	enddo
473	# endif
474	# else
475	# ifdef NS_COM_PERIODIC
476	# define J_RANGE JstrV,Jend
477	# else
478	# define J_RANGE Jstr,JendR
479	# endif
480	do k=1,N ! Wall: free-slip (gamma2=+1)
481	do j=J_RANGE ! ==== no-slip (gamma2=-1)
482	v(Iend+1,j,k,nnew)=gamma2*v(Iend,j,k,nnew)
483	# ifdef MASKING
484	& *vmask(Iend+1,j)
485	# endif
486	enddo
487	enddo
488	# undef J_RANGE
489	# endif
490	endif !<-- EASTERN_EDGE
491	# endif /* !EW_COM_PERIODIC */
492
493	! Corners between adjacent open boundaries
494	! ======= ======= ======== ==== ==========
495
496	# if defined OBC_COM_SOUTH && defined OBC_COM_WEST
497	if (WESTERN_EDGE .and. SOUTHERN_EDGE) then
498	do k=1,N
499	v(Istr-1,Jstr,k,nnew)=0.5*( v(Istr-1,Jstr+1,k,nnew)
500	& +v(Istr ,Jstr ,k,nnew))
501	# ifdef MASKING
502	& *vmask(Istr-1,Jstr)
503	# endif
504	enddo
505	endif
506	# endif
507	# if defined OBC_COM_SOUTH && defined OBC_COM_EAST
508	if (EASTERN_EDGE .and. SOUTHERN_EDGE) then
509	do k=1,N
510	v(Iend+1,Jstr,k,nnew)=0.5*( v(Iend+1,Jstr+1,k,nnew)
511	& +v(Iend ,Jstr ,k,nnew))
512	# ifdef MASKING
513	& *vmask(Iend+1,Jstr)
514	# endif
515	enddo
516	endif
517	# endif
518	# if defined OBC_COM_NORTH && defined OBC_COM_WEST
519	if (WESTERN_EDGE .and. NORTHERN_EDGE) then
520	do k=1,N
521	v(Istr-1,Jend+1,k,nnew)=0.5*( v(Istr-1,Jend,k,nnew)
522	& +v(Istr,Jend+1,k,nnew))
523	# ifdef MASKING
524	& *vmask(Istr-1,Jend+1)
525	# endif
526	enddo
527	endif
528	# endif
529	# if defined OBC_COM_NORTH && defined OBC_COM_EAST
530	if (EASTERN_EDGE .and. NORTHERN_EDGE) then
531	do k=1,N
532	v(Iend+1,Jend+1,k,nnew)=0.5*( v(Iend+1,Jend,k,nnew)
533	& +v(Iend,Jend+1,k,nnew))
534	# ifdef MASKING
535	& *vmask(Iend+1,Jend+1)
536	# endif
537	enddo
538	endif
539	# endif
540	return
541	end
542	#else
543	subroutine v3dbc_empty
544	end
545	#endif /* SOLVE3D */
546	#ifndef CHILD
547	# define CHILD
548	# ifdef AGRIF
549	# include "v3dbc.F"
550	# endif
551	# undef CHILD
552	#endif /* !CHILD */
553

Note: See TracBrowser for help on using the repository browser.

Context Navigation

source: trunk/v3dbc.F @ 2

Download in other formats: