56
57
61#if defined ADJUST_BOUNDARY || defined ADJUST_STFLUX || \
62 defined adjust_wstress
64#endif
65#ifdef DISTRIBUTE
66
68#endif
69
70
71
72 integer, intent(in) :: ng, tile, model
73 integer, intent(in) :: LBi, UBi, LBj, UBj, LBij, UBij
74
75#ifdef ASSUMED_SHAPE
76# ifdef MASKING
77 real(r8), intent(in) :: rmask(LBi:,LBj:)
78 real(r8), intent(in) :: umask(LBi:,LBj:)
79 real(r8), intent(in) :: vmask(LBi:,LBj:)
80# endif
81# ifdef ADJUST_BOUNDARY
82# ifdef SOLVE3D
83 real(r8), intent(in) :: s1_t_obc(LBij:,:,:,:,:)
84 real(r8), intent(in) :: s2_t_obc(LBij:,:,:,:,:)
85 real(r8), intent(in) :: s1_u_obc(LBij:,:,:,:)
86 real(r8), intent(in) :: s2_u_obc(LBij:,:,:,:)
87 real(r8), intent(in) :: s1_v_obc(LBij:,:,:,:)
88 real(r8), intent(in) :: s2_v_obc(LBij:,:,:,:)
89
90 real(r8), intent(inout) :: s3_t_obc(LBij:,:,:,:,:)
91 real(r8), intent(inout) :: s3_u_obc(LBij:,:,:,:)
92 real(r8), intent(inout) :: s3_v_obc(LBij:,:,:,:)
93# endif
94 real(r8), intent(in) :: s1_ubar_obc(LBij:,:,:)
95 real(r8), intent(in) :: s2_ubar_obc(LBij:,:,:)
96 real(r8), intent(in) :: s1_vbar_obc(LBij:,:,:)
97 real(r8), intent(in) :: s2_vbar_obc(LBij:,:,:)
98 real(r8), intent(in) :: s1_zeta_obc(LBij:,:,:)
99 real(r8), intent(in) :: s2_zeta_obc(LBij:,:,:)
100
101 real(r8), intent(inout) :: s3_ubar_obc(LBij:,:,:)
102 real(r8), intent(inout) :: s3_vbar_obc(LBij:,:,:)
103 real(r8), intent(inout) :: s3_zeta_obc(LBij:,:,:)
104# endif
105# ifdef ADJUST_WSTRESS
106 real(r8), intent(in) :: s1_sustr(LBi:,LBj:,:)
107 real(r8), intent(in) :: s2_sustr(LBi:,LBj:,:)
108 real(r8), intent(in) :: s1_svstr(LBi:,LBj:,:)
109 real(r8), intent(in) :: s2_svstr(LBi:,LBj:,:)
110
111 real(r8), intent(inout) :: s3_sustr(LBi:,LBj:,:)
112 real(r8), intent(inout) :: s3_svstr(LBi:,LBj:,:)
113# endif
114# ifdef SOLVE3D
115# ifdef ADJUST_STFLUX
116 real(r8), intent(in) :: s1_tflux(LBi:,LBj:,:,:)
117 real(r8), intent(in) :: s2_tflux(LBi:,LBj:,:,:)
118
119 real(r8), intent(inout) :: s3_tflux(LBi:,LBj:,:,:)
120# endif
121 real(r8), intent(in) :: s1_t(LBi:,LBj:,:,:)
122 real(r8), intent(in) :: s2_t(LBi:,LBj:,:,:)
123 real(r8), intent(in) :: s1_u(LBi:,LBj:,:)
124 real(r8), intent(in) :: s2_u(LBi:,LBj:,:)
125 real(r8), intent(in) :: s1_v(LBi:,LBj:,:)
126 real(r8), intent(in) :: s2_v(LBi:,LBj:,:)
127
128 real(r8), intent(inout) :: s3_t(LBi:,LBj:,:,:)
129 real(r8), intent(inout) :: s3_u(LBi:,LBj:,:)
130 real(r8), intent(inout) :: s3_v(LBi:,LBj:,:)
131# else
132 real(r8), intent(in) :: s1_ubar(LBi:,LBj:)
133 real(r8), intent(in) :: s2_ubar(LBi:,LBj:)
134 real(r8), intent(in) :: s1_vbar(LBi:,LBj:)
135 real(r8), intent(in) :: s2_vbar(LBi:,LBj:)
136
137 real(r8), intent(inout) :: s3_ubar(LBi:,LBj:)
138 real(r8), intent(inout) :: s3_vbar(LBi:,LBj:)
139# endif
140 real(r8), intent(in) :: s1_zeta(LBi:,LBj:)
141 real(r8), intent(in) :: s2_zeta(LBi:,LBj:)
142
143 real(r8), intent(inout) :: s3_zeta(LBi:,LBj:)
144
145#else
146
147# ifdef MASKING
148 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
149 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
150 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
151# endif
152
153# ifdef ADJUST_BOUNDARY
154# ifdef SOLVE3D
155 real(r8), intent(in) :: s1_t_obc(LBij:UBij,N(ng),4, &
156 & Nbrec(ng),NT(ng))
157 real(r8), intent(in) :: s2_t_obc(LBij:UBij,N(ng),4, &
158 & Nbrec(ng),NT(ng))
159 real(r8), intent(in) :: s1_u_obc(LBij:UBij,N(ng),4,Nbrec(ng))
160 real(r8), intent(in) :: s2_u_obc(LBij:UBij,N(ng),4,Nbrec(ng))
161 real(r8), intent(in) :: s1_v_obc(LBij:UBij,N(ng),4,Nbrec(ng))
162 real(r8), intent(in) :: s2_v_obc(LBij:UBij,N(ng),4,Nbrec(ng))
163
164 real(r8), intent(inout) :: s3_t_obc(LBij:UBij,N(ng),4, &
165 & Nbrec(ng),NT(ng))
166 real(r8), intent(inout) :: s3_u_obc(LBij:UBij,N(ng),4,Nbrec(ng))
167 real(r8), intent(inout) :: s3_v_obc(LBij:UBij,N(ng),4,Nbrec(ng))
168# endif
169 real(r8), intent(in) :: s1_ubar_obc(LBij:UBij,4,Nbrec(ng))
170 real(r8), intent(in) :: s2_ubar_obc(LBij:UBij,4,Nbrec(ng))
171 real(r8), intent(in) :: s1_vbar_obc(LBij:UBij,4,Nbrec(ng))
172 real(r8), intent(in) :: s2_vbar_obc(LBij:UBij,4,Nbrec(ng))
173 real(r8), intent(in) :: s1_zeta_obc(LBij:UBij,4,Nbrec(ng))
174 real(r8), intent(in) :: s2_zeta_obc(LBij:UBij,4,Nbrec(ng))
175
176 real(r8), intent(inout) :: s3_ubar_obc(LBij:UBij,4,Nbrec(ng))
177 real(r8), intent(inout) :: s3_vbar_obc(LBij:UBij,4,Nbrec(ng))
178 real(r8), intent(inout) :: s3_zeta_obc(LBij:UBij,4,Nbrec(ng))
179# endif
180# ifdef ADJUST_WSTRESS
181 real(r8), intent(in) :: s1_sustr(LBi:UBi,LBj:UBj,Nfrec(ng))
182 real(r8), intent(in) :: s2_sustr(LBi:UBi,LBj:UBj,Nfrec(ng))
183 real(r8), intent(in) :: s1_svstr(LBi:UBi,LBj:UBj,Nfrec(ng))
184 real(r8), intent(in) :: s2_svstr(LBi:UBi,LBj:UBj,Nfrec(ng))
185
186 real(r8), intent(inout) :: s3_sustr(LBi:UBi,LBj:UBj,Nfrec(ng))
187 real(r8), intent(inout) :: s3_svstr(LBi:UBi,LBj:UBj,Nfrec(ng))
188# endif
189# ifdef SOLVE3D
190# ifdef ADJUST_STFLUX
191 real(r8), intent(in) :: s1_tflux(LBi:UBi,LBj:UBj,Nfrec(ng),NT(ng))
192 real(r8), intent(in) :: s2_tflux(LBi:UBi,LBj:UBj,Nfrec(ng),NT(ng))
193
194 real(r8), intent(inout) :: s3_tflux(LBi:UBi,LBj:UBj, &
195 & Nfrec(ng),NT(ng))
196# endif
197 real(r8), intent(in) :: s1_t(LBi:UBi,LBj:UBj,N(ng),NT(ng))
198 real(r8), intent(in) :: s2_t(LBi:UBi,LBj:UBj,N(ng),NT(ng))
199 real(r8), intent(in) :: s1_u(LBi:UBi,LBj:UBj,N(ng))
200 real(r8), intent(in) :: s2_u(LBi:UBi,LBj:UBj,N(ng))
201 real(r8), intent(in) :: s1_v(LBi:UBi,LBj:UBj,N(ng))
202 real(r8), intent(in) :: s2_v(LBi:UBi,LBj:UBj,N(ng))
203
204 real(r8), intent(inout) :: s3_t(LBi:UBi,LBj:UBj,N(ng),NT(ng))
205 real(r8), intent(inout) :: s3_u(LBi:UBi,LBj:UBj,N(ng))
206 real(r8), intent(inout) :: s3_v(LBi:UBi,LBj:UBj,N(ng))
207# else
208 real(r8), intent(in) :: s1_ubar(LBi:UBi,LBj:UBj)
209 real(r8), intent(in) :: s2_ubar(LBi:UBi,LBj:UBj)
210 real(r8), intent(in) :: s1_vbar(LBi:UBi,LBj:UBj)
211 real(r8), intent(in) :: s2_vbar(LBi:UBi,LBj:UBj)
212
213 real(r8), intent(inout) :: s3_ubar(LBi:UBi,LBj:UBj)
214 real(r8), intent(inout) :: s3_vbar(LBi:UBi,LBj:UBj)
215# endif
216 real(r8), intent(in) :: s1_zeta(LBi:UBi,LBj:UBj)
217 real(r8), intent(in) :: s2_zeta(LBi:UBi,LBj:UBj)
218
219 real(r8), intent(inout) :: s3_zeta(LBi:UBi,LBj:UBj)
220#endif
221
222
223
224 integer :: NSUB, i, j, k
225 integer :: ir, it
226
227 real(r8) :: cff
228
229#include "set_bounds.h"
230
231
232
233
234
235
236
237
238 DO j=jstrt,jendt
239 DO i=istrt,iendt
240 cff=s1_zeta(i,j)*s2_zeta(i,j)
241#ifdef MASKING
242 cff=cff*rmask(i,j)
243#endif
244 s3_zeta(i,j)=cff
245 END DO
246 END DO
247
248#ifdef ADJUST_BOUNDARY
249
250
251
255 &
domain(ng)%Western_Edge(tile))
THEN
256 DO j=jstr,jend
257 cff=s1_zeta_obc(j,
iwest,ir)*s2_zeta_obc(j,
iwest,ir)
258# ifdef MASKING
259 cff=cff*rmask(istr-1,j)
260# endif
261 s3_zeta_obc(j,
iwest,ir)=cff
262 END DO
263 END IF
265 &
domain(ng)%Eastern_Edge(tile))
THEN
266 DO j=jstr,jend
267 cff=s1_zeta_obc(j,
ieast,ir)* &
268 & s2_zeta_obc(j,
ieast,ir)
269# ifdef MASKING
270 cff=cff*rmask(iend+1,j)
271# endif
272 s3_zeta_obc(j,
ieast,ir)=cff
273 END DO
274 END IF
276 &
domain(ng)%Southern_Edge(tile))
THEN
277 DO i=istr,iend
278 cff=s1_zeta_obc(i,
isouth,ir)* &
279 & s2_zeta_obc(i,
isouth,ir)
280# ifdef MASKING
281 cff=cff*rmask(i,jstr-1)
282# endif
283 s3_zeta_obc(i,
isouth,ir)=cff
284 END DO
285 END IF
287 &
domain(ng)%Northern_Edge(tile))
THEN
288 DO i=istr,iend
289 cff=s1_zeta_obc(i,
inorth,ir)* &
290 & s2_zeta_obc(i,
inorth,ir)
291# ifdef MASKING
292 cff=cff*rmask(i,jend+1)
293# endif
294 s3_zeta_obc(i,
inorth,ir)=cff
295 END DO
296 END IF
297 END DO
298 END IF
299#endif
300
301#ifndef SOLVE3D
302
303
304
305 DO j=jstrt,jendt
306 DO i=istrp,iendt
307 cff=s1_ubar(i,j)*s2_ubar(i,j)
308# ifdef MASKING
309 cff=cff*umask(i,j)
310# endif
311 s3_ubar(i,j)=cff
312 END DO
313 END DO
314#endif
315
316#ifdef ADJUST_BOUNDARY
317
318
319
323 &
domain(ng)%Western_Edge(tile))
THEN
324 DO j=jstr,jend
325 cff=s1_ubar_obc(j,
iwest,ir)* &
326 & s2_ubar_obc(j,
iwest,ir)
327# ifdef MASKING
328 cff=cff*umask(istr,j)
329# endif
330 s3_ubar_obc(j,
iwest,ir)=cff
331 END DO
332 END IF
334 &
domain(ng)%Eastern_Edge(tile))
THEN
335 DO j=jstr,jend
336 cff=s1_ubar_obc(j,
ieast,ir)* &
337 & s2_ubar_obc(j,
ieast,ir)
338# ifdef MASKING
339 cff=cff*umask(iend+1,j)
340# endif
341 s3_ubar_obc(j,
ieast,ir)=cff
342 END DO
343 END IF
345 &
domain(ng)%Southern_Edge(tile))
THEN
346 DO i=istru,iend
347 cff=s1_ubar_obc(i,
isouth,ir)* &
348 & s2_ubar_obc(i,
isouth,ir)
349# ifdef MASKING
350 cff=cff*umask(i,jstr-1)
351# endif
352 s3_ubar_obc(i,
isouth,ir)=cff
353 END DO
354 END IF
356 &
domain(ng)%Northern_Edge(tile))
THEN
357 DO i=istru,iend
358 cff=s1_ubar_obc(i,
inorth,ir)* &
359 & s2_ubar_obc(i,
inorth,ir)
360# ifdef MASKING
361 cff=cff*umask(i,jend+1)
362# endif
363 s3_ubar_obc(i,
inorth,ir)=cff
364 END DO
365 END IF
366 END DO
367 END IF
368#endif
369
370#ifndef SOLVE3D
371
372
373
374 DO j=jstrp,jendt
375 DO i=istrt,iendt
376 cff=s1_vbar(i,j)*s2_vbar(i,j)
377# ifdef MASKING
378 cff=cff*vmask(i,j)
379# endif
380 s3_vbar(i,j)=cff
381 END DO
382 END DO
383#endif
384
385#ifdef ADJUST_BOUNDARY
386
387
388
392 &
domain(ng)%Western_Edge(tile))
THEN
393 DO j=jstrv,jend
394 cff=s1_vbar_obc(j,
iwest,ir)* &
395 & s2_vbar_obc(j,
iwest,ir)
396# ifdef MASKING
397 cff=cff*vmask(istr-1,j)
398# endif
399 s3_vbar_obc(j,
iwest,ir)=cff
400 END DO
401 END IF
403 &
domain(ng)%Eastern_Edge(tile))
THEN
404 DO j=jstrv,jend
405 cff=s1_vbar_obc(j,
ieast,ir)* &
406 & s2_vbar_obc(j,
ieast,ir)
407# ifdef MASKING
408 cff=cff*vmask(iend+1,j)
409# endif
410 s3_vbar_obc(j,
ieast,ir)=cff
411 END DO
412 END IF
414 &
domain(ng)%Southern_Edge(tile))
THEN
415 DO i=istr,iend
416 cff=s1_vbar_obc(i,
isouth,ir)* &
417 & s2_vbar_obc(i,
isouth,ir)
418# ifdef MASKING
419 cff=cff*vmask(i,jstr)
420# endif
421 s3_vbar_obc(i,
isouth,ir)=cff
422 END DO
423 END IF
425 &
domain(ng)%Northern_Edge(tile))
THEN
426 DO i=istr,iend
427 cff=s1_vbar_obc(i,
inorth,ir)* &
428 & s2_vbar_obc(i,
inorth,ir)
429# ifdef MASKING
430 cff=cff*vmask(i,jend+1)
431# endif
432 s3_vbar_obc(i,
inorth,ir)=cff
433 END DO
434 END IF
435 END DO
436 END IF
437#endif
438
439#ifdef ADJUST_WSTRESS
440
441
442
444 DO j=jstrt,jendt
445 DO i=istrp,iendt
446 cff=s1_sustr(i,j,ir)*s2_sustr(i,j,ir)
447# ifdef MASKING
448 cff=cff*umask(i,j)
449# endif
450 s3_sustr(i,j,ir)=cff
451 END DO
452 END DO
453 DO j=jstrp,jendt
454 DO i=istrt,iendt
455 cff=s1_svstr(i,j,ir)*s2_svstr(i,j,ir)
456# ifdef MASKING
457 cff=cff*vmask(i,j)
458# endif
459 s3_svstr(i,j,ir)=cff
460 END DO
461 END DO
462 END DO
463#endif
464
465#ifdef SOLVE3D
466
467
468
470 DO j=jstrt,jendt
471 DO i=istrp,iendt
472 cff=s1_u(i,j,k)*s2_u(i,j,k)
473# ifdef MASKING
474 cff=cff*umask(i,j)
475# endif
476 s3_u(i,j,k)=cff
477 END DO
478 END DO
479 END DO
480
481# ifdef ADJUST_BOUNDARY
482
483
484
488 &
domain(ng)%Western_Edge(tile))
THEN
490 DO j=jstr,jend
491 cff=s1_u_obc(j,k,
iwest,ir)* &
492 & s2_u_obc(j,k,
iwest,ir)
493# ifdef MASKING
494 cff=cff*umask(istr,j)
495# endif
496 s3_u_obc(j,k,
iwest,ir)=cff
497 END DO
498 END DO
499 END IF
501 &
domain(ng)%Eastern_Edge(tile))
THEN
503 DO j=jstr,jend
504 cff=s1_u_obc(j,k,
ieast,ir)* &
505 & s2_u_obc(j,k,
ieast,ir)
506# ifdef MASKING
507 cff=cff*umask(iend+1,j)
508# endif
509 s3_u_obc(j,k,
ieast,ir)=cff
510 END DO
511 END DO
512 END IF
514 &
domain(ng)%Southern_Edge(tile))
THEN
516 DO i=istru,iend
517 cff=s1_u_obc(i,k,
isouth,ir)* &
519# ifdef MASKING
520 cff=cff*umask(i,jstr-1)
521# endif
522 s3_u_obc(i,k,
isouth,ir)=cff
523 END DO
524 END DO
525 END IF
527 &
domain(ng)%Northern_Edge(tile))
THEN
529 DO i=istru,iend
530 cff=s1_u_obc(i,k,
inorth,ir)* &
532# ifdef MASKING
533 cff=cff*umask(i,jend+1)
534# endif
535 s3_u_obc(i,k,
inorth,ir)=cff
536 END DO
537 END DO
538 END IF
539 END DO
540 END IF
541# endif
542
543
544
546 DO j=jstrp,jendt
547 DO i=istrt,iendt
548 cff=s1_v(i,j,k)*s2_v(i,j,k)
549# ifdef MASKING
550 cff=cff*vmask(i,j)
551# endif
552 s3_v(i,j,k)=cff
553 END DO
554 END DO
555 END DO
556
557# ifdef ADJUST_BOUNDARY
558
559
560
564 &
domain(ng)%Western_Edge(tile))
THEN
566 DO j=jstrv,jend
567 cff=s1_v_obc(j,k,
iwest,ir)* &
568 & s2_v_obc(j,k,
iwest,ir)
569# ifdef MASKING
570 cff=cff*vmask(istr-1,j)
571# endif
572 s3_v_obc(j,k,
iwest,ir)=cff
573 END DO
574 END DO
575 END IF
577 &
domain(ng)%Eastern_Edge(tile))
THEN
579 DO j=jstrv,jend
580 cff=s1_v_obc(j,k,
ieast,ir)* &
581 & s2_v_obc(j,k,
ieast,ir)
582# ifdef MASKING
583 cff=cff*vmask(iend+1,j)
584# endif
585 s3_v_obc(j,k,
ieast,ir)=cff
586 END DO
587 END DO
588 END IF
590 &
domain(ng)%Southern_Edge(tile))
THEN
592 DO i=istr,iend
593 cff=s1_v_obc(i,k,
isouth,ir)* &
595# ifdef MASKING
596 cff=cff*vmask(i,jstr)
597# endif
598 s3_v_obc(i,k,
isouth,ir)=cff
599 END DO
600 END DO
601 END IF
603 &
domain(ng)%Northern_Edge(tile))
THEN
605 DO i=istr,iend
606 cff=s1_v_obc(i,k,
inorth,ir)* &
608# ifdef MASKING
609 cff=cff*vmask(i,jend+1)
610# endif
611 s3_v_obc(i,k,
inorth,ir)=cff
612 END DO
613 END DO
614 END IF
615 END DO
616 END IF
617# endif
618
619
620
623 DO j=jstrt,jendt
624 DO i=istrt,iendt
625 cff=s1_t(i,j,k,it)*s2_t(i,j,k,it)
626# ifdef MASKING
627 cff=cff*rmask(i,j)
628# endif
629 s3_t(i,j,k,it)=cff
630 END DO
631 END DO
632 END DO
633 END DO
634
635# ifdef ADJUST_BOUNDARY
636
637
638
643 &
domain(ng)%Western_Edge(tile))
THEN
645 DO j=jstr,jend
646 cff=s1_t_obc(j,k,
iwest,ir,it)* &
647 & s2_t_obc(j,k,
iwest,ir,it)
648# ifdef MASKING
649 cff=cff*rmask(istr-1,j)
650# endif
651 s3_t_obc(j,k,
iwest,ir,it)=cff
652 END DO
653 END DO
654 END IF
656 &
domain(ng)%Eastern_Edge(tile))
THEN
658 DO j=jstr,jend
659 cff=s1_t_obc(j,k,
ieast,ir,it)* &
660 & s2_t_obc(j,k,
ieast,ir,it)
661# ifdef MASKING
662 cff=cff*rmask(iend+1,j)
663# endif
664 s3_t_obc(j,k,
ieast,ir,it)=cff
665 END DO
666 END DO
667 END IF
669 &
domain(ng)%Southern_Edge(tile))
THEN
671 DO i=istr,iend
672 cff=s1_t_obc(i,k,
isouth,ir,it)* &
673 & s2_t_obc(i,k,
isouth,ir,it)
674# ifdef MASKING
675 cff=cff*rmask(i,jstr-1)
676# endif
677 s3_t_obc(i,k,
isouth,ir,it)=cff
678 END DO
679 END DO
680 END IF
682 &
domain(ng)%Northern_Edge(tile))
THEN
684 DO i=istr,iend
685 cff=s1_t_obc(i,k,
inorth,ir,it)* &
686 & s2_t_obc(i,k,
inorth,ir,it)
687# ifdef MASKING
688 cff=cff*rmask(i,jend+1)
689# endif
690 s3_t_obc(i,k,
inorth,ir,it)=cff
691 END DO
692 END DO
693 END IF
694 END DO
695 END IF
696 END DO
697# endif
698
699# ifdef ADJUST_STFLUX
700
701
702
706 DO j=jstrt,jendt
707 DO i=istrt,iendt
708 cff=s1_tflux(i,j,ir,it)*s2_tflux(i,j,ir,it)
709# ifdef MASKING
710 cff=cff*rmask(i,j)
711# endif
712 s3_tflux(i,j,ir,it)=cff
713 END DO
714 END DO
715 END DO
716 END IF
717 END DO
718# endif
719
720#endif
721
722 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