ROMS
Loading...
Searching...
No Matches
rp_prsgrd_mod Module Reference

Functions/Subroutines

subroutine, public rp_prsgrd (ng, tile)
 
subroutine rp_prsgrd31_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, om_v, on_u, hz, tl_hz, z_r, tl_z_r, z_w, tl_z_w, rho, tl_rho, eq_tide, tl_eq_tide, pair, tl_ru, tl_rv)
 
subroutine rp_prsgrd32_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, umask, vmask, om_v, on_u, hz, tl_hz, z_r, tl_z_r, z_w, tl_z_w, rho, tl_rho, eq_tide, tl_eq_tide, pair, tl_ru, tl_rv)
 
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)
 

Function/Subroutine Documentation

◆ rp_prsgrd()

subroutine public rp_prsgrd_mod::rp_prsgrd ( integer, intent(in) ng,
integer, intent(in) tile )

Definition at line 35 of file rp_prsgrd31.h.

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, irpm, 23, __line__, myfile)
62#endif
63 CALL rp_prsgrd31_tile (ng, tile, &
64 & lbi, ubi, lbj, ubj, &
65 & imins, imaxs, jmins, jmaxs, &
66 & nrhs(ng), &
67 & grid(ng) % om_v, &
68 & grid(ng) % on_u, &
69 & grid(ng) % Hz, &
70 & grid(ng) % tl_Hz, &
71 & grid(ng) % z_r, &
72 & grid(ng) % tl_z_r, &
73 & grid(ng) % z_w, &
74 & grid(ng) % tl_z_w, &
75 & ocean(ng) % rho, &
76 & ocean(ng) % tl_rho, &
77#ifdef TIDE_GENERATING_FORCES
78 & ocean(ng) % eq_tide, &
79 & ocean(ng) % tl_eq_tide, &
80#endif
81#ifdef ATM_PRESS
82 & forces(ng) % Pair, &
83#endif
84#ifdef DIAGNOSTICS_UV
85!! & DIAGS(ng) % DiaRU, &
86!! & DIAGS(ng) % DiaRV, &
87#endif
88 & ocean(ng) % tl_ru, &
89 & ocean(ng) % tl_rv)
90#ifdef PROFILE
91 CALL wclock_off (ng, irpm, 23, __line__, myfile)
92#endif
93!
94 RETURN
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
integer, dimension(:), allocatable nrhs
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

References mod_forces::forces, mod_grid::grid, mod_param::irpm, mod_stepping::nrhs, mod_ocean::ocean, rp_prsgrd31_tile(), wclock_off(), and wclock_on().

Referenced by rp_rhs3d_mod::rp_rhs3d().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ rp_prsgrd31_tile()

subroutine rp_prsgrd_mod::rp_prsgrd31_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) imins,
integer, intent(in) imaxs,
integer, intent(in) jmins,
integer, intent(in) jmaxs,
integer, intent(in) nrhs,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) om_v,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) on_u,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) hz,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) tl_hz,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) z_r,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) tl_z_r,
real(r8), dimension(lbi:ubi,lbj:ubj,0:n(ng)), intent(in) z_w,
real(r8), dimension(lbi:ubi,lbj:ubj,0:n(ng)), intent(in) tl_z_w,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) rho,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) tl_rho,
real(r8), dimension(lbi:,lbj:), intent(in) eq_tide,
real(r8), dimension(lbi:,lbj:), intent(in) tl_eq_tide,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) pair,
real(r8), dimension(lbi:ubi,lbj:ubj,0:n(ng),2), intent(inout) tl_ru,
real(r8), dimension(lbi:ubi,lbj:ubj,0:n(ng),2), intent(inout) tl_rv )
private

Definition at line 98 of file rp_prsgrd31.h.

117!***********************************************************************
118!
119 USE mod_param
120 USE mod_scalars
121!
122! Imported variable declarations.
123!
124 integer, intent(in) :: ng, tile
125 integer, intent(in) :: LBi, UBi, LBj, UBj
126 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
127 integer, intent(in) :: nrhs
128
129#ifdef ASSUMED_SHAPE
130 real(r8), intent(in) :: om_v(LBi:,LBj:)
131 real(r8), intent(in) :: on_u(LBi:,LBj:)
132 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
133 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
134 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
135 real(r8), intent(in) :: rho(LBi:,LBj:,:)
136
137 real(r8), intent(in) :: tl_Hz(LBi:,LBj:,:)
138 real(r8), intent(in) :: tl_z_r(LBi:,LBj:,:)
139 real(r8), intent(in) :: tl_z_w(LBi:,LBj:,0:)
140 real(r8), intent(in) :: tl_rho(LBi:,LBj:,:)
141
142# ifdef TIDE_GENERATING_FORCES
143 real(r8), intent(in) :: eq_tide(LBi:,LBj:)
144 real(r8), intent(in) :: tl_eq_tide(LBi:,LBj:)
145# endif
146# ifdef ATM_PRESS
147 real(r8), intent(in) :: Pair(LBi:,LBj:)
148# endif
149# ifdef DIAGNOSTICS_UV
150!! real(r8), intent(inout) :: DiaRU(LBi:,LBj:,:,:,:)
151!! real(r8), intent(inout) :: DiaRV(LBi:,LBj:,:,:,:)
152# endif
153 real(r8), intent(inout) :: tl_ru(LBi:,LBj:,0:,:)
154 real(r8), intent(inout) :: tl_rv(LBi:,LBj:,0:,:)
155#else
156 real(r8), intent(in) :: om_v(LBi:UBi,LBj:UBj)
157 real(r8), intent(in) :: on_u(LBi:UBi,LBj:UBj)
158 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
159 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
160 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
161 real(r8), intent(in) :: rho(LBi:UBi,LBj:UBj,N(ng))
162
163 real(r8), intent(in) :: tl_Hz(LBi:UBi,LBj:UBj,N(ng))
164 real(r8), intent(in) :: tl_z_r(LBi:UBi,LBj:UBj,N(ng))
165 real(r8), intent(in) :: tl_z_w(LBi:UBi,LBj:UBj,0:N(ng))
166 real(r8), intent(in) :: tl_rho(LBi:UBi,LBj:UBj,N(ng))
167
168# ifdef TIDE_GENERATING_FORCES
169 real(r8), intent(in) :: eq_tide(LBi:,LBj:)
170 real(r8), intent(in) :: tl_eq_tide(LBi:,LBj:)
171# endif
172# ifdef ATM_PRESS
173 real(r8), intent(in) :: Pair(LBi:UBi,LBj:UBj)
174# endif
175# ifdef DIAGNOSTICS_UV
176!! real(r8), intent(inout) :: DiaRU(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
177!! real(r8), intent(inout) :: DiaRV(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
178# endif
179 real(r8), intent(inout) :: tl_ru(LBi:UBi,LBj:UBj,0:N(ng),2)
180 real(r8), intent(inout) :: tl_rv(LBi:UBi,LBj:UBj,0:N(ng),2)
181#endif
182!
183! Local variable declarations.
184!
185 integer :: i, j, k
186
187 real(r8) :: fac, fac1, fac2, fac3
188 real(r8) :: cff1, cff2, cff3, cff4
189 real(r8) :: tl_cff1, tl_cff2, tl_cff3, tl_cff4
190#ifdef WJ_GRADP
191 real(r8) :: gamma, tl_gamma
192#endif
193
194 real(r8), dimension(IminS:ImaxS) :: phie
195 real(r8), dimension(IminS:ImaxS) :: phix
196
197 real(r8), dimension(IminS:ImaxS) :: tl_phie
198 real(r8), dimension(IminS:ImaxS) :: tl_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 tl_cff1=tl_z_w(i ,j,n(ng))-tl_eq_tide(i ,j)- &
221 & tl_z_r(i ,j,n(ng))+ &
222 & tl_z_w(i-1,j,n(ng))-tl_eq_tide(i-1,j)- &
223 & tl_z_r(i-1,j,n(ng))
224#else
225 cff1=z_w(i ,j,n(ng))-z_r(i ,j,n(ng))+ &
226 & z_w(i-1,j,n(ng))-z_r(i-1,j,n(ng))
227 tl_cff1=tl_z_w(i ,j,n(ng))-tl_z_r(i ,j,n(ng))+ &
228 & tl_z_w(i-1,j,n(ng))-tl_z_r(i-1,j,n(ng))
229#endif
230 phix(i)=fac1*(rho(i,j,n(ng))-rho(i-1,j,n(ng)))*cff1
231 tl_phix(i)=fac1* &
232 & ((tl_rho(i,j,n(ng))-tl_rho(i-1,j,n(ng)))*cff1+ &
233 & (rho(i,j,n(ng))-rho(i-1,j,n(ng)))*tl_cff1)- &
234#ifdef TL_IOMS
235 & phix(i)
236#endif
237#ifdef ATM_PRESS
238 phix(i)=phix(i)+fac*(pair(i,j)-pair(i-1,j))
239#endif
240#ifdef RHO_SURF
241 phix(i)=phix(i)+ &
242 & (fac2+fac1*(rho(i,j,n(ng))+rho(i-1,j,n(ng))))* &
243 & (z_w(i,j,n(ng))-z_w(i-1,j,n(ng)))
244 tl_phix(i)=tl_phix(i)+ &
245 & (fac1*(tl_rho(i,j,n(ng))+tl_rho(i-1,j,n(ng))))* &
246 & (z_w(i,j,n(ng))-z_w(i-1,j,n(ng)))+ &
247 & (fac2+fac1*(rho(i,j,n(ng))+rho(i-1,j,n(ng))))* &
248 & (tl_z_w(i,j,n(ng))-tl_z_w(i-1,j,n(ng)))- &
249# ifdef TL_IOMS
250 & fac1*(rho(i,j,n(ng))+rho(i-1,j,n(ng)))* &
251 & (z_w(i,j,n(ng))-z_w(i-1,j,n(ng)))
252# endif
253#endif
254!^ ru(i,j,N(ng),nrhs)=-0.5_r8*(Hz(i,j,N(ng))+Hz(i-1,j,N(ng)))* &
255!^ & phix(i)*on_u(i,j)
256!^
257 tl_ru(i,j,n(ng),nrhs)=-0.5_r8*on_u(i,j)* &
258 & ((tl_hz(i ,j,n(ng))+ &
259 & tl_hz(i-1,j,n(ng)))*phix(i)+ &
260 & (hz(i ,j,n(ng))+ &
261 & hz(i-1,j,n(ng)))*tl_phix(i))+ &
262#ifdef TL_IOMS
263 & 0.5_r8*on_u(i,j)* &
264 & (hz(i ,j,n(ng))+ &
265 & hz(i-1,j,n(ng)))*phix(i)
266#endif
267#ifdef DIAGNOSTICS_UV
268!! DiaRU(i,j,N(ng),nrhs,M3pgrd)=ru(i,j,N(ng),nrhs)
269#endif
270 END DO
271!
272! Compute interior baroclinic pressure gradient. Differentiate and
273! then vertically integrate.
274!
275 DO k=n(ng)-1,1,-1
276 DO i=istru,iend
277#ifdef WJ_GRADP
278 cff1=1.0_r8/((z_r(i ,j,k+1)-z_r(i ,j,k))* &
279 & (z_r(i-1,j,k+1)-z_r(i-1,j,k)))
280 tl_cff1=-cff1*cff1*((tl_z_r(i ,j,k+1)-tl_z_r(i ,j,k))* &
281 & (z_r(i-1,j,k+1)-z_r(i-1,j,k))+ &
282 & (z_r(i ,j,k+1)-z_r(i ,j,k))* &
283 & (tl_z_r(i-1,j,k+1)-tl_z_r(i-1,j,k)))+ &
284# ifdef TL_IOMS
285 & 3.0_r8*cff1
286# endif
287 cff2=z_r(i ,j,k )-z_r(i-1,j,k )+ &
288 & z_r(i ,j,k+1)-z_r(i-1,j,k+1)
289 tl_cff2=tl_z_r(i ,j,k )-tl_z_r(i-1,j,k )+ &
290 & tl_z_r(i ,j,k+1)-tl_z_r(i-1,j,k+1)
291 cff3=z_r(i ,j,k+1)-z_r(i ,j,k )- &
292 & z_r(i-1,j,k+1)+z_r(i-1,j,k )
293 tl_cff3=tl_z_r(i ,j,k+1)-tl_z_r(i ,j,k )- &
294 & tl_z_r(i-1,j,k+1)+tl_z_r(i-1,j,k )
295 gamma=0.125_r8*cff1*cff2*cff3
296 tl_gamma=0.125_r8*(tl_cff1*cff2*cff3+ &
297 & cff1*(tl_cff2*cff3+ &
298 & cff2*tl_cff3))- &
299# ifdef TL_IOMS
300 & 2.0_r8*gamma
301# endif
302
303 cff1=(1.0_r8+gamma)*(rho(i,j,k+1)-rho(i-1,j,k+1))+ &
304 & (1.0_r8-gamma)*(rho(i,j,k )-rho(i-1,j,k ))
305 tl_cff1=tl_gamma*(rho(i,j,k+1)-rho(i-1,j,k+1)- &
306 & rho(i,j,k )+rho(i-1,j,k ))+ &
307 & (1.0_r8+gamma)*(tl_rho(i ,j,k+1)- &
308 & tl_rho(i-1,j,k+1))+ &
309 & (1.0_r8-gamma)*(tl_rho(i ,j,k )- &
310 & tl_rho(i-1,j,k ))- &
311# ifdef TL_IOMS
312 & gamma*((rho(i,j,k+1)-rho(i-1,j,k+1))- &
313 & (rho(i,j,k )-rho(i-1,j,k )))
314# endif
315 cff2=rho(i,j,k+1)+rho(i-1,j,k+1)- &
316 & rho(i,j,k )-rho(i-1,j,k )
317 tl_cff2=tl_rho(i,j,k+1)+tl_rho(i-1,j,k+1)- &
318 & tl_rho(i,j,k )-tl_rho(i-1,j,k )
319 cff3=z_r(i,j,k+1)+z_r(i-1,j,k+1)- &
320 & z_r(i,j,k )-z_r(i-1,j,k )
321 tl_cff3=tl_z_r(i,j,k+1)+tl_z_r(i-1,j,k+1)- &
322 & tl_z_r(i,j,k )-tl_z_r(i-1,j,k )
323 cff4=(1.0_r8+gamma)*(z_r(i,j,k+1)-z_r(i-1,j,k+1))+ &
324 & (1.0_r8-gamma)*(z_r(i,j,k )-z_r(i-1,j,k ))
325 tl_cff4=tl_gamma*(z_r(i,j,k+1)-z_r(i-1,j,k+1)- &
326 & z_r(i,j,k )+z_r(i-1,j,k ))+ &
327 & (1.0_r8+gamma)*(tl_z_r(i ,j,k+1)- &
328 & tl_z_r(i-1,j,k+1))+ &
329 & (1.0_r8-gamma)*(tl_z_r(i ,j,k )- &
330 & tl_z_r(i-1,j,k ))- &
331# ifdef TL_IOMS
332 & gamma*((z_r(i,j,k+1)-z_r(i-1,j,k+1))- &
333 & (z_r(i,j,k )-z_r(i-1,j,k )))
334# endif
335 phix(i)=phix(i)+ &
336 & fac3*(cff1*cff3-cff2*cff4)
337 tl_phix(i)=tl_phix(i)+ &
338 & fac3*(tl_cff1*cff3+ &
339 & cff1*tl_cff3- &
340 & tl_cff2*cff4- &
341 & cff2*tl_cff4)- &
342# ifdef TL_IOMS
343 & fac3*(cff1*cff3- &
344 & cff2*cff4)
345# endif
346#else
347 cff1=rho(i,j,k+1)-rho(i-1,j,k+1)+ &
348 & rho(i,j,k )-rho(i-1,j,k )
349 cff2=rho(i,j,k+1)+rho(i-1,j,k+1)- &
350 & rho(i,j,k )-rho(i-1,j,k )
351 tl_cff1=tl_rho(i,j,k+1)-tl_rho(i-1,j,k+1)+ &
352 & tl_rho(i,j,k )-tl_rho(i-1,j,k )
353 tl_cff2=tl_rho(i,j,k+1)+tl_rho(i-1,j,k+1)- &
354 & tl_rho(i,j,k )-tl_rho(i-1,j,k )
355 cff3=z_r(i,j,k+1)+z_r(i-1,j,k+1)- &
356 & z_r(i,j,k )-z_r(i-1,j,k )
357 cff4=z_r(i,j,k+1)-z_r(i-1,j,k+1)+ &
358 & z_r(i,j,k )-z_r(i-1,j,k )
359 tl_cff3=tl_z_r(i,j,k+1)+tl_z_r(i-1,j,k+1)- &
360 & tl_z_r(i,j,k )-tl_z_r(i-1,j,k )
361 tl_cff4=tl_z_r(i,j,k+1)-tl_z_r(i-1,j,k+1)+ &
362 & tl_z_r(i,j,k )-tl_z_r(i-1,j,k )
363 phix(i)=phix(i)+ &
364 & fac3*(cff1*cff3-cff2*cff4)
365 tl_phix(i)=tl_phix(i)+ &
366 & fac3*(tl_cff1*cff3+ &
367 & cff1*tl_cff3- &
368 & tl_cff2*cff4- &
369 & cff2*tl_cff4)- &
370# ifdef TL_IOMS
371 & fac3*(cff1*cff3- &
372 & cff2*cff4)
373# endif
374#endif
375!^ ru(i,j,k,nrhs)=-0.5_r8*(Hz(i,j,k)+Hz(i-1,j,k))* &
376!^ & phix(i)*on_u(i,j)
377!^
378 tl_ru(i,j,k,nrhs)=-0.5_r8*on_u(i,j)* &
379 & ((tl_hz(i,j,k)+tl_hz(i-1,j,k))* &
380 & phix(i)+ &
381 & (hz(i,j,k)+hz(i-1,j,k))* &
382 & tl_phix(i))+ &
383#ifdef TL_IOMS
384 & 0.5_r8*on_u(i,j)* &
385 & (hz(i,j,k)+hz(i-1,j,k))*phix(i)
386#endif
387#ifdef DIAGNOSTICS_UV
388!! DiaRU(i,j,k,nrhs,M3pgrd)=ru(i,j,k,nrhs)
389#endif
390 END DO
391 END DO
392!
393!-----------------------------------------------------------------------
394! Calculate pressure gradient in the ETA-direction (m4/s2).
395!-----------------------------------------------------------------------
396!
397! Compute surface baroclinic pressure gradient.
398!
399 IF (j.ge.jstrv) THEN
400 DO i=istr,iend
401 cff1=z_w(i,j ,n(ng))-z_r(i,j ,n(ng))+ &
402 & z_w(i,j-1,n(ng))-z_r(i,j-1,n(ng))
403 tl_cff1=tl_z_w(i,j ,n(ng))-tl_z_r(i,j ,n(ng))+ &
404 & tl_z_w(i,j-1,n(ng))-tl_z_r(i,j-1,n(ng))
405 phie(i)=fac1*(rho(i,j,n(ng))-rho(i,j-1,n(ng)))*cff1
406 tl_phie(i)=fac1* &
407 & ((tl_rho(i,j,n(ng))-tl_rho(i,j-1,n(ng)))*cff1+ &
408 & (rho(i,j,n(ng))-rho(i,j-1,n(ng)))*tl_cff1)- &
409#ifdef TL_IOMS
410 & phie(i)
411#endif
412#ifdef ATM_PRESS
413 phie(i)=phie(i)+fac*(pair(i,j)-pair(i,j-1))
414#endif
415#ifdef RHO_SURF
416 phie(i)=phie(i)+ &
417 & (fac2+fac1*(rho(i,j,n(ng))+rho(i,j-1,n(ng))))* &
418 & (z_w(i,j,n(ng))-z_w(i,j-1,n(ng)))
419 tl_phie(i)=tl_phie(i)+ &
420 & (fac1*(tl_rho(i,j,n(ng))+tl_rho(i,j-1,n(ng))))* &
421 & (z_w(i,j,n(ng))-z_w(i,j-1,n(ng)))+ &
422 & (fac2+fac1*(rho(i,j,n(ng))+rho(i,j-1,n(ng))))* &
423 & (tl_z_w(i,j,n(ng))-tl_z_w(i,j-1,n(ng)))- &
424# ifdef TL_IOMS
425 & fac1*(rho(i,j,n(ng))+rho(i,j-1,n(ng)))* &
426 & (z_w(i,j,n(ng))-z_w(i,j-1,n(ng)))
427# endif
428#endif
429!^ rv(i,j,N(ng),nrhs)=-0.5_r8*(Hz(i,j,N(ng))+Hz(i,j-1,N(ng)))* &
430!^ & phie(i)*om_v(i,j)
431!^
432 tl_rv(i,j,n(ng),nrhs)=-0.5_r8*om_v(i,j)* &
433 & ((tl_hz(i,j ,n(ng))+ &
434 & tl_hz(i,j-1,n(ng)))*phie(i)+ &
435 & (hz(i,j ,n(ng))+ &
436 & hz(i,j-1,n(ng)))*tl_phie(i))+ &
437#ifdef TL_IOMS
438 & 0.5_r8*om_v(i,j)* &
439 & (hz(i,j ,n(ng))+ &
440 & hz(i,j-1,n(ng)))*phie(i)
441#endif
442#ifdef DIAGNOSTICS_UV
443!! DiaRV(i,j,N(ng),nrhs,M3pgrd)=rv(i,j,N(ng),nrhs)
444#endif
445 END DO
446!
447! Compute interior baroclinic pressure gradient. Differentiate and
448! then vertically integrate.
449!
450 DO k=n(ng)-1,1,-1
451 DO i=istr,iend
452#ifdef WJ_GRADP
453 cff1=1.0_r8/((z_r(i,j ,k+1)-z_r(i,j ,k))* &
454 & (z_r(i,j-1,k+1)-z_r(i,j-1,k)))
455 tl_cff1=-cff1*cff1*((tl_z_r(i,j ,k+1)-tl_z_r(i,j ,k))* &
456 & (z_r(i,j-1,k+1)-z_r(i,j-1,k))+ &
457 & (z_r(i,j ,k+1)-z_r(i,j ,k))* &
458 & (tl_z_r(i,j-1,k+1)-tl_z_r(i,j-1,k)))+ &
459# ifdef TL_IOMS
460 & 3.0_r8*cff1
461# endif
462 cff2=z_r(i,j ,k )-z_r(i,j-1,k )+ &
463 & z_r(i,j ,k+1)-z_r(i,j-1,k+1)
464 tl_cff2=tl_z_r(i,j ,k )-tl_z_r(i,j-1,k )+ &
465 & tl_z_r(i,j ,k+1)-tl_z_r(i,j-1,k+1)
466 cff3=z_r(i,j ,k+1)-z_r(i,j ,k )- &
467 & z_r(i,j-1,k+1)+z_r(i,j-1,k )
468 tl_cff3=tl_z_r(i,j ,k+1)-tl_z_r(i,j ,k )- &
469 & tl_z_r(i,j-1,k+1)+tl_z_r(i,j-1,k )
470 gamma=0.125_r8*cff1*cff2*cff3
471 tl_gamma=0.125_r8*(tl_cff1*cff2*cff3+ &
472 & cff1*(tl_cff2*cff3+ &
473 & cff2*tl_cff3))- &
474# ifdef TL_IOMS
475 & 2.0_r8*gamma
476# endif
477
478 cff1=(1.0_r8+gamma)*(rho(i,j,k+1)-rho(i,j-1,k+1))+ &
479 & (1.0_r8-gamma)*(rho(i,j,k )-rho(i,j-1,k ))
480 tl_cff1=tl_gamma*(rho(i,j,k+1)-rho(i,j-1,k+1)- &
481 & rho(i,j,k )+rho(i,j-1,k ))+ &
482 & (1.0_r8+gamma)*(tl_rho(i,j ,k+1)- &
483 & tl_rho(i,j-1,k+1))+ &
484 & (1.0_r8-gamma)*(tl_rho(i,j ,k )- &
485 & tl_rho(i,j-1,k ))- &
486# ifdef TL_IOMS
487 & gamma*((rho(i,j,k+1)-rho(i,j-1,k+1))- &
488 & (rho(i,j,k )-rho(i,j-1,k )))
489# endif
490 cff2=rho(i,j,k+1)+rho(i,j-1,k+1)- &
491 & rho(i,j,k )-rho(i,j-1,k )
492 tl_cff2=tl_rho(i,j,k+1)+tl_rho(i,j-1,k+1)- &
493 & tl_rho(i,j,k )-tl_rho(i,j-1,k )
494 cff3=z_r(i,j,k+1)+z_r(i,j-1,k+1)- &
495 & z_r(i,j,k )-z_r(i,j-1,k )
496 tl_cff3=tl_z_r(i,j,k+1)+tl_z_r(i,j-1,k+1)- &
497 & tl_z_r(i,j,k )-tl_z_r(i,j-1,k )
498 cff4=(1.0_r8+gamma)*(z_r(i,j,k+1)-z_r(i,j-1,k+1))+ &
499 & (1.0_r8-gamma)*(z_r(i,j,k )-z_r(i,j-1,k ))
500 tl_cff4=tl_gamma*(z_r(i,j,k+1)-z_r(i,j-1,k+1)- &
501 & z_r(i,j,k )+z_r(i,j-1,k ))+ &
502 & (1.0_r8+gamma)*(tl_z_r(i,j ,k+1)- &
503 & tl_z_r(i,j-1,k+1))+ &
504 & (1.0_r8-gamma)*(tl_z_r(i,j ,k )- &
505 & tl_z_r(i,j-1,k ))- &
506# ifdef TL_IOMS
507 & gamma*((z_r(i,j,k+1)-z_r(i,j-1,k+1))- &
508 & (z_r(i,j,k )-z_r(i,j-1,k )))
509# endif
510 phie(i)=phie(i)+ &
511 & fac3*(cff1*cff3-cff2*cff4)
512 tl_phie(i)=tl_phie(i)+ &
513 & fac3*(tl_cff1*cff3+ &
514 & cff1*tl_cff3- &
515 & tl_cff2*cff4- &
516 & cff2*tl_cff4)- &
517# ifdef TL_IOMS
518 & fac3*(cff1*cff3- &
519 & cff2*cff4)
520# endif
521#else
522 cff1=rho(i,j,k+1)-rho(i,j-1,k+1)+ &
523 & rho(i,j,k )-rho(i,j-1,k )
524 cff2=rho(i,j,k+1)+rho(i,j-1,k+1)- &
525 & rho(i,j,k )-rho(i,j-1,k )
526 tl_cff1=tl_rho(i,j,k+1)-tl_rho(i,j-1,k+1)+ &
527 & tl_rho(i,j,k )-tl_rho(i,j-1,k )
528 tl_cff2=tl_rho(i,j,k+1)+tl_rho(i,j-1,k+1)- &
529 & tl_rho(i,j,k )-tl_rho(i,j-1,k )
530 cff3=z_r(i,j,k+1)+z_r(i,j-1,k+1)- &
531 & z_r(i,j,k )-z_r(i,j-1,k )
532 cff4=z_r(i,j,k+1)-z_r(i,j-1,k+1)+ &
533 & z_r(i,j,k )-z_r(i,j-1,k )
534 tl_cff3=tl_z_r(i,j,k+1)+tl_z_r(i,j-1,k+1)- &
535 & tl_z_r(i,j,k )-tl_z_r(i,j-1,k )
536 tl_cff4=tl_z_r(i,j,k+1)-tl_z_r(i,j-1,k+1)+ &
537 & tl_z_r(i,j,k )-tl_z_r(i,j-1,k )
538 phie(i)=phie(i)+ &
539 & fac3*(cff1*cff3-cff2*cff4)
540 tl_phie(i)=tl_phie(i)+ &
541 & fac3*(tl_cff1*cff3+ &
542 & cff1*tl_cff3- &
543 & tl_cff2*cff4- &
544 & cff2*tl_cff4)- &
545# ifdef TL_IOMS
546 & fac3*(cff1*cff3- &
547 & cff2*cff4)
548# endif
549#endif
550!^ rv(i,j,k,nrhs)=-0.5_r8*(Hz(i,j,k)+Hz(i,j-1,k))* &
551!^ & phie(i)*om_v(i,j)
552!^
553 tl_rv(i,j,k,nrhs)=-0.5_r8*om_v(i,j)* &
554 & ((tl_hz(i,j,k)+tl_hz(i,j-1,k))* &
555 & phie(i)+ &
556 & (hz(i,j,k)+hz(i,j-1,k))* &
557 & tl_phie(i))+ &
558#ifdef TL_IOMS
559 & 0.5_r8*om_v(i,j)* &
560 & (hz(i,j,k)+hz(i,j-1,k))*phie(i)
561#endif
562#ifdef DIAGNOSTICS_UV
563!! DiaRV(i,j,k,nrhs,M3pgrd)=rv(i,j,k,nrhs)
564#endif
565 END DO
566 END DO
567 END IF
568 END DO
569!
570 RETURN
integer, dimension(:), allocatable n
Definition mod_param.F:479
real(dp) g
real(dp) rho0

References mod_scalars::g, and mod_scalars::rho0.

Referenced by rp_prsgrd().

Here is the caller graph for this function:

◆ rp_prsgrd32_tile()

subroutine rp_prsgrd_mod::rp_prsgrd32_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) imins,
integer, intent(in) imaxs,
integer, intent(in) jmins,
integer, intent(in) jmaxs,
integer, intent(in) nrhs,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) umask,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) vmask,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) om_v,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) on_u,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) hz,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) tl_hz,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) z_r,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) tl_z_r,
real(r8), dimension(lbi:ubi,lbj:ubj,0:n(ng)), intent(in) z_w,
real(r8), dimension(lbi:ubi,lbj:ubj,0:n(ng)), intent(in) tl_z_w,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) rho,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) tl_rho,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) eq_tide,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) tl_eq_tide,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) pair,
real(r8), dimension(lbi:ubi,lbj:ubj,0:n(ng),2), intent(inout) tl_ru,
real(r8), dimension(lbi:ubi,lbj:ubj,0:n(ng),2), intent(inout) tl_rv )
private

Definition at line 107 of file rp_prsgrd32.h.

129!***********************************************************************
130!
131 USE mod_param
132 USE mod_scalars
133!
134! Imported variable declarations.
135!
136 integer, intent(in) :: ng, tile
137 integer, intent(in) :: LBi, UBi, LBj, UBj
138 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
139 integer, intent(in) :: nrhs
140
141#ifdef ASSUMED_SHAPE
142# ifdef MASKING
143 real(r8), intent(in) :: umask(LBi:,LBj:)
144 real(r8), intent(in) :: vmask(LBi:,LBj:)
145# endif
146 real(r8), intent(in) :: om_v(LBi:,LBj:)
147 real(r8), intent(in) :: on_u(LBi:,LBj:)
148 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
149 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
150 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
151 real(r8), intent(in) :: rho(LBi:,LBj:,:)
152
153 real(r8), intent(in) :: tl_Hz(LBi:,LBj:,:)
154 real(r8), intent(in) :: tl_z_r(LBi:,LBj:,:)
155 real(r8), intent(in) :: tl_z_w(LBi:,LBj:,0:)
156 real(r8), intent(in) :: tl_rho(LBi:,LBj:,:)
157# ifdef TIDE_GENERATING_FORCES
158 real(r8), intent(in) :: eq_tide(LBi:,LBj:)
159 real(r8), intent(in) :: tl_eq_tide(LBi:,LBj:)
160# endif
161# ifdef ATM_PRESS
162 real(r8), intent(in) :: Pair(LBi:,LBj:)
163# endif
164# ifdef DIAGNOSTICS_UV
165!! real(r8), intent(inout) :: DiaRU(LBi:,LBj:,:,:,:)
166!! real(r8), intent(inout) :: DiaRV(LBi:,LBj:,:,:,:)
167# endif
168 real(r8), intent(inout) :: tl_ru(LBi:,LBj:,0:,:)
169 real(r8), intent(inout) :: tl_rv(LBi:,LBj:,0:,:)
170#else
171# ifdef MASKING
172 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
173 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
174# endif
175 real(r8), intent(in) :: om_v(LBi:UBi,LBj:UBj)
176 real(r8), intent(in) :: on_u(LBi:UBi,LBj:UBj)
177 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
178 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
179 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
180 real(r8), intent(in) :: rho(LBi:UBi,LBj:UBj,N(ng))
181
182 real(r8), intent(in) :: tl_Hz(LBi:UBi,LBj:UBj,N(ng))
183 real(r8), intent(in) :: tl_z_r(LBi:UBi,LBj:UBj,N(ng))
184 real(r8), intent(in) :: tl_z_w(LBi:UBi,LBj:UBj,0:N(ng))
185 real(r8), intent(in) :: tl_rho(LBi:UBi,LBj:UBj,N(ng))
186# ifdef TIDE_GENERATING_FORCES
187 real(r8), intent(in) :: eq_tide(LBi:UBi,LBj:UBj)
188 real(r8), intent(in) :: tl_eq_tide(LBi:UBi,LBj:UBj)
189# endif
190# ifdef ATM_PRESS
191 real(r8), intent(in) :: Pair(LBi:UBi,LBj:UBj)
192# endif
193# ifdef DIAGNOSTICS_UV
194!! real(r8), intent(inout) :: DiaRU(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
195!! real(r8), intent(inout) :: DiaRV(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
196# endif
197 real(r8), intent(inout) :: tl_ru(LBi:UBi,LBj:UBj,0:N(ng),2)
198 real(r8), intent(inout) :: tl_rv(LBi:UBi,LBj:UBj,0:N(ng),2)
199#endif
200!
201! Local variable declarations.
202!
203 integer :: i, j, k
204
205 real(r8), parameter :: OneFifth = 0.2_r8
206 real(r8), parameter :: OneTwelfth = 1.0_r8/12.0_r8
207 real(r8), parameter :: eps = 1.0e-10_r8
208
209 real(r8) :: GRho, GRho0, HalfGRho
210 real(r8) :: cff, cff1, cff2
211 real(r8) :: tl_cff, tl_cff1, tl_cff2
212#ifdef ATM_PRESS
213 real(r8) :: OneAtm, fac
214#endif
215 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,N(ng)) :: P
216
217 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,N(ng)) :: tl_P
218
219 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: dR
220 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: dR1
221 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: dZ
222 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: dZ1
223
224 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: tl_dR
225 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: tl_dZ
226
227 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FC
228 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: aux
229 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: dRx
230 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: dZx
231
232 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_FC
233 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_aux
234 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_dRx
235 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_dZx
236
237#include "set_bounds.h"
238!
239!-----------------------------------------------------------------------
240! Preliminary step (same for XI- and ETA-components:
241!-----------------------------------------------------------------------
242!
243 grho=g/rho0
244 grho0=1000.0_r8*grho
245 halfgrho=0.5_r8*grho
246#ifdef ATM_PRESS
247 oneatm=1013.25_r8 ! 1 atm = 1013.25 mb
248 fac=100.0_r8/rho0
249#endif
250!
251 DO j=jstrv-1,jend
252 DO k=1,n(ng)-1
253 DO i=istru-1,iend
254 dr(i,k)=rho(i,j,k+1)-rho(i,j,k)
255 tl_dr(i,k)=tl_rho(i,j,k+1)-tl_rho(i,j,k)
256 dz(i,k)=z_r(i,j,k+1)-z_r(i,j,k)
257 tl_dz(i,k)=tl_z_r(i,j,k+1)-tl_z_r(i,j,k)
258 END DO
259 END DO
260 DO i=istru-1,iend
261 dr(i,n(ng))=dr(i,n(ng)-1)
262 tl_dr(i,n(ng))=tl_dr(i,n(ng)-1)
263 dz(i,n(ng))=dz(i,n(ng)-1)
264 tl_dz(i,n(ng))=tl_dz(i,n(ng)-1)
265 dr(i,0)=dr(i,1)
266 tl_dr(i,0)=tl_dr(i,1)
267 dz(i,0)=dz(i,1)
268 tl_dz(i,0)=tl_dz(i,1)
269 END DO
270 DO k=n(ng),1,-1
271 DO i=istru-1,iend
272 cff=2.0_r8*dr(i,k)*dr(i,k-1)
273 tl_cff=2.0_r8*(tl_dr(i,k)*dr(i,k-1)+ &
274 & dr(i,k)*tl_dr(i,k-1))- &
275#ifdef TL_IOMS
276 & cff
277#endif
278 dr1(i,k)=dr(i,k)
279 IF (cff.gt.eps) THEN
280 dr(i,k)=cff/(dr(i,k)+dr(i,k-1))
281 tl_dr(i,k)=(tl_cff-dr(i,k)*(tl_dr(i,k)+tl_dr(i,k-1)))/ &
282 & (dr1(i,k)+dr(i,k-1))+ &
283#ifdef TL_IOMS
284 & dr(i,k)
285#endif
286 ELSE
287 dr(i,k)=0.0_r8
288 tl_dr(i,k)=0.0_r8
289 END IF
290 dz1(i,k)=dz(i,k)
291 dz(i,k)=2.0_r8*dz(i,k)*dz(i,k-1)/(dz(i,k)+dz(i,k-1))
292 tl_dz(i,k)=(2.0_r8*(tl_dz(i,k)*dz(i,k-1)+ &
293 & dz1(i,k)*tl_dz(i,k-1))- &
294 & dz(i,k)*(tl_dz(i,k)+tl_dz(i,k-1)))/ &
295 & (dz1(i,k)+dz(i,k-1))
296 END DO
297 END DO
298 DO i=istru-1,iend
299 cff1=1.0_r8/(z_r(i,j,n(ng))-z_r(i,j,n(ng)-1))
300 tl_cff1=-cff1*cff1*(tl_z_r(i,j,n(ng))-tl_z_r(i,j,n(ng)-1))+ &
301#ifdef TL_IOMS
302 & 2.0_r8*cff1
303#endif
304 cff2=0.5_r8*(rho(i,j,n(ng))-rho(i,j,n(ng)-1))* &
305 & (z_w(i,j,n(ng))-z_r(i,j,n(ng)))*cff1
306 tl_cff2=0.5_r8*((tl_rho(i,j,n(ng))-tl_rho(i,j,n(ng)-1))* &
307 & (z_w(i,j,n(ng))-z_r(i,j,n(ng)))*cff1+ &
308 & (rho(i,j,n(ng))-rho(i,j,n(ng)-1))* &
309 & ((tl_z_w(i,j,n(ng))-tl_z_r(i,j,n(ng)))*cff1+ &
310 & (z_w(i,j,n(ng))-z_r(i,j,n(ng)))*tl_cff1))- &
311#ifdef TL_IOMS
312 & 2.0_r8*cff2
313#endif
314 p(i,j,n(ng))=g*z_w(i,j,n(ng))+ &
315#ifdef ATM_PRESS
316 & fac*(pair(i,j)-oneatm)+ &
317#endif
318 & grho*(rho(i,j,n(ng))+cff2)* &
319 & (z_w(i,j,n(ng))-z_r(i,j,n(ng)))
320 tl_p(i,j,n(ng))=g*tl_z_w(i,j,n(ng))+ &
321 & grho*((tl_rho(i,j,n(ng))+tl_cff2)* &
322 & (z_w(i,j,n(ng))-z_r(i,j,n(ng)))+ &
323 & (rho(i,j,n(ng))+cff2)* &
324 & (tl_z_w(i,j,n(ng))-tl_z_r(i,j,n(ng))))- &
325#ifdef TL_IOMS
326 & grho*(rho(i,j,n(ng))+cff2)* &
327 & (z_w(i,j,n(ng))-z_r(i,j,n(ng)))
328#endif
329#ifdef TIDE_GENERATING_FORCES
330 p(i,j,n(ng))=p(i,j,n(ng))-g*eq_tide(i,j)
331 tl_p(i,j,n(ng))=tl_p(i,j,n(ng))-g*tl_eq_tide(i,j)
332#endif
333 END DO
334 DO k=n(ng)-1,1,-1
335 DO i=istru-1,iend
336 cff=halfgrho*((rho(i,j,k+1)+rho(i,j,k))* &
337 & (z_r(i,j,k+1)-z_r(i,j,k))- &
338 & onefifth* &
339 & ((dr(i,k+1)-dr(i,k))* &
340 & (z_r(i,j,k+1)-z_r(i,j,k)- &
341 & onetwelfth* &
342 & (dz(i,k+1)+dz(i,k)))- &
343 & (dz(i,k+1)-dz(i,k))* &
344 & (rho(i,j,k+1)-rho(i,j,k)- &
345 & onetwelfth* &
346 & (dr(i,k+1)+dr(i,k)))))
347 tl_cff=halfgrho*((tl_rho(i,j,k+1)+tl_rho(i,j,k))* &
348 & (z_r(i,j,k+1)-z_r(i,j,k))+ &
349 & (rho(i,j,k+1)+rho(i,j,k))* &
350 & (tl_z_r(i,j,k+1)-tl_z_r(i,j,k))- &
351 & onefifth* &
352 & ((tl_dr(i,k+1)-tl_dr(i,k))* &
353 & (z_r(i,j,k+1)-z_r(i,j,k)- &
354 & onetwelfth* &
355 & (dz(i,k+1)+dz(i,k)))+ &
356 & (dr(i,k+1)-dr(i,k))* &
357 & (tl_z_r(i,j,k+1)-tl_z_r(i,j,k)- &
358 & onetwelfth* &
359 & (tl_dz(i,k+1)+tl_dz(i,k)))- &
360 & (tl_dz(i,k+1)-tl_dz(i,k))* &
361 & (rho(i,j,k+1)-rho(i,j,k)- &
362 & onetwelfth* &
363 & (dr(i,k+1)+dr(i,k)))- &
364 & (dz(i,k+1)-dz(i,k))* &
365 & (tl_rho(i,j,k+1)-tl_rho(i,j,k)- &
366 & onetwelfth* &
367 & (tl_dr(i,k+1)+tl_dr(i,k)))))- &
368#ifdef TL_IOMS
369 & cff
370#endif
371 p(i,j,k)=p(i,j,k+1)+cff
372 tl_p(i,j,k)=tl_p(i,j,k+1)+tl_cff
373 END DO
374 END DO
375 END DO
376!
377!-----------------------------------------------------------------------
378! Compute XI-component pressure gradient term.
379!-----------------------------------------------------------------------
380!
381 DO k=n(ng),1,-1
382 DO j=jstr,jend
383 DO i=istru-1,iend+1
384 aux(i,j)=z_r(i,j,k)-z_r(i-1,j,k)
385 tl_aux(i,j)=tl_z_r(i,j,k)-tl_z_r(i-1,j,k)
386#ifdef MASKING
387 aux(i,j)=aux(i,j)*umask(i,j)
388 tl_aux(i,j)=tl_aux(i,j)*umask(i,j)
389#endif
390 fc(i,j)=rho(i,j,k)-rho(i-1,j,k)
391 tl_fc(i,j)=tl_rho(i,j,k)-tl_rho(i-1,j,k)
392#ifdef MASKING
393 fc(i,j)=fc(i,j)*umask(i,j)
394 tl_fc(i,j)=tl_fc(i,j)*umask(i,j)
395#endif
396 END DO
397 END DO
398!
399 DO j=jstr,jend
400 DO i=istru-1,iend
401 cff=2.0_r8*aux(i,j)*aux(i+1,j)
402 tl_cff=2.0_r8*(tl_aux(i,j)*aux(i+1,j)+ &
403 & aux(i,j)*tl_aux(i+1,j))- &
404#ifdef TL_IOMS
405 & cff
406#endif
407 IF (cff.gt.eps) THEN
408 cff1=1.0_r8/(aux(i,j)+aux(i+1,j))
409 tl_cff1=-cff1*cff1*(tl_aux(i,j)+tl_aux(i+1,j))+ &
410#ifdef TL_IOMS
411 & 2.0_r8*cff1
412#endif
413 dzx(i,j)=cff*cff1
414 tl_dzx(i,j)=tl_cff*cff1+cff*tl_cff1- &
415#ifdef TL_IOMS
416 & dzx(i,j)
417#endif
418 ELSE
419 dzx(i,j)=0.0_r8
420 tl_dzx(i,j)=0.0_r8
421 END IF
422 cff1=2.0_r8*fc(i,j)*fc(i+1,j)
423 tl_cff1=2.0_r8*(tl_fc(i,j)*fc(i+1,j)+ &
424 & fc(i,j)*tl_fc(i+1,j))- &
425#ifdef TL_IOMS
426 & cff1
427#endif
428 IF (cff1.gt.eps) THEN
429 cff2=1.0_r8/(fc(i,j)+fc(i+1,j))
430 tl_cff2=-cff2*cff2*(tl_fc(i,j)+tl_fc(i+1,j))+ &
431#ifdef TL_IOMS
432 & 2.0_r8*cff2
433#endif
434 drx(i,j)=cff1*cff2
435 tl_drx(i,j)=tl_cff1*cff2+cff1*tl_cff2- &
436#ifdef TL_IOMS
437 & drx(i,j)
438#endif
439 ELSE
440 drx(i,j)=0.0_r8
441 tl_drx(i,j)=0.0_r8
442 END IF
443 END DO
444 END DO
445!
446 DO j=jstr,jend
447 DO i=istru,iend
448!^ ru(i,j,k,nrhs)=on_u(i,j)*0.5_r8* &
449!^ & (Hz(i,j,k)+Hz(i-1,j,k))* &
450!^ & (P(i-1,j,k)-P(i,j,k)- &
451!^ & HalfGRho* &
452!^ & ((rho(i,j,k)+rho(i-1,j,k))* &
453!^ & (z_r(i,j,k)-z_r(i-1,j,k))- &
454!^ & OneFifth* &
455!^ & ((dRx(i,j)-dRx(i-1,j))* &
456!^ & (z_r(i,j,k)-z_r(i-1,j,k)- &
457!^ & OneTwelfth* &
458!^ & (dZx(i,j)+dZx(i-1,j)))- &
459!^ & (dZx(i,j)-dZx(i-1,j))* &
460!^ & (rho(i,j,k)-rho(i-1,j,k)- &
461!^ & OneTwelfth* &
462!^ & (dRx(i,j)+dRx(i-1,j))))))
463!^
464 tl_ru(i,j,k,nrhs)=on_u(i,j)*0.5_r8* &
465 & ((tl_hz(i,j,k)+tl_hz(i-1,j,k))* &
466 & (p(i-1,j,k)-p(i,j,k)- &
467 & halfgrho* &
468 & ((rho(i,j,k)+rho(i-1,j,k))* &
469 & (z_r(i,j,k)-z_r(i-1,j,k))- &
470 & onefifth* &
471 & ((drx(i,j)-drx(i-1,j))* &
472 & (z_r(i,j,k)-z_r(i-1,j,k)- &
473 & onetwelfth* &
474 & (dzx(i,j)+dzx(i-1,j)))- &
475 & (dzx(i,j)-dzx(i-1,j))* &
476 & (rho(i,j,k)-rho(i-1,j,k)- &
477 & onetwelfth* &
478 & (drx(i,j)+drx(i-1,j))))))+ &
479 & (hz(i,j,k)+hz(i-1,j,k))* &
480 & (tl_p(i-1,j,k)-tl_p(i,j,k)- &
481 & halfgrho* &
482 & ((tl_rho(i,j,k)+tl_rho(i-1,j,k))* &
483 & (z_r(i,j,k)-z_r(i-1,j,k))+ &
484 & (rho(i,j,k)+rho(i-1,j,k))* &
485 & (tl_z_r(i,j,k)-tl_z_r(i-1,j,k))- &
486 & onefifth* &
487 & ((tl_drx(i,j)-tl_drx(i-1,j))* &
488 & (z_r(i,j,k)-z_r(i-1,j,k)- &
489 & onetwelfth* &
490 & (dzx(i,j)+dzx(i-1,j)))+ &
491 & (drx(i,j)-drx(i-1,j))* &
492 & (tl_z_r(i,j,k)-tl_z_r(i-1,j,k)- &
493 & onetwelfth* &
494 & (tl_dzx(i,j)+tl_dzx(i-1,j)))- &
495 & (tl_dzx(i,j)-tl_dzx(i-1,j))* &
496 & (rho(i,j,k)-rho(i-1,j,k)- &
497 & onetwelfth* &
498 & (drx(i,j)+drx(i-1,j)))- &
499 & (dzx(i,j)-dzx(i-1,j))* &
500 & (tl_rho(i,j,k)-tl_rho(i-1,j,k)- &
501 & onetwelfth* &
502 & (tl_drx(i,j)+tl_drx(i-1,j)))))))- &
503#ifdef TL_IOMS
504 & on_u(i,j)*0.5_r8* &
505 & (hz(i,j,k)+hz(i-1,j,k))* &
506 & (p(i-1,j,k)-p(i,j,k)- &
507 & 2.0_r8*halfgrho* &
508 & ((rho(i,j,k)+rho(i-1,j,k))* &
509 & (z_r(i,j,k)-z_r(i-1,j,k))- &
510 & onefifth* &
511 & ((drx(i,j)-drx(i-1,j))* &
512 & (z_r(i,j,k)-z_r(i-1,j,k)- &
513 & onetwelfth* &
514 & (dzx(i,j)+dzx(i-1,j)))- &
515 & (dzx(i,j)-dzx(i-1,j))* &
516 & (rho(i,j,k)-rho(i-1,j,k)- &
517 & onetwelfth* &
518 & (drx(i,j)+drx(i-1,j))))))
519#endif
520#ifdef DIAGNOSTICS_UV
521!! DiaRU(i,j,k,nrhs,M3pgrd)=ru(i,j,k,nrhs)
522#endif
523 END DO
524 END DO
525 END DO
526!
527!-----------------------------------------------------------------------
528! ETA-component pressure gradient term.
529!-----------------------------------------------------------------------
530!
531 DO k=n(ng),1,-1
532 DO j=jstrv-1,jend+1
533 DO i=istr,iend
534 aux(i,j)=z_r(i,j,k)-z_r(i,j-1,k)
535 tl_aux(i,j)=tl_z_r(i,j,k)-tl_z_r(i,j-1,k)
536#ifdef MASKING
537 aux(i,j)=aux(i,j)*vmask(i,j)
538 tl_aux(i,j)=tl_aux(i,j)*vmask(i,j)
539#endif
540 fc(i,j)=rho(i,j,k)-rho(i,j-1,k)
541 tl_fc(i,j)=tl_rho(i,j,k)-tl_rho(i,j-1,k)
542#ifdef MASKING
543 fc(i,j)=fc(i,j)*vmask(i,j)
544 tl_fc(i,j)=tl_fc(i,j)*vmask(i,j)
545#endif
546 END DO
547 END DO
548!
549 DO j=jstrv-1,jend
550 DO i=istr,iend
551 cff=2.0_r8*aux(i,j)*aux(i,j+1)
552 tl_cff=2.0_r8*(tl_aux(i,j)*aux(i,j+1)+ &
553 & aux(i,j)*tl_aux(i,j+1))- &
554#ifdef TL_IOMS
555 & cff
556#endif
557 IF (cff.gt.eps) THEN
558 cff1=1.0_r8/(aux(i,j)+aux(i,j+1))
559 tl_cff1=-cff1*cff1*(tl_aux(i,j)+tl_aux(i,j+1))+ &
560#ifdef TL_IOMS
561 & 2.0_r8*cff1
562#endif
563 dzx(i,j)=cff*cff1
564 tl_dzx(i,j)=tl_cff*cff1+cff*tl_cff1- &
565#ifdef TL_IOMS
566 & dzx(i,j)
567#endif
568 ELSE
569 dzx(i,j)=0.0_r8
570 tl_dzx(i,j)=0.0_r8
571 END IF
572 cff1=2.0_r8*fc(i,j)*fc(i,j+1)
573 tl_cff1=2.0_r8*(tl_fc(i,j)*fc(i,j+1)+ &
574 & fc(i,j)*tl_fc(i,j+1))- &
575#ifdef TL_IOMS
576 & cff1
577#endif
578 IF (cff1.gt.eps) THEN
579 cff2=1.0_r8/(fc(i,j)+fc(i,j+1))
580 tl_cff2=-cff2*cff2*(tl_fc(i,j)+tl_fc(i,j+1))+ &
581#ifdef TL_IOMS
582 & 2.0_r8*cff2
583#endif
584 drx(i,j)=cff1*cff2
585 tl_drx(i,j)=tl_cff1*cff2+cff1*tl_cff2- &
586#ifdef TL_IOMS
587 & drx(i,j)
588#endif
589 ELSE
590 drx(i,j)=0.0_r8
591 tl_drx(i,j)=0.0_r8
592 END IF
593 END DO
594 END DO
595!
596 DO j=jstrv,jend
597 DO i=istr,iend
598!^ rv(i,j,k,nrhs)=om_v(i,j)*0.5_r8* &
599!^ & (Hz(i,j,k)+Hz(i,j-1,k))* &
600!^ & (P(i,j-1,k)-P(i,j,k)- &
601!^ & HalfGRho* &
602!^ & ((rho(i,j,k)+rho(i,j-1,k))* &
603!^ & (z_r(i,j,k)-z_r(i,j-1,k))- &
604!^ & OneFifth* &
605!^ & ((dRx(i,j)-dRx(i,j-1))* &
606!^ & (z_r(i,j,k)-z_r(i,j-1,k)- &
607!^ & OneTwelfth* &
608!^ & (dZx(i,j)+dZx(i,j-1)))- &
609!^ & (dZx(i,j)-dZx(i,j-1))* &
610!^ & (rho(i,j,k)-rho(i,j-1,k)- &
611!^ & OneTwelfth* &
612!^ & (dRx(i,j)+dRx(i,j-1))))))
613!^
614 tl_rv(i,j,k,nrhs)=om_v(i,j)*0.5_r8* &
615 & ((tl_hz(i,j,k)+tl_hz(i,j-1,k))* &
616 & (p(i,j-1,k)-p(i,j,k)- &
617 & halfgrho* &
618 & ((rho(i,j,k)+rho(i,j-1,k))* &
619 & (z_r(i,j,k)-z_r(i,j-1,k))- &
620 & onefifth* &
621 & ((drx(i,j)-drx(i,j-1))* &
622 & (z_r(i,j,k)-z_r(i,j-1,k)- &
623 & onetwelfth* &
624 & (dzx(i,j)+dzx(i,j-1)))- &
625 & (dzx(i,j)-dzx(i,j-1))* &
626 & (rho(i,j,k)-rho(i,j-1,k)- &
627 & onetwelfth* &
628 & (drx(i,j)+drx(i,j-1))))))+ &
629 & (hz(i,j,k)+hz(i,j-1,k))* &
630 & (tl_p(i,j-1,k)-tl_p(i,j,k)- &
631 & halfgrho* &
632 & ((tl_rho(i,j,k)+tl_rho(i,j-1,k))* &
633 & (z_r(i,j,k)-z_r(i,j-1,k))+ &
634 & (rho(i,j,k)+rho(i,j-1,k))* &
635 & (tl_z_r(i,j,k)-tl_z_r(i,j-1,k))- &
636 & onefifth* &
637 & ((tl_drx(i,j)-tl_drx(i,j-1))* &
638 & (z_r(i,j,k)-z_r(i,j-1,k)- &
639 & onetwelfth* &
640 & (dzx(i,j)+dzx(i,j-1)))+ &
641 & (drx(i,j)-drx(i,j-1))* &
642 & (tl_z_r(i,j,k)-tl_z_r(i,j-1,k)- &
643 & onetwelfth* &
644 & (tl_dzx(i,j)+tl_dzx(i,j-1)))- &
645 & (tl_dzx(i,j)-tl_dzx(i,j-1))* &
646 & (rho(i,j,k)-rho(i,j-1,k)- &
647 & onetwelfth* &
648 & (drx(i,j)+drx(i,j-1)))- &
649 & (dzx(i,j)-dzx(i,j-1))* &
650 & (tl_rho(i,j,k)-tl_rho(i,j-1,k)- &
651 & onetwelfth* &
652 & (tl_drx(i,j)+tl_drx(i,j-1)))))))- &
653#ifdef TL_IOMS
654 & om_v(i,j)*0.5_r8* &
655 & (hz(i,j,k)+hz(i,j-1,k))* &
656 & (p(i,j-1,k)-p(i,j,k)- &
657 & 2.0_r8*halfgrho* &
658 & ((rho(i,j,k)+rho(i,j-1,k))* &
659 & (z_r(i,j,k)-z_r(i,j-1,k))- &
660 & onefifth* &
661 & ((drx(i,j)-drx(i,j-1))* &
662 & (z_r(i,j,k)-z_r(i,j-1,k)- &
663 & onetwelfth* &
664 & (dzx(i,j)+dzx(i,j-1)))- &
665 & (dzx(i,j)-dzx(i,j-1))* &
666 & (rho(i,j,k)-rho(i,j-1,k)- &
667 & onetwelfth* &
668 & (drx(i,j)+drx(i,j-1))))))
669#endif
670#ifdef DIAGNOSTICS_UV
671!! DiaRV(i,j,k,nrhs,M3pgrd)=rv(i,j,k,nrhs)
672#endif
673 END DO
674 END DO
675 END DO
676!
677 RETURN

References mod_scalars::g, and mod_scalars::rho0.

◆ rp_prsgrd40_tile()

subroutine rp_prsgrd_mod::rp_prsgrd40_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) imins,
integer, intent(in) imaxs,
integer, intent(in) jmins,
integer, intent(in) jmaxs,
integer, intent(in) nrhs,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) om_v,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) on_u,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) hz,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) tl_hz,
real(r8), dimension(lbi:ubi,lbj:ubj,0:n(ng)), intent(in) z_w,
real(r8), dimension(lbi:ubi,lbj:ubj,0:n(ng)), intent(in) tl_z_w,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) rho,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) tl_rho,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) eq_tide,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) tl_eq_tide,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) pair,
real(r8), dimension(lbi:ubi,lbj:ubj,0:n(ng),2), intent(inout) tl_ru,
real(r8), dimension(lbi:ubi,lbj:ubj,0:n(ng),2), intent(inout) tl_rv )
private

Definition at line 94 of file rp_prsgrd40.h.

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

References mod_scalars::g, and mod_scalars::rho0.