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