ROMS
Loading...
Searching...
No Matches
ad_t3dmix2_geo.h
Go to the documentation of this file.
2!
3!git $Id$
4!================================================== Hernan G. Arango ===
5! Copyright (c) 2002-2025 The ROMS Group Andrew M. Moore !
6! Licensed under a MIT/X style license !
7! See License_ROMS.md !
8!=======================================================================
9! !
10! This subroutine computes adjoint horizontal harmonic mixing of !
11! tracers along geopotential surfaces. !
12! !
13! BASIC STATE variables needed: diff2, Hz, t, z_r !
14! !
15!=======================================================================
16!
17 implicit none
18!
19 PRIVATE
20 PUBLIC ad_t3dmix2
21!
22 CONTAINS
23!
24!***********************************************************************
25 SUBROUTINE ad_t3dmix2 (ng, tile)
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
92 END SUBROUTINE ad_t3dmix2
93!
94!***********************************************************************
95 SUBROUTINE ad_t3dmix2_geo_tile (ng, tile, &
96 & LBi, UBi, LBj, UBj, &
97 & IminS, ImaxS, JminS, JmaxS, &
98 & nrhs, nstp, nnew, &
99#ifdef MASKING
100 & umask, vmask, &
101#endif
102#ifdef WET_DRY_NOT_YET
103 & umask_wet, vmask_wet, &
104#endif
105 & om_v, on_u, pm, pn, &
106 & Hz, ad_Hz, &
107 & z_r, ad_z_r, &
108#ifdef DIFF_3DCOEF
109 & diff3d_r, &
110#else
111 & diff2, &
112#endif
113#ifdef TS_MIX_CLIMA
114 & tclm, &
115#endif
116#ifdef DIAGNOSTICS_TS
117!! & DiaTwrk, &
118#endif
119 & t, ad_t)
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
962 END SUBROUTINE ad_t3dmix2_geo_tile
963
964 END MODULE ad_t3dmix2_mod
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)
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
real(dp), dimension(:), allocatable dt
logical, dimension(:,:), allocatable ltracerclm
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