57
58
60#if defined ADJUST_BOUNDARY || defined ADJUST_STFLUX || \
61 defined adjust_wstress
64#endif
65
66
67
68 integer, intent(in) :: ng, tile
69 integer, intent(in) :: LBi, UBi, LBj, UBj, LBij, UBij
70 integer, intent(in) :: Linp, Lout
71
72#ifdef ASSUMED_SHAPE
73# ifdef ADJUST_BOUNDARY
74# ifdef SOLVE3D
75 real(r8), intent(in) :: s2_t_obc(LBij:,:,:,:,:,:)
76 real(r8), intent(in) :: s2_u_obc(LBij:,:,:,:,:)
77 real(r8), intent(in) :: s2_v_obc(LBij:,:,:,:,:)
78# endif
79 real(r8), intent(in) :: s2_ubar_obc(LBij:,:,:,:)
80 real(r8), intent(in) :: s2_vbar_obc(LBij:,:,:,:)
81 real(r8), intent(in) :: s2_zeta_obc(LBij:,:,:,:)
82# endif
83# ifdef ADJUST_WSTRESS
84 real(r8), intent(in) :: s2_sustr(LBi:,LBj:,:,:)
85 real(r8), intent(in) :: s2_svstr(LBi:,LBj:,:,:)
86# endif
87# ifdef SOLVE3D
88# ifdef ADJUST_STFLUX
89 real(r8), intent(in) :: s2_tflux(LBi:,LBj:,:,:,:)
90# endif
91 real(r8), intent(in) :: s2_t(LBi:,LBj:,:,:,:)
92 real(r8), intent(in) :: s2_u(LBi:,LBj:,:,:)
93 real(r8), intent(in) :: s2_v(LBi:,LBj:,:,:)
94# if defined WEAK_CONSTRAINT && defined TIME_CONV
95 real(r8), intent(in) :: s2_ubar(LBi:,LBj:,:)
96 real(r8), intent(in) :: s2_vbar(LBi:,LBj:,:)
97# endif
98# else
99 real(r8), intent(in) :: s2_ubar(LBi:,LBj:,:)
100 real(r8), intent(in) :: s2_vbar(LBi:,LBj:,:)
101# endif
102 real(r8), intent(in) :: s2_zeta(LBi:,LBj:,:)
103
104# ifdef ADJUST_BOUNDARY
105# ifdef SOLVE3D
106 real(r8), intent(inout) :: s1_t_obc(LBij:,:,:,:,:,:)
107 real(r8), intent(inout) :: s1_u_obc(LBij:,:,:,:,:)
108 real(r8), intent(inout) :: s1_v_obc(LBij:,:,:,:,:)
109# endif
110 real(r8), intent(inout) :: s1_ubar_obc(LBij:,:,:,:)
111 real(r8), intent(inout) :: s1_vbar_obc(LBij:,:,:,:)
112 real(r8), intent(inout) :: s1_zeta_obc(LBij:,:,:,:)
113# endif
114# ifdef ADJUST_WSTRESS
115 real(r8), intent(inout) :: s1_sustr(LBi:,LBj:,:,:)
116 real(r8), intent(inout) :: s1_svstr(LBi:,LBj:,:,:)
117# endif
118# ifdef SOLVE3D
119# ifdef ADJUST_STFLUX
120 real(r8), intent(inout) :: s1_tflux(LBi:,LBj:,:,:,:)
121# endif
122 real(r8), intent(inout) :: s1_t(LBi:,LBj:,:,:,:)
123 real(r8), intent(inout) :: s1_u(LBi:,LBj:,:,:)
124 real(r8), intent(inout) :: s1_v(LBi:,LBj:,:,:)
125# if defined WEAK_CONSTRAINT && defined TIME_CONV
126 real(r8), intent(inout) :: s1_ubar(LBi:,LBj:,:)
127 real(r8), intent(inout) :: s1_vbar(LBi:,LBj:,:)
128# endif
129# else
130 real(r8), intent(inout) :: s1_ubar(LBi:,LBj:,:)
131 real(r8), intent(inout) :: s1_vbar(LBi:,LBj:,:)
132# endif
133 real(r8), intent(inout) :: s1_zeta(LBi:,LBj:,:)
134
135#else
136
137# ifdef ADJUST_BOUNDARY
138# ifdef SOLVE3D
139 real(r8), intent(in) :: s2_t_obc(LBij:UBij,N(ng),4, &
140 & Nbrec(ng),2,NT(ng))
141 real(r8), intent(in) :: s2_u_obc(LBij:UBij,N(ng),4,Nbrec(ng),2)
142 real(r8), intent(in) :: s2_v_obc(LBij:UBij,N(ng),4,Nbrec(ng),2)
143# endif
144 real(r8), intent(in) :: s2_ubar_obc(LBij:UBij,4,Nbrec(ng),2)
145 real(r8), intent(in) :: s2_vbar_obc(LBij:UBij,4,Nbrec(ng),2)
146 real(r8), intent(in) :: s2_zeta_obc(LBij:UBij,4,Nbrec(ng),2)
147# endif
148# ifdef ADJUST_WSTRESS
149 real(r8), intent(in) :: s2_sustr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
150 real(r8), intent(in) :: s2_svstr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
151# endif
152# ifdef SOLVE3D
153# ifdef ADJUST_STFLUX
154 real(r8), intent(in) :: s2_tflux(LBi:UBi,LBj:UBj, &
155 & Nfrec(ng),2,NT(ng))
156# endif
157 real(r8), intent(in) :: s2_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
158 real(r8), intent(in) :: s2_u(LBi:UBi,LBj:UBj,N(ng),2)
159 real(r8), intent(in) :: s2_v(LBi:UBi,LBj:UBj,N(ng),2)
160# if defined WEAK_CONSTRAINT && defined TIME_CONV
161 real(r8), intent(in) :: s2_ubar(LBi:UBi,LBj:UBj,:)
162 real(r8), intent(in) :: s2_vbar(LBi:UBi,LBj:UBj,:)
163# endif
164# else
165 real(r8), intent(in) :: s2_ubar(LBi:UBi,LBj:UBj,:)
166 real(r8), intent(in) :: s2_vbar(LBi:UBi,LBj:UBj,:)
167# endif
168 real(r8), intent(in) :: s2_zeta(LBi:UBi,LBj:UBj,:)
169
170# ifdef ADJUST_BOUNDARY
171# ifdef SOLVE3D
172 real(r8), intent(inout) :: s1_t_obc(LBij:UBij,N(ng),4, &
173 & Nbrec(ng),2,NT(ng))
174 real(r8), intent(inout) :: s1_u_obc(LBij:UBij,N(ng),4,Nbrec(ng),2)
175 real(r8), intent(inout) :: s1_v_obc(LBij:UBij,N(ng),4,Nbrec(ng),2)
176# endif
177 real(r8), intent(inout) :: s1_ubar_obc(LBij:UBij,4,Nbrec(ng),2)
178 real(r8), intent(inout) :: s1_vbar_obc(LBij:UBij,4,Nbrec(ng),2)
179 real(r8), intent(inout) :: s1_zeta_obc(LBij:UBij,4,Nbrec(ng),2)
180# endif
181# ifdef ADJUST_WSTRESS
182 real(r8), intent(inout) :: s1_sustr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
183 real(r8), intent(inout) :: s1_svstr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
184# endif
185# ifdef SOLVE3D
186# ifdef ADJUST_STFLUX
187 real(r8), intent(inout) :: s1_tflux(LBi:UBi,LBj:UBj, &
188 & Nfrec(ng),2,NT(ng))
189# endif
190 real(r8), intent(inout) :: s1_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
191 real(r8), intent(inout) :: s1_u(LBi:UBi,LBj:UBj,N(ng),2)
192 real(r8), intent(inout) :: s1_v(LBi:UBi,LBj:UBj,N(ng),2)
193# if defined WEAK_CONSTRAINT && defined TIME_CONV
194 real(r8), intent(inout) :: s1_ubar(LBi:UBi,LBj:UBj,:)
195 real(r8), intent(inout) :: s1_vbar(LBi:UBi,LBj:UBj,:)
196# endif
197# else
198 real(r8), intent(inout) :: s1_ubar(LBi:UBi,LBj:UBj,:)
199 real(r8), intent(inout) :: s1_vbar(LBi:UBi,LBj:UBj,:)
200# endif
201 real(r8), intent(inout) :: s1_zeta(LBi:UBi,LBj:UBj,:)
202#endif
203
204
205
206 integer :: i, j, k
207 integer :: ib, ir, it
208
209#include "set_bounds.h"
210
211
212
213
214
215
216
217
218
219 DO j=jstrt,jendt
220 DO i=istrt,iendt
221 s1_zeta(i,j,lout)=s2_zeta(i,j,linp)
222 END DO
223 END DO
224
225#ifdef ADJUST_BOUNDARY
226
227
228
232 &
domain(ng)%Western_Edge(tile))
THEN
234 DO j=jstr,jend
235 s1_zeta_obc(j,ib,ir,lout)=s2_zeta_obc(j,ib,ir,linp)
236 END DO
237 END IF
239 &
domain(ng)%Eastern_Edge(tile))
THEN
241 DO j=jstr,jend
242 s1_zeta_obc(j,ib,ir,lout)=s2_zeta_obc(j,ib,ir,linp)
243 END DO
244 END IF
246 &
domain(ng)%Southern_Edge(tile))
THEN
248 DO i=istr,iend
249 s1_zeta_obc(i,ib,ir,lout)=s2_zeta_obc(i,ib,ir,linp)
250 END DO
251 END IF
253 &
domain(ng)%Northern_Edge(tile))
THEN
255 DO i=istr,iend
256 s1_zeta_obc(i,ib,ir,lout)=s2_zeta_obc(i,ib,ir,linp)
257 END DO
258 END IF
259 END DO
260 END IF
261#endif
262
263#if !defined SOLVE3D || (defined WEAK_CONSTRAINT && defined TIME_CONV)
264
265
266
267 DO j=jstrt,jendt
268 DO i=istrp,iendt
269 s1_ubar(i,j,lout)=s2_ubar(i,j,linp)
270 END DO
271 END DO
272#endif
273
274#ifdef ADJUST_BOUNDARY
275
276
277
281 &
domain(ng)%Western_Edge(tile))
THEN
283 DO j=jstr,jend
284 s1_ubar_obc(j,ib,ir,lout)=s2_ubar_obc(j,ib,ir,linp)
285 END DO
286 END IF
288 &
domain(ng)%Eastern_Edge(tile))
THEN
290 DO j=jstr,jend
291 s1_ubar_obc(j,ib,ir,lout)=s2_ubar_obc(j,ib,ir,linp)
292 END DO
293 END IF
295 &
domain(ng)%Southern_Edge(tile))
THEN
297 DO i=istru,iend
298 s1_ubar_obc(i,ib,ir,lout)=s2_ubar_obc(i,ib,ir,linp)
299 END DO
300 END IF
302 &
domain(ng)%Northern_Edge(tile))
THEN
304 DO i=istru,iend
305 s1_ubar_obc(i,ib,ir,lout)=s2_ubar_obc(i,ib,ir,linp)
306 END DO
307 END IF
308 END DO
309 END IF
310#endif
311
312#if !defined SOLVE3D || (defined WEAK_CONSTRAINT && defined TIME_CONV)
313
314
315
316 DO j=jstrp,jendt
317 DO i=istrt,iendt
318 s1_vbar(i,j,lout)=s2_vbar(i,j,linp)
319 END DO
320 END DO
321#endif
322
323#ifdef ADJUST_BOUNDARY
324
325
326
330 &
domain(ng)%Western_Edge(tile))
THEN
332 DO j=jstrv,jend
333 s1_vbar_obc(j,ib,ir,lout)=s2_vbar_obc(j,ib,ir,linp)
334 END DO
335 END IF
337 &
domain(ng)%Eastern_Edge(tile))
THEN
339 DO j=jstrv,jend
340 s1_vbar_obc(j,ib,ir,lout)=s2_vbar_obc(j,ib,ir,linp)
341 END DO
342 END IF
344 &
domain(ng)%Southern_Edge(tile))
THEN
346 DO i=istr,iend
347 s1_vbar_obc(i,ib,ir,lout)=s2_vbar_obc(i,ib,ir,linp)
348 END DO
349 END IF
351 &
domain(ng)%Northern_Edge(tile))
THEN
353 DO i=istr,iend
354 s1_vbar_obc(i,ib,ir,lout)=s2_vbar_obc(i,ib,ir,linp)
355 END DO
356 END IF
357 END DO
358 END IF
359#endif
360
361#ifdef ADJUST_WSTRESS
362
363
364
366 DO j=jstrt,jendt
367 DO i=istrp,iendt
368 s1_sustr(i,j,ir,lout)=s2_sustr(i,j,ir,linp)
369 END DO
370 END DO
371 DO j=jstrp,jendt
372 DO i=istrt,iendt
373 s1_svstr(i,j,ir,lout)=s2_svstr(i,j,ir,linp)
374 END DO
375 END DO
376 END DO
377#endif
378
379#ifdef SOLVE3D
380
381
382
384 DO j=jstrt,jendt
385 DO i=istrp,iendt
386 s1_u(i,j,k,lout)=s2_u(i,j,k,linp)
387 END DO
388 END DO
389 END DO
390
391# ifdef ADJUST_BOUNDARY
392
393
394
398 &
domain(ng)%Western_Edge(tile))
THEN
401 DO j=jstr,jend
402 s1_u_obc(j,k,ib,ir,lout)=s2_u_obc(j,k,ib,ir,linp)
403 END DO
404 END DO
405 END IF
407 &
domain(ng)%Eastern_Edge(tile))
THEN
410 DO j=jstr,jend
411 s1_u_obc(j,k,ib,ir,lout)=s2_u_obc(j,k,ib,ir,linp)
412 END DO
413 END DO
414 END IF
416 &
domain(ng)%Southern_Edge(tile))
THEN
419 DO i=istru,iend
420 s1_u_obc(i,k,ib,ir,lout)=s2_u_obc(i,k,ib,ir,linp)
421 END DO
422 END DO
423 END IF
425 &
domain(ng)%Northern_Edge(tile))
THEN
428 DO i=istru,iend
429 s1_u_obc(i,k,ib,ir,lout)=s2_u_obc(i,k,ib,ir,linp)
430 END DO
431 END DO
432 END IF
433 END DO
434 END IF
435# endif
436
437
438
440 DO j=jstrp,jendt
441 DO i=istrt,iendt
442 s1_v(i,j,k,lout)=s2_v(i,j,k,linp)
443 END DO
444 END DO
445 END DO
446
447# ifdef ADJUST_BOUNDARY
448
449
450
454 &
domain(ng)%Western_Edge(tile))
THEN
457 DO j=jstrv,jend
458 s1_v_obc(j,k,ib,ir,lout)=s2_v_obc(j,k,ib,ir,linp)
459 END DO
460 END DO
461 END IF
463 &
domain(ng)%Eastern_Edge(tile))
THEN
466 DO j=jstrv,jend
467 s1_v_obc(j,k,ib,ir,lout)=s2_v_obc(j,k,ib,ir,linp)
468 END DO
469 END DO
470 END IF
472 &
domain(ng)%Southern_Edge(tile))
THEN
475 DO i=istr,iend
476 s1_v_obc(i,k,ib,ir,lout)=s2_v_obc(i,k,ib,ir,linp)
477 END DO
478 END DO
479 END IF
481 &
domain(ng)%Northern_Edge(tile))
THEN
484 DO i=istr,iend
485 s1_v_obc(i,k,ib,ir,lout)=s2_v_obc(i,k,ib,ir,linp)
486 END DO
487 END DO
488 END IF
489 END DO
490 END IF
491# endif
492
493
494
497 DO j=jstrt,jendt
498 DO i=istrt,iendt
499 s1_t(i,j,k,lout,it)=s2_t(i,j,k,linp,it)
500 END DO
501 END DO
502 END DO
503 END DO
504
505# ifdef ADJUST_BOUNDARY
506
507
508
513 &
domain(ng)%Western_Edge(tile))
THEN
516 DO j=jstr,jend
517 s1_t_obc(j,k,ib,ir,lout,it)= &
518 & s2_t_obc(j,k,ib,ir,linp,it)
519 END DO
520 END DO
521 END IF
523 &
domain(ng)%Eastern_Edge(tile))
THEN
526 DO j=jstr,jend
527 s1_t_obc(j,k,ib,ir,lout,it)= &
528 & s2_t_obc(j,k,ib,ir,linp,it)
529 END DO
530 END DO
531 END IF
533 &
domain(ng)%Southern_Edge(tile))
THEN
536 DO i=istr,iend
537 s1_t_obc(i,k,ib,ir,lout,it)= &
538 & s2_t_obc(i,k,ib,ir,linp,it)
539 END DO
540 END DO
541 END IF
543 &
domain(ng)%Northern_Edge(tile))
THEN
546 DO i=istr,iend
547 s1_t_obc(i,k,ib,ir,lout,it)= &
548 & s2_t_obc(i,k,ib,ir,linp,it)
549 END DO
550 END DO
551 END IF
552 END DO
553 END IF
554 END DO
555# endif
556
557# ifdef ADJUST_STFLUX
558
559
560
564 DO j=jstrt,jendt
565 DO i=istrt,iendt
566 s1_tflux(i,j,ir,lout,it)=s2_tflux(i,j,ir,linp,it)
567 END DO
568 END DO
569 END DO
570 END IF
571 END DO
572# endif
573
574#endif
575
576 RETURN
integer, dimension(:), allocatable istvar
integer, dimension(:), allocatable n
type(t_domain), dimension(:), allocatable domain
integer, dimension(:), allocatable nt
logical, dimension(:,:,:), allocatable lobc
logical, dimension(:,:), allocatable lstflux
integer, dimension(:), allocatable nfrec
integer, parameter isouth
integer, parameter inorth
integer, dimension(:), allocatable nbrec