ROMS
Loading...
Searching...
No Matches
state_dotprod.F
Go to the documentation of this file.
1#include "cppdefs.h"
3!
4!git $Id$
5!================================================== Hernan G. Arango ===
6! Copyright (c) 2002-2025 The ROMS Group !
7! Licensed under a MIT/X style license !
8! See License_ROMS.md !
9!=======================================================================
10! !
11! This routine computes the dot product between two model states: !
12! !
13! DotProd(0:NstateVars) = < s1, s2 > !
14! !
15! where !
16! !
17! DotProd(0) All state variable dot product !
18#ifdef ADJUST_WSTRESS
19! DotProd(isUbar) Surface U-momentum stress contribution !
20! DotProd(isVbar) Surface V-momentum stress contribution !
21#endif
22#ifdef SOLVE3D
23# ifdef ADJUST_STFLUX
24! DotProd(isTsur(:)) Surface Tracer flux contribution !
25# endif
26! DotProd(isUvel) 3D U-momentum contribution !
27! DotProd(isVvel) 3D V-momentum contribution !
28! DotProd(isTvar(:)) Tracer-type variables contribution !
29#else
30! DotProd(isUbar) 2D U-momentum contribution !
31! DotProd(isVbar) 2D V-momentum contribution !
32#endif
33! DotProd(isFsur) Free-surface contribution !
34! !
35#ifdef ADJUST_BOUNDARY
36! !
37! Notice that the state variables are processed over the full grid !
38! even when the adjustment of open boundaries is activated. This !
39! is harmless because the S1 and S2 states are originated from the !
40! tangent linear and adjoit models and currently have as zero !
41! value at the boundary edges when data is imposed (clamped). !
42! !
43#endif
44!=======================================================================
45!
46 implicit none
47!
48 PUBLIC :: state_dotprod
49!
50 CONTAINS
51!
52!***********************************************************************
53 SUBROUTINE state_dotprod (ng, tile, model, &
54 & LBi, UBi, LBj, UBj, LBij, UBij, &
55 & NstateVars, DotProd, &
56#ifdef MASKING
57 & rmask, umask, vmask, &
58#endif
59#ifdef ADJUST_BOUNDARY
60# ifdef SOLVE3D
61 & s1_t_obc, s2_t_obc, &
62 & s1_u_obc, s2_u_obc, &
63 & s1_v_obc, s2_v_obc, &
64# endif
65 & s1_ubar_obc, s2_ubar_obc, &
66 & s1_vbar_obc, s2_vbar_obc, &
67 & s1_zeta_obc, s2_zeta_obc, &
68#endif
69#ifdef ADJUST_WSTRESS
70 & s1_sustr, s2_sustr, &
71 & s1_svstr, s2_svstr, &
72#endif
73#ifdef SOLVE3D
74# ifdef ADJUST_STFLUX
75 & s1_tflux, s2_tflux, &
76# endif
77 & s1_t, s2_t, &
78 & s1_u, s2_u, &
79 & s1_v, s2_v, &
80#else
81 & s1_ubar, s2_ubar, &
82 & s1_vbar, s2_vbar, &
83#endif
84 & s1_zeta, s2_zeta)
85!***********************************************************************
86!
87 USE mod_param
88 USE mod_parallel
89 USE mod_ncparam
90#if defined ADJUST_BOUNDARY || defined ADJUST_STFLUX || \
91 defined adjust_wstress
92 USE mod_scalars
93#endif
94#ifdef DISTRIBUTE
95!
96 USE distribute_mod, ONLY : mp_reduce
97#endif
98!
99! Imported variable declarations.
100!
101 integer, intent(in) :: ng, tile, model
102 integer, intent(in) :: lbi, ubi, lbj, ubj, lbij, ubij
103 integer, intent(in) :: nstatevars
104!
105#ifdef ASSUMED_SHAPE
106# ifdef MASKING
107 real(r8), intent(in) :: rmask(lbi:,lbj:)
108 real(r8), intent(in) :: umask(lbi:,lbj:)
109 real(r8), intent(in) :: vmask(lbi:,lbj:)
110# endif
111# ifdef ADJUST_BOUNDARY
112# ifdef SOLVE3D
113 real(r8), intent(in) :: s1_t_obc(lbij:,:,:,:,:)
114 real(r8), intent(in) :: s2_t_obc(lbij:,:,:,:,:)
115 real(r8), intent(in) :: s1_u_obc(lbij:,:,:,:)
116 real(r8), intent(in) :: s2_u_obc(lbij:,:,:,:)
117 real(r8), intent(in) :: s1_v_obc(lbij:,:,:,:)
118 real(r8), intent(in) :: s2_v_obc(lbij:,:,:,:)
119# endif
120 real(r8), intent(in) :: s1_ubar_obc(lbij:,:,:)
121 real(r8), intent(in) :: s2_ubar_obc(lbij:,:,:)
122 real(r8), intent(in) :: s1_vbar_obc(lbij:,:,:)
123 real(r8), intent(in) :: s2_vbar_obc(lbij:,:,:)
124 real(r8), intent(in) :: s1_zeta_obc(lbij:,:,:)
125 real(r8), intent(in) :: s2_zeta_obc(lbij:,:,:)
126# endif
127# ifdef ADJUST_WSTRESS
128 real(r8), intent(in) :: s1_sustr(lbi:,lbj:,:)
129 real(r8), intent(in) :: s2_sustr(lbi:,lbj:,:)
130 real(r8), intent(in) :: s1_svstr(lbi:,lbj:,:)
131 real(r8), intent(in) :: s2_svstr(lbi:,lbj:,:)
132# endif
133# ifdef SOLVE3D
134# ifdef ADJUST_STFLUX
135 real(r8), intent(in) :: s1_tflux(lbi:,lbj:,:,:)
136 real(r8), intent(in) :: s2_tflux(lbi:,lbj:,:,:)
137# endif
138 real(r8), intent(in) :: s1_t(lbi:,lbj:,:,:)
139 real(r8), intent(in) :: s2_t(lbi:,lbj:,:,:)
140 real(r8), intent(in) :: s1_u(lbi:,lbj:,:)
141 real(r8), intent(in) :: s2_u(lbi:,lbj:,:)
142 real(r8), intent(in) :: s1_v(lbi:,lbj:,:)
143 real(r8), intent(in) :: s2_v(lbi:,lbj:,:)
144# else
145 real(r8), intent(in) :: s1_ubar(lbi:,lbj:)
146 real(r8), intent(in) :: s2_ubar(lbi:,lbj:)
147 real(r8), intent(in) :: s1_vbar(lbi:,lbj:)
148 real(r8), intent(in) :: s2_vbar(lbi:,lbj:)
149# endif
150 real(r8), intent(in) :: s1_zeta(lbi:,lbj:)
151 real(r8), intent(in) :: s2_zeta(lbi:,lbj:)
152
153#else
154
155# ifdef MASKING
156 real(r8), intent(in) :: rmask(lbi:ubi,lbj:ubj)
157 real(r8), intent(in) :: umask(lbi:ubi,lbj:ubj)
158 real(r8), intent(in) :: vmask(lbi:ubi,lbj:ubj)
159# endif
160
161# ifdef ADJUST_BOUNDARY
162# ifdef SOLVE3D
163 real(r8), intent(in) :: s1_t_obc(lbij:ubij,n(ng),4, &
164 & Nbrec(ng),NT(ng))
165 real(r8), intent(in) :: s2_t_obc(lbij:ubij,n(ng),4, &
166 & Nbrec(ng),NT(ng))
167 real(r8), intent(in) :: s1_u_obc(lbij:ubij,n(ng),4,nbrec(ng))
168 real(r8), intent(in) :: s2_u_obc(lbij:ubij,n(ng),4,nbrec(ng))
169 real(r8), intent(in) :: s1_v_obc(lbij:ubij,n(ng),4,nbrec(ng))
170 real(r8), intent(in) :: s2_v_obc(lbij:ubij,n(ng),4,nbrec(ng))
171# endif
172 real(r8), intent(in) :: s1_ubar_obc(lbij:ubij,4,nbrec(ng))
173 real(r8), intent(in) :: s2_ubar_obc(lbij:ubij,4,nbrec(ng))
174 real(r8), intent(in) :: s1_vbar_obc(lbij:ubij,4,nbrec(ng))
175 real(r8), intent(in) :: s2_vbar_obc(lbij:ubij,4,nbrec(ng))
176 real(r8), intent(in) :: s1_zeta_obc(lbij:ubij,4,nbrec(ng))
177 real(r8), intent(in) :: s2_zeta_obc(lbij:ubij,4,nbrec(ng))
178# endif
179# ifdef ADJUST_WSTRESS
180 real(r8), intent(in) :: s1_sustr(lbi:ubi,lbj:ubj,nfrec(ng))
181 real(r8), intent(in) :: s2_sustr(lbi:ubi,lbj:ubj,nfrec(ng))
182 real(r8), intent(in) :: s1_svstr(lbi:ubi,lbj:ubj,nfrec(ng))
183 real(r8), intent(in) :: s2_svstr(lbi:ubi,lbj:ubj,nfrec(ng))
184# endif
185# ifdef SOLVE3D
186# ifdef ADJUST_STFLUX
187 real(r8), intent(in) :: s1_tflux(lbi:ubi,lbj:ubj,nfrec(ng),nt(ng))
188 real(r8), intent(in) :: s2_tflux(lbi:ubi,lbj:ubj,nfrec(ng),nt(ng))
189# endif
190 real(r8), intent(in) :: s1_t(lbi:ubi,lbj:ubj,n(ng),nt(ng))
191 real(r8), intent(in) :: s2_t(lbi:ubi,lbj:ubj,n(ng),nt(ng))
192 real(r8), intent(in) :: s1_u(lbi:ubi,lbj:ubj,n(ng))
193 real(r8), intent(in) :: s2_u(lbi:ubi,lbj:ubj,n(ng))
194 real(r8), intent(in) :: s1_v(lbi:ubi,lbj:ubj,n(ng))
195 real(r8), intent(in) :: s2_v(lbi:ubi,lbj:ubj,n(ng))
196# else
197 real(r8), intent(in) :: s1_ubar(lbi:ubi,lbj:ubj)
198 real(r8), intent(in) :: s2_ubar(lbi:ubi,lbj:ubj)
199 real(r8), intent(in) :: s1_vbar(lbi:ubi,lbj:ubj)
200 real(r8), intent(in) :: s2_vbar(lbi:ubi,lbj:ubj)
201# endif
202 real(r8), intent(in) :: s1_zeta(lbi:ubi,lbj:ubj)
203 real(r8), intent(in) :: s2_zeta(lbi:ubi,lbj:ubj)
204#endif
205!
206 real(r8), intent(out), dimension(0:NstateVars) :: dotprod
207!
208! Local variable declarations.
209!
210 integer :: nsub, i, j, k
211 integer :: ir, it
212
213 real(r8) :: cff
214 real(r8), dimension(0:NstateVars) :: my_dotprod
215#ifdef DISTRIBUTE
216 character (len=3), dimension(0:NstateVars) :: op_handle
217#endif
218
219#include "set_bounds.h"
220!
221!-----------------------------------------------------------------------
222! Compute dot product between S1 and S2 model state trajectories.
223!-----------------------------------------------------------------------
224!
225 DO i=0,nstatevars
226 my_dotprod(i)=0.0_r8
227 END DO
228!
229! Free-surface.
230!
231 DO j=jstrt,jendt
232 DO i=istrt,iendt
233 cff=s1_zeta(i,j)*s2_zeta(i,j)
234#ifdef MASKING
235 cff=cff*rmask(i,j)
236#endif
237 my_dotprod(0)=my_dotprod(0)+cff
238 my_dotprod(isfsur)=my_dotprod(isfsur)+cff
239 END DO
240 END DO
241
242#ifdef ADJUST_BOUNDARY
243!
244! Free-surface open boundaries.
245!
246 IF (any(lobc(:,isfsur,ng))) THEN
247 DO ir=1,nbrec(ng)
248 IF ((lobc(iwest,isfsur,ng)).and. &
249 & domain(ng)%Western_Edge(tile)) THEN
250 DO j=jstr,jend
251 cff=s1_zeta_obc(j,iwest,ir)* &
252 & s2_zeta_obc(j,iwest,ir)
253# ifdef MASKING
254 cff=cff*rmask(istr-1,j)
255# endif
256 my_dotprod(0)=my_dotprod(0)+cff
257 my_dotprod(isfsur)=my_dotprod(isfsur)+cff
258 END DO
259 END IF
260 IF ((lobc(ieast,isfsur,ng)).and. &
261 & domain(ng)%Eastern_Edge(tile)) THEN
262 DO j=jstr,jend
263 cff=s1_zeta_obc(j,ieast,ir)* &
264 & s2_zeta_obc(j,ieast,ir)
265# ifdef MASKING
266 cff=cff*rmask(iend+1,j)
267# endif
268 my_dotprod(0)=my_dotprod(0)+cff
269 my_dotprod(isfsur)=my_dotprod(isfsur)+cff
270 END DO
271 END IF
272 IF ((lobc(isouth,isfsur,ng)).and. &
273 & domain(ng)%Southern_Edge(tile)) THEN
274 DO i=istr,iend
275 cff=s1_zeta_obc(i,isouth,ir)* &
276 & s2_zeta_obc(i,isouth,ir)
277# ifdef MASKING
278 cff=cff*rmask(i,jstr-1)
279# endif
280 my_dotprod(0)=my_dotprod(0)+cff
281 my_dotprod(isfsur)=my_dotprod(isfsur)+cff
282 END DO
283 END IF
284 IF ((lobc(inorth,isfsur,ng)).and. &
285 & domain(ng)%Northern_Edge(tile)) THEN
286 DO i=istr,iend
287 cff=s1_zeta_obc(i,inorth,ir)* &
288 & s2_zeta_obc(i,inorth,ir)
289# ifdef MASKING
290 cff=cff*rmask(i,jend+1)
291# endif
292 my_dotprod(0)=my_dotprod(0)+cff
293 my_dotprod(isfsur)=my_dotprod(isfsur)+cff
294 END DO
295 END IF
296 END DO
297 END IF
298#endif
299
300#ifndef SOLVE3D
301!
302! 2D U-momentum component.
303!
304 DO j=jstrt,jendt
305 DO i=istrp,iendt
306 cff=s1_ubar(i,j)*s2_ubar(i,j)
307# ifdef MASKING
308 cff=cff*umask(i,j)
309# endif
310 my_dotprod(0)=my_dotprod(0)+cff
311 my_dotprod(isubar)=my_dotprod(isubar)+cff
312 END DO
313 END DO
314#endif
315
316#ifdef ADJUST_BOUNDARY
317!
318! 2D U-momentum open boundaries.
319!
320 IF (any(lobc(:,isubar,ng))) THEN
321 DO ir=1,nbrec(ng)
322 IF ((lobc(iwest,isubar,ng)).and. &
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 my_dotprod(0)=my_dotprod(0)+cff
331 my_dotprod(isubar)=my_dotprod(isubar)+cff
332 END DO
333 END IF
334 IF ((lobc(ieast,isubar,ng)).and. &
335 & domain(ng)%Eastern_Edge(tile)) THEN
336 DO j=jstr,jend
337 cff=s1_ubar_obc(j,ieast,ir)* &
338 & s2_ubar_obc(j,ieast,ir)
339# ifdef MASKING
340 cff=cff*umask(iend+1,j)
341# endif
342 my_dotprod(0)=my_dotprod(0)+cff
343 my_dotprod(isubar)=my_dotprod(isubar)+cff
344 END DO
345 END IF
346 IF ((lobc(isouth,isubar,ng)).and. &
347 & domain(ng)%Southern_Edge(tile)) THEN
348 DO i=istru,iend
349 cff=s1_ubar_obc(i,isouth,ir)* &
350 & s2_ubar_obc(i,isouth,ir)
351# ifdef MASKING
352 cff=cff*umask(i,jstr-1)
353# endif
354 my_dotprod(0)=my_dotprod(0)+cff
355 my_dotprod(isubar)=my_dotprod(isubar)+cff
356 END DO
357 END IF
358 IF ((lobc(inorth,isubar,ng)).and. &
359 & domain(ng)%Northern_Edge(tile)) THEN
360 DO i=istru,iend
361 cff=s1_ubar_obc(i,inorth,ir)* &
362 & s2_ubar_obc(i,inorth,ir)
363# ifdef MASKING
364 cff=cff*umask(i,jend+1)
365# endif
366 my_dotprod(0)=my_dotprod(0)+cff
367 my_dotprod(isubar)=my_dotprod(isubar)+cff
368 END DO
369 END IF
370 END DO
371 END IF
372#endif
373
374#ifndef SOLVE3D
375!
376! 2D V-momentum component.
377!
378 DO j=jstrp,jendt
379 DO i=istrt,iendt
380 cff=s1_vbar(i,j)*s2_vbar(i,j)
381# ifdef MASKING
382 cff=cff*vmask(i,j)
383# endif
384 my_dotprod(0)=my_dotprod(0)+cff
385 my_dotprod(isvbar)=my_dotprod(isvbar)+cff
386 END DO
387 END DO
388#endif
389
390#ifdef ADJUST_BOUNDARY
391!
392! 2D V-momentum open boundaries.
393!
394 IF (any(lobc(:,isvbar,ng))) THEN
395 DO ir=1,nbrec(ng)
396 IF ((lobc(iwest,isvbar,ng)).and. &
397 & domain(ng)%Western_Edge(tile)) THEN
398 DO j=jstrv,jend
399 cff=s1_vbar_obc(j,iwest,ir)* &
400 & s2_vbar_obc(j,iwest,ir)
401# ifdef MASKING
402 cff=cff*vmask(istr-1,j)
403# endif
404 my_dotprod(0)=my_dotprod(0)+cff
405 my_dotprod(isvbar)=my_dotprod(isvbar)+cff
406 END DO
407 END IF
408 IF ((lobc(ieast,isvbar,ng)).and. &
409 & domain(ng)%Eastern_Edge(tile)) THEN
410 DO j=jstrv,jend
411 cff=s1_vbar_obc(j,ieast,ir)* &
412 & s2_vbar_obc(j,ieast,ir)
413# ifdef MASKING
414 cff=cff*vmask(iend+1,j)
415# endif
416 my_dotprod(0)=my_dotprod(0)+cff
417 my_dotprod(isvbar)=my_dotprod(isvbar)+cff
418 END DO
419 END IF
420 IF ((lobc(isouth,isvbar,ng)).and. &
421 & domain(ng)%Southern_Edge(tile)) THEN
422 DO i=istr,iend
423 cff=s1_vbar_obc(i,isouth,ir)* &
424 & s2_vbar_obc(i,isouth,ir)
425# ifdef MASKING
426 cff=cff*vmask(i,jstr)
427# endif
428 my_dotprod(0)=my_dotprod(0)+cff
429 my_dotprod(isvbar)=my_dotprod(isvbar)+cff
430 END DO
431 END IF
432 IF ((lobc(inorth,isvbar,ng)).and. &
433 & domain(ng)%Northern_Edge(tile)) THEN
434 DO i=istr,iend
435 cff=s1_vbar_obc(i,inorth,ir)* &
436 & s2_vbar_obc(i,inorth,ir)
437# ifdef MASKING
438 cff=cff*vmask(i,jend+1)
439# endif
440 my_dotprod(0)=my_dotprod(0)+cff
441 my_dotprod(isvbar)=my_dotprod(isvbar)+cff
442 END DO
443 END IF
444 END DO
445 END IF
446#endif
447
448#ifdef ADJUST_WSTRESS
449!
450! Surface momentum stress.
451!
452 DO ir=1,nfrec(ng)
453 DO j=jstrt,jendt
454 DO i=istrp,iendt
455 cff=s1_sustr(i,j,ir)*s2_sustr(i,j,ir)
456# ifdef MASKING
457 cff=cff*umask(i,j)
458# endif
459 my_dotprod(0)=my_dotprod(0)+cff
460 my_dotprod(isustr)=my_dotprod(isustr)+cff
461 END DO
462 END DO
463 DO j=jstrp,jendt
464 DO i=istrt,iendt
465 cff=s1_svstr(i,j,ir)*s2_svstr(i,j,ir)
466# ifdef MASKING
467 cff=cff*vmask(i,j)
468# endif
469 my_dotprod(0)=my_dotprod(0)+cff
470 my_dotprod(isvstr)=my_dotprod(isvstr)+cff
471 END DO
472 END DO
473 END DO
474#endif
475
476#ifdef SOLVE3D
477!
478! 3D U-momentum component.
479!
480 DO k=1,n(ng)
481 DO j=jstrt,jendt
482 DO i=istrp,iendt
483 cff=s1_u(i,j,k)*s2_u(i,j,k)
484# ifdef MASKING
485 cff=cff*umask(i,j)
486# endif
487 my_dotprod(0)=my_dotprod(0)+cff
488 my_dotprod(isuvel)=my_dotprod(isuvel)+cff
489 END DO
490 END DO
491 END DO
492
493# ifdef ADJUST_BOUNDARY
494!
495! 3D U-momentum open boundaries.
496!
497 IF (any(lobc(:,isuvel,ng))) THEN
498 DO ir=1,nbrec(ng)
499 IF ((lobc(iwest,isuvel,ng)).and. &
500 & domain(ng)%Western_Edge(tile)) THEN
501 DO k=1,n(ng)
502 DO j=jstr,jend
503 cff=s1_u_obc(j,k,iwest,ir)* &
504 & s2_u_obc(j,k,iwest,ir)
505# ifdef MASKING
506 cff=cff*umask(istr,j)
507# endif
508 my_dotprod(0)=my_dotprod(0)+cff
509 my_dotprod(isuvel)=my_dotprod(isuvel)+cff
510 END DO
511 END DO
512 END IF
513 IF ((lobc(ieast,isuvel,ng)).and. &
514 & domain(ng)%Eastern_Edge(tile)) THEN
515 DO k=1,n(ng)
516 DO j=jstr,jend
517 cff=s1_u_obc(j,k,ieast,ir)* &
518 & s2_u_obc(j,k,ieast,ir)
519# ifdef MASKING
520 cff=cff*umask(iend+1,j)
521# endif
522 my_dotprod(0)=my_dotprod(0)+cff
523 my_dotprod(isuvel)=my_dotprod(isuvel)+cff
524 END DO
525 END DO
526 END IF
527 IF ((lobc(isouth,isuvel,ng)).and. &
528 & domain(ng)%Southern_Edge(tile)) THEN
529 DO k=1,n(ng)
530 DO i=istru,iend
531 cff=s1_u_obc(i,k,isouth,ir)* &
532 & s2_u_obc(i,k,isouth,ir)
533# ifdef MASKING
534 cff=cff*umask(i,jstr-1)
535# endif
536 my_dotprod(0)=my_dotprod(0)+cff
537 my_dotprod(isuvel)=my_dotprod(isuvel)+cff
538 END DO
539 END DO
540 END IF
541 IF ((lobc(inorth,isuvel,ng)).and. &
542 & domain(ng)%Northern_Edge(tile)) THEN
543 DO k=1,n(ng)
544 DO i=istru,iend
545 cff=s1_u_obc(i,k,inorth,ir)* &
546 & s2_u_obc(i,k,inorth,ir)
547# ifdef MASKING
548 cff=cff*umask(i,jend+1)
549# endif
550 my_dotprod(0)=my_dotprod(0)+cff
551 my_dotprod(isuvel)=my_dotprod(isuvel)+cff
552 END DO
553 END DO
554 END IF
555 END DO
556 END IF
557# endif
558!
559! 3D V-momentum component.
560!
561 DO k=1,n(ng)
562 DO j=jstrp,jendt
563 DO i=istrt,iendt
564 cff=s1_v(i,j,k)*s2_v(i,j,k)
565# ifdef MASKING
566 cff=cff*vmask(i,j)
567# endif
568 my_dotprod(0)=my_dotprod(0)+cff
569 my_dotprod(isvvel)=my_dotprod(isvvel)+cff
570 END DO
571 END DO
572 END DO
573
574# ifdef ADJUST_BOUNDARY
575!
576! 3D V-momentum open boundaries.
577!
578 IF (any(lobc(:,isvvel,ng))) THEN
579 DO ir=1,nbrec(ng)
580 IF ((lobc(iwest,isvvel,ng)).and. &
581 & domain(ng)%Western_Edge(tile)) THEN
582 DO k=1,n(ng)
583 DO j=jstrv,jend
584 cff=s1_v_obc(j,k,iwest,ir)* &
585 & s2_v_obc(j,k,iwest,ir)
586# ifdef MASKING
587 cff=cff*vmask(istr-1,j)
588# endif
589 my_dotprod(0)=my_dotprod(0)+cff
590 my_dotprod(isvvel)=my_dotprod(isvvel)+cff
591 END DO
592 END DO
593 END IF
594 IF ((lobc(ieast,isvvel,ng)).and. &
595 & domain(ng)%Eastern_Edge(tile)) THEN
596 DO k=1,n(ng)
597 DO j=jstrv,jend
598 cff=s1_v_obc(j,k,ieast,ir)* &
599 & s2_v_obc(j,k,ieast,ir)
600# ifdef MASKING
601 cff=cff*vmask(iend+1,j)
602# endif
603 my_dotprod(0)=my_dotprod(0)+cff
604 my_dotprod(isvvel)=my_dotprod(isvvel)+cff
605 END DO
606 END DO
607 END IF
608 IF ((lobc(isouth,isvvel,ng)).and. &
609 & domain(ng)%Southern_Edge(tile)) THEN
610 DO k=1,n(ng)
611 DO i=istr,iend
612 cff=s1_v_obc(i,k,isouth,ir)* &
613 & s2_v_obc(i,k,isouth,ir)
614# ifdef MASKING
615 cff=cff*vmask(i,jstr)
616# endif
617 my_dotprod(0)=my_dotprod(0)+cff
618 my_dotprod(isvvel)=my_dotprod(isvvel)+cff
619 END DO
620 END DO
621 END IF
622 IF ((lobc(inorth,isvvel,ng)).and. &
623 & domain(ng)%Northern_Edge(tile)) THEN
624 DO k=1,n(ng)
625 DO i=istr,iend
626 cff=s1_v_obc(i,k,inorth,ir)* &
627 & s2_v_obc(i,k,inorth,ir)
628# ifdef MASKING
629 cff=cff*vmask(i,jend+1)
630# endif
631 my_dotprod(0)=my_dotprod(0)+cff
632 my_dotprod(isvvel)=my_dotprod(isvvel)+cff
633 END DO
634 END DO
635 END IF
636 END DO
637 END IF
638# endif
639!
640! Tracers.
641!
642 DO it=1,nt(ng)
643 DO k=1,n(ng)
644 DO j=jstrt,jendt
645 DO i=istrt,iendt
646 cff=s1_t(i,j,k,it)*s2_t(i,j,k,it)
647# ifdef MASKING
648 cff=cff*rmask(i,j)
649# endif
650 my_dotprod(0)=my_dotprod(0)+cff
651 my_dotprod(istvar(it))=my_dotprod(istvar(it))+cff
652 END DO
653 END DO
654 END DO
655 END DO
656
657# ifdef ADJUST_BOUNDARY
658!
659! Tracers open boundaries.
660!
661 DO it=1,nt(ng)
662 IF (any(lobc(:,istvar(it),ng))) THEN
663 DO ir=1,nbrec(ng)
664 IF ((lobc(iwest,istvar(it),ng)).and. &
665 & domain(ng)%Western_Edge(tile)) THEN
666 DO k=1,n(ng)
667 DO j=jstr,jend
668 cff=s1_t_obc(j,k,iwest,ir,it)* &
669 & s2_t_obc(j,k,iwest,ir,it)
670# ifdef MASKING
671 cff=cff*rmask(istr-1,j)
672# endif
673 my_dotprod(0)=my_dotprod(0)+cff
674 my_dotprod(istvar(it))=my_dotprod(istvar(it))+cff
675 END DO
676 END DO
677 END IF
678 IF ((lobc(ieast,istvar(it),ng)).and. &
679 & domain(ng)%Eastern_Edge(tile)) THEN
680 DO k=1,n(ng)
681 DO j=jstr,jend
682 cff=s1_t_obc(j,k,ieast,ir,it)* &
683 & s2_t_obc(j,k,ieast,ir,it)
684# ifdef MASKING
685 cff=cff*rmask(iend+1,j)
686# endif
687 my_dotprod(0)=my_dotprod(0)+cff
688 my_dotprod(istvar(it))=my_dotprod(istvar(it))+cff
689 END DO
690 END DO
691 END IF
692 IF ((lobc(isouth,istvar(it),ng)).and. &
693 & domain(ng)%Southern_Edge(tile)) THEN
694 DO k=1,n(ng)
695 DO i=istr,iend
696 cff=s1_t_obc(i,k,isouth,ir,it)* &
697 & s2_t_obc(i,k,isouth,ir,it)
698# ifdef MASKING
699 cff=cff*rmask(i,jstr-1)
700# endif
701 my_dotprod(0)=my_dotprod(0)+cff
702 my_dotprod(istvar(it))=my_dotprod(istvar(it))+cff
703 END DO
704 END DO
705 END IF
706 IF ((lobc(inorth,istvar(it),ng)).and. &
707 & domain(ng)%Northern_Edge(tile)) THEN
708 DO k=1,n(ng)
709 DO i=istr,iend
710 cff=s1_t_obc(i,k,inorth,ir,it)* &
711 & s2_t_obc(i,k,inorth,ir,it)
712# ifdef MASKING
713 cff=cff*rmask(i,jend+1)
714# endif
715 my_dotprod(0)=my_dotprod(0)+cff
716 my_dotprod(istvar(it))=my_dotprod(istvar(it))+cff
717 END DO
718 END DO
719 END IF
720 END DO
721 END IF
722 END DO
723# endif
724
725# ifdef ADJUST_STFLUX
726!
727! Surface tracers flux.
728!
729 DO it=1,nt(ng)
730 IF (lstflux(it,ng)) THEN
731 DO ir=1,nfrec(ng)
732 DO j=jstrt,jendt
733 DO i=istrt,iendt
734 cff=s1_tflux(i,j,ir,it)*s2_tflux(i,j,ir,it)
735# ifdef MASKING
736 cff=cff*rmask(i,j)
737# endif
738 my_dotprod(0)=my_dotprod(0)+cff
739 my_dotprod(istsur(it))=my_dotprod(istsur(it))+cff
740 END DO
741 END DO
742 END DO
743 END IF
744 END DO
745# endif
746
747#endif
748!
749!-----------------------------------------------------------------------
750! Perform parallel global reduction operations.
751!-----------------------------------------------------------------------
752!
753#ifdef DISTRIBUTE
754 nsub=1 ! distributed-memory
755#else
756 IF (domain(ng)%SouthWest_Corner(tile).and. &
757 & domain(ng)%NorthEast_Corner(tile)) THEN
758 nsub=1 ! non-tiled application
759 ELSE
760 nsub=ntilex(ng)*ntilee(ng) ! tiled application
761 END IF
762#endif
763!$OMP CRITICAL (DOT_PROD)
764 IF (tile_count.eq.0) THEN
765 DO i=0,nstatevars
766 dotprod(i)=0.0_r8
767 END DO
768 END IF
769 DO i=0,nstatevars
770 dotprod(i)=dotprod(i)+my_dotprod(i)
771 END DO
773 IF (tile_count.eq.nsub) THEN
774 tile_count=0
775#ifdef DISTRIBUTE
776 DO i=0,nstatevars
777 op_handle(i)='SUM'
778 END DO
779 CALL mp_reduce (ng, model, nstatevars+1, dotprod(0:), &
780 & op_handle(0:))
781#endif
782 END IF
783!$OMP END CRITICAL (DOT_PROD)
784!
785 RETURN
786 END SUBROUTINE state_dotprod
787
788 END MODULE state_dotprod_mod
integer isvvel
integer isvbar
integer isvstr
integer, dimension(:), allocatable istvar
integer isuvel
integer isfsur
integer isustr
integer isubar
integer, dimension(:), allocatable istsur
integer tile_count
integer, dimension(:), allocatable n
Definition mod_param.F:479
integer, dimension(:), allocatable ntilex
Definition mod_param.F:685
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer, dimension(:), allocatable ntilee
Definition mod_param.F:686
integer, dimension(:), allocatable nt
Definition mod_param.F:489
logical, dimension(:,:,:), allocatable lobc
integer, parameter iwest
logical, dimension(:,:), allocatable lstflux
integer, dimension(:), allocatable nfrec
integer, parameter isouth
integer, parameter ieast
integer, parameter inorth
integer, dimension(:), allocatable nbrec
subroutine, public state_dotprod(ng, tile, model, lbi, ubi, lbj, ubj, lbij, ubij, nstatevars, dotprod, rmask, umask, vmask, s1_t_obc, s2_t_obc, s1_u_obc, s2_u_obc, s1_v_obc, s2_v_obc, s1_ubar_obc, s2_ubar_obc, s1_vbar_obc, s2_vbar_obc, s1_zeta_obc, s2_zeta_obc, s1_sustr, s2_sustr, s1_svstr, s2_svstr, s1_tflux, s2_tflux, s1_t, s2_t, s1_u, s2_u, s1_v, s2_v, s1_zeta, s2_zeta)