ROMS
Loading...
Searching...
No Matches
rp_prsgrd40.h
Go to the documentation of this file.
1 MODULE rp_prsgrd_mod
2!
3!git $Id$
4!================================================== Hernan G. Arango ===
5! Copyright (c) 2002-2025 The ROMS Group Andrew M. Moore !
6! Licensed under a MIT/X style license !
7! See License_ROMS.md !
8!=======================================================================
9! !
10! This routine evaluates the representers tangent linear baroclinic, !
11! hydrostatic pressure gradient term using the finite-volume pressure !
12! Jacobian scheme of Lin (1997). !
13! !
14! The pressure gradient terms (m4/s2) are loaded into right-hand-side !
15! arrays "tl_ru" and "tl_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 :: rp_prsgrd
29!
30 CONTAINS
31!
32!***********************************************************************
33 SUBROUTINE rp_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, irpm, 23, __line__, myfile)
60#endif
61 CALL rp_prsgrd40_tile (ng, tile, &
62 & lbi, ubi, lbj, ubj, &
63 & imins, imaxs, jmins, jmaxs, &
64 & nrhs(ng), &
65 & grid(ng) % om_v, &
66 & grid(ng) % on_u, &
67 & grid(ng) % Hz, &
68 & grid(ng) % tl_Hz, &
69 & grid(ng) % z_w, &
70 & grid(ng) % tl_z_w, &
71 & ocean(ng) % rho, &
72 & ocean(ng) % tl_rho, &
73#ifdef TIDE_GENERATING_FORCES
74 & ocean(ng) % eq_tide, &
75 & ocean(ng) % tl_eq_tide, &
76#endif
77#ifdef ATM_PRESS
78 & forces(ng) % Pair, &
79#endif
80#ifdef DIAGNOSTICS_UV
81!! & DIAGS(ng) % DiaRU, &
82!! & DIAGS(ng) % DiaRV, &
83#endif
84 & ocean(ng) % tl_ru, &
85 & ocean(ng) % tl_rv)
86#ifdef PROFILE
87 CALL wclock_off (ng, irpm, 23, __line__, myfile)
88#endif
89!
90 RETURN
91 END SUBROUTINE rp_prsgrd
92!
93!***********************************************************************
94 SUBROUTINE rp_prsgrd40_tile (ng, tile, &
95 & LBi, UBi, LBj, UBj, &
96 & IminS, ImaxS, JminS, JmaxS, &
97 & nrhs, &
98 & om_v, on_u, &
99 & Hz, tl_Hz, &
100 & z_w, tl_z_w, &
101 & rho, tl_rho, &
102#ifdef TIDE_GENERATING_FORCES
103 & eq_tide, tl_eq_tide, &
104#endif
105#ifdef ATM_PRESS
106 & Pair, &
107#endif
108#ifdef DIAGNOSTICS_UV
109!! & DiaRU, DiaRV, &
110#endif
111 & tl_ru, tl_rv)
112!***********************************************************************
113!
114 USE mod_param
115 USE mod_scalars
116!
117! Imported variable declarations.
118!
119 integer, intent(in) :: ng, tile
120 integer, intent(in) :: LBi, UBi, LBj, UBj
121 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
122 integer, intent(in) :: nrhs
123
124#ifdef ASSUMED_SHAPE
125 real(r8), intent(in) :: om_v(LBi:,LBj:)
126 real(r8), intent(in) :: on_u(LBi:,LBj:)
127 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
128 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
129 real(r8), intent(in) :: rho(LBi:,LBj:,:)
130
131 real(r8), intent(in) :: tl_Hz(LBi:,LBj:,:)
132 real(r8), intent(in) :: tl_z_w(LBi:,LBj:,0:)
133 real(r8), intent(in) :: tl_rho(LBi:,LBj:,:)
134# ifdef TIDE_GENERATING_FORCES
135 real(r8), intent(in) :: eq_tide(LBi:,LBj:)
136 real(r8), intent(in) :: tl_eq_tide(LBi:,LBj:)
137# endif
138# ifdef ATM_PRESS
139 real(r8), intent(in) :: Pair(LBi:,LBj:)
140# endif
141# ifdef DIAGNOSTICS_UV
142!! real(r8), intent(inout) :: DiaRU(LBi:,LBj:,:,:,:)
143!! real(r8), intent(inout) :: DiaRV(LBi:,LBj:,:,:,:)
144# endif
145 real(r8), intent(inout) :: tl_ru(LBi:,LBj:,0:,:)
146 real(r8), intent(inout) :: tl_rv(LBi:,LBj:,0:,:)
147#else
148 real(r8), intent(in) :: om_v(LBi:UBi,LBj:UBj)
149 real(r8), intent(in) :: on_u(LBi:UBi,LBj:UBj)
150 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
151 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
152 real(r8), intent(in) :: rho(LBi:UBi,LBj:UBj,N(ng))
153
154 real(r8), intent(in) :: tl_Hz(LBi:UBi,LBj:UBj,N(ng))
155 real(r8), intent(in) :: tl_z_w(LBi:UBi,LBj:UBj,0:N(ng))
156 real(r8), intent(in) :: tl_rho(LBi:UBi,LBj:UBj,N(ng))
157# ifdef TIDE_GENERATING_FORCES
158 real(r8), intent(in) :: eq_tide(LBi:UBi,LBj:UBj)
159 real(r8), intent(in) :: tl_eq_tide(LBi:UBi,LBj:UBj)
160# endif
161# ifdef ATM_PRESS
162 real(r8), intent(in) :: Pair(LBi:UBi,LBj:UBj)
163# endif
164# ifdef DIAGNOSTICS_UV
165!! real(r8), intent(inout) :: DiaRU(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
166!! real(r8), intent(inout) :: DiaRV(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
167# endif
168 real(r8), intent(inout) :: tl_ru(LBi:UBi,LBj:UBj,0:N(ng),2)
169 real(r8), intent(inout) :: tl_rv(LBi:UBi,LBj:UBj,0:N(ng),2)
170#endif
171!
172! Local variable declarations.
173!
174 integer :: i, j, k
175
176 real(r8) :: cff, cff1, dh
177 real(r8) :: tl_dh
178#ifdef ATM_PRESS
179 real(r8) :: OneAtm, fac
180#endif
181 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: tl_FC
182
183 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,N(ng)) :: tl_FX
184
185 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,0:N(ng)) :: P
186 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,0:N(ng)) :: tl_P
187
188#include "set_bounds.h"
189!
190!-----------------------------------------------------------------------
191! Finite Volume pressure gradient algorithm (Lin, 1997).
192!-----------------------------------------------------------------------
193!
194! Compute pressure and its vertical integral. Initialize pressure at
195! the free-surface as zero.
196!
197#ifdef ATM_PRESS
198 oneatm=1013.25_r8 ! 1 atm = 1013.25 mb
199 fac=100.0_r8/g
200#endif
201 j_loop : DO j=jstrv-1,jend
202 DO i=istru-1,iend
203 p(i,j,n(ng))=0.0_r8
204#ifdef ATM_PRESS
205 p(i,j,n(ng))=p(i,j,n(ng))+fac*(pair(i,j)-oneatm)
206#endif
207 tl_p(i,j,n(ng))=0.0_r8
208#ifdef TIDE_GENERATING_FORCES
209 p(i,j,n(ng))=p(i,j,n(ng))-g*eq_tide(i,j)
210 tl_p(i,j,n(ng))=tl_p(i,j,n(ng))--g*tl_eq_tide(i,j)
211#endif
212 END DO
213 DO k=n(ng),1,-1
214 DO i=istru-1,iend
215 p(i,j,k-1)=p(i,j,k)+ &
216 & hz(i,j,k)*rho(i,j,k)
217 tl_p(i,j,k-1)=tl_p(i,j,k)+ &
218 & tl_hz(i,j,k)*rho(i,j,k)+ &
219 & hz(i,j,k)*tl_rho(i,j,k)
220!^ FX(i,j,k)=0.5_r8*Hz(i,j,k)*(P(i,j,k)+P(i,j,k-1))
221!^
222 tl_fx(i,j,k)=0.5_r8* &
223 & (tl_hz(i,j,k)*(p(i,j,k)+p(i,j,k-1))+ &
224 & hz(i,j,k)*(tl_p(i,j,k)+tl_p(i,j,k-1)))
225 END DO
226 END DO
227!
228! Calculate pressure gradient in the XI-direction (m4/s2).
229!
230 IF (j.ge.jstr) THEN
231 DO i=istru,iend
232!^ FC(i,N(ng))=0.0_r8
233!^
234 tl_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 tl_dh=tl_z_w(i,j,k-1)-tl_z_w(i-1,j,k-1)
242!^ FC(i,k-1)=0.5_r8*dh*(P(i,j,k-1)+P(i-1,j,k-1))
243!^
244 tl_fc(i,k-1)=0.5_r8* &
245 & (tl_dh*(p(i,j,k-1)+p(i-1,j,k-1))+ &
246 & dh*(tl_p(i,j,k-1)+tl_p(i-1,j,k-1)))
247!^ ru(i,j,k,nrhs)=(cff*(Hz(i-1,j,k)+ &
248!^ & Hz(i ,j,k))* &
249!^ & (z_w(i-1,j,N(ng))- &
250!^ & z_w(i ,j,N(ng)))+ &
251!^ & cff1*(FX(i-1,j,k)- &
252!^ & FX(i ,j,k)+ &
253!^ & FC(i,k )- &
254!^ & FC(i,k-1)))*on_u(i,j)
255!^
256 tl_ru(i,j,k,nrhs)=(cff*((tl_hz(i-1,j,k)+ &
257 & tl_hz(i ,j,k))* &
258 & (z_w(i-1,j,n(ng))- &
259 & z_w(i ,j,n(ng)))+ &
260 & (hz(i-1,j,k)+ &
261 & hz(i ,j,k))* &
262 & (tl_z_w(i-1,j,n(ng))- &
263 & tl_z_w(i ,j,n(ng))))+ &
264 & cff1*(tl_fx(i-1,j,k)- &
265 & tl_fx(i ,j,k)+ &
266 & tl_fc(i,k )- &
267 & tl_fc(i,k-1)))*on_u(i,j)
268#ifdef DIAGNOSTICS_UV
269!! DiaRU(i,j,k,nrhs,M3pgrd)=ru(i,j,k,nrhs)
270#endif
271 END DO
272 END DO
273 END IF
274!
275! Calculate pressure gradient in the ETA-direction (m4/s2).
276!
277 IF (j.ge.jstrv) THEN
278 DO i=istr,iend
279!^ FC(i,N(ng))=0.0_r8
280!^
281 tl_fc(i,n(ng))=0.0_r8
282 END DO
283 cff=0.5_r8*g
284 cff1=g/rho0
285 DO k=n(ng),1,-1
286 DO i=istr,iend
287 dh=z_w(i,j,k-1)-z_w(i,j-1,k-1)
288 tl_dh=tl_z_w(i,j,k-1)-tl_z_w(i,j-1,k-1)
289!^ FC(i,k-1)=0.5_r8*dh*(P(i,j,k-1)+P(i,j-1,k-1))
290!^
291 tl_fc(i,k-1)=0.5_r8* &
292 & (tl_dh*(p(i,j,k-1)+p(i,j-1,k-1))+ &
293 & dh*(tl_p(i,j,k-1)+tl_p(i,j-1,k-1))
294!^ rv(i,j,k,nrhs)=(cff*(Hz(i,j-1,k)+ &
295!^ & Hz(i,j ,k))* &
296!^ & (z_w(i,j-1,N(ng))- &
297!^ & z_w(i,j ,N(ng)))+ &
298!^ & cff1*(FX(i,j-1,k)- &
299!^ & FX(i,j ,k)+ &
300!^ & FC(i,k )- &
301!^ & FC(i,k-1)))*om_v(i,j)
302!^
303 tl_rv(i,j,k,nrhs)=(cff*((tl_hz(i,j-1,k)+ &
304 & tl_hz(i,j ,k))* &
305 & (z_w(i,j-1,n(ng))- &
306 & z_w(i,j ,n(ng)))+ &
307 & (hz(i,j-1,k)+ &
308 & hz(i,j ,k))* &
309 & (tl_z_w(i,j-1,n(ng))- &
310 & tl_z_w(i,j ,n(ng))))+ &
311 & cff1*(tl_fx(i,j-1,k)- &
312 & tl_fx(i,j ,k)+ &
313 & tl_fc(i,k )- &
314 & tl_fc(i,k-1)))*om_v(i,j)
315#ifdef DIAGNOSTICS_UV
316!! DiaRV(i,j,k,nrhs,M3pgrd)=rv(i,j,k,nrhs)
317#endif
318 END DO
319 END DO
320 END IF
321 END DO j_loop
322!
323 RETURN
324 END SUBROUTINE rp_prsgrd40_tile
325
326 END MODULE rp_prsgrd_mod
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 irpm
Definition mod_param.F:664
real(dp) g
real(dp) rho0
integer, dimension(:), allocatable nrhs
subroutine, public rp_prsgrd(ng, tile)
Definition rp_prsgrd31.h:36
subroutine rp_prsgrd40_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, om_v, on_u, hz, tl_hz, z_w, tl_z_w, rho, tl_rho, eq_tide, tl_eq_tide, pair, tl_ru, tl_rv)
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