ROMS
Loading...
Searching...
No Matches
back_cost.F
Go to the documentation of this file.
1#include "cppdefs.h"
2 MODULE back_cost_mod
3#ifdef BACKGROUND
4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group Andrew M. Moore !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This routine computes the background cost function (Jb) as the !
13! misfit (squared difference) between model and background state: !
14! !
15! Jb = 1/2 * transpose(X - Xb) * B^(-1) * (X - Xb) !
16! !
17! where !
18! !
19! Xb : background or reference state (first guess) !
20! X : nonlinear model initial state !
21! B : background-error covariance !
22! !
23! In incremental 4DVAR, the initial conditions estimate is done in !
24! in minimization space by defining the following transformation: !
25! !
26! deltaV = B^(-1/2) deltaX (v-space) !
27! or !
28! deltaX = B^(1/2) deltaV (x-space) !
29! !
30! where !
31! !
32! B = tanspose{B^(1/2)} B^(1/2) !
33! !
34! Then, the background cost function becomes: !
35! !
36! Jb = 1/2 * transpose{(X - Xb)} * B^(-1) * (X - Xb) x-space !
37! !
38! or !
39! !
40! Jb = 1/2 * transpose(deltaX) * B^(-1) * deltaX x-space !
41! !
42! or !
43! !
44! Jb = 1/2 * transpose(deltaV) * deltaV v-space !
45! !
46! Therefore, in v-space the background cost function gradient is: !
47! !
48! GRADv(Jb) = deltaV !
49! !
50! Notice that initially, Jb is zero since the model is initialized !
51! with the background state (X=Xb) and the minimization increment, !
52! deltaX is zero. !
53! !
54! !
55! WARNING: !
56! ------- !
57! !
58! The background cost function term is computed in v-space. Recall !
59! that in the inner loop the increment vector deltaX (x-space) and !
60! deltaV (v-space) are written into ITL(ng)%name NetCDF file at !
61! records 1 and 2, respectively. !
62! !
63!=======================================================================
64!
65 implicit none
66
67 CONTAINS
68!
69!***********************************************************************
70 SUBROUTINE back_cost (ng, tile, Lsum)
71!***********************************************************************
72!
73 USE mod_param
74# ifdef ADJUST_BOUNDARY
75 USE mod_boundary
76# endif
77# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
78 USE mod_forces
79# endif
80# ifdef MASKING
81 USE mod_grid
82# endif
83 USE mod_ocean
84!
85! Imported variable declarations.
86!
87 integer, intent(in) :: ng, tile, Lsum
88!
89! Local variable declarations.
90!
91# include "tile.h"
92!
93 CALL back_cost_tile (ng, tile, &
94 & lbi, ubi, lbj, ubj, lbij, ubij, &
95 & imins, imaxs, jmins, jmaxs, &
96 & lsum, &
97# ifdef MASKING
98 & grid(ng) % rmask, &
99 & grid(ng) % umask, &
100 & grid(ng) % vmask, &
101# endif
102# ifdef ADJUST_BOUNDARY
103# ifdef SOLVE3D
104 & boundary(ng) % tl_t_obc, &
105 & boundary(ng) % tl_u_obc, &
106 & boundary(ng) % tl_v_obc, &
107# endif
108 & boundary(ng) % tl_ubar_obc, &
109 & boundary(ng) % tl_vbar_obc, &
110 & boundary(ng) % tl_zeta_obc, &
111# endif
112# ifdef ADJUST_WSTRESS
113 & forces(ng) % tl_ustr, &
114 & forces(ng) % tl_vstr, &
115# endif
116# ifdef SOLVE3D
117# ifdef ADJUST_STFLUX
118 & forces(ng) % tl_tflux, &
119# endif
120 & ocean(ng) % tl_t, &
121 & ocean(ng) % tl_u, &
122 & ocean(ng) % tl_v, &
123# else
124 & ocean(ng) % tl_ubar, &
125 & ocean(ng) % tl_vbar, &
126# endif
127 & ocean(ng) % tl_zeta)
128
129 RETURN
130 END SUBROUTINE back_cost
131!
132!***********************************************************************
133 SUBROUTINE back_cost_tile (ng, tile, &
134 & LBi, UBi, LBj, UBj, LBij, UBij, &
135 & IminS, ImaxS, JminS, JmaxS, &
136 & Lsum, &
137# ifdef MASKING
138 & rmask, umask, vmask, &
139# endif
140# ifdef ADJUST_BOUNDARY
141# ifdef SOLVE3D
142 & tl_t_obc, tl_u_obc, tl_v_obc, &
143# endif
144 & tl_ubar_obc, tl_vbar_obc, &
145 & tl_zeta_obc, &
146# endif
147# ifdef ADJUST_WSTRESS
148 & tl_ustr, &
149 & tl_vstr, &
150# endif
151# ifdef SOLVE3D
152# ifdef ADJUST_STFLUX
153 & tl_tflux, &
154# endif
155 & tl_t, &
156 & tl_u, &
157 & tl_v, &
158# else
159 & tl_ubar, &
160 & tl_vbar, &
161# endif
162 & tl_zeta)
163!***********************************************************************
164!
165 USE mod_param
166 USE mod_parallel
167 USE mod_fourdvar
168 USE mod_ncparam
169 USE mod_scalars
170
171# ifdef DISTRIBUTE
172!
173 USE distribute_mod, ONLY : mp_reduce
174# endif
175!
176! Imported variable declarations.
177!
178 integer, intent(in) :: ng, tile
179 integer, intent(in) :: LBi, UBi, LBj, UBj, LBij, UBij
180 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
181 integer, intent(in) :: Lsum
182!
183# ifdef ASSUMED_SHAPE
184# ifdef MASKING
185 real(r8), intent(in) :: rmask(LBi:,LBj:)
186 real(r8), intent(in) :: umask(LBi:,LBj:)
187 real(r8), intent(in) :: vmask(LBi:,LBj:)
188# endif
189# ifdef ADJUST_BOUNDARY
190# ifdef SOLVE3D
191 real(r8), intent(in) :: tl_t_obc(LBij:,:,:,:,:,:)
192 real(r8), intent(in) :: tl_u_obc(LBij:,:,:,:,:)
193 real(r8), intent(in) :: tl_v_obc(LBij:,:,:,:,:)
194# endif
195 real(r8), intent(in) :: tl_ubar_obc(LBij:,:,:,:)
196 real(r8), intent(in) :: tl_vbar_obc(LBij:,:,:,:)
197 real(r8), intent(in) :: tl_zeta_obc(LBij:,:,:,:)
198# endif
199# ifdef ADJUST_WSTRESS
200 real(r8), intent(in) :: tl_ustr(LBi:,LBj:,:,:)
201 real(r8), intent(in) :: tl_vstr(LBi:,LBj:,:,:)
202# endif
203# ifdef SOLVE3D
204# ifdef ADJUST_STFLUX
205 real(r8), intent(in) :: tl_tflux(LBi:,LBj:,:,:,:)
206# endif
207 real(r8), intent(in) :: tl_t(LBi:,LBj:,:,:,:)
208 real(r8), intent(in) :: tl_u(LBi:,LBj:,:,:)
209 real(r8), intent(in) :: tl_v(LBi:,LBj:,:,:)
210# else
211 real(r8), intent(in) :: tl_ubar(LBi:,LBj:,:)
212 real(r8), intent(in) :: tl_vbar(LBi:,LBj:,:)
213# endif
214 real(r8), intent(in) :: tl_zeta(LBi:,LBj:,:)
215
216# else
217
218# ifdef MASKING
219 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
220 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
221 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
222# endif
223# ifdef ADJUST_BOUNDARY
224# ifdef SOLVE3D
225 real(r8), intent(in) :: tl_t_obc(LBij:UBij,N(ng),4, &
226 & Nbrec(ng),2,NT(ng))
227 real(r8), intent(in) :: tl_u_obc(LBij:UBij,N(ng),4,Nbrec(ng),2)
228 real(r8), intent(in) :: tl_v_obc(LBij:UBij,N(ng),4,Nbrec(ng),2)
229# endif
230 real(r8), intent(in) :: tl_ubar_obc(LBij:UBij,4,Nbrec(ng),2)
231 real(r8), intent(in) :: tl_vbar_obc(LBij:UBij,4,Nbrec(ng),2)
232 real(r8), intent(in) :: tl_zeta_obc(LBij:UBij,4,Nbrec(ng),2)
233# endif
234# ifdef ADJUST_WSTRESS
235 real(r8), intent(in) :: tl_ustr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
236 real(r8), intent(in) :: tl_vstr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
237# endif
238# ifdef SOLVE3D
239# ifdef ADJUST_STFLUX
240 real(r8), intent(in) :: tl_tflux(LBi:UBi,LBj:UBj, &
241 & Nfrec(ng),2,NT(ng))
242# endif
243 real(r8), intent(in) :: tl_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
244 real(r8), intent(in) :: tl_u(LBi:UBi,LBj:UBj,N(ng),2)
245 real(r8), intent(in) :: tl_v(LBi:UBi,LBj:UBj,N(ng),2)
246# else
247 real(r8), intent(in) :: tl_ubar(LBi:UBi,LBj:UBj,:)
248 real(r8), intent(in) :: tl_vbar(LBi:UBi,LBj:UBj,:)
249# endif
250 real(r8), intent(in) :: tl_zeta(LBi:UBi,LBj:UBj,:)
251# endif
252!
253! Local variable declarations.
254!
255 integer :: NSUB, i, ib, ir, j, k
256# ifdef SOLVE3D
257 integer :: itrc
258# endif
259
260 real(r8), dimension(0:NstateVar(ng)) :: my_BackCost
261
262 real(r8) :: cff1, cff2
263
264# ifdef DISTRIBUTE
265 character (len=3), dimension(0:NstateVar(ng)) :: op_handle
266# endif
267
268# include "set_bounds.h"
269!
270!----------------------------------------------------------------------
271! Compute v-space background cost function (Jb) as misfit between
272! model and background states at initial time of the assimilation
273! window using the sum of all the previous outer-loop increments
274! (Lsum).
275!
276! Initially, the misfit innovation matrix (X-Xb) is zero. As the
277! assimilation algorithm iterates, Jb becomes greater than zero.
278!----------------------------------------------------------------------
279!
280 DO i=0,nstatevar(ng)
281 my_backcost(i)=0.0_r8
282 END DO
283!
284! Free-surface contribution.
285!
286 DO j=jstr,jend
287 DO i=istr,iend
288 cff1=tl_zeta(i,j,lsum)
289# ifdef MASKING
290 cff1=cff1*rmask(i,j)
291# endif
292 cff2=0.5_r8*cff1*cff1
293 my_backcost(0)=my_backcost(0)+cff2
294 my_backcost(isfsur)=my_backcost(isfsur)+cff2
295 END DO
296 END DO
297
298# ifdef ADJUST_BOUNDARY
299!
300! Free-surface open boundaries.
301!
302 IF (any(lobc(:,isfsur,ng))) THEN
303 DO ir=1,nbrec(ng)
304 IF ((lobc(iwest,isfsur,ng)).and. &
305 & domain(ng)%Western_Edge(tile)) THEN
306 ib=iwest
307 DO j=jstr,jend
308 cff1=tl_zeta_obc(j,ib,ir,lsum)
309# ifdef MASKING
310 cff1=cff1*rmask(istr-1,j)
311# endif
312 cff2=0.5_r8*cff1*cff1
313 my_backcost(0)=my_backcost(0)+cff2
314 my_backcost(isfsur)=my_backcost(isfsur)+cff2
315 END DO
316 END IF
317 IF ((lobc(ieast,isfsur,ng)).and. &
318 & domain(ng)%Eastern_Edge(tile)) THEN
319 ib=ieast
320 DO j=jstr,jend
321 cff1=tl_zeta_obc(j,ib,ir,lsum)
322# ifdef MASKING
323 cff1=cff1*rmask(iend+1,j)
324# endif
325 cff2=0.5_r8*cff1*cff1
326 my_backcost(0)=my_backcost(0)+cff2
327 my_backcost(isfsur)=my_backcost(isfsur)+cff2
328 END DO
329 END IF
330 IF ((lobc(isouth,isfsur,ng)).and. &
331 & domain(ng)%Southern_Edge(tile)) THEN
332 ib=isouth
333 DO i=istr,iend
334 cff1=tl_zeta_obc(i,ib,ir,lsum)
335# ifdef MASKING
336 cff1=cff1*rmask(i,jstr-1)
337# endif
338 cff2=0.5_r8*cff1*cff1
339 my_backcost(0)=my_backcost(0)+cff2
340 my_backcost(isfsur)=my_backcost(isfsur)+cff2
341 END DO
342 END IF
343 IF ((lobc(inorth,isfsur,ng)).and. &
344 & domain(ng)%Northern_Edge(tile)) THEN
345 ib=inorth
346 DO i=istr,iend
347 cff1=tl_zeta_obc(i,ib,ir,lsum)
348# ifdef MASKING
349 cff1=cff1*rmask(i,jend+1)
350# endif
351 cff2=0.5_r8*cff1*cff1
352 my_backcost(0)=my_backcost(0)+cff2
353 my_backcost(isfsur)=my_backcost(isfsur)+cff2
354 END DO
355 END IF
356 END DO
357 END IF
358# endif
359
360# ifndef SOLVE3D
361!
362! 2D U-momentum contribution.
363!
364 DO j=jstr,jend
365 DO i=istru,iend
366 cff1=tl_ubar(i,j,lsum)
367# ifdef MASKING
368 cff1=cff1*umask(i,j)
369# endif
370 cff2=0.5_r8*cff1*cff1
371 my_backcost(0)=my_backcost(0)+cff2
372 my_backcost(isubar)=my_backcost(isubar)+cff2
373 END DO
374 END DO
375# endif
376
377# ifdef ADJUST_BOUNDARY
378!
379! 2D U-momentum open boundaries.
380!
381 IF (any(lobc(:,isubar,ng))) THEN
382 DO ir=1,nbrec(ng)
383 IF ((lobc(iwest,isubar,ng)).and. &
384 & domain(ng)%Western_Edge(tile)) THEN
385 ib=iwest
386 DO j=jstr,jend
387 cff1=tl_ubar_obc(j,ib,ir,lsum)
388# ifdef MASKING
389 cff1=cff1*umask(istr,j)
390# endif
391 cff2=0.5_r8*cff1*cff1
392 my_backcost(0)=my_backcost(0)+cff2
393 my_backcost(isubar)=my_backcost(isubar)+cff2
394 END DO
395 END IF
396 IF ((lobc(ieast,isubar,ng)).and. &
397 & domain(ng)%Eastern_Edge(tile)) THEN
398 ib=ieast
399 DO j=jstr,jend
400 cff1=tl_ubar_obc(j,ib,ir,lsum)
401# ifdef MASKING
402 cff1=cff1*umask(iend+1,j)
403# endif
404 cff2=0.5_r8*cff1*cff1
405 my_backcost(0)=my_backcost(0)+cff2
406 my_backcost(isubar)=my_backcost(isubar)+cff2
407 END DO
408 END IF
409 IF ((lobc(isouth,isubar,ng)).and. &
410 & domain(ng)%Southern_Edge(tile)) THEN
411 ib=isouth
412 DO i=istru,iend
413 cff1=tl_ubar_obc(i,ib,ir,lsum)
414# ifdef MASKING
415 cff1=cff1*umask(i,jstr-1)
416# endif
417 cff2=0.5_r8*cff1*cff1
418 my_backcost(0)=my_backcost(0)+cff2
419 my_backcost(isubar)=my_backcost(isubar)+cff2
420 END DO
421 END IF
422 IF ((lobc(inorth,isubar,ng)).and. &
423 & domain(ng)%Northern_Edge(tile)) THEN
424 ib=inorth
425 DO i=istru,iend
426 cff1=tl_ubar_obc(i,ib,ir,lsum)
427# ifdef MASKING
428 cff1=cff1*umask(i,jend+1)
429# endif
430 cff2=0.5_r8*cff1*cff1
431 my_backcost(0)=my_backcost(0)+cff2
432 my_backcost(isubar)=my_backcost(isubar)+cff2
433 END DO
434 END IF
435 END DO
436 END IF
437# endif
438
439# ifndef SOLVE3D
440!
441! 2D V-momentum contribution.
442!
443 DO j=jstrv,jend
444 DO i=istr,iend
445 cff1=tl_vbar(i,j,lsum)
446# ifdef MASKING
447 cff1=cff1*vmask(i,j)
448# endif
449 cff2=0.5_r8*cff1*cff1
450 my_backcost(0)=my_backcost(0)+cff2
451 my_backcost(isvbar)=my_backcost(isvbar)+cff2
452 END DO
453 END DO
454# endif
455
456# ifdef ADJUST_BOUNDARY
457!
458! 2D V-momentum open boundaries.
459!
460 IF (any(lobc(:,isvbar,ng))) THEN
461 DO ir=1,nbrec(ng)
462 IF ((lobc(iwest,isvbar,ng)).and. &
463 & domain(ng)%Western_Edge(tile)) THEN
464 ib=iwest
465 DO j=jstrv,jend
466 cff1=tl_vbar_obc(j,ib,ir,lsum)
467# ifdef MASKING
468 cff1=cff1*vmask(istr-1,j)
469# endif
470 cff2=0.5_r8*cff1*cff1
471 my_backcost(0)=my_backcost(0)+cff2
472 my_backcost(isvbar)=my_backcost(isvbar)+cff2
473 END DO
474 END IF
475 IF ((lobc(ieast,isvbar,ng)).and. &
476 & domain(ng)%Eastern_Edge(tile)) THEN
477 ib=ieast
478 DO j=jstrv,jend
479 cff1=tl_vbar_obc(j,ib,ir,lsum)
480# ifdef MASKING
481 cff1=cff1*vmask(iend+1,j)
482# endif
483 cff2=0.5_r8*cff1*cff1
484 my_backcost(0)=my_backcost(0)+cff2
485 my_backcost(isvbar)=my_backcost(isvbar)+cff2
486 END DO
487 END IF
488 IF ((lobc(isouth,isvbar,ng)).and. &
489 & domain(ng)%Southern_Edge(tile)) THEN
490 ib=isouth
491 DO i=istr,iend
492 cff1=tl_vbar_obc(i,ib,ir,lsum)
493# ifdef MASKING
494 cff1=cff1*vmask(i,jstr)
495# endif
496 cff2=0.5_r8*cff1*cff1
497 my_backcost(0)=my_backcost(0)+cff2
498 my_backcost(isvbar)=my_backcost(isvbar)+cff2
499 END DO
500 END IF
501 IF ((lobc(inorth,isvbar,ng)).and. &
502 & domain(ng)%Northern_Edge(tile)) THEN
503 ib=inorth
504 DO i=istr,iend
505 cff1=tl_vbar_obc(i,ib,ir,lsum)
506# ifdef MASKING
507 cff1=cff1*vmask(i,jend+1)
508# endif
509 cff2=0.5_r8*cff1*cff1
510 my_backcost(0)=my_backcost(0)+cff2
511 my_backcost(isvbar)=my_backcost(isvbar)+cff2
512 END DO
513 END IF
514 END DO
515 END IF
516# endif
517
518# ifdef ADJUST_WSTRESS
519!
520! Surface momentum stress contribution.
521!
522 DO ir=1,nfrec(ng)
523 DO j=jstrr,jendr
524 DO i=istr,iendr
525 cff1=tl_ustr(i,j,ir,lsum)
526# ifdef MASKING
527 cff1=cff1*umask(i,j)
528# endif
529 cff2=0.5_r8*cff1*cff1
530 my_backcost(0)=my_backcost(0)+cff2
531 my_backcost(isustr)=my_backcost(isustr)+cff2
532 END DO
533 END DO
534 DO j=jstr,jendr
535 DO i=istrr,iendr
536 cff1=tl_vstr(i,j,ir,lsum)
537# ifdef MASKING
538 cff1=cff1*vmask(i,j)
539# endif
540 cff2=0.5_r8*cff1*cff1
541 my_backcost(0)=my_backcost(0)+cff2
542 my_backcost(isvstr)=my_backcost(isvstr)+cff2
543 END DO
544 END DO
545 END DO
546# endif
547
548# ifdef SOLVE3D
549!
550! 3D U-momentum contribution.
551!
552 DO k=1,n(ng)
553 DO j=jstr,jend
554 DO i=istru,iend
555 cff1=tl_u(i,j,k,lsum)
556# ifdef MASKING
557 cff1=cff1*umask(i,j)
558# endif
559 cff2=0.5_r8*cff1*cff1
560 my_backcost(0)=my_backcost(0)+cff2
561 my_backcost(isuvel)=my_backcost(isuvel)+cff2
562 END DO
563 END DO
564 END DO
565
566# ifdef ADJUST_BOUNDARY
567!
568! 3D U-momentum open boundaries.
569!
570 IF (any(lobc(:,isuvel,ng))) THEN
571 DO ir=1,nbrec(ng)
572 IF ((lobc(iwest,isuvel,ng)).and. &
573 & domain(ng)%Western_Edge(tile)) THEN
574 ib=iwest
575 DO k=1,n(ng)
576 DO j=jstr,jend
577 cff1=tl_u_obc(j,k,ib,ir,lsum)
578# ifdef MASKING
579 cff1=cff1*umask(istr,j)
580# endif
581 cff2=0.5_r8*cff1*cff1
582 my_backcost(0)=my_backcost(0)+cff2
583 my_backcost(isuvel)=my_backcost(isuvel)+cff2
584 END DO
585 END DO
586 END IF
587 IF ((lobc(ieast,isuvel,ng)).and. &
588 & domain(ng)%Eastern_Edge(tile)) THEN
589 ib=ieast
590 DO k=1,n(ng)
591 DO j=jstr,jend
592 cff1=tl_u_obc(j,k,ib,ir,lsum)
593# ifdef MASKING
594 cff1=cff1*umask(iend+1,j)
595# endif
596 cff2=0.5_r8*cff1*cff1
597 my_backcost(0)=my_backcost(0)+cff2
598 my_backcost(isuvel)=my_backcost(isuvel)+cff2
599 END DO
600 END DO
601 END IF
602 IF ((lobc(isouth,isuvel,ng)).and. &
603 & domain(ng)%Southern_Edge(tile)) THEN
604 ib=isouth
605 DO k=1,n(ng)
606 DO i=istru,iend
607 cff1=tl_u_obc(i,k,ib,ir,lsum)
608# ifdef MASKING
609 cff1=cff1*umask(i,jstr-1)
610# endif
611 cff2=0.5_r8*cff1*cff1
612 my_backcost(0)=my_backcost(0)+cff2
613 my_backcost(isuvel)=my_backcost(isuvel)+cff2
614 END DO
615 END DO
616 END IF
617 IF ((lobc(inorth,isuvel,ng)).and. &
618 & domain(ng)%Northern_Edge(tile)) THEN
619 ib=inorth
620 DO k=1,n(ng)
621 DO i=istru,iend
622 cff1=tl_u_obc(i,k,ib,ir,lsum)
623# ifdef MASKING
624 cff1=cff1*umask(i,jend+1)
625# endif
626 cff2=0.5_r8*cff1*cff1
627 my_backcost(0)=my_backcost(0)+cff2
628 my_backcost(isuvel)=my_backcost(isuvel)+cff2
629 END DO
630 END DO
631 END IF
632 END DO
633 END IF
634# endif
635!
636! 3D V-momentum contribution.
637!
638 DO k=1,n(ng)
639 DO j=jstrv,jend
640 DO i=istr,iend
641 cff1=tl_v(i,j,k,lsum)
642# ifdef MASKING
643 cff1=cff1*vmask(i,j)
644# endif
645 cff2=0.5_r8*cff1*cff1
646 my_backcost(0)=my_backcost(0)+cff2
647 my_backcost(isvvel)=my_backcost(isvvel)+cff2
648 END DO
649 END DO
650 END DO
651
652# ifdef ADJUST_BOUNDARY
653!
654! 3D V-momentum open boundaries.
655!
656 IF (any(lobc(:,isvvel,ng))) THEN
657 DO ir=1,nbrec(ng)
658 IF ((lobc(iwest,isvvel,ng)).and. &
659 & domain(ng)%Western_Edge(tile)) THEN
660 ib=iwest
661 DO k=1,n(ng)
662 DO j=jstrv,jend
663 cff1=tl_v_obc(j,k,ib,ir,lsum)
664# ifdef MASKING
665 cff1=cff1*vmask(istr-1,j)
666# endif
667 cff2=0.5_r8*cff1*cff1
668 my_backcost(0)=my_backcost(0)+cff2
669 my_backcost(isvvel)=my_backcost(isvvel)+cff2
670 END DO
671 END DO
672 END IF
673 IF ((lobc(ieast,isvvel,ng)).and. &
674 & domain(ng)%Eastern_Edge(tile)) THEN
675 ib=ieast
676 DO k=1,n(ng)
677 DO j=jstrv,jend
678 cff1=tl_v_obc(j,k,ib,ir,lsum)
679# ifdef MASKING
680 cff1=cff1*vmask(iend+1,j)
681# endif
682 cff2=0.5_r8*cff1*cff1
683 my_backcost(0)=my_backcost(0)+cff2
684 my_backcost(isvvel)=my_backcost(isvvel)+cff2
685 END DO
686 END DO
687 END IF
688 IF ((lobc(isouth,isvvel,ng)).and. &
689 & domain(ng)%Southern_Edge(tile)) THEN
690 ib=isouth
691 DO k=1,n(ng)
692 DO i=istr,iend
693 cff1=tl_v_obc(i,k,ib,ir,lsum)
694# ifdef MASKING
695 cff1=cff1*vmask(i,jstr)
696# endif
697 cff2=0.5_r8*cff1*cff1
698 my_backcost(0)=my_backcost(0)+cff2
699 my_backcost(isvvel)=my_backcost(isvvel)+cff2
700 END DO
701 END DO
702 END IF
703 IF ((lobc(inorth,isvvel,ng)).and. &
704 & domain(ng)%Northern_Edge(tile)) THEN
705 ib=inorth
706 DO k=1,n(ng)
707 DO i=istr,iend
708 cff1=tl_v_obc(i,k,ib,ir,lsum)
709# ifdef MASKING
710 cff1=cff1*vmask(i,jend+1)
711# endif
712 cff2=0.5_r8*cff1*cff1
713 my_backcost(0)=my_backcost(0)+cff2
714 my_backcost(isvvel)=my_backcost(isvvel)+cff2
715 END DO
716 END DO
717 END IF
718 END DO
719 END IF
720# endif
721!
722! Tracers contribution.
723!
724 DO itrc=1,nt(ng)
725 DO k=1,n(ng)
726 DO j=jstr,jend
727 DO i=istr,iend
728 cff1=tl_t(i,j,k,lsum,itrc)
729# ifdef MASKING
730 cff1=cff1*rmask(i,j)
731# endif
732 cff2=0.5_r8*cff1*cff1
733 my_backcost(0)=my_backcost(0)+cff2
734 my_backcost(istvar(itrc))=my_backcost(istvar(itrc))+cff2
735 END DO
736 END DO
737 END DO
738 END DO
739
740# ifdef ADJUST_BOUNDARY
741!
742! Tracers open boundaries.
743!
744 DO itrc=1,nt(ng)
745 IF (any(lobc(:,istvar(itrc),ng))) THEN
746 DO ir=1,nbrec(ng)
747 IF ((lobc(iwest,istvar(itrc),ng)).and. &
748 & domain(ng)%Western_Edge(tile)) THEN
749 ib=iwest
750 DO k=1,n(ng)
751 DO j=jstr,jend
752 cff1=tl_t_obc(j,k,ib,ir,lsum,itrc)
753# ifdef MASKING
754 cff1=cff1*rmask(istr-1,j)
755# endif
756 cff2=0.5_r8*cff1*cff1
757 my_backcost(0)=my_backcost(0)+cff2
758 my_backcost(istvar(itrc))=my_backcost(istvar(itrc))+ &
759 & cff2
760 END DO
761 END DO
762 END IF
763 IF ((lobc(ieast,istvar(itrc),ng)).and. &
764 & domain(ng)%Eastern_Edge(tile)) THEN
765 ib=ieast
766 DO k=1,n(ng)
767 DO j=jstr,jend
768 cff1=tl_t_obc(j,k,ib,ir,lsum,itrc)
769# ifdef MASKING
770 cff1=cff1*rmask(iend+1,j)
771# endif
772 cff2=0.5_r8*cff1*cff1
773 my_backcost(0)=my_backcost(0)+cff2
774 my_backcost(istvar(itrc))=my_backcost(istvar(itrc))+ &
775 & cff2
776 END DO
777 END DO
778 END IF
779 IF ((lobc(isouth,istvar(itrc),ng)).and. &
780 & domain(ng)%Southern_Edge(tile)) THEN
781 ib=isouth
782 DO k=1,n(ng)
783 DO i=istr,iend
784 cff1=tl_t_obc(i,k,ib,ir,lsum,itrc)
785# ifdef MASKING
786 cff1=cff1*rmask(i,jstr-1)
787# endif
788 cff2=0.5_r8*cff1*cff1
789 my_backcost(0)=my_backcost(0)+cff2
790 my_backcost(istvar(itrc))=my_backcost(istvar(itrc))+ &
791 & cff2
792 END DO
793 END DO
794 END IF
795 IF ((lobc(inorth,istvar(itrc),ng)).and. &
796 & domain(ng)%Northern_Edge(tile)) THEN
797 ib=inorth
798 DO k=1,n(ng)
799 DO i=istr,iend
800 cff1=tl_t_obc(i,k,ib,ir,lsum,itrc)
801# ifdef MASKING
802 cff1=cff1*rmask(i,jend+1)
803# endif
804 cff2=0.5_r8*cff1*cff1
805 my_backcost(0)=my_backcost(0)+cff2
806 my_backcost(istvar(itrc))=my_backcost(istvar(itrc))+ &
807 & cff2
808 END DO
809 END DO
810 END IF
811 END DO
812 END IF
813 END DO
814# endif
815
816# ifdef ADJUST_STFLUX
817!
818! Surface tracers flux contribution.
819!
820 DO itrc=1,nt(ng)
821 IF (lstflux(itrc,ng)) THEN
822 DO ir=1,nfrec(ng)
823 DO j=jstrr,jendr
824 DO i=istrr,iendr
825 cff1=tl_tflux(i,j,ir,lsum,itrc)
826# ifdef MASKING
827 cff1=cff1*rmask(i,j)
828# endif
829 cff2=0.5_r8*cff1*cff1
830 my_backcost(0)=my_backcost(0)+cff2
831 my_backcost(istsur(itrc))=my_backcost(istsur(itrc))+cff2
832 END DO
833 END DO
834 END DO
835 END IF
836 END DO
837# endif
838# endif
839!
840!-----------------------------------------------------------------------
841! Compute global background cost function.
842!-----------------------------------------------------------------------
843!
844# ifdef DISTRIBUTE
845 nsub=1 ! distributed-memory
846# else
847 IF (domain(ng)%SouthWest_Corner(tile).and. &
848 & domain(ng)%NorthEast_Corner(tile)) THEN
849 nsub=1 ! non-tiled application
850 ELSE
851 nsub=ntilex(ng)*ntilee(ng) ! tiled application
852 END IF
853# endif
854!$OMP CRITICAL (BACKCOST)
855 IF (tile_count.eq.0) THEN
856 DO i=0,nstatevar(ng)
857 fourdvar(ng)%BackCost(i)=my_backcost(i)
858 END DO
859 ELSE
860 DO i=0,nstatevar(ng)
861 fourdvar(ng)%BackCost(i)=fourdvar(ng)%BackCost(i)+ &
862 & my_backcost(i)
863 END DO
864 END IF
865 tile_count=tile_count+1
866 IF (tile_count.eq.nsub) THEN
867 tile_count=0
868# ifdef DISTRIBUTE
869 DO i=0,nstatevar(ng)
870 op_handle(i)='SUM'
871 END DO
872 CALL mp_reduce (ng, inlm, nstatevar(ng)+1, &
873 & fourdvar(ng)%BackCost(0:), op_handle(0:))
874# endif
875 END IF
876!$OMP END CRITICAL (BACKCOST)
877
878 RETURN
879 END SUBROUTINE back_cost_tile
880#endif
881 END MODULE back_cost_mod
distribute_mod::mp_reduce
Definition distribute.F:113
back_cost_mod
Definition back_cost.F:2
back_cost_mod::back_cost
subroutine back_cost(ng, tile, lsum)
Definition back_cost.F:71
back_cost_mod::back_cost_tile
subroutine back_cost_tile(ng, tile, lbi, ubi, lbj, ubj, lbij, ubij, imins, imaxs, jmins, jmaxs, lsum, rmask, umask, vmask, tl_t_obc, tl_u_obc, tl_v_obc, tl_ubar_obc, tl_vbar_obc, tl_zeta_obc, tl_ustr, tl_vstr, tl_tflux, tl_t, tl_u, tl_v, tl_zeta)
Definition back_cost.F:163
distribute_mod
Definition distribute.F:3
mod_boundary
Definition mod_boundary.F:2
mod_boundary::boundary
type(t_boundary), dimension(:), allocatable boundary
Definition mod_boundary.F:746
mod_forces
Definition mod_forces.F:2
mod_forces::forces
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
mod_fourdvar
Definition mod_fourdvar.F:2
mod_fourdvar::fourdvar
type(t_fourdvar), dimension(:), allocatable fourdvar
Definition mod_fourdvar.F:211
mod_fourdvar::nstatevar
integer, dimension(:), allocatable nstatevar
Definition mod_fourdvar.F:362
mod_grid
Definition mod_grid.F:2
mod_grid::grid
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
mod_ncparam
Definition mod_ncparam.F:2
mod_ncparam::isvvel
integer isvvel
Definition mod_ncparam.F:507
mod_ncparam::isvbar
integer isvbar
Definition mod_ncparam.F:505
mod_ncparam::isvstr
integer isvstr
Definition mod_ncparam.F:511
mod_ncparam::istvar
integer, dimension(:), allocatable istvar
Definition mod_ncparam.F:520
mod_ncparam::isuvel
integer isuvel
Definition mod_ncparam.F:506
mod_ncparam::isfsur
integer isfsur
Definition mod_ncparam.F:503
mod_ncparam::isustr
integer isustr
Definition mod_ncparam.F:510
mod_ncparam::isubar
integer isubar
Definition mod_ncparam.F:504
mod_ncparam::istsur
integer, dimension(:), allocatable istsur
Definition mod_ncparam.F:519
mod_ocean
Definition mod_ocean.F:2
mod_ocean::ocean
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
mod_parallel
Definition mod_parallel.F:2
mod_parallel::tile_count
integer tile_count
Definition mod_parallel.F:89
mod_param
Definition mod_param.F:2
mod_param::inlm
integer, parameter inlm
Definition mod_param.F:662
mod_param::ntilex
integer, dimension(:), allocatable ntilex
Definition mod_param.F:685
mod_param::domain
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
mod_param::ntilee
integer, dimension(:), allocatable ntilee
Definition mod_param.F:686
mod_scalars
Definition mod_scalars.F:2
mod_scalars::lobc
logical, dimension(:,:,:), allocatable lobc
Definition mod_scalars.F:1127
mod_scalars::iwest
integer, parameter iwest
Definition mod_scalars.F:1432
mod_scalars::lstflux
logical, dimension(:,:), allocatable lstflux
Definition mod_scalars.F:1156
mod_scalars::isouth
integer, parameter isouth
Definition mod_scalars.F:1433
mod_scalars::ieast
integer, parameter ieast
Definition mod_scalars.F:1434
mod_scalars::inorth
integer, parameter inorth
Definition mod_scalars.F:1435