ROMS
Loading...
Searching...
No Matches
prsgrd31.h
Go to the documentation of this file.
1 MODULE prsgrd_mod
2!
3!git $Id$
4!=======================================================================
5! Copyright (c) 2002-2025 The ROMS Group !
6! Licensed under a MIT/X style license !
7! See License_ROMS.md Hernan G. Arango !
8!========================================== Alexander F. Shchepetkin ===
9! !
10! This subroutine evaluates the baroclinic hydrostatic pressure !
11! gradient term using the STANDARD density Jacobian or WEIGHTED !
12! density Jacobian scheme of Song (1998). Both of these approaches !
13! compute horizontal differences of density before of the vertical !
14! integration. !
15! !
16! The pressure gradient terms (m4/s2) are loaded into right-hand- !
17! side arrays "ru" and "rv". !
18! !
19! Reference: !
20! !
21! Song, Y.T., 1998: A general pressure gradient formulation for !
22! numerical ocean models. Part I: Scheme design and diagnostic !
23! analysis, Monthly Weather Rev., 126, 3213-3230. !
24! !
25!=======================================================================
26!
27 implicit none
28!
29 PRIVATE
30 PUBLIC :: prsgrd
31!
32 CONTAINS
33!
34!***********************************************************************
35 SUBROUTINE prsgrd (ng, tile)
36!***********************************************************************
37!
38 USE mod_param
39#ifdef DIAGNOSTICS
40 USE mod_diags
41#endif
42#ifdef ATM_PRESS
43 USE mod_forces
44#endif
45 USE mod_grid
46 USE mod_ocean
47 USE mod_stepping
48!
49! Imported variable declarations.
50!
51 integer, intent(in) :: ng, tile
52!
53! Local variable declarations.
54!
55 character (len=*), parameter :: myfile = &
56 & __FILE__
57!
58#include "tile.h"
59!
60#ifdef PROFILE
61 CALL wclock_on (ng, inlm, 23, __line__, myfile)
62#endif
63 CALL prsgrd31_tile (ng, tile, &
64 & lbi, ubi, lbj, ubj, &
65 & imins, imaxs, jmins, jmaxs, &
66 & nrhs(ng), &
67#ifdef WET_DRY
68 & grid(ng)%umask_wet, &
69 & grid(ng)%vmask_wet, &
70#endif
71 & grid(ng) % Hz, &
72 & grid(ng) % om_v, &
73 & grid(ng) % on_u, &
74 & grid(ng) % z_r, &
75 & grid(ng) % z_w, &
76 & ocean(ng) % rho, &
77#ifdef TIDE_GENERATING_FORCES
78 & ocean(ng) % eq_tide, &
79#endif
80#ifdef WEC_VF
81 & ocean(ng) % zetat, &
82#endif
83#ifdef ATM_PRESS
84 & forces(ng) % Pair, &
85#endif
86#ifdef DIAGNOSTICS_UV
87 & diags(ng) % DiaRU, &
88 & diags(ng) % DiaRV, &
89#endif
90 & ocean(ng) % ru, &
91 & ocean(ng) % rv)
92#ifdef PROFILE
93 CALL wclock_off (ng, inlm, 23, __line__, myfile)
94#endif
95!
96 RETURN
97 END SUBROUTINE prsgrd
98!
99!***********************************************************************
100 SUBROUTINE prsgrd31_tile (ng, tile, &
101 & LBi, UBi, LBj, UBj, &
102 & IminS, ImaxS, JminS, JmaxS, &
103 & nrhs, &
104#ifdef WET_DRY
105 & umask_wet, vmask_wet, &
106#endif
107 & Hz, om_v, on_u, z_r, z_w, &
108 & rho, &
109#ifdef TIDE_GENERATING_FORCES
110 & eq_tide, &
111#endif
112#ifdef WEC_VF
113 & zetat, &
114#endif
115#ifdef ATM_PRESS
116 & Pair, &
117#endif
118#ifdef DIAGNOSTICS_UV
119 & DiaRU, DiaRV, &
120#endif
121 & ru, rv)
122!***********************************************************************
123!
124 USE mod_param
125 USE mod_scalars
126!
127! Imported variable declarations.
128!
129 integer, intent(in) :: ng, tile
130 integer, intent(in) :: LBi, UBi, LBj, UBj
131 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
132 integer, intent(in) :: nrhs
133
134#ifdef ASSUMED_SHAPE
135# ifdef WET_DRY
136 real(r8), intent(in) :: umask_wet(LBi:,LBj:)
137 real(r8), intent(in) :: vmask_wet(LBi:,LBj:)
138# endif
139 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
140 real(r8), intent(in) :: om_v(LBi:,LBj:)
141 real(r8), intent(in) :: on_u(LBi:,LBj:)
142 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
143 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
144 real(r8), intent(in) :: rho(LBi:,LBj:,:)
145# ifdef TIDE_GENERATING_FORCES
146 real(r8), intent(in) :: eq_tide(LBi:,LBj:)
147# endif
148# ifdef WEC_VF
149 real(r8), intent(in) :: zetat(LBi:,LBj:)
150# endif
151# ifdef ATM_PRESS
152 real(r8), intent(in) :: Pair(LBi:,LBj:)
153# endif
154# ifdef DIAGNOSTICS_UV
155 real(r8), intent(inout) :: DiaRU(LBi:,LBj:,:,:,:)
156 real(r8), intent(inout) :: DiaRV(LBi:,LBj:,:,:,:)
157# endif
158 real(r8), intent(inout) :: ru(LBi:,LBj:,0:,:)
159 real(r8), intent(inout) :: rv(LBi:,LBj:,0:,:)
160#else
161# ifdef WET_DRY
162 real(r8), intent(in) :: umask_wet(LBi:UBi,LBj:UBj)
163 real(r8), intent(in) :: vmask_wet(LBi:UBi,LBj:UBj)
164# endif
165 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
166 real(r8), intent(in) :: om_v(LBi:UBi,LBj:UBj)
167 real(r8), intent(in) :: on_u(LBi:UBi,LBj:UBj)
168 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
169 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
170 real(r8), intent(in) :: rho(LBi:UBi,LBj:UBj,N(ng))
171# ifdef TIDE_GENERATING_FORCES
172 real(r8), intent(in) :: eq_tide(LBi:UBi,LBj:UBj)
173# endif
174# ifdef WEC_VF
175 real(r8), intent(in) :: zetat(LBi:UBi,LBj:UBj)
176# endif
177# ifdef ATM_PRESS
178 real(r8), intent(in) :: Pair(LBi:UBi,LBj:UBj)
179# endif
180# ifdef DIAGNOSTICS_UV
181 real(r8), intent(inout) :: DiaRU(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
182 real(r8), intent(inout) :: DiaRV(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
183# endif
184 real(r8), intent(inout) :: ru(LBi:UBi,LBj:UBj,0:N(ng),2)
185 real(r8), intent(inout) :: rv(LBi:UBi,LBj:UBj,0:N(ng),2)
186#endif
187!
188! Local variable declarations.
189!
190 integer :: i, j, k
191 real(r8) :: fac, fac1, fac2, fac3
192 real(r8) :: cff1, cff2, cff3, cff4
193#ifdef WJ_GRADP
194 real(r8) :: gamma
195#endif
196
197 real(r8), dimension(IminS:ImaxS) :: phie
198 real(r8), dimension(IminS:ImaxS) :: phix
199
200#include "set_bounds.h"
201!
202!-----------------------------------------------------------------------
203! Calculate pressure gradient in the XI-direction (m4/s2).
204!-----------------------------------------------------------------------
205!
206! Compute surface baroclinic pressure gradient.
207!
208#ifdef ATM_PRESS
209 fac=100.0_r8/rho0
210#endif
211 fac1=0.5_r8*g/rho0
212 fac2=1000.0_r8*g/rho0
213 fac3=0.25_r8*g/rho0
214
215 DO j=jstr,jend
216 DO i=istru,iend
217#ifdef TIDE_GENERATING_FORCES
218 cff1=z_w(i ,j,n(ng))-eq_tide(i ,j)-z_r(i ,j,n(ng))+ &
219 & z_w(i-1,j,n(ng))-eq_tide(i-1,j)-z_r(i-1,j,n(ng))
220#else
221 cff1=z_w(i ,j,n(ng))-z_r(i ,j,n(ng))+ &
222 & z_w(i-1,j,n(ng))-z_r(i-1,j,n(ng))
223#endif
224 phix(i)=fac1*(rho(i,j,n(ng))-rho(i-1,j,n(ng)))*cff1
225#ifdef WEC_VF
226 phix(i)=phix(i)+zetat(i,j)-zetat(i-1,j)
227#endif
228#ifdef ATM_PRESS
229 phix(i)=phix(i)+fac*(pair(i,j)-pair(i-1,j))
230#endif
231#ifdef RHO_SURF
232 phix(i)=phix(i)+ &
233 & (fac2+fac1*(rho(i,j,n(ng))+rho(i-1,j,n(ng))))* &
234 & (z_w(i,j,n(ng))-z_w(i-1,j,n(ng)))
235#endif
236 ru(i,j,n(ng),nrhs)=-0.5_r8*(hz(i,j,n(ng))+hz(i-1,j,n(ng)))* &
237 & phix(i)*on_u(i,j)
238#ifdef WET_DRY
239 ru(i,j,n(ng),nrhs)=ru(i,j,n(ng),nrhs)*umask_wet(i,j)
240#endif
241#ifdef DIAGNOSTICS_UV
242 diaru(i,j,n(ng),nrhs,m3pgrd)=ru(i,j,n(ng),nrhs)
243#endif
244 END DO
245!
246! Compute interior baroclinic pressure gradient. Differentiate and
247! then vertically integrate.
248!
249 DO k=n(ng)-1,1,-1
250 DO i=istru,iend
251#ifdef WJ_GRADP
252 cff1=1.0_r8/((z_r(i ,j,k+1)-z_r(i ,j,k))* &
253 & (z_r(i-1,j,k+1)-z_r(i-1,j,k)))
254 cff2=z_r(i ,j,k )-z_r(i-1,j,k )+ &
255 & z_r(i ,j,k+1)-z_r(i-1,j,k+1)
256 cff3=z_r(i ,j,k+1)-z_r(i ,j,k )- &
257 & z_r(i-1,j,k+1)+z_r(i-1,j,k )
258 gamma=0.125_r8*cff1*cff2*cff3
259
260 cff1=(1.0_r8+gamma)*(rho(i,j,k+1)-rho(i-1,j,k+1))+ &
261 & (1.0_r8-gamma)*(rho(i,j,k )-rho(i-1,j,k ))
262 cff2=rho(i,j,k+1)+rho(i-1,j,k+1)- &
263 & rho(i,j,k )-rho(i-1,j,k )
264 cff3=z_r(i,j,k+1)+z_r(i-1,j,k+1)- &
265 & z_r(i,j,k )-z_r(i-1,j,k )
266 cff4=(1.0_r8+gamma)*(z_r(i,j,k+1)-z_r(i-1,j,k+1))+ &
267 & (1.0_r8-gamma)*(z_r(i,j,k )-z_r(i-1,j,k ))
268 phix(i)=phix(i)+ &
269 & fac3*(cff1*cff3-cff2*cff4)
270#else
271 cff1=rho(i,j,k+1)-rho(i-1,j,k+1)+ &
272 & rho(i,j,k )-rho(i-1,j,k )
273 cff2=rho(i,j,k+1)+rho(i-1,j,k+1)- &
274 & rho(i,j,k )-rho(i-1,j,k )
275 cff3=z_r(i,j,k+1)+z_r(i-1,j,k+1)- &
276 & z_r(i,j,k )-z_r(i-1,j,k )
277 cff4=z_r(i,j,k+1)-z_r(i-1,j,k+1)+ &
278 & z_r(i,j,k )-z_r(i-1,j,k )
279 phix(i)=phix(i)+ &
280 & fac3*(cff1*cff3-cff2*cff4)
281#endif
282 ru(i,j,k,nrhs)=-0.5_r8*(hz(i,j,k)+hz(i-1,j,k))* &
283 & phix(i)*on_u(i,j)
284#ifdef WET_DRY
285 ru(i,j,k,nrhs)=ru(i,j,k,nrhs)*umask_wet(i,j)
286#endif
287#ifdef DIAGNOSTICS_UV
288 diaru(i,j,k,nrhs,m3pgrd)=ru(i,j,k,nrhs)
289#endif
290 END DO
291 END DO
292!
293!-----------------------------------------------------------------------
294! Calculate pressure gradient in the ETA-direction (m4/s2).
295!-----------------------------------------------------------------------
296!
297! Compute surface baroclinic pressure gradient.
298!
299 IF (j.ge.jstrv) THEN
300 DO i=istr,iend
301#ifdef TIDE_GENERATING_FORCES
302 cff1=z_w(i,j ,n(ng))-eq_tide(i,j )-z_r(i,j ,n(ng))+ &
303 & z_w(i,j-1,n(ng))-eq_tide(i,j-1)-z_r(i,j-1,n(ng))
304#else
305 cff1=z_w(i,j ,n(ng))-z_r(i,j ,n(ng))+ &
306 & z_w(i,j-1,n(ng))-z_r(i,j-1,n(ng))
307#endif
308 phie(i)=fac1*(rho(i,j,n(ng))-rho(i,j-1,n(ng)))*cff1
309#ifdef WEC_VF
310 phie(i)=phie(i)+zetat(i,j)-zetat(i,j-1)
311#endif
312#ifdef ATM_PRESS
313 phie(i)=phie(i)+fac*(pair(i,j)-pair(i,j-1))
314#endif
315#ifdef RHO_SURF
316 phie(i)=phie(i)+ &
317 & (fac2+fac1*(rho(i,j,n(ng))+rho(i,j-1,n(ng))))* &
318 & (z_w(i,j,n(ng))-z_w(i,j-1,n(ng)))
319#endif
320 rv(i,j,n(ng),nrhs)=-0.5_r8*(hz(i,j,n(ng))+hz(i,j-1,n(ng)))* &
321 & phie(i)*om_v(i,j)
322#ifdef WET_DRY
323 rv(i,j,n(ng),nrhs)=rv(i,j,n(ng),nrhs)*vmask_wet(i,j)
324#endif
325#ifdef DIAGNOSTICS_UV
326 diarv(i,j,n(ng),nrhs,m3pgrd)=rv(i,j,n(ng),nrhs)
327#endif
328 END DO
329!
330! Compute interior baroclinic pressure gradient. Differentiate and
331! then vertically integrate.
332!
333 DO k=n(ng)-1,1,-1
334 DO i=istr,iend
335#ifdef WJ_GRADP
336 cff1=1.0_r8/((z_r(i,j ,k+1)-z_r(i,j ,k))* &
337 & (z_r(i,j-1,k+1)-z_r(i,j-1,k)))
338 cff2=z_r(i,j ,k )-z_r(i,j-1,k )+ &
339 & z_r(i,j ,k+1)-z_r(i,j-1,k+1)
340 cff3=z_r(i,j ,k+1)-z_r(i,j ,k )- &
341 & z_r(i,j-1,k+1)+z_r(i,j-1,k )
342 gamma=0.125_r8*cff1*cff2*cff3
343
344 cff1=(1.0_r8+gamma)*(rho(i,j,k+1)-rho(i,j-1,k+1))+ &
345 & (1.0_r8-gamma)*(rho(i,j,k )-rho(i,j-1,k ))
346 cff2=rho(i,j,k+1)+rho(i,j-1,k+1)- &
347 & rho(i,j,k )-rho(i,j-1,k )
348 cff3=z_r(i,j,k+1)+z_r(i,j-1,k+1)- &
349 & z_r(i,j,k )-z_r(i,j-1,k )
350 cff4=(1.0_r8+gamma)*(z_r(i,j,k+1)-z_r(i,j-1,k+1))+ &
351 & (1.0_r8-gamma)*(z_r(i,j,k )-z_r(i,j-1,k ))
352 phie(i)=phie(i)+ &
353 & fac3*(cff1*cff3-cff2*cff4)
354#else
355 cff1=rho(i,j,k+1)-rho(i,j-1,k+1)+ &
356 & rho(i,j,k )-rho(i,j-1,k )
357 cff2=rho(i,j,k+1)+rho(i,j-1,k+1)- &
358 & rho(i,j,k )-rho(i,j-1,k )
359 cff3=z_r(i,j,k+1)+z_r(i,j-1,k+1)- &
360 & z_r(i,j,k )-z_r(i,j-1,k )
361 cff4=z_r(i,j,k+1)-z_r(i,j-1,k+1)+ &
362 & z_r(i,j,k )-z_r(i,j-1,k )
363 phie(i)=phie(i)+ &
364 & fac3*(cff1*cff3-cff2*cff4)
365#endif
366 rv(i,j,k,nrhs)=-0.5_r8*(hz(i,j,k)+hz(i,j-1,k))* &
367 & phie(i)*om_v(i,j)
368#ifdef WET_DRY
369 rv(i,j,k,nrhs)=rv(i,j,k,nrhs)*vmask_wet(i,j)
370#endif
371#ifdef DIAGNOSTICS_UV
372 diarv(i,j,k,nrhs,m3pgrd)=rv(i,j,k,nrhs)
373#endif
374 END DO
375 END DO
376 END IF
377 END DO
378!
379 RETURN
380 END SUBROUTINE prsgrd31_tile
381
382 END MODULE prsgrd_mod
mod_diags
Definition mod_diags.F:2
mod_diags::diags
type(t_diags), dimension(:), allocatable diags
Definition mod_diags.F:100
mod_forces
Definition mod_forces.F:2
mod_forces::forces
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
mod_grid
Definition mod_grid.F:2
mod_grid::grid
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
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_scalars
Definition mod_scalars.F:2
mod_scalars::g
real(dp) g
Definition mod_scalars.F:464
mod_scalars::rho0
real(dp) rho0
Definition mod_scalars.F:787
mod_scalars::m3pgrd
integer m3pgrd
Definition mod_scalars.F:190
mod_stepping
Definition mod_stepping.F:2
mod_stepping::nrhs
integer, dimension(:), allocatable nrhs
Definition mod_stepping.F:70
prsgrd_mod
Definition prsgrd31.h:1
prsgrd_mod::prsgrd31_tile
subroutine prsgrd31_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, umask_wet, vmask_wet, hz, om_v, on_u, z_r, z_w, rho, eq_tide, zetat, pair, diaru, diarv, ru, rv)
Definition prsgrd31.h:122
prsgrd_mod::prsgrd
subroutine, public prsgrd(ng, tile)
Definition prsgrd31.h:36
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