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

Functions/Subroutines

subroutine, public ad_prsgrd (ng, tile)
 
subroutine ad_prsgrd31_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, om_v, on_u, hz, ad_hz, z_r, ad_z_r, z_w, ad_z_w, rho, ad_rho, eq_tide, ad_eq_tide, pair, ad_ru, ad_rv)
 
subroutine ad_prsgrd32_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, umask, vmask, om_v, on_u, hz, ad_hz, z_r, ad_z_r, z_w, ad_z_w, rho, ad_rho, eq_tide, ad_eq_tide, pair, ad_ru, ad_rv)
 
subroutine ad_prsgrd40_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, om_v, on_u, hz, ad_hz, z_w, ad_z_w, rho, ad_rho, eq_tide, ad_eq_tide, pair, ad_ru, ad_rv)
 

Function/Subroutine Documentation

◆ ad_prsgrd()

subroutine public ad_prsgrd_mod::ad_prsgrd ( integer, intent(in) ng,
integer, intent(in) tile )

Definition at line 36 of file ad_prsgrd31.h.

37!***********************************************************************
38!
39 USE mod_param
40#ifdef DIAGNOSTICS
41!! USE mod_diags
42#endif
43#ifdef ATM_PRESS
44 USE mod_forces
45#endif
46 USE mod_grid
47 USE mod_ocean
48 USE mod_stepping
49!
50! Imported variable declarations.
51!
52 integer, intent(in) :: ng, tile
53!
54! Local variable declarations.
55!
56 character (len=*), parameter :: MyFile = &
57 & __FILE__
58!
59#include "tile.h"
60!
61#ifdef PROFILE
62 CALL wclock_on (ng, iadm, 23, __line__, myfile)
63#endif
64 CALL ad_prsgrd31_tile (ng, tile, &
65 & lbi, ubi, lbj, ubj, &
66 & imins, imaxs, jmins, jmaxs, &
67 & nrhs(ng), &
68 & grid(ng) % om_v, &
69 & grid(ng) % on_u, &
70 & grid(ng) % Hz, &
71 & grid(ng) % ad_Hz, &
72 & grid(ng) % z_r, &
73 & grid(ng) % ad_z_r, &
74 & grid(ng) % z_w, &
75 & grid(ng) % ad_z_w, &
76 & ocean(ng) % rho, &
77 & ocean(ng) % ad_rho, &
78#ifdef TIDE_GENERATING_FORCES
79 & ocean(ng) % eq_tide, &
80 & ocean(ng) % ad_eq_tide, &
81#endif
82#ifdef ATM_PRESS
83 & forces(ng) % Pair, &
84#endif
85#ifdef DIAGNOSTICS_UV
86!! & DIAGS(ng) % DiaRU, &
87!! & DIAGS(ng) % DiaRV, &
88#endif
89 & ocean(ng) % ad_ru, &
90 & ocean(ng) % ad_rv)
91#ifdef PROFILE
92 CALL wclock_off (ng, iadm, 23, __line__, myfile)
93#endif
94!
95 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 iadm
Definition mod_param.F:665
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 ad_prsgrd31_tile(), mod_forces::forces, mod_grid::grid, mod_param::iadm, mod_stepping::nrhs, mod_ocean::ocean, wclock_off(), and wclock_on().

Referenced by ad_rhs3d_mod::ad_rhs3d().

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

◆ ad_prsgrd31_tile()

subroutine ad_prsgrd_mod::ad_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(inout) ad_hz,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) z_r,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(inout) ad_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(inout) ad_z_w,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) rho,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(inout) ad_rho,
real(r8), dimension(lbi:,lbj:), intent(in) eq_tide,
real(r8), dimension(lbi:,lbj:), intent(out) ad_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) ad_ru,
real(r8), dimension(lbi:ubi,lbj:ubj,0:n(ng),2), intent(inout) ad_rv )
private

Definition at line 99 of file ad_prsgrd31.h.

118!***********************************************************************
119!
120 USE mod_param
121 USE mod_scalars
122!
123! Imported variable declarations.
124!
125 integer, intent(in) :: ng, tile
126 integer, intent(in) :: LBi, UBi, LBj, UBj
127 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
128 integer, intent(in) :: nrhs
129
130#ifdef ASSUMED_SHAPE
131 real(r8), intent(in) :: om_v(LBi:,LBj:)
132 real(r8), intent(in) :: on_u(LBi:,LBj:)
133 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
134 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
135 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
136 real(r8), intent(in) :: rho(LBi:,LBj:,:)
137# ifdef ATM_PRESS
138 real(r8), intent(in) :: Pair(LBi:,LBj:)
139# endif
140# ifdef TIDE_GENERATING_FORCES
141 real(r8), intent(in) :: eq_tide(LBi:,LBj:)
142 real(r8), intent(inout) :: ad_eq_tide(LBi:,LBj:)
143# endif
144# ifdef DIAGNOSTICS_UV
145!! real(r8), intent(inout) :: DiaRU(LBi:,LBj:,:,:,:)
146!! real(r8), intent(inout) :: DiaRV(LBi:,LBj:,:,:,:)
147# endif
148 real(r8), intent(inout) :: ad_Hz(LBi:,LBj:,:)
149 real(r8), intent(inout) :: ad_z_r(LBi:,LBj:,:)
150 real(r8), intent(inout) :: ad_z_w(LBi:,LBj:,0:)
151 real(r8), intent(inout) :: ad_rho(LBi:,LBj:,:)
152 real(r8), intent(inout) :: ad_ru(LBi:,LBj:,0:,:)
153 real(r8), intent(inout) :: ad_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# ifdef ATM_PRESS
162 real(r8), intent(in) :: Pair(LBi:UBi,LBj:UBj)
163# endif
164# ifdef TIDE_GENERATING_FORCES
165 real(r8), intent(in) :: eq_tide(LBi:,LBj:)
166 real(r8), intent(out) :: ad_eq_tide(LBi:,LBj:)
167# endif
168# ifdef DIAGNOSTICS_UV
169!! real(r8), intent(inout) :: DiaRU(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
170!! real(r8), intent(inout) :: DiaRV(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
171# endif
172 real(r8), intent(inout) :: ad_Hz(LBi:UBi,LBj:UBj,N(ng))
173 real(r8), intent(inout) :: ad_z_r(LBi:UBi,LBj:UBj,N(ng))
174 real(r8), intent(inout) :: ad_z_w(LBi:UBi,LBj:UBj,0:N(ng))
175 real(r8), intent(inout) :: ad_rho(LBi:UBi,LBj:UBj,N(ng))
176 real(r8), intent(inout) :: ad_ru(LBi:UBi,LBj:UBj,0:N(ng),2)
177 real(r8), intent(inout) :: ad_rv(LBi:UBi,LBj:UBj,0:N(ng),2)
178#endif
179!
180! Local variable declarations.
181!
182 integer :: i, j, k, kk
183
184 real(r8) :: fac, fac1, fac2, fac3
185 real(r8) :: cff1, cff2, cff3, cff4
186 real(r8) :: adfac, adfac1, adfac2
187 real(r8) :: ad_cff1, ad_cff2, ad_cff3, ad_cff4
188#ifdef WJ_GRADP
189 real(r8) :: gamma, ad_gamma
190#endif
191
192 real(r8), dimension(IminS:ImaxS) :: phie
193 real(r8), dimension(IminS:ImaxS) :: phix
194
195 real(r8), dimension(IminS:ImaxS) :: ad_phie
196 real(r8), dimension(IminS:ImaxS) :: ad_phix
197
198#include "set_bounds.h"
199!
200!-----------------------------------------------------------------------
201! Initialize adjoint private variables.
202!-----------------------------------------------------------------------
203!
204 ad_cff1=0.0_r8
205 ad_cff2=0.0_r8
206 ad_cff3=0.0_r8
207 ad_cff4=0.0_r8
208
209#ifdef WJ_GRADP
210 ad_gamma=0.0_r8
211#endif
212 DO i=imins,imaxs
213 ad_phie(i)=0.0_r8
214 ad_phix(i)=0.0_r8
215 END DO
216!
217!-----------------------------------------------------------------------
218! Calculate adjoint pressure gradient in the ETA-direction (m4/s2).
219!-----------------------------------------------------------------------
220!
221#ifdef ATM_PRESS
222 fac=100.0_r8/rho0
223#endif
224 fac1=0.5_r8*g/rho0
225 fac2=1000.0_r8*g/rho0
226 fac3=0.25_r8*g/rho0
227
228 j_loop : DO j=jstr,jend
229 IF (j.ge.jstrv) THEN
230!
231! Compute appropriate BASIC STATE "phie". Notice that a reverse
232! vertical integration of "phie" is carried out over kk-index.
233!
234!^ DO k=N(ng)-1,1,-1
235!^
236 DO k=1,n(ng)-1
237 DO i=istr,iend
238#ifdef TIDE_GENERATING_FORCES
239 cff1=z_w(i,j ,n(ng))-eq_tide(i,j )-z_r(i,j ,n(ng))+ &
240 & z_w(i,j-1,n(ng))-eq_tide(i,j-1)-z_r(i,j-1,n(ng))
241#else
242 cff1=z_w(i,j ,n(ng))-z_r(i,j ,n(ng))+ &
243 & z_w(i,j-1,n(ng))-z_r(i,j-1,n(ng))
244#endif
245 phie(i)=fac1*(rho(i,j,n(ng))-rho(i,j-1,n(ng)))*cff1
246#ifdef ATM_PRESS
247 phie(i)=phie(i)+fac*(pair(i,j)-pair(i,j-1))
248#endif
249#ifdef RHO_SURF
250 phie(i)=phie(i)+ &
251 & (fac2+fac1*(rho(i,j,n(ng))+rho(i,j-1,n(ng))))* &
252 & (z_w(i,j,n(ng))-z_w(i,j-1,n(ng)))
253#endif
254 END DO
255 DO kk=n(ng)-1,k,-1
256 DO i=istr,iend
257#ifdef WJ_GRADP
258 cff1=1.0_r8/((z_r(i,j ,kk+1)-z_r(i,j ,kk))* &
259 & (z_r(i,j-1,kk+1)-z_r(i,j-1,kk)))
260 cff2=z_r(i,j ,kk )-z_r(i,j-1,kk )+ &
261 & z_r(i,j ,kk+1)-z_r(i,j-1,kk+1)
262 cff3=z_r(i,j ,kk+1)-z_r(i,j ,kk )- &
263 & z_r(i,j-1,kk+1)+z_r(i,j-1,kk )
264 gamma=0.125_r8*cff1*cff2*cff3
265
266 cff1=(1.0_r8+gamma)*(rho(i,j,kk+1)-rho(i,j-1,kk+1))+ &
267 & (1.0_r8-gamma)*(rho(i,j,kk )-rho(i,j-1,kk ))
268 cff2=rho(i,j,kk+1)+rho(i,j-1,kk+1)- &
269 & rho(i,j,kk )-rho(i,j-1,kk )
270 cff3=z_r(i,j,kk+1)+z_r(i,j-1,kk+1)- &
271 & z_r(i,j,kk )-z_r(i,j-1,kk )
272 cff4=(1.0_r8+gamma)*(z_r(i,j,kk+1)-z_r(i,j-1,kk+1))+ &
273 & (1.0_r8-gamma)*(z_r(i,j,kk )-z_r(i,j-1,kk ))
274 phie(i)=phie(i)+ &
275 & fac3*(cff1*cff3-cff2*cff4)
276#else
277 cff1=rho(i,j,kk+1)-rho(i,j-1,kk+1)+ &
278 & rho(i,j,kk )-rho(i,j-1,kk )
279 cff2=rho(i,j,kk+1)+rho(i,j-1,kk+1)- &
280 & rho(i,j,kk )-rho(i,j-1,kk )
281 cff3=z_r(i,j,kk+1)+z_r(i,j-1,kk+1)- &
282 & z_r(i,j,kk )-z_r(i,j-1,kk )
283 cff4=z_r(i,j,kk+1)-z_r(i,j-1,kk+1)+ &
284 & z_r(i,j,kk )-z_r(i,j-1,kk )
285 phie(i)=phie(i)+ &
286 & fac3*(cff1*cff3-cff2*cff4)
287#endif
288 END DO
289 END DO
290!
291! Compute interior adjoint baroclinic pressure gradient. Differentiate
292! and then vertically integrate.
293!
294 DO i=istr,iend
295#ifdef DIAGNOSTICS_UV
296!! DiaRV(i,j,k,nrhs,M3pgrd)=rv(i,j,k,nrhs)
297#endif
298!^ tl_rv(i,j,k,nrhs)=-0.5_r8*om_v(i,j)* &
299!^ & ((tl_Hz(i,j,k)+tl_Hz(i,j-1,k))* &
300!^ & phie(i)+ &
301!^ & (Hz(i,j,k)+Hz(i,j-1,k))* &
302!^ & tl_phie(i))
303!^
304 adfac=-0.5_r8*om_v(i,j)*ad_rv(i,j,k,nrhs)
305 adfac1=adfac*phie(i)
306 ad_phie(i)=ad_phie(i)+ &
307 & (hz(i,j,k)+hz(i,j-1,k))*adfac
308 ad_hz(i,j-1,k)=ad_hz(i,j-1,k)+adfac1
309 ad_hz(i,j ,k)=ad_hz(i,j ,k)+adfac1
310 ad_rv(i,j,k,nrhs)=0.0_r8
311#ifdef WJ_GRADP
312 cff1=1.0_r8/((z_r(i,j ,k+1)-z_r(i,j ,k))* &
313 & (z_r(i,j-1,k+1)-z_r(i,j-1,k)))
314 cff2=z_r(i,j ,k )-z_r(i,j-1,k )+ &
315 & z_r(i,j ,k+1)-z_r(i,j-1,k+1)
316 cff3=z_r(i,j ,k+1)-z_r(i,j ,k )- &
317 & z_r(i,j-1,k+1)+z_r(i,j-1,k )
318 gamma=0.125_r8*cff1*cff2*cff3
319
320 cff1=(1.0_r8+gamma)*(rho(i,j,k+1)-rho(i,j-1,k+1))+ &
321 & (1.0_r8-gamma)*(rho(i,j,k )-rho(i,j-1,k ))
322 cff2=rho(i,j,k+1)+rho(i,j-1,k+1)- &
323 & rho(i,j,k )-rho(i,j-1,k )
324 cff3=z_r(i,j,k+1)+z_r(i,j-1,k+1)- &
325 & z_r(i,j,k )-z_r(i,j-1,k )
326 cff4=(1.0_r8+gamma)*(z_r(i,j,k+1)-z_r(i,j-1,k+1))+ &
327 & (1.0_r8-gamma)*(z_r(i,j,k )-z_r(i,j-1,k ))
328!^ tl_phie(i)=tl_phie(i)+ &
329!^ & fac3*(tl_cff1*cff3+ &
330!^ & cff1*tl_cff3- &
331!^ & tl_cff2*cff4- &
332!^ & cff2*tl_cff4)
333!^
334 adfac=fac3*ad_phie(i)
335 ad_cff1=ad_cff1+cff3*adfac
336 ad_cff2=ad_cff2-cff4*adfac
337 ad_cff3=ad_cff3+cff1*adfac
338 ad_cff4=ad_cff4-cff2*adfac
339!^ tl_cff4=tl_gamma*(z_r(i,j,k+1)-z_r(i,j-1,k+1)- &
340!^ & z_r(i,j,k )+z_r(i,j-1,k ))+ &
341!^ & (1.0_r8+gamma)*(tl_z_r(i,j ,k+1)- &
342!^ & tl_z_r(i,j-1,k+1))+ &
343!^ & (1.0_r8-gamma)*(tl_z_r(i,j ,k )- &
344!^ & tl_z_r(i,j-1,k ))
345!^ tl_cff3=tl_z_r(i,j,k+1)+tl_z_r(i,j-1,k+1)- &
346!^ & tl_z_r(i,j,k )-tl_z_r(i,j-1,k )
347!^
348 adfac1=(1.0_r8+gamma)*ad_cff4
349 adfac2=(1.0_r8-gamma)*ad_cff4
350 ad_z_r(i,j-1,k )=ad_z_r(i,j-1,k )-adfac2-ad_cff3
351 ad_z_r(i,j ,k )=ad_z_r(i,j ,k )+adfac2-ad_cff3
352 ad_z_r(i,j-1,k+1)=ad_z_r(i,j-1,k+1)-adfac1+ad_cff3
353 ad_z_r(i,j ,k+1)=ad_z_r(i,j ,k+1)+adfac1+ad_cff3
354 ad_gamma=ad_gamma+ &
355 & (z_r(i,j,k+1)-z_r(i,j-1,k+1)- &
356 & z_r(i,j,k )+z_r(i,j-1,k ))*ad_cff4
357 ad_cff4=0.0_r8
358 ad_cff3=0.0_r8
359!^ tl_cff2=tl_rho(i,j,k+1)+tl_rho(i,j-1,k+1)- &
360!^ & tl_rho(i,j,k )-tl_rho(i,j-1,k )
361!^ tl_cff1=tl_gamma*(rho(i,j,k+1)-rho(i,j-1,k+1)- &
362!^ & rho(i,j,k )+rho(i,j-1,k ))+ &
363!^ & (1.0_r8+gamma)*(tl_rho(i,j ,k+1)- &
364!^ & tl_rho(i,j-1,k+1))+ &
365!^ & (1.0_r8-gamma)*(tl_rho(i,j ,k )- &
366!^ & tl_rho(i,j-1,k ))
367!^
368 adfac1=(1.0_r8+gamma)*ad_cff1
369 adfac2=(1.0_r8-gamma)*ad_cff1
370 ad_rho(i,j-1,k )=ad_rho(i,j-1,k )-adfac2-ad_cff2
371 ad_rho(i,j ,k )=ad_rho(i,j ,k )+adfac2-ad_cff2
372 ad_rho(i,j-1,k+1)=ad_rho(i,j-1,k+1)-adfac1+ad_cff2
373 ad_rho(i,j ,k+1)=ad_rho(i,j ,k+1)+adfac1+ad_cff2
374 ad_gamma=ad_gamma+ &
375 & (rho(i,j,k+1)-rho(i,j-1,k+1)- &
376 & rho(i,j,k )+rho(i,j-1,k ))*ad_cff1
377 ad_cff2=0.0_r8
378 ad_cff1=0.0_r8
379!
380 cff1=1.0_r8/((z_r(i,j ,k+1)-z_r(i,j ,k))* &
381 & (z_r(i,j-1,k+1)-z_r(i,j-1,k)))
382 cff2=z_r(i,j ,k )-z_r(i,j-1,k )+ &
383 & z_r(i,j ,k+1)-z_r(i,j-1,k+1)
384 cff3=z_r(i,j ,k+1)-z_r(i,j ,k )- &
385 & z_r(i,j-1,k+1)+z_r(i,j-1,k )
386
387!^ tl_gamma=0.125_r8*(tl_cff1*cff2*cff3+ &
388!^ & cff1*(tl_cff2*cff3+ &
389!^ & cff2*tl_cff3))
390!^
391 adfac=0.125_r8*ad_gamma
392 adfac1=adfac*cff1
393 ad_cff3=ad_cff3+cff2*adfac1
394 ad_cff2=ad_cff2+cff3*adfac1
395 ad_cff1=ad_cff1+cff2*cff3*adfac
396 ad_gamma=0.0_r8
397!^ tl_cff3=tl_z_r(i,j ,k+1)-tl_z_r(i,j ,k )- &
398!^ & tl_z_r(i,j-1,k+1)+tl_z_r(i,j-1,k )
399!^ tl_cff2=tl_z_r(i,j ,k )-tl_z_r(i,j-1,k )+ &
400!^ & tl_z_r(i,j ,k+1)-tl_z_r(i,j-1,k+1)
401!^
402 ad_z_r(i,j-1,k )=ad_z_r(i,j-1,k )-ad_cff2+ad_cff3
403 ad_z_r(i,j ,k )=ad_z_r(i,j ,k )+ad_cff2-ad_cff3
404 ad_z_r(i,j-1,k+1)=ad_z_r(i,j-1,k+1)-ad_cff2-ad_cff3
405 ad_z_r(i,j ,k+1)=ad_z_r(i,j ,k+1)+ad_cff2+ad_cff3
406 ad_cff3=0.0_r8
407 ad_cff2=0.0_r8
408!^ tl_cff1=-cff1*cff1*((tl_z_r(i,j ,k+1)-tl_z_r(i,j ,k))* &
409!^ & (z_r(i,j-1,k+1)-z_r(i,j-1,k))+ &
410!^ & (z_r(i,j ,k+1)-z_r(i,j ,k))* &
411!^ & (tl_z_r(i,j-1,k+1)-tl_z_r(i,j-1,k)))
412!^
413 adfac=-cff1*cff1*ad_cff1
414 adfac1=adfac*(z_r(i,j-1,k+1)-z_r(i,j-1,k))
415 adfac2=adfac*(z_r(i,j ,k+1)-z_r(i,j ,k))
416 ad_z_r(i,j-1,k )=ad_z_r(i,j-1,k )-adfac2
417 ad_z_r(i,j ,k )=ad_z_r(i,j ,k )-adfac1
418 ad_z_r(i,j-1,k+1)=ad_z_r(i,j-1,k+1)+adfac2
419 ad_z_r(i,j ,k+1)=ad_z_r(i,j ,k+1)+adfac1
420 ad_cff1=0.0_r8
421#else
422!
423 cff1=rho(i,j,k+1)-rho(i,j-1,k+1)+ &
424 & rho(i,j,k )-rho(i,j-1,k )
425 cff2=rho(i,j,k+1)+rho(i,j-1,k+1)- &
426 & rho(i,j,k )-rho(i,j-1,k )
427 cff3=z_r(i,j,k+1)+z_r(i,j-1,k+1)- &
428 & z_r(i,j,k )-z_r(i,j-1,k )
429 cff4=z_r(i,j,k+1)-z_r(i,j-1,k+1)+ &
430 & z_r(i,j,k )-z_r(i,j-1,k )
431
432!^ tl_phie(i)=tl_phie(i)+ &
433!^ & fac3*(tl_cff1*cff3+ &
434!^ & cff1*tl_cff3- &
435!^ & tl_cff2*cff4- &
436!^ & cff2*tl_cff4)
437!^
438 adfac=fac3*ad_phie(i)
439 ad_cff1=ad_cff1+cff3*adfac
440 ad_cff2=ad_cff2-cff4*adfac
441 ad_cff3=ad_cff3+cff1*adfac
442 ad_cff4=ad_cff4-cff2*adfac
443!^ tl_cff4=tl_z_r(i,j,k+1)-tl_z_r(i,j-1,k+1)+ &
444!^ & tl_z_r(i,j,k )-tl_z_r(i,j-1,k )
445!^ tl_cff3=tl_z_r(i,j,k+1)+tl_z_r(i,j-1,k+1)- &
446!^ & tl_z_r(i,j,k )-tl_z_r(i,j-1,k )
447!^
448 ad_z_r(i,j-1,k )=ad_z_r(i,j-1,k )-ad_cff3-ad_cff4
449 ad_z_r(i,j ,k )=ad_z_r(i,j ,k )-ad_cff3+ad_cff4
450 ad_z_r(i,j-1,k+1)=ad_z_r(i,j-1,k+1)+ad_cff3-ad_cff4
451 ad_z_r(i,j ,k+1)=ad_z_r(i,j ,k+1)+ad_cff3+ad_cff4
452 ad_cff4=0.0_r8
453 ad_cff3=0.0_r8
454!^ tl_cff2=tl_rho(i,j,k+1)+tl_rho(i,j-1,k+1)- &
455!^ & tl_rho(i,j,k )-tl_rho(i,j-1,k )
456!^ tl_cff1=tl_rho(i,j,k+1)-tl_rho(i,j-1,k+1)+ &
457!^ & tl_rho(i,j,k )-tl_rho(i,j-1,k )
458!^
459 ad_rho(i,j-1,k )=ad_rho(i,j-1,k )-ad_cff1-ad_cff2
460 ad_rho(i,j ,k )=ad_rho(i,j ,k )+ad_cff1-ad_cff2
461 ad_rho(i,j-1,k+1)=ad_rho(i,j-1,k+1)-ad_cff1+ad_cff2
462 ad_rho(i,j ,k+1)=ad_rho(i,j ,k+1)+ad_cff1+ad_cff2
463 ad_cff2=0.0_r8
464 ad_cff1=0.0_r8
465#endif
466 END DO
467 END DO
468!
469! Compute surface adjoint baroclinic pressure gradient.
470!
471 DO i=istr,iend
472#ifdef TIDE_GENERATING_FORCES
473 cff1=z_w(i,j ,n(ng))-eq_tide(i,j )-z_r(i,j ,n(ng))+ &
474 & z_w(i,j-1,n(ng))-eq_tide(i,j-1)-z_r(i,j-1,n(ng))
475#else
476 cff1=z_w(i,j ,n(ng))-z_r(i,j ,n(ng))+ &
477 & z_w(i,j-1,n(ng))-z_r(i,j-1,n(ng))
478#endif
479 phie(i)=fac1*(rho(i,j,n(ng))-rho(i,j-1,n(ng)))*cff1
480#ifdef ATM_PRESS
481 phie(i)=phie(i)+fac*(pair(i,j)-pair(i,j-1))
482#endif
483#ifdef RHO_SURF
484 phie(i)=phie(i)+ &
485 & (fac2+fac1*(rho(i,j,n(ng))+rho(i,j-1,n(ng))))* &
486 & (z_w(i,j,n(ng))-z_w(i,j-1,n(ng)))
487#endif
488# ifdef DIAGNOSTICS_UV
489!! DiaRV(i,j,N(ng),nrhs,M3pgrd)=rv(i,j,N(ng),nrhs)
490# endif
491!^ tl_rv(i,j,N(ng),nrhs)=-0.5_r8*om_v(i,j)* &
492!^ & ((tl_Hz(i,j ,N(ng))+ &
493!^ & tl_Hz(i,j-1,N(ng)))*phie(i)+ &
494!^ & (Hz(i,j ,N(ng))+ &
495!^ & Hz(i,j-1,N(ng)))*tl_phie(i))
496!^
497 adfac=-0.5_r8*om_v(i,j)*ad_rv(i,j,n(ng),nrhs)
498 adfac1=adfac*phie(i)
499 ad_phie(i)=ad_phie(i)+(hz(i,j ,n(ng))+ &
500 & hz(i,j-1,n(ng)))*adfac
501 ad_hz(i,j-1,n(ng))=ad_hz(i,j-1,n(ng))+adfac1
502 ad_hz(i,j ,n(ng))=ad_hz(i,j ,n(ng))+adfac1
503 ad_rv(i,j,n(ng),nrhs)=0.0
504#ifdef RHO_SURF
505!^ tl_phie(i)=tl_phie(i)+ &
506!^ & (fac1*(tl_rho(i,j,N(ng))+tl_rho(i,j-1,N(ng))))* &
507!^ & (z_w(i,j,N(ng))-z_w(i,j-1,N(ng)))+ &
508!^ & (fac2+fac1*(rho(i,j,N(ng))+rho(i,j-1,N(ng))))* &
509!^ & (tl_z_w(i,j,N(ng))-tl_z_w(i,j-1,N(ng)))
510!^
511 adfac1=fac1*(z_w(i,j,n(ng))-z_w(i,j-1,n(ng)))* &
512 & ad_phie(i)
513 adfac2=(fac2+fac1*(rho(i,j,n(ng))+rho(i,j-1,n(ng))))* &
514 & ad_phie(i)
515 ad_rho(i,j-1,n(ng))=ad_rho(i,j-1,n(ng))+adfac1
516 ad_rho(i,j ,n(ng))=ad_rho(i,j ,n(ng))+adfac1
517 ad_z_w(i,j-1,n(ng))=ad_z_w(i,j-1,n(ng))-adfac2
518 ad_z_w(i,j ,n(ng))=ad_z_w(i,j ,n(ng))+adfac2
519#endif
520!^ tl_phie(i)=fac1* &
521!^ & ((tl_rho(i,j,N(ng))-tl_rho(i,j-1,N(ng)))*cff1+ &
522!^ & (rho(i,j,N(ng))-rho(i,j-1,N(ng)))*tl_cff1)
523!^
524 adfac=fac1*ad_phie(i)
525 adfac1=adfac*cff1
526 ad_rho(i,j-1,n(ng))=ad_rho(i,j-1,n(ng))-adfac1
527 ad_rho(i,j ,n(ng))=ad_rho(i,j ,n(ng))+adfac1
528 ad_cff1=ad_cff1+ &
529 & (rho(i,j,n(ng))-rho(i,j-1,n(ng)))*adfac
530 ad_phie(i)=0.0_r8
531#ifdef TIDE_GENERATING_FORCES
532!^ tl_cff1=tl_z_w(i,j ,N(ng))-tl_eq_tide(i,j )- &
533!^ & tl_z_r(i,j ,N(ng))+ &
534!^ & tl_z_w(i,j-1,N(ng))-tl_eq_tide(i,j-1)- &
535!^ & tl_z_r(i,j-1,N(ng))
536!^
537 ad_eq_tide(i,j-1)=ad_eq_tide(i,j-1)-ad_cff1
538 ad_eq_tide(i,j )=ad_eq_tide(i,j )-ad_cff1
539 ad_z_r(i,j-1,n(ng))=ad_z_r(i,j-1,n(ng))-ad_cff1
540 ad_z_r(i,j ,n(ng))=ad_z_r(i,j ,n(ng))-ad_cff1
541 ad_z_w(i,j-1,n(ng))=ad_z_w(i,j-1,n(ng))+ad_cff1
542 ad_z_w(i,j ,n(ng))=ad_z_w(i,j ,n(ng))+ad_cff1
543 ad_cff1=0.0_r8
544#else
545!^ tl_cff1=tl_z_w(i,j ,N(ng))-tl_z_r(i,j ,N(ng))+ &
546!^ & tl_z_w(i,j-1,N(ng))-tl_z_r(i,j-1,N(ng))
547!^
548 ad_z_r(i,j-1,n(ng))=ad_z_r(i,j-1,n(ng))-ad_cff1
549 ad_z_r(i,j ,n(ng))=ad_z_r(i,j ,n(ng))-ad_cff1
550 ad_z_w(i,j-1,n(ng))=ad_z_w(i,j-1,n(ng))+ad_cff1
551 ad_z_w(i,j ,n(ng))=ad_z_w(i,j ,n(ng))+ad_cff1
552 ad_cff1=0.0_r8
553#endif
554 END DO
555 END IF
556!
557!-----------------------------------------------------------------------
558! Calculate adjoint pressure gradient in the XI-direction (m4/s2).
559!-----------------------------------------------------------------------
560!
561! Compute appropriate BASIC STATE "phix". Notice that a reverse
562! vertical integration of "phix" is carried out over kk-index.
563!
564!^ DO k=N(ng)-1,1,-1
565!^
566 DO k=1,n(ng)-1
567 DO i=istru,iend
568#ifdef TIDE_GENERATING_FORCES
569 cff1=z_w(i ,j,n(ng))-eq_tide(i ,j)-z_r(i ,j,n(ng))+ &
570 & z_w(i-1,j,n(ng))-eq_tide(i-1,j)-z_r(i-1,j,n(ng))
571#else
572 cff1=z_w(i ,j,n(ng))-z_r(i ,j,n(ng))+ &
573 & z_w(i-1,j,n(ng))-z_r(i-1,j,n(ng))
574#endif
575 phix(i)=fac1*(rho(i,j,n(ng))-rho(i-1,j,n(ng)))*cff1
576#ifdef ATM_PRESS
577 phix(i)=phix(i)+fac*(pair(i,j)-pair(i-1,j))
578#endif
579#ifdef RHO_SURF
580 phix(i)=phix(i)+ &
581 & (fac2+fac1*(rho(i,j,n(ng))+rho(i-1,j,n(ng))))* &
582 & (z_w(i,j,n(ng))-z_w(i-1,j,n(ng)))
583#endif
584 END DO
585 DO kk=n(ng)-1,k,-1
586 DO i=istru,iend
587#ifdef WJ_GRADP
588 cff1=1.0_r8/((z_r(i ,j,kk+1)-z_r(i ,j,kk))* &
589 & (z_r(i-1,j,kk+1)-z_r(i-1,j,kk)))
590 cff2=z_r(i ,j,kk )-z_r(i-1,j,kk )+ &
591 & z_r(i ,j,kk+1)-z_r(i-1,j,k+1)
592 cff3=z_r(i ,j,kk+1)-z_r(i ,j,kk )- &
593 & z_r(i-1,j,kk+1)+z_r(i-1,j,kk )
594 gamma=0.125_r8*cff1*cff2*cff3
595
596 cff1=(1.0_r8+gamma)*(rho(i,j,kk+1)-rho(i-1,j,kk+1))+ &
597 & (1.0_r8-gamma)*(rho(i,j,kk )-rho(i-1,j,kk ))
598 cff2=rho(i,j,kk+1)+rho(i-1,j,kk+1)- &
599 & rho(i,j,kk )-rho(i-1,j,kk )
600 cff3=z_r(i,j,kk+1)+z_r(i-1,j,kk+1)- &
601 & z_r(i,j,kk )-z_r(i-1,j,kk )
602 cff4=(1.0_r8+gamma)*(z_r(i,j,kk+1)-z_r(i-1,j,kk+1))+ &
603 & (1.0_r8-gamma)*(z_r(i,j,kk )-z_r(i-1,j,kk ))
604 phix(i)=phix(i)+ &
605 & fac3*(cff1*cff3-cff2*cff4)
606#else
607 cff1=rho(i,j,kk+1)-rho(i-1,j,kk+1)+ &
608 & rho(i,j,kk )-rho(i-1,j,kk )
609 cff2=rho(i,j,kk+1)+rho(i-1,j,kk+1)- &
610 & rho(i,j,kk )-rho(i-1,j,kk )
611 cff3=z_r(i,j,kk+1)+z_r(i-1,j,kk+1)- &
612 & z_r(i,j,kk )-z_r(i-1,j,kk )
613 cff4=z_r(i,j,kk+1)-z_r(i-1,j,kk+1)+ &
614 & z_r(i,j,kk )-z_r(i-1,j,kk )
615 phix(i)=phix(i)+ &
616 & fac3*(cff1*cff3-cff2*cff4)
617#endif
618 END DO
619 END DO
620!
621! Compute interior adjoint baroclinic pressure gradient. Differentiate
622! and then vertically integrate.
623!
624 DO i=istru,iend
625# ifdef DIAGNOSTICS_UV
626!! DiaRU(i,j,k,nrhs,M3pgrd)=ru(i,j,k,nrhs)
627# endif
628!^ tl_ru(i,j,k,nrhs)=-0.5_r8*on_u(i,j)* &
629!^ & ((tl_Hz(i,j,k)+tl_Hz(i-1,j,k))* &
630!^ & phix(i)+ &
631!^ & (Hz(i,j,k)+Hz(i-1,j,k))* &
632!^ & tl_phix(i))
633!^
634 adfac=-0.5_r8*on_u(i,j)*ad_ru(i,j,k,nrhs)
635 adfac1=adfac*phix(i)
636 ad_phix(i)=ad_phix(i)+ &
637 & (hz(i,j,k)+hz(i-1,j,k))*adfac
638 ad_hz(i-1,j,k)=ad_hz(i-1,j,k)+adfac1
639 ad_hz(i ,j,k)=ad_hz(i ,j,k)+adfac1
640 ad_ru(i,j,k,nrhs)=0.0
641#ifdef WJ_GRADP
642 cff1=1.0_r8/((z_r(i ,j,k+1)-z_r(i ,j,k))* &
643 & (z_r(i-1,j,k+1)-z_r(i-1,j,k)))
644 cff2=z_r(i ,j,k )-z_r(i-1,j,k )+ &
645 & z_r(i ,j,k+1)-z_r(i-1,j,k+1)
646 cff3=z_r(i ,j,k+1)-z_r(i ,j,k )- &
647 & z_r(i-1,j,k+1)+z_r(i-1,j,k )
648 gamma=0.125_r8*cff1*cff2*cff3
649
650 cff1=(1.0_r8+gamma)*(rho(i,j,k+1)-rho(i-1,j,k+1))+ &
651 & (1.0_r8-gamma)*(rho(i,j,k )-rho(i-1,j,k ))
652 cff2=rho(i,j,k+1)+rho(i-1,j,k+1)- &
653 & rho(i,j,k )-rho(i-1,j,k )
654 cff3=z_r(i,j,k+1)+z_r(i-1,j,k+1)- &
655 & z_r(i,j,k )-z_r(i-1,j,k )
656 cff4=(1.0_r8+gamma)*(z_r(i,j,k+1)-z_r(i-1,j,k+1))+ &
657 & (1.0_r8-gamma)*(z_r(i,j,k )-z_r(i-1,j,k ))
658!^ tl_phix(i)=tl_phix(i)+ &
659!^ & fac3*(tl_cff1*cff3+ &
660!^ & cff1*tl_cff3- &
661!^ & tl_cff2*cff4- &
662!^ & cff2*tl_cff4)
663!^
664 adfac=fac3*ad_phix(i)
665 ad_cff1=ad_cff1+cff3*adfac
666 ad_cff2=ad_cff2-cff4*adfac
667 ad_cff3=ad_cff3+cff1*adfac
668 ad_cff4=ad_cff4-cff2*adfac
669!^ tl_cff4=tl_gamma*(z_r(i,j,k+1)-z_r(i-1,j,k+1)- &
670!^ & z_r(i,j,k )+z_r(i-1,j,k ))+ &
671!^ & (1.0_r8+gamma)*(tl_z_r(i ,j,k+1)- &
672!^ & tl_z_r(i-1,j,k+1))+ &
673!^ & (1.0_r8-gamma)*(tl_z_r(i ,j,k )- &
674!^ & tl_z_r(i-1,j,k ))
675!^ tl_cff3=tl_z_r(i,j,k+1)+tl_z_r(i-1,j,k+1)- &
676!^ & tl_z_r(i,j,k )-tl_z_r(i-1,j,k )
677!^
678 adfac1=(1.0_r8+gamma)*ad_cff4
679 adfac2=(1.0_r8-gamma)*ad_cff4
680 ad_z_r(i-1,j,k )=ad_z_r(i-1,j,k )-adfac2-ad_cff3
681 ad_z_r(i ,j,k )=ad_z_r(i ,j,k )+adfac2-ad_cff3
682 ad_z_r(i-1,j,k+1)=ad_z_r(i-1,j,k+1)-adfac1+ad_cff3
683 ad_z_r(i ,j,k+1)=ad_z_r(i ,j,k+1)+adfac1+ad_cff3
684 ad_gamma=ad_gamma+ &
685 & (z_r(i,j,k+1)-z_r(i-1,j,k+1)- &
686 & z_r(i,j,k )+z_r(i-1,j,k ))*ad_cff4
687 ad_cff4=0.0_r8
688 ad_cff3=0.0_r8
689!^ & tl_cff2=tl_rho(i,j,k+1)+tl_rho(i-1,j,k+1)- &
690!^ & tl_rho(i,j,k )-tl_rho(i-1,j,k )
691!^ tl_cff1=tl_gamma*(rho(i,j,k+1)-rho(i-1,j,k+1)- &
692!^ & rho(i,j,k )+rho(i-1,j,k ))+ &
693!^ & (1.0_r8+gamma)*(tl_rho(i ,j,k+1)- &
694!^ & tl_rho(i-1,j,k+1))+ &
695!^ & (1.0_r8-gamma)*(tl_rho(i ,j,k )- &
696!^ & tl_rho(i-1,j,k ))
697!^
698 adfac1=(1.0_r8+gamma)*ad_cff1
699 adfac2=(1.0_r8-gamma)*ad_cff1
700 ad_rho(i-1,j,k )=ad_rho(i-1,j,k )-adfac2-ad_cff2
701 ad_rho(i ,j,k )=ad_rho(i ,j,k )+adfac2-ad_cff2
702 ad_rho(i-1,j,k+1)=ad_rho(i-1,j,k+1)-adfac1+ad_cff2
703 ad_rho(i ,j,k+1)=ad_rho(i ,j,k+1)+adfac1+ad_cff2
704 ad_gamma=ad_gamma+ &
705 & (rho(i,j,k+1)-rho(i-1,j,k+1)- &
706 & rho(i,j,k )+rho(i-1,j,k ))*ad_cff1
707 ad_cff2=0.0_r8
708 ad_cff1=0.0_r8
709!
710 cff1=1.0_r8/((z_r(i ,j,k+1)-z_r(i ,j,k))* &
711 & (z_r(i-1,j,k+1)-z_r(i-1,j,k)))
712 cff2=z_r(i ,j,k )-z_r(i-1,j,k )+ &
713 & z_r(i ,j,k+1)-z_r(i-1,j,k+1)
714 cff3=z_r(i ,j,k+1)-z_r(i ,j,k )- &
715 & z_r(i-1,j,k+1)+z_r(i-1,j,k )
716
717!^ tl_gamma=0.125_r8*(tl_cff1*cff2*cff3+ &
718!^ & cff1*(tl_cff2*cff3+ &
719!^ & cff2*tl_cff3))
720!^
721 adfac=0.125_r8*ad_gamma
722 adfac1=adfac*cff1
723 ad_cff3=ad_cff3+cff2*adfac1
724 ad_cff2=ad_cff2+cff3*adfac1
725 ad_cff1=ad_cff1+cff2*cff3*adfac
726 ad_gamma=0.0_r8
727!^ & tl_cff3=tl_z_r(i ,j,k+1)-tl_z_r(i ,j,k )- &
728!^ & tl_z_r(i-1,j,k+1)+tl_z_r(i-1,j,k )
729!^ tl_cff2=tl_z_r(i ,j,k )-tl_z_r(i-1,j,k )+ &
730!^ & tl_z_r(i ,j,k+1)-tl_z_r(i-1,j,k+1)
731!^
732 ad_z_r(i-1,j,k )=ad_z_r(i-1,j,k )-ad_cff2+ad_cff3
733 ad_z_r(i ,j,k )=ad_z_r(i ,j,k )+ad_cff2-ad_cff3
734 ad_z_r(i-1,j,k+1)=ad_z_r(i-1,j,k+1)-ad_cff2-ad_cff3
735 ad_z_r(i ,j,k+1)=ad_z_r(i ,j,k+1)+ad_cff2+ad_cff3
736 ad_cff3=0.0
737 ad_cff2=0.0
738!^ tl_cff1=-cff1*cff1*((tl_z_r(i ,j,k+1)-tl_z_r(i ,j,k))* &
739!^ & (z_r(i-1,j,k+1)-z_r(i-1,j,k))+ &
740!^ & (z_r(i ,j,k+1)-z_r(i ,j,k))* &
741!^ & (tl_z_r(i-1,j,k+1)-tl_z_r(i-1,j,k)))
742!^
743 adfac=-cff1*cff1*ad_cff1
744 adfac1=adfac*(z_r(i-1,j,k+1)-z_r(i-1,j,k))
745 adfac2=adfac*(z_r(i ,j,k+1)-z_r(i ,j,k))
746 ad_z_r(i-1,j,k )=ad_z_r(i-1,j,k )-adfac2
747 ad_z_r(i ,j,k )=ad_z_r(i ,j,k )-adfac1
748 ad_z_r(i-1,j,k+1)=ad_z_r(i-1,j,k+1)+adfac2
749 ad_z_r(i ,j,k+1)=ad_z_r(i ,j,k+1)+adfac1
750 ad_cff1=0.0_r8
751#else
752 cff1=rho(i,j,k+1)-rho(i-1,j,k+1)+ &
753 & rho(i,j,k )-rho(i-1,j,k )
754 cff2=rho(i,j,k+1)+rho(i-1,j,k+1)- &
755 & rho(i,j,k )-rho(i-1,j,k )
756 cff3=z_r(i,j,k+1)+z_r(i-1,j,k+1)- &
757 & z_r(i,j,k )-z_r(i-1,j,k )
758 cff4=z_r(i,j,k+1)-z_r(i-1,j,k+1)+ &
759 & z_r(i,j,k )-z_r(i-1,j,k )
760!^ tl_phix(i)=tl_phix(i)+ &
761!^ & fac3*(tl_cff1*cff3+ &
762!^ & cff1*tl_cff3- &
763!^ & tl_cff2*cff4- &
764!^ & cff2*tl_cff4)
765!^
766 adfac=fac3*ad_phix(i)
767 ad_cff1=ad_cff1+cff3*adfac
768 ad_cff2=ad_cff2-cff4*adfac
769 ad_cff3=ad_cff3+cff1*adfac
770 ad_cff4=ad_cff4-cff2*adfac
771!^ tl_cff4=tl_z_r(i,j,k+1)-tl_z_r(i-1,j,k+1)+ &
772!^ & tl_z_r(i,j,k )-tl_z_r(i-1,j,k )
773!^ tl_cff3=tl_z_r(i,j,k+1)+tl_z_r(i-1,j,k+1)- &
774!^ & tl_z_r(i,j,k )-tl_z_r(i-1,j,k )
775!^
776 ad_z_r(i-1,j,k )=ad_z_r(i-1,j,k )-ad_cff3-ad_cff4
777 ad_z_r(i ,j,k )=ad_z_r(i ,j,k )-ad_cff3+ad_cff4
778 ad_z_r(i-1,j,k+1)=ad_z_r(i-1,j,k+1)+ad_cff3-ad_cff4
779 ad_z_r(i ,j,k+1)=ad_z_r(i ,j,k+1)+ad_cff3+ad_cff4
780 ad_cff4=0.0_r8
781 ad_cff3=0.0_r8
782!^ tl_cff1=tl_rho(i,j,k+1)-tl_rho(i-1,j,k+1)+ &
783!^ & tl_rho(i,j,k )-tl_rho(i-1,j,k )
784!^ tl_cff2=tl_rho(i,j,k+1)+tl_rho(i-1,j,k+1)- &
785!^ & tl_rho(i,j,k )-tl_rho(i-1,j,k )
786!^
787 ad_rho(i-1,j,k )=ad_rho(i-1,j,k )-ad_cff2-ad_cff1
788 ad_rho(i ,j,k )=ad_rho(i ,j,k )-ad_cff2+ad_cff1
789 ad_rho(i-1,j,k+1)=ad_rho(i-1,j,k+1)+ad_cff2-ad_cff1
790 ad_rho(i ,j,k+1)=ad_rho(i ,j,k+1)+ad_cff2+ad_cff1
791 ad_cff2=0.0_r8
792 ad_cff1=0.0_r8
793#endif
794 END DO
795 END DO
796!
797! Compute surface adjoint baroclinic pressure gradient.
798!
799 DO i=istru,iend
800#ifdef TIDE_GENERATING_FORCES
801 cff1=z_w(i ,j,n(ng))-eq_tide(i ,j)-z_r(i ,j,n(ng))+ &
802 & z_w(i-1,j,n(ng))-eq_tide(i-1,j)-z_r(i-1,j,n(ng))
803#else
804 cff1=z_w(i ,j,n(ng))-z_r(i ,j,n(ng))+ &
805 & z_w(i-1,j,n(ng))-z_r(i-1,j,n(ng))
806#endif
807 phix(i)=fac1*(rho(i,j,n(ng))-rho(i-1,j,n(ng)))*cff1
808#ifdef ATM_PRESS
809 phix(i)=phix(i)+fac*(pair(i,j)-pair(i-1,j))
810#endif
811#ifdef RHO_SURF
812 phix(i)=phix(i)+ &
813 & (fac2+fac1*(rho(i,j,n(ng))+rho(i-1,j,n(ng))))* &
814 & (z_w(i,j,n(ng))-z_w(i-1,j,n(ng)))
815#endif
816#ifdef DIAGNOSTICS_UV
817!! DiaRU(i,j,N(ng),nrhs,M3pgrd)=ru(i,j,N(ng),nrhs)
818#endif
819!^ tl_ru(i,j,N(ng),nrhs)=-0.5_r8*on_u(i,j)* &
820!^ & ((tl_Hz(i ,j,N(ng))+ &
821!^ & tl_Hz(i-1,j,N(ng)))*phix(i)+ &
822!^ & (Hz(i ,j,N(ng))+ &
823!^ & Hz(i-1,j,N(ng)))*tl_phix(i))
824!^
825 adfac=-0.5_r8*on_u(i,j)*ad_ru(i,j,n(ng),nrhs)
826 adfac1=adfac*phix(i)
827 ad_phix(i)=ad_phix(i)+(hz(i ,j,n(ng))+ &
828 & hz(i-1,j,n(ng)))*adfac
829 ad_hz(i-1,j,n(ng))=ad_hz(i-1,j,n(ng))+adfac1
830 ad_hz(i ,j,n(ng))=ad_hz(i ,j,n(ng))+adfac1
831 ad_ru(i,j,n(ng),nrhs)=0.0_r8
832#ifdef RHO_SURF
833!^ tl_phix(i)=tl_phix(i)+ &
834!^ & (fac1*(tl_rho(i,j,N(ng))+tl_rho(i-1,j,N(ng))))* &
835!^ & (z_w(i,j,N(ng))-z_w(i-1,j,N(ng)))+ &
836!^ & (fac2+fac1*(rho(i,j,N(ng))+rho(i-1,j,N(ng))))*
837!^ & (tl_z_w(i,j,N(ng))-tl_z_w(i-1,j,N(ng)))
838!^
839 adfac1=fac1*(z_w(i,j,n(ng))-z_w(i-1,j,n(ng)))* &
840 & ad_phix(i)
841 adfac2=(fac2+fac1*(rho(i,j,n(ng))+rho(i-1,j,n(ng))))* &
842 & ad_phix(i)
843 ad_rho(i-1,j,n(ng))=ad_rho(i-1,j,n(ng))+adfac1
844 ad_rho(i ,j,n(ng))=ad_rho(i ,j,n(ng))+adfac1
845 ad_z_w(i-1,j,n(ng))=ad_z_w(i-1,j,n(ng))-adfac2
846 ad_z_w(i ,j,n(ng))=ad_z_w(i ,j,n(ng))+adfac2
847#endif
848!^ tl_phix(i)=fac1* &
849!^ & ((tl_rho(i,j,N(ng))-tl_rho(i-1,j,N(ng)))*cff1+ &
850!^ & (rho(i,j,N(ng))-rho(i-1,j,N(ng)))*tl_cff1)
851!^
852 adfac=fac1*ad_phix(i)
853 adfac1=adfac*cff1
854 ad_rho(i-1,j,n(ng))=ad_rho(i-1,j,n(ng))-adfac1
855 ad_rho(i ,j,n(ng))=ad_rho(i ,j,n(ng))+adfac1
856 ad_cff1=ad_cff1+ &
857 & (rho(i,j,n(ng))-rho(i-1,j,n(ng)))*adfac
858 ad_phix(i)=0.0
859#ifdef TIDE_GENERATING_FORCES
860!^ tl_cff1=tl_z_w(i ,j,N(ng))-tl_eq_tide(i ,j)- &
861!^ & tl_z_r(i ,j,N(ng))+ &
862!^ & tl_z_w(i-1,j,N(ng))-tl_eq_tide(i-1,j)- &
863!^ & tl_z_r(i-1,j,N(ng))
864!^
865 ad_eq_tide(i-1,j)=ad_eq_tide(i-1,j)-ad_cff1
866 ad_eq_tide(i ,j)=ad_eq_tide(i ,j)-ad_cff1
867 ad_z_r(i-1,j,n(ng))=ad_z_r(i-1,j,n(ng))-ad_cff1
868 ad_z_r(i ,j,n(ng))=ad_z_r(i ,j,n(ng))-ad_cff1
869 ad_z_w(i-1,j,n(ng))=ad_z_w(i-1,j,n(ng))+ad_cff1
870 ad_z_w(i ,j,n(ng))=ad_z_w(i ,j,n(ng))+ad_cff1
871 ad_cff1=0.0_r8
872#else
873!^ tl_cff1=tl_z_w(i ,j,N(ng))-tl_z_r(i ,j,N(ng))+ &
874!^ & tl_z_w(i-1,j,N(ng))-tl_z_r(i-1,j,N(ng))
875!^
876 ad_z_r(i-1,j,n(ng))=ad_z_r(i-1,j,n(ng))-ad_cff1
877 ad_z_r(i ,j,n(ng))=ad_z_r(i ,j,n(ng))-ad_cff1
878 ad_z_w(i-1,j,n(ng))=ad_z_w(i-1,j,n(ng))+ad_cff1
879 ad_z_w(i ,j,n(ng))=ad_z_w(i ,j,n(ng))+ad_cff1
880 ad_cff1=0.0_r8
881#endif
882 END DO
883 END DO j_loop
884!
885 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 ad_prsgrd().

Here is the caller graph for this function:

◆ ad_prsgrd32_tile()

subroutine ad_prsgrd_mod::ad_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(inout) ad_hz,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) z_r,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(inout) ad_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(inout) ad_z_w,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) rho,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(inout) ad_rho,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) eq_tide,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(inout) ad_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) ad_ru,
real(r8), dimension(lbi:ubi,lbj:ubj,0:n(ng),2), intent(inout) ad_rv )
private

Definition at line 107 of file ad_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# ifdef ATM_PRESS
153 real(r8), intent(in) :: Pair(LBi:,LBj:)
154# endif
155# ifdef TIDE_GENERATING_FORCES
156 real(r8), intent(in) :: eq_tide(LBi:,LBj:)
157 real(r8), intent(inout) :: ad_eq_tide(LBi:,LBj:)
158# endif
159# ifdef DIAGNOSTICS_UV
160!! real(r8), intent(inout) :: DiaRU(LBi:,LBj:,:,:,:)
161!! real(r8), intent(inout) :: DiaRV(LBi:,LBj:,:,:,:)
162# endif
163 real(r8), intent(inout) :: ad_Hz(LBi:,LBj:,:)
164 real(r8), intent(inout) :: ad_z_r(LBi:,LBj:,:)
165 real(r8), intent(inout) :: ad_z_w(LBi:,LBj:,0:)
166 real(r8), intent(inout) :: ad_rho(LBi:,LBj:,:)
167 real(r8), intent(inout) :: ad_ru(LBi:,LBj:,0:,:)
168 real(r8), intent(inout) :: ad_rv(LBi:,LBj:,0:,:)
169#else
170# ifdef MASKING
171 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
172 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
173# endif
174 real(r8), intent(in) :: om_v(LBi:UBi,LBj:UBj)
175 real(r8), intent(in) :: on_u(LBi:UBi,LBj:UBj)
176 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
177 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
178 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
179 real(r8), intent(in) :: rho(LBi:UBi,LBj:UBj,N(ng))
180# ifdef ATM_PRESS
181 real(r8), intent(in) :: Pair(LBi:UBi,LBj:UBj)
182# endif
183# ifdef TIDE_GENERATING_FORCES
184 real(r8), intent(in) :: eq_tide(LBi:UBi,LBj:UBj)
185 real(r8), intent(inout) :: ad_eq_tide(LBi:UBi,LBj:UBj)
186# endif
187# ifdef DIAGNOSTICS_UV
188!! real(r8), intent(inout) :: DiaRU(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
189!! real(r8), intent(inout) :: DiaRV(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
190# endif
191 real(r8), intent(inout) :: ad_Hz(LBi:UBi,LBj:UBj,N(ng))
192 real(r8), intent(inout) :: ad_z_r(LBi:UBi,LBj:UBj,N(ng))
193 real(r8), intent(inout) :: ad_z_w(LBi:UBi,LBj:UBj,0:N(ng))
194 real(r8), intent(inout) :: ad_rho(LBi:UBi,LBj:UBj,N(ng))
195 real(r8), intent(inout) :: ad_ru(LBi:UBi,LBj:UBj,0:N(ng),2)
196 real(r8), intent(inout) :: ad_rv(LBi:UBi,LBj:UBj,0:N(ng),2)
197#endif
198!
199! Local variable declarations.
200!
201 integer :: i, j, k
202
203 real(r8), parameter :: OneFifth = 0.2_r8
204 real(r8), parameter :: OneTwelfth = 1.0_r8/12.0_r8
205 real(r8), parameter :: eps = 1.0e-10_r8
206
207 real(r8) :: GRho, GRho0, HalfGRho
208 real(r8) :: cff, cff1, cff2
209#ifdef ATM_PRESS
210 real(r8) :: OneAtm, fac
211#endif
212 real(r8) :: ad_cff, ad_cff1, ad_cff2, adfac
213 real(r8) :: adfac1, adfac2, adfac3, adfac4, adfac5, adfac6
214 real(r8) :: adfac7, adfac8, adfac9, adfac10, adfac11, adfac12
215
216 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,N(ng)) :: P
217
218 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,N(ng)) :: ad_P
219
220 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: dR
221 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: dR1
222 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: dZ
223 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: dZ1
224
225 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: ad_dR
226 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: ad_dZ
227
228 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FC
229 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: aux
230 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: dRx
231 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: dZx
232
233 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_FC
234 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_aux
235 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_dRx
236 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_dZx
237
238#include "set_bounds.h"
239!
240!-----------------------------------------------------------------------
241! Initialize adjoint private variables.
242!-----------------------------------------------------------------------
243!
244 ad_cff=0.0_r8
245 ad_cff1=0.0_r8
246 ad_cff2=0.0_r8
247 DO j=jmins,jmaxs
248 DO i=imins,imaxs
249 ad_fc(i,j)=0.0_r8
250 ad_aux(i,j)=0.0_r8
251 ad_drx(i,j)=0.0_r8
252 ad_dzx(i,j)=0.0_r8
253 END DO
254 DO k=1,n(ng)
255 DO i=imins,imaxs
256 ad_p(i,j,k)=0.0_r8
257 END DO
258 END DO
259 END DO
260 DO k=0,n(ng)
261 DO i=imins,imaxs
262 ad_dr(i,k)=0.0_r8
263 ad_dz(i,k)=0.0_r8
264 END DO
265 END DO
266!
267!-----------------------------------------------------------------------
268! Preliminary step (same for XI- and ETA-components):
269!-----------------------------------------------------------------------
270!
271! Compute BASIC STATE dynamic pressure, P.
272!
273 grho=g/rho0
274 grho0=1000.0_r8*grho
275 halfgrho=0.5_r8*grho
276#ifdef ATM_PRESS
277 oneatm=1013.25_r8 ! 1 atm = 1013.25 mb
278 fac=100.0_r8/rho0
279#endif
280!
281 DO j=jstrv-1,jend
282 DO k=1,n(ng)-1
283 DO i=istru-1,iend
284 dr(i,k)=rho(i,j,k+1)-rho(i,j,k)
285 dz(i,k)=z_r(i,j,k+1)-z_r(i,j,k)
286 END DO
287 END DO
288 DO i=istru-1,iend
289 dr(i,n(ng))=dr(i,n(ng)-1)
290 dz(i,n(ng))=dz(i,n(ng)-1)
291 dr(i,0)=dr(i,1)
292 dz(i,0)=dz(i,1)
293 END DO
294 DO k=n(ng),1,-1
295 DO i=istru-1,iend
296 cff=2.0_r8*dr(i,k)*dr(i,k-1)
297 IF (cff.gt.eps) THEN
298 dr(i,k)=cff/(dr(i,k)+dr(i,k-1))
299 ELSE
300 dr(i,k)=0.0_r8
301 END IF
302 dz(i,k)=2.0_r8*dz(i,k)*dz(i,k-1)/(dz(i,k)+dz(i,k-1))
303 END DO
304 END DO
305 DO i=istru-1,iend
306 cff1=1.0_r8/(z_r(i,j,n(ng))-z_r(i,j,n(ng)-1))
307 cff2=0.5_r8*(rho(i,j,n(ng))-rho(i,j,n(ng)-1))* &
308 & (z_w(i,j,n(ng))-z_r(i,j,n(ng)))*cff1
309 p(i,j,n(ng))=g*z_w(i,j,n(ng))+ &
310#ifdef ATM_PRESS
311 & fac*(pair(i,j)-oneatm)+ &
312#endif
313 & grho*(rho(i,j,n(ng))+cff2)* &
314 & (z_w(i,j,n(ng))-z_r(i,j,n(ng)))
315#ifdef TIDE_GENERATING_FORCES
316 p(i,j,n(ng))=p(i,j,n(ng))-g*eq_tide(i,j)
317#endif
318 END DO
319 DO k=n(ng)-1,1,-1
320 DO i=istru-1,iend
321 p(i,j,k)=p(i,j,k+1)+ &
322 & halfgrho*((rho(i,j,k+1)+rho(i,j,k))* &
323 & (z_r(i,j,k+1)-z_r(i,j,k))- &
324 & onefifth* &
325 & ((dr(i,k+1)-dr(i,k))* &
326 & (z_r(i,j,k+1)-z_r(i,j,k)- &
327 & onetwelfth* &
328 & (dz(i,k+1)+dz(i,k)))- &
329 & (dz(i,k+1)-dz(i,k))* &
330 & (rho(i,j,k+1)-rho(i,j,k)- &
331 & onetwelfth* &
332 & (dr(i,k+1)+dr(i,k)))))
333 END DO
334 END DO
335 END DO
336!
337!-----------------------------------------------------------------------
338! Adjoint ETA-component pressure gradient term.
339!-----------------------------------------------------------------------
340!
341 DO k=1,n(ng)
342!
343! Compute BASIC STATE variables.
344!
345 DO j=jstrv-1,jend+1
346 DO i=istr,iend
347 aux(i,j)=z_r(i,j,k)-z_r(i,j-1,k)
348# ifdef MASKING
349 aux(i,j)=aux(i,j)*vmask(i,j)
350# endif
351 fc(i,j)=rho(i,j,k)-rho(i,j-1,k)
352# ifdef MASKING
353 fc(i,j)=fc(i,j)*vmask(i,j)
354# endif
355 END DO
356 END DO
357!
358 DO j=jstrv-1,jend
359 DO i=istr,iend
360 cff=2.0_r8*aux(i,j)*aux(i,j+1)
361 IF (cff.gt.eps) THEN
362 cff1=1.0_r8/(aux(i,j)+aux(i,j+1))
363 dzx(i,j)=cff*cff1
364 ELSE
365 dzx(i,j)=0.0_r8
366 END IF
367 cff1=2.0_r8*fc(i,j)*fc(i,j+1)
368 IF (cff1.gt.eps) THEN
369 cff2=1.0_r8/(fc(i,j)+fc(i,j+1))
370 drx(i,j)=cff1*cff2
371 ELSE
372 drx(i,j)=0.0_r8
373 END IF
374 END DO
375 END DO
376!
377 DO j=jstrv,jend
378 DO i=istr,iend
379#ifdef DIAGNOSTICS_UV
380!! DiaRV(i,j,k,nrhs,M3pgrd)=rv(i,j,k,nrhs)
381#endif
382!^ tl_rv(i,j,k,nrhs)=om_v(i,j)*0.5_r8* &
383!^ & ((tl_Hz(i,j,k)+tl_Hz(i,j-1,k))* &
384!^ & (P(i,j-1,k)-P(i,j,k)- &
385!^ & HalfGRho* &
386!^ & ((rho(i,j,k)+rho(i,j-1,k))* &
387!^ & (z_r(i,j,k)-z_r(i,j-1,k))- &
388!^ & OneFifth* &
389!^ & ((dRx(i,j)-dRx(i,j-1))* &
390!^ & (z_r(i,j,k)-z_r(i,j-1,k)- &
391!^ & OneTwelfth* &
392!^ & (dZx(i,j)+dZx(i,j-1)))- &
393!^ & (dZx(i,j)-dZx(i,j-1))* &
394!^ & (rho(i,j,k)-rho(i,j-1,k)- &
395!^ & OneTwelfth* &
396!^ & (dRx(i,j)+dRx(i,j-1))))))+ &
397!^ & (Hz(i,j,k)+Hz(i,j-1,k))* &
398!^ & (tl_P(i,j-1,k)-tl_P(i,j,k)- &
399!^ & HalfGRho* &
400!^ & ((tl_rho(i,j,k)+tl_rho(i,j-1,k))* &
401!^ & (z_r(i,j,k)-z_r(i,j-1,k))+ &
402!^ & (rho(i,j,k)+rho(i,j-1,k))* &
403!^ & (tl_z_r(i,j,k)-tl_z_r(i,j-1,k))- &
404!^ & OneFifth* &
405!^ & ((tl_dRx(i,j)-tl_dRx(i,j-1))* &
406!^ & (z_r(i,j,k)-z_r(i,j-1,k)- &
407!^ & OneTwelfth* &
408!^ & (dZx(i,j)+dZx(i,j-1)))+ &
409!^ & (dRx(i,j)-dRx(i,j-1))* &
410!^ & (tl_z_r(i,j,k)-tl_z_r(i,j-1,k)- &
411!^ & OneTwelfth* &
412!^ & (tl_dZx(i,j)+tl_dZx(i,j-1)))- &
413!^ & (tl_dZx(i,j)-tl_dZx(i,j-1))* &
414!^ & (rho(i,j,k)-rho(i,j-1,k)- &
415!^ & OneTwelfth* &
416!^ & (dRx(i,j)+dRx(i,j-1)))- &
417!^ & (dZx(i,j)-dZx(i,j-1))* &
418!^ & (tl_rho(i,j,k)-tl_rho(i,j-1,k)- &
419!^ & OneTwelfth* &
420!^ & (tl_dRx(i,j)+tl_dRx(i,j-1)))))))
421!^
422 adfac=om_v(i,j)*0.5_r8*ad_rv(i,j,k,nrhs)
423 adfac1=adfac*(p(i,j-1,k)-p(i,j,k)- &
424 & halfgrho* &
425 & ((rho(i,j,k)+rho(i,j-1,k))* &
426 & (z_r(i,j,k)-z_r(i,j-1,k))- &
427 & onefifth* &
428 & ((drx(i,j)-drx(i,j-1))* &
429 & (z_r(i,j,k)-z_r(i,j-1,k)- &
430 & onetwelfth* &
431 & (dzx(i,j)+dzx(i,j-1)))- &
432 & (dzx(i,j)-dzx(i,j-1))* &
433 & (rho(i,j,k)-rho(i,j-1,k)- &
434 & onetwelfth* &
435 & (drx(i,j)+drx(i,j-1))))))
436 adfac2=adfac*(hz(i,j,k)+hz(i,j-1,k))
437 adfac3=adfac2*halfgrho
438 adfac4=adfac3*(z_r(i,j,k)-z_r(i,j-1,k))
439 adfac5=adfac3*(rho(i,j,k)+rho(i,j-1,k))
440 adfac6=adfac3*onefifth
441 adfac7=adfac6*(z_r(i,j,k)-z_r(i,j-1,k)- &
442 & onetwelfth*(dzx(i,j)+dzx(i,j-1)))
443 adfac8=adfac6*(drx(i,j)-drx(i,j-1))
444 adfac9=adfac8*onetwelfth
445 adfac10=adfac6*(rho(i,j,k)-rho(i,j-1,k)- &
446 & onetwelfth*(drx(i,j)+drx(i,j-1)))
447 adfac11=adfac6*(dzx(i,j)-dzx(i,j-1))
448 adfac12=adfac11*onetwelfth
449 ad_hz(i,j-1,k)=ad_hz(i,j-1,k)+adfac1
450 ad_hz(i,j ,k)=ad_hz(i,j ,k)+adfac1
451 ad_p(i,j-1,k)=ad_p(i,j-1,k)+adfac2
452 ad_p(i,j ,k)=ad_p(i,j ,k)-adfac2
453 ad_rho(i,j-1,k)=ad_rho(i,j-1,k)-adfac4+adfac11
454 ad_rho(i,j ,k)=ad_rho(i,j ,k)-adfac4-adfac11
455 ad_z_r(i,j-1,k)=ad_z_r(i,j-1,k)+adfac5-adfac8
456 ad_z_r(i,j ,k)=ad_z_r(i,j ,k)-adfac5+adfac8
457 ad_drx(i,j-1)=ad_drx(i,j-1)-adfac7+adfac12
458 ad_drx(i,j )=ad_drx(i,j )+adfac7+adfac12
459 ad_dzx(i,j-1)=ad_dzx(i,j-1)-adfac9+adfac10
460 ad_dzx(i,j )=ad_dzx(i,j )-adfac9-adfac10
461 ad_rv(i,j,k,nrhs)=0.0_r8
462 END DO
463 END DO
464!
465 DO j=jstrv-1,jend
466 DO i=istr,iend
467 cff1=2.0_r8*fc(i,j)*fc(i,j+1)
468 IF (cff1.gt.eps) THEN
469 cff2=1.0_r8/(fc(i,j)+fc(i,j+1))
470!^ tl_dRx(i,j)=tl_cff1*cff2+cff1*tl_cff2
471!^
472 ad_cff1=ad_cff1+cff2*ad_drx(i,j)
473 ad_cff2=ad_cff2+cff1*ad_drx(i,j)
474 ad_drx(i,j)=0.0_r8
475!^ tl_cff2=-cff2*cff2*(tl_FC(i,j)+tl_FC(i,j+1))
476!^
477 adfac=-cff2*cff2*ad_cff2
478 ad_fc(i,j )=ad_fc(i,j )+adfac
479 ad_fc(i,j+1)=ad_fc(i,j+1)+adfac
480 ad_cff2=0.0_r8
481 ELSE
482!^ tl_dRx(i,j)=0.0_r8
483!^
484 ad_drx(i,j)=0.0_r8
485 END IF
486!^ tl_cff1=2.0_r8*(tl_FC(i,j)*FC(i,j+1)+ &
487!^ & FC(i,j)*tl_FC(i,j+1))
488!^
489 adfac=2.0_r8*ad_cff1
490 ad_fc(i,j )=ad_fc(i,j )+fc(i,j+1)*adfac
491 ad_fc(i,j+1)=ad_fc(i,j+1)+fc(i,j )*adfac
492 ad_cff1=0.0_r8
493
494 cff=2.0_r8*aux(i,j)*aux(i,j+1)
495 IF (cff.gt.eps) THEN
496 cff1=1.0_r8/(aux(i,j)+aux(i,j+1))
497!^ tl_dZx(i,j)=tl_cff*cff1+cff*tl_cff1
498!^
499 ad_cff=ad_cff+cff1*ad_dzx(i,j)
500 ad_cff1=ad_cff1+cff*ad_dzx(i,j)
501 ad_dzx(i,j)=0.0_r8
502!^ tl_cff1=-cff1*cff1*(tl_aux(i,j)+tl_aux(i,j+1))
503!^
504 adfac=-cff1*cff1*ad_cff1
505 ad_aux(i,j )=ad_aux(i,j )+adfac
506 ad_aux(i,j+1)=ad_aux(i,j+1)+adfac
507 ad_cff1=0.0_r8
508 ELSE
509!^ tl_dZx(i,j)=0.0_r8
510!^
511 ad_dzx(i,j)=0.0_r8
512 END IF
513!^ tl_cff=2.0_r8*(tl_aux(i,j)*aux(i,j+1)+ &
514!^ & aux(i,j)*tl_aux(i,j+1))
515!^
516 adfac=2.0_r8*ad_cff
517 ad_aux(i,j )=ad_aux(i,j )+aux(i,j+1)*adfac
518 ad_aux(i,j+1)=ad_aux(i,j+1)+aux(i,j )*adfac
519 ad_cff=0.0_r8
520 END DO
521 END DO
522 DO j=jstrv-1,jend+1
523 DO i=istr,iend
524#ifdef MASKING
525!^ tl_FC(i,j)=tl_FC(i,j)*vmask(i,j)
526!^
527 ad_fc(i,j)=ad_fc(i,j)*vmask(i,j)
528#endif
529!^ tl_FC(i,j)=tl_rho(i,j,k)-tl_rho(i,j-1,k)
530!^
531 ad_rho(i,j-1,k)=ad_rho(i,j-1,k)-ad_fc(i,j)
532 ad_rho(i,j, k)=ad_rho(i,j ,k)+ad_fc(i,j)
533 ad_fc(i,j)=0.0_r8
534#ifdef MASKING
535!^ tl_aux(i,j)=tl_aux(i,j)*vmask(i,j)
536!^
537 ad_aux(i,j)=ad_aux(i,j)*vmask(i,j)
538#endif
539!^ tl_aux(i,j)=tl_z_r(i,j,k)-tl_z_r(i,j-1,k)
540!^
541 ad_z_r(i,j-1,k)=ad_z_r(i,j-1,k)-ad_aux(i,j)
542 ad_z_r(i,j ,k)=ad_z_r(i,j ,k)+ad_aux(i,j)
543 ad_aux(i,j)=0.0_r8
544 END DO
545 END DO
546 END DO
547!
548!-----------------------------------------------------------------------
549! Compute adjoint XI-component pressure gradient term.
550!-----------------------------------------------------------------------
551!
552 DO k=1,n(ng)
553!
554! Compute BASIC STATE variables.
555!
556 DO j=jstr,jend
557 DO i=istru-1,iend+1
558 aux(i,j)=z_r(i,j,k)-z_r(i-1,j,k)
559# ifdef MASKING
560 aux(i,j)=aux(i,j)*umask(i,j)
561# endif
562 fc(i,j)=rho(i,j,k)-rho(i-1,j,k)
563# ifdef MASKING
564 fc(i,j)=fc(i,j)*umask(i,j)
565# endif
566 END DO
567 END DO
568!
569 DO j=jstr,jend
570 DO i=istru-1,iend
571 cff=2.0_r8*aux(i,j)*aux(i+1,j)
572 IF (cff.gt.eps) THEN
573 cff1=1.0_r8/(aux(i,j)+aux(i+1,j))
574 dzx(i,j)=cff*cff1
575 ELSE
576 dzx(i,j)=0.0_r8
577 END IF
578 cff1=2.0_r8*fc(i,j)*fc(i+1,j)
579 IF (cff1.gt.eps) THEN
580 cff2=1.0_r8/(fc(i,j)+fc(i+1,j))
581 drx(i,j)=cff1*cff2
582 ELSE
583 drx(i,j)=0.0_r8
584 END IF
585 END DO
586 END DO
587!
588 DO j=jstr,jend
589 DO i=istru,iend
590#ifdef DIAGNOSTICS_UV
591!! DiaRU(i,j,k,nrhs,M3pgrd)=ru(i,j,k,nrhs)
592#endif
593!^ tl_ru(i,j,k,nrhs)=on_u(i,j)*0.5_r8* &
594!^ & ((tl_Hz(i,j,k)+tl_Hz(i-1,j,k))* &
595!^ & (P(i-1,j,k)-P(i,j,k)- &
596!^ & HalfGRho* &
597!^ & ((rho(i,j,k)+rho(i-1,j,k))* &
598!^ & (z_r(i,j,k)-z_r(i-1,j,k))- &
599!^ & OneFifth* &
600!^ & ((dRx(i,j)-dRx(i-1,j))* &
601!^ & (z_r(i,j,k)-z_r(i-1,j,k)- &
602!^ & OneTwelfth* &
603!^ & (dZx(i,j)+dZx(i-1,j)))- &
604!^ & (dZx(i,j)-dZx(i-1,j))* &
605!^ & (rho(i,j,k)-rho(i-1,j,k)- &
606!^ & OneTwelfth* &
607!^ & (dRx(i,j)+dRx(i-1,j))))))+ &
608!^ & (Hz(i,j,k)+Hz(i-1,j,k))* &
609!^ & (tl_P(i-1,j,k)-tl_P(i,j,k)- &
610!^ & HalfGRho* &
611!^ & ((tl_rho(i,j,k)+tl_rho(i-1,j,k))* &
612!^ & (z_r(i,j,k)-z_r(i-1,j,k))+ &
613!^ & (rho(i,j,k)+rho(i-1,j,k))* &
614!^ & (tl_z_r(i,j,k)-tl_z_r(i-1,j,k))- &
615!^ & OneFifth* &
616!^ & ((tl_dRx(i,j)-tl_dRx(i-1,j))* &
617!^ & (z_r(i,j,k)-z_r(i-1,j,k)- &
618!^ & OneTwelfth* &
619!^ & (dZx(i,j)+dZx(i-1,j)))+ &
620!^ & (dRx(i,j)-dRx(i-1,j))* &
621!^ & (tl_z_r(i,j,k)-tl_z_r(i-1,j,k)- &
622!^ & OneTwelfth* &
623!^ & (tl_dZx(i,j)+tl_dZx(i-1,j)))- &
624!^ & (tl_dZx(i,j)-tl_dZx(i-1,j))* &
625!^ & (rho(i,j,k)-rho(i-1,j,k)- &
626!^ & OneTwelfth* &
627!^ & (dRx(i,j)+dRx(i-1,j)))- &
628!^ & (dZx(i,j)-dZx(i-1,j))* &
629!^ & (tl_rho(i,j,k)-tl_rho(i-1,j,k)- &
630!^ & OneTwelfth* &
631!^ & (tl_dRx(i,j)+tl_dRx(i-1,j)))))))
632!^
633 adfac=on_u(i,j)*0.5_r8*ad_ru(i,j,k,nrhs)
634 adfac1=adfac*(p(i-1,j,k)-p(i,j,k)- &
635 & halfgrho* &
636 & ((rho(i,j,k)+rho(i-1,j,k))* &
637 & (z_r(i,j,k)-z_r(i-1,j,k))- &
638 & onefifth* &
639 & ((drx(i,j)-drx(i-1,j))* &
640 & (z_r(i,j,k)-z_r(i-1,j,k)- &
641 & onetwelfth* &
642 & (dzx(i,j)+dzx(i-1,j)))- &
643 & (dzx(i,j)-dzx(i-1,j))* &
644 & (rho(i,j,k)-rho(i-1,j,k)- &
645 & onetwelfth* &
646 & (drx(i,j)+drx(i-1,j))))))
647 adfac2=adfac*(hz(i,j,k)+hz(i-1,j,k))
648 adfac3=adfac2*halfgrho
649 adfac4=adfac3*(z_r(i,j,k)-z_r(i-1,j,k))
650 adfac5=adfac3*(rho(i,j,k)+rho(i-1,j,k))
651 adfac6=adfac3*onefifth
652 adfac7=adfac6*(z_r(i,j,k)-z_r(i-1,j,k)- &
653 & onetwelfth*(dzx(i,j)+dzx(i-1,j)))
654 adfac8=adfac6*(drx(i,j)-drx(i-1,j))
655 adfac9=adfac8*onetwelfth
656 adfac10=adfac6*(rho(i,j,k)-rho(i-1,j,k)- &
657 & onetwelfth*(drx(i,j)+drx(i-1,j)))
658 adfac11=adfac6*(dzx(i,j)-dzx(i-1,j))
659 adfac12=adfac11*onetwelfth
660 ad_hz(i-1,j,k)=ad_hz(i-1,j,k)+adfac1
661 ad_hz(i ,j,k)=ad_hz(i ,j,k)+adfac1
662 ad_p(i-1,j,k)=ad_p(i-1,j,k)+adfac2
663 ad_p(i ,j,k)=ad_p(i ,j,k)-adfac2
664 ad_rho(i-1,j,k)=ad_rho(i-1,j,k)-adfac4+adfac11
665 ad_rho(i ,j,k)=ad_rho(i ,j,k)-adfac4-adfac11
666 ad_z_r(i-1,j,k)=ad_z_r(i-1,j,k)+adfac5-adfac8
667 ad_z_r(i ,j,k)=ad_z_r(i ,j,k)-adfac5+adfac8
668 ad_drx(i-1,j)=ad_drx(i-1,j)-adfac7+adfac12
669 ad_drx(i ,j)=ad_drx(i ,j)+adfac7+adfac12
670 ad_dzx(i-1,j)=ad_dzx(i-1,j)-adfac9+adfac10
671 ad_dzx(i ,j)=ad_dzx(i ,j)-adfac9-adfac10
672 ad_ru(i,j,k,nrhs)=0.0_r8
673 END DO
674 END DO
675!
676 DO j=jstr,jend
677 DO i=istru-1,iend
678 cff1=2.0_r8*fc(i,j)*fc(i+1,j)
679 IF (cff1.gt.eps) THEN
680 cff2=1.0_r8/(fc(i,j)+fc(i+1,j))
681!^ tl_dRx(i,j)=tl_cff1*cff2+cff1*tl_cff2
682!^
683 ad_cff1=ad_cff1+cff2*ad_drx(i,j)
684 ad_cff2=ad_cff2+cff1*ad_drx(i,j)
685 ad_drx(i,j)=0.0_r8
686!^ tl_cff2=-cff2*cff2*(tl_FC(i,j)+tl_FC(i+1,j))
687!^
688 adfac=-cff2*cff2*ad_cff2
689 ad_fc(i ,j)=ad_fc(i ,j)+adfac
690 ad_fc(i+1,j)=ad_fc(i+1,j)+adfac
691 ad_cff2=0.0_r8
692 ELSE
693!^ tl_dRx(i,j)=0.0_r8
694!^
695 ad_drx(i,j)=0.0_r8
696 END IF
697!^ tl_cff1=2.0_r8*(tl_FC(i,j)*FC(i+1,j)+ &
698!^ & FC(i,j)*tl_FC(i+1,j)
699!^
700 adfac=2.0_r8*ad_cff1
701 ad_fc(i ,j)=ad_fc(i ,j)+fc(i+1,j)*adfac
702 ad_fc(i+1,j)=ad_fc(i+1,j)+fc(i ,j)*adfac
703 ad_cff1=0.0_r8
704
705 cff=2.0_r8*aux(i,j)*aux(i+1,j)
706 IF (cff.gt.eps) THEN
707 cff1=1.0_r8/(aux(i,j)+aux(i+1,j))
708!^ tl_dZx(i,j)=tl_cff*cff1+cff*tl_cff1
709!^
710 ad_cff=ad_cff+cff1*ad_dzx(i,j)
711 ad_cff1=ad_cff1+cff*ad_dzx(i,j)
712 ad_dzx(i,j)=0.0_r8
713!^ tl_cff1=-cff1*cff1*(tl_aux(i,j)+tl_aux(i+1,j))
714!^
715 adfac=-cff1*cff1*ad_cff1
716 ad_aux(i ,j)=ad_aux(i ,j)+adfac
717 ad_aux(i+1,j)=ad_aux(i+1,j)+adfac
718 ad_cff1=0.0_r8
719 ELSE
720!^ tl_dZx(i,j)=0.0_r8
721!^
722 ad_dzx(i,j)=0.0_r8
723 END IF
724!^ tl_cff=2.0_r8*(tl_aux(i,j)*aux(i+1,j)+ &
725!^ & aux(i,j)*tl_aux(i+1,j)
726!^
727 adfac=2.0_r8*ad_cff
728 ad_aux(i ,j)=ad_aux(i ,j)+aux(i+1,j)*adfac
729 ad_aux(i+1,j)=ad_aux(i+1,j)+aux(i ,j)*adfac
730 ad_cff=0.0_r8
731 END DO
732 END DO
733 DO j=jstr,jend
734 DO i=istru-1,iend+1
735#ifdef MASKING
736!^ tl_FC(i,j)=tl_FC(i,j)*umask(i,j)
737!^
738 ad_fc(i,j)=ad_fc(i,j)*umask(i,j)
739#endif
740!^ tl_FC(i,j)=tl_rho(i,j,k)-tl_rho(i-1,j,k)
741!^
742 ad_rho(i-1,j,k)=ad_rho(i-1,j,k)-ad_fc(i,j)
743 ad_rho(i ,j,k)=ad_rho(i ,j,k)+ad_fc(i,j)
744 ad_fc(i,j)=0.0_r8
745#ifdef MASKING
746!^ tl_aux(i,j)=tl_aux(i,j)*umask(i,j)
747!^
748 ad_aux(i,j)=ad_aux(i,j)*umask(i,j)
749#endif
750!^ tl_aux(i,j)=tl_z_r(i,j,k)-tl_z_r(i-1,j,k)
751!^
752 ad_z_r(i-1,j,k)=ad_z_r(i-1,j,k)-ad_aux(i,j)
753 ad_z_r(i ,j,k)=ad_z_r(i ,j,k)+ad_aux(i,j)
754 ad_aux(i,j)=0.0_r8
755 END DO
756 END DO
757 END DO
758!
759!-----------------------------------------------------------------------
760! Adjoint of Preliminary step (same for XI- and ETA-components):
761!-----------------------------------------------------------------------
762!
763 DO j=jstrv-1,jend
764 DO k=1,n(ng)-1
765 DO i=istru-1,iend
766 dr(i,k)=rho(i,j,k+1)-rho(i,j,k)
767 dz(i,k)=z_r(i,j,k+1)-z_r(i,j,k)
768 dr1(i,k)=dr(i,k)
769 dz1(i,k)=dz(i,k)
770 END DO
771 END DO
772 DO i=istru-1,iend
773 dr(i,n(ng))=dr(i,n(ng)-1)
774 dr1(i,n(ng))=dr(i,n(ng))
775 dz(i,n(ng))=dz(i,n(ng)-1)
776 dz1(i,n(ng))=dz(i,n(ng))
777 dr(i,0)=dr(i,1)
778 dr1(i,0)=dr(i,0)
779 dz(i,0)=dz(i,1)
780 dz1(i,0)=dz(i,0)
781 END DO
782 DO k=n(ng),1,-1
783 DO i=istru-1,iend
784 cff=2.0_r8*dr(i,k)*dr(i,k-1)
785 IF (cff.gt.eps) THEN
786 dr(i,k)=cff/(dr(i,k)+dr(i,k-1))
787 ELSE
788 dr(i,k)=0.0_r8
789 END IF
790 dz(i,k)=2.0_r8*dz(i,k)*dz(i,k-1)/(dz(i,k)+dz(i,k-1))
791 END DO
792 END DO
793!
794 DO k=1,n(ng)-1
795 DO i=istru-1,iend
796!^ tl_P(i,j,k)=tl_P(i,j,k+1)+tl_cff
797!^
798 ad_cff=ad_cff+ad_p(i,j,k)
799!^ tl_cff=HalfGRho*((tl_rho(i,j,k+1)+tl_rho(i,j,k))* &
800!^ & (z_r(i,j,k+1)-z_r(i,j,k))+ &
801!^ & (rho(i,j,k+1)+rho(i,j,k))* &
802!^ & (tl_z_r(i,j,k+1)-tl_z_r(i,j,k))- &
803!^ & OneFifth* &
804!^ & ((tl_dR(i,k+1)-tl_dR(i,k))* &
805!^ & (z_r(i,j,k+1)-z_r(i,j,k)- &
806!^ & OneTwelfth* &
807!^ & (dZ(i,k+1)+dZ(i,k)))+ &
808!^ & (dR(i,k+1)-dR(i,k))* &
809!^ & (tl_z_r(i,j,k+1)-tl_z_r(i,j,k)- &
810!^ & OneTwelfth* &
811!^ & (tl_dZ(i,k+1)+tl_dZ(i,k)))- &
812!^ & (tl_dZ(i,k+1)-tl_dZ(i,k))* &
813!^ & (rho(i,j,k+1)-rho(i,j,k)- &
814!^ & OneTwelfth* &
815!^ & (dR(i,k+1)+dR(i,k)))- &
816!^ & (dZ(i,k+1)-dZ(i,k))* &
817!^ & (tl_rho(i,j,k+1)-tl_rho(i,j,k)- &
818!^ & OneTwelfth* &
819!^ & (tl_dR(i,k+1)+tl_dR(i,k)))))
820!^
821 adfac=halfgrho*ad_cff
822 adfac1=adfac*(z_r(i,j,k+1)-z_r(i,j,k))
823 adfac2=adfac*(rho(i,j,k+1)+rho(i,j,k))
824 adfac3=adfac*onefifth
825 adfac4=adfac3*(z_r(i,j,k+1)-z_r(i,j,k)- &
826 & onetwelfth*(dz(i,k+1)+dz(i,k)))
827 adfac5=adfac3*(dr(i,k+1)-dr(i,k))
828 adfac6=adfac5*onetwelfth
829 adfac7=adfac3*(rho(i,j,k+1)-rho(i,j,k)- &
830 & onetwelfth*(dr(i,k+1)+dr(i,k)))
831 adfac8=adfac3*(dz(i,k+1)-dz(i,k))
832 adfac9=adfac8*onetwelfth
833 ad_p(i,j,k+1)=ad_p(i,j,k+1)+ad_p(i,j,k)
834 ad_rho(i,j,k )=ad_rho(i,j,k )+adfac1-adfac8
835 ad_rho(i,j,k+1)=ad_rho(i,j,k+1)+adfac1+adfac8
836 ad_z_r(i,j,k )=ad_z_r(i,j,k )-adfac2+adfac5
837 ad_z_r(i,j,k+1)=ad_z_r(i,j,k+1)+adfac2-adfac5
838 ad_dr(i,k )=ad_dr(i,k )+adfac4-adfac9
839 ad_dr(i,k+1)=ad_dr(i,k+1)-adfac4-adfac9
840 ad_dz(i,k )=ad_dz(i,k )+adfac6-adfac7
841 ad_dz(i,k+1)=ad_dz(i,k+1)+adfac6+adfac7
842 ad_cff=0.0_r8
843 END DO
844 END DO
845 DO i=istru-1,iend
846 cff1=1.0_r8/(z_r(i,j,n(ng))-z_r(i,j,n(ng)-1))
847 cff2=0.5_r8*(rho(i,j,n(ng))-rho(i,j,n(ng)-1))* &
848 & (z_w(i,j,n(ng))-z_r(i,j,n(ng)))*cff1
849#ifdef TIDE_GENERATING_FORCES
850!^ tl_P(i,j,N(ng))=tl_P(i,j,N(ng))-g*tl_eq_tide(i,j)
851!^
852 ad_eq_tide(i,j)=ad_eq_tide(i,j)-g*ad_p(i,j,n(ng))
853#endif
854!^ tl_P(i,j,N(ng))=g*tl_z_w(i,j,N(ng))+ &
855!^ & GRho*((tl_rho(i,j,N(ng))+tl_cff2)* &
856!^ & (z_w(i,j,N(ng))-z_r(i,j,N(ng)))+ &
857!^ & (rho(i,j,N(ng))+cff2)* &
858!^ & (tl_z_w(i,j,N(ng))-tl_z_r(i,j,N(ng))))
859!^
860 adfac=grho*ad_p(i,j,n(ng))
861 adfac1=adfac*(z_w(i,j,n(ng))-z_r(i,j,n(ng)))
862 adfac2=adfac*(rho(i,j,n(ng))+cff2)
863 ad_z_r(i,j,n(ng))=ad_z_r(i,j,n(ng))-adfac2
864 ad_z_w(i,j,n(ng))=ad_z_w(i,j,n(ng))+adfac2+ &
865 & g*ad_p(i,j,n(ng))
866 ad_rho(i,j,n(ng))=ad_rho(i,j,n(ng))+adfac1
867 ad_cff2=ad_cff2+adfac1
868 ad_p(i,j,n(ng))=0.0_r8
869!^ tl_cff2=0.5_r8*((tl_rho(i,j,N(ng))-tl_rho(i,j,N(ng)-1))* &
870!^ & (z_w(i,j,N(ng))-z_r(i,j,N(ng)))*cff1+ &
871!^ & (rho(i,j,N(ng))-rho(i,j,N(ng)-1))* &
872!^ & ((tl_z_w(i,j,N(ng))-tl_z_r(i,j,N(ng)))*cff1+ &
873!^ & (z_w(i,j,N(ng))-z_r(i,j,N(ng)))*tl_cff1))
874!^
875 adfac=0.5_r8*ad_cff2
876 adfac1=adfac*(z_w(i,j,n(ng))-z_r(i,j,n(ng)))*cff1
877 adfac2=adfac*(rho(i,j,n(ng))-rho(i,j,n(ng)-1))
878 adfac3=adfac2*cff1
879 ad_rho(i,j,n(ng)-1)=ad_rho(i,j,n(ng)-1)-adfac1
880 ad_rho(i,j,n(ng) )=ad_rho(i,j,n(ng) )+adfac1
881 ad_z_r(i,j,n(ng))=ad_z_r(i,j,n(ng))-adfac3
882 ad_z_w(i,j,n(ng))=ad_z_w(i,j,n(ng))+adfac3
883 ad_cff1=ad_cff1+(z_w(i,j,n(ng))-z_r(i,j,n(ng)))*adfac2
884 ad_cff2=0.0_r8
885!^ tl_cff1=-cff1*cff1*(tl_z_r(i,j,N(ng))-tl_z_r(i,j,N(ng)-1))
886!^
887 adfac=-cff1*cff1*ad_cff1
888 ad_z_r(i,j,n(ng)-1)=ad_z_r(i,j,n(ng)-1)-adfac
889 ad_z_r(i,j,n(ng) )=ad_z_r(i,j,n(ng) )+adfac
890 ad_cff1=0.0_r8
891 END DO
892!
893! The forward code has recursive statements, so we need to dR1(ik) and
894! dZ1(i,k) for BOTH terms in denominator of adjoint code.
895!
896 DO k=1,n(ng)
897 DO i=istru-1,iend
898!^ tl_dZ(i,k)=(2.0_r8*(tl_dZ(i,k)*dZ1(i,k-1)+ &
899!^ & dZ1(i,k)*tl_dZ(i,k-1))-
900!^ & dZ(i,k)*(tl_dZ(i,k)+tl_dZ(i,k-1)))/
901!^ & (dZ1(i,k)+dZ1(i,k-1))
902!^ recursive
903 adfac=ad_dz(i,k)/(dz1(i,k)+dz1(i,k-1))
904 adfac1=adfac*2.0_r8
905 adfac2=adfac*dz(i,k)
906 ad_dz(i,k-1)=ad_dz(i,k-1)+dz1(i,k)*adfac1-adfac2
907 ad_dz(i,k )=dz1(i,k-1)*adfac1-adfac2
908
909 cff=2.0_r8*dr1(i,k)*dr1(i,k-1)
910 IF (cff.gt.eps) THEN
911!^ tl_dR(i,k)=(tl_cff-dR(i,k)*(tl_dR(i,k)+tl_dR(i,k-1)))/ &
912!^ & (dR1(i,k)+dR1(i,k-1))
913!^
914 adfac=ad_dr(i,k)/(dr1(i,k)+dr1(i,k-1))
915 adfac1=adfac*dr(i,k)
916 ad_dr(i,k-1)=ad_dr(i,k-1)-adfac1
917 ad_dr(i,k )=-adfac1
918 ad_cff=ad_cff+adfac
919 ELSE
920!^ tl_dR(i,k)=0.0_r8
921!^
922 ad_dr(i,k)=0.0_r8
923 END IF
924!^ tl_cff=2.0_r8*(tl_dR(i,k)*dR1(i,k-1)+ &
925!^ & dR1(i,k)*tl_dR(i,k-1))
926!^
927 adfac=2.0_r8*ad_cff
928 ad_dr(i,k-1)=ad_dr(i,k-1)+dr1(i,k )*adfac
929 ad_dr(i,k )=ad_dr(i,k )+dr1(i,k-1)*adfac
930 ad_cff=0.0_r8
931 END DO
932 END DO
933 DO i=istru-1,iend
934!^ tl_dZ(i,0)=tl_dZ(i,1)
935!^
936 ad_dz(i,1)=ad_dz(i,1)+ad_dz(i,0)
937 ad_dz(i,0)=0.0_r8
938!^ tl_dR(i,0)=tl_dR(i,1)
939!^
940 ad_dr(i,1)=ad_dr(i,1)+ad_dr(i,0)
941 ad_dr(i,0)=0.0_r8
942!^ tl_dZ(i,N(ng))=tl_dZ(i,N(ng)-1)
943!^
944 ad_dz(i,n(ng)-1)=ad_dz(i,n(ng)-1)+ad_dz(i,n(ng))
945 ad_dz(i,n(ng))=0.0_r8
946
947!^ tl_dR(i,N(ng))=tl_dR(i,N(ng)-1)
948!^
949 ad_dr(i,n(ng)-1)=ad_dr(i,n(ng)-1)+ad_dr(i,n(ng))
950 ad_dr(i,n(ng))=0.0_r8
951 END DO
952 DO k=n(ng)-1,1,-1
953 DO i=istru-1,iend
954!^ tl_dZ(i,k)=tl_z_r(i,j,k+1)-tl_z_r(i,j,k)
955!^
956 ad_z_r(i,j,k )=ad_z_r(i,j,k )-ad_dz(i,k)
957 ad_z_r(i,j,k+1)=ad_z_r(i,j,k+1)+ad_dz(i,k)
958 ad_dz(i,k)=0.0_r8
959
960!^ tl_dR(i,k)=tl_rho(i,j,k+1)-tl_rho(i,j,k)
961!^
962 ad_rho(i,j,k )=ad_rho(i,j,k )-ad_dr(i,k)
963 ad_rho(i,j,k+1)=ad_rho(i,j,k+1)+ad_dr(i,k)
964 ad_dr(i,k)=0.0_r8
965 END DO
966 END DO
967 END DO
968!
969 RETURN

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

◆ ad_prsgrd40_tile()

subroutine ad_prsgrd_mod::ad_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(inout) ad_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(inout) ad_z_w,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) rho,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(inout) ad_rho,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) eq_tide,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(inout) ad_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) ad_ru,
real(r8), dimension(lbi:ubi,lbj:ubj,0:n(ng),2), intent(inout) ad_rv )
private

Definition at line 94 of file ad_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# ifdef ATM_PRESS
131 real(r8), intent(in) :: Pair(LBi:,LBj:)
132# endif
133# ifdef TIDE_GENERATING_FORCES
134 real(r8), intent(in) :: eq_tide(LBi:,LBj:)
135 real(r8), intent(inout) :: ad_eq_tide(LBi:,LBj:)
136# endif
137# ifdef DIAGNOSTICS_UV
138!! real(r8), intent(inout) :: DiaRU(LBi:,LBj:,:,:,:)
139!! real(r8), intent(inout) :: DiaRV(LBi:,LBj:,:,:,:)
140# endif
141 real(r8), intent(inout) :: ad_Hz(LBi:,LBj:,:)
142 real(r8), intent(inout) :: ad_z_w(LBi:,LBj:,0:)
143 real(r8), intent(inout) :: ad_rho(LBi:,LBj:,:)
144 real(r8), intent(inout) :: ad_ru(LBi:,LBj:,0:,:)
145 real(r8), intent(inout) :: ad_rv(LBi:,LBj:,0:,:)
146#else
147 real(r8), intent(in) :: om_v(LBi:UBi,LBj:UBj)
148 real(r8), intent(in) :: on_u(LBi:UBi,LBj:UBj)
149 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
150 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
151 real(r8), intent(in) :: rho(LBi:UBi,LBj:UBj,N(ng))
152# ifdef ATM_PRESS
153 real(r8), intent(in) :: Pair(LBi:UBi,LBj:UBj)
154# endif
155# ifdef TIDE_GENERATING_FORCES
156 real(r8), intent(in) :: eq_tide(LBi:UBi,LBj:UBj)
157 real(r8), intent(inout) :: ad_eq_tide(LBi:UBi,LBj:UBj)
158# endif
159# ifdef DIAGNOSTICS_UV
160!! real(r8), intent(inout) :: DiaRU(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
161!! real(r8), intent(inout) :: DiaRV(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
162# endif
163 real(r8), intent(inout) :: ad_Hz(LBi:UBi,LBj:UBj,N(ng))
164 real(r8), intent(inout) :: ad_z_w(LBi:UBi,LBj:UBj,0:N(ng))
165 real(r8), intent(inout) :: ad_rho(LBi:UBi,LBj:UBj,N(ng))
166 real(r8), intent(inout) :: ad_ru(LBi:UBi,LBj:UBj,0:N(ng),2)
167 real(r8), intent(inout) :: ad_rv(LBi:UBi,LBj:UBj,0:N(ng),2)
168#endif
169!
170! Local variable declarations.
171!
172 integer :: i, j, k
173
174 real(r8) :: cff, cff1, dh
175 real(r8) :: ad_dh
176#ifdef ATM_PRESS
177 real(r8) :: OneAtm, fac
178#endif
179 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: ad_FC
180
181 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,N(ng)) :: ad_FX
182
183 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,0:N(ng)) :: P
184 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,0:N(ng)) :: ad_P
185
186#include "set_bounds.h"
187!
188!-----------------------------------------------------------------------
189! Initialize adjoint private variables.
190!-----------------------------------------------------------------------
191!
192 ad_dh=0.0_r8
193 DO j=jmins,jmaxs
194 DO i=imins,imaxs
195 ad_fx=(i,j)=0.0_r8
196 END DO
197 DO k=0,n(ng)
198 DO i=imins,imaxs
199 ad_p(i,j,k)=0.0_r8
200 END DO
201 END DO
202 END DO
203 DO k=0,n(ng)
204 DO i=imins,imaxs
205 ad_fc=(i,k)=0.0_r8
206 END DO
207 END DO
208!
209!-----------------------------------------------------------------------
210! Finite Volume pressure gradient algorithm (Lin, 1997).
211!-----------------------------------------------------------------------
212!
213! Compute BASIC STATE dynamic pressure, P.
214!
215#ifdef ATM_PRESS
216 oneatm=1013.25_r8 ! 1 atm = 1013.25 mb
217 fac=100.0_r8/g
218#endif
219 DO j=jstrv-1,jend
220 DO i=istru-1,iend
221 p(i,j,n(ng))=0.0_r8
222#ifdef ATM_PRESS
223 p(i,j,n(ng))=p(i,j,n(ng))+fac*(pair(i,j)-oneatm)
224#endif
225#ifdef TIDE_GENERATING_FORCES
226 p(i,j,n(ng))=p(i,j,n(ng))-g*eq_tide(i,j)
227#endif
228 END DO
229 DO k=n(ng),1,-1
230 DO i=istru-1,iend
231 p(i,j,k-1)=p(i,j,k)+hz(i,j,k)*rho(i,j,k)
232 END DO
233 END DO
234 END DO
235!
236! Calculate adjoint of pressure gradient in the ETA-direction (m4/s2).
237!
238 j_loop: DO j=jend,jstrv-1,-1
239 IF (j.ge.jstrv) THEN
240 cff=0.5_r8*g
241 cff1=g/rho0
242 DO k=1,n(ng)
243 DO i=istr,iend
244 dh=z_w(i,j,k-1)-z_w(i,j-1,k-1)
245#ifdef DIAGNOSTICS_UV
246!! DiaRV(i,j,k,nrhs,M3pgrd)=rv(i,j,k,nrhs)
247#endif
248!^ tl_rv(i,j,k,nrhs)=(cff*((tl_Hz(i,j-1,k)+ &
249!^ & tl_Hz(i,j ,k))* &
250!^ & (z_w(i,j-1,N(ng))- &
251!^ & z_w(i,j ,N(ng)))+ &
252!^ & (Hz(i,j-1,k)+ &
253!^ & Hz(i,j ,k))* &
254!^ & (tl_z_w(i,j-1,N(ng))- &
255!^ & tl_z_w(i,j ,N(ng))))+ &
256!^ & cff1*(tl_FX(i,j-1,k)- &
257!^ & tl_FX(i,j ,k)+ &
258!^ & tl_FC(i,k )- &
259!^ & tl_FC(i,k-1)))*om_v(i,j)
260!^
261 adfac=om_v(i,j)*ad_rv(i,j,k,nrhs)
262 adfac1=adfac*cff
263 adfac2=adfac1*(z_w(i,j-1,n(ng))- &
264 & z_w(i,j ,n(ng)))
265 adfac3=adfac1*(hz(i,j-1,k)+ &
266 & hz(i,j ,k))
267 adfac4=adfac*cff1
268 ad_hz(i,j-1,k)=ad_hz(i,j-1,k)+adfac2
269 ad_hz(i,j ,k)=ad_hz(i,j ,k)+adfac2
270 ad_z_w(i,j-1,n(ng))=ad_z_w(i,j-1,n(ng))+adfac3
271 ad_z_w(i,j ,n(ng))=ad_z_w(i,j ,n(ng))-adfac3
272 ad_fx(i,j-1,k)=ad_fx(i,j-1,k)+adfac4
273 ad_fx(i,j ,k)=ad_fx(i,j ,k)-adfac4
274 ad_fc(i,k-1)=ad_fc(i,k-1)-adfac4
275 ad_fc(i,k )=ad_fc(i,k )+adfac4
276 ad_rv(i,j,k,nrhs)=0.0_r8
277!^ tl_FC(i,k-1)=0.5_r8* &
278!^ & (tl_dh*(P(i,j,k-1)+P(i,j-1,k-1))+ &
279!^ & dh*(tl_P(i,j,k-1)+tl_P(i,j-1,k-1))
280!^
281 adfac=0.5_r8*ad_fc(i,k-1)
282 adfac1=adfac*dh
283 ad_dh=ad_dh+(p(i,j,k-1)+p(i,j-1,k-1))*adfac
284 ad_p(i,j-1,k-1)=ad_p(i,j-1,k-1)+adfac1
285 ad_p(i,j ,k-1)=ad_p(i,j ,k-1)+adfac1
286 ad_fc(i,k-1)=0.0_r8
287!^ tl_dh=tl_z_w(i,j,k-1)-tl_z_w(i,j-1,k-1)
288!^
289 ad_z_w(i,j-1,k-1)=ad_z_w(i,j-1,k-1)-ad_dh
290 ad_z_w(i,j ,k-1)=ad_z_w(i,j ,k-1)+ad_dh
291 ad_dh=0.0_r8
292 END DO
293 END DO
294 DO i=istr,iend
295!^ tl_FC(i,N(ng))=0.0_r8
296!^
297 ad_fc(i,n(ng))=0.0_r8
298 END DO
299 END IF
300!
301! Calculate the adjoint of pressure gradient in the XI-direction (m4/s2).
302!
303 IF (j.ge.jstr) THEN
304 cff=0.5_r8*g
305 cff1=g/rho0
306 DO k=1,n(ng)
307 DO i=iend,istru,-1
308 dh=z_w(i,j,k-1)-z_w(i-1,j,k-1)
309#ifdef DIAGNOSTICS_UV
310!! DiaRU(i,j,k,nrhs,M3pgrd)=ru(i,j,k,nrhs)
311#endif
312!^ tl_ru(i,j,k,nrhs)=(cff*((tl_Hz(i-1,j,k)+ &
313!^ & tl_Hz(i ,j,k))* &
314!^ & (z_w(i-1,j,N(ng))- &
315!^ & z_w(i ,j,N(ng)))+ &
316!^ & (Hz(i-1,j,k)+ &
317!^ & Hz(i ,j,k))* &
318!^ & (tl_z_w(i-1,j,N(ng))- &
319!^ & tl_z_w(i ,j,N(ng))))+ &
320!^ & cff1*(tl_FX(i-1,j,k)- &
321!^ & tl_FX(i ,j,k)+ &
322!^ & tl_FC(i,k )- &
323!^ & tl_FC(i,k-1)))*on_u(i,j)
324!^
325 adfac=on_u(i,j)*tl_ru(i,j,k,nrhs)
326 adfac1=adfac*cff
327 adfac2=adfac1*(z_w(i-1,j,n(ng))- &
328 & z_w(i ,j,n(ng)))
329 adfac3=adfac1*(hz(i-1,j,k)+ &
330 & hz(i ,j,k))
331 adfac4=adfac*cff1
332 ad_hz(i-1,j,k)=ad_hz(i-1,j,k)+adfac2
333 ad_hz(i ,j,k)=ad_hz(i ,j,k)+adfac2
334 ad_z_w(i-1,j,n(ng))=ad_z_w(i-1,j,n(ng))+adfac3
335 ad_z_w(i ,j,n(ng))=ad_z_w(i ,j,n(ng))-adfac3
336 ad_fx(i-1,j,k)=ad_fx(i-1,j,k)+adfac4
337 ad_fx(i ,j,k)=ad_fx(i ,j,k)-adfac4
338 ad_fc(i,k )=ad_fc(i,k )+adfac4
339 ad_fc(i,k-1)=ad_fc(i,k-1)-adfac4
340 ad_ru(i,j,k,nrhs)=0.0_r8
341!^ tl_FC(i,k-1)=0.5_r8* &
342!^ & (tl_dh*(P(i,j,k-1)+P(i-1,j,k-1))+ &
343!^ & dh*(tl_P(i,j,k-1)+tl_P(i-1,j,k-1)))
344!^
345 adfac=0.5_r8*tl_fc(i,k-1)
346 adfac1=adfac*dh
347 ad_dh=ad_dh+(p(i,j,k-1)+p(i-1,j,k-1))*adfac
348 ad_p(i-1,j,k-1)=ad_p(i-1,j,k-1)+adfac1
349 ad_p(i ,j,k-1)=ad_p(i ,j,k-1)+adfac1
350 ad_fc(i,k-1)=0.0_r8
351!^ tl_dh=tl_z_w(i,j,k-1)-tl_z_w(i-1,j,k-1)
352!^
353 ad_z_w(i-1,j,k-1)=ad_z_w(i-1,j,k-1)-ad_dh
354 ad_z_w(i ,j,k-1)=ad_z_w(i ,j,k-1)+ad_dh
355 ad_dh=0.0_r8
356 END DO
357 END DO
358 DO i=istru,iend
359!^ tl_FC(i,N(ng))=0.0_r8
360 ad_fc(i,n(ng))=0.0_r8
361 END DO
362 END IF
363!
364! Compute adjoint pressure and its vertical integral.
365!
366 DO k=1,n(ng)
367 DO i=istru-1,iend
368!^ tl_FX(i,j,k)=0.5_r8* &
369!^ & (tl_Hz(i,j,k)*(P(i,j,k)+P(i,j,k-1))+ &
370!^ & Hz(i,j,k)*(tl_P(i,j,k)+tl_P(i,j,k-1)))
371!^
372 adfac=0.5_r8*ad_fx(i,j,k)
373 adfac1=adfac*hz(i,j,k)
374 ad_hz(i,j,k)=ad_hz(i,j,k)+(p(i,j,k)+p(i,j,k-1))*adfac
375 ad_p(i,j,k-1)=ad_p(i,j,k-1)+adfac1
376 ad_p(i,j,k )=ad_p(i,j,k )+adfac1
377 ad_fx(i,j,k)=0.0_r8
378!^ tl_P(i,j,k-1)=tl_P(i,j,k)+ &
379!^ & tl_Hz(i,j,k)*rho(i,j,k)+ &
380!^ & Hz(i,j,k)*tl_rho(i,j,k)
381!^
382 ad_p(i,j,k)=ad_p(i,j,k)+ad_p(i,j,k-1)
383 ad_hz(i,j,k)=ad_hz(i,j,k)+rho(i,j,k)*ad_p(i,j,k-1)
384 ad_rho(i,j,k)=ad_rho(i,j,k)+hz(i,j,k)*ad_p(i,j,k-1)
385 END DO
386 END DO
387 DO i=istru-1,iend
388#ifdef TIDE_GENERATING_FORCES
389!^ tl_P(i,j,N(ng))=tl_P(i,j,N(ng))-g*tl_eq_tide(i,j)
390!^
391 ad_eq_tide(i,j)=ad_eq_tide(i,j)-g*ad_p(i,j,n(ng))
392#endif
393!^ tl_P(i,j,N(ng))=0.0_r8
394!^
395 ad_p(i,j,n(ng))=0.0_r8
396 END DO
397 END DO j_loop
398!
399 RETURN

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