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

Functions/Subroutines

subroutine, public tl_prsgrd (ng, tile)
 
subroutine tl_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 tl_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 tl_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

◆ tl_prsgrd()

subroutine public tl_prsgrd_mod::tl_prsgrd ( integer, intent(in) ng,
integer, intent(in) tile )

Definition at line 34 of file tl_prsgrd31.h.

35!***********************************************************************
36!
37 USE mod_param
38#ifdef DIAGNOSTICS
39!! USE mod_diags
40#endif
41#ifdef ATM_PRESS
42 USE mod_forces
43#endif
44 USE mod_grid
45 USE mod_ocean
46 USE mod_stepping
47!
48! Imported variable declarations.
49!
50 integer, intent(in) :: ng, tile
51!
52! Local variable declarations.
53!
54 character (len=*), parameter :: MyFile = &
55 & __FILE__
56!
57#include "tile.h"
58!
59#ifdef PROFILE
60 CALL wclock_on (ng, itlm, 23, __line__, myfile)
61#endif
62 CALL tl_prsgrd31_tile (ng, tile, &
63 & lbi, ubi, lbj, ubj, &
64 & imins, imaxs, jmins, jmaxs, &
65 & nrhs(ng), &
66 & grid(ng) % om_v, &
67 & grid(ng) % on_u, &
68 & grid(ng) % Hz, &
69 & grid(ng) % tl_Hz, &
70 & grid(ng) % z_r, &
71 & grid(ng) % tl_z_r, &
72 & grid(ng) % z_w, &
73 & grid(ng) % tl_z_w, &
74 & ocean(ng) % rho, &
75 & ocean(ng) % tl_rho, &
76#ifdef TIDE_GENERATING_FORCES
77 & ocean(ng) % eq_tide, &
78 & ocean(ng) % tl_eq_tide, &
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) % tl_ru, &
88 & ocean(ng) % tl_rv)
89#ifdef PROFILE
90 CALL wclock_off (ng, itlm, 23, __line__, myfile)
91#endif
92!
93 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 itlm
Definition mod_param.F:663
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::itlm, mod_stepping::nrhs, mod_ocean::ocean, tl_prsgrd31_tile(), wclock_off(), and wclock_on().

Referenced by tl_rhs3d_mod::tl_rhs3d().

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

◆ tl_prsgrd31_tile()

subroutine tl_prsgrd_mod::tl_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 97 of file tl_prsgrd31.h.

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

Here is the caller graph for this function:

◆ tl_prsgrd32_tile()

subroutine tl_prsgrd_mod::tl_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 tl_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 dr1(i,k)=dr(i,k)
276 IF (cff.gt.eps) THEN
277 dr(i,k)=cff/(dr(i,k)+dr(i,k-1))
278 tl_dr(i,k)=(tl_cff-dr(i,k)*(tl_dr(i,k)+tl_dr(i,k-1)))/ &
279 & (dr1(i,k)+dr(i,k-1))
280 ELSE
281 dr(i,k)=0.0_r8
282 tl_dr(i,k)=0.0_r8
283 END IF
284 dz1(i,k)=dz(i,k)
285 dz(i,k)=2.0_r8*dz(i,k)*dz(i,k-1)/(dz(i,k)+dz(i,k-1))
286 tl_dz(i,k)=(2.0_r8*(tl_dz(i,k)*dz(i,k-1)+ &
287 & dz1(i,k)*tl_dz(i,k-1))- &
288 & dz(i,k)*(tl_dz(i,k)+tl_dz(i,k-1)))/ &
289 & (dz1(i,k)+dz(i,k-1))
290 END DO
291 END DO
292 DO i=istru-1,iend
293 cff1=1.0_r8/(z_r(i,j,n(ng))-z_r(i,j,n(ng)-1))
294 tl_cff1=-cff1*cff1*(tl_z_r(i,j,n(ng))-tl_z_r(i,j,n(ng)-1))
295 cff2=0.5_r8*(rho(i,j,n(ng))-rho(i,j,n(ng)-1))* &
296 & (z_w(i,j,n(ng))-z_r(i,j,n(ng)))*cff1
297 tl_cff2=0.5_r8*((tl_rho(i,j,n(ng))-tl_rho(i,j,n(ng)-1))* &
298 & (z_w(i,j,n(ng))-z_r(i,j,n(ng)))*cff1+ &
299 & (rho(i,j,n(ng))-rho(i,j,n(ng)-1))* &
300 & ((tl_z_w(i,j,n(ng))-tl_z_r(i,j,n(ng)))*cff1+ &
301 & (z_w(i,j,n(ng))-z_r(i,j,n(ng)))*tl_cff1))
302 p(i,j,n(ng))=g*z_w(i,j,n(ng))+ &
303#ifdef ATM_PRESS
304 & fac*(pair(i,j)-oneatm)+ &
305#endif
306 & grho*(rho(i,j,n(ng))+cff2)* &
307 & (z_w(i,j,n(ng))-z_r(i,j,n(ng)))
308 tl_p(i,j,n(ng))=g*tl_z_w(i,j,n(ng))+ &
309 & grho*((tl_rho(i,j,n(ng))+tl_cff2)* &
310 & (z_w(i,j,n(ng))-z_r(i,j,n(ng)))+ &
311 & (rho(i,j,n(ng))+cff2)* &
312 & (tl_z_w(i,j,n(ng))-tl_z_r(i,j,n(ng))))
313#ifdef TIDE_GENERATING_FORCES
314 p(i,j,n(ng))=p(i,j,n(ng))-g*eq_tide(i,j)
315 tl_p(i,j,n(ng))=tl_p(i,j,n(ng))-g*tl_eq_tide(i,j)
316#endif
317 END DO
318 DO k=n(ng)-1,1,-1
319 DO i=istru-1,iend
320 cff=halfgrho*((rho(i,j,k+1)+rho(i,j,k))* &
321 & (z_r(i,j,k+1)-z_r(i,j,k))- &
322 & onefifth* &
323 & ((dr(i,k+1)-dr(i,k))* &
324 & (z_r(i,j,k+1)-z_r(i,j,k)- &
325 & onetwelfth* &
326 & (dz(i,k+1)+dz(i,k)))- &
327 & (dz(i,k+1)-dz(i,k))* &
328 & (rho(i,j,k+1)-rho(i,j,k)- &
329 & onetwelfth* &
330 & (dr(i,k+1)+dr(i,k)))))
331 tl_cff=halfgrho*((tl_rho(i,j,k+1)+tl_rho(i,j,k))* &
332 & (z_r(i,j,k+1)-z_r(i,j,k))+ &
333 & (rho(i,j,k+1)+rho(i,j,k))* &
334 & (tl_z_r(i,j,k+1)-tl_z_r(i,j,k))- &
335 & onefifth* &
336 & ((tl_dr(i,k+1)-tl_dr(i,k))* &
337 & (z_r(i,j,k+1)-z_r(i,j,k)- &
338 & onetwelfth* &
339 & (dz(i,k+1)+dz(i,k)))+ &
340 & (dr(i,k+1)-dr(i,k))* &
341 & (tl_z_r(i,j,k+1)-tl_z_r(i,j,k)- &
342 & onetwelfth* &
343 & (tl_dz(i,k+1)+tl_dz(i,k)))- &
344 & (tl_dz(i,k+1)-tl_dz(i,k))* &
345 & (rho(i,j,k+1)-rho(i,j,k)- &
346 & onetwelfth* &
347 & (dr(i,k+1)+dr(i,k)))- &
348 & (dz(i,k+1)-dz(i,k))* &
349 & (tl_rho(i,j,k+1)-tl_rho(i,j,k)- &
350 & onetwelfth* &
351 & (tl_dr(i,k+1)+tl_dr(i,k)))))
352 p(i,j,k)=p(i,j,k+1)+cff
353 tl_p(i,j,k)=tl_p(i,j,k+1)+tl_cff
354 END DO
355 END DO
356 END DO
357!
358!-----------------------------------------------------------------------
359! Compute XI-component pressure gradient term.
360!-----------------------------------------------------------------------
361!
362 DO k=n(ng),1,-1
363 DO j=jstr,jend
364 DO i=istru-1,iend+1
365 aux(i,j)=z_r(i,j,k)-z_r(i-1,j,k)
366 tl_aux(i,j)=tl_z_r(i,j,k)-tl_z_r(i-1,j,k)
367#ifdef MASKING
368 aux(i,j)=aux(i,j)*umask(i,j)
369 tl_aux(i,j)=tl_aux(i,j)*umask(i,j)
370#endif
371 fc(i,j)=rho(i,j,k)-rho(i-1,j,k)
372 tl_fc(i,j)=tl_rho(i,j,k)-tl_rho(i-1,j,k)
373#ifdef MASKING
374 fc(i,j)=fc(i,j)*umask(i,j)
375 tl_fc(i,j)=tl_fc(i,j)*umask(i,j)
376#endif
377 END DO
378 END DO
379!
380 DO j=jstr,jend
381 DO i=istru-1,iend
382 cff=2.0_r8*aux(i,j)*aux(i+1,j)
383 tl_cff=2.0_r8*(tl_aux(i,j)*aux(i+1,j)+ &
384 & aux(i,j)*tl_aux(i+1,j))
385 IF (cff.gt.eps) THEN
386 cff1=1.0_r8/(aux(i,j)+aux(i+1,j))
387 tl_cff1=-cff1*cff1*(tl_aux(i,j)+tl_aux(i+1,j))
388 dzx(i,j)=cff*cff1
389 tl_dzx(i,j)=tl_cff*cff1+cff*tl_cff1
390 ELSE
391 dzx(i,j)=0.0_r8
392 tl_dzx(i,j)=0.0_r8
393 END IF
394 cff1=2.0_r8*fc(i,j)*fc(i+1,j)
395 tl_cff1=2.0_r8*(tl_fc(i,j)*fc(i+1,j)+ &
396 & fc(i,j)*tl_fc(i+1,j))
397 IF (cff1.gt.eps) THEN
398 cff2=1.0_r8/(fc(i,j)+fc(i+1,j))
399 tl_cff2=-cff2*cff2*(tl_fc(i,j)+tl_fc(i+1,j))
400 drx(i,j)=cff1*cff2
401 tl_drx(i,j)=tl_cff1*cff2+cff1*tl_cff2
402 ELSE
403 drx(i,j)=0.0_r8
404 tl_drx(i,j)=0.0_r8
405 END IF
406 END DO
407 END DO
408!
409 DO j=jstr,jend
410 DO i=istru,iend
411!^ ru(i,j,k,nrhs)=on_u(i,j)*0.5_r8* &
412!^ & (Hz(i,j,k)+Hz(i-1,j,k))* &
413!^ & (P(i-1,j,k)-P(i,j,k)- &
414!^ & HalfGRho* &
415!^ & ((rho(i,j,k)+rho(i-1,j,k))* &
416!^ & (z_r(i,j,k)-z_r(i-1,j,k))- &
417!^ & OneFifth* &
418!^ & ((dRx(i,j)-dRx(i-1,j))* &
419!^ & (z_r(i,j,k)-z_r(i-1,j,k)- &
420!^ & OneTwelfth* &
421!^ & (dZx(i,j)+dZx(i-1,j)))- &
422!^ & (dZx(i,j)-dZx(i-1,j))* &
423!^ & (rho(i,j,k)-rho(i-1,j,k)- &
424!^ & OneTwelfth* &
425!^ & (dRx(i,j)+dRx(i-1,j))))))
426!^
427 tl_ru(i,j,k,nrhs)=on_u(i,j)*0.5_r8* &
428 & ((tl_hz(i,j,k)+tl_hz(i-1,j,k))* &
429 & (p(i-1,j,k)-p(i,j,k)- &
430 & halfgrho* &
431 & ((rho(i,j,k)+rho(i-1,j,k))* &
432 & (z_r(i,j,k)-z_r(i-1,j,k))- &
433 & onefifth* &
434 & ((drx(i,j)-drx(i-1,j))* &
435 & (z_r(i,j,k)-z_r(i-1,j,k)- &
436 & onetwelfth* &
437 & (dzx(i,j)+dzx(i-1,j)))- &
438 & (dzx(i,j)-dzx(i-1,j))* &
439 & (rho(i,j,k)-rho(i-1,j,k)- &
440 & onetwelfth* &
441 & (drx(i,j)+drx(i-1,j))))))+ &
442 & (hz(i,j,k)+hz(i-1,j,k))* &
443 & (tl_p(i-1,j,k)-tl_p(i,j,k)- &
444 & halfgrho* &
445 & ((tl_rho(i,j,k)+tl_rho(i-1,j,k))* &
446 & (z_r(i,j,k)-z_r(i-1,j,k))+ &
447 & (rho(i,j,k)+rho(i-1,j,k))* &
448 & (tl_z_r(i,j,k)-tl_z_r(i-1,j,k))- &
449 & onefifth* &
450 & ((tl_drx(i,j)-tl_drx(i-1,j))* &
451 & (z_r(i,j,k)-z_r(i-1,j,k)- &
452 & onetwelfth* &
453 & (dzx(i,j)+dzx(i-1,j)))+ &
454 & (drx(i,j)-drx(i-1,j))* &
455 & (tl_z_r(i,j,k)-tl_z_r(i-1,j,k)- &
456 & onetwelfth* &
457 & (tl_dzx(i,j)+tl_dzx(i-1,j)))- &
458 & (tl_dzx(i,j)-tl_dzx(i-1,j))* &
459 & (rho(i,j,k)-rho(i-1,j,k)- &
460 & onetwelfth* &
461 & (drx(i,j)+drx(i-1,j)))- &
462 & (dzx(i,j)-dzx(i-1,j))* &
463 & (tl_rho(i,j,k)-tl_rho(i-1,j,k)- &
464 & onetwelfth* &
465 & (tl_drx(i,j)+tl_drx(i-1,j)))))))
466#ifdef DIAGNOSTICS_UV
467!! DiaRU(i,j,k,nrhs,M3pgrd)=ru(i,j,k,nrhs)
468#endif
469 END DO
470 END DO
471 END DO
472!
473!-----------------------------------------------------------------------
474! ETA-component pressure gradient term.
475!-----------------------------------------------------------------------
476!
477 DO k=n(ng),1,-1
478 DO j=jstrv-1,jend+1
479 DO i=istr,iend
480 aux(i,j)=z_r(i,j,k)-z_r(i,j-1,k)
481 tl_aux(i,j)=tl_z_r(i,j,k)-tl_z_r(i,j-1,k)
482#ifdef MASKING
483 aux(i,j)=aux(i,j)*vmask(i,j)
484 tl_aux(i,j)=tl_aux(i,j)*vmask(i,j)
485#endif
486 fc(i,j)=rho(i,j,k)-rho(i,j-1,k)
487 tl_fc(i,j)=tl_rho(i,j,k)-tl_rho(i,j-1,k)
488#ifdef MASKING
489 fc(i,j)=fc(i,j)*vmask(i,j)
490 tl_fc(i,j)=tl_fc(i,j)*vmask(i,j)
491#endif
492 END DO
493 END DO
494!
495 DO j=jstrv-1,jend
496 DO i=istr,iend
497 cff=2.0_r8*aux(i,j)*aux(i,j+1)
498 tl_cff=2.0_r8*(tl_aux(i,j)*aux(i,j+1)+ &
499 & aux(i,j)*tl_aux(i,j+1))
500 IF (cff.gt.eps) THEN
501 cff1=1.0_r8/(aux(i,j)+aux(i,j+1))
502 tl_cff1=-cff1*cff1*(tl_aux(i,j)+tl_aux(i,j+1))
503 dzx(i,j)=cff*cff1
504 tl_dzx(i,j)=tl_cff*cff1+cff*tl_cff1
505 ELSE
506 dzx(i,j)=0.0_r8
507 tl_dzx(i,j)=0.0_r8
508 END IF
509 cff1=2.0_r8*fc(i,j)*fc(i,j+1)
510 tl_cff1=2.0_r8*(tl_fc(i,j)*fc(i,j+1)+ &
511 & fc(i,j)*tl_fc(i,j+1))
512 IF (cff1.gt.eps) THEN
513 cff2=1.0_r8/(fc(i,j)+fc(i,j+1))
514 tl_cff2=-cff2*cff2*(tl_fc(i,j)+tl_fc(i,j+1))
515 drx(i,j)=cff1*cff2
516 tl_drx(i,j)=tl_cff1*cff2+cff1*tl_cff2
517 ELSE
518 drx(i,j)=0.0_r8
519 tl_drx(i,j)=0.0_r8
520 END IF
521 END DO
522 END DO
523!
524 DO j=jstrv,jend
525 DO i=istr,iend
526!^ rv(i,j,k,nrhs)=om_v(i,j)*0.5_r8* &
527!^ & (Hz(i,j,k)+Hz(i,j-1,k))* &
528!^ & (P(i,j-1,k)-P(i,j,k)- &
529!^ & HalfGRho* &
530!^ & ((rho(i,j,k)+rho(i,j-1,k))* &
531!^ & (z_r(i,j,k)-z_r(i,j-1,k))- &
532!^ & OneFifth* &
533!^ & ((dRx(i,j)-dRx(i,j-1))* &
534!^ & (z_r(i,j,k)-z_r(i,j-1,k)- &
535!^ & OneTwelfth* &
536!^ & (dZx(i,j)+dZx(i,j-1)))- &
537!^ & (dZx(i,j)-dZx(i,j-1))* &
538!^ & (rho(i,j,k)-rho(i,j-1,k)- &
539!^ & OneTwelfth* &
540!^ & (dRx(i,j)+dRx(i,j-1))))))
541!^
542 tl_rv(i,j,k,nrhs)=om_v(i,j)*0.5_r8* &
543 & ((tl_hz(i,j,k)+tl_hz(i,j-1,k))* &
544 & (p(i,j-1,k)-p(i,j,k)- &
545 & halfgrho* &
546 & ((rho(i,j,k)+rho(i,j-1,k))* &
547 & (z_r(i,j,k)-z_r(i,j-1,k))- &
548 & onefifth* &
549 & ((drx(i,j)-drx(i,j-1))* &
550 & (z_r(i,j,k)-z_r(i,j-1,k)- &
551 & onetwelfth* &
552 & (dzx(i,j)+dzx(i,j-1)))- &
553 & (dzx(i,j)-dzx(i,j-1))* &
554 & (rho(i,j,k)-rho(i,j-1,k)- &
555 & onetwelfth* &
556 & (drx(i,j)+drx(i,j-1))))))+ &
557 & (hz(i,j,k)+hz(i,j-1,k))* &
558 & (tl_p(i,j-1,k)-tl_p(i,j,k)- &
559 & halfgrho* &
560 & ((tl_rho(i,j,k)+tl_rho(i,j-1,k))* &
561 & (z_r(i,j,k)-z_r(i,j-1,k))+ &
562 & (rho(i,j,k)+rho(i,j-1,k))* &
563 & (tl_z_r(i,j,k)-tl_z_r(i,j-1,k))- &
564 & onefifth* &
565 & ((tl_drx(i,j)-tl_drx(i,j-1))* &
566 & (z_r(i,j,k)-z_r(i,j-1,k)- &
567 & onetwelfth* &
568 & (dzx(i,j)+dzx(i,j-1)))+ &
569 & (drx(i,j)-drx(i,j-1))* &
570 & (tl_z_r(i,j,k)-tl_z_r(i,j-1,k)- &
571 & onetwelfth* &
572 & (tl_dzx(i,j)+tl_dzx(i,j-1)))- &
573 & (tl_dzx(i,j)-tl_dzx(i,j-1))* &
574 & (rho(i,j,k)-rho(i,j-1,k)- &
575 & onetwelfth* &
576 & (drx(i,j)+drx(i,j-1)))- &
577 & (dzx(i,j)-dzx(i,j-1))* &
578 & (tl_rho(i,j,k)-tl_rho(i,j-1,k)- &
579 & onetwelfth* &
580 & (tl_drx(i,j)+tl_drx(i,j-1)))))))
581#ifdef DIAGNOSTICS_UV
582!! DiaRV(i,j,k,nrhs,M3pgrd)=rv(i,j,k,nrhs)
583#endif
584 END DO
585 END DO
586 END DO
587!
588 RETURN

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

◆ tl_prsgrd40_tile()

subroutine tl_prsgrd_mod::tl_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 tl_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.