ROMS
Loading...
Searching...
No Matches
ad_t3dmix2_iso.h
Go to the documentation of this file.
1 MODULE ad_t3dmix2_mod
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 isopycnic surfaces. !
12! !
13! BASIC STATE variables needed: diff2, Hz, rho, 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, 26, __line__, myfile)
53#endif
54 CALL ad_t3dmix2_iso_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 & grid(ng) % om_v, &
63 & grid(ng) % on_u, &
64 & grid(ng) % pm, &
65 & grid(ng) % pn, &
66 & grid(ng) % Hz, &
67 & grid(ng) % ad_Hz, &
68 & grid(ng) % z_r, &
69 & grid(ng) % ad_z_r, &
70 & mixing(ng) % diff2, &
71 & ocean(ng) % pden, &
72 & ocean(ng) % ad_pden, &
73#ifdef TS_MIX_CLIMA
74 & clima(ng) % tclm, &
75#endif
76#ifdef DIAGNOSTICS_TS
77!! & DIAGS(ng) % DiaTwrk, &
78#endif
79 & ocean(ng) % t, &
80 & ocean(ng) % ad_t)
81#ifdef PROFILE
82 CALL wclock_off (ng, iadm, 26, __line__, myfile)
83#endif
84!
85 RETURN
86 END SUBROUTINE ad_t3dmix2
87!
88!***********************************************************************
89 SUBROUTINE ad_t3dmix2_iso_tile (ng, tile, &
90 & LBi, UBi, LBj, UBj, &
91 & IminS, ImaxS, JminS, JmaxS, &
92 & nrhs, nstp, nnew, &
93#ifdef MASKING
94 & umask, vmask, &
95#endif
96 & om_v, on_u, pm, pn, &
97 & Hz, ad_Hz, &
98 & z_r, ad_z_r, &
99 & diff2, &
100 & pden, ad_pden, &
101#ifdef TS_MIX_CLIMA
102 & tclm, &
103#endif
104#ifdef DIAGNOSTICS_TS
105!! & DiaTwrk, &
106#endif
107 & t, ad_t)
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
949 END SUBROUTINE ad_t3dmix2_iso_tile
950
951 END MODULE ad_t3dmix2_mod
subroutine, public ad_t3dmix2(ng, tile)
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)
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
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