ROMS
Loading...
Searching...
No Matches
conv_2d.F
Go to the documentation of this file.
1#include "cppdefs.h"
2
3 MODULE conv_2d_mod
4
5#if defined NONLINEAR && defined FOUR_DVAR
6!
7!git $Id$
8!================================================== Hernan G. Arango ===
9! Copyright (c) 2002-2025 The ROMS Group Andrew M. Moore !
10! Licensed under a MIT/X style license !
11! See License_ROMS.md !
12!=======================================================================
13! !
14! These routines applies the background error covariance to data !
15! assimilation fields via the space convolution of the diffusion !
16! equation (filter) for 2D state variables. The diffusion filter !
17! is solved using an explicit (inefficient) algorithm. !
18! !
19! For Gaussian (bell-shaped) correlations, the space convolution !
20! of the diffusion operator is an efficient way to estimate the !
21! finite domain error covariances. !
22! !
23! On Input: !
24! !
25! ng Nested grid number. !
26! model Calling model identifier. !
27! Istr Starting tile index in the I-direction. !
28! Iend Ending tile index in the I-direction. !
29! Jstr Starting tile index in the J-direction. !
30! Jend Ending tile index in the J-direction. !
31! LBi I-dimension Lower bound. !
32! UBi I-dimension Upper bound. !
33! LBj J-dimension Lower bound. !
34! UBj J-dimension Upper bound. !
35! Nghost Number of ghost points. !
36! NHsteps Number of horizontal diffusion integration steps. !
37! DTsizeH Horizontal diffusion pseudo time-step size. !
38! Kh Horizontal diffusion coefficients. !
39! A 2D state variable to diffuse. !
40! !
41! On Output: !
42! !
43! A Diffused 2D state variable. !
44! !
45! Routines: !
46! !
47! conv_r2d_tile Nonlinear 2D convolution at RHO-points !
48! conv_u2d_tile Nonlinear 2D convolution at U-points !
49! conv_v2d_tile Nonlinear 2D convolution at V-points !
50! !
51!=======================================================================
52!
53 implicit none
54!
55 PUBLIC
56!
57 CONTAINS
58!
59!***********************************************************************
60 SUBROUTINE conv_r2d_tile (ng, tile, model, &
61 & LBi, UBi, LBj, UBj, &
62 & IminS, ImaxS, JminS, JmaxS, &
63 & Nghost, NHsteps, DTsizeH, &
64 & Kh, &
65 & pm, pn, pmon_u, pnom_v, &
66# ifdef MASKING
67 & rmask, umask, vmask, &
68# endif
69 & A)
70!***********************************************************************
71!
72 USE mod_param
73 USE mod_scalars
74!
75 USE bc_2d_mod, ONLY: dabc_r2d_tile
76# ifdef DISTRIBUTE
77 USE mp_exchange_mod, ONLY : mp_exchange2d
78# endif
79!
80! Imported variable declarations.
81!
82 integer, intent(in) :: ng, tile, model
83 integer, intent(in) :: LBi, UBi, LBj, UBj
84 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
85 integer, intent(in) :: Nghost, NHsteps
86
87 real(r8), intent(in) :: DTsizeH
88!
89# ifdef ASSUMED_SHAPE
90 real(r8), intent(in) :: pm(LBi:,LBj:)
91 real(r8), intent(in) :: pn(LBi:,LBj:)
92 real(r8), intent(in) :: pmon_u(LBi:,LBj:)
93 real(r8), intent(in) :: pnom_v(LBi:,LBj:)
94# ifdef MASKING
95 real(r8), intent(in) :: rmask(LBi:,LBj:)
96 real(r8), intent(in) :: umask(LBi:,LBj:)
97 real(r8), intent(in) :: vmask(LBi:,LBj:)
98# endif
99 real(r8), intent(in) :: Kh(LBi:,LBj:)
100 real(r8), intent(inout) :: A(LBi:,LBj:)
101# else
102 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
103 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
104 real(r8), intent(in) :: pmon_u(LBi:UBi,LBj:UBj)
105 real(r8), intent(in) :: pnom_v(LBi:UBi,LBj:UBj)
106# ifdef MASKING
107 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
108 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
109 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
110# endif
111 real(r8), intent(in) :: Kh(LBi:UBi,LBj:UBj)
112 real(r8), intent(inout) :: A(LBi:UBi,LBj:UBj)
113# endif
114!
115! Local variable declarations.
116!
117 integer :: Nnew, Nold, Nsav, i, j, step
118
119 real(r8), dimension(LBi:UBi,LBj:UBj,2) :: Awrk
120
121 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FE
122 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FX
123 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Hfac
124
125# include "set_bounds.h"
126!
127!-----------------------------------------------------------------------
128! Space convolution of the diffusion equation for a 2D state variable
129! at RHO-points.
130!-----------------------------------------------------------------------
131!
132! Compute metrics factor.
133!
134 DO j=jstr,jend
135 DO i=istr,iend
136 hfac(i,j)=dtsizeh*pm(i,j)*pn(i,j)
137 END DO
138 END DO
139!
140! Set integration indices and initial conditions.
141!
142 nold=1
143 nnew=2
144 CALL dabc_r2d_tile (ng, tile, &
145 & lbi, ubi, lbj, ubj, &
146 & a)
147# ifdef DISTRIBUTE
148 CALL mp_exchange2d (ng, tile, model, 1, &
149 & lbi, ubi, lbj, ubj, &
150 & nghost, &
151 & ewperiodic(ng), nsperiodic(ng), &
152 & a)
153# endif
154 DO j=jstr-1,jend+1
155 DO i=istr-1,iend+1
156 awrk(i,j,nold)=a(i,j)
157 END DO
158 END DO
159!
160!-----------------------------------------------------------------------
161! Integrate horizontal diffusion terms.
162!-----------------------------------------------------------------------
163!
164 DO step=1,nhsteps
165!
166! Compute XI- and ETA-components of diffusive flux.
167!
168 DO j=jstr,jend
169 DO i=istr,iend+1
170 fx(i,j)=pmon_u(i,j)*0.5_r8*(kh(i-1,j)+kh(i,j))* &
171 & (awrk(i,j,nold)-awrk(i-1,j,nold))
172# ifdef MASKING
173 fx(i,j)=fx(i,j)*umask(i,j)
174# endif
175 END DO
176 END DO
177 DO j=jstr,jend+1
178 DO i=istr,iend
179 fe(i,j)=pnom_v(i,j)*0.5_r8*(kh(i,j-1)+kh(i,j))* &
180 & (awrk(i,j,nold)-awrk(i,j-1,nold))
181# ifdef MASKING
182 fe(i,j)=fe(i,j)*vmask(i,j)
183# endif
184 END DO
185 END DO
186!
187! Time-step horizontal diffusion terms.
188!
189 DO j=jstr,jend
190 DO i=istr,iend
191 awrk(i,j,nnew)=awrk(i,j,nold)+ &
192 & hfac(i,j)* &
193 & (fx(i+1,j)-fx(i,j)+ &
194 & fe(i,j+1)-fe(i,j))
195 END DO
196 END DO
197!
198! Apply boundary conditions. If applicable, exchange boundary data.
199!
200 CALL dabc_r2d_tile (ng, tile, &
201 & lbi, ubi, lbj, ubj, &
202 & awrk(:,:,nnew))
203# ifdef DISTRIBUTE
204 CALL mp_exchange2d (ng, tile, model, 1, &
205 & lbi, ubi, lbj, ubj, &
206 & nghost, &
207 & ewperiodic(ng), nsperiodic(ng), &
208 & awrk(:,:,nnew))
209# endif
210!
211! Update integration indices.
212!
213 nsav=nold
214 nold=nnew
215 nnew=nsav
216 END DO
217!
218!-----------------------------------------------------------------------
219! Load convolved solution.
220!-----------------------------------------------------------------------
221!
222 DO j=jstr,jend
223 DO i=istr,iend
224 a(i,j)=awrk(i,j,nold)
225 END DO
226 END DO
227 CALL dabc_r2d_tile (ng, tile, &
228 & lbi, ubi, lbj, ubj, &
229 & a)
230# ifdef DISTRIBUTE
231 CALL mp_exchange2d (ng, tile, model, 1, &
232 & lbi, ubi, lbj, ubj, &
233 & nghost, &
234 & ewperiodic(ng), nsperiodic(ng), &
235 & a)
236# endif
237
238 RETURN
239 END SUBROUTINE conv_r2d_tile
240!
241!***********************************************************************
242 SUBROUTINE conv_u2d_tile (ng, tile, model, &
243 & LBi, UBi, LBj, UBj, &
244 & IminS, ImaxS, JminS, JmaxS, &
245 & Nghost, NHsteps, DTsizeH, &
246 & Kh, &
247 & pm, pn, pmon_r, pnom_p, &
248# ifdef MASKING
249 & umask, pmask, &
250# endif
251 & A)
252!***********************************************************************
253!
254 USE mod_param
255 USE mod_scalars
256!
257 USE bc_2d_mod, ONLY: dabc_u2d_tile
258# ifdef DISTRIBUTE
259 USE mp_exchange_mod, ONLY : mp_exchange2d
260# endif
261!
262! Imported variable declarations.
263!
264 integer, intent(in) :: ng, tile, model
265 integer, intent(in) :: LBi, UBi, LBj, UBj
266 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
267 integer, intent(in) :: Nghost, NHsteps
268
269 real(r8), intent(in) :: DTsizeH
270!
271# ifdef ASSUMED_SHAPE
272 real(r8), intent(in) :: pm(LBi:,LBj:)
273 real(r8), intent(in) :: pn(LBi:,LBj:)
274 real(r8), intent(in) :: pmon_r(LBi:,LBj:)
275 real(r8), intent(in) :: pnom_p(LBi:,LBj:)
276# ifdef MASKING
277 real(r8), intent(in) :: umask(LBi:,LBj:)
278 real(r8), intent(in) :: pmask(LBi:,LBj:)
279# endif
280 real(r8), intent(in) :: Kh(LBi:,LBj:)
281 real(r8), intent(inout) :: A(LBi:,LBj:)
282# else
283 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
284 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
285 real(r8), intent(in) :: pmon_r(LBi:UBi,LBj:UBj)
286 real(r8), intent(in) :: pnom_p(LBi:UBi,LBj:UBj)
287# ifdef MASKING
288 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
289 real(r8), intent(in) :: pmask(LBi:UBi,LBj:UBj)
290# endif
291 real(r8), intent(in) :: Kh(LBi:UBi,LBj:UBj)
292 real(r8), intent(inout) :: A(LBi:UBi,LBj:UBj)
293# endif
294!
295! Local variable declarations.
296!
297 integer :: Nnew, Nold, Nsav, i, j, step
298
299 real(r8) :: cff
300
301 real(r8), dimension(LBi:UBi,LBj:UBj,2) :: Awrk
302
303 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FE
304 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FX
305 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Hfac
306
307# include "set_bounds.h"
308!
309!-----------------------------------------------------------------------
310! Space convolution of the diffusion equation for a 2D state variable
311! at U-points.
312!-----------------------------------------------------------------------
313!
314! Compute metrics factor.
315!
316 cff=dtsizeh*0.25_r8
317 DO j=jstr,jend
318 DO i=istru,iend
319 hfac(i,j)=cff*(pm(i-1,j)+pm(i,j))*(pn(i-1,j)+pn(i,j))
320 END DO
321 END DO
322!
323! Set integration indices and initial conditions.
324!
325 nold=1
326 nnew=2
327 CALL dabc_u2d_tile (ng, tile, &
328 & lbi, ubi, lbj, ubj, &
329 & a)
330# ifdef DISTRIBUTE
331 CALL mp_exchange2d (ng, tile, model, 1, &
332 & lbi, ubi, lbj, ubj, &
333 & nghost, &
334 & ewperiodic(ng), nsperiodic(ng), &
335 & a)
336# endif
337 DO j=jstr-1,jend+1
338 DO i=istru-1,iend+1
339 awrk(i,j,nold)=a(i,j)
340 END DO
341 END DO
342!
343!-----------------------------------------------------------------------
344! Integrate horizontal diffusion terms.
345!-----------------------------------------------------------------------
346!
347 DO step=1,nhsteps
348!
349! Compute XI- and ETA-components of diffusive flux.
350!
351 DO j=jstr,jend
352 DO i=istru-1,iend
353 fx(i,j)=pmon_r(i,j)*kh(i,j)* &
354 & (awrk(i+1,j,nold)-awrk(i,j,nold))
355 END DO
356 END DO
357 DO j=jstr,jend+1
358 DO i=istru,iend
359 fe(i,j)=pnom_p(i,j)*0.25_r8*(kh(i-1,j )+kh(i,j )+ &
360 & kh(i-1,j-1)+kh(i,j-1))* &
361 & (awrk(i,j,nold)-awrk(i,j-1,nold))
362# ifdef MASKING
363 fe(i,j)=fe(i,j)*pmask(i,j)
364# endif
365 END DO
366 END DO
367!
368! Time-step horizontal diffusion terms.
369!
370 DO j=jstr,jend
371 DO i=istru,iend
372 awrk(i,j,nnew)=awrk(i,j,nold)+ &
373 & hfac(i,j)* &
374 & (fx(i,j)-fx(i-1,j)+ &
375 & fe(i,j+1)-fe(i,j))
376 END DO
377 END DO
378!
379! Apply boundary conditions. If applicable, exchange boundary data.
380!
381 CALL dabc_u2d_tile (ng, tile, &
382 & lbi, ubi, lbj, ubj, &
383 & awrk(:,:,nnew))
384# ifdef DISTRIBUTE
385 CALL mp_exchange2d (ng, tile, model, 1, &
386 & lbi, ubi, lbj, ubj, &
387 & nghost, &
388 & ewperiodic(ng), nsperiodic(ng), &
389 & awrk(:,:,nnew))
390# endif
391!
392! Update integration indices.
393!
394 nsav=nold
395 nold=nnew
396 nnew=nsav
397 END DO
398!
399!-----------------------------------------------------------------------
400! Load convolved solution.
401!-----------------------------------------------------------------------
402!
403 DO j=jstr,jend
404 DO i=istru,iend
405 a(i,j)=awrk(i,j,nold)
406 END DO
407 END DO
408 CALL dabc_u2d_tile (ng, tile, &
409 & lbi, ubi, lbj, ubj, &
410 & a)
411# ifdef DISTRIBUTE
412 CALL mp_exchange2d (ng, tile, model, 1, &
413 & lbi, ubi, lbj, ubj, &
414 & nghost, &
415 & ewperiodic(ng), nsperiodic(ng), &
416 & a)
417# endif
418
419 RETURN
420 END SUBROUTINE conv_u2d_tile
421!
422!***********************************************************************
423 SUBROUTINE conv_v2d_tile (ng, tile, model, &
424 & LBi, UBi, LBj, UBj, &
425 & IminS, ImaxS, JminS, JmaxS, &
426 & Nghost, NHsteps, DTsizeH, &
427 & Kh, &
428 & pm, pn, pmon_p, pnom_r, &
429# ifdef MASKING
430 & vmask, pmask, &
431# endif
432 & A)
433!***********************************************************************
434!
435 USE mod_param
436 USE mod_scalars
437!
438 USE bc_2d_mod, ONLY: dabc_v2d_tile
439# ifdef DISTRIBUTE
440 USE mp_exchange_mod, ONLY : mp_exchange2d
441# endif
442!
443! Imported variable declarations.
444!
445 integer, intent(in) :: ng, tile, model
446 integer, intent(in) :: LBi, UBi, LBj, UBj
447 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
448 integer, intent(in) :: Nghost, NHsteps
449
450 real(r8), intent(in) :: DTsizeH
451!
452# ifdef ASSUMED_SHAPE
453 real(r8), intent(in) :: pm(LBi:,LBj:)
454 real(r8), intent(in) :: pn(LBi:,LBj:)
455 real(r8), intent(in) :: pmon_p(LBi:,LBj:)
456 real(r8), intent(in) :: pnom_r(LBi:,LBj:)
457# ifdef MASKING
458 real(r8), intent(in) :: vmask(LBi:,LBj:)
459 real(r8), intent(in) :: pmask(LBi:,LBj:)
460# endif
461 real(r8), intent(in) :: Kh(LBi:,LBj:)
462 real(r8), intent(inout) :: A(LBi:,LBj:)
463# else
464 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
465 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
466 real(r8), intent(in) :: pmon_p(LBi:UBi,LBj:UBj)
467 real(r8), intent(in) :: pnom_r(LBi:UBi,LBj:UBj)
468# ifdef MASKING
469 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
470 real(r8), intent(in) :: pmask(LBi:UBi,LBj:UBj)
471# endif
472 real(r8), intent(in) :: Kh(LBi:UBi,LBj:UBj)
473 real(r8), intent(inout) :: A(LBi:UBi,LBj:UBj)
474# endif
475!
476! Local variable declarations.
477!
478 integer :: Nnew, Nold, Nsav, i, j, step
479
480 real(r8) :: cff
481
482 real(r8), dimension(LBi:UBi,LBj:UBj,2) :: Awrk
483
484 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FE
485 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FX
486 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Hfac
487
488# include "set_bounds.h"
489!
490!-----------------------------------------------------------------------
491! Space convolution of the diffusion equation for a 2D state variable
492! at V-points.
493!-----------------------------------------------------------------------
494!
495! Compute metrics factor.
496!
497 cff=dtsizeh*0.25_r8
498 DO j=jstrv,jend
499 DO i=istr,iend
500 hfac(i,j)=cff*(pm(i,j-1)+pm(i,j))*(pn(i,j-1)+pn(i,j))
501 END DO
502 END DO
503!
504! Set integration indices and initial conditions.
505!
506 nold=1
507 nnew=2
508 CALL dabc_v2d_tile (ng, tile, &
509 & lbi, ubi, lbj, ubj, &
510 & a)
511# ifdef DISTRIBUTE
512 CALL mp_exchange2d (ng, tile, model, 1, &
513 & lbi, ubi, lbj, ubj, &
514 & nghost, &
515 & ewperiodic(ng), nsperiodic(ng), &
516 & a)
517# endif
518 DO j=jstrv-1,jend+1
519 DO i=istr-1,iend+1
520 awrk(i,j,nold)=a(i,j)
521 END DO
522 END DO
523!
524!-----------------------------------------------------------------------
525! Integrate horizontal diffusion terms.
526!-----------------------------------------------------------------------
527!
528 DO step=1,nhsteps
529!
530! Compute XI- and ETA-components of diffusive flux.
531!
532 DO j=jstrv,jend
533 DO i=istr,iend+1
534 fx(i,j)=pmon_p(i,j)*0.25_r8*(kh(i-1,j )+kh(i,j )+ &
535 & kh(i-1,j-1)+kh(i,j-1))* &
536 & (awrk(i,j,nold)-awrk(i-1,j,nold))
537# ifdef MASKING
538 fx(i,j)=fx(i,j)*pmask(i,j)
539# endif
540 END DO
541 END DO
542 DO j=jstrv-1,jend
543 DO i=istr,iend
544 fe(i,j)=pnom_r(i,j)*kh(i,j)* &
545 & (awrk(i,j+1,nold)-awrk(i,j,nold))
546 END DO
547 END DO
548!
549! Time-step horizontal diffusion terms.
550!
551 DO j=jstrv,jend
552 DO i=istr,iend
553 awrk(i,j,nnew)=awrk(i,j,nold)+ &
554 & hfac(i,j)* &
555 & (fx(i+1,j)-fx(i,j)+ &
556 & fe(i,j)-fe(i,j-1))
557 END DO
558 END DO
559!
560! Apply boundary conditions. If applicable, exchange boundary data.
561!
562 CALL dabc_v2d_tile (ng, tile, &
563 & lbi, ubi, lbj, ubj, &
564 & awrk(:,:,nnew))
565# ifdef DISTRIBUTE
566 CALL mp_exchange2d (ng, tile, model, 1, &
567 & lbi, ubi, lbj, ubj, &
568 & nghost, &
569 & ewperiodic(ng), nsperiodic(ng), &
570 & awrk(:,:,nnew))
571# endif
572!
573! Update integration indices.
574!
575 nsav=nold
576 nold=nnew
577 nnew=nsav
578 END DO
579!
580!-----------------------------------------------------------------------
581! Load convolved solution.
582!-----------------------------------------------------------------------
583!
584 DO j=jstrv,jend
585 DO i=istr,iend
586 a(i,j)=awrk(i,j,nold)
587 END DO
588 END DO
589 CALL dabc_v2d_tile (ng, tile, &
590 & lbi, ubi, lbj, ubj, &
591 & a)
592# ifdef DISTRIBUTE
593 CALL mp_exchange2d (ng, tile, model, 1, &
594 & lbi, ubi, lbj, ubj, &
595 & nghost, &
596 & ewperiodic(ng), nsperiodic(ng), &
597 & a)
598# endif
599
600 RETURN
601 END SUBROUTINE conv_v2d_tile
602#endif
603 END MODULE conv_2d_mod
bc_2d_mod
Definition bc_2d.F:2
bc_2d_mod::dabc_v2d_tile
subroutine dabc_v2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
Definition bc_2d.F:771
bc_2d_mod::dabc_u2d_tile
subroutine dabc_u2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
Definition bc_2d.F:646
bc_2d_mod::dabc_r2d_tile
subroutine dabc_r2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
Definition bc_2d.F:523
conv_2d_mod
Definition conv_2d.F:3
conv_2d_mod::conv_u2d_tile
subroutine 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, a)
Definition conv_2d.F:252
conv_2d_mod::conv_v2d_tile
subroutine 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, a)
Definition conv_2d.F:433
conv_2d_mod::conv_r2d_tile
subroutine 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, a)
Definition conv_2d.F:70
mod_param
Definition mod_param.F:2
mod_scalars
Definition mod_scalars.F:2
mod_scalars::ewperiodic
logical, dimension(:), allocatable ewperiodic
Definition mod_scalars.F:510
mod_scalars::nsperiodic
logical, dimension(:), allocatable nsperiodic
Definition mod_scalars.F:511
mp_exchange_mod
Definition mp_exchange.F:2
mp_exchange_mod::mp_exchange2d
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)
Definition mp_exchange.F:294