ROMS
Loading...
Searching...
No Matches
tl_convolution.F
Go to the documentation of this file.
1#include "cppdefs.h"
3
4#if defined TANGENT && defined FOUR_DVAR
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 performs a spatial convolution of the tangeny state !
14! solution to model the background error correlations, C, using !
15! a generalized diffusion operator. This allows the observational !
16! information to spread spatially in 4DVAR data assimilation. !
17! !
18! The background error covariance is defined as: !
19! !
20! B = S C S !
21! !
22! C = C^(1/2) C^(T/2) !
23! !
24! C^(1/2) = G L^(1/2) W^(-1/2) TLM !
25! C^(T/2) = W^(-1/2) L^(T/2) G ADM !
26! !
27! where !
28! !
29! B : background-error covariance matrix !
30! S : diagonal matrix of background-error standard deviations !
31! C : symmetric matrix of background-error correlations !
32! G : normalization coefficients matrix used to ensure that the !
33! diagonal variances of C are equal to unity. !
34! L : tangent linear and adjoint diffusion operators !
35! W : diagonal matrix of local area or volume metrics used to !
36! convert L into a symmetric matrix: LW^(-1). !
37! !
38! Here, T/2 denote the transpose if a squared-root factor. !
39! !
40! This routine is used to provide a better preconditioning of the !
41! minimization problem, which is expressed as a function of a new !
42! state vector, v, given by: !
43! !
44! v = B^(-1/2) delta_x (v-space) !
45! or !
46! delta_x = B^(1/2) v !
47! !
48! where !
49! !
50! B = B^(T/2) B^(1/2) !
51! !
52! Therefore, the cost function, J, gradient becomes: !
53! !
54! GRAD_v(J) = v + {B^(T/2)} GRAD_x(J) !
55! !
56! In incremental 4DVAR, these spatial convolutions constitutes a !
57! smoothing action on the correlation operator and they are used !
58! to transform between model space to minimization space and vice !
59! versa: !
60! !
61! ad_convolution compute GRAD_v(J) from GRAD_x(J) !
62! tl_convolution compute x from v !
63! !
64! The minimization of of J in the descent algorithm is in v-space. !
65! !
66! Reference: !
67! !
68! Weaver, A. and P. Courtier, 2001: Correlation modeling on the !
69! sphere using a generalized diffusion equation, Q.J.R. Meteo. !
70! Soc, 127, 1815-1846. !
71! !
72!======================================================================!
73!
74 USE mod_kinds
75!
76 implicit none
77!
78 PRIVATE
79 PUBLIC :: tl_convolution
80!
81 CONTAINS
82!
83!***********************************************************************
84 SUBROUTINE tl_convolution (ng, tile, Linp, Lweak, ifac)
85!***********************************************************************
86!
87 USE mod_param
88# ifdef ADJUST_BOUNDARY
89 USE mod_boundary
90# endif
91# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
92 USE mod_forces
93# endif
94 USE mod_grid
95 USE mod_mixing
96 USE mod_ocean
97# if defined SEDIMENT && defined SED_MORPH && defined SOLVE3D
98 USE mod_sedbed
99# endif
100 USE mod_stepping
101!
102! Imported variable declarations.
103!
104 logical, intent(in) :: lweak
105
106 integer, intent(in) :: ng, tile, linp, ifac
107!
108! Local variable declarations.
109!
110# include "tile.h"
111!
112 CALL tl_convolution_tile (ng, tile, &
113 & lbi, ubi, lbj, ubj, lbij, ubij, &
114 & imins, imaxs, jmins, jmaxs, &
115 & nstp(ng), nnew(ng), linp, lweak, ifac, &
116 & grid(ng) % pm, &
117 & grid(ng) % om_p, &
118 & grid(ng) % om_r, &
119 & grid(ng) % om_u, &
120 & grid(ng) % om_v, &
121 & grid(ng) % pn, &
122 & grid(ng) % on_p, &
123 & grid(ng) % on_r, &
124 & grid(ng) % on_u, &
125 & grid(ng) % on_v, &
126 & grid(ng) % pmon_p, &
127 & grid(ng) % pmon_r, &
128 & grid(ng) % pmon_u, &
129 & grid(ng) % pnom_p, &
130 & grid(ng) % pnom_r, &
131 & grid(ng) % pnom_v, &
132# ifdef MASKING
133 & grid(ng) % rmask, &
134 & grid(ng) % pmask, &
135 & grid(ng) % umask, &
136 & grid(ng) % vmask, &
137# endif
138# ifdef SOLVE3D
139 & grid(ng) % h, &
140# ifdef ICESHELF
141 & grid(ng) % zice, &
142# endif
143# if defined SEDIMENT && defined SED_MORPH
144 & sedbed(ng) % bed_thick, &
145# endif
146 & grid(ng) % Hz, &
147 & grid(ng) % z_r, &
148 & grid(ng) % z_w, &
149# endif
150 & mixing(ng) % Kh, &
151# ifdef SOLVE3D
152 & mixing(ng) % Kv, &
153# endif
154# ifdef ADJUST_BOUNDARY
155# ifdef SOLVE3D
156 & boundary(ng) % b_t_obc, &
157 & boundary(ng) % b_u_obc, &
158 & boundary(ng) % b_v_obc, &
159# endif
160 & boundary(ng) % b_ubar_obc, &
161 & boundary(ng) % b_vbar_obc, &
162 & boundary(ng) % b_zeta_obc, &
163# endif
164# ifdef ADJUST_WSTRESS
165 & forces(ng) % b_sustr, &
166 & forces(ng) % b_svstr, &
167# endif
168# if defined ADJUST_STFLUX && defined SOLVE3D
169 & forces(ng) % b_stflx, &
170# endif
171# ifdef SOLVE3D
172 & ocean(ng) % b_t, &
173 & ocean(ng) % b_u, &
174 & ocean(ng) % b_v, &
175# endif
176 & ocean(ng) % b_zeta, &
177 & ocean(ng) % b_ubar, &
178 & ocean(ng) % b_vbar, &
179# ifdef ADJUST_BOUNDARY
180# ifdef SOLVE3D
181 & boundary(ng) % tl_t_obc, &
182 & boundary(ng) % tl_u_obc, &
183 & boundary(ng) % tl_v_obc, &
184# endif
185 & boundary(ng) % tl_ubar_obc, &
186 & boundary(ng) % tl_vbar_obc, &
187 & boundary(ng) % tl_zeta_obc, &
188# endif
189# ifdef ADJUST_WSTRESS
190 & forces(ng) % tl_ustr, &
191 & forces(ng) % tl_vstr, &
192# endif
193# if defined ADJUST_STFLUX && defined SOLVE3D
194 & forces(ng) % tl_tflux, &
195# endif
196# ifdef SOLVE3D
197 & ocean(ng) % tl_t, &
198 & ocean(ng) % tl_u, &
199 & ocean(ng) % tl_v, &
200# endif
201 & ocean(ng) % tl_ubar, &
202 & ocean(ng) % tl_vbar, &
203 & ocean(ng) % tl_zeta)
204
205 RETURN
206 END SUBROUTINE tl_convolution
207!
208!***********************************************************************
209 SUBROUTINE tl_convolution_tile (ng, tile, &
210 & LBi, UBi, LBj, UBj, LBij, UBij, &
211 & IminS, ImaxS, JminS, JmaxS, &
212 & nstp, nnew, Linp, Lweak, ifac, &
213 & pm, om_p, om_r, om_u, om_v, &
214 & pn, on_p, on_r, on_u, on_v, &
215 & pmon_p, pmon_r, pmon_u, &
216 & pnom_p, pnom_r, pnom_v, &
217# ifdef MASKING
218 & rmask, pmask, umask, vmask, &
219# endif
220# ifdef SOLVE3D
221 & h, &
222# ifdef ICESHELF
223 & zice, &
224# endif
225# if defined SEDIMENT && defined SED_MORPH
226 & bed_thick, &
227# endif
228 & Hz, z_r, z_w, &
229# endif
230 & Kh, &
231# ifdef SOLVE3D
232 & Kv, &
233# endif
234# ifdef ADJUST_BOUNDARY
235# ifdef SOLVE3D
236 & VnormRobc, VnormUobc, VnormVobc, &
237# endif
238 & HnormRobc, HnormUobc, HnormVobc, &
239# endif
240# ifdef ADJUST_WSTRESS
241 & HnormSUS, HnormSVS, &
242# endif
243# if defined ADJUST_STFLUX && defined SOLVE3D
244 & HnormSTF, &
245# endif
246# ifdef SOLVE3D
247 & VnormR, VnormU, VnormV, &
248# endif
249 & HnormR, HnormU, HnormV, &
250# ifdef ADJUST_BOUNDARY
251# ifdef SOLVE3D
252 & tl_t_obc, tl_u_obc, tl_v_obc, &
253# endif
254 & tl_ubar_obc, tl_vbar_obc, &
255 & tl_zeta_obc, &
256# endif
257# ifdef ADJUST_WSTRESS
258 & tl_ustr, tl_vstr, &
259# endif
260# if defined ADJUST_STFLUX && defined SOLVE3D
261 & tl_tflux, &
262# endif
263# ifdef SOLVE3D
264 & tl_t, tl_u, tl_v, &
265# endif
266 & tl_ubar, tl_vbar, &
267 & tl_zeta)
268!***********************************************************************
269!
270 USE mod_param
271 USE mod_fourdvar
272 USE mod_ncparam
273 USE mod_scalars
274!
276# ifdef SOLVE3D
278# endif
279# ifdef ADJUST_BOUNDARY
281# ifdef SOLVE3D
283# endif
284# endif
285# ifdef DISTRIBUTE
287# endif
288 USE set_depth_mod
289!
290! Imported variable declarations.
291!
292 logical, intent(in) :: Lweak
293
294 integer, intent(in) :: ng, tile
295 integer, intent(in) :: LBi, UBi, LBj, UBj, LBij, UBij
296 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
297 integer, intent(in) :: nstp, nnew, Linp, ifac
298!
299# ifdef ASSUMED_SHAPE
300 real(r8), intent(in) :: pm(LBi:,LBj:)
301 real(r8), intent(in) :: om_p(LBi:,LBj:)
302 real(r8), intent(in) :: om_r(LBi:,LBj:)
303 real(r8), intent(in) :: om_u(LBi:,LBj:)
304 real(r8), intent(in) :: om_v(LBi:,LBj:)
305 real(r8), intent(in) :: pn(LBi:,LBj:)
306 real(r8), intent(in) :: on_p(LBi:,LBj:)
307 real(r8), intent(in) :: on_r(LBi:,LBj:)
308 real(r8), intent(in) :: on_u(LBi:,LBj:)
309 real(r8), intent(in) :: on_v(LBi:,LBj:)
310 real(r8), intent(in) :: pmon_p(LBi:,LBj:)
311 real(r8), intent(in) :: pmon_r(LBi:,LBj:)
312 real(r8), intent(in) :: pmon_u(LBi:,LBj:)
313 real(r8), intent(in) :: pnom_p(LBi:,LBj:)
314 real(r8), intent(in) :: pnom_r(LBi:,LBj:)
315 real(r8), intent(in) :: pnom_v(LBi:,LBj:)
316# ifdef MASKING
317 real(r8), intent(in) :: rmask(LBi:,LBj:)
318 real(r8), intent(in) :: pmask(LBi:,LBj:)
319 real(r8), intent(in) :: umask(LBi:,LBj:)
320 real(r8), intent(in) :: vmask(LBi:,LBj:)
321# endif
322 real(r8), intent(in) :: Kh(LBi:,LBj:)
323# ifdef SOLVE3D
324 real(r8), intent(in) :: Kv(LBi:,LBj:,0:)
325# ifdef ICESHELF
326 real(r8), intent(in) :: zice(LBi:,LBj:)
327# endif
328# if defined SEDIMENT && defined SED_MORPH
329 real(r8), intent(inout):: bed_thick(LBi:,LBj:,:)
330# endif
331 real(r8), intent(inout) :: h(LBi:,LBj:)
332# endif
333# ifdef ADJUST_BOUNDARY
334# ifdef SOLVE3D
335 real(r8), intent (in) :: VnormRobc(LBij:,:,:,:)
336 real(r8), intent (in) :: VnormUobc(LBij:,:,:)
337 real(r8), intent (in) :: VnormVobc(LBij:,:,:)
338# endif
339 real(r8), intent (in) :: HnormRobc(LBij:,:)
340 real(r8), intent (in) :: HnormUobc(LBij:,:)
341 real(r8), intent (in) :: HnormVobc(LBij:,:)
342# endif
343# ifdef ADJUST_WSTRESS
344 real(r8), intent(in) :: HnormSUS(LBi:,LBj:)
345 real(r8), intent(in) :: HnormSVS(LBi:,LBj:)
346# endif
347# if defined ADJUST_STFLUX && defined SOLVE3D
348 real(r8), intent(in) :: HnormSTF(LBi:,LBj:,:)
349# endif
350# ifdef SOLVE3D
351 real(r8), intent(in) :: VnormR(LBi:,LBj:,:,:,:)
352 real(r8), intent(in) :: VnormU(LBi:,LBj:,:,:)
353 real(r8), intent(in) :: VnormV(LBi:,LBj:,:,:)
354# endif
355 real(r8), intent(in) :: HnormR(LBi:,LBj:,:)
356 real(r8), intent(in) :: HnormU(LBi:,LBj:,:)
357 real(r8), intent(in) :: HnormV(LBi:,LBj:,:)
358# ifdef ADJUST_BOUNDARY
359# ifdef SOLVE3D
360 real(r8), intent(inout) :: tl_t_obc(LBij:,:,:,:,:,:)
361 real(r8), intent(inout) :: tl_u_obc(LBij:,:,:,:,:)
362 real(r8), intent(inout) :: tl_v_obc(LBij:,:,:,:,:)
363# endif
364 real(r8), intent(inout) :: tl_ubar_obc(LBij:,:,:,:)
365 real(r8), intent(inout) :: tl_vbar_obc(LBij:,:,:,:)
366 real(r8), intent(inout) :: tl_zeta_obc(LBij:,:,:,:)
367# endif
368# ifdef ADJUST_WSTRESS
369 real(r8), intent(inout) :: tl_ustr(LBi:,LBj:,:,:)
370 real(r8), intent(inout) :: tl_vstr(LBi:,LBj:,:,:)
371# endif
372# if defined ADJUST_STFLUX && defined SOLVE3D
373 real(r8), intent(inout) :: tl_tflux(LBi:,LBj:,:,:,:)
374# endif
375# ifdef SOLVE3D
376 real(r8), intent(inout) :: tl_t(LBi:,LBj:,:,:,:)
377 real(r8), intent(inout) :: tl_u(LBi:,LBj:,:,:)
378 real(r8), intent(inout) :: tl_v(LBi:,LBj:,:,:)
379# endif
380 real(r8), intent(inout) :: tl_ubar(LBi:,LBj:,:)
381 real(r8), intent(inout) :: tl_vbar(LBi:,LBj:,:)
382 real(r8), intent(inout) :: tl_zeta(LBi:,LBj:,:)
383# ifdef SOLVE3D
384 real(r8), intent(out) :: Hz(LBi:,LBj:,:)
385 real(r8), intent(out) :: z_r(LBi:,LBj:,:)
386 real(r8), intent(out) :: z_w(LBi:,LBj:,0:)
387# endif
388# else
389 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
390 real(r8), intent(in) :: om_p(LBi:UBi,LBj:UBj)
391 real(r8), intent(in) :: om_r(LBi:UBi,LBj:UBj)
392 real(r8), intent(in) :: om_u(LBi:UBi,LBj:UBj)
393 real(r8), intent(in) :: om_v(LBi:UBi,LBj:UBj)
394 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
395 real(r8), intent(in) :: on_p(LBi:UBi,LBj:UBj)
396 real(r8), intent(in) :: on_r(LBi:UBi,LBj:UBj)
397 real(r8), intent(in) :: on_u(LBi:UBi,LBj:UBj)
398 real(r8), intent(in) :: on_v(LBi:UBi,LBj:UBj)
399 real(r8), intent(in) :: pmon_p(LBi:UBi,LBj:UBj)
400 real(r8), intent(in) :: pmon_r(LBi:UBi,LBj:UBj)
401 real(r8), intent(in) :: pmon_u(LBi:UBi,LBj:UBj)
402 real(r8), intent(in) :: pnom_p(LBi:UBi,LBj:UBj)
403 real(r8), intent(in) :: pnom_r(LBi:UBi,LBj:UBj)
404 real(r8), intent(in) :: pnom_v(LBi:UBi,LBj:UBj)
405# ifdef MASKING
406 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
407 real(r8), intent(in) :: pmask(LBi:UBi,LBj:UBj)
408 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
409 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
410# endif
411 real(r8), intent(in) :: Kh(LBi:UBi,LBj:UBj)
412# ifdef SOLVE3D
413 real(r8), intent(in) :: Kv(LBi:UBi,LBj:UBj,0:N(ng))
414# ifdef ICESHELF
415 real(r8), intent(in) :: zice(LBi:UBi,LBj:UBj)
416# endif
417# if defined SEDIMENT && defined SED_MORPH
418 real(r8), intent(inout):: bed_thick(LBi:UBi,LBj:UBj,:)
419# endif
420 real(r8), intent(inout) :: h(LBi:UBi,LBj:UBj)
421# endif
422# ifdef ADJUST_BOUNDARY
423# ifdef SOLVE3D
424 real(r8), intent (in) :: VnormRobc(LBij:UBij,N(ng),4,NT(ng))
425 real(r8), intent (in) :: VnormUobc(LBij:UBij,N(ng),4)
426 real(r8), intent (in) :: VnormVobc(LBij:UBij,N(ng),4)
427# endif
428 real(r8), intent (in) :: HnormRobc(LBij:UBij,4)
429 real(r8), intent (in) :: HnormUobc(LBij:UBij,4)
430 real(r8), intent (in) :: HnormVobc(LBij:UBij,4)
431# endif
432# ifdef ADJUST_WSTRESS
433 real(r8), intent(in) :: HnormSUS(LBi:UBi,LBj:UBj)
434 real(r8), intent(in) :: HnormSVS(LBi:UBi,LBj:UBj)
435# endif
436# if defined ADJUST_STFLUX && defined SOLVE3D
437 real(r8), intent(in) :: HnormSTF(LBi:UBi,LBj:UBj,NT(ng))
438# endif
439# ifdef SOLVE3D
440 real(r8), intent(in) :: VnormR(LBi:UBi,LBj:UBj,N(ng),NSA,NT(ng))
441 real(r8), intent(in) :: VnormU(LBi:UBi,LBj:UBj,NSA,N(ng))
442 real(r8), intent(in) :: VnormV(LBi:UBi,LBj:UBj,NSA,N(ng))
443# endif
444 real(r8), intent(in) :: HnormR(LBi:UBi,LBj:UBj,NSA)
445 real(r8), intent(in) :: HnormU(LBi:UBi,LBj:UBj,NSA)
446 real(r8), intent(in) :: HnormV(LBi:UBi,LBj:UBj,NSA)
447# ifdef ADJUST_BOUNDARY
448# ifdef SOLVE3D
449 real(r8), intent(inout) :: tl_t_obc(LBij:UBij,N(ng),4, &
450 & Nbrec(ng),2,NT(ng))
451 real(r8), intent(inout) :: tl_u_obc(LBij:UBij,N(ng),4,Nbrec(ng),2)
452 real(r8), intent(inout) :: tl_v_obc(LBij:UBij,N(ng),4,Nbrec(ng),2)
453# endif
454 real(r8), intent(inout) :: tl_ubar_obc(LBij:UBij,4,Nbrec(ng),2)
455 real(r8), intent(inout) :: tl_vbar_obc(LBij:UBij,4,Nbrec(ng),2)
456 real(r8), intent(inout) :: tl_zeta_obc(LBij:UBij,4,Nbrec(ng),2)
457# endif
458# ifdef ADJUST_WSTRESS
459 real(r8), intent(inout) :: tl_ustr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
460 real(r8), intent(inout) :: tl_vstr(LBi:UBi,LBj:UBj,Nfrec(ng),2)
461# endif
462# if defined ADJUST_STFLUX && defined SOLVE3D
463 real(r8), intent(inout) :: tl_tflux(LBi:UBi,LBj:UBj, &
464 & Nfrec(ng),2,NT(ng))
465# endif
466# ifdef SOLVE3D
467 real(r8), intent(inout) :: tl_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
468 real(r8), intent(inout) :: tl_u(LBi:UBi,LBj:UBj,N(ng),2)
469 real(r8), intent(inout) :: tl_v(LBi:UBi,LBj:UBj,N(ng),2)
470# endif
471 real(r8), intent(inout) :: tl_ubar(LBi:UBi,LBj:UBj,:)
472 real(r8), intent(inout) :: tl_vbar(LBi:UBi,LBj:UBj,:)
473 real(r8), intent(inout) :: tl_zeta(LBi:UBi,LBj:UBj,:)
474# ifdef SOLVE3D
475 real(r8), intent(out) :: Hz(LBi:UBi,LBj:UBj,N(ng))
476 real(r8), intent(out) :: z_r(LBi:UBi,LBj:UBj,N(ng))
477 real(r8), intent(out) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
478# endif
479# endif
480!
481! Local variable declarations.
482!
483# ifdef ADJUST_BOUNDARY
484 logical, dimension(4) :: Lconvolve
485# endif
486 integer :: i, ib, ir, is, it, j, k, rec
487 real(r8) :: cff
488# ifdef SOLVE3D
489 real(r8) :: fac
490# endif
491# ifdef SOLVE3D
492 real(r8), dimension(LBi:UBi,LBj:UBj) :: work
493# endif
494!
495# include "set_bounds.h"
496!
497!-----------------------------------------------------------------------
498! Determine error covariance normalization factors to use.
499!-----------------------------------------------------------------------
500!
501 IF (lweak) THEN
502 rec=2 ! weak constraint
503 ELSE
504 rec=1 ! strong constraint
505 END IF
506
507# ifdef ADJUST_BOUNDARY
508!
509! Set switch to convolve boundary segments by the appropriate
510! tiles.
511!
512 lconvolve(iwest )=domain(ng)%Western_Edge (tile)
513 lconvolve(ieast )=domain(ng)%Eastern_Edge (tile)
514 lconvolve(isouth)=domain(ng)%Southern_Edge(tile)
515 lconvolve(inorth)=domain(ng)%Northern_Edge(tile)
516# endif
517
518# ifdef SOLVE3D
519!
520!-----------------------------------------------------------------------
521! Compute time invariant depths (use zero free-surface).
522!-----------------------------------------------------------------------
523!
524 DO i=lbi,ubi
525 DO j=lbj,ubj
526 work(i,j)=0.0_r8
527 END DO
528 END DO
529
530 CALL set_depth_tile (ng, tile, itlm, &
531 & lbi, ubi, lbj, ubj, &
532 & imins, imaxs, jmins, jmaxs, &
533 & nstp, nnew, &
534 & h, &
535# ifdef ICESHELF
536 & zice, &
537# endif
538# if defined SEDIMENT && defined SED_MORPH
539 & bed_thick, &
540# endif
541 & work, &
542 & hz, z_r, z_w)
543# endif
544!
545!-----------------------------------------------------------------------
546! Multiply tangent linear state by the inverse squared root of its
547! associated area (2D) or volume (3D).
548!-----------------------------------------------------------------------
549!
550! Tangent linear free-surface.
551!
552 DO j=jstrt,jendt
553 DO i=istrt,iendt
554 tl_zeta(i,j,linp)=tl_zeta(i,j,linp)/ &
555 & sqrt(om_r(i,j)*on_r(i,j))
556 END DO
557 END DO
558!
559! Tangent linear 2D momentum.
560!
561 DO j=jstrt,jendt
562 DO i=istrp,iendt
563 tl_ubar(i,j,linp)=tl_ubar(i,j,linp)/ &
564 & sqrt(om_u(i,j)*on_u(i,j))
565 END DO
566 END DO
567 DO j=jstrp,jendt
568 DO i=istrt,iendt
569 tl_vbar(i,j,linp)=tl_vbar(i,j,linp)/ &
570 & sqrt(om_v(i,j)*on_v(i,j))
571 END DO
572 END DO
573# ifdef DISTRIBUTE
574 CALL mp_exchange2d (ng, tile, itlm, 3, &
575 & lbi, ubi, lbj, ubj, &
576 & nghostpoints, &
577 & ewperiodic(ng), nsperiodic(ng), &
578 & tl_zeta(:,:,linp), &
579 & tl_ubar(:,:,linp), &
580 & tl_vbar(:,:,linp))
581# endif
582# ifdef SOLVE3D
583!
584! Tangent linear 3D momentum.
585!
586 DO j=jstrt,jendt
587 DO i=istrp,iendt
588 cff=om_u(i,j)*on_u(i,j)*0.5_r8
589 DO k=1,n(ng)
590 tl_u(i,j,k,linp)=tl_u(i,j,k,linp)/ &
591 & sqrt(cff*(hz(i-1,j,k)+hz(i,j,k)))
592 END DO
593 END DO
594 END DO
595 DO j=jstrp,jendt
596 DO i=istrt,iendt
597 cff=om_v(i,j)*on_v(i,j)*0.5_r8
598 DO k=1,n(ng)
599 tl_v(i,j,k,linp)=tl_v(i,j,k,linp)/ &
600 & sqrt(cff*(hz(i,j-1,k)+hz(i,j,k)))
601 END DO
602 END DO
603 END DO
604# ifdef DISTRIBUTE
605 CALL mp_exchange3d (ng, tile, itlm, 2, &
606 & lbi, ubi, lbj, ubj, 1, n(ng), &
607 & nghostpoints, &
608 & ewperiodic(ng), nsperiodic(ng), &
609 & tl_u(:,:,:,linp), &
610 & tl_v(:,:,:,linp))
611# endif
612!
613! Tangent linear tracers.
614!
615 DO j=jstrt,jendt
616 DO i=istrt,iendt
617 cff=om_r(i,j)*on_r(i,j)
618 DO k=1,n(ng)
619 fac=1.0_r8/sqrt(cff*hz(i,j,k))
620 DO it=1,nt(ng)
621 tl_t(i,j,k,linp,it)=fac*tl_t(i,j,k,linp,it)
622 END DO
623 END DO
624 END DO
625 END DO
626# ifdef DISTRIBUTE
627 CALL mp_exchange4d (ng, tile, itlm, 1, &
628 & lbi, ubi, lbj, ubj, 1, n(ng), 1, nt(ng), &
629 & nghostpoints, &
630 & ewperiodic(ng), nsperiodic(ng), &
631 & tl_t(:,:,:,linp,:))
632# endif
633# endif
634
635# ifdef ADJUST_BOUNDARY
636!
637! Tangent linear free-surface open boundaries.
638!
639 DO ir=1,nbrec(ng)
640 DO ib=1,4
641 IF (.not.lweak.and.lobc(ib,isfsur,ng)) THEN
642 IF (lconvolve(ib)) THEN
643 IF ((ib.eq.iwest).or.(ib.eq.ieast)) THEN
644 i=bounds(ng)%edge(ib,r2dvar)
645 DO j=jstrt,jendt
646 tl_zeta_obc(j,ib,ir,linp)=tl_zeta_obc(j,ib,ir,linp)/ &
647 & sqrt(on_r(i,j))
648 END DO
649 ELSE IF ((ib.eq.isouth).or.(ib.eq.inorth)) THEN
650 j=bounds(ng)%edge(ib,r2dvar)
651 DO i=istrt,iendt
652 tl_zeta_obc(i,ib,ir,linp)=tl_zeta_obc(i,ib,ir,linp)/ &
653 & sqrt(om_r(i,j))
654
655 END DO
656 END IF
657 END IF
658# ifdef DISTRIBUTE
659 CALL mp_exchange2d_bry (ng, tile, itlm, 1, ib, &
660 & lbij, ubij, &
661 & nghostpoints, &
662 & ewperiodic(ng), nsperiodic(ng), &
663 & tl_zeta_obc(:,ib,ir,linp))
664# endif
665 END IF
666 END DO
667 END DO
668!
669! Tangent linear 2D U-momentum open boundaries.
670!
671 DO ir=1,nbrec(ng)
672 DO ib=1,4
673 IF (.not.lweak.and.lobc(ib,isubar,ng)) THEN
674 IF (lconvolve(ib)) THEN
675 IF ((ib.eq.iwest).or.(ib.eq.ieast)) THEN
676 i=bounds(ng)%edge(ib,u2dvar)
677 DO j=jstrt,jendt
678 tl_ubar_obc(j,ib,ir,linp)=tl_ubar_obc(j,ib,ir,linp)/ &
679 & sqrt(on_u(i,j))
680 END DO
681 ELSE IF ((ib.eq.isouth).or.(ib.eq.inorth)) THEN
682 j=bounds(ng)%edge(ib,u2dvar)
683 DO i=istrp,iendt
684 tl_ubar_obc(i,ib,ir,linp)=tl_ubar_obc(i,ib,ir,linp)/ &
685 & sqrt(om_u(i,j))
686 END DO
687 END IF
688 END IF
689# ifdef DISTRIBUTE
690 CALL mp_exchange2d_bry (ng, tile, itlm, 1, ib, &
691 & lbij, ubij, &
692 & nghostpoints, &
693 & ewperiodic(ng), nsperiodic(ng), &
694 & tl_ubar_obc(:,ib,ir,linp))
695# endif
696 END IF
697 END DO
698 END DO
699!
700! Tangent linear 2D V-momentum open boundaries.
701!
702 DO ir=1,nbrec(ng)
703 DO ib=1,4
704 IF (.not.lweak.and.lobc(ib,isvbar,ng)) THEN
705 IF (lconvolve(ib)) THEN
706 IF ((ib.eq.iwest).or.(ib.eq.ieast)) THEN
707 i=bounds(ng)%edge(ib,v2dvar)
708 DO j=jstrp,jendt
709 tl_vbar_obc(j,ib,ir,linp)=tl_vbar_obc(j,ib,ir,linp)/ &
710 & sqrt(on_v(i,j))
711 END DO
712 ELSE IF ((ib.eq.isouth).or.(ib.eq.inorth)) THEN
713 j=bounds(ng)%edge(ib,v2dvar)
714 DO i=istrt,iendt
715 tl_vbar_obc(i,ib,ir,linp)=tl_vbar_obc(i,ib,ir,linp)/ &
716 & sqrt(om_v(i,j))
717 END DO
718 END IF
719 END IF
720# ifdef DISTRIBUTE
721 CALL mp_exchange2d_bry (ng, tile, itlm, 1, ib, &
722 & lbij, ubij, &
723 & nghostpoints, &
724 & ewperiodic(ng), nsperiodic(ng), &
725 & tl_vbar_obc(:,ib,ir,linp))
726# endif
727 END IF
728 END DO
729 END DO
730
731# ifdef SOLVE3D
732!
733! Tangent linear 3D U-momentum open boundaries.
734!
735 DO ir=1,nbrec(ng)
736 DO ib=1,4
737 IF (.not.lweak.and.lobc(ib,isuvel,ng)) THEN
738 IF (lconvolve(ib)) THEN
739 IF ((ib.eq.iwest).or.(ib.eq.ieast)) THEN
740 i=bounds(ng)%edge(ib,u2dvar)
741 DO j=jstrt,jendt
742 cff=on_u(i,j)*0.5_r8
743 DO k=1,n(ng)
744 tl_u_obc(j,k,ib,ir,linp)=tl_u_obc(j,k,ib,ir,linp)/ &
745 & sqrt(cff* &
746 & (hz(i-1,j,k)+ &
747 & hz(i ,j,k)))
748 END DO
749 END DO
750 ELSE IF ((ib.eq.isouth).or.(ib.eq.inorth)) THEN
751 j=bounds(ng)%edge(ib,u2dvar)
752 DO i=istrp,iendt
753 cff=om_u(i,j)*0.5_r8
754 DO k=1,n(ng)
755 tl_u_obc(i,k,ib,ir,linp)=tl_u_obc(i,k,ib,ir,linp)/ &
756 & sqrt(cff* &
757 & (hz(i-1,j,k)+ &
758 & hz(i ,j,k)))
759 END DO
760 END DO
761 END IF
762 END IF
763# ifdef DISTRIBUTE
764 CALL mp_exchange3d_bry (ng, tile, itlm, 1, ib, &
765 & lbij, ubij, 1, n(ng), &
766 & nghostpoints, &
767 & ewperiodic(ng), nsperiodic(ng), &
768 & tl_u_obc(:,:,ib,ir,linp))
769# endif
770 END IF
771 END DO
772 END DO
773!
774! Tangent linear 3D V-momentum open boundaries.
775!
776 DO ir=1,nbrec(ng)
777 DO ib=1,4
778 IF (.not.lweak.and.lobc(ib,isvvel,ng)) THEN
779 IF (lconvolve(ib)) THEN
780 IF ((ib.eq.iwest).or.(ib.eq.ieast)) THEN
781 i=bounds(ng)%edge(ib,v2dvar)
782 DO j=jstrp,jendt
783 cff=on_v(i,j)*0.5_r8
784 DO k=1,n(ng)
785 tl_v_obc(j,k,ib,ir,linp)=tl_v_obc(j,k,ib,ir,linp)/ &
786 & sqrt(cff* &
787 & (hz(i,j-1,k)+ &
788 & hz(i,j ,k)))
789 END DO
790 END DO
791 ELSE IF ((ib.eq.isouth).or.(ib.eq.inorth)) THEN
792 j=bounds(ng)%edge(ib,v2dvar)
793 DO i=istrt,iendt
794 cff=om_v(i,j)*0.5_r8
795 DO k=1,n(ng)
796 tl_v_obc(i,k,ib,ir,linp)=tl_v_obc(i,k,ib,ir,linp)/ &
797 & sqrt(cff* &
798 & (hz(i,j-1,k)+ &
799 & hz(i,j ,k)))
800 END DO
801 END DO
802 END IF
803 END IF
804# ifdef DISTRIBUTE
805 CALL mp_exchange3d_bry (ng, tile, itlm, 1, ib, &
806 & lbij, ubij, 1, n(ng), &
807 & nghostpoints, &
808 & ewperiodic(ng), nsperiodic(ng), &
809 & tl_v_obc(:,:,ib,ir,linp))
810# endif
811 END IF
812 END DO
813 END DO
814!
815! Tangent linear tracers open boundaries.
816!
817 DO it=1,nt(ng)
818 DO ir=1,nbrec(ng)
819 DO ib=1,4
820 IF (.not.lweak.and.lobc(ib,istvar(it),ng)) THEN
821 IF (lconvolve(ib)) THEN
822 IF ((ib.eq.iwest).or.(ib.eq.ieast)) THEN
823 i=bounds(ng)%edge(ib,r2dvar)
824 DO j=jstrt,jendt
825 cff=on_r(i,j)
826 DO k=1,n(ng)
827 tl_t_obc(j,k,ib,ir,linp,it)= &
828 & tl_t_obc(j,k,ib,ir,linp,it)/ &
829 & sqrt(cff*hz(i,j,k))
830 END DO
831 END DO
832 ELSE IF ((ib.eq.isouth).or.(ib.eq.inorth)) THEN
833 j=bounds(ng)%edge(ib,r2dvar)
834 DO i=istrt,iendt
835 cff=om_r(i,j)
836 DO k=1,n(ng)
837 tl_t_obc(i,k,ib,ir,linp,it)= &
838 & tl_t_obc(i,k,ib,ir,linp,it)/ &
839 & sqrt(cff*hz(i,j,k))
840 END DO
841 END DO
842 END IF
843 END IF
844# ifdef DISTRIBUTE
845 CALL mp_exchange3d_bry (ng, tile, itlm, 1, ib, &
846 & lbij, ubij, 1, n(ng), &
847 & nghostpoints, &
848 & ewperiodic(ng), nsperiodic(ng), &
849 & tl_t_obc(:,:,ib,ir,linp,it))
850# endif
851 END IF
852 END DO
853 END DO
854 END DO
855# endif
856# endif
857
858# ifdef ADJUST_WSTRESS
859!
860! Tangent linear surface momentum stress.
861!
862 IF (.not.lweak) THEN
863 DO ir=1,nfrec(ng)
864 DO j=jstrt,jendt
865 DO i=istrp,iendt
866 tl_ustr(i,j,ir,linp)=tl_ustr(i,j,ir,linp)/ &
867 & sqrt(om_u(i,j)*on_u(i,j))
868 END DO
869 END DO
870 DO j=jstrp,jendt
871 DO i=istrt,iendt
872 tl_vstr(i,j,ir,linp)=tl_vstr(i,j,ir,linp)/ &
873 & sqrt(om_v(i,j)*on_v(i,j))
874 END DO
875 END DO
876 END DO
877# ifdef DISTRIBUTE
878 CALL mp_exchange3d (ng, tile, itlm, 2, &
879 & lbi, ubi, lbj, ubj, 1, nfrec(ng), &
880 & nghostpoints, &
881 & ewperiodic(ng), nsperiodic(ng), &
882 & tl_ustr(:,:,:,linp), &
883 & tl_vstr(:,:,:,linp))
884# endif
885 END IF
886# endif
887# if defined ADJUST_STFLUX && defined SOLVE3D
888!
889! Tangent linear surface tracers flux.
890!
891 IF (.not.lweak) THEN
892 DO j=jstrt,jendt
893 DO i=istrt,iendt
894 fac=1.0_r8/sqrt(om_r(i,j)*on_r(i,j))
895 DO it=1,nt(ng)
896 IF (lstflux(it,ng)) THEN
897 DO ir=1,nfrec(ng)
898 tl_tflux(i,j,ir,linp,it)=fac*tl_tflux(i,j,ir,linp,it)
899 END DO
900 END IF
901 END DO
902 END DO
903 END DO
904# ifdef DISTRIBUTE
905 DO it=1,nt(ng)
906 IF (lstflux(it,ng)) THEN
907 CALL mp_exchange3d (ng, tile, itlm, 1, &
908 & lbi, ubi, lbj, ubj, 1, nfrec(ng), &
909 & nghostpoints, &
910 & ewperiodic(ng), nsperiodic(ng), &
911 & tl_tflux(:,:,:,linp,it))
912 END IF
913 END DO
914# endif
915 END IF
916# endif
917!
918!-----------------------------------------------------------------------
919! Initial conditions and model error convariance: Convolve tangent
920! linear state vector with a generalized diffusion equation to filter
921! solution with specified horizontal scales. Convert from minimization
922! space (v-space) to model space.
923!-----------------------------------------------------------------------
924!
925! Tangent linear free-surface.
926!
927 CALL tl_conv_r2d_tile (ng, tile, itlm, &
928 & lbi, ubi, lbj, ubj, &
929 & imins, imaxs, jmins, jmaxs, &
930 & nghostpoints, &
931 & nhsteps(rec,isfsur)/ifac, &
932 & dtsizeh(rec,isfsur), &
933 & kh, &
934 & pm, pn, pmon_u, pnom_v, &
935# ifdef MASKING
936 & rmask, umask, vmask, &
937# endif
938 & tl_zeta(:,:,linp))
939!
940! Tangent linear 2D momentum.
941!
942 CALL tl_conv_u2d_tile (ng, tile, itlm, &
943 & lbi, ubi, lbj, ubj, &
944 & imins, imaxs, jmins, jmaxs, &
945 & nghostpoints, &
946 & nhsteps(rec,isubar)/ifac, &
947 & dtsizeh(rec,isubar), &
948 & kh, &
949 & pm, pn, pmon_r, pnom_p, &
950# ifdef MASKING
951 & umask, pmask, &
952# endif
953 & tl_ubar(:,:,linp))
954
955 CALL tl_conv_v2d_tile (ng, tile, itlm, &
956 & lbi, ubi, lbj, ubj, &
957 & imins, imaxs, jmins, jmaxs, &
958 & nghostpoints, &
959 & nhsteps(rec,isvbar)/ifac, &
960 & dtsizeh(rec,isvbar), &
961 & kh, &
962 & pm, pn, pmon_p, pnom_r, &
963# ifdef MASKING
964 & vmask, pmask, &
965# endif
966 & tl_vbar(:,:,linp))
967# ifdef SOLVE3D
968!
969! Tangent linear 3D momentum.
970!
971 CALL tl_conv_u3d_tile (ng, tile, itlm, &
972 & lbi, ubi, lbj, ubj, 1, n(ng), &
973 & imins, imaxs, jmins, jmaxs, &
974 & nghostpoints, &
975 & nhsteps(rec,isuvel)/ifac, &
976 & nvsteps(rec,isuvel)/ifac, &
977 & dtsizeh(rec,isuvel), &
978 & dtsizev(rec,isuvel), &
979 & kh, kv, &
980 & pm, pn, &
981# ifdef GEOPOTENTIAL_HCONV
982 & on_r, om_p, &
983# else
984 & pmon_r, pnom_p, &
985# endif
986# ifdef MASKING
987# ifdef GEOPOTENTIAL_HCONV
988 & pmask, rmask, umask, vmask, &
989# else
990 & umask, pmask, &
991# endif
992# endif
993 & hz, z_r, &
994 & tl_u(:,:,:,linp))
995
996 CALL tl_conv_v3d_tile (ng, tile, itlm, &
997 & lbi, ubi, lbj, ubj, 1, n(ng), &
998 & imins, imaxs, jmins, jmaxs, &
999 & nghostpoints, &
1000 & nhsteps(rec,isvvel)/ifac, &
1001 & nvsteps(rec,isvvel)/ifac, &
1002 & dtsizeh(rec,isvvel), &
1003 & dtsizev(rec,isvvel), &
1004 & kh, kv, &
1005 & pm, pn, &
1006# ifdef GEOPOTENTIAL_HCONV
1007 & on_p, om_r, &
1008# else
1009 & pmon_p, pnom_r, &
1010# endif
1011# ifdef MASKING
1012# ifdef GEOPOTENTIAL_HCONV
1013 & pmask, rmask, umask, vmask, &
1014# else
1015 & vmask, pmask, &
1016# endif
1017# endif
1018 & hz, z_r, &
1019 & tl_v(:,:,:,linp))
1020!
1021! Tangent linear tracers.
1022!
1023 DO it=1,nt(ng)
1024 is=istvar(it)
1025 CALL tl_conv_r3d_tile (ng, tile, itlm, &
1026 & lbi, ubi, lbj, ubj, 1, n(ng), &
1027 & imins, imaxs, jmins, jmaxs, &
1028 & nghostpoints, &
1029 & nhsteps(rec,is)/ifac, &
1030 & nvsteps(rec,is)/ifac, &
1031 & dtsizeh(rec,is), &
1032 & dtsizev(rec,is), &
1033 & kh, kv, &
1034 & pm, pn, &
1035# ifdef GEOPOTENTIAL_HCONV
1036 & on_u, om_v, &
1037# else
1038 & pmon_u, pnom_v, &
1039# endif
1040# ifdef MASKING
1041 & rmask, umask, vmask, &
1042# endif
1043 & hz, z_r, &
1044 & tl_t(:,:,:,linp,it))
1045 END DO
1046# endif
1047
1048# ifdef ADJUST_BOUNDARY
1049!
1050!-----------------------------------------------------------------------
1051! Open boundaries error convariance: Convolve tangent linear boundary
1052! edges with a generalized diffusion equation to filter solution with
1053! specified horizontal scales. Convert from minimization space
1054! (v-space) to model space.
1055!-----------------------------------------------------------------------
1056!
1057! Tangent linear free-surface open boundaries.
1058!
1059 DO ir=1,nbrec(ng)
1060 DO ib=1,4
1061 IF (.not.lweak.and.lobc(ib,isfsur,ng)) THEN
1062 CALL tl_conv_r2d_bry_tile (ng, tile, itlm, ib, &
1063 & bounds(ng)%edge(:,r2dvar), &
1064 & lbij, ubij, &
1065 & lbi, ubi, lbj, ubj, &
1066 & imins, imaxs, jmins, jmaxs, &
1067 & nghostpoints, &
1068 & nhstepsb(ib,isfsur)/ifac, &
1069 & dtsizehb(ib,isfsur), &
1070 & kh, &
1071 & pm, pn, pmon_u, pnom_v, &
1072# ifdef MASKING
1073 & rmask, umask, vmask, &
1074# endif
1075 & tl_zeta_obc(:,ib,ir,linp))
1076 END IF
1077 END DO
1078 END DO
1079!
1080! Tangent linear 2D U-momentum open boundaries.
1081!
1082 DO ir=1,nbrec(ng)
1083 DO ib=1,4
1084 IF (.not.lweak.and.lobc(ib,isubar,ng)) THEN
1085 CALL tl_conv_u2d_bry_tile (ng, tile, itlm, ib, &
1086 & bounds(ng)%edge(:,u2dvar), &
1087 & lbij, ubij, &
1088 & lbi, ubi, lbj, ubj, &
1089 & imins, imaxs, jmins, jmaxs, &
1090 & nghostpoints, &
1091 & nhstepsb(ib,isubar)/ifac, &
1092 & dtsizehb(ib,isubar), &
1093 & kh, &
1094 & pm, pn, pmon_r, pnom_p, &
1095# ifdef MASKING
1096 & umask, pmask, &
1097# endif
1098 & tl_ubar_obc(:,ib,ir,linp))
1099 END IF
1100 END DO
1101 END DO
1102!
1103! Tangent linear 2D V-momentum open boundaries.
1104!
1105 DO ir=1,nbrec(ng)
1106 DO ib=1,4
1107 IF (.not.lweak.and.lobc(ib,isvbar,ng)) THEN
1108 CALL tl_conv_v2d_bry_tile (ng, tile, itlm, ib, &
1109 & bounds(ng)%edge(:,v2dvar), &
1110 & lbij, ubij, &
1111 & lbi, ubi, lbj, ubj, &
1112 & imins, imaxs, jmins, jmaxs, &
1113 & nghostpoints, &
1114 & nhstepsb(ib,isvbar)/ifac, &
1115 & dtsizehb(ib,isvbar), &
1116 & kh, &
1117 & pm, pn, pmon_p, pnom_r, &
1118# ifdef MASKING
1119 & vmask, pmask, &
1120# endif
1121 & tl_vbar_obc(:,ib,ir,linp))
1122 END IF
1123 END DO
1124 END DO
1125
1126# ifdef SOLVE3D
1127!
1128! Tangent linear 3D U-momentum open boundaries.
1129!
1130 DO ir=1,nbrec(ng)
1131 DO ib=1,4
1132 IF (.not.lweak.and.lobc(ib,isuvel,ng)) THEN
1133 CALL tl_conv_u3d_bry_tile (ng, tile, itlm, ib, &
1134 & bounds(ng)%edge(:,u2dvar), &
1135 & lbij, ubij, &
1136 & lbi, ubi, lbj, ubj, 1, n(ng), &
1137 & imins, imaxs, jmins, jmaxs, &
1138 & nghostpoints, &
1139 & nhstepsb(ib,isuvel)/ifac, &
1140 & nvstepsb(ib,isuvel)/ifac, &
1141 & dtsizehb(ib,isuvel), &
1142 & dtsizevb(ib,isuvel), &
1143 & kh, kv, &
1144 & pm, pn, pmon_r, pnom_p, &
1145# ifdef MASKING
1146 & umask, pmask, &
1147# endif
1148 & hz, z_r, &
1149 & tl_u_obc(:,:,ib,ir,linp))
1150 END IF
1151 END DO
1152 END DO
1153!
1154! Tangent linear 3D V-momentum open boundaries.
1155!
1156 DO ir=1,nbrec(ng)
1157 DO ib=1,4
1158 IF (.not.lweak.and.lobc(ib,isvvel,ng)) THEN
1159 CALL tl_conv_v3d_bry_tile (ng, tile, itlm, ib, &
1160 & bounds(ng)%edge(:,v2dvar), &
1161 & lbij, ubij, &
1162 & lbi, ubi, lbj, ubj, 1, n(ng), &
1163 & imins, imaxs, jmins, jmaxs, &
1164 & nghostpoints, &
1165 & nhstepsb(ib,isvvel)/ifac, &
1166 & nvstepsb(ib,isvvel)/ifac, &
1167 & dtsizehb(ib,isvvel), &
1168 & dtsizevb(ib,isvvel), &
1169 & kh, kv, &
1170 & pm, pn, pmon_p, pnom_r, &
1171# ifdef MASKING
1172 & vmask, pmask, &
1173# endif
1174 & hz, z_r, &
1175 & tl_v_obc(:,:,ib,ir,linp))
1176 END IF
1177 END DO
1178 END DO
1179!
1180! Tangent linear tracers open boundaries.
1181!
1182 DO it=1,nt(ng)
1183 is=istvar(it)
1184 DO ir=1,nbrec(ng)
1185 DO ib=1,4
1186 IF (.not.lweak.and.lobc(ib,is,ng)) THEN
1187 CALL tl_conv_r3d_bry_tile (ng, tile, itlm, ib, &
1188 & bounds(ng)%edge(:,r2dvar), &
1189 & lbij, ubij, &
1190 & lbi, ubi, lbj, ubj, 1, n(ng), &
1191 & imins, imaxs, jmins, jmaxs, &
1192 & nghostpoints, &
1193 & nhstepsb(ib,is)/ifac, &
1194 & nvstepsb(ib,is)/ifac, &
1195 & dtsizehb(ib,is), &
1196 & dtsizevb(ib,is), &
1197 & kh, kv, &
1198 & pm, pn, pmon_u, pnom_v, &
1199# ifdef MASKING
1200 & rmask, umask, vmask, &
1201# endif
1202 & hz, z_r, &
1203 & tl_t_obc(:,:,ib,ir,linp,it))
1204 END IF
1205 END DO
1206 END DO
1207 END DO
1208# endif
1209# endif
1210
1211# if defined ADJUST_WSTRESS || defined ADJUST_STFLUX
1212!
1213!-----------------------------------------------------------------------
1214! Surface forcing error convariance: Convolve tangent linear state
1215! vector with a generalized diffusion equation to filter solution with
1216! specified horizontal scales. Convert from minimization space
1217! (v-space) to model space.
1218!-----------------------------------------------------------------------
1219
1220# ifdef ADJUST_WSTRESS
1221!
1222! Tangent linear surface momentum stress.
1223!
1224 IF (.not.lweak) THEN
1225 DO k=1,nfrec(ng)
1226 CALL tl_conv_u2d_tile (ng, tile, itlm, &
1227 & lbi, ubi, lbj, ubj, &
1228 & imins, imaxs, jmins, jmaxs, &
1229 & nghostpoints, &
1230 & nhsteps(rec,isustr)/ifac, &
1231 & dtsizeh(rec,isustr), &
1232 & kh, &
1233 & pm, pn, pmon_r, pnom_p, &
1234# ifdef MASKING
1235 & umask, pmask, &
1236# endif
1237 & tl_ustr(:,:,k,linp))
1238
1239 CALL tl_conv_v2d_tile (ng, tile, itlm, &
1240 & lbi, ubi, lbj, ubj, &
1241 & imins, imaxs, jmins, jmaxs, &
1242 & nghostpoints, &
1243 & nhsteps(rec,isvstr)/ifac, &
1244 & dtsizeh(rec,isvstr), &
1245 & kh, &
1246 & pm, pn, pmon_p, pnom_r, &
1247# ifdef MASKING
1248 & vmask, pmask, &
1249# endif
1250 & tl_vstr(:,:,k,linp))
1251 END DO
1252 END IF
1253# endif
1254# if defined ADJUST_STFLUX && defined SOLVE3D
1255!
1256! Tangent linear surface tracers flux.
1257!
1258 IF (.not.lweak) THEN
1259 DO it=1,nt(ng)
1260 IF (lstflux(it,ng)) THEN
1261 is=istsur(it)
1262 DO k=1,nfrec(ng)
1263 CALL tl_conv_r2d_tile (ng, tile, itlm, &
1264 & lbi, ubi, lbj, ubj, &
1265 & imins, imaxs, jmins, jmaxs, &
1266 & nghostpoints, &
1267 & nhsteps(rec,is)/ifac, &
1268 & dtsizeh(rec,is), &
1269 & kh, &
1270 & pm, pn, pmon_u, pnom_v, &
1271# ifdef MASKING
1272 & rmask, umask, vmask, &
1273# endif
1274 & tl_tflux(:,:,k,linp,it))
1275 END DO
1276 END IF
1277 END DO
1278 END IF
1279# endif
1280# endif
1281!
1282!-----------------------------------------------------------------------
1283! Multiply convolved tangent linear state by its corresponding
1284! normalization factor.
1285!-----------------------------------------------------------------------
1286!
1287! Tangent linear free-surface.
1288!
1289 DO j=jstrt,jendt
1290 DO i=istrt,iendt
1291 tl_zeta(i,j,linp)=tl_zeta(i,j,linp)*hnormr(i,j,rec)
1292 END DO
1293 END DO
1294!
1295! Tangent linear 2D momentum.
1296!
1297 DO j=jstrt,jendt
1298 DO i=istrp,iendt
1299 tl_ubar(i,j,linp)=tl_ubar(i,j,linp)*hnormu(i,j,rec)
1300 END DO
1301 END DO
1302 DO j=jstrp,jendt
1303 DO i=istrt,iendt
1304 tl_vbar(i,j,linp)=tl_vbar(i,j,linp)*hnormv(i,j,rec)
1305 END DO
1306 END DO
1307# ifdef DISTRIBUTE
1308 CALL mp_exchange2d (ng, tile, itlm, 3, &
1309 & lbi, ubi, lbj, ubj, &
1310 & nghostpoints, &
1311 & ewperiodic(ng), nsperiodic(ng), &
1312 & tl_zeta(:,:,linp), &
1313 & tl_ubar(:,:,linp), &
1314 & tl_vbar(:,:,linp))
1315# endif
1316# ifdef SOLVE3D
1317!
1318! Tangent linear 3D momentum.
1319!
1320 DO k=1,n(ng)
1321 DO j=jstrt,jendt
1322 DO i=istrp,iendt
1323 tl_u(i,j,k,linp)=tl_u(i,j,k,linp)*vnormu(i,j,k,rec)
1324 END DO
1325 END DO
1326 DO j=jstrp,jendt
1327 DO i=istrt,iendt
1328 tl_v(i,j,k,linp)=tl_v(i,j,k,linp)*vnormv(i,j,k,rec)
1329 END DO
1330 END DO
1331 END DO
1332# ifdef DISTRIBUTE
1333 CALL mp_exchange3d (ng, tile, itlm, 2, &
1334 & lbi, ubi, lbj, ubj, 1, n(ng), &
1335 & nghostpoints, &
1336 & ewperiodic(ng), nsperiodic(ng), &
1337 & tl_u(:,:,:,linp), &
1338 & tl_v(:,:,:,linp))
1339# endif
1340!
1341! Tangent linear tracers.
1342!
1343 DO it=1,nt(ng)
1344 DO k=1,n(ng)
1345 DO j=jstrt,jendt
1346 DO i=istrt,iendt
1347 tl_t(i,j,k,linp,it)=tl_t(i,j,k,linp,it)* &
1348 & vnormr(i,j,k,rec,it)
1349 END DO
1350 END DO
1351 END DO
1352 END DO
1353# ifdef DISTRIBUTE
1354 CALL mp_exchange4d (ng, tile, itlm, 1, &
1355 & lbi, ubi, lbj, ubj, 1, n(ng), 1, nt(ng), &
1356 & nghostpoints, &
1357 & ewperiodic(ng), nsperiodic(ng), &
1358 & tl_t(:,:,:,linp,:))
1359# endif
1360# endif
1361
1362# ifdef ADJUST_BOUNDARY
1363!
1364! Tangent linear free-surface open boundaries.
1365!
1366 DO ir=1,nbrec(ng)
1367 DO ib=1,4
1368 IF (.not.lweak.and.lobc(ib,isfsur,ng)) THEN
1369 IF (lconvolve(ib)) THEN
1370 IF ((ib.eq.iwest).or.(ib.eq.ieast)) THEN
1371 DO j=jstrt,jendt
1372 tl_zeta_obc(j,ib,ir,linp)=tl_zeta_obc(j,ib,ir,linp)* &
1373 & hnormrobc(j,ib)
1374 END DO
1375 ELSE IF ((ib.eq.isouth).or.(ib.eq.inorth)) THEN
1376 DO i=istrt,iendt
1377 tl_zeta_obc(i,ib,ir,linp)=tl_zeta_obc(i,ib,ir,linp)* &
1378 & hnormrobc(i,ib)
1379 END DO
1380 END IF
1381 END IF
1382# ifdef DISTRIBUTE
1383 CALL mp_exchange2d_bry (ng, tile, itlm, 1, ib, &
1384 & lbij, ubij, &
1385 & nghostpoints, &
1386 & ewperiodic(ng), nsperiodic(ng), &
1387 & tl_zeta_obc(:,ib,ir,linp))
1388# endif
1389 END IF
1390 END DO
1391 END DO
1392!
1393! Tangent linear 2D U-momentum open boundaries.
1394!
1395 DO ir=1,nbrec(ng)
1396 DO ib=1,4
1397 IF (.not.lweak.and.lobc(ib,isubar,ng)) THEN
1398 IF (lconvolve(ib)) THEN
1399 IF ((ib.eq.iwest).or.(ib.eq.ieast)) THEN
1400 DO j=jstrt,jendt
1401 tl_ubar_obc(j,ib,ir,linp)=tl_ubar_obc(j,ib,ir,linp)* &
1402 & hnormuobc(j,ib)
1403 END DO
1404 ELSE IF ((ib.eq.isouth).or.(ib.eq.inorth)) THEN
1405 DO i=istrp,iendt
1406 tl_ubar_obc(i,ib,ir,linp)=tl_ubar_obc(i,ib,ir,linp)* &
1407 & hnormuobc(i,ib)
1408 END DO
1409 END IF
1410 END IF
1411# ifdef DISTRIBUTE
1412 CALL mp_exchange2d_bry (ng, tile, itlm, 1, ib, &
1413 & lbij, ubij, &
1414 & nghostpoints, &
1415 & ewperiodic(ng), nsperiodic(ng), &
1416 & tl_ubar_obc(:,ib,ir,linp))
1417# endif
1418 END IF
1419 END DO
1420 END DO
1421!
1422! Tangent linear 2D V-momentum open boundaries.
1423!
1424 DO ir=1,nbrec(ng)
1425 DO ib=1,4
1426 IF (.not.lweak.and.lobc(ib,isvbar,ng)) THEN
1427 IF (lconvolve(ib)) THEN
1428 IF ((ib.eq.iwest).or.(ib.eq.ieast)) THEN
1429 DO j=jstrp,jendt
1430 tl_vbar_obc(j,ib,ir,linp)=tl_vbar_obc(j,ib,ir,linp)* &
1431 & hnormvobc(j,ib)
1432 END DO
1433 ELSE IF ((ib.eq.isouth).or.(ib.eq.inorth)) THEN
1434 DO i=istrt,iendt
1435 tl_vbar_obc(i,ib,ir,linp)=tl_vbar_obc(i,ib,ir,linp)* &
1436 & hnormvobc(i,ib)
1437 END DO
1438 END IF
1439 END IF
1440# ifdef DISTRIBUTE
1441 CALL mp_exchange2d_bry (ng, tile, itlm, 1, ib, &
1442 & lbij, ubij, &
1443 & nghostpoints, &
1444 & ewperiodic(ng), nsperiodic(ng), &
1445 & tl_vbar_obc(:,ib,ir,linp))
1446# endif
1447 END IF
1448 END DO
1449 END DO
1450
1451# ifdef SOLVE3D
1452!
1453! Tangent linear 3D U-momentum open boundaries.
1454!
1455 DO ir=1,nbrec(ng)
1456 DO ib=1,4
1457 IF (.not.lweak.and.lobc(ib,isuvel,ng)) THEN
1458 IF (lconvolve(ib)) THEN
1459 IF ((ib.eq.iwest).or.(ib.eq.ieast)) THEN
1460 DO k=1,n(ng)
1461 DO j=jstrt,jendt
1462 tl_u_obc(j,k,ib,ir,linp)=tl_u_obc(j,k,ib,ir,linp)* &
1463 & vnormuobc(j,k,ib)
1464 END DO
1465 END DO
1466 ELSE IF ((ib.eq.isouth).or.(ib.eq.inorth)) THEN
1467 DO k=1,n(ng)
1468 DO i=istrp,iendt
1469 tl_u_obc(i,k,ib,ir,linp)=tl_u_obc(i,k,ib,ir,linp)* &
1470 & vnormuobc(i,k,ib)
1471 END DO
1472 END DO
1473 END IF
1474 END IF
1475# ifdef DISTRIBUTE
1476 CALL mp_exchange3d_bry (ng, tile, itlm, 1, ib, &
1477 & lbij, ubij, 1, n(ng), &
1478 & nghostpoints, &
1479 & ewperiodic(ng), nsperiodic(ng), &
1480 & tl_u_obc(:,:,ib,ir,linp))
1481# endif
1482 END IF
1483 END DO
1484 END DO
1485!
1486! Tangent linear 3D V-momentum open boundaries.
1487!
1488 DO ir=1,nbrec(ng)
1489 DO ib=1,4
1490 IF (.not.lweak.and.lobc(ib,isvvel,ng)) THEN
1491 IF (lconvolve(ib)) THEN
1492 IF ((ib.eq.iwest).or.(ib.eq.ieast)) THEN
1493 DO k=1,n(ng)
1494 DO j=jstrp,jendt
1495 tl_v_obc(j,k,ib,ir,linp)=tl_v_obc(j,k,ib,ir,linp)* &
1496 & vnormvobc(j,k,ib)
1497 END DO
1498 END DO
1499 ELSE IF ((ib.eq.isouth).or.(ib.eq.inorth)) THEN
1500 DO k=1,n(ng)
1501 DO i=istrt,iendt
1502 tl_v_obc(i,k,ib,ir,linp)=tl_v_obc(i,k,ib,ir,linp)* &
1503 & vnormvobc(i,k,ib)
1504 END DO
1505 END DO
1506 END IF
1507 END IF
1508# ifdef DISTRIBUTE
1509 CALL mp_exchange3d_bry (ng, tile, itlm, 1, ib, &
1510 & lbij, ubij, 1, n(ng), &
1511 & nghostpoints, &
1512 & ewperiodic(ng), nsperiodic(ng), &
1513 & tl_v_obc(:,:,ib,ir,linp))
1514# endif
1515 END IF
1516 END DO
1517 END DO
1518!
1519! Tangent linear tracers open boundaries.
1520!
1521 DO it=1,nt(ng)
1522 DO ir=1,nbrec(ng)
1523 DO ib=1,4
1524 IF (.not.lweak.and.lobc(ib,istvar(it),ng)) THEN
1525 IF (lconvolve(ib)) THEN
1526 IF ((ib.eq.iwest).or.(ib.eq.ieast)) THEN
1527 DO k=1,n(ng)
1528 DO j=jstrt,jendt
1529 tl_t_obc(j,k,ib,ir,linp,it)= &
1530 & tl_t_obc(j,k,ib,ir,linp,it)* &
1531 & vnormrobc(j,k,ib,it)
1532 END DO
1533 END DO
1534 ELSE IF ((ib.eq.isouth).or.(ib.eq.inorth)) THEN
1535 DO k=1,n(ng)
1536 DO i=istrt,iendt
1537 tl_t_obc(i,k,ib,ir,linp,it)= &
1538 & tl_t_obc(i,k,ib,ir,linp,it)* &
1539 & vnormrobc(i,k,ib,it)
1540 END DO
1541 END DO
1542 END IF
1543 END IF
1544# ifdef DISTRIBUTE
1545 CALL mp_exchange3d_bry (ng, tile, itlm, 1, ib, &
1546 & lbij, ubij, 1, n(ng), &
1547 & nghostpoints, &
1548 & ewperiodic(ng), nsperiodic(ng), &
1549 & tl_t_obc(:,:,ib,ir,linp,it))
1550# endif
1551 END IF
1552 END DO
1553 END DO
1554 END DO
1555# endif
1556# endif
1557
1558# ifdef ADJUST_WSTRESS
1559!
1560! Tangent linear surface momentum stress.
1561!
1562 IF (.not.lweak) THEN
1563 DO k=1,nfrec(ng)
1564 DO j=jstrt,jendt
1565 DO i=istrp,iendt
1566 tl_ustr(i,j,k,linp)=tl_ustr(i,j,k,linp)*hnormsus(i,j)
1567 END DO
1568 END DO
1569 DO j=jstrp,jendt
1570 DO i=istrt,iendt
1571 tl_vstr(i,j,k,linp)=tl_vstr(i,j,k,linp)*hnormsvs(i,j)
1572 END DO
1573 END DO
1574 END DO
1575# ifdef DISTRIBUTE
1576 CALL mp_exchange3d (ng, tile, itlm, 2, &
1577 & lbi, ubi, lbj, ubj, 1, nfrec(ng), &
1578 & nghostpoints, &
1579 & ewperiodic(ng), nsperiodic(ng), &
1580 & tl_ustr(:,:,:,linp), &
1581 & tl_vstr(:,:,:,linp))
1582# endif
1583 END IF
1584# endif
1585
1586# if defined ADJUST_STFLUX && defined SOLVE3D
1587!
1588! Tangent linear surface tracers flux.
1589!
1590 IF (.not.lweak) THEN
1591 DO it=1,nt(ng)
1592 IF (lstflux(it,ng)) THEN
1593 DO k=1,nfrec(ng)
1594 DO j=jstrt,jendt
1595 DO i=istrt,iendt
1596 tl_tflux(i,j,k,linp,it)=tl_tflux(i,j,k,linp,it)* &
1597 & hnormstf(i,j,it)
1598 END DO
1599 END DO
1600 END DO
1601# ifdef DISTRIBUTE
1602 CALL mp_exchange3d (ng, tile, itlm, 1, &
1603 & lbi, ubi, lbj, ubj, 1, nfrec(ng), &
1604 & nghostpoints, &
1605 & ewperiodic(ng), nsperiodic(ng), &
1606 & tl_tflux(:,:,:,linp,it))
1607# endif
1608 END IF
1609 END DO
1610 END IF
1611# endif
1612
1613 RETURN
1614 END SUBROUTINE tl_convolution_tile
1615#endif
1616 END MODULE tl_convolution_mod
type(t_boundary), dimension(:), allocatable boundary
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
integer, dimension(:,:), allocatable nhsteps
integer, dimension(:,:), allocatable nvsteps
real(r8), dimension(:,:), allocatable dtsizehb
real(r8), dimension(:,:), allocatable dtsizeh
integer, dimension(:,:), allocatable nvstepsb
real(r8), dimension(:,:), allocatable dtsizevb
real(r8), dimension(:,:), allocatable dtsizev
integer, dimension(:,:), allocatable nhstepsb
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_mixing), dimension(:), allocatable mixing
Definition mod_mixing.F:399
integer isvvel
integer isvbar
integer isvstr
integer, dimension(:), allocatable istvar
integer isuvel
integer isfsur
integer isustr
integer isubar
integer, dimension(:), allocatable istsur
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
type(t_bounds), dimension(:), allocatable bounds
Definition mod_param.F:232
integer nghostpoints
Definition mod_param.F:710
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer, parameter u2dvar
Definition mod_param.F:718
integer, parameter itlm
Definition mod_param.F:663
integer, parameter r2dvar
Definition mod_param.F:717
integer, parameter v2dvar
Definition mod_param.F:719
logical, dimension(:,:,:), allocatable lobc
logical, dimension(:), allocatable ewperiodic
integer, parameter iwest
logical, dimension(:), allocatable nsperiodic
logical, dimension(:,:), allocatable lstflux
integer, parameter isouth
integer, parameter ieast
integer, parameter inorth
type(t_sedbed), dimension(:), allocatable sedbed
Definition sedbed_mod.h:157
integer, dimension(:), allocatable nnew
integer, dimension(:), allocatable nstp
subroutine mp_exchange2d_bry(ng, tile, model, nvar, boundary, lbij, ubij, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine mp_exchange4d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, lbt, ubt, nghost, ew_periodic, ns_periodic, a, b, c)
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine mp_exchange3d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine mp_exchange3d_bry(ng, tile, model, nvar, boundary, lbij, ubij, lbk, ubk, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine, public set_depth_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nstp, nnew, h, zice, zt_avg1, hz, z_r, z_w)
Definition set_depth.F:86
subroutine tl_conv_u2d_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nghost, nhsteps, dtsizeh, kh, pm, pn, pmon_r, pnom_p, umask, pmask, tl_a)
Definition tl_conv_2d.F:301
subroutine tl_conv_r2d_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nghost, nhsteps, dtsizeh, kh, pm, pn, pmon_u, pnom_v, rmask, umask, vmask, tl_a)
Definition tl_conv_2d.F:70
subroutine tl_conv_v2d_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nghost, nhsteps, dtsizeh, kh, pm, pn, pmon_p, pnom_r, vmask, pmask, tl_a)
Definition tl_conv_2d.F:530
subroutine tl_conv_u3d_tile(ng, tile, model, lbi, ubi, lbj, ubj, lbk, ubk, imins, imaxs, jmins, jmaxs, nghost, nhsteps, nvsteps, dtsizeh, dtsizev, kh, kv, pm, pn, on_r, om_p, pmask, rmask, umask, vmask, hz, z_r, tl_a)
Definition tl_conv_3d.F:847
subroutine tl_conv_v3d_tile(ng, tile, model, lbi, ubi, lbj, ubj, lbk, ubk, imins, imaxs, jmins, jmaxs, nghost, nhsteps, nvsteps, dtsizeh, dtsizev, kh, kv, pm, pn, on_p, om_r, pmask, rmask, umask, vmask, hz, z_r, tl_a)
subroutine tl_conv_r3d_tile(ng, tile, model, lbi, ubi, lbj, ubj, lbk, ubk, imins, imaxs, jmins, jmaxs, nghost, nhsteps, nvsteps, dtsizeh, dtsizev, kh, kv, pm, pn, on_u, om_v, rmask, umask, vmask, hz, z_r, tl_a)
Definition tl_conv_3d.F:85
subroutine tl_conv_v2d_bry_tile(ng, tile, model, boundary, edge, lbij, ubij, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nghost, nhsteps, dtsizeh, kh, pm, pn, pmon_p, pnom_r, vmask, pmask, tl_a)
subroutine tl_conv_u2d_bry_tile(ng, tile, model, boundary, edge, lbij, ubij, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nghost, nhsteps, dtsizeh, kh, pm, pn, pmon_r, pnom_p, umask, pmask, tl_a)
subroutine tl_conv_r2d_bry_tile(ng, tile, model, boundary, edge, lbij, ubij, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nghost, nhsteps, dtsizeh, kh, pm, pn, pmon_u, pnom_v, rmask, umask, vmask, tl_a)
subroutine tl_conv_v3d_bry_tile(ng, tile, model, boundary, edge, lbij, ubij, lbi, ubi, lbj, ubj, lbk, ubk, imins, imaxs, jmins, jmaxs, nghost, nhsteps, nvsteps, dtsizeh, dtsizev, kh, kv, pm, pn, pmon_p, pnom_r, vmask, pmask, hz, z_r, tl_a)
subroutine tl_conv_u3d_bry_tile(ng, tile, model, boundary, edge, lbij, ubij, lbi, ubi, lbj, ubj, lbk, ubk, imins, imaxs, jmins, jmaxs, nghost, nhsteps, nvsteps, dtsizeh, dtsizev, kh, kv, pm, pn, pmon_r, pnom_p, umask, pmask, hz, z_r, tl_a)
subroutine tl_conv_r3d_bry_tile(ng, tile, model, boundary, edge, lbij, ubij, lbi, ubi, lbj, ubj, lbk, ubk, imins, imaxs, jmins, jmaxs, nghost, nhsteps, nvsteps, dtsizeh, dtsizev, kh, kv, pm, pn, pmon_u, pnom_v, rmask, umask, vmask, hz, z_r, tl_a)
subroutine tl_convolution_tile(ng, tile, lbi, ubi, lbj, ubj, lbij, ubij, imins, imaxs, jmins, jmaxs, nstp, nnew, linp, lweak, ifac, pm, om_p, om_r, om_u, om_v, pn, on_p, on_r, on_u, on_v, pmon_p, pmon_r, pmon_u, pnom_p, pnom_r, pnom_v, rmask, pmask, umask, vmask, h, zice, bed_thick, hz, z_r, z_w, kh, kv, vnormrobc, vnormuobc, vnormvobc, hnormrobc, hnormuobc, hnormvobc, hnormsus, hnormsvs, hnormstf, vnormr, vnormu, vnormv, hnormr, hnormu, hnormv, 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_ubar, tl_vbar, tl_zeta)
subroutine, public tl_convolution(ng, tile, linp, lweak, ifac)