ROMS
Loading...
Searching...
No Matches
tl_prsgrd31.h
Go to the documentation of this file.
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 evalutes the tangent linear baroclinic hydrostatic !
11! pressure gradient term using the standard and weighted Jacobians !
12! scheme of Song and Wright (1997). Note that horizontal gradients !
13! are computed before of the vertical integration. !
14! !
15! The pressure gradient terms (m4/s2) are loaded into right-hand- !
16! side arrays "tl_ru" and "tl_rv". !
17! !
18! Reference: !
19! !
20! Song, Y.T., 1998: A general pressure gradient formulation for !
21! numerical ocean models. Part I: Scheme design and diagnostic !
22! analysis, Monthly Weather Rev., 126, 3213-3230. !
23! !
24!=======================================================================
25!
26 implicit none
27!
28 PRIVATE
29 PUBLIC :: tl_prsgrd
30!
31 CONTAINS
32!
33!***********************************************************************
34 SUBROUTINE tl_prsgrd (ng, tile)
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
94 END SUBROUTINE tl_prsgrd
95!
96!***********************************************************************
97 SUBROUTINE tl_prsgrd31_tile (ng, tile, &
98 & LBi, UBi, LBj, UBj, &
99 & IminS, ImaxS, JminS, JmaxS, &
100 & nrhs, &
101 & om_v, on_u, &
102 & Hz, tl_Hz, &
103 & z_r, tl_z_r, &
104 & z_w, tl_z_w, &
105 & rho, tl_rho, &
106#ifdef TIDE_GENERATING_FORCES
107 & eq_tide, tl_eq_tide, &
108#endif
109#ifdef ATM_PRESS
110 & Pair, &
111#endif
112#ifdef DIAGNOSTICS_UV
113!! & DiaRU, DiaRV, &
114#endif
115 & tl_ru, tl_rv)
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
503 END SUBROUTINE tl_prsgrd31_tile
504
505 END MODULE tl_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 itlm
Definition mod_param.F:663
real(dp) g
real(dp) rho0
integer, dimension(:), allocatable nrhs
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, public tl_prsgrd(ng, tile)
Definition tl_prsgrd31.h:35
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