ROMS
Loading...
Searching...
No Matches
vorticity.F
Go to the documentation of this file.
1#include "cppdefs.h"
2 MODULE vorticity_mod
3
4#if defined AVERAGES || \
5 (defined ad_averages && defined adjoint) || \
6 (defined rp_averages && defined tl_ioms) || \
7 (defined tl_averages && defined tangent)
8!
9!git $Id$
10!=======================================================================
11! Copyright (c) 2002-2025 The ROMS Group !
12! Licensed under a MIT/X style license !
13! See License_ROMS.md !
14!================================================== Hernan G. Arango ===
15! !
16! This routine computes relative (s-1) and potential (m-1 s-1) !
17! vorticity for an adiabatic Boussinesq fluid where the potential !
18! density is conserved: !
19! !
20! pvor = 1/rho0 dot_product(avor, grad(pden)) !
21! !
22! where "avor" is the absolute (relative plus planetary) vorticity !
23! and "pden" is the potential density (a conserved quantity). !
24! !
25! avor = rvor + f !
26! !
27! In curvilinear coordinates, the vertical component of relative !
28! vorticity and potential vorticity are: !
29! are: !
30! !
31! rvor = mn * [d(v/n)/d(xi) - d(u/m)/d(eta)] !
32! !
33! pvor = mn/rho0 * [f/mn + !
34! !
35! d(v/n)/d(xi) - d(u/m)/d(eta)] * d(pden)/d(z) + !
36! !
37! 1/rho0 * [1/n d(pden)/d(eta) d(u)/d(z) - !
38! !
39! 1/m d(pden)/d(xi) d(v)/d(z)] !
40! !
41! In addition, the vertically integrated (shallow water) relative !
42! and potential vorticity are computed. !
43! !
44! The relative and potential vorticity is discretized at horizontal !
45! PSI-points and vertical RHO-points. !
46! !
47!=======================================================================
48!
49 implicit none
50!
51 PRIVATE
52 PUBLIC :: vorticity
53 PUBLIC :: vorticity_tile
54!
55 CONTAINS
56!
57!***********************************************************************
58 SUBROUTINE vorticity (ng, tile)
59!***********************************************************************
60!
61 USE mod_param
62 USE mod_average
63 USE mod_grid
64 USE mod_ocean
65 USE mod_stepping
66!
67! Imported variable declarations.
68!
69 integer, intent(in) :: ng, tile
70!
71! Local variable declarations.
72!
73 character (len=*), parameter :: myfile = &
74 & __FILE__
75!
76# include "tile.h"
77!
78# ifdef PROFILE
79 CALL wclock_on (ng, inlm, 5, __line__, myfile)
80# endif
81 CALL vorticity_tile (ng, tile, &
82 & lbi, ubi, lbj, ubj, &
83 & imins, imaxs, jmins, jmaxs, &
84# ifdef SOLVE3D
85 & kout, nout, &
86# else
87 & kout, &
88# endif
89# ifdef MASKING
90 & grid(ng) % pmask, &
91 & grid(ng) % umask, &
92 & grid(ng) % vmask, &
93# endif
94 & grid(ng) % f, &
95 & grid(ng) % h, &
96 & grid(ng) % om_u, &
97 & grid(ng) % on_v, &
98 & grid(ng) % pm, &
99 & grid(ng) % pn, &
100# ifdef SOLVE3D
101 & grid(ng) % z_r, &
102 & ocean(ng) % pden, &
103 & ocean(ng) % u, &
104 & ocean(ng) % v, &
105# endif
106 & ocean(ng) % ubar, &
107 & ocean(ng) % vbar, &
108 & ocean(ng) % zeta, &
109# ifdef SOLVE3D
110 & average(ng) % avgpvor3d, &
111 & average(ng) % avgrvor3d, &
112# endif
113 & average(ng) % avgpvor2d, &
114 & average(ng) % avgrvor2d)
115
116# ifdef PROFILE
117 CALL wclock_off (ng, inlm, 5, __line__, myfile)
118# endif
119!
120 RETURN
121 END SUBROUTINE vorticity
122!
123!***********************************************************************
124 SUBROUTINE vorticity_tile (ng, tile, &
125 & LBi, UBi, LBj, UBj, &
126 & IminS, ImaxS, JminS, JmaxS, &
127# ifdef SOLVE3D
128 & kout, nout, &
129# else
130 & kout, &
131# endif
132# ifdef MASKING
133 & pmask, umask, vmask, &
134# endif
135 & f, h, om_u, on_v, pm, pn, &
136# ifdef SOLVE3D
137 & z_r, pden, u, v, &
138# endif
139 & ubar, vbar, zeta, &
140# ifdef SOLVE3D
141 & pvor, rvor, &
142# endif
143 & pvor_bar, rvor_bar)
144!***********************************************************************
145!
146 USE mod_param
147 USE mod_ncparam
148 USE mod_scalars
149!
150 USE exchange_2d_mod, ONLY : exchange_p2d_tile
151# ifdef SOLVE3D
152 USE exchange_3d_mod, ONLY : exchange_p3d_tile
153# endif
154# ifdef DISTRIBUTE
155 USE mp_exchange_mod, ONLY : mp_exchange2d
156# ifdef SOLVE3D
157 USE mp_exchange_mod, ONLY : mp_exchange3d
158# endif
159# endif
160!
161! Imported variable declarations.
162!
163 integer, intent(in) :: ng, tile
164 integer, intent(in) :: lbi, ubi, lbj, ubj
165 integer, intent(in) :: imins, imaxs, jmins, jmaxs
166# ifdef SOLVE3D
167 integer, intent(in) :: kout, nout
168# else
169 integer, intent(in) :: kout
170# endif
171
172# ifdef ASSUMED_SHAPE
173# ifdef MASKING
174 real(r8), intent(in) :: pmask(lbi:,lbj:)
175 real(r8), intent(in) :: umask(lbi:,lbj:)
176 real(r8), intent(in) :: vmask(lbi:,lbj:)
177# endif
178 real(r8), intent(in) :: f(lbi:,lbj:)
179 real(r8), intent(in) :: h(lbi:,lbj:)
180 real(r8), intent(in) :: om_u(lbi:,lbj:)
181 real(r8), intent(in) :: on_v(lbi:,lbj:)
182 real(r8), intent(in) :: pm(lbi:,lbj:)
183 real(r8), intent(in) :: pn(lbi:,lbj:)
184# ifdef SOLVE3D
185 real(r8), intent(in) :: z_r(lbi:,lbj:,:)
186 real(r8), intent(in) :: pden(lbi:,lbj:,:)
187 real(r8), intent(in) :: u(lbi:,lbj:,:,:)
188 real(r8), intent(in) :: v(lbi:,lbj:,:,:)
189# endif
190 real(r8), intent(in) :: ubar(lbi:,lbj:,:)
191 real(r8), intent(in) :: vbar(lbi:,lbj:,:)
192 real(r8), intent(in) :: zeta(lbi:,lbj:,:)
193
194 real(r8), intent(out) :: pvor_bar(lbi:,lbj:)
195 real(r8), intent(out) :: rvor_bar(lbi:,lbj:)
196# ifdef SOLVE3D
197 real(r8), intent(out) :: pvor(lbi:,lbj:,:)
198 real(r8), intent(out) :: rvor(lbi:,lbj:,:)
199# endif
200
201# else
202
203# ifdef MASKING
204 real(r8), intent(in) :: pmask(lbi:ubi,lbj:ubj)
205 real(r8), intent(in) :: umask(lbi:ubi,lbj:ubj)
206 real(r8), intent(in) :: vmask(lbi:ubi,lbj:ubj)
207# endif
208 real(r8), intent(in) :: f(lbi:ubi,lbj:ubj)
209 real(r8), intent(in) :: h(lbi:ubi,lbj:ubj)
210 real(r8), intent(in) :: om_u(lbi:ubi,lbj:ubj)
211 real(r8), intent(in) :: on_v(lbi:ubi,lbj:ubj)
212 real(r8), intent(in) :: pm(lbi:ubi,lbj:ubj)
213 real(r8), intent(in) :: pn(lbi:ubi,lbj:ubj)
214# ifdef SOLVE3D
215 real(r8), intent(in) :: z_r(lbi:ubi,lbj:ubj,n(ng))
216 real(r8), intent(in) :: pden(lbi:ubi,lbj:ubj,n(ng))
217 real(r8), intent(in) :: u(lbi:ubi,lbj:ubj,n(ng),2)
218 real(r8), intent(in) :: v(lbi:ubi,lbj:ubj,n(ng),2)
219# endif
220 real(r8), intent(in) :: ubar(lbi:ubi,lbj:ubj,:)
221 real(r8), intent(in) :: vbar(lbi:ubi,lbj:ubj,:)
222 real(r8), intent(in) :: zeta(lbi:ubi,lbj:ubj,:)
223
224 real(r8), intent(out) :: pvor_bar(lbi:ubi,lbj:ubj)
225 real(r8), intent(out) :: rvor_bar(lbi:ubi,lbj:ubj)
226# ifdef SOLVE3D
227 real(r8), intent(out) :: pvor(lbi:ubi,lbj:ubj,n(ng))
228 real(r8), intent(out) :: rvor(lbi:ubi,lbj:ubj,n(ng))
229# endif
230# endif
231!
232! Local variable declarations.
233!
234 integer :: i, j
235# ifdef SOLVE3D
236 integer :: k, k1, k2
237# endif
238 real(r8) :: cff
239 real(r8) :: dvdx_p, dude_p, fomn_p
240# ifdef SOLVE3D
241 real(r8) :: drde_pr, drdx_pr, drdz_pr, dudz_pr, dvdz_pr
242 real(r8) :: orho0
243
244 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: drde
245 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: drdx
246 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: dude
247 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: dvdx
248
249 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: drdz
250 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dudz
251 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dvdz
252# endif
253
254# include "set_bounds.h"
255
256# ifdef SOLVE3D
257!
258!-----------------------------------------------------------------------
259! Compute 3D relative and potential vorticity.
260!-----------------------------------------------------------------------
261!
262! Compute horizontal and vertical gradients. Notice the recursive
263! blocking sequence for vertical placement of the gradients is:
264!
265! dRdz,dUdz,dVdz(:,:,k1) k-1/2 W-points
266! dRdz,dUdz,dVdz(:,:,k2) k+1/2 W-points
267!
268 orho0=1.0_r8/rho0
269
270 k2=1
271 k_loop : DO k=0,n(ng)
272 k1=k2
273 k2=3-k1
274 IF (k.gt.0) THEN
275 DO j=jstr-1,jendr
276 DO i=istr,iendr
277 cff=0.5_r8*(pm(i,j)+pm(i-1,j))
278# ifdef MASKING
279 cff=cff*umask(i,j)
280# endif
281 drdx(i,j)=cff*(pden(i ,j,k)- &
282 & pden(i-1,j,k))
283 END DO
284 END DO
285 DO j=jstr,jendr
286 DO i=istr-1,iendr
287 cff=0.5_r8*(pn(i,j)+pn(i,j-1))
288# ifdef MASKING
289 cff=cff*vmask(i,j)
290# endif
291 drde(i,j)=cff*(pden(i,j ,k)- &
292 & pden(i,j-1,k))
293 END DO
294 END DO
295 DO j=jstr,jendr
296 DO i=istr,iendr
297 dude(i,j)=om_u(i,j )*u(i,j ,k,nout)- &
298 & om_u(i,j-1)*u(i,j-1,k,nout)
299# ifdef MASKING
300 dude(i,j)=dude(i,j)*pmask(i,j)
301# endif
302 dvdx(i,j)=on_v(i ,j)*v(i ,j,k,nout)- &
303 & on_v(i-1,j)*v(i-1,j,k,nout)
304# ifdef MASKING
305 dvdx(i,j)=dvdx(i,j)*pmask(i,j)
306# endif
307 END DO
308 END DO
309 END IF
310 IF ((k.eq.0).or.(k.eq.n(ng))) THEN
311 DO j=jstr-1,jendr
312 DO i=istr-1,iendr
313 drdz(i,j,k2)=0.0_r8
314 END DO
315 END DO
316 DO j=jstr-1,jendr
317 DO i=istr,iendr
318 dudz(i,j,k2)=0.0_r8
319 END DO
320 END DO
321 DO j=jstr,jendr
322 DO i=istr-1,iendr
323 dvdz(i,j,k2)=0.0_r8
324 END DO
325 END DO
326 ELSE
327 DO j=jstr-1,jendr
328 DO i=istr-1,iendr
329 cff=1.0_r8/(z_r(i,j,k+1)-z_r(i,j,k))
330 drdz(i,j,k2)=cff*(pden(i,j,k+1)- &
331 & pden(i,j,k ))
332 END DO
333 END DO
334 DO j=jstr-1,jendr
335 DO i=istr,iendr
336 cff=1.0_r8/(0.5_r8*(z_r(i-1,j,k+1)-z_r(i-1,j,k)+ &
337 & z_r(i ,j,k+1)-z_r(i ,j,k)))
338 dudz(i,j,k2)=cff*(u(i,j,k+1,nout)- &
339 & u(i,j,k ,nout))
340 END DO
341 END DO
342 DO j=jstr,jendr
343 DO i=istr-1,iendr
344 cff=1.0_r8/(0.5_r8*(z_r(i,j-1,k+1)-z_r(i,j-1,k)+ &
345 & z_r(i,j ,k+1)-z_r(i,j ,k)))
346 dvdz(i,j,k2)=cff*(v(i,j,k+1,nout)- &
347 & v(i,j,k ,nout))
348 END DO
349 END DO
350 END IF
351!
352! Compute relative vorticity (second-1) and potential vorticity
353! (meter-1 second-1) at horizontal PSI-points and vertical RHO-points.
354!
355 IF (k.gt.0) THEN
356 DO j=jstr,jendr
357 DO i=istr,iendr
358 cff=0.0625_r8* &
359 & (pm(i-1,j-1)+pm(i-1,j)+pm(i,j-1)+pm(i,j))* &
360 & (pn(i-1,j-1)+pn(i-1,j)+pn(i,j-1)+pn(i,j))
361 fomn_p=0.25_r8*(f(i-1,j-1)+f(i-1,j)+f(i,j-1)+f(i,j))/cff
362 drde_pr=drde(i-1,j )+drde(i,j)
363 drdx_pr=drdx(i ,j-1)+drdx(i,j)
364 drdz_pr=0.125_r8*(drdz(i-1,j-1,k1)+drdz(i-1,j-1,k2)+ &
365 & drdz(i ,j-1,k1)+drdz(i ,j-1,k2)+ &
366 & drdz(i-1,j ,k1)+drdz(i-1,j ,k2)+ &
367 & drdz(i ,j ,k1)+drdz(i ,j ,k2))
368 dudz_pr=dudz(i ,j-1,k1)+dudz(i ,j-1,k2)+ &
369 & dudz(i ,j ,k1)+dudz(i ,j ,k2)
370 dvdz_pr=dvdz(i-1,j ,k1)+dvdz(i-1,j ,k2)+ &
371 & dvdz(i ,j ,k1)+dvdz(i ,j ,k2)
372 rvor(i,j,k)=cff*(dvdx(i,j)-dude(i,j))
373 pvor(i,j,k)=orho0* &
374 & (cff*drdz_pr*(fomn_p+ &
375 & dvdx(i,j)-dude(i,j))+ &
376 & 0.125_r8*(dudz_pr*drde_pr-dvdz_pr*drdx_pr))
377 END DO
378 END DO
379 END IF
380 END DO k_loop
381!
382! Exchange boundary data.
383!
384 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
385 CALL exchange_p3d_tile (ng, tile, &
386 & lbi, ubi, lbj, ubj, 1, n(ng), &
387 & pvor)
388 CALL exchange_p3d_tile (ng, tile, &
389 & lbi, ubi, lbj, ubj, 1, n(ng), &
390 & rvor)
391 END IF
392
393# ifdef DISTRIBUTE
394 CALL mp_exchange3d (ng, tile, inlm, 2, &
395 & lbi, ubi, lbj, ubj, 1, n(ng), &
396 & nghostpoints, ewperiodic(ng), nsperiodic(ng), &
397 & pvor, &
398 & rvor)
399# endif
400# endif
401!
402!-----------------------------------------------------------------------
403! Compute 2D relative and potential vorticity.
404!-----------------------------------------------------------------------
405!
406! Compute vertically-integrated relative vorticity (second-1) and
407! potential vorticity (meter-1 second-1) at PSI-points.
408!
409 DO j=jstr,jendr
410 DO i=istr,iendr
411 cff=0.0625_r8* &
412 & (pm(i-1,j-1)+pm(i-1,j)+pm(i,j-1)+pm(i,j))* &
413 & (pn(i-1,j-1)+pn(i-1,j)+pn(i,j-1)+pn(i,j))
414 fomn_p=0.25_r8*(f(i-1,j-1)+f(i-1,j)+f(i,j-1)+f(i,j))/cff
415 cff=pm(i,j)*pn(i,j)
416 dvdx_p=on_v(i ,j)*vbar(i ,j,kout)- &
417 & on_v(i-1,j)*vbar(i-1,j,kout)
418# ifdef MASKING
419 dvdx_p=dvdx_p*pmask(i,j)
420# endif
421 dude_p=om_u(i,j )*ubar(i,j ,kout)- &
422 & om_u(i,j-1)*ubar(i,j-1,kout)
423# ifdef MASKING
424 dude_p=dude_p*pmask(i,j)
425# endif
426 rvor_bar(i,j)=cff*(dvdx_p-dude_p)
427 pvor_bar(i,j)=cff*((fomn_p+dvdx_p-dude_p)/ &
428 & (h(i,j)+zeta(i,j,kout)))
429 END DO
430 END DO
431!
432! Exchange boundary information.
433!
434 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
435 CALL exchange_p2d_tile (ng, tile, &
436 & lbi, ubi, lbj, ubj, &
437 & pvor_bar)
438 CALL exchange_p2d_tile (ng, tile, &
439 & lbi, ubi, lbj, ubj, &
440 & rvor_bar)
441 END IF
442
443# ifdef DISTRIBUTE
444 CALL mp_exchange2d (ng, tile, inlm, 2, &
445 & lbi, ubi, lbj, ubj, &
446 & nghostpoints, ewperiodic(ng), nsperiodic(ng), &
447 & pvor_bar, &
448 & rvor_bar)
449# endif
450!
451 RETURN
452 END SUBROUTINE vorticity_tile
453#endif
454 END MODULE vorticity_mod
exchange_2d_mod
Definition exchange_2d.F:2
exchange_2d_mod::exchange_p2d_tile
subroutine exchange_p2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
Definition exchange_2d.F:66
exchange_3d_mod
Definition exchange_3d.F:2
exchange_3d_mod::exchange_p3d_tile
subroutine exchange_p3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
Definition exchange_3d.F:70
mod_average
Definition mod_average.F:2
mod_average::average
type(t_average), dimension(:), allocatable average
Definition mod_average.F:351
mod_grid
Definition mod_grid.F:2
mod_grid::grid
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
mod_ncparam
Definition mod_ncparam.F:2
mod_ocean
Definition mod_ocean.F:2
mod_ocean::ocean
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
mod_param
Definition mod_param.F:2
mod_param::inlm
integer, parameter inlm
Definition mod_param.F:662
mod_param::n
integer, dimension(:), allocatable n
Definition mod_param.F:479
mod_param::nghostpoints
integer nghostpoints
Definition mod_param.F:710
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
mod_scalars::rho0
real(dp) rho0
Definition mod_scalars.F:787
mod_stepping
Definition mod_stepping.F:2
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
mp_exchange_mod::mp_exchange3d
subroutine mp_exchange3d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, nghost, ew_periodic, ns_periodic, a, b, c, d)
Definition mp_exchange.F:2037
vorticity_mod
Definition vorticity.F:2
vorticity_mod::vorticity
subroutine, public vorticity(ng, tile)
Definition vorticity.F:59
vorticity_mod::vorticity_tile
subroutine, public vorticity_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, kout, nout, pmask, umask, vmask, f, h, om_u, on_v, pm, pn, z_r, pden, u, v, ubar, vbar, zeta, pvor, rvor, pvor_bar, rvor_bar)
Definition vorticity.F:144
wclock_off
recursive subroutine wclock_off(ng, model, region, line, routine)
Definition timers.F:148
wclock_on
recursive subroutine wclock_on(ng, model, region, line, routine)
Definition timers.F:3