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

Functions/Subroutines

subroutine, public ad_t3dmix2 (ng, tile)
 
subroutine ad_t3dmix2_geo_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, nstp, nnew, umask, vmask, umask_wet, vmask_wet, om_v, on_u, pm, pn, hz, ad_hz, z_r, ad_z_r, diff3d_r, diff2, tclm, t, ad_t)
 
subroutine ad_t3dmix2_iso_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, nstp, nnew, umask, vmask, om_v, on_u, pm, pn, hz, ad_hz, z_r, ad_z_r, diff2, pden, ad_pden, tclm, t, ad_t)
 
subroutine ad_t3dmix2_s_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, nstp, nnew, umask, vmask, umask_wet, vmask_wet, hz, ad_hz, pmon_u, pnom_v, pm, pn, diff3d_r, diff2, tclm, t, ad_t)
 

Function/Subroutine Documentation

◆ ad_t3dmix2()

subroutine public ad_t3dmix2_mod::ad_t3dmix2 ( integer, intent(in) ng,
integer, intent(in) tile )

Definition at line 25 of file ad_t3dmix2_geo.h.

26!***********************************************************************
27!
28 USE mod_param
29#ifdef TS_MIX_CLIMA
30 USE mod_clima
31#endif
32#ifdef DIAGNOSTICS_TS
33!! USE mod_diags
34#endif
35 USE mod_grid
36 USE mod_mixing
37 USE mod_ocean
38 USE mod_stepping
39!
40! Imported variable declarations.
41!
42 integer, intent(in) :: ng, tile
43!
44! Local variable declarations.
45!
46 character (len=*), parameter :: MyFile = &
47 & __FILE__
48!
49#include "tile.h"
50!
51#ifdef PROFILE
52 CALL wclock_on (ng, iadm, 25, __line__, myfile)
53#endif
54 CALL ad_t3dmix2_geo_tile (ng, tile, &
55 & lbi, ubi, lbj, ubj, &
56 & imins, imaxs, jmins, jmaxs, &
57 & nrhs(ng), nstp(ng), nnew(ng), &
58#ifdef MASKING
59 & grid(ng) % umask, &
60 & grid(ng) % vmask, &
61#endif
62#ifdef WET_DRY_NOT_YET
63 & grid(ng) % umask_wet, &
64 & grid(ng) % vmask_wet, &
65#endif
66 & grid(ng) % om_v, &
67 & grid(ng) % on_u, &
68 & grid(ng) % pm, &
69 & grid(ng) % pn, &
70 & grid(ng) % Hz, &
71 & grid(ng) % ad_Hz, &
72 & grid(ng) % z_r, &
73 & grid(ng) % ad_z_r, &
74#ifdef DIFF_3DCOEF
75 & mixing(ng) % diff3d_r, &
76#else
77 & mixing(ng) % diff2, &
78#endif
79#ifdef TS_MIX_CLIMA
80 & clima(ng) % tclm, &
81#endif
82#ifdef DIAGNOSTICS_TS
83!! & DIAGS(ng) % DiaTwrk, &
84#endif
85 & ocean(ng) % t, &
86 & ocean(ng) % ad_t)
87#ifdef PROFILE
88 CALL wclock_off (ng, iadm, 25, __line__, myfile)
89#endif
90!
91 RETURN
type(t_clima), dimension(:), allocatable clima
Definition mod_clima.F:153
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_mixing), dimension(:), allocatable mixing
Definition mod_mixing.F:399
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
integer, parameter iadm
Definition mod_param.F:665
integer, dimension(:), allocatable nrhs
integer, dimension(:), allocatable nnew
integer, dimension(:), allocatable nstp
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_t3dmix2_geo_tile(), mod_clima::clima, mod_grid::grid, mod_param::iadm, mod_mixing::mixing, mod_stepping::nnew, mod_stepping::nrhs, mod_stepping::nstp, 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_t3dmix2_geo_tile()

subroutine ad_t3dmix2_mod::ad_t3dmix2_geo_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,
integer, intent(in) nstp,
integer, intent(in) nnew,
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) umask_wet,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) vmask_wet,
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), intent(in) pm,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) pn,
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,n(ng)), intent(in) diff3d_r,
real(r8), dimension(lbi:ubi,lbj:ubj,nt(ng)), intent(in) diff2,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),nt(ng)), intent(in) tclm,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),3,nt(ng)), intent(in) t,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),3,nt(ng)), intent(inout) ad_t )
private

Definition at line 95 of file ad_t3dmix2_geo.h.

120!***********************************************************************
121!
122 USE mod_param
123 USE mod_scalars
124!
125! Imported variable declarations.
126!
127 integer, intent(in) :: ng, tile
128 integer, intent(in) :: LBi, UBi, LBj, UBj
129 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
130 integer, intent(in) :: nrhs, nstp, nnew
131
132#ifdef ASSUMED_SHAPE
133# ifdef MASKING
134 real(r8), intent(in) :: umask(LBi:,LBj:)
135 real(r8), intent(in) :: vmask(LBi:,LBj:)
136# endif
137# ifdef WET_DRY_NOT_YET
138 real(r8), intent(in) :: umask_wet(LBi:,LBj:)
139 real(r8), intent(in) :: vmask_wet(LBi:,LBj:)
140# endif
141# ifdef DIFF_3DCOEF
142 real(r8), intent(in) :: diff3d_r(LBi:,LBj:,:)
143# else
144 real(r8), intent(in) :: diff2(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) :: pm(LBi:,LBj:)
149 real(r8), intent(in) :: pn(LBi:,LBj:)
150 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
151 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
152 real(r8), intent(in) :: t(LBi:,LBj:,:,:,:)
153# ifdef TS_MIX_CLIMA
154 real(r8), intent(in) :: tclm(LBi:,LBj:,:,:)
155# endif
156# ifdef DIAGNOSTICS_TS
157!! real(r8), intent(inout) :: DiaTwrk(LBi:,LBj:,:,:,:)
158# endif
159 real(r8), intent(inout) :: ad_Hz(LBi:,LBj:,:)
160 real(r8), intent(inout) :: ad_z_r(LBi:,LBj:,:)
161 real(r8), intent(inout) :: ad_t(LBi:,LBj:,:,:,:)
162#else
163# ifdef MASKING
164 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
165 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
166# endif
167# ifdef WET_DRY
168 real(r8), intent(in) :: umask_wet(LBi:UBi,LBj:UBj)
169 real(r8), intent(in) :: vmask_wet(LBi:UBi,LBj:UBj)
170# endif
171# ifdef DIFF_3DCOEF
172 real(r8), intent(in) :: diff3d_r(LBi:UBi,LBj:UBj,N(ng))
173# else
174 real(r8), intent(in) :: diff2(LBi:UBi,LBj:UBj,NT(ng))
175# endif
176 real(r8), intent(in) :: om_v(LBi:UBi,LBj:UBj)
177 real(r8), intent(in) :: on_u(LBi:UBi,LBj:UBj)
178 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
179 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
180 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
181 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
182 real(r8), intent(in) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
183# ifdef TS_MIX_CLIMA
184 real(r8), intent(in) :: tclm(LBi:UBi,LBj:UBj,N(ng),NT(ng))
185# endif
186# ifdef DIAGNOSTICS_TS
187!! real(r8), intent(inout) :: DiaTwrk(LBi:UBi,LBj:UBj,N(ng),NT(ng), &
188!! & NDT)
189# endif
190 real(r8), intent(inout) :: ad_Hz(LBi:UBi,LBj:UBj,N(ng))
191 real(r8), intent(inout) :: ad_z_r(LBi:UBi,LBj:UBj,N(ng))
192 real(r8), intent(inout) :: ad_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
193#endif
194!
195! Local variable declarations.
196!
197 integer :: i, itrc, j, k, kk, kt, k1, k1b, k2, k2b
198
199 real(r8) :: cff, cff1, cff2, cff3, cff4
200 real(r8) :: ad_cff, ad_cff1, ad_cff2, ad_cff3, ad_cff4
201 real(r8) :: adfac, adfac1, adfac2, adfac3, adfac4
202
203 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_FE
204 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_FX
205
206 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dTdz
207 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dTdx
208 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dTde
209 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dZdx
210 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dZde
211
212 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: ad_FS
213 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: ad_dTdz
214 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: ad_dTdx
215 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: ad_dTde
216 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: ad_dZdx
217 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: ad_dZde
218
219#include "set_bounds.h"
220!
221!-----------------------------------------------------------------------
222! Initialize adjoint private variables.
223!-----------------------------------------------------------------------
224!
225 ad_cff=0.0_r8
226 ad_cff1=0.0_r8
227 ad_cff2=0.0_r8
228 ad_cff3=0.0_r8
229 ad_cff4=0.0_r8
230
231 ad_fe(imins:imaxs,jmins:jmaxs)=0.0_r8
232 ad_fx(imins:imaxs,jmins:jmaxs)=0.0_r8
233
234 ad_fs(imins:imaxs,jmins:jmaxs,1:2)=0.0_r8
235
236 ad_dtdz(imins:imaxs,jmins:jmaxs,1:2)=0.0_r8
237 ad_dtdx(imins:imaxs,jmins:jmaxs,1:2)=0.0_r8
238 ad_dtde(imins:imaxs,jmins:jmaxs,1:2)=0.0_r8
239 ad_dzdx(imins:imaxs,jmins:jmaxs,1:2)=0.0_r8
240 ad_dzde(imins:imaxs,jmins:jmaxs,1:2)=0.0_r8
241!
242!----------------------------------------------------------------------
243! Compute horizontal harmonic diffusion along geopotential surfaces.
244!----------------------------------------------------------------------
245!
246! Compute horizontal and vertical gradients. Notice the recursive
247! blocking sequence. The vertical placement of the gradients is:
248!
249! dTdx,dTde(:,:,k1) k rho-points
250! dTdx,dTde(:,:,k2) k+1 rho-points
251! FS,dTdz(:,:,k1) k-1/2 W-points
252! FS,dTdz(:,:,k2) k+1/2 W-points
253
254#ifdef TS_MIX_STABILITY
255!
256! In order to increase stability, the biharmonic operator is applied
257! as: 3/4 t(:,:,:,nrhs,:) + 1/4 t(:,:,:,nstp,:).
258#endif
259!
260! Compute adjoint of starting values of k1 and k2.
261!
262 t_loop : DO itrc=1,nt(ng)
263 k1=2
264 k2=1
265 DO k=0,n(ng)
266!!
267!! Note: The following code is equivalent to
268!!
269!! kt=k1
270!! k1=k2
271!! k2=kt
272!!
273!! We use the adjoint of above code.
274!!
275 k1=k2
276 k2=3-k1
277 END DO
278 k_loop : DO k=n(ng),0,-1
279!
280! Compute required BASIC STATE fields. Need to look forward in
281! recursive kk index.
282!
283 k2b=1
284 DO kk=0,k
285 k1b=k2b
286 k2b=3-k1b
287!
288! Compute components of the rotated tracer flux (T m3/s) along
289! geopotential surfaces (required BASIC STATE fields).
290!
291 IF (kk.lt.n(ng)) THEN
292 DO j=jstr,jend
293 DO i=istr,iend+1
294 cff=0.5_r8*(pm(i,j)+pm(i-1,j))
295#ifdef MASKING
296 cff=cff*umask(i,j)
297#endif
298#ifdef WET_DRY_NOT_YET
299 cff=cff*umask_wet(i,j)
300#endif
301 dzdx(i,j,k2b)=cff*(z_r(i ,j,kk+1)- &
302 & z_r(i-1,j,kk+1))
303#if defined TS_MIX_STABILITY
304 dtdx(i,j,k2b)=cff*(0.75_r8*(t(i ,j,kk+1,nrhs,itrc)- &
305 & t(i-1,j,kk+1,nrhs,itrc))+ &
306 & 0.25_r8*(t(i ,j,kk+1,nstp,itrc)- &
307 & t(i-1,j,kk+1,nstp,itrc)))
308#elif defined TS_MIX_CLIMA
309 IF (ltracerclm(itrc,ng)) THEN
310 dtdx(i,j,k2b)=cff*((t(i ,j,kk+1,nrhs,itrc)- &
311 & tclm(i ,j,kk+1,itrc))- &
312 & (t(i-1,j,kk+1,nrhs,itrc)- &
313 & tclm(i-1,j,kk+1,itrc)))
314 ELSE
315 dtdx(i,j,k2b)=cff*(t(i ,j,kk+1,nrhs,itrc)- &
316 & t(i-1,j,kk+1,nrhs,itrc))
317 END IF
318#else
319 dtdx(i,j,k2b)=cff*(t(i ,j,kk+1,nrhs,itrc)- &
320 & t(i-1,j,kk+1,nrhs,itrc))
321#endif
322 END DO
323 END DO
324 IF (kk.eq.0) THEN
325 DO j=jstr,jend
326 DO i=istr,iend+1
327 dzdx(i,j,k1b)=0.0_r8
328 dtdx(i,j,k1b)=0.0_r8
329 END DO
330 END DO
331 END IF
332 DO j=jstr,jend+1
333 DO i=istr,iend
334 cff=0.5_r8*(pn(i,j)+pn(i,j-1))
335#ifdef MASKING
336 cff=cff*vmask(i,j)
337#endif
338#ifdef WET_DRY_NOT_YET
339 cff=cff*vmask_wet(i,j)
340#endif
341 dzde(i,j,k2b)=cff*(z_r(i,j ,kk+1)- &
342 & z_r(i,j-1,kk+1))
343#if defined TS_MIX_STABILITY
344 dtde(i,j,k2b)=cff*(0.75_r8*(t(i,j ,kk+1,nrhs,itrc)- &
345 & t(i,j-1,kk+1,nrhs,itrc))+ &
346 & 0.25_r8*(t(i,j ,kk+1,nstp,itrc)- &
347 & t(i,j-1,kk+1,nstp,itrc)))
348#elif defined TS_MIX_CLIMA
349 IF (ltracerclm(itrc,ng)) THEN
350 dtde(i,j,k2b)=cff*((t(i,j ,kk+1,nrhs,itrc)- &
351 & tclm(i,j ,kk+1,itrc))- &
352 & (t(i,j-1,kk+1,nrhs,itrc)- &
353 & tclm(i,j-1,kk+1,itrc)))
354 ELSE
355 dtde(i,j,k2b)=cff*(t(i,j ,kk+1,nrhs,itrc)- &
356 & t(i,j-1,kk+1,nrhs,itrc))
357 END IF
358#else
359 dtde(i,j,k2b)=cff*(t(i,j ,kk+1,nrhs,itrc)- &
360 & t(i,j-1,kk+1,nrhs,itrc))
361#endif
362 END DO
363 END DO
364 IF (kk.eq.0) THEN
365 DO j=jstr,jend+1
366 DO i=istr,iend
367 dzde(i,j,k1b)=0.0_r8
368 dtde(i,j,k1b)=0.0_r8
369 END DO
370 END DO
371 END IF
372 END IF
373 IF ((kk.eq.0).or.(kk.eq.n(ng))) THEN
374 DO j=jstr-1,jend+1
375 DO i=istr-1,iend+1
376 dtdz(i,j,k2b)=0.0_r8
377 END DO
378 END DO
379 IF (kk.eq.0) THEN
380 DO j=jstr-1,jend+1
381 DO i=istr-1,iend+1
382 dtdz(i,j,k1b)=0.0_r8
383 END DO
384 END DO
385 END IF
386 ELSE
387 DO j=jstr-1,jend+1
388 DO i=istr-1,iend+1
389 cff=1.0_r8/(z_r(i,j,kk+1)-z_r(i,j,kk))
390#if defined TS_MIX_STABILITY
391 dtdz(i,j,k2b)=cff*(0.75_r8*(t(i,j,kk+1,nrhs,itrc)- &
392 & t(i,j,kk ,nrhs,itrc))+ &
393 & 0.25_r8*(t(i,j,kk+1,nstp,itrc)- &
394 & t(i,j,kk ,nstp,itrc)))
395#elif defined TS_MIX_CLIMA
396 IF (ltracerclm(itrc,ng)) THEN
397 dtdz(i,j,k2b)=cff*((t(i,j,kk+1,nrhs,itrc)- &
398 & tclm(i,j,kk+1,itrc))- &
399 & (t(i,j,kk ,nrhs,itrc)- &
400 & tclm(i,j,kk ,itrc)))
401 ELSE
402 dtdz(i,j,k2b)=cff*(t(i,j,kk+1,nrhs,itrc)- &
403 & t(i,j,kk ,nrhs,itrc))
404 END IF
405#else
406 dtdz(i,j,k2b)=cff*(t(i,j,kk+1,nrhs,itrc)- &
407 & t(i,j,kk ,nrhs,itrc))
408#endif
409 END DO
410 END DO
411 END IF
412 END DO
413!
414 IF (k.gt.0) THEN
415!
416! Time-step adjoint harmonic, geopotential diffusion term.
417!
418 DO j=jstr,jend
419 DO i=istr,iend
420#ifdef DIAGNOSTICS_TS
421!! DiaTwrk(i,j,k,itrc,iThdif)=cff
422#endif
423!^ tl_t(i,j,k,nnew,itrc)=tl_t(i,j,k,nnew,itrc)+tl_cff
424!^
425 ad_cff=ad_cff+ad_t(i,j,k,nnew,itrc)
426!^ tl_cff=dt(ng)*pm(i,j)*pn(i,j)* &
427!^ & (tl_FX(i+1,j)-tl_FX(i,j)+ &
428!^ & tl_FE(i,j+1)-tl_FE(i,j))+ &
429!^ & dt(ng)*(tl_FS(i,j,k2)-tl_FS(i,j,k1))
430!^
431 adfac=dt(ng)*ad_cff
432 adfac1=adfac*pm(i,j)*pn(i,j)
433 ad_fs(i,j,k2)=ad_fs(i,j,k2)+adfac
434 ad_fs(i,j,k1)=ad_fs(i,j,k1)-adfac
435 ad_fe(i,j )=ad_fe(i,j )-adfac1
436 ad_fe(i,j+1)=ad_fe(i,j+1)+adfac1
437 ad_fx(i ,j)=ad_fx(i ,j)-adfac1
438 ad_fx(i+1,j)=ad_fx(i+1,j)+adfac1
439 ad_cff=0.0_r8
440 END DO
441 END DO
442!
443! Compute components of the rotated tracer flux (T m4/s) along
444! geopotential surfaces.
445!
446 IF (k.lt.n(ng)) THEN
447 DO j=jstr,jend
448 DO i=istr,iend
449#ifdef DIFF_3DCOEF
450 cff=0.5_r8*diff3d_r(i,j,k)
451#else
452 cff=0.5_r8*diff2(i,j,itrc)
453#endif
454 cff1=min(dzde(i,j ,k1),0.0_r8)
455 cff2=min(dzde(i,j+1,k2),0.0_r8)
456 cff3=max(dzde(i,j ,k2),0.0_r8)
457 cff4=max(dzde(i,j+1,k1),0.0_r8)
458!^ tl_FS(i,j,k2)=tl_FS(i,j,k2)+ &
459!^ & cff* &
460!^ & (tl_cff1*(cff1*dTdz(i,j,k2)- &
461!^ & dTde(i,j ,k1))+ &
462!^ & tl_cff2*(cff2*dTdz(i,j,k2)- &
463!^ & dTde(i,j+1,k2))+ &
464!^ & tl_cff3*(cff3*dTdz(i,j,k2)- &
465!^ & dTde(i,j ,k2))+ &
466!^ & tl_cff4*(cff4*dTdz(i,j,k2)- &
467!^ & dTde(i,j+1,k1))+ &
468!^ & cff1*(tl_cff1*dTdz(i,j,k2)+ &
469!^ & cff1*tl_dTdz(i,j,k2)- &
470!^ & tl_dTde(i,j ,k1))+ &
471!^ & cff2*(tl_cff2*dTdz(i,j,k2)+ &
472!^ & cff2*tl_dTdz(i,j,k2)- &
473!^ & tl_dTde(i,j+1,k2))+ &
474!^ & cff3*(tl_cff3*dTdz(i,j,k2)+ &
475!^ & cff3*tl_dTdz(i,j,k2)- &
476!^ & tl_dTde(i,j ,k2))+ &
477!^ & cff4*(tl_cff4*dTdz(i,j,k2)+ &
478!^ & cff4*tl_dTdz(i,j,k2)- &
479!^ & tl_dTde(i,j+1,k1)))
480!^
481 adfac=cff*ad_fs(i,j,k2)
482 ad_cff1=ad_cff1+ &
483 & (2.0_r8*cff1*dtdz(i,j,k2)-dtde(i,j ,k1))* &
484 & adfac
485 ad_cff2=ad_cff2+ &
486 & (2.0_r8*cff2*dtdz(i,j,k2)-dtde(i,j+1,k2))* &
487 & adfac
488 ad_cff3=ad_cff3+ &
489 & (2.0_r8*cff3*dtdz(i,j,k2)-dtde(i,j ,k2))* &
490 & adfac
491 ad_cff4=ad_cff4+ &
492 & (2.0_r8*cff4*dtdz(i,j,k2)-dtde(i,j+1,k1))* &
493 & adfac
494 ad_dtdz(i,j,k2)=ad_dtdz(i,j,k2)+ &
495 & (cff1*cff1+ &
496 & cff2*cff2+ &
497 & cff3*cff3+ &
498 & cff4*cff4)*adfac
499 ad_dtde(i,j ,k1)=ad_dtde(i,j ,k1)-cff1*adfac
500 ad_dtde(i,j+1,k2)=ad_dtde(i,j+1,k2)-cff2*adfac
501 ad_dtde(i,j ,k2)=ad_dtde(i,j ,k2)-cff3*adfac
502 ad_dtde(i,j+1,k1)=ad_dtde(i,j+1,k1)-cff4*adfac
503!^ tl_cff4=(0.5_r8+SIGN(0.5_r8, dZde(i,j+1,k1)))* &
504!^ & tl_dZde(i,j+1,k1)
505!^
506 ad_dzde(i,j+1,k1)=ad_dzde(i,j+1,k1)+ &
507 & (0.5_r8+sign(0.5_r8, &
508 & dzde(i,j+1,k1)))* &
509 & ad_cff4
510 ad_cff4=0.0_r8
511!^ tl_cff3=(0.5_r8+SIGN(0.5_r8, dZde(i,j ,k2)))* &
512!^ & tl_dZde(i,j ,k2)
513!^
514 ad_dzde(i,j ,k2)=ad_dzde(i,j ,k2)+ &
515 & (0.5_r8+sign(0.5_r8, &
516 & dzde(i,j ,k2)))* &
517 & ad_cff3
518 ad_cff3=0.0_r8
519!^ tl_cff2=(0.5_r8+SIGN(0.5_r8,-dZde(i,j+1,k2)))* &
520!^ & tl_dZde(i,j+1,k2)
521!^
522 ad_dzde(i,j+1,k2)=ad_dzde(i,j+1,k2)+ &
523 & (0.5_r8+sign(0.5_r8, &
524 & -dzde(i,j+1,k2)))* &
525 & ad_cff2
526 ad_cff2=0.0_r8
527!^ tl_cff1=(0.5_r8+SIGN(0.5_r8,-dZde(i,j ,k1)))* &
528!^ & tl_dZde(i,j ,k1)
529!^
530 ad_dzde(i,j ,k1)=ad_dzde(i,j ,k1)+ &
531 & (0.5_r8+sign(0.5_r8, &
532 & -dzde(i,j ,k1)))* &
533 & ad_cff1
534 ad_cff1=0.0_r8
535!
536 cff1=min(dzdx(i ,j,k1),0.0_r8)
537 cff2=min(dzdx(i+1,j,k2),0.0_r8)
538 cff3=max(dzdx(i ,j,k2),0.0_r8)
539 cff4=max(dzdx(i+1,j,k1),0.0_r8)
540!^ tl_FS(i,j,k2)=cff* &
541!^ & (tl_cff1*(cff1*dTdz(i,j,k2)- &
542!^ & dTdx(i ,j,k1))+ &
543!^ & tl_cff2*(cff2*dTdz(i,j,k2)- &
544!^ & dTdx(i+1,j,k2))+ &
545!^ & tl_cff3*(cff3*dTdz(i,j,k2)- &
546!^ & dTdx(i ,j,k2))+ &
547!^ & tl_cff4*(cff4*dTdz(i,j,k2)- &
548!^ & dTdx(i+1,j,k1))+ &
549!^ & cff1*(tl_cff1*dTdz(i,j,k2)+ &
550!^ & cff1*tl_dTdz(i,j,k2)- &
551!^ & tl_dTdx(i ,j,k1))+ &
552!^ & cff2*(tl_cff2*dTdz(i,j,k2)+ &
553!^ & cff2*tl_dTdz(i,j,k2)- &
554!^ & tl_dTdx(i+1,j,k2))+ &
555!^ & cff3*(tl_cff3*dTdz(i,j,k2)+ &
556!^ & cff3*tl_dTdz(i,j,k2)- &
557!^ & tl_dTdx(i ,j,k2))+ &
558!^ & cff4*(tl_cff4*dTdz(i,j,k2)+ &
559!^ & cff4*tl_dTdz(i,j,k2)- &
560!^ & tl_dTdx(i+1,j,k1)))
561!^
562 ad_cff1=ad_cff1+ &
563 & (2.0_r8*cff1*dtdz(i,j,k2)-dtdx(i ,j,k1))* &
564 & adfac
565 ad_cff2=ad_cff2+ &
566 & (2.0_r8*cff2*dtdz(i,j,k2)-dtdx(i+1,j,k2))* &
567 & adfac
568 ad_cff3=ad_cff3+ &
569 & (2.0_r8*cff3*dtdz(i,j,k2)-dtdx(i ,j,k2))* &
570 & adfac
571 ad_cff4=ad_cff4+ &
572 & (2.0_r8*cff4*dtdz(i,j,k2)-dtdx(i+1,j,k1))* &
573 & adfac
574 ad_dtdz(i,j,k2)=ad_dtdz(i,j,k2)+ &
575 & (cff1*cff1+ &
576 & cff2*cff2+ &
577 & cff3*cff3+ &
578 & cff4*cff4)*adfac
579 ad_dtdx(i ,j,k1)=ad_dtdx(i ,j,k1)-cff1*adfac
580 ad_dtdx(i+1,j,k2)=ad_dtdx(i+1,j,k2)-cff2*adfac
581 ad_dtdx(i ,j,k2)=ad_dtdx(i ,j,k2)-cff3*adfac
582 ad_dtdx(i+1,j,k1)=ad_dtdx(i+1,j,k1)-cff4*adfac
583 ad_fs(i,j,k2)=0.0_r8
584!^ tl_cff4=(0.5_r8+SIGN(0.5_r8, dZdx(i+1,j,k1)))* &
585!^ & tl_dZdx(i+1,j,k1)
586!^
587 ad_dzdx(i+1,j,k1)=ad_dzdx(i+1,j,k1)+ &
588 & (0.5_r8+sign(0.5_r8, &
589 & dzdx(i+1,j,k1)))* &
590 & ad_cff4
591 ad_cff4=0.0_r8
592!^ tl_cff3=(0.5_r8+SIGN(0.5_r8, dZdx(i ,j,k2)))* &
593!^ & tl_dZdx(i ,j,k2)
594!^
595 ad_dzdx(i ,j,k2)=ad_dzdx(i ,j,k2)+ &
596 & (0.5_r8+sign(0.5_r8, &
597 & dzdx(i ,j,k2)))* &
598 & ad_cff3
599 ad_cff3=0.0_r8
600!^ tl_cff2=(0.5_r8+SIGN(0.5_r8,-dZdx(i+1,j,k2)))* &
601!^ & tl_dZdx(i+1,j,k2)
602!^
603 ad_dzdx(i+1,j,k2)=ad_dzdx(i+1,j,k2)+ &
604 & (0.5_r8+sign(0.5_r8, &
605 & -dzdx(i+1,j,k2)))* &
606 & ad_cff2
607 ad_cff2=0.0_r8
608!^ tl_cff1=(0.5_r8+SIGN(0.5_r8,-dZdx(i ,j,k1)))* &
609!^ & tl_dZdx(i ,j,k1)
610!^
611 ad_dzdx(i ,j,k1)=ad_dzdx(i ,j,k1)+ &
612 & (0.5_r8+sign(0.5_r8, &
613 & -dzdx(i ,j,k1)))* &
614 & ad_cff1
615 ad_cff1=0.0_r8
616 END DO
617 END DO
618 END IF
619 DO j=jstr,jend+1
620 DO i=istr,iend
621#ifdef DIFF_3DCOEF
622 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i,j-1,k))* &
623 & om_v(i,j)
624#else
625 cff=0.25_r8*(diff2(i,j,itrc)+diff2(i,j-1,itrc))* &
626 & om_v(i,j)
627#endif
628!^ tl_FE(i,j)=cff* &
629!^ & (((tl_Hz(i,j,k)+tl_Hz(i,j-1,k))* &
630!^ & (dTde(i,j,k1)- &
631!^ & 0.5_r8*(MIN(dZde(i,j,k1),0.0_r8)* &
632!^ & (dTdz(i,j-1,k1)+ &
633!^ & dTdz(i,j ,k2))+ &
634!^ & MAX(dZde(i,j,k1),0.0_r8)* &
635!^ & (dTdz(i,j-1,k2)+ &
636!^ & dTdz(i,j ,k1)))))+ &
637!^ & ((Hz(i,j,k)+Hz(i,j-1,k))* &
638!^ & (tl_dTde(i,j,k1)- &
639!^ & 0.5_r8*(MIN(dZde(i,j,k1),0.0_r8)* &
640!^ & (tl_dTdz(i,j-1,k1)+ &
641!^ & tl_dTdz(i,j ,k2))+ &
642!^ & MAX(dZde(i,j,k1),0.0_r8)* &
643!^ & (tl_dTdz(i,j-1,k2)+ &
644!^ & tl_dTdz(i,j ,k1)))- &
645!^ & 0.5_r8*((0.5_r8+ &
646!^ & SIGN(0.5_r8,-dZde(i,j,k1)))* &
647!^ & tl_dZde(i,j,k1)* &
648!^ & (dTdz(i,j-1,k1)+dTdz(i,j,k2))+ &
649!^ & (0.5_r8+ &
650!^ & SIGN(0.5_r8, dZde(i,j,k1)))* &
651!^ & tl_dZde(i,j,k1)* &
652!^ & (dTdz(i,j-1,k2)+dTdz(i,j,k1))))))
653!^
654 adfac=cff*ad_fe(i,j)
655 adfac1=adfac*(dtde(i,j,k1)- &
656 & 0.5_r8*(min(dzde(i,j,k1),0.0_r8)* &
657 & (dtdz(i,j-1,k1)+ &
658 & dtdz(i,j ,k2))+ &
659 & max(dzde(i,j,k1),0.0_r8)* &
660 & (dtdz(i,j-1,k2)+ &
661 & dtdz(i,j ,k1))))
662 adfac2=adfac*(hz(i,j,k)+hz(i,j-1,k))
663 adfac3=adfac2*0.5_r8*min(dzde(i,j,k1),0.0_r8)
664 adfac4=adfac2*0.5_r8*max(dzde(i,j,k1),0.0_r8)
665 ad_hz(i,j-1,k)=ad_hz(i,j-1,k)+adfac1
666 ad_hz(i,j ,k)=ad_hz(i,j ,k)+adfac1
667 ad_dtde(i,j,k1)=ad_dtde(i,j,k1)+adfac2
668 ad_dtdz(i,j-1,k1)=ad_dtdz(i,j-1,k1)-adfac3
669 ad_dtdz(i,j ,k2)=ad_dtdz(i,j ,k2)-adfac3
670 ad_dtdz(i,j-1,k2)=ad_dtdz(i,j-1,k2)-adfac4
671 ad_dtdz(i,j ,k1)=ad_dtdz(i,j ,k1)-adfac4
672 ad_dzde(i,j,k1)=ad_dzde(i,j,k1)- &
673 & adfac2*0.5_r8* &
674 & ((0.5_r8+sign(0.5_r8,-dzde(i,j,k1)))* &
675 & (dtdz(i,j-1,k1)+dtdz(i,j,k2))+ &
676 & (0.5_r8+sign(0.5_r8, dzde(i,j,k1)))* &
677 & (dtdz(i,j-1,k2)+dtdz(i,j,k1)))
678 ad_fe(i,j)=0.0_r8
679 END DO
680 END DO
681 DO j=jstr,jend
682 DO i=istr,iend+1
683#ifdef DIFF_3DCOEF
684 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i-1,j,k))* &
685 & on_u(i,j)
686#else
687 cff=0.25_r8*(diff2(i,j,itrc)+diff2(i-1,j,itrc))* &
688 & on_u(i,j)
689#endif
690!^ tl_FX(i,j)=cff* &
691!^ & (((tl_Hz(i,j,k)+tl_Hz(i-1,j,k))* &
692!^ & (dTdx(i,j,k1)- &
693!^ & 0.5_r8*(MIN(dZdx(i,j,k1),0.0_r8)* &
694!^ & (dTdz(i-1,j,k1)+ &
695!^ & dTdz(i ,j,k2))+ &
696!^ & MAX(dZdx(i,j,k1),0.0_r8)* &
697!^ & (dTdz(i-1,j,k2)+ &
698!^ & dTdz(i ,j,k1)))))+ &
699!^ & ((Hz(i,j,k)+Hz(i-1,j,k))* &
700!^ & (tl_dTdx(i,j,k1)- &
701!^ & 0.5_r8*(MIN(dZdx(i,j,k1),0.0_r8)* &
702!^ & (tl_dTdz(i-1,j,k1)+ &
703!^ & tl_dTdz(i ,j,k2))+ &
704!^ & MAX(dZdx(i,j,k1),0.0_r8)* &
705!^ & (tl_dTdz(i-1,j,k2)+ &
706!^ & tl_dTdz(i ,j,k1)))- &
707!^ & 0.5_r8*((0.5_r8+ &
708!^ & SIGN(0.5_r8,-dZdx(i,j,k1)))* &
709!^ & tl_dZdx(i,j,k1)* &
710!^ & (dTdz(i-1,j,k1)+dTdz(i,j,k2))+ &
711!^ & (0.5_r8+ &
712!^ & SIGN(0.5_r8, dZdx(i,j,k1)))* &
713!^ & tl_dZdx(i,j,k1)* &
714!^ & (dTdz(i-1,j,k2)+dTdz(i,j,k1))))))
715!^
716 adfac=cff*ad_fx(i,j)
717 adfac1=adfac*(dtdx(i,j,k1)- &
718 & 0.5_r8*(min(dzdx(i,j,k1),0.0_r8)* &
719 & (dtdz(i-1,j,k1)+ &
720 & dtdz(i ,j,k2))+ &
721 & max(dzdx(i,j,k1),0.0_r8)* &
722 & (dtdz(i-1,j,k2)+ &
723 & dtdz(i ,j,k1))))
724 adfac2=adfac*(hz(i,j,k)+hz(i-1,j,k))
725 adfac3=adfac2*0.5_r8*min(dzdx(i,j,k1),0.0_r8)
726 adfac4=adfac2*0.5_r8*max(dzdx(i,j,k1),0.0_r8)
727 ad_hz(i-1,j,k)=ad_hz(i-1,j,k)+adfac1
728 ad_hz(i ,j,k)=ad_hz(i ,j,k)+adfac1
729 ad_dtdx(i,j,k1)=ad_dtdx(i,j,k1)+adfac2
730 ad_dtdz(i-1,j,k1)=ad_dtdz(i-1,j,k1)-adfac3
731 ad_dtdz(i ,j,k2)=ad_dtdz(i ,j,k2)-adfac3
732 ad_dtdz(i-1,j,k2)=ad_dtdz(i-1,j,k2)-adfac4
733 ad_dtdz(i ,j,k1)=ad_dtdz(i ,j,k1)-adfac4
734 ad_dzdx(i,j,k1)=ad_dzdx(i,j,k1)- &
735 & adfac2*0.5_r8* &
736 & ((0.5_r8+sign(0.5_r8,-dzdx(i,j,k1)))* &
737 & (dtdz(i-1,j,k1)+dtdz(i,j,k2))+ &
738 & (0.5_r8+sign(0.5_r8, dzdx(i,j,k1)))* &
739 & (dtdz(i-1,j,k2)+dtdz(i,j,k1)))
740 ad_fx(i,j)=0.0_r8
741 END DO
742 END DO
743 END IF
744 IF ((k.eq.0).or.(k.eq.n(ng))) THEN
745 DO j=jstr-1,jend+1
746 DO i=istr-1,iend+1
747!^ tl_FS(i,j,k2)=0.0_r8
748!^
749 ad_fs(i,j,k2)=0.0_r8
750!^ tl_dTdz(i,j,k2)=0.0_r8
751!^
752 ad_dtdz(i,j,k2)=0.0_r8
753 END DO
754 END DO
755 ELSE
756 DO j=jstr-1,jend+1
757 DO i=istr-1,iend+1
758 cff=1.0_r8/(z_r(i,j,k+1)-z_r(i,j,k))
759#if defined TS_MIX_STABILITY
760!^ tl_dTdz(i,j,k2)=tl_cff*(0.75_r8*(t(i,j,k+1,nrhs,itrc)- &
761!^ & t(i,j,k ,nrhs,itrc))+ &
762!^ & 0.25_r8*(t(i,j,k+1,nstp,itrc)- &
763!^ & t(i,j,k ,nstp,itrc)))+&
764!^ & cff*(0.75_r8*(tl_t(i,j,k+1,nrhs,itrc)- &
765!^ & tl_t(i,j,k ,nrhs,itrc))+ &
766!^ & 0.25_r8*(tl_t(i,j,k+1,nstp,itrc)- &
767!^ & tl_t(i,j,k ,nstp,itrc)))
768!^
769 adfac=cff*ad_dtdz(i,j,k2)
770 adfac1=adfac*0.75_r8
771 adfac2=adfac*0.25_r8
772 ad_t(i,j,k ,nrhs,itrc)=ad_t(i,j,k ,nrhs,itrc)-adfac1
773 ad_t(i,j,k+1,nrhs,itrc)=ad_t(i,j,k+1,nrhs,itrc)+adfac1
774 ad_t(i,j,k ,nstp,itrc)=ad_t(i,j,k ,nstp,itrc)-adfac2
775 ad_t(i,j,k+1,nstp,itrc)=ad_t(i,j,k+1,nstp,itrc)+adfac2
776 ad_cff=ad_cff+(0.75_r8*(t(i,j,k+1,nrhs,itrc)- &
777 & t(i,j,k ,nrhs,itrc))+ &
778 & 0.25_r8*(t(i,j,k+1,nstp,itrc)- &
779 & t(i,j,k ,nstp,itrc)))* &
780 & ad_dtdz(i,j,k2)
781 ad_dtdz(i,j,k2)=0.0_r8
782#elif defined TS_MIX_CLIMA
783 IF (ltracerclm(itrc,ng)) THEN
784!^ tl_dTdz(i,j,k2)=tl_cff*((t(i,j,k+1,nrhs,itrc)- &
785!^ & tclm(i,j,k+1,itrc))- &
786!^ & (t(i,j,k ,nrhs,itrc)- &
787!^ & tclm(i,j,k ,itrc)))+ &
788!^ & cff*(tl_t(i,j,k+1,nrhs,itrc)- &
789!^ & tl_t(i,j,k ,nrhs,itrc))
790!^
791 adfac=cff*ad_dtdz(i,j,k2)
792 ad_t(i,j,k ,nrhs,itrc)=ad_t(i,j,k ,nrhs,itrc)-adfac
793 ad_t(i,j,k+1,nrhs,itrc)=ad_t(i,j,k+1,nrhs,itrc)+adfac
794 ad_cff=ad_cff+((t(i,j,k+1,nrhs,itrc)- &
795 & tclm(i,j,k+1,itrc))- &
796 & (t(i,j,k ,nrhs,itrc)- &
797 & tclm(i,j,k ,itrc)))*ad_dtdz(i,j,k2)
798 ad_dtdz(i,j,k2)=0.0_r8
799 ELSE
800!^ tl_dTdz(i,j,k2)=tl_cff*(t(i,j,k+1,nrhs,itrc)- &
801!^ & t(i,j,k ,nrhs,itrc))+ &
802!^ & cff*(tl_t(i,j,k+1,nrhs,itrc)- &
803!^ & tl_t(i,j,k ,nrhs,itrc))
804!^
805 adfac=cff*ad_dtdz(i,j,k2)
806 ad_t(i,j,k ,nrhs,itrc)=ad_t(i,j,k ,nrhs,itrc)-adfac
807 ad_t(i,j,k+1,nrhs,itrc)=ad_t(i,j,k+1,nrhs,itrc)+adfac
808 ad_cff=ad_cff+(t(i,j,k+1,nrhs,itrc)- &
809 & t(i,j,k ,nrhs,itrc))*ad_dtdz(i,j,k2)
810 ad_dtdz(i,j,k2)=0.0_r8
811 END IF
812#else
813!^ tl_dTdz(i,j,k2)=tl_cff*(t(i,j,k+1,nrhs,itrc)- &
814!^ & t(i,j,k ,nrhs,itrc))+ &
815!^ & cff*(tl_t(i,j,k+1,nrhs,itrc)- &
816!^ & tl_t(i,j,k ,nrhs,itrc))
817!^
818 adfac=cff*ad_dtdz(i,j,k2)
819 ad_t(i,j,k ,nrhs,itrc)=ad_t(i,j,k ,nrhs,itrc)-adfac
820 ad_t(i,j,k+1,nrhs,itrc)=ad_t(i,j,k+1,nrhs,itrc)+adfac
821 ad_cff=ad_cff+(t(i,j,k+1,nrhs,itrc)- &
822 & t(i,j,k ,nrhs,itrc))*ad_dtdz(i,j,k2)
823 ad_dtdz(i,j,k2)=0.0_r8
824#endif
825!^ tl_cff=-cff*cff*(tl_z_r(i,j,k+1)-tl_z_r(i,j,k))
826!^
827 adfac=cff*cff*ad_cff
828 ad_z_r(i,j,k )=ad_z_r(i,j,k )+adfac
829 ad_z_r(i,j,k+1)=ad_z_r(i,j,k+1)-adfac
830 ad_cff=0.0_r8
831 END DO
832 END DO
833 END IF
834 IF (k.lt.n(ng)) THEN
835 DO j=jstr,jend+1
836 DO i=istr,iend
837 cff=0.5_r8*(pn(i,j)+pn(i,j-1))
838#ifdef MASKING
839 cff=cff*vmask(i,j)
840#endif
841#ifdef WET_DRY_NOT_YET
842 cff=cff*vmask_wet(i,j)
843#endif
844#if defined TS_MIX_STABILITY
845!^ tl_dTde(i,j,k2)=cff* &
846!^ & (0.75_r8*(tl_t(i,j ,k+1,nrhs,itrc)- &
847!^ & tl_t(i,j-1,k+1,nrhs,itrc))+ &
848!^ & 0.25_r8*(tl_t(i,j ,k+1,nstp,itrc)- &
849!^ & tl_t(i,j-1,k+1,nstp,itrc)))
850!^
851 adfac=cff*ad_dtde(i,j,k2)
852 adfac1=adfac*0.75_r8
853 adfac2=adfac*0.25_r8
854 ad_t(i,j-1,k+1,nrhs,itrc)=ad_t(i,j-1,k+1,nrhs,itrc)- &
855 & adfac1
856 ad_t(i,j ,k+1,nrhs,itrc)=ad_t(i,j ,k+1,nrhs,itrc)+ &
857 & adfac1
858 ad_t(i,j-1,k+1,nstp,itrc)=ad_t(i,j-1,k+1,nstp,itrc)- &
859 & adfac2
860 ad_t(i,j ,k+1,nstp,itrc)=ad_t(i,j ,k+1,nstp,itrc)+ &
861 & adfac2
862 ad_dtde(i,j,k2)=0.0_r8
863#elif defined TS_MIX_CLIMA
864!^ tl_dTde(i,j,k2)=cff*(tl_t(i,j ,k+1,nrhs,itrc)- &
865!^ & tl_t(i,j-1,k+1,nrhs,itrc))
866!^
867 adfac=cff*ad_dtde(i,j,k2)
868 ad_t(i,j-1,k+1,nrhs,itrc)=ad_t(i,j-1,k+1,nrhs,itrc)- &
869 & adfac
870 ad_t(i,j ,k+1,nrhs,itrc)=ad_t(i,j ,k+1,nrhs,itrc)+ &
871 & adfac
872 ad_dtde(i,j,k2)=0.0_r8
873#else
874!^ tl_dTde(i,j,k2)=cff*(tl_t(i,j ,k+1,nrhs,itrc)- &
875!^ & tl_t(i,j-1,k+1,nrhs,itrc))
876!^
877 adfac=cff*ad_dtde(i,j,k2)
878 ad_t(i,j-1,k+1,nrhs,itrc)=ad_t(i,j-1,k+1,nrhs,itrc)- &
879 & adfac
880 ad_t(i,j ,k+1,nrhs,itrc)=ad_t(i,j ,k+1,nrhs,itrc)+ &
881 & adfac
882 ad_dtde(i,j,k2)=0.0_r8
883#endif
884!^ tl_dZde(i,j,k2)=cff*(tl_z_r(i,j ,k+1)- &
885!^ & tl_z_r(i,j-1,k+1))
886!^
887 adfac=cff*ad_dzde(i,j,k2)
888 ad_z_r(i,j-1,k+1)=ad_z_r(i,j-1,k+1)-adfac
889 ad_z_r(i,j ,k+1)=ad_z_r(i,j ,k+1)+adfac
890 ad_dzde(i,j,k2)=0.0_r8
891 END DO
892 END DO
893 DO j=jstr,jend
894 DO i=istr,iend+1
895 cff=0.5_r8*(pm(i,j)+pm(i-1,j))
896#ifdef MASKING
897 cff=cff*umask(i,j)
898#endif
899#ifdef WET_DRY_NOT_YET
900 cff=cff*umask_wet(i,j)
901#endif
902#if defined TS_MIX_STABILITY
903!^ tl_dTdx(i,j,k2)=cff* &
904!^ & (0.75_r8*(tl_t(i ,j,k+1,nrhs,itrc)- &
905!^ & tl_t(i-1,j,k+1,nrhs,itrc))+ &
906!^ & 0.25_r8*(tl_t(i ,j,k+1,nstp,itrc)- &
907!^ & tl_t(i-1,j,k+1,nstp,itrc)))
908!^
909 adfac=cff*ad_dtdx(i,j,k2)
910 adfac1=adfac*0.75_r8
911 adfac2=adfac*0.25_r8
912 ad_t(i-1,j,k+1,nrhs,itrc)=ad_t(i-1,j,k+1,nrhs,itrc)- &
913 & adfac1
914 ad_t(i ,j,k+1,nrhs,itrc)=ad_t(i ,j,k+1,nrhs,itrc)+ &
915 & adfac1
916 ad_t(i-1,j,k+1,nstp,itrc)=ad_t(i-1,j,k+1,nstp,itrc)- &
917 & adfac1
918 ad_t(i ,j,k+1,nstp,itrc)=ad_t(i ,j,k+1,nstp,itrc)+ &
919 & adfac1
920 ad_dtdx(i,j,k2)=0.0_r8
921#elif defined TS_MIX_CLIMA
922!^ tl_dTdx(i,j,k2)=cff*(tl_t(i ,j,k+1,nrhs,itrc)- &
923!^ & tl_t(i-1,j,k+1,nrhs,itrc))
924!^
925 adfac=cff*ad_dtdx(i,j,k2)
926 ad_t(i-1,j,k+1,nrhs,itrc)=ad_t(i-1,j,k+1,nrhs,itrc)- &
927 & adfac
928 ad_t(i ,j,k+1,nrhs,itrc)=ad_t(i ,j,k+1,nrhs,itrc)+ &
929 & adfac
930 ad_dtdx(i,j,k2)=0.0_r8
931#else
932!^ tl_dTdx(i,j,k2)=cff*(tl_t(i ,j,k+1,nrhs,itrc)- &
933!^ & tl_t(i-1,j,k+1,nrhs,itrc))
934!^
935 adfac=cff*ad_dtdx(i,j,k2)
936 ad_t(i-1,j,k+1,nrhs,itrc)=ad_t(i-1,j,k+1,nrhs,itrc)- &
937 & adfac
938 ad_t(i ,j,k+1,nrhs,itrc)=ad_t(i ,j,k+1,nrhs,itrc)+ &
939 & adfac
940 ad_dtdx(i,j,k2)=0.0_r8
941#endif
942!^ tl_dZdx(i,j,k2)=cff*(tl_z_r(i ,j,k+1)- &
943!^ & tl_z_r(i-1,j,k+1))
944!^
945 adfac=cff*ad_dzdx(i,j,k2)
946 ad_z_r(i-1,j,k+1)=ad_z_r(i-1,j,k+1)-adfac
947 ad_z_r(i ,j,k+1)=ad_z_r(i ,j,k+1)+adfac
948 ad_dzdx(i,j,k2)=0.0_r8
949 END DO
950 END DO
951 END IF
952!
953! Compute new storage recursive indices.
954!
955 kt=k2
956 k2=k1
957 k1=kt
958 END DO k_loop
959 END DO t_loop
960!
961 RETURN
integer, dimension(:), allocatable n
Definition mod_param.F:479
integer, dimension(:), allocatable nt
Definition mod_param.F:489
real(dp), dimension(:), allocatable dt
logical, dimension(:,:), allocatable ltracerclm

References mod_scalars::dt, and mod_scalars::ltracerclm.

Referenced by ad_t3dmix2().

Here is the caller graph for this function:

◆ ad_t3dmix2_iso_tile()

subroutine ad_t3dmix2_mod::ad_t3dmix2_iso_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,
integer, intent(in) nstp,
integer, intent(in) nnew,
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), intent(in) pm,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) pn,
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,nt(ng)), intent(in) diff2,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) pden,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(inout) ad_pden,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),nt(ng)), intent(in) tclm,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),3,nt(ng)), intent(in) t,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),3,nt(ng)), intent(inout) ad_t )
private

Definition at line 89 of file ad_t3dmix2_iso.h.

108!***********************************************************************
109!
110 USE mod_param
111 USE mod_scalars
112!
113! Imported variable declarations.
114!
115 integer, intent(in) :: ng, tile
116 integer, intent(in) :: LBi, UBi, LBj, UBj
117 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
118 integer, intent(in) :: nrhs, nstp, nnew
119
120#ifdef ASSUMED_SHAPE
121# ifdef MASKING
122 real(r8), intent(in) :: umask(LBi:,LBj:)
123 real(r8), intent(in) :: vmask(LBi:,LBj:)
124# endif
125 real(r8), intent(in) :: diff2(LBi:,LBj:,:)
126 real(r8), intent(in) :: om_v(LBi:,LBj:)
127 real(r8), intent(in) :: on_u(LBi:,LBj:)
128 real(r8), intent(in) :: pm(LBi:,LBj:)
129 real(r8), intent(in) :: pn(LBi:,LBj:)
130 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
131 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
132 real(r8), intent(in) :: pden(LBi:,LBj:,:)
133 real(r8), intent(in) :: t(LBi:,LBj:,:,:,:)
134# ifdef TS_MIX_CLIMA
135 real(r8), intent(in) :: tclm(LBi:,LBj:,:,:)
136# endif
137# ifdef DIAGNOSTICS_TS
138 real(r8), intent(inout) :: DiaTwrk(LBi:,LBj:,:,:,:)
139# endif
140 real(r8), intent(inout) :: ad_Hz(LBi:,LBj:,:)
141 real(r8), intent(inout) :: ad_z_r(LBi:,LBj:,:)
142 real(r8), intent(inout) :: ad_pden(LBi:,LBj:,:)
143 real(r8), intent(inout) :: ad_t(LBi:,LBj:,:,:,:)
144#else
145# ifdef MASKING
146 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
147 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
148# endif
149 real(r8), intent(in) :: diff2(LBi:UBi,LBj:UBj,NT(ng))
150 real(r8), intent(in) :: om_v(LBi:UBi,LBj:UBj)
151 real(r8), intent(in) :: on_u(LBi:UBi,LBj:UBj)
152 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
153 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
154 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
155 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
156 real(r8), intent(in) :: pden(LBi:UBi,LBj:UBj,N(ng))
157 real(r8), intent(in) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
158# ifdef TS_MIX_CLIMA
159 real(r8), intent(in) :: tclm(LBi:UBi,LBj:UBj,N(ng),NT(ng))
160# endif
161# ifdef DIAGNOSTICS_TS
162!! real(r8), intent(inout) :: DiaTwrk(LBi:UBi,LBj:UBj,N(ng),NT(ng), &
163!! & NDT)
164# endif
165 real(r8), intent(inout) :: ad_Hz(LBi:UBi,LBj:UBj,N(ng))
166 real(r8), intent(inout) :: ad_z_r(LBi:UBi,LBj:UBj,N(ng))
167 real(r8), intent(inout) :: ad_pden(LBi:UBi,LBj:UBj,N(ng))
168 real(r8), intent(inout) :: ad_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
169#endif
170!
171! Local variable declarations.
172!
173 integer :: i, itrc, j, k, kk, kt, k1, k1b, k2, k2b
174
175 real(r8), parameter :: eps = 0.5_r8
176 real(r8), parameter :: small = 1.0e-14_r8
177 real(r8), parameter :: slope_max = 0.0001_r8
178 real(r8), parameter :: strat_min = 0.1_r8
179
180 real(r8) :: cff, cff1, cff2, cff3, cff4
181 real(r8) :: ad_cff, ad_cff1, ad_cff2, ad_cff3, ad_cff4
182 real(r8) :: adfac, adfac1, adfac2, adfac3, adfac4
183
184 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_FE
185 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_FX
186
187 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: FS
188 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dRde
189 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dRdx
190 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dTde
191 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dTdr
192 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dTdx
193
194 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: ad_FS
195 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: ad_dRde
196 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: ad_dRdx
197 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: ad_dTde
198 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: ad_dTdr
199 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: ad_dTdx
200
201#include "set_bounds.h"
202!
203!-----------------------------------------------------------------------
204! Initialize adjoint private variables.
205!-----------------------------------------------------------------------
206!
207 ad_cff=0.0_r8
208 ad_cff1=0.0_r8
209 ad_cff2=0.0_r8
210 ad_cff3=0.0_r8
211 ad_cff4=0.0_r8
212
213 ad_fe(imins:imaxs,jmins:jmaxs)=0.0_r8
214 ad_fx(imins:imaxs,jmins:jmaxs)=0.0_r8
215
216 ad_fs(imins:imaxs,jmins:jmaxs,1:2)=0.0_r8
217
218 ad_drde(imins:imaxs,jmins:jmaxs,1:2)=0.0_r8
219 ad_drdx(imins:imaxs,jmins:jmaxs,1:2)=0.0_r8
220 ad_dtde(imins:imaxs,jmins:jmaxs,1:2)=0.0_r8
221 ad_dtdr(imins:imaxs,jmins:jmaxs,1:2)=0.0_r8
222 ad_dtdx(imins:imaxs,jmins:jmaxs,1:2)=0.0_r8
223!
224!----------------------------------------------------------------------
225! Compute horizontal harmonic diffusion along isopycnic surfaces.
226!----------------------------------------------------------------------
227!
228! Compute horizontal and density gradients. Notice the recursive
229! blocking sequence. The vertical placement of the gradients is:
230!
231! dTdx,dTde(:,:,k1) k rho-points
232! dTdx,dTde(:,:,k2) k+1 rho-points
233! FS,dTdr(:,:,k1) k-1/2 W-points
234! FS,dTdr(:,:,k2) k+1/2 W-points
235!
236! Compute adjoint of starting values of k1 and k2.
237!
238 t_loop : DO itrc=1,nt(ng)
239 k1=2
240 k2=1
241 DO k=0,n(ng)
242!!
243!! Note: The following code is equivalent to
244!!
245!! kt=k1
246!! k1=k2
247!! k2=kt
248!!
249!! We use the adjoint of above code.
250!!
251 k1=k2
252 k2=3-k1
253 END DO
254 k_loop : DO k=n(ng),0,-1
255!
256! Compute required BASIC STATE fields. Need to look forward in
257! recursive kk index.
258!
259 k2b=1
260 DO kk=0,k
261 k1b=k2b
262 k2b=3-k1b
263!
264! Compute components of the rotated tracer flux (T m3/s) along
265! isopycnic surfaces (required BASIC STATE fields).
266!
267 IF (kk.lt.n(ng)) THEN
268 DO j=jstr,jend
269 DO i=istr,iend+1
270 cff=0.5_r8*(pm(i,j)+pm(i-1,j))
271#ifdef MASKING
272 cff=cff*umask(i,j)
273#endif
274 drdx(i,j,k2b)=cff*(pden(i ,j,kk+1)- &
275 & pden(i-1,j,kk+1))
276#ifdef TS_MIX_CLIMA
277 dtdx(i,j,k2b)=cff*((t(i ,j,k+1,nrhs,itrc)- &
278 & tclm(i ,j,k+1,itrc))- &
279 & (t(i-1,j,k+1,nrhs,itrc)- &
280 & tclm(i-1,j,k+1,itrc)))
281#else
282 dtdx(i,j,k2b)=cff*(t(i ,j,kk+1,nrhs,itrc)- &
283 & t(i-1,j,kk+1,nrhs,itrc))
284#endif
285 END DO
286 END DO
287 IF (kk.eq.0) THEN
288 DO j=jstr,jend
289 DO i=istr,iend+1
290 drdx(i,j,k1b)=0.0_r8
291 dtdx(i,j,k1b)=0.0_r8
292 END DO
293 END DO
294 END IF
295 DO j=jstr,jend+1
296 DO i=istr,iend
297 cff=0.5_r8*(pn(i,j)+pn(i,j-1))
298#ifdef MASKING
299 cff=cff*vmask(i,j)
300#endif
301 drde(i,j,k2b)=cff*(pden(i,j ,kk+1)- &
302 & pden(i,j-1,kk+1))
303#ifdef TS_MIX_CLIMA
304 dtde(i,j,k2b)=cff*((t(i,j ,k+1,nrhs,itrc)- &
305 & tclm(i,j ,k+1,itrc))- &
306 & (t(i,j-1,k+1,nrhs,itrc)- &
307 & tclm(i,j-1,k+1,itrc)))
308#else
309 dtde(i,j,k2b)=cff*(t(i,j ,kk+1,nrhs,itrc)- &
310 & t(i,j-1,kk+1,nrhs,itrc))
311#endif
312 END DO
313 END DO
314 IF (kk.eq.0) THEN
315 DO j=jstr,jend+1
316 DO i=istr,iend
317 drde(i,j,k1b)=0.0_r8
318 dtde(i,j,k1b)=0.0_r8
319 END DO
320 END DO
321 END IF
322 END IF
323 IF ((kk.eq.0).or.(kk.eq.n(ng))) THEN
324 DO j=jstr-1,jend+1
325 DO i=istr-1,iend+1
326 dtdr(i,j,k2b)=0.0_r8
327 fs(i,j,k2b)=0.0_r8
328 END DO
329 END DO
330 IF (kk.eq.0) THEN
331 DO j=jstr-1,jend+1
332 DO i=istr-1,iend+1
333 dtdr(i,j,k1b)=0.0_r8
334 fs(i,j,k1b)=0.0_r8
335 END DO
336 END DO
337 END IF
338 ELSE
339 DO j=jstr-1,jend+1
340 DO i=istr-1,iend+1
341#if defined TS_MIX_MAX_SLOPE
342 cff1=sqrt(drdx(i,j,k2b)**2+drdx(i+1,j,k2b)**2+ &
343 & drdx(i,j,k1b)**2+drdx(i+1,j,k1b)**2+ &
344 & drde(i,j,k2b)**2+drde(i,j+1,k2b)**2+ &
345 & drde(i,j,k1b)**2+drde(i,j+1,k1b)**2)
346 cff2=0.25_r8*slope_max* &
347 & (z_r(i,j,kk+1)-z_r(i,j,kk))*cff1
348 cff3=max(pden(i,j,kk)-pden(i,j,kk+1),small)
349 cff4=max(cff2,cff3)
350 cff=-1.0_r8/cff4
351#elif defined TS_MIX_MIN_STRAT
352 cff1=max(pden(i,j,kk)-pden(i,j,kk+1), &
353 & strat_min*(z_r(i,j,kk+1)-z_r(i,j,kk)))
354 cff=-1.0_r8/cff1
355#else
356 cff1=max(pden(i,j,kk)-pden(i,j,kk+1),eps)
357 cff=-1.0_r8/cff1
358#endif
359#ifdef TS_MIX_CLIMA
360 dtdr(i,j,k2b)=cff*((t(i,j,k+1,nrhs,itrc)- &
361 & tclm(i,j,k+1,itrc))- &
362 & (t(i,j,k ,nrhs,itrc)- &
363 & tclm(i,j,k ,itrc)))
364#else
365 dtdr(i,j,k2b)=cff*(t(i,j,kk+1,nrhs,itrc)- &
366 & t(i,j,kk ,nrhs,itrc))
367#endif
368 fs(i,j,k2b)=cff*(z_r(i,j,kk+1)-z_r(i,j,kk))
369 END DO
370 END DO
371 END IF
372 END DO
373!
374 IF (k.gt.0) THEN
375!
376! Time-step harmonic, isopycnic diffusion term.
377!
378 DO j=jstr,jend
379 DO i=istr,iend
380#ifdef DIAGNOSTICS_TS
381!! DiaTwrk(i,j,k,itrc,iThdif)=cff
382#endif
383!^ tl_t(i,j,k,nnew,itrc)=tl_t(i,j,k,nnew,itrc)+tl_cff
384!^
385 ad_cff=ad_cff+ad_t(i,j,k,nnew,itrc)
386!^ tl_cff=dt(ng)*pm(i,j)*pn(i,j)* &
387!^ & (tl_FX(i+1,j)-tl_FX(i,j)+ &
388!^ & tl_FE(i,j+1)-tl_FE(i,j))+ &
389!^ & dt(ng)*(tl_FS(i,j,k2)-tl_FS(i,j,k1))
390!^
391 adfac=dt(ng)*ad_cff
392 adfac1=adfac*pm(i,j)*pn(i,j)
393 ad_fs(i,j,k2)=ad_fs(i,j,k2)+adfac
394 ad_fs(i,j,k1)=ad_fs(i,j,k1)-adfac
395 ad_fe(i,j )=ad_fe(i,j )-adfac1
396 ad_fe(i,j+1)=ad_fe(i,j+1)+adfac1
397 ad_fx(i ,j)=ad_fx(i ,j)-adfac1
398 ad_fx(i+1,j)=ad_fx(i+1,j)+adfac1
399 ad_cff=0.0_r8
400 END DO
401 END DO
402!
403! Compute components of the rotated tracer flux (T m4/s) along
404! isopycnic surfaces.
405!
406 IF (k.lt.n(ng)) THEN
407 DO j=jstr,jend
408 DO i=istr,iend
409 cff1=max(drdx(i ,j,k1),0.0_r8)
410 cff2=max(drdx(i+1,j,k2),0.0_r8)
411 cff3=min(drdx(i ,j,k2),0.0_r8)
412 cff4=min(drdx(i+1,j,k1),0.0_r8)
413 cff=cff1*(cff1*dtdr(i,j,k2)-dtdx(i ,j,k1))+ &
414 & cff2*(cff2*dtdr(i,j,k2)-dtdx(i+1,j,k2))+ &
415 & cff3*(cff3*dtdr(i,j,k2)-dtdx(i ,j,k2))+ &
416 & cff4*(cff4*dtdr(i,j,k2)-dtdx(i+1,j,k1))
417 cff1=max(drde(i,j ,k1),0.0_r8)
418 cff2=max(drde(i,j+1,k2),0.0_r8)
419 cff3=min(drde(i,j ,k2),0.0_r8)
420 cff4=min(drde(i,j+1,k1),0.0_r8)
421 cff=cff+ &
422 & cff1*(cff1*dtdr(i,j,k2)-dtde(i,j ,k1))+ &
423 & cff2*(cff2*dtdr(i,j,k2)-dtde(i,j+1,k2))+ &
424 & cff3*(cff3*dtdr(i,j,k2)-dtde(i,j ,k2))+ &
425 & cff4*(cff4*dtdr(i,j,k2)-dtde(i,j+1,k1))
426!^ tl_FS(i,j,k2)=0.5_r8*diff2(i,j,itrc)* &
427!^ & (tl_cff*FS(i,j,k2)+ &
428!^ & cff*tl_FS(i,j,k2))
429!^
430 adfac=0.5_r8*diff2(i,j,itrc)*ad_fs(i,j,k2)
431 ad_cff=ad_cff+adfac*fs(i,j,k2)
432 ad_fs(i,j,k2)=cff*adfac
433!^ tl_cff=tl_cff+ &
434!^ & tl_cff1*(cff1*dTdr(i,j,k2)- &
435!^ & dTde(i,j ,k1))+ &
436!^ & tl_cff2*(cff2*dTdr(i,j,k2)- &
437!^ & dTde(i,j+1,k2))+ &
438!^ & tl_cff3*(cff3*dTdr(i,j,k2)- &
439!^ & dTde(i,j ,k2))+ &
440!^ & tl_cff4*(cff4*dTdr(i,j,k2)- &
441!^ & dTde(i,j+1,k1))+ &
442!^ & cff1*(tl_cff1*dTdr(i,j,k2)+ &
443!^ & cff1*tl_dTdr(i,j,k2)- &
444!^ & tl_dTde(i,j ,k1))+ &
445!^ & cff2*(tl_cff2*dTdr(i,j,k2)+ &
446!^ & cff2*tl_dTdr(i,j,k2)- &
447!^ & tl_dTde(i,j+1,k2))+ &
448!^ & cff3*(tl_cff3*dTdr(i,j,k2)+ &
449!^ & cff3*tl_dTdr(i,j,k2)- &
450!^ & tl_dTde(i,j ,k2))+ &
451!^ & cff4*(tl_cff4*dTdr(i,j,k2)+ &
452!^ & cff4*tl_dTdr(i,j,k2)- &
453!^ & tl_dTde(i,j+1,k1))
454!^
455 ad_cff1=ad_cff1+ &
456 & (2.0_r8*cff1*dtdr(i,j,k2)-dtde(i,j ,k1))* &
457 & ad_cff
458 ad_cff2=ad_cff2+ &
459 & (2.0_r8*cff2*dtdr(i,j,k2)-dtde(i,j+1,k2))* &
460 & ad_cff
461 ad_cff3=ad_cff3+ &
462 & (2.0_r8*cff3*dtdr(i,j,k2)-dtde(i,j ,k2))* &
463 & ad_cff
464 ad_cff4=ad_cff4+ &
465 & (2.0_r8*cff4*dtdr(i,j,k2)-dtde(i,j+1,k1))* &
466 & ad_cff
467 ad_dtdr(i,j,k2)=ad_dtdr(i,j,k2)+ &
468 & (cff1*cff1+ &
469 & cff2*cff2+ &
470 & cff3*cff3+ &
471 & cff4*cff4)*ad_cff
472 ad_dtde(i,j ,k1)=ad_dtde(i,j ,k1)-cff1*ad_cff
473 ad_dtde(i,j+1,k2)=ad_dtde(i,j+1,k2)-cff2*ad_cff
474 ad_dtde(i,j ,k2)=ad_dtde(i,j ,k2)-cff3*ad_cff
475 ad_dtde(i,j+1,k1)=ad_dtde(i,j+1,k1)-cff4*ad_cff
476!^ tl_cff4=(0.5_r8+SIGN(0.5_r8,-dRde(i,j+1,k1)))* &
477!^ & tl_dRde(i,j+1,k1)
478!^
479 ad_drde(i,j+1,k1)=ad_drde(i,j+1,k1)+ &
480 & (0.5_r8+sign(0.5_r8, &
481 & -drde(i,j+1,k1)))* &
482 & ad_cff4
483 ad_cff4=0.0_r8
484!^ tl_cff3=(0.5_r8+SIGN(0.5_r8,-dRde(i,j ,k2)))* &
485!^ & tl_dRde(i,j ,k2)
486!^
487 ad_drde(i,j ,k2)=ad_drde(i,j ,k2)+ &
488 & (0.5_r8+sign(0.5_r8, &
489 & -drde(i,j ,k2)))* &
490 & ad_cff3
491 ad_cff3=0.0_r8
492!^ tl_cff2=(0.5_r8+SIGN(0.5_r8, dRde(i,j+1,k2)))* &
493!^ & tl_dRde(i,j+1,k2)
494!^
495 ad_drde(i,j+1,k2)=ad_drde(i,j+1,k2)+ &
496 & (0.5_r8+sign(0.5_r8, &
497 & drde(i,j+1,k2)))* &
498 & ad_cff2
499 ad_cff2=0.0_r8
500!^ tl_cff1=(0.5_r8+SIGN(0.5_r8, dRde(i,j ,k1)))* &
501!^ & tl_dRde(i,j ,k1)
502!^
503 ad_drde(i ,j,k1)=ad_drde(i ,j,k1)+ &
504 & (0.5_r8+sign(0.5_r8, &
505 & drde(i ,j,k1)))* &
506 & ad_cff1
507 ad_cff1=0.0_r8
508!
509 cff1=max(drdx(i ,j,k1),0.0_r8)
510 cff2=max(drdx(i+1,j,k2),0.0_r8)
511 cff3=min(drdx(i ,j,k2),0.0_r8)
512 cff4=min(drdx(i+1,j,k1),0.0_r8)
513!^ tl_cff=tl_cff1*(cff1*dTdr(i ,j,k2)- &
514!^ & dTdx(i ,j,k1))+ &
515!^ & tl_cff2*(cff2*dTdr(i,j,k2)- &
516!^ & dTdx(i+1,j,k2))+ &
517!^ & tl_cff3*(cff3*dTdr(i,j,k2)- &
518!^ & dTdx(i ,j,k2))+ &
519!^ & tl_cff4*(cff4*dTdr(i,j,k2)- &
520!^ & dTdx(i+1,j,k1))+ &
521!^ & cff1*(tl_cff1*dTdr(i,j,k2)+ &
522!^ & cff1*tl_dTdr(i,j,k2)- &
523!^ & tl_dTdx(i ,j,k1))+ &
524!^ & cff2*(tl_cff2*dTdr(i,j,k2)+ &
525!^ & cff2*tl_dTdr(i,j,k2)- &
526!^ & tl_dTdx(i+1,j,k2))+ &
527!^ & cff3*(tl_cff3*dTdr(i,j,k2)+ &
528!^ & cff3*tl_dTdr(i,j,k2)- &
529!^ & tl_dTdx(i ,j,k2))+ &
530!^ & cff4*(tl_cff4*dTdr(i,j,k2)+ &
531!^ & cff4*tl_dTdr(i,j,k2)- &
532!^ & tl_dTdx(i+1,j,k1))
533!^
534 ad_cff1=ad_cff1+ &
535 & (2.0_r8*cff1*dtdr(i,j,k2)-dtdx(i ,j,k1))* &
536 & ad_cff
537 ad_cff2=ad_cff2+ &
538 & (2.0_r8*cff2*dtdr(i,j,k2)-dtdx(i+1,j,k2))* &
539 & ad_cff
540 ad_cff3=ad_cff3+ &
541 & (2.0_r8*cff3*dtdr(i,j,k2)-dtdx(i ,j,k2))* &
542 & ad_cff
543 ad_cff4=ad_cff4+ &
544 & (2.0_r8*cff4*dtdr(i,j,k2)-dtdx(i+1,j,k1))* &
545 & ad_cff
546 ad_dtdr(i,j,k2)=ad_dtdr(i,j,k2)+ &
547 & (cff1*cff1+ &
548 & cff2*cff2+ &
549 & cff3*cff3+ &
550 & cff4*cff4)*ad_cff
551 ad_dtdx(i ,j,k1)=ad_dtdx(i ,j,k1)-cff1*ad_cff
552 ad_dtdx(i+1,j,k2)=ad_dtdx(i+1,j,k2)-cff2*ad_cff
553 ad_dtdx(i ,j,k2)=ad_dtdx(i ,j,k2)-cff3*ad_cff
554 ad_dtdx(i+1,j,k1)=ad_dtdx(i+1,j,k1)-cff4*ad_cff
555 ad_cff=0.0_r8
556!^ tl_cff4=(0.5_r8+SIGN(0.5_r8,-dRdx(i+1,j,k1)))* &
557!^ & tl_dRdx(i+1,j,k1)
558!^
559 ad_drdx(i+1,j,k1)=ad_drdx(i+1,j,k1)+ &
560 & (0.5_r8+sign(0.5_r8, &
561 & -drdx(i+1,j,k1)))* &
562 & ad_cff4
563 ad_cff4=0.0_r8
564!^ tl_cff3=(0.5_r8+SIGN(0.5_r8,-dRdx(i ,j,k2)))* &
565!^ & tl_dRdx(i ,j,k2)
566!^
567 ad_drdx(i ,j,k2)=ad_drdx(i ,j,k2)+ &
568 & (0.5_r8+sign(0.5_r8, &
569 & -drdx(i ,j,k2)))* &
570 & ad_cff3
571 ad_cff3=0.0_r8
572!^ tl_cff2=(0.5_r8+SIGN(0.5_r8, dRdx(i+1,j,k2)))* &
573!^ & tl_dRdx(i+1,j,k2)
574!^
575 ad_drdx(i+1,j,k2)=ad_drdx(i+1,j,k2)+ &
576 & (0.5_r8+sign(0.5_r8, &
577 & drdx(i+1,j,k2)))* &
578 & ad_cff2
579 ad_cff2=0.0_r8
580!^ tl_cff1=(0.5_r8+SIGN(0.5_r8, dRdx(i ,j,k1)))* &
581!^ & tl_dRdx(i ,j,k1)
582!^
583 ad_drdx(i ,j,k1)=ad_drdx(i ,j,k1)+ &
584 & (0.5_r8+sign(0.5_r8, &
585 & drdx(i ,j,k1)))* &
586 & ad_cff1
587 ad_cff1=0.0_r8
588 END DO
589 END DO
590 END IF
591 DO j=jstr,jend+1
592 DO i=istr,iend
593 cff=0.25_r8*(diff2(i,j,itrc)+diff2(i,j-1,itrc))* &
594 & om_v(i,j)
595!^ tl_FE(i,j)=cff* &
596!^ & (((tl_Hz(i,j,k)+tl_Hz(i,j-1,k))* &
597!^ & (dTde(i,j,k1)- &
598!^ & 0.5_r8*(MAX(dRde(i,j,k1),0.0_r8)* &
599!^ & (dTdr(i,j-1,k1)+ &
600!^ & dTdr(i,j ,k2))+ &
601!^ & MIN(dRde(i,j,k1),0.0_r8)* &
602!^ & (dTdr(i,j-1,k2)+ &
603!^ & dTdr(i,j ,k1)))))+ &
604!^ & ((Hz(i,j,k)+Hz(i,j-1,k))* &
605!^ & (tl_dTde(i,j,k1)- &
606!^ & 0.5_r8*(MAX(dRde(i,j,k1),0.0_r8)* &
607!^ & (tl_dTdr(i,j-1,k1)+ &
608!^ & tl_dTdr(i,j ,k2))+ &
609!^ & MIN(dRde(i,j,k1),0.0_r8)* &
610!^ & (tl_dTdr(i,j-1,k2)+ &
611!^ & tl_dTdr(i,j ,k1)))- &
612!^ & 0.5_r8*((0.5_r8+ &
613!^ & SIGN(0.5_r8, dRde(i,j,k1)))* &
614!^ & tl_dRde(i,j,k1)* &
615!^ & (dTdr(i,j-1,k1)+dTdr(i,j,k2))+ &
616!^ & (0.5_r8+ &
617!^ & SIGN(0.5_r8,-dRde(i,j,k1)))* &
618!^ & tl_dRde(i,j,k1)* &
619!^ & (dTdr(i,j-1,k2)+dTdr(i,j,k1))))))
620!^
621 adfac=cff*ad_fe(i,j)
622 adfac1=adfac*(dtde(i,j,k1)- &
623 & 0.5_r8*(max(drde(i,j,k1),0.0_r8)* &
624 & (dtdr(i,j-1,k1)+ &
625 & dtdr(i,j ,k2))+ &
626 & min(drde(i,j,k1),0.0_r8)* &
627 & (dtdr(i,j-1,k2)+ &
628 & dtdr(i,j ,k1))))
629 adfac2=adfac*(hz(i,j,k)+hz(i,j-1,k))
630 adfac3=adfac2*0.5_r8*max(drde(i,j,k1),0.0_r8)
631 adfac4=adfac2*0.5_r8*min(drde(i,j,k1),0.0_r8)
632 ad_hz(i,j-1,k)=ad_hz(i,j-1,k)+adfac1
633 ad_hz(i,j ,k)=ad_hz(i,j ,k)+adfac1
634 ad_dtde(i,j,k1)=ad_dtde(i,j,k1)+adfac2
635 ad_dtdr(i,j-1,k1)=ad_dtdr(i,j-1,k1)-adfac3
636 ad_dtdr(i,j ,k2)=ad_dtdr(i,j ,k2)-adfac3
637 ad_dtdr(i,j-1,k2)=ad_dtdr(i,j-1,k2)-adfac4
638 ad_dtdr(i,j ,k1)=ad_dtdr(i,j ,k1)-adfac4
639 ad_drde(i,j,k1)=ad_drde(i,j,k1)- &
640 & 0.5_r8* &
641 & ((0.5_r8+sign(0.5_r8, drde(i,j,k1)))* &
642 & (dtdr(i,j-1,k1)+dtdr(i,j,k2))+ &
643 & (0.5_r8+sign(0.5_r8,-drde(i,j,k1)))* &
644 & (dtdr(i,j-1,k2)+dtdr(i,j,k1)))*adfac2
645 ad_fe(i,j)=0.0_r8
646 END DO
647 END DO
648 DO j=jstr,jend
649 DO i=istr,iend+1
650 cff=0.25_r8*(diff2(i,j,itrc)+diff2(i-1,j,itrc))* &
651 & on_u(i,j)
652!^ tl_FX(i,j)=cff* &
653!^ & (((tl_Hz(i,j,k)+tl_Hz(i-1,j,k))* &
654!^ & (dTdx(i,j,k1)- &
655!^ & 0.5_r8*(MAX(dRdx(i,j,k1),0.0_r8)* &
656!^ & (dTdr(i-1,j,k1)+ &
657!^ & dTdr(i ,j,k2))+ &
658!^ & MIN(dRdx(i,j,k1),0.0_r8)* &
659!^ & (dTdr(i-1,j,k2)+ &
660!^ & dTdr(i ,j,k1)))))+ &
661!^ & ((Hz(i,j,k)+Hz(i-1,j,k))* &
662!^ & (tl_dTdx(i,j,k1)- &
663!^ & 0.5_r8*(MAX(dRdx(i,j,k1),0.0_r8)* &
664!^ & (tl_dTdr(i-1,j,k1)+ &
665!^ & tl_dTdr(i ,j,k2))+ &
666!^ & MIN(dRdx(i,j,k1),0.0_r8)* &
667!^ & (tl_dTdr(i-1,j,k2)+ &
668!^ & tl_dTdr(i ,j,k1)))- &
669!^ & 0.5_r8*((0.5_r8+ &
670!^ & SIGN(0.5_r8, dRdx(i,j,k1)))* &
671!^ & tl_dRdx(i,j,k1)* &
672!^ & (dTdr(i-1,j,k1)+dTdr(i,j,k2))+ &
673!^ & (0.5_r8+ &
674!^ & SIGN(0.5_r8,-dRdx(i,j,k1)))* &
675!^ & tl_dRdx(i,j,k1)* &
676!^ & (dTdr(i-1,j,k2)+dTdr(i,j,k1))))))
677!^
678 adfac=cff*ad_fx(i,j)
679 adfac1=adfac*(dtdx(i,j,k1)- &
680 & 0.5_r8*(max(drdx(i,j,k1),0.0_r8)* &
681 & (dtdr(i-1,j,k1)+ &
682 & dtdr(i ,j,k2))+ &
683 & min(drdx(i,j,k1),0.0_r8)* &
684 & (dtdr(i-1,j,k2)+ &
685 & dtdr(i ,j,k1))))
686 adfac2=adfac*(hz(i,j,k)+hz(i-1,j,k))
687 adfac3=adfac2*0.5_r8*max(drdx(i,j,k1),0.0_r8)
688 adfac4=adfac2*0.5_r8*min(drdx(i,j,k1),0.0_r8)
689 ad_hz(i-1,j,k)=ad_hz(i-1,j,k)+adfac1
690 ad_hz(i ,j,k)=ad_hz(i ,j,k)+adfac1
691 ad_dtdx(i,j,k1)=ad_dtdx(i,j,k1)+adfac2
692 ad_dtdr(i-1,j,k1)=ad_dtdr(i-1,j,k1)-adfac3
693 ad_dtdr(i ,j,k2)=ad_dtdr(i ,j,k2)-adfac3
694 ad_dtdr(i-1,j,k2)=ad_dtdr(i-1,j,k2)-adfac4
695 ad_dtdr(i ,j,k1)=ad_dtdr(i ,j,k1)-adfac4
696 ad_drdx(i,j,k1)=ad_drdx(i,j,k1)- &
697 & 0.5_r8* &
698 & ((0.5_r8+sign(0.5_r8, drdx(i,j,k1)))* &
699 & (dtdr(i-1,j,k1)+dtdr(i,j,k2))+ &
700 & (0.5_r8+sign(0.5_r8,-drdx(i,j,k1)))* &
701 & (dtdr(i-1,j,k2)+dtdr(i,j,k1)))*adfac2
702 ad_fx(i,j)=0.0_r8
703 END DO
704 END DO
705 END IF
706 IF ((k.eq.0).or.(k.eq.n(ng))) THEN
707 DO j=jstr-1,jend+1
708 DO i=istr-1,iend+1
709!^ tl_FS(i,j,k2)=0.0_r8
710!^
711 ad_fs(i,j,k2)=0.0_r8
712!^ tl_dTdr(i,j,k2)=0.0_r8
713!^
714 ad_dtdr(i,j,k2)=0.0_r8
715 END DO
716 END DO
717 ELSE
718 DO j=jstr-1,jend+1
719 DO i=istr-1,iend+1
720#if defined TS_MIX_MAX_SLOPE
721 cff1=sqrt(drdx(i,j,k2)**2+drdx(i+1,j,k2)**2+ &
722 & drdx(i,j,k1)**2+drdx(i+1,j,k1)**2+ &
723 & drde(i,j,k2)**2+drde(i,j+1,k2)**2+ &
724 & drde(i,j,k1)**2+drde(i,j+1,k1)**2)
725 cff2=0.25_r8*slope_max* &
726 & (z_r(i,j,k+1)-z_r(i,j,k))*cff1
727 cff3=max(pden(i,j,k)-pden(i,j,k+1),small)
728 cff4=max(cff2,cff3)
729 cff=-1.0_r8/cff4
730#elif defined TS_MIX_MIN_STRAT
731 cff1=max(pden(i,j,k)-pden(i,j,k+1), &
732 & strat_min*(z_r(i,j,k+1)-z_r(i,j,k)))
733 cff=-1.0_r8/cff1
734#else
735 cff1=max(pden(i,j,k)-pden(i,j,k+1),eps)
736 cff=-1.0_r8/cff1
737#endif
738!^ tl_FS(i,j,k2)=tl_cff*(z_r(i,j,k+1)-z_r(i,j,k))+ &
739!^ & cff*(tl_z_r(i,j,k+1)-tl_z_r(i,j,k))
740!^
741 adfac=cff*ad_fs(i,j,k2)
742 ad_z_r(i,j,k )=ad_z_r(i,j,k )-adfac
743 ad_z_r(i,j,k+1)=ad_z_r(i,j,k+1)+adfac
744 ad_cff=ad_cff+(z_r(i,j,k+1)- &
745 & z_r(i,j,k ))*ad_fs(i,j,k2)
746 ad_fs(i,j,k2)=0.0_r8
747#ifdef TS_MIX_CLIMA
748!^ tl_dTdr(i,j,k2)=tl_cff*((t(i,j,k+1,nrhs,itrc)- &
749!^ & tclm(i,j,k+1,itrc))- &
750!^ & (t(i,j,k ,nrhs,itrc)- &
751!^ & tclm(i,j,k ,itrc)))+ &
752!^ & cff*(tl_t(i,j,k+1,nrhs,itrc)- &
753!^ & tl_t(i,j,k ,nrhs,itrc))
754#else
755!^ tl_dTdr(i,j,k2)=tl_cff*(t(i,j,k+1,nrhs,itrc)- &
756!^ & t(i,j,k ,nrhs,itrc))+ &
757!^ & cff*(tl_t(i,j,k+1,nrhs,itrc)- &
758!^ & tl_t(i,j,k ,nrhs,itrc))
759#endif
760!^
761 adfac=cff*ad_dtdr(i,j,k2)
762 ad_t(i,j,k ,nrhs,itrc)=ad_t(i,j,k ,nrhs,itrc)-adfac
763 ad_t(i,j,k+1,nrhs,itrc)=ad_t(i,j,k+1,nrhs,itrc)+adfac
764#ifdef TS_MIX_CLIMA
765 ad_cff=ad_cff+((t(i,j,k+1,nrhs,itrc)- &
766 & tclm(i,j,k+1,itrc))- &
767 & (t(i,j,k ,nrhs,itrc)- &
768 & tclm(i,j,k ,itrc)))*ad_dtdr(i,j,k2)
769#else
770 ad_cff=ad_cff+(t(i,j,k+1,nrhs,itrc)- &
771 & t(i,j,k ,nrhs,itrc))*ad_dtdr(i,j,k2)
772#endif
773 ad_dtdr(i,j,k2)=0.0_r8
774#if defined TS_MIX_MAX_SLOPE
775!^ tl_cff=cff*cff*tl_cff4
776!^
777 ad_cff4=ad_cff4+cff*cff*ad_cff
778 ad_cff=0.0_r8
779!^ tl_cff4=(0.5_r8+SIGN(0.5_r8,cff2-cff3))*tl_cff2+ &
780!^ & (0.5_r8-SIGN(0.5_r8,cff2-cff3))*tl_cff3
781!^
782 ad_cff3=ad_cff3+ &
783 & (0.5_r8-sign(0.5_r8,cff2-cff3))*ad_cff4
784 ad_cff2=ad_cff2+ &
785 & (0.5_r8+sign(0.5_r8,cff2-cff3))*ad_cff4
786 ad_cff4=0.0_r8
787!^ tl_cff3=(0.5_r8+SIGN(0.5_r8,pden(i,j,k)-pden(i,j,k+1)- &
788!^ & small))* &
789!^ & (tl_pden(i,j,k)-tl_pden(i,j,k+1))
790!^
791 adfac=(0.5_r8+sign(0.5_r8,pden(i,j,k)-pden(i,j,k+1)- &
792 & small))*ad_cff3
793 ad_pden(i,j,k )=ad_pden(i,j,k )+adfac
794 ad_pden(i,j,k+1)=ad_pden(i,j,k+1)-adfac
795 ad_cff3=0.0_r8
796!^ tl_cff2=0.25_r8*slope_max* &
797!^ & ((tl_z_r(i,j,k+1)-tl_z_r(i,j,k))*cff1+ &
798!^ & (z_r(i,j,k+1)-z_r(i,j,k))*tl_cff1)
799!^
800 adfac=0.25_r8*slope_max*ad_cff2
801 adfac1=adfac*cff1
802 ad_cff1=ad_cff1+(z_r(i,j,k+1)-z_r(i,j,k))*adfac
803 ad_z_r(i,j,k )=ad_z_r(i,j,k )-adfac1
804 ad_z_r(i,j,k+1)=ad_z_r(i,j,k+1)+adfac1
805 ad_cff2=0.0_r8
806 IF (cff1.ne.0.0_r8) THEN
807!^ tl_cff1=(dRdx(i ,j,k2)*tl_dRdx(i ,j,k2)+ &
808!^ & dRdx(i+1,j,k2)*tl_dRdx(i+1,j,k2)+ &
809!^ & dRdx(i ,j,k1)*tl_dRdx(i ,j,k1)+ &
810!^ & dRdx(i+1,j,k1)*tl_dRdx(i+1,j,k1)+ &
811!^ & dRde(i,j ,k2)*tl_dRde(i,j ,k2)+ &
812!^ & dRde(i,j+1,k2)*tl_dRde(i,j+1,k2)+ &
813!^ & dRde(i,j ,k1)*tl_dRde(i,j ,k1)+ &
814!^ & dRde(i,j+1,k1)*tl_dRde(i,j+1,k1))/cff1
815!^
816 adfac=ad_cff1/cff1
817 ad_drdx(i ,j,k1)=ad_drdx(i ,j,k1)+ &
818 & drdx(i ,j,k1)*adfac
819 ad_drdx(i+1,j,k1)=ad_drdx(i+1,j,k1)+ &
820 & drdx(i+1,j,k1)*adfac
821 ad_drdx(i ,j,k2)=ad_drdx(i ,j,k2)+ &
822 & drdx(i ,j,k2)*adfac
823 ad_drdx(i+1,j,k2)=ad_drdx(i+1,j,k2)+ &
824 & drdx(i+1,j,k2)*adfac
825 ad_drde(i,j ,k2)=ad_drde(i,j ,k2)+ &
826 & drde(i,j ,k2)*adfac
827 ad_drde(i,j+1,k2)=ad_drde(i,j+1,k2)+ &
828 & drde(i,j+1,k2)*adfac
829 ad_drde(i,j ,k1)=ad_drde(i,j ,k1)+ &
830 & drde(i,j ,k1)*adfac
831 ad_drde(i,j+1,k1)=ad_drde(i,j+1,k1)+ &
832 & drde(i,j+1,k1)*adfac
833 ad_cff1=0.0_r8
834 ELSE
835!^ tl_cff1=0.0_r8
836!^
837 ad_cff1=0.0_r8
838 END IF
839#elif defined TS_MIX_MIN_STRAT
840!^ tl_cff=cff*cff*tl_cff1
841!^
842 ad_cff1=ad_cff1+cff*cff*ad_cff
843 ad_cff=0.0_r8
844!^ tl_cff1=(0.5_r8+SIGN(0.5_r8, &
845!^ & pden(i,j,k)-pden(i,j,k+1)- &
846!^ & strat_min*(z_r(i,j,k+1)- &
847!^ & z_r(i,j,k ))))* &
848!^ & (tl_pden(i,j,k)-tl_pden(i,j,k+1))+ &
849!^ & (0.5_r8-SIGN(0.5_r8, &
850!^ & pden(i,j,k)-pden(i,j,k+1)- &
851!^ & strat_min*(z_r(i,j,k+1)- &
852!^ & z_r(i,j,k ))))* &
853!^ & (strat_min*(tl_z_r(i,j,k+1)-tl_z_r(i,j,k )))
854!^
855 adfac1=(0.5_r8+sign(0.5_r8, &
856 & pden(i,j,k)-pden(i,j,k+1)- &
857 & strat_min*(z_r(i,j,k+1)- &
858 & z_r(i,j,k ))))* &
859 & ad_cff1
860 adfac2=(0.5_r8-sign(0.5_r8, &
861 & pden(i,j,k)-pden(i,j,k+1)- &
862 & strat_min*(z_r(i,j,k+1)- &
863 & z_r(i,j,k ))))* &
864 & strat_min*ad_cff1
865 ad_pden(i,j,k )=ad_pden(i,j,k )+adfac1
866 ad_pden(i,j,k+1)=ad_pden(i,j,k+1)-adfac1
867 ad_z_r(i,j,k )=ad_z_r(i,j,k )-adfac2
868 ad_z_r(i,j,k+1)=ad_z_r(i,j,k+1)+adfac2
869 ad_cff1=0.0_r8
870#else
871!^ tl_cff=cff*cff*tl_cff1
872!^
873 ad_cff1=ad_cff1+cff*cff*ad_cff
874 ad_cff=0.0_r8
875!^ tl_cff1=(0.5_r8+SIGN(0.5_r8, &
876!^ & pden(i,j,k)-pden(i,j,k+1)-eps))* &
877!^ & (tl_pden(i,j,k)-tl_pden(i,j,k+1))
878!^
879 adfac=(0.5_r8+sign(0.5_r8, &
880 & pden(i,j,k)-pden(i,j,k+1)-eps))* &
881 & ad_cff1
882 ad_pden(i,j,k )=ad_pden(i,j,k )+adfac
883 ad_pden(i,j,k+1)=ad_pden(i,j,k+1)-adfac
884 ad_cff1=0.0_r8
885#endif
886 END DO
887 END DO
888 END IF
889 IF (k.lt.n(ng)) THEN
890 DO j=jstr,jend+1
891 DO i=istr,iend
892 cff=0.5_r8*(pn(i,j)+pn(i,j-1))
893# ifdef MASKING
894 cff=cff*vmask(i,j)
895# endif
896!^ tl_dTde(i,j,k2)=cff*(tl_t(i,j ,k+1,nrhs,itrc)- &
897!^ & tl_t(i,j-1,k+1,nrhs,itrc))
898!^
899 adfac=cff*ad_dtde(i,j,k2)
900 ad_t(i,j-1,k+1,nrhs,itrc)=ad_t(i,j-1,k+1,nrhs,itrc)- &
901 & adfac
902 ad_t(i,j ,k+1,nrhs,itrc)=ad_t(i,j ,k+1,nrhs,itrc)+ &
903 & adfac
904 ad_dtde(i,j,k2)=0.0_r8
905!^ tl_dRde(i,j,k2)=cff*(tl_pden(i,j ,k+1)- &
906!^ & tl_pden(i,j-1,k+1))
907!^
908 adfac=cff*ad_drde(i,j,k2)
909 ad_pden(i,j-1,k+1)=ad_pden(i,j-1,k+1)-adfac
910 ad_pden(i,j ,k+1)=ad_pden(i,j ,k+1)+adfac
911 ad_drde(i,j,k2)=0.0_r8
912 END DO
913 END DO
914 DO j=jstr,jend
915 DO i=istr,iend+1
916 cff=0.5_r8*(pm(i,j)+pm(i-1,j))
917# ifdef MASKING
918 cff=cff*umask(i,j)
919# endif
920!^ tl_dTdx(i,j,k2)=cff*(tl_t(i ,j,k+1,nrhs,itrc)- &
921!^ & tl_t(i-1,j,k+1,nrhs,itrc))
922!^
923 adfac=cff*ad_dtdx(i,j,k2)
924 ad_t(i-1,j,k+1,nrhs,itrc)=ad_t(i-1,j,k+1,nrhs,itrc)- &
925 & adfac
926 ad_t(i ,j,k+1,nrhs,itrc)=ad_t(i ,j,k+1,nrhs,itrc)+ &
927 & adfac
928 ad_dtdx(i,j,k2)=0.0_r8
929!^ tl_dRdx(i,j,k2)=cff*(tl_pden(i ,j,k+1)- &
930!^ & tl_pden(i-1,j,k+1))
931!^
932 adfac=cff*ad_drdx(i,j,k2)
933 ad_pden(i-1,j,k+1)=ad_pden(i-1,j,k+1)-adfac
934 ad_pden(i ,j,k+1)=ad_pden(i ,j,k+1)+adfac
935 ad_drdx(i,j,k2)=0.0_r8
936 END DO
937 END DO
938 END IF
939!
940! Compute new storage recursive indices.
941!
942 kt=k2
943 k2=k1
944 k1=kt
945 END DO k_loop
946 END DO t_loop
947!
948 RETURN

References mod_scalars::dt.

◆ ad_t3dmix2_s_tile()

subroutine ad_t3dmix2_mod::ad_t3dmix2_s_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,
integer, intent(in) nstp,
integer, intent(in) nnew,
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) umask_wet,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) vmask_wet,
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), intent(in) pmon_u,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) pnom_v,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) pm,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) pn,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) diff3d_r,
real(r8), dimension(lbi:ubi,lbj:ubj,nt(ng)), intent(in) diff2,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),nt(ng)), intent(in) tclm,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),3,nt(ng)), intent(in) t,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),3,nt(ng)), intent(inout) ad_t )
private

Definition at line 93 of file ad_t3dmix2_s.h.

117!***********************************************************************
118!
119 USE mod_param
120 USE mod_scalars
121!
122! Imported variable declarations.
123!
124 integer, intent(in) :: ng, tile
125 integer, intent(in) :: LBi, UBi, LBj, UBj
126 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
127 integer, intent(in) :: nrhs, nstp, nnew
128
129#ifdef ASSUMED_SHAPE
130# ifdef MASKING
131 real(r8), intent(in) :: umask(LBi:,LBj:)
132 real(r8), intent(in) :: vmask(LBi:,LBj:)
133# endif
134# ifdef WET_DRY_NOT_YET
135 real(r8), intent(in) :: umask_wet(LBi:,LBj:)
136 real(r8), intent(in) :: vmask_wet(LBi:,LBj:)
137# endif
138# ifdef DIFF_3DCOEF
139 real(r8), intent(in) :: diff3d_r(LBi:,LBj:,:)
140# else
141 real(r8), intent(in) :: diff2(LBi:,LBj:,:)
142# endif
143 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
144 real(r8), intent(in) :: pmon_u(LBi:,LBj:)
145 real(r8), intent(in) :: pnom_v(LBi:,LBj:)
146 real(r8), intent(in) :: pm(LBi:,LBj:)
147 real(r8), intent(in) :: pn(LBi:,LBj:)
148# ifdef TS_MIX_CLIMA
149 real(r8), intent(in) :: tclm(LBi:,LBj:,:,:)
150# endif
151# ifdef DIAGNOSTICS_TS
152!! real(r8), intent(inout) :: DiaTwrk(LBi:,LBj:,:,:,:)
153# endif
154 real(r8), intent(in) :: t(LBi:,LBj:,:,:,:)
155
156 real(r8), intent(inout) :: ad_Hz(LBi:,LBj:,:)
157 real(r8), intent(inout) :: ad_t(LBi:,LBj:,:,:,:)
158#else
159# ifdef MASKING
160 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
161 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
162# endif
163# ifdef WET_DRY_NOT_YET
164 real(r8), intent(in) :: umask_wet(LBi:UBi,LBj:UBj)
165 real(r8), intent(in) :: vmask_wet(LBi:UBi,LBj:UBj)
166# endif
167# ifdef DIFF_3DCOEF
168 real(r8), intent(in) :: diff3d_r(LBi:UBi,LBj:UBj,N(ng))
169# else
170 real(r8), intent(in) :: diff2(LBi:UBi,LBj:UBj,NT(ng))
171# endif
172 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
173 real(r8), intent(in) :: pmon_u(LBi:UBi,LBj:UBj)
174 real(r8), intent(in) :: pnom_v(LBi:UBi,LBj:UBj)
175 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
176 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
177# ifdef TS_MIX_CLIMA
178 real(r8), intent(in) :: tclm(LBi:UBi,LBj:UBj,N(ng),NT(ng))
179# endif
180# ifdef DIAGNOSTICS_TS
181!! real(r8), intent(inout) :: DiaTwrk(LBi:UBi,LBj:UBj,N(ng),NT(ng), &
182!! & NDT)
183# endif
184 real(r8), intent(in) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
185
186 real(r8), intent(inout) :: ad_Hz(LBi:UBi,LBj:UBj,N(ng))
187 real(r8), intent(inout) :: ad_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
188#endif
189!
190! Local variable declarations.
191!
192 integer :: i, itrc, j, k
193
194 real(r8) :: cff, cff1
195 real(r8) :: adfac, adfac1, adfac2, adfac3, adfac4, ad_cff, ad_cff1
196
197 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_FE
198 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_FX
199
200#include "set_bounds.h"
201!
202!-----------------------------------------------------------------------
203! Initialize adjoint private variables.
204!-----------------------------------------------------------------------
205!
206 ad_cff=0.0_r8
207 ad_cff1=0.0_r8
208
209 ad_fe(imins:imaxs,jmins:jmaxs)=0.0_r8
210 ad_fx(imins:imaxs,jmins:jmaxs)=0.0_r8
211!
212!----------------------------------------------------------------------
213! Compute adjoint horizontal harmonic diffusion along constant
214! S-surfaces.
215!----------------------------------------------------------------------
216!
217 DO itrc=1,nt(ng)
218 DO k=1,n(ng)
219!
220! Time-step harmonic, S-surfaces diffusion term (m Tunits).
221!
222 DO j=jstr,jend
223 DO i=istr,iend
224#ifdef DIAGNOSTICS_TS
225!! DiaTwrk(i,j,k,itrc,iThdif)=cff
226#endif
227!^ tl_t(i,j,k,nnew,itrc)=tl_t(i,j,k,nnew,itrc)+tl_cff
228!^
229 ad_cff=ad_cff+ad_t(i,j,k,nnew,itrc)
230!^ tl_cff=dt(ng)*pm(i,j)*pn(i,j)* &
231!^ & (tl_FX(i+1,j)-tl_FX(i,j)+ &
232!^ & tl_FE(i,j+1)-tl_FE(i,j))
233!^
234 adfac=dt(ng)*pm(i,j)*pn(i,j)*ad_cff
235 ad_fx(i ,j)=ad_fx(i ,j)-adfac
236 ad_fx(i+1,j)=ad_fx(i+1,j)+adfac
237 ad_fe(i,j )=ad_fe(i,j )-adfac
238 ad_fe(i,j+1)=ad_fe(i,j+1)+adfac
239 ad_cff=0.0_r8
240 END DO
241 END DO
242!
243! Compute XI- and ETA-components of diffusive tracer flux (T m3/s).
244!
245 DO j=jstr,jend+1
246 DO i=istr,iend
247#ifdef DIFF_3DCOEF
248 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i,j-1,k))* &
249 & pnom_v(i,j)
250#else
251 cff=0.25_r8*(diff2(i,j,itrc)+diff2(i,j-1,itrc))* &
252 & pnom_v(i,j)
253#endif
254#ifdef WET_DRY_NOT_YET
255 fe(i,j)=fe(i,j)*vmask_wet(i,j)
256#endif
257#ifdef MASKING
258!^ tl_FE(i,j)=tl_FE(i,j)*vmask(i,j)
259!^
260 ad_fe(i,j)=ad_fe(i,j)*vmask(i,j)
261#endif
262#if defined TS_MIX_STABILITY
263!^ tl_FE(i,j)=cff* &
264!^ & ((tl_Hz(i,j,k)+tl_Hz(i,j-1,k))* &
265!^ & (0.75_r8*(t(i,j ,k,nrhs,itrc)- &
266!^ & t(i,j-1,k,nrhs,itrc))+ &
267!^ & 0.25_r8*(t(i,j ,k,nstp,itrc)- &
268!^ & t(i,j-1,k,nstp,itrc)))+ &
269!^ & (Hz(i,j,k)+Hz(i,j-1,k))* &
270!^ & (0.75_r8*(tl_t(i,j ,k,nrhs,itrc)- &
271!^ & tl_t(i,j-1,k,nrhs,itrc))+ &
272!^ & 0.25_r8*(tl_t(i,j ,k,nstp,itrc)- &
273!^ & tl_t(i,j-1,k,nstp,itrc))))
274!^
275 adfac=cff*ad_fe(i,j)
276 adfac1=adfac*(0.75_r8*(t(i,j ,k,nrhs,itrc)- &
277 & t(i,j-1,k,nrhs,itrc))+ &
278 & 0.25_r8*(t(i,j ,k,nstp,itrc)- &
279 & t(i,j-1,k,nstp,itrc)))
280 adfac2=adfac*(hz(i,j,k)+hz(i,j-1,k))
281 adfac3=adfac2*0.75_r8
282 adfac4=adfac2*0.25_r8
283 ad_hz(i,j-1,k)=ad_hz(i,j-1,k)+adfac1
284 ad_hz(i,j ,k)=ad_hz(i,j ,k)+adfac1
285 ad_t(i,j-1,k,nrhs,itrc)=ad_t(i,j-1,k,nrhs,itrc)-adfac3
286 ad_t(i,j ,k,nrhs,itrc)=ad_t(i,j ,k,nrhs,itrc)+adfac3
287 ad_t(i,j-1,k,nstp,itrc)=ad_t(i,j-1,k,nstp,itrc)-adfac4
288 ad_t(i,j ,k,nstp,itrc)=ad_t(i,j ,k,nstp,itrc)+adfac4
289 ad_fe(i,j)=0.0_r8
290#elif defined TS_MIX_CLIMA
291 IF (ltracerclm(itrc,ng)) THEN
292!^ tl_FE(i,j)=cff* &
293!^ & ((tl_Hz(i,j,k)+tl_Hz(i,j-1,k))* &
294!^ & ((t(i,j ,k,nrhs,itrc)- &
295!^ & tclm(i,j ,k,itrc))- &
296!^ & (t(i,j-1,k,nrhs,itrc)- &
297!^ & tclm(i,j-1,k,itrc)))+ &
298!^ & (Hz(i,j,k)+Hz(i,j-1,k))* &
299!^ & (tl_t(i,j ,k,nrhs,itrc)- &
300!^ & tl_t(i,j-1,k,nrhs,itrc)))
301!^
302 adfac=cff*ad_fe(i,j)
303 adfac1=adfac*((t(i,j ,k,nrhs,itrc)- &
304 & tclm(i,j ,k,itrc))- &
305 & (t(i,j-1,k,nrhs,itrc)- &
306 & tclm(i,j-1,k,itrc)))
307 adfac2=adfac*(hz(i,j,k)+hz(i,j-1,k))
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_t(i,j-1,k,nrhs,itrc)=ad_t(i,j-1,k,nrhs,itrc)-adfac2
311 ad_t(i,j ,k,nrhs,itrc)=ad_t(i,j ,k,nrhs,itrc)+adfac2
312 ad_fe(i,j)=0.0_r8
313 ELSE
314!^ tl_FE(i,j)=cff* &
315!^ & ((tl_Hz(i,j,k)+tl_Hz(i,j-1,k))* &
316!^ & (t(i,j ,k,nrhs,itrc)- &
317!^ & t(i,j-1,k,nrhs,itrc))+ &
318!^ & (Hz(i,j,k)+Hz(i,j-1,k))* &
319!^ & (tl_t(i,j ,k,nrhs,itrc)- &
320!^ & tl_t(i,j-1,k,nrhs,itrc)))
321!^
322 adfac=cff*ad_fe(i,j)
323 adfac1=adfac*(t(i,j,k,nrhs,itrc)-t(i,j-1,k,nrhs,itrc))
324 adfac2=adfac*(hz(i,j,k)+hz(i,j-1,k))
325 ad_hz(i,j-1,k)=ad_hz(i,j-1,k)+adfac1
326 ad_hz(i,j ,k)=ad_hz(i,j ,k)+adfac1
327 ad_t(i,j-1,k,nrhs,itrc)=ad_t(i,j-1,k,nrhs,itrc)-adfac2
328 ad_t(i,j ,k,nrhs,itrc)=ad_t(i,j ,k,nrhs,itrc)+adfac2
329 ad_fe(i,j)=0.0_r8
330 END IF
331#else
332!^ tl_FE(i,j)=cff* &
333!^ & ((tl_Hz(i,j,k)+tl_Hz(i,j-1,k))* &
334!^ & (t(i,j ,k,nrhs,itrc)- &
335!^ & t(i,j-1,k,nrhs,itrc))+ &
336!^ & (Hz(i,j,k)+Hz(i,j-1,k))* &
337!^ & (tl_t(i,j ,k,nrhs,itrc)- &
338!^ & tl_t(i,j-1,k,nrhs,itrc)))
339!^
340 adfac=cff*ad_fe(i,j)
341 adfac1=adfac*(t(i,j,k,nrhs,itrc)-t(i,j-1,k,nrhs,itrc))
342 adfac2=adfac*(hz(i,j,k)+hz(i,j-1,k))
343 ad_hz(i,j-1,k)=ad_hz(i,j-1,k)+adfac1
344 ad_hz(i,j ,k)=ad_hz(i,j ,k)+adfac1
345 ad_t(i,j-1,k,nrhs,itrc)=ad_t(i,j-1,k,nrhs,itrc)-adfac2
346 ad_t(i,j ,k,nrhs,itrc)=ad_t(i,j ,k,nrhs,itrc)+adfac2
347 ad_fe(i,j)=0.0_r8
348#endif
349 END DO
350 END DO
351 DO j=jstr,jend
352 DO i=istr,iend+1
353#ifdef DIFF_3DCOEF
354 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i-1,j,k))* &
355 & pmon_u(i,j)
356#else
357 cff=0.25_r8*(diff2(i,j,itrc)+diff2(i-1,j,itrc))* &
358 & pmon_u(i,j)
359#endif
360#ifdef WET_DRY_NOT_YET
361 fx(i,j)=fx(i,j)*umask_wet(i,j)
362#endif
363#ifdef MASKING
364!^ tl_FX(i,j)=tl_FX(i,j)*umask(i,j)
365!^
366 ad_fx(i,j)=ad_fx(i,j)*umask(i,j)
367#endif
368#if defined TS_MIX_STABILITY
369!^ tl_FX(i,j)=cff* &
370!^ & ((tl_Hz(i,j,k)+tl_Hz(i-1,j,k))* &
371!^ & (0.75_r8*(t(i ,j,k,nrhs,itrc)- &
372!^ & t(i-1,j,k,nrhs,itrc))+ &
373!^ & 0.25_r8*(t(i ,j,k,nstp,itrc)- &
374!^ & t(i-1,j,k,nstp,itrc)))+ &
375!^ & (Hz(i,j,k)+Hz(i-1,j,k))* &
376!^ & (0.75_r8*(tl_t(i ,j,k,nrhs,itrc)- &
377!^ & tl_t(i-1,j,k,nrhs,itrc))+ &
378!^ & 0.25_r8*(tl_t(i ,j,k,nstp,itrc)- &
379!^ & tl_t(i-1,j,k,nstp,itrc))))
380!^
381 adfac=cff*ad_fx(i,j)
382 adfac1=adfac*(0.75_r8*(t(i ,j,k,nrhs,itrc)- &
383 & t(i-1,j,k,nrhs,itrc))+ &
384 & 0.25_r8*(t(i ,j,k,nstp,itrc)- &
385 & t(i-1,j,k,nstp,itrc)))
386 adfac2=adfac*(hz(i,j,k)+hz(i-1,j,k))
387 adfac3=adfac2*0.75_r8
388 adfac4=adfac2*0.25_r8
389 ad_hz(i-1,j,k)=ad_hz(i-1,j,k)+adfac1
390 ad_hz(i ,j,k)=ad_hz(i ,j,k)+adfac1
391 ad_t(i-1,j,k,nrhs,itrc)=ad_t(i-1,j,k,nrhs,itrc)-adfac3
392 ad_t(i ,j,k,nrhs,itrc)=ad_t(i ,j,k,nrhs,itrc)+adfac3
393 ad_t(i-1,j,k,nstp,itrc)=ad_t(i-1,j,k,nstp,itrc)-adfac4
394 ad_t(i ,j,k,nstp,itrc)=ad_t(i ,j,k,nstp,itrc)+adfac4
395 ad_fx(i,j)=0.0_r8
396#elif defined TS_MIX_CLIMA
397 IF (ltracerclm(itrc,ng)) THEN
398!^ tl_FX(i,j)=cff* &
399!^ & ((tl_Hz(i,j,k)+tl_Hz(i-1,j,k))* &
400!^ & ((t(i ,j,k,nrhs,itrc)- &
401!^ & tclm(i ,j,k,itrc))- &
402!^ & (t(i-1,j,k,nrhs,itrc)- &
403!^ & tclm(i-1,j,k,itrc)))+ &
404!^ & (Hz(i,j,k)+Hz(i-1,j,k))* &
405!^ & (tl_t(i ,j,k,nrhs,itrc)- &
406!^ & tl_t(i-1,j,k,nrhs,itrc)))
407!^
408 adfac=cff*ad_fx(i,j)
409 adfac1=adfac*((t(i ,j,k,nrhs,itrc)- &
410 & tclm(i ,j,k,itrc))- &
411 & (t(i-1,j,k,nrhs,itrc)- &
412 & tclm(i-1,j,k,itrc)))
413 adfac2=adfac*(hz(i,j,k)+hz(i-1,j,k))
414 ad_hz(i-1,j,k)=ad_hz(i-1,j,k)+adfac1
415 ad_hz(i ,j,k)=ad_hz(i ,j,k)+adfac1
416 ad_t(i-1,j,k,nrhs,itrc)=ad_t(i-1,j,k,nrhs,itrc)-adfac2
417 ad_t(i ,j,k,nrhs,itrc)=ad_t(i ,j,k,nrhs,itrc)+adfac2
418 ad_fx(i,j)=0.0_r8
419 ELSE
420!^ tl_FX(i,j)=cff* &
421!^ & ((tl_Hz(i,j,k)+tl_Hz(i-1,j,k))* &
422!^ & (t(i ,j,k,nrhs,itrc)- &
423!^ & t(i-1,j,k,nrhs,itrc))+ &
424!^ & (Hz(i,j,k)+Hz(i-1,j,k))* &
425!^ & (tl_t(i ,j,k,nrhs,itrc)- &
426!^ & tl_t(i-1,j,k,nrhs,itrc)))
427!^
428 adfac=cff*ad_fx(i,j)
429 adfac1=adfac*(t(i,j,k,nrhs,itrc)-t(i-1,j,k,nrhs,itrc))
430 adfac2=adfac*(hz(i,j,k)+hz(i-1,j,k))
431 ad_hz(i-1,j,k)=ad_hz(i-1,j,k)+adfac1
432 ad_hz(i ,j,k)=ad_hz(i ,j,k)+adfac1
433 ad_t(i-1,j,k,nrhs,itrc)=ad_t(i-1,j,k,nrhs,itrc)-adfac2
434 ad_t(i ,j,k,nrhs,itrc)=ad_t(i ,j,k,nrhs,itrc)+adfac2
435 ad_fx(i,j)=0.0_r8
436 END IF
437#else
438!^ tl_FX(i,j)=cff* &
439!^ & ((tl_Hz(i,j,k)+tl_Hz(i-1,j,k))* &
440!^ & (t(i ,j,k,nrhs,itrc)- &
441!^ & t(i-1,j,k,nrhs,itrc))+ &
442!^ & (Hz(i,j,k)+Hz(i-1,j,k))* &
443!^ & (tl_t(i ,j,k,nrhs,itrc)- &
444!^ & tl_t(i-1,j,k,nrhs,itrc)))
445!^
446 adfac=cff*ad_fx(i,j)
447 adfac1=adfac*(t(i,j,k,nrhs,itrc)-t(i-1,j,k,nrhs,itrc))
448 adfac2=adfac*(hz(i,j,k)+hz(i-1,j,k))
449 ad_hz(i-1,j,k)=ad_hz(i-1,j,k)+adfac1
450 ad_hz(i ,j,k)=ad_hz(i ,j,k)+adfac1
451 ad_t(i-1,j,k,nrhs,itrc)=ad_t(i-1,j,k,nrhs,itrc)-adfac2
452 ad_t(i ,j,k,nrhs,itrc)=ad_t(i ,j,k,nrhs,itrc)+adfac2
453 ad_fx(i,j)=0.0_r8
454#endif
455 END DO
456 END DO
457 END DO
458 END DO
459!
460 RETURN

References mod_scalars::dt, and mod_scalars::ltracerclm.