ROMS
Loading...
Searching...
No Matches
prsgrd40.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 nonlinear baroclinic, hydrostatic !
11! pressure gradient term using the finite-volume pressure Jacobian !
12! scheme of Lin (1997). !
13! !
14! The pressure gradient terms (m4/s2) are loaded into right-hand- !
15! side arrays "ru" and "rv". !
16! !
17! Reference: !
18! !
19! Lin, Shian-Jiann, 1997: A finite volume integration method !
20! for computing pressure gradient force in general vertical !
21! coordinates, Q. J. R. Meteorol. Soc., 123, 1749-1762. !
22! !
23!=======================================================================
24!
25 implicit none
26!
27 PRIVATE
28 PUBLIC :: prsgrd
29!
30 CONTAINS
31!
32!***********************************************************************
33 SUBROUTINE prsgrd (ng, tile)
34!***********************************************************************
35!
36 USE mod_param
37#ifdef DIAGNOSTICS
38 USE mod_diags
39#endif
40#ifdef ATM_PRESS
41 USE mod_forces
42#endif
43 USE mod_grid
44 USE mod_ocean
45 USE mod_stepping
46!
47! Imported variable declarations.
48!
49 integer, intent(in) :: ng, tile
50!
51! Local variable declarations.
52!
53 character (len=*), parameter :: MyFile = &
54 & __FILE__
55!
56#include "tile.h"
57!
58#ifdef PROFILE
59 CALL wclock_on (ng, inlm, 23, __line__, myfile)
60#endif
61 CALL prsgrd40_tile (ng, tile, &
62 & lbi, ubi, lbj, ubj, &
63 & imins, imaxs, jmins, jmaxs, &
64 & nrhs(ng), &
65#ifdef WET_DRY
66 & grid(ng)%umask_wet, &
67 & grid(ng)%vmask_wet, &
68#endif
69 & grid(ng) % om_v, &
70 & grid(ng) % on_u, &
71 & grid(ng) % Hz, &
72 & grid(ng) % z_w, &
73 & ocean(ng) % rho, &
74#ifdef TIDE_GENERATING_FORCES
75 & ocean(ng) % eq_tide, &
76#endif
77#ifdef WEC_VF
78 & ocean(ng) % zetat, &
79#endif
80#ifdef ATM_PRESS
81 & forces(ng) % Pair, &
82#endif
83#ifdef DIAGNOSTICS_UV
84 & diags(ng) % DiaRU, &
85 & diags(ng) % DiaRV, &
86#endif
87 & ocean(ng) % ru, &
88 & ocean(ng) % rv)
89#ifdef PROFILE
90 CALL wclock_off (ng, inlm, 23, __line__, myfile)
91#endif
92!
93 RETURN
94 END SUBROUTINE prsgrd
95!
96!***********************************************************************
97 SUBROUTINE prsgrd40_tile (ng, tile, &
98 & LBi, UBi, LBj, UBj, &
99 & IminS, ImaxS, JminS, JmaxS, &
100 & nrhs, &
101#ifdef WET_DRY
102 & umask_wet, vmask_wet, &
103#endif
104 & om_v, on_u, &
105 & Hz, z_w, &
106 & rho, &
107#ifdef TIDE_GENERATING_FORCES
108 & eq_tide, &
109#endif
110#ifdef WEC_VF
111 & zetat, &
112#endif
113#ifdef ATM_PRESS
114 & Pair, &
115#endif
116#ifdef DIAGNOSTICS_UV
117 & DiaRU, DiaRV, &
118#endif
119 & ru, rv)
120!***********************************************************************
121!
122 USE mod_param
123 USE mod_scalars
124!
125! Imported variable declarations.
126!
127 integer, intent(in) :: ng, tile
128 integer, intent(in) :: LBi, UBi, LBj, UBj
129 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
130 integer, intent(in) :: nrhs
131
132#ifdef ASSUMED_SHAPE
133 real(r8), intent(in) :: om_v(LBi:,LBj:)
134 real(r8), intent(in) :: on_u(LBi:,LBj:)
135 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
136 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
137 real(r8), intent(in) :: rho(LBi:,LBj:,:)
138# ifdef TIDE_GENERATING_FORCES
139 real(r8), intent(in) :: eq_tide(LBi:,LBj:)
140# endif
141# ifdef WEC_VF
142 real(r8), intent(in) :: zetat(LBi:,LBj:)
143# endif
144# ifdef ATM_PRESS
145 real(r8), intent(in) :: Pair(LBi:,LBj:)
146# endif
147# ifdef DIAGNOSTICS_UV
148 real(r8), intent(inout) :: DiaRU(LBi:,LBj:,:,:,:)
149 real(r8), intent(inout) :: DiaRV(LBi:,LBj:,:,:,:)
150# endif
151 real(r8), intent(inout) :: ru(LBi:,LBj:,0:,:)
152 real(r8), intent(inout) :: rv(LBi:,LBj:,0:,:)
153#else
154 real(r8), intent(in) :: om_v(LBi:UBi,LBj:UBj)
155 real(r8), intent(in) :: on_u(LBi:UBi,LBj:UBj)
156 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
157 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
158 real(r8), intent(in) :: rho(LBi:UBi,LBj:UBj,N(ng))
159# ifdef TIDE_GENERATING_FORCES
160 real(r8), intent(in) :: eq_tide(LBi:UBi,LBj:UBj)
161# endif
162# ifdef WEC_VF
163 real(r8), intent(in) :: zetat(LBi:UBi,LBj:UBj)
164# endif
165# ifdef ATM_PRESS
166 real(r8), intent(in) :: Pair(LBi:UBi,LBj:UBj)
167# endif
168# ifdef DIAGNOSTICS_UV
169 real(r8), intent(inout) :: DiaRU(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
170 real(r8), intent(inout) :: DiaRV(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
171# endif
172 real(r8), intent(inout) :: ru(LBi:UBi,LBj:UBj,0:N(ng),2)
173 real(r8), intent(inout) :: rv(LBi:UBi,LBj:UBj,0:N(ng),2)
174#endif
175!
176! Local variable declarations.
177!
178 integer :: i, j, k
179
180 real(r8) :: cff, cff1, dh
181#ifdef ATM_PRESS
182 real(r8) :: OneAtm, fac
183#endif
184#ifdef WEC_VF
185 real(r8) :: fac1
186#endif
187 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: FC
188
189 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,N(ng)) :: FX
190
191 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,0:N(ng)) :: P
192
193#include "set_bounds.h"
194!
195!-----------------------------------------------------------------------
196! Finite Volume pressure gradient algorithm (Lin, 1997).
197!-----------------------------------------------------------------------
198!
199! Compute pressure and its vertical integral. Initialize pressure at
200! the free-surface as zero.
201!
202#ifdef ATM_PRESS
203 oneatm=1013.25_r8 ! 1 atm = 1013.25 mb
204 fac=100.0_r8/g
205#endif
206#ifdef WEC_VF
207 fac1=rho0/g
208#endif
209 j_loop : DO j=jstrv-1,jend
210 DO i=istru-1,iend
211 p(i,j,n(ng))=0.0_r8
212#ifdef WEC_VF
213 p(i,j,n(ng))=p(i,j,n(ng))+fac1*zetat(i,j)
214#endif
215#ifdef ATM_PRESS
216 p(i,j,n(ng))=p(i,j,n(ng))+fac*(pair(i,j)-oneatm)
217#endif
218#ifdef TIDE_GENERATING_FORCES
219 p(i,j,n(ng))=p(i,j,n(ng))-g*eq_tide(i,j)
220#endif
221 END DO
222 DO k=n(ng),1,-1
223 DO i=istru-1,iend
224 p(i,j,k-1)=p(i,j,k)+ &
225 & hz(i,j,k)*rho(i,j,k)
226 fx(i,j,k)=0.5_r8*hz(i,j,k)*(p(i,j,k)+p(i,j,k-1))
227 END DO
228 END DO
229!
230! Calculate pressure gradient in the XI-direction (m4/s2).
231!
232 IF (j.ge.jstr) THEN
233 DO i=istru,iend
234 fc(i,n(ng))=0.0_r8
235 END DO
236 cff=0.5_r8*g
237 cff1=g/rho0
238 DO k=n(ng),1,-1
239 DO i=istru,iend
240 dh=z_w(i,j,k-1)-z_w(i-1,j,k-1)
241 fc(i,k-1)=0.5_r8*dh*(p(i,j,k-1)+p(i-1,j,k-1))
242 ru(i,j,k,nrhs)=(cff*(hz(i-1,j,k)+ &
243 & hz(i ,j,k))* &
244 & (z_w(i-1,j,n(ng))- &
245 & z_w(i ,j,n(ng)))+ &
246 & cff1*(fx(i-1,j,k)- &
247 & fx(i ,j,k)+ &
248 & fc(i,k )- &
249 & fc(i,k-1)))*on_u(i,j)
250#ifdef WET_DRY
251 ru(i,j,k,nrhs)=ru(i,j,k,nrhs)*umask_wet(i,j)
252#endif
253#ifdef DIAGNOSTICS_UV
254 diaru(i,j,k,nrhs,m3pgrd)=ru(i,j,k,nrhs)
255#endif
256 END DO
257 END DO
258 END IF
259!
260! Calculate pressure gradient in the ETA-direction (m4/s2).
261!
262 IF (j.ge.jstrv) THEN
263 DO i=istr,iend
264 fc(i,n(ng))=0.0_r8
265 END DO
266 cff=0.5_r8*g
267 cff1=g/rho0
268 DO k=n(ng),1,-1
269 DO i=istr,iend
270 dh=z_w(i,j,k-1)-z_w(i,j-1,k-1)
271 fc(i,k-1)=0.5_r8*dh*(p(i,j,k-1)+p(i,j-1,k-1))
272 rv(i,j,k,nrhs)=(cff*(hz(i,j-1,k)+ &
273 & hz(i,j ,k))* &
274 & (z_w(i,j-1,n(ng))- &
275 & z_w(i,j ,n(ng)))+ &
276 & cff1*(fx(i,j-1,k)- &
277 & fx(i,j ,k)+ &
278 & fc(i,k )- &
279 & fc(i,k-1)))*om_v(i,j)
280#ifdef WET_DRY
281 rv(i,j,k,nrhs)=rv(i,j,k,nrhs)*vmask_wet(i,j)
282#endif
283#ifdef DIAGNOSTICS_UV
284 diarv(i,j,k,nrhs,m3pgrd)=rv(i,j,k,nrhs)
285#endif
286 END DO
287 END DO
288 END IF
289 END DO j_loop
290!
291 RETURN
292 END SUBROUTINE prsgrd40_tile
293
294 END MODULE prsgrd_mod
type(t_diags), dimension(:), allocatable diags
Definition mod_diags.F:100
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
integer, parameter inlm
Definition mod_param.F:662
real(dp) g
real(dp) rho0
integer m3pgrd
integer, dimension(:), allocatable nrhs
subroutine prsgrd40_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, umask_wet, vmask_wet, om_v, on_u, hz, z_w, rho, eq_tide, zetat, pair, diaru, diarv, ru, rv)
Definition prsgrd40.h:120
subroutine, public prsgrd(ng, tile)
Definition prsgrd31.h:36
recursive subroutine wclock_off(ng, model, region, line, routine)
Definition timers.F:148
recursive subroutine wclock_on(ng, model, region, line, routine)
Definition timers.F:3