ROMS
Loading...
Searching...
No Matches
ad_pre_step3d.F
Go to the documentation of this file.
1#include "cppdefs.h"
2 MODULE ad_pre_step3d_mod
3
4#if defined ADJOINT && defined SOLVE3D
5# define NEUMANN
6!
7!git $Id$
8!================================================== Hernan G. Arango ===
9! Copyright (c) 2002-2025 The ROMS Group Andrew M. Moore !
10! Licensed under a MIT/X style license !
11! See License_ROMS.md !
12!=======================================================================
13! !
14! This subroutine initialize computations for new time step of the !
15! adjoint 3D primitive variables. !
16! !
17! Both n-1 and n-2 time-step contributions of the Adams/Bashforth !
18! scheme are added here to u and v at time index "nnew", since the !
19! right-hand-side arrays ru and rv at n-2 will be overwriten in !
20! subsequent calls to routines within the 3D engine. !
21! !
22! It also computes the time "n" vertical viscosity and diffusion !
23! contributions of the Crank-Nicholson implicit scheme because the !
24! thicknesses "Hz" will be overwriten at the end of the 2D engine !
25! (barotropic mode) computations. !
26! !
27! The actual time step will be carried out in routines "step3d_uv" !
28! and "step3d_t". !
29! !
30! BASIC STATE variables needed: AKt, AKv, Huon, Hvom, Hz, Tsrc, W, !
31! bustr, bvstr, ghats, srflx, sustr, !
32! svstr, t, z_r, z_w !
33! !
34!=======================================================================
35!
36 implicit none
37!
38 PRIVATE
39 PUBLIC :: ad_pre_step3d
40!
41 CONTAINS
42!
43!***********************************************************************
44 SUBROUTINE ad_pre_step3d (ng, tile)
45!***********************************************************************
46!
47 USE mod_param
48# ifdef DIAGNOSTICS
49!! USE mod_diags
50# endif
51 USE mod_forces
52 USE mod_grid
53 USE mod_mixing
54 USE mod_ocean
55 USE mod_stepping
56!
57! Imported variable declarations.
58!
59 integer, intent(in) :: ng, tile
60!
61! Local variable declarations.
62!
63 character (len=*), parameter :: myfile = &
64 & __FILE__
65!
66# include "tile.h"
67!
68# ifdef PROFILE
69 CALL wclock_on (ng, iadm, 22, __line__, myfile)
70# endif
71 CALL ad_pre_step3d_tile (ng, tile, &
72 & lbi, ubi, lbj, ubj, &
73 & imins, imaxs, jmins, jmaxs, &
74# ifdef FOUR_DVAR
75 & kstp(ng), knew(ng), &
76# endif
77 & nrhs(ng), nstp(ng), nnew(ng), &
78# ifdef MASKING
79 & grid(ng) % rmask, &
80 & grid(ng) % umask, &
81 & grid(ng) % vmask, &
82# if defined SOLAR_SOURCE && defined WET_DRY
83 & grid(ng) % rmask_wet, &
84# endif
85# endif
86# ifdef FOUR_DVAR
87 & grid(ng) % om_v, &
88 & grid(ng) % on_u, &
89# endif
90 & grid(ng) % pm, &
91 & grid(ng) % pn, &
92 & grid(ng) % Hz, &
93 & grid(ng) % ad_Hz, &
94 & grid(ng) % Huon, &
95 & grid(ng) % ad_Huon, &
96 & grid(ng) % Hvom, &
97 & grid(ng) % ad_Hvom, &
98 & grid(ng) % z_r, &
99 & grid(ng) % ad_z_r, &
100 & grid(ng) % z_w, &
101 & grid(ng) % ad_z_w, &
102 & forces(ng) % ad_btflx, &
103 & forces(ng) % ad_bustr, &
104 & forces(ng) % ad_bvstr, &
105 & forces(ng) % ad_stflx, &
106 & forces(ng) % ad_sustr, &
107 & forces(ng) % ad_svstr, &
108# ifdef SOLAR_SOURCE
109 & forces(ng) % srflx, &
110# endif
111 & mixing(ng) % Akt, &
112 & mixing(ng) % ad_Akt, &
113 & mixing(ng) % Akv, &
114 & mixing(ng) % ad_Akv, &
115# ifdef LMD_NONLOCAL_NOT_YET
116 & mixing(ng) % ad_ghats, &
117# endif
118# ifdef FOUR_DVAR
119 & ocean(ng) % ad_ubar, &
120 & ocean(ng) % ad_vbar, &
121# endif
122 & ocean(ng) % W, &
123 & ocean(ng) % ad_W, &
124 & ocean(ng) % ad_ru, &
125 & ocean(ng) % ad_rv, &
126# ifdef DIAGNOSTICS_TS
127!! & DIAGS(ng) % DiaTwrk, &
128# endif
129# ifdef DIAGNOSTICS_UV
130!! & DIAGS(ng) % DiaU3wrk, &
131!! & DIAGS(ng) % DiaV3wrk, &
132!! & DIAGS(ng) % DiaRU, &
133!! & DIAGS(ng) % DiaRV, &
134# endif
135 & ocean(ng) % t, &
136 & ocean(ng) % ad_t, &
137 & ocean(ng) % u, &
138 & ocean(ng) % ad_u, &
139 & ocean(ng) % v, &
140 & ocean(ng) % ad_v)
141# ifdef PROFILE
142 CALL wclock_off (ng, iadm, 22, __line__, myfile)
143# endif
144!
145 RETURN
146 END SUBROUTINE ad_pre_step3d
147!
148!***********************************************************************
149 SUBROUTINE ad_pre_step3d_tile (ng, tile, &
150 & LBi, UBi, LBj, UBj, &
151 & IminS, ImaxS, JminS, JmaxS, &
152# ifdef FOUR_DVAR
153 & kstp, knew, &
154# endif
155 & nrhs, nstp, nnew, &
156# ifdef MASKING
157 & rmask, umask, vmask, &
158# if defined SOLAR_SOURCE && defined WET_DRY
159 & rmask_wet, &
160# endif
161# endif
162# ifdef FOUR_DVAR
163 & om_v, on_u, &
164# endif
165 & pm, pn, &
166 & Hz, ad_Hz, &
167 & Huon, ad_Huon, &
168 & Hvom, ad_Hvom, &
169 & z_r, ad_z_r, &
170 & z_w, ad_z_w, &
171 & ad_btflx, ad_bustr, ad_bvstr, &
172 & ad_stflx, ad_sustr, ad_svstr, &
173# ifdef SOLAR_SOURCE
174 & srflx, &
175# endif
176 & Akt, ad_Akt, &
177 & Akv, ad_Akv, &
178# ifdef LMD_NONLOCAL_NOT_YET
179 & ad_ghats, &
180# endif
181# ifdef FOUR_DVAR
182 & ad_ubar, ad_vbar, &
183# endif
184 & W, ad_W, &
185 & ad_ru, ad_rv, &
186# ifdef DIAGNOSTICS_TS
187!! & DiaTwrk, &
188# endif
189# ifdef DIAGNOSTICS_UV
190!! & DiaU3wrk, DiaV3wrk, &
191!! & DiaRU, DiaRV, &
192# endif
193 & t, ad_t, &
194 & u, ad_u, &
195 & v, ad_v)
196!***********************************************************************
197!
198 USE mod_param
199 USE mod_scalars
200 USE mod_sources
201!
202 USE ad_exchange_3d_mod, ONLY : ad_exchange_r3d_tile
203# ifdef DISTRIBUTE
204 USE mp_exchange_mod, ONLY : ad_mp_exchange4d
205# endif
206 USE ad_t3dbc_mod, ONLY : ad_t3dbc_tile
207!
208! Imported variable declarations.
209!
210 integer, intent(in) :: ng, tile
211 integer, intent(in) :: LBi, UBi, LBj, UBj
212 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
213# ifdef FOUR_DVAR
214 integer, intent(in) :: kstp, knew
215# endif
216 integer, intent(in) :: nrhs, nstp, nnew
217!
218# ifdef ASSUMED_SHAPE
219# ifdef MASKING
220 real(r8), intent(in) :: rmask(LBi:,LBj:)
221 real(r8), intent(in) :: umask(LBi:,LBj:)
222 real(r8), intent(in) :: vmask(LBi:,LBj:)
223# if defined SOLAR_SOURCE && defined WET_DRY
224 real(r8), intent(in) :: rmask_wet(LBi:,LBj:)
225# endif
226# endif
227# ifdef FOUR_DVAR
228 real(r8), intent(in) :: om_v(LBi:,LBj:)
229 real(r8), intent(in) :: on_u(LBi:,LBj:)
230# endif
231 real(r8), intent(in) :: pm(LBi:,LBj:)
232 real(r8), intent(in) :: pn(LBi:,LBj:)
233 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
234 real(r8), intent(in) :: Huon(LBi:,LBj:,:)
235 real(r8), intent(in) :: Hvom(LBi:,LBj:,:)
236 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
237 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
238# ifdef SOLAR_SOURCE
239 real(r8), intent(in) :: srflx(LBi:,LBj:)
240# endif
241# ifdef SUN
242 real(r8), intent(in) :: Akt(LBi:UBi,LBj:UBj,0:N(ng),NAT)
243# else
244 real(r8), intent(in) :: Akt(LBi:,LBj:,0:,:)
245# endif
246 real(r8), intent(in) :: Akv(LBi:,LBj:,0:)
247# ifdef FOUR_DVAR
248 real(r8), intent(in) :: ad_ubar(LBi:,LBj:,:)
249 real(r8), intent(in) :: ad_vbar(LBi:,LBj:,:)
250# endif
251 real(r8), intent(in) :: W(LBi:,LBj:,0:)
252# ifdef SUN
253 real(r8), intent(in) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
254# else
255 real(r8), intent(in) :: t(LBi:,LBj:,:,:,:)
256# endif
257 real(r8), intent(in) :: u(LBi:,LBj:,:,:)
258 real(r8), intent(in) :: v(LBi:,LBj:,:,:)
259
260# ifdef DIAGNOSTICS_TS
261!! real(r8), intent(inout) :: DiaTwrk(LBi:,LBj:,:,:,:)
262# endif
263# ifdef DIAGNOSTICS_UV
264!! real(r8), intent(inout) :: DiaU3wrk(LBi:,LBj:,:,:)
265!! real(r8), intent(inout) :: DiaV3wrk(LBi:,LBj:,:,:)
266!! real(r8), intent(inout) :: DiaRU(LBi:,LBj:,:,:,:)
267!! real(r8), intent(inout) :: DiaRV(LBi:,LBj:,:,:,:)
268# endif
269 real(r8), intent(inout) :: ad_Hz(LBi:,LBj:,:)
270 real(r8), intent(inout) :: ad_Huon(LBi:,LBj:,:)
271 real(r8), intent(inout) :: ad_Hvom(LBi:,LBj:,:)
272 real(r8), intent(inout) :: ad_z_r(LBi:,LBj:,:)
273 real(r8), intent(inout) :: ad_z_w(LBi:,LBj:,0:)
274 real(r8), intent(inout) :: ad_btflx(LBi:,LBj:,:)
275 real(r8), intent(inout) :: ad_bustr(LBi:,LBj:)
276 real(r8), intent(inout) :: ad_bvstr(LBi:,LBj:)
277 real(r8), intent(inout) :: ad_stflx(LBi:,LBj:,:)
278 real(r8), intent(inout) :: ad_sustr(LBi:,LBj:)
279 real(r8), intent(inout) :: ad_svstr(LBi:,LBj:)
280 real(r8), intent(inout) :: ad_ru(LBi:,LBj:,0:,:)
281 real(r8), intent(inout) :: ad_rv(LBi:,LBj:,0:,:)
282# ifdef LMD_NONLOCAL_NOT_YET
283 real(r8), intent(inout) :: ad_ghats(LBi:,LBj:,0:,:)
284# endif
285# ifdef SUN
286 real(r8), intent(inout) :: ad_Akt(LBi:UBi,LBj:UBj,0:N(ng),NAT)
287# else
288 real(r8), intent(inout) :: ad_Akt(LBi:,LBj:,0:,:)
289# endif
290 real(r8), intent(inout) :: ad_Akv(LBi:,LBj:,0:)
291 real(r8), intent(inout) :: ad_W(LBi:,LBj:,0:)
292# ifdef SUN
293 real(r8), intent(inout) :: ad_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
294# else
295 real(r8), intent(inout) :: ad_t(LBi:,LBj:,:,:,:)
296# endif
297 real(r8), intent(inout) :: ad_u(LBi:,LBj:,:,:)
298 real(r8), intent(inout) :: ad_v(LBi:,LBj:,:,:)
299
300# else
301
302# ifdef MASKING
303 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
304 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
305 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
306# if defined SOLAR_SOURCE && defined WET_DRY
307 real(r8), intent(in) :: rmask_wet(LBi:UBi,LBj:UBj)
308# endif
309# endif
310# ifdef FOUR_DVAR
311 real(r8), intent(in) :: om_v(LBi:UBi,LBj:UBj)
312 real(r8), intent(in) :: on_u(LBi:UBi,LBj:UBj)
313# endif
314 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
315 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
316 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
317 real(r8), intent(in) :: Huon(LBi:UBi,LBj:UBj,N(ng))
318 real(r8), intent(in) :: Hvom(LBi:UBi,LBj:UBj,N(ng))
319 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
320 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
321# ifdef SOLAR_SOURCE
322 real(r8), intent(in) :: srflx(LBi:UBi,LBj:UBj)
323# endif
324 real(r8), intent(in) :: Akt(LBi:UBi,LBj:UBj,0:N(ng),NAT)
325 real(r8), intent(in) :: Akv(LBi:UBi,LBj:UBj,0:N(ng))
326# ifdef FOUR_DVAR
327 real(r8), intent(in) :: ad_ubar(LBi:UBi,LBj:UBj,:)
328 real(r8), intent(in) :: ad_vbar(LBi:UBi,LBj:UBj,:)
329# endif
330 real(r8), intent(in) :: W(LBi:UBi,LBj:UBj,0:N(ng))
331 real(r8), intent(in) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
332 real(r8), intent(in) :: u(LBi:UBi,LBj:UBj,N(ng),2)
333 real(r8), intent(in) :: v(LBi:UBi,LBj:UBj,N(ng),2)
334
335# ifdef DIAGNOSTICS_TS
336!! real(r8), intent(inout) :: DiaTwrk(LBi:UBi,LBj:UBj,N(ng),NT(ng), &
337!! & NDT)
338# endif
339# ifdef DIAGNOSTICS_UV
340!! real(r8), intent(inout) :: DiaU3wrk(LBi:UBi,LBj:UBj,N(ng),NDM3d)
341!! real(r8), intent(inout) :: DiaV3wrk(LBi:UBi,LBj:UBj,N(ng),NDM3d)
342!! real(r8), intent(inout) :: DiaRU(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
343!! real(r8), intent(inout) :: DiaRV(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
344# endif
345 real(r8), intent(inout) :: ad_Hz(LBi:UBi,LBj:UBj,N(ng))
346 real(r8), intent(inout) :: ad_Huon(LBi:UBi,LBj:UBj,N(ng))
347 real(r8), intent(inout) :: ad_Hvom(LBi:UBi,LBj:UBj,N(ng))
348 real(r8), intent(inout) :: ad_z_r(LBi:UBi,LBj:UBj,N(ng))
349 real(r8), intent(inout) :: ad_z_w(LBi:UBi,LBj:UBj,0:N(ng))
350 real(r8), intent(inout) :: ad_btflx(LBi:UBi,LBj:UBj,NT(ng))
351 real(r8), intent(inout) :: ad_bustr(LBi:UBi,LBj:UBj)
352 real(r8), intent(inout) :: ad_bvstr(LBi:UBi,LBj:UBj)
353 real(r8), intent(inout) :: ad_stflx(LBi:UBi,LBj:UBj,NT(ng))
354 real(r8), intent(inout) :: ad_sustr(LBi:UBi,LBj:UBj)
355 real(r8), intent(inout) :: ad_svstr(LBi:UBi,LBj:UBj)
356 real(r8), intent(inout) :: ad_ru(LBi:UBi,LBj:UBj,0:N(ng),2)
357 real(r8), intent(inout) :: ad_rv(LBi:UBi,LBj:UBj,0:N(ng),2)
358# ifdef LMD_NONLOCAL_NOT_YET
359 real(r8), intent(inout) :: ad_ghats(LBi:UBi,LBj:UBj,0:N(ng),NAT)
360# endif
361 real(r8), intent(inout) :: ad_Akt(LBi:UBi,LBj:UBj,0:N(ng),NAT)
362 real(r8), intent(inout) :: ad_Akv(LBi:UBi,LBj:UBj,0:N(ng))
363 real(r8), intent(inout) :: ad_W(LBi:UBi,LBj:UBj,0:N(ng))
364 real(r8), intent(inout) :: ad_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
365 real(r8), intent(inout) :: ad_u(LBi:UBi,LBj:UBj,N(ng),2)
366 real(r8), intent(inout) :: ad_v(LBi:UBi,LBj:UBj,N(ng),2)
367# endif
368!
369! Local variable declarations.
370!
371 integer :: Isrc, Jsrc
372 integer :: i, ic, indx, is, itrc, j, k, ltrc
373# if defined AGE_MEAN && defined T_PASSIVE
374 integer :: iage
375# endif
376# if defined DIAGNOSTICS_TS || defined DIAGNOSTICS_UV
377 integer :: idiag
378# endif
379 real(r8), parameter :: eps = 1.0e-16_r8
380
381 real(r8) :: cff, cff1, cff2, cff3, cff4
382 real(r8) :: ad_cff, ad_cff1, ad_cff2, ad_cff3, ad_cff4
383 real(r8) :: adfac, adfac1, adfac2, adfac3
384 real(r8) :: Gamma
385
386 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: CF
387 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: DC
388 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: FC
389
390 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: ad_CF
391 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: ad_DC
392 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: ad_FC
393
394# ifdef SOLAR_SOURCE
395 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,0:N(ng)) :: ad_swdk
396# endif
397
398 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FE
399 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FX
400 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: curv
401 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: grad
402
403 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_FE
404 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_FX
405 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_curv
406 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_grad
407
408# include "set_bounds.h"
409!
410! In the tangent model, we have:
411!
412!^ ! nnew=3-nstp
413!^ indx=3-nrhs ! nrhs=nstp
414!^
415 indx=nnew
416!
417!-----------------------------------------------------------------------
418! Initialize adjoint private variables.
419!-----------------------------------------------------------------------
420!
421 ad_cff=0.0_r8
422 ad_cff1=0.0_r8
423 ad_cff2=0.0_r8
424 ad_cff3=0.0_r8
425 ad_cff4=0.0_r8
426 DO j=jmins,jmaxs
427 DO i=imins,imaxs
428 ad_fe(i,j)=0.0_r8
429 ad_fx(i,j)=0.0_r8
430 ad_curv(i,j)=0.0_r8
431 ad_grad(i,j)=0.0_r8
432 END DO
433# ifdef SOLAR_SOURCE
434 DO k=0,n(ng)
435 DO i=imins,imaxs
436 ad_swdk(i,j,k)=0.0_r8
437 END DO
438 END DO
439# endif
440 END DO
441 DO k=0,n(ng)
442 DO i=imins,imaxs
443 ad_cf(i,k)=0.0_r8
444 ad_dc(i,k)=0.0_r8
445 ad_fc(i,k)=0.0_r8
446 END DO
447 END DO
448
449# ifndef TS_FIXED
450!
451!=======================================================================
452! Apply tracers lateral boundary conditions.
453!=======================================================================
454!
455# ifdef DISTRIBUTE
456!^ CALL mp_exchange4d (ng, tile, iTLM, 1, &
457!^ & LBi, UBi, LBj, UBj, 1, N(ng), 1, NT(ng), &
458!^ & NghostPoints, &
459!^ & EWperiodic(ng), NSperiodic(ng), &
460!^ & tl_t(:,:,:,3,:))
461!^
462 CALL ad_mp_exchange4d (ng, tile, iadm, 1, &
463 & lbi, ubi, lbj, ubj, 1, n(ng), 1, nt(ng), &
464 & nghostpoints, &
465 & ewperiodic(ng), nsperiodic(ng), &
466 & ad_t(:,:,:,3,:))
467!
468# endif
469!
470 ic=0
471 DO itrc=1,nt(ng)
472 IF (ltracerclm(itrc,ng).and.lnudgetclm(itrc,ng)) THEN
473 ic=ic+1
474 END IF
475 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
476!^ CALL exchange_r3d_tile (ng, tile, &
477!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
478!^ & tl_t(:,:,:,3,itrc))
479!^
480 CALL ad_exchange_r3d_tile (ng, tile, &
481 & lbi, ubi, lbj, ubj, 1, n(ng), &
482 & ad_t(:,:,:,3,itrc))
483 END IF
484!^ CALL tl_t3dbc_tile (ng, tile, itrc, ic, &
485!^ & LBi, UBi, LBj, UBj, N(ng), NT(ng), &
486!^ & IminS, ImaxS, JminS, JmaxS, &
487!^ & nstp, 3, &
488!^ & tl_t)
489!^
490 CALL ad_t3dbc_tile (ng, tile, itrc, ic, &
491 & lbi, ubi, lbj, ubj, n(ng), nt(ng), &
492 & imins, imaxs, jmins, jmaxs, &
493 & nstp, 3, &
494 & ad_t)
495 END DO
496# endif
497!
498!=======================================================================
499! 3D adjoint momentum equation in the ETA-direction.
500!=======================================================================
501!
502 j_loop2 : DO j=jstr,jend
503 IF (j.ge.jstrv) THEN
504!
505! Compute new V-momentum (m m/s).
506!
507 cff=dt(ng)*0.25_r8
508 DO i=istr,iend
509 dc(i,0)=cff*(pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))
510 END DO
511 IF (iic(ng).eq.ntfirst(ng)) THEN
512 DO k=1,n(ng)
513 DO i=istr,iend
514# ifdef DIAGNOSTICS_UV
515!! DiaV3wrk(i,j,k,M3rate)=cff1
516!! DiaV3wrk(i,j,k,M3vvis)=cff2
517!! DO idiag=1,M3pgrd
518!! DiaV3wrk(i,j,k,idiag)=0.0_r8
519!! END DO
520# endif
521!^ tl_v(i,j,k,nnew)=tl_cff1+tl_cff2
522!^
523 ad_cff1=ad_cff1+ad_v(i,j,k,nnew)
524 ad_cff2=ad_cff2+ad_v(i,j,k,nnew)
525 ad_v(i,j,k,nnew)=0.0_r8
526!^ tl_cff2=tl_FC(i,k)-tl_FC(i,k-1)
527!^
528 ad_fc(i,k-1)=ad_fc(i,k-1)-ad_cff2
529 ad_fc(i,k )=ad_fc(i,k )+ad_cff2
530 ad_cff2=0.0_r8
531!^ tl_cff1=0.5_r8*(tl_v(i,j,k,nstp)* &
532!^ & (Hz(i,j,k)+Hz(i,j-1,k))+ &
533!^ & v(i,j,k,nstp)* &
534!^ & (tl_Hz(i,j,k)+tl_Hz(i,j-1,k)))
535!^
536 adfac=0.5_r8*ad_cff1
537 adfac1=adfac*v(i,j,k,nstp)
538 ad_v(i,j,k,nstp)=ad_v(i,j,k,nstp)+ &
539 & (hz(i,j,k)+hz(i,j-1,k))*adfac
540 ad_hz(i,j-1,k)=ad_hz(i,j-1,k)+adfac1
541 ad_hz(i,j ,k)=ad_hz(i,j ,k)+adfac1
542 ad_cff1=0.0_r8
543 END DO
544 END DO
545 ELSE IF (iic(ng).eq.(ntfirst(ng)+1)) THEN
546 DO k=1,n(ng)
547 DO i=istr,iend
548 cff3=0.5_r8*dc(i,0)
549# ifdef DIAGNOSTICS_UV
550!! DiaV3wrk(i,j,k,M3rate)=cff1
551# ifdef BODYFORCE
552!! DiaV3wrk(i,j,k,M3vvis)=DiaV3wrk(i,j,k,M3vvis)- &
553!! & cff3*DiaRV(i,j,k,indx,M3vvis)
554# endif
555!! DiaV3wrk(i,j,k,M3vvis)=cff2
556!! DO idiag=1,M3pgrd
557!! DiaV3wrk(i,j,k,idiag)=-cff3*DiaRV(i,j,k,indx,idiag)
558!! END DO
559# endif
560!^ tl_v(i,j,k,nnew)=tl_cff1- &
561!^ & cff3*tl_rv(i,j,k,indx)+ &
562!^ & tl_cff2
563!
564 ad_rv(i,j,k,indx)=ad_rv(i,j,k,indx)- &
565 & cff3*ad_v(i,j,k,nnew)
566 ad_cff1=ad_cff1+ad_v(i,j,k,nnew)
567 ad_cff2=ad_cff2+ad_v(i,j,k,nnew)
568 ad_v(i,j,k,nnew)=0.0_r8
569!^ tl_cff2=tl_FC(i,k)-tl_FC(i,k-1)
570!^
571 ad_fc(i,k-1)=ad_fc(i,k-1)-ad_cff2
572 ad_fc(i,k )=ad_fc(i,k )+ad_cff2
573 ad_cff2=0.0_r8
574!^ tl_cff1=0.5_r8*(tl_v(i,j,k,nstp)* &
575!^ & (Hz(i,j,k)+Hz(i,j-1,k))+ &
576!^ & v(i,j,k,nstp)* &
577!^ & (tl_Hz(i,j,k)+tl_Hz(i,j-1,k)))
578!^
579 adfac=0.5_r8*ad_cff1
580 adfac1=adfac*v(i,j,k,nstp)
581 ad_hz(i,j-1,k)=ad_hz(i,j-1,k)+adfac1
582 ad_hz(i,j ,k)=ad_hz(i,j ,k)+adfac1
583 ad_v(i,j,k,nstp)=ad_v(i,j,k,nstp)+ &
584 & (hz(i,j,k)+hz(i,j-1,k))*adfac
585 ad_cff1=0.0_r8
586 END DO
587 END DO
588 ELSE
589 cff1= 5.0_r8/12.0_r8
590 cff2=16.0_r8/12.0_r8
591 DO k=1,n(ng)
592 DO i=istr,iend
593# ifdef DIAGNOSTICS_UV
594!! DiaV3wrk(i,j,k,M3rate)=cff3
595# ifdef BODYFORCE
596!! DiaV3wrk(i,j,k,M3vvis)=DiaV3wrk(i,j,k,M3vvis)+ &
597!! & DC(i,0)* &
598!! & (cff1*DiaRV(i,j,k,nrhs,M3vvis)- &
599!! & cff2*DiaRV(i,j,k,indx,M3vvis))
600# endif
601!! DiaV3wrk(i,j,k,M3vvis)=cff4
602!! DO idiag=1,M3pgrd
603!! DiaV3wrk(i,j,k,idiag)=DC(i,0)* &
604!! & (cff1*DiaRV(i,j,k,nrhs,idiag)- &
605!! & cff2*DiaRV(i,j,k,indx,idiag))
606!! END DO
607# endif
608!^ tl_v(i,j,k,nnew)=tl_cff3+ &
609!^ & DC(i,0)*(cff1*tl_rv(i,j,k,nrhs)- &
610!^ & cff2*tl_rv(i,j,k,indx))+ &
611!^ & tl_cff4
612!^
613 adfac=dc(i,0)*ad_v(i,j,k,nnew)
614 ad_rv(i,j,k,nrhs)=ad_rv(i,j,k,nrhs)+cff1*adfac
615 ad_rv(i,j,k,indx)=ad_rv(i,j,k,indx)-cff2*adfac
616 ad_cff3=ad_cff3+ad_v(i,j,k,nnew)
617 ad_cff4=ad_cff4+ad_v(i,j,k,nnew)
618 ad_v(i,j,k,nnew)=0.0_r8
619!^ tl_cff4=tl_FC(i,k)-tl_FC(i,k-1)
620!^
621 ad_fc(i,k-1)=ad_fc(i,k-1)-ad_cff4
622 ad_fc(i,k )=ad_fc(i,k )+ad_cff4
623 ad_cff4=0.0_r8
624!^ tl_cff3=0.5_r8*(tl_v(i,j,k,nstp)* &
625!^ & (Hz(i,j,k)+Hz(i,j-1,k))+ &
626!^ & v(i,j,k,nstp)* &
627!^ & (tl_Hz(i,j,k)+tl_Hz(i,j-1,k)))
628!^
629 adfac=0.5_r8*ad_cff3
630 adfac1=adfac*v(i,j,k,nstp)
631 ad_hz(i,j-1,k)=ad_hz(i,j-1,k)+adfac1
632 ad_hz(i,j ,k)=ad_hz(i,j ,k)+adfac1
633 ad_v(i,j,k,nstp)=ad_v(i,j,k,nstp)+ &
634 & (hz(i,j,k)+hz(i,j-1,k))*adfac
635 ad_cff3=0.0_r8
636 END DO
637 END DO
638 END IF
639!
640! Apply bottom and surface stresses, if so is prescribed.
641!
642 DO i=istr,iend
643# ifdef BODYFORCE
644!^ tl_FC(i,N(ng))=0.0_r8
645!^
646 ad_fc(i,n(ng))=0.0_r8
647!^ tl_FC(i,0)=0.0_r8
648!^
649 ad_fc(i,0)=0.0_r8
650# else
651!^ tl_FC(i,N(ng))=dt(ng)*tl_svstr(i,j)
652!^
653 ad_svstr(i,j)=ad_svstr(i,j)+dt(ng)*ad_fc(i,n(ng))
654 ad_fc(i,n(ng))=0.0_r8
655!^ tl_FC(i,0)=dt(ng)*tl_bvstr(i,j)
656!^
657 ad_bvstr(i,j)=ad_bvstr(i,j)+dt(ng)*ad_fc(i,0)
658 ad_fc(i,0)=0.0_r8
659# endif
660 END DO
661!
662! Compute adjoint V-component viscous vertical momentum fluxes "FC"
663! at current time-step n, and at horizontal V-points and vertical
664! W-points.
665!
666 cff3=dt(ng)*(1.0_r8-lambda)
667 DO k=1,n(ng)-1
668 DO i=istr,iend
669 cff=1.0_r8/(z_r(i,j,k+1)+z_r(i,j-1,k+1)- &
670 & z_r(i,j,k )-z_r(i,j-1,k ))
671!^ tl_FC(i,k)=cff3* &
672!^ & (cff*((tl_v(i,j,k+1,nstp)-tl_v(i,j,k,nstp))* &
673!^ & (Akv(i,j,k)+Akv(i,j-1,k))+ &
674!^ & (v(i,j,k+1,nstp)-v(i,j,k,nstp))* &
675!^ & (tl_Akv(i,j,k)+tl_Akv(i,j-1,k)))+ &
676!^ & tl_cff*(v(i,j,k+1,nstp)-v(i,j,k,nstp))* &
677!^ & (Akv(i,j,k)+Akv(i,j-1,k)))
678!^
679 adfac=cff3*ad_fc(i,k)
680 adfac1=adfac*cff
681 adfac2=adfac1*(akv(i,j,k)+akv(i,j-1,k))
682 adfac3=adfac1*(v(i,j,k+1,nstp)-v(i,j,k,nstp))
683 ad_v(i,j,k ,nstp)=ad_v(i,j,k ,nstp)-adfac2
684 ad_v(i,j,k+1,nstp)=ad_v(i,j,k+1,nstp)+adfac2
685 ad_akv(i,j-1,k)=ad_akv(i,j-1,k)+adfac3
686 ad_akv(i,j ,k)=ad_akv(i,j ,k)+adfac3
687 ad_cff=ad_cff+ &
688 & (v(i,j,k+1,nstp)-v(i,j,k,nstp))* &
689 & (akv(i,j,k)+akv(i,j-1,k))*adfac
690 ad_fc(i,k)=0.0_r8
691!^ tl_cff=-cff*cff*(tl_z_r(i,j,k+1)+tl_z_r(i,j-1,k+1)- &
692!^ & tl_z_r(i,j,k )-tl_z_r(i,j-1,k ))
693!^
694 adfac=-cff*cff*ad_cff
695 ad_z_r(i,j-1,k )=ad_z_r(i,j-1,k )-adfac
696 ad_z_r(i,j ,k )=ad_z_r(i,j ,k )-adfac
697 ad_z_r(i,j-1,k+1)=ad_z_r(i,j-1,k+1)+adfac
698 ad_z_r(i,j ,k+1)=ad_z_r(i,j ,k+1)+adfac
699 ad_cff=0.0_r8
700 END DO
701 END DO
702 END IF
703!
704!=====================================================================
705! 3D adjoint momentum equation in the XI-direction.
706!=====================================================================
707!
708! Compute new U-momentum (m m/s).
709!
710 cff=dt(ng)*0.25_r8
711 DO i=istru,iend
712 dc(i,0)=cff*(pm(i,j)+pm(i-1,j))*(pn(i,j)+pn(i-1,j))
713 END DO
714 IF (iic(ng).eq.ntfirst(ng)) THEN
715 DO k=1,n(ng)
716 DO i=istru,iend
717# ifdef DIAGNOSTICS_UV
718!! DiaU3wrk(i,j,k,M3rate)=cff1
719!! DiaU3wrk(i,j,k,M3vvis)=cff2
720!! DO idiag=1,M3pgrd
721!! DiaU3wrk(i,j,k,idiag)=0.0_r8
722!! END DO
723# endif
724!^ tl_u(i,j,k,nnew)=tl_cff1+tl_cff2
725!^
726 ad_cff1=ad_cff1+ad_u(i,j,k,nnew)
727 ad_cff2=ad_cff2+ad_u(i,j,k,nnew)
728 ad_u(i,j,k,nnew)=0.0_r8
729!^ tl_cff2=tl_FC(i,k)-tl_FC(i,k-1)
730!^
731 ad_fc(i,k-1)=ad_fc(i,k-1)-ad_cff2
732 ad_fc(i,k )=ad_fc(i,k )+ad_cff2
733 ad_cff2=0.0_r8
734!^ tl_cff1=0.5_r8*(tl_u(i,j,k,nstp)* &
735!^ & (Hz(i,j,k)+Hz(i-1,j,k))+ &
736!^ & u(i,j,k,nstp)* &
737!^ & (tl_Hz(i,j,k)+tl_Hz(i-1,j,k)))
738!^
739 adfac=0.5_r8*ad_cff1
740 adfac1=adfac*u(i,j,k,nstp)
741 ad_hz(i-1,j,k)=ad_hz(i-1,j,k)+adfac1
742 ad_hz(i ,j,k)=ad_hz(i ,j,k)+adfac1
743 ad_u(i,j,k,nstp)=ad_u(i,j,k,nstp)+ &
744 & (hz(i,j,k)+hz(i-1,j,k))*adfac
745 ad_cff1=0.0_r8
746 END DO
747 END DO
748 ELSE IF (iic(ng).eq.(ntfirst(ng)+1)) THEN
749 DO k=1,n(ng)
750 DO i=istru,iend
751 cff3=0.5_r8*dc(i,0)
752# ifdef DIAGNOSTICS_UV
753!! DiaU3wrk(i,j,k,M3rate)=cff1
754# ifdef BODYFORCE
755!! DiaU3wrk(i,j,k,M3vvis)=DiaU3wrk(i,j,k,M3vvis)- &
756!! & cff3*DiaRU(i,j,k,indx,M3vvis)
757# endif
758!! DO idiag=1,M3pgrd
759!! DiaU3wrk(i,j,k,idiag)=-cff3*DiaRU(i,j,k,indx,idiag)
760!! END DO
761!! DiaU3wrk(i,j,k,M3vvis)=cff2
762# endif
763!^ tl_u(i,j,k,nnew)=tl_cff1- &
764!^ & cff3*tl_ru(i,j,k,indx)+ &
765!^ & tl_cff2
766!^
767 ad_ru(i,j,k,indx)=ad_ru(i,j,k,indx)- &
768 & cff3*ad_u(i,j,k,nnew)
769 ad_cff1=ad_cff1+ad_u(i,j,k,nnew)
770 ad_cff2=ad_cff2+ad_u(i,j,k,nnew)
771 ad_u(i,j,k,nnew)=0.0_r8
772!^ tl_cff2=tl_FC(i,k)-tl_FC(i,k-1)
773!^
774 ad_fc(i,k-1)=ad_fc(i,k-1)-ad_cff2
775 ad_fc(i,k )=ad_fc(i,k )+ad_cff2
776 ad_cff2=0.0_r8
777!^ tl_cff1=0.5_r8*(tl_u(i,j,k,nstp)* &
778!^ & (Hz(i,j,k)+Hz(i-1,j,k))+ &
779!^ & u(i,j,k,nstp)* &
780!^ & (tl_Hz(i,j,k)+tl_Hz(i-1,j,k)))
781!^
782 adfac=0.5_r8*ad_cff1
783 adfac1=adfac*u(i,j,k,nstp)
784 ad_hz(i-1,j,k)=ad_hz(i-1,j,k)+adfac1
785 ad_hz(i ,j,k)=ad_hz(i ,j,k)+adfac1
786 ad_u(i,j,k,nstp)=ad_u(i,j,k,nstp)+ &
787 & (hz(i,j,k)+hz(i-1,j,k))*adfac
788 ad_cff1=0.0_r8
789 END DO
790 END DO
791 ELSE
792 cff1= 5.0_r8/12.0_r8
793 cff2=16.0_r8/12.0_r8
794 DO k=1,n(ng)
795 DO i=istru,iend
796# ifdef DIAGNOSTICS_UV
797!! DiaU3wrk(i,j,k,M3rate)=cff3
798# ifdef BODYFORCE
799!! DiaU3wrk(i,j,k,M3vvis)=DiaU3wrk(i,j,k,M3vvis)+ &
800!! & DC(i,0)* &
801!! & (cff1*DiaRU(i,j,k,nrhs,M3vvis)- &
802!! & cff2*DiaRU(i,j,k,indx,M3vvis))
803# endif
804!! DiaU3wrk(i,j,k,M3vvis)=cff4
805!! DO idiag=1,M3pgrd
806!! DiaU3wrk(i,j,k,idiag)=DC(i,0)* &
807!! & (cff1*DiaRU(i,j,k,nrhs,idiag)- &
808!! & cff2*DiaRU(i,j,k,indx,idiag))
809!! END DO
810# endif
811!^ tl_u(i,j,k,nnew)=tl_cff3+ &
812!^ & DC(i,0)*(cff1*tl_ru(i,j,k,nrhs)- &
813!^ & cff2*tl_ru(i,j,k,indx))+ &
814!^ & tl_cff4
815!^
816 adfac=dc(i,0)*ad_u(i,j,k,nnew)
817 ad_ru(i,j,k,nrhs)=ad_ru(i,j,k,nrhs)+cff1*adfac
818 ad_ru(i,j,k,indx)=ad_ru(i,j,k,indx)-cff2*adfac
819 ad_cff3=ad_cff3+ad_u(i,j,k,nnew)
820 ad_cff4=ad_cff4+ad_u(i,j,k,nnew)
821 ad_u(i,j,k,nnew)=0.0_r8
822!^ tl_cff4=tl_FC(i,k)-tl_FC(i,k-1)
823!^
824 ad_fc(i,k-1)=ad_fc(i,k-1)-ad_cff4
825 ad_fc(i,k )=ad_fc(i,k )+ad_cff4
826 ad_cff4=0.0_r8
827!^ tl_cff3=0.5_r8*(tl_u(i,j,k,nstp)* &
828!^ & (Hz(i,j,k)+Hz(i-1,j,k))+ &
829!^ & u(i,j,k,nstp)* &
830!^ & (tl_Hz(i,j,k)+tl_Hz(i-1,j,k)))
831!^
832 adfac=0.5_r8*ad_cff3
833 adfac1=adfac*u(i,j,k,nstp)
834 ad_hz(i-1,j,k)=ad_hz(i-1,j,k)+adfac1
835 ad_hz(i ,j,k)=ad_hz(i ,j,k)+adfac1
836 ad_u(i,j,k,nstp)=ad_u(i,j,k,nstp)+ &
837 & (hz(i,j,k)+hz(i-1,j,k))*adfac
838 ad_cff3=0.0_r8
839 END DO
840 END DO
841 END IF
842!
843! Apply bottom and surface stresses, if so is prescribed.
844!
845 DO i=istru,iend
846# ifdef BODYFORCE
847!^ tl_FC(i,N(ng))=0.0_r8
848!^
849 ad_fc(i,n(ng))=0.0_r8
850!^ tl_FC(i,0)=0.0_r8
851!^
852 ad_fc(i,0)=0.0_r8
853# else
854!^ tl_FC(i,N(ng))=dt(ng)*tl_sustr(i,j)
855!^
856 ad_sustr(i,j)=ad_sustr(i,j)+dt(ng)*ad_fc(i,n(ng))
857 ad_fc(i,n(ng))=0.0_r8
858!^ tl_FC(i,0)=dt(ng)*tl_bustr(i,j)
859!^
860 ad_bustr(i,j)=ad_bustr(i,j)+dt(ng)*ad_fc(i,0)
861 ad_fc(i,0)=0.0_r8
862# endif
863 END DO
864!
865! Compute adjoint U-component viscous vertical momentum fluxes "FC"
866! at current time-step n, and at horizontal U-points and vertical
867! W-points.
868!
869 cff3=dt(ng)*(1.0_r8-lambda)
870 DO k=1,n(ng)-1
871 DO i=istru,iend
872 cff=1.0/(z_r(i,j,k+1)+z_r(i-1,j,k+1)- &
873 & z_r(i,j,k )-z_r(i-1,j,k ))
874!^ tl_FC(i,k)=cff3* &
875!^ & (cff*((tl_u(i,j,k+1,nstp)-tl_u(i,j,k,nstp))* &
876!^ & (Akv(i,j,k)+Akv(i-1,j,k))+ &
877!^ & (u(i,j,k+1,nstp)-u(i,j,k,nstp))* &
878!^ & (tl_Akv(i,j,k)+tl_Akv(i-1,j,k)))+ &
879!^ & tl_cff*(u(i,j,k+1,nstp)-u(i,j,k,nstp))* &
880!^ & (Akv(i,j,k)+Akv(i-1,j,k)))
881!^
882 adfac=cff3*ad_fc(i,k)
883 adfac1=adfac*cff
884 adfac2=adfac1*(akv(i,j,k)+akv(i-1,j,k))
885 adfac3=adfac1*(u(i,j,k+1,nstp)-u(i,j,k,nstp))
886 ad_u(i,j,k ,nstp)=ad_u(i,j,k ,nstp)-adfac2
887 ad_u(i,j,k+1,nstp)=ad_u(i,j,k+1,nstp)+adfac2
888 ad_akv(i-1,j,k)=ad_akv(i-1,j,k)+adfac3
889 ad_akv(i ,j,k)=ad_akv(i ,j,k)+adfac3
890 ad_cff=ad_cff+ &
891 & (u(i,j,k+1,nstp)-u(i,j,k,nstp))* &
892 & (akv(i,j,k)+akv(i-1,j,k))*adfac
893 ad_fc(i,k)=0.0_r8
894!^ tl_cff=-cff*cff*(tl_z_r(i,j,k+1)+tl_z_r(i-1,j,k+1)- &
895!^ & tl_z_r(i,j,k )-tl_z_r(i-1,j,k ))
896!^
897 adfac=-cff*cff*ad_cff
898 ad_z_r(i-1,j,k )=ad_z_r(i-1,j,k )-adfac
899 ad_z_r(i ,j,k )=ad_z_r(i ,j,k )-adfac
900 ad_z_r(i-1,j,k+1)=ad_z_r(i-1,j,k+1)+adfac
901 ad_z_r(i ,j,k+1)=ad_z_r(i ,j,k+1)+adfac
902 ad_cff=0.0_r8
903 END DO
904 END DO
905 END DO j_loop2
906
907# ifndef TS_FIXED
908!
909!=======================================================================
910! Adjoint tracer equation(s).
911!=======================================================================
912!
913!-----------------------------------------------------------------------
914! Start computation of tracers at n+1 time-step, t(i,j,k,nnew,itrc).
915!-----------------------------------------------------------------------
916!
917! Compute vertical diffusive fluxes "FC" of the tracer fields at
918! current time step n, and at horizontal RHO-points and vertical
919! W-points.
920!
921 DO j=jstr,jend
922 DO itrc=1,nt(ng)
923!
924! Compute new tracer field (m Tunits).
925!
926 DO k=1,n(ng)
927 DO i=istr,iend
928# ifdef DIAGNOSTICS_TS
929!! DiaTwrk(i,j,k,itrc,iTvdif)=cff2
930!! DiaTwrk(i,j,k,itrc,iTrate)=cff1
931# endif
932!^ tl_t(i,j,k,nnew,itrc)=tl_cff1+tl_cff2
933!^
934 ad_cff1=ad_cff1+ad_t(i,j,k,nnew,itrc)
935 ad_cff2=ad_cff2+ad_t(i,j,k,nnew,itrc)
936 ad_t(i,j,k,nnew,itrc)=0.0_r8
937!^ tl_cff2=tl_FC(i,k)-tl_FC(i,k-1)
938!^
939 ad_fc(i,k-1)=ad_fc(i,k-1)-ad_cff2
940 ad_fc(i,k )=ad_fc(i,k )+ad_cff2
941 ad_cff2=0.0_r8
942!^ tl_cff1=tl_Hz(i,j,k)*t(i,j,k,nstp,itrc)+ &
943!^ & Hz(i,j,k)*tl_t(i,j,k,nstp,itrc)
944!^
945 ad_t(i,j,k,nstp,itrc)=ad_t(i,j,k,nstp,itrc)+ &
946 & hz(i,j,k)*ad_cff1
947 ad_hz(i,j,k)=ad_hz(i,j,k)+t(i,j,k,nstp,itrc)*ad_cff1
948 ad_cff1=0.0_r8
949 END DO
950 END DO
951!
952! Apply bottom and surface tracer flux conditions.
953!
954 DO i=istr,iend
955!^ tl_FC(i,N(ng))=dt(ng)*tl_stflx(i,j,itrc)
956!^
957 ad_stflx(i,j,itrc)=ad_stflx(i,j,itrc)+dt(ng)*ad_fc(i,n(ng))
958 ad_fc(i,n(ng))=0.0_r8
959!^ tl_FC(i,0)=dt(ng)*tl_btflx(i,j,itrc)
960!^
961 ad_btflx(i,j,itrc)=ad_btflx(i,j,itrc)+dt(ng)*ad_fc(i,0)
962 ad_fc(i,0)=0.0_r8
963 END DO
964
965# ifdef SOLAR_SOURCE
966!
967! Add in incoming solar radiation at interior W-points using decay
968! decay penetration function based on Jerlov water type.
969!
970 IF (itrc.eq.itemp) THEN
971 DO k=1,n(ng)-1
972 DO i=istr,iend
973!^ tl_FC(i,k)=tl_FC(i,k)+ &
974!^ & dt(ng)*srflx(i,j)* &
975# ifdef WET_DRY_NOT_YET
976!^ & rmask_wet(i,j)* &
977# endif
978!^ & tl_swdk(i,j,k)
979!^
980 ad_swdk(i,j,k)=ad_swdk(i,j,k)+ &
981 & dt(ng)*srflx(i,j)* &
982# ifdef WET_DRY_NOT_YET
983 & rmask_wet(i,j)* &
984# endif
985 & ad_fc(i,k)
986 END DO
987 END DO
988 END IF
989# endif
990# ifdef LMD_NONLOCAL_NOT_YET
991!
992! Add in the nonlocal transport flux for unstable (convective)
993! forcing conditions into matrix FC when using the Large et al.
994! KPP scheme. The nonlocal transport is only applied to active
995! tracers.
996!
997 IF (itrc.le.nat) THEN
998 DO k=1,n(ng)-1
999 DO i=istr,iend
1000!^ tl_FC(i,k)=tl_FC(i,k)- &
1001!^ & dt(ng)*(tl_Akt(i,j,k,itrc)* &
1002!^ & ghats(i,j,k,itrc)+ &
1003!^ & Akt(i,j,k,itrc)* &
1004!^ & tl_ghats(i,j,k,itrc))
1005!^
1006 adfac=dt(ng)*ad_fc(i,k)
1007 ad_ghats(i,j,k,itrc)=ad_ghats(i,j,k,itrc)- &
1008 & akt(i,j,k,itrc)*adfac
1009 ad_akt(i,j,k,itrc)=ad_akt(i,j,k,itrc)- &
1010 & ghats(i,j,k,itrc)*adfac
1011 END DO
1012 END DO
1013 END IF
1014# endif
1015!
1016! Compute adjoint vertical diffusive fluxes "FC" of the tracer fields
1017! at current time step n, and at horizontal RHO-points and vertical
1018! W-points. Notice that the vertical diffusion coefficients for
1019! passive tracers is the same as that for salinity (ltrc=NAT).
1020!
1021 cff3=dt(ng)*(1.0_r8-lambda)
1022 ltrc=min(nat,itrc)
1023 DO k=1,n(ng)-1
1024 DO i=istr,iend
1025 cff=1.0_r8/(z_r(i,j,k+1)-z_r(i,j,k))
1026!^ tl_FC(i,k)=cff3* &
1027!^ & (cff*(tl_Akt(i,j,k,ltrc)* &
1028!^ & (t(i,j,k+1,nstp,itrc)- &
1029!^ & t(i,j,k ,nstp,itrc))+ &
1030!^ & Akt(i,j,k,ltrc)* &
1031!^ & (tl_t(i,j,k+1,nstp,itrc)- &
1032!^ & tl_t(i,j,k ,nstp,itrc)))+ &
1033!^ & tl_cff*(Akt(i,j,k,ltrc)* &
1034!^ & (t(i,j,k+1,nstp,itrc)- &
1035!^ & t(i,j,k,nstp,itrc))))
1036!^
1037 adfac=cff3*ad_fc(i,k)
1038 adfac1=adfac*cff
1039 adfac2=adfac1*akt(i,j,k,ltrc)
1040 ad_akt(i,j,k,ltrc)=ad_akt(i,j,k,ltrc)+ &
1041 & (t(i,j,k+1,nstp,itrc)- &
1042 & t(i,j,k ,nstp,itrc))*adfac1
1043 ad_t(i,j,k ,nstp,itrc)=ad_t(i,j,k ,nstp,itrc)-adfac2
1044 ad_t(i,j,k+1,nstp,itrc)=ad_t(i,j,k+1,nstp,itrc)+adfac2
1045 ad_cff=ad_cff+ &
1046 & (akt(i,j,k,ltrc)* &
1047 & (t(i,j,k+1,nstp,itrc)- &
1048 & t(i,j,k ,nstp,itrc)))*adfac
1049 ad_fc(i,k)=0.0_r8
1050!^ tl_cff=-cff*cff*(tl_z_r(i,j,k+1)-tl_z_r(i,j,k))
1051!^
1052 adfac=-cff*cff*ad_cff
1053 ad_z_r(i,j,k )=ad_z_r(i,j,k )-adfac
1054 ad_z_r(i,j,k+1)=ad_z_r(i,j,k+1)+adfac
1055 ad_cff=0.0_r8
1056 END DO
1057 END DO
1058 END DO
1059 END DO
1060!
1061!------------------------------------------------------------------------
1062! Time-step adjoint vertical advection of tracers. Impose artificial
1063! continuity equation.
1064!------------------------------------------------------------------------
1065!
1066 j_loop1 : DO j=jstr,jend
1067 t_loop2 : DO itrc=1,nt(ng)
1068!
1069 vadv_flux1 : IF (ad_vadvection(itrc,ng)%SPLINES) THEN
1070!
1071! Build BASIC conservative parabolic splines for the vertical
1072! derivatives "FC" of the tracer. Then, the interfacial "FC" values
1073! are converted to vertical advective flux.
1074!
1075 DO i=istr,iend
1076# ifdef NEUMANN
1077 fc(i,0)=1.5_r8*t(i,j,1,nstp,itrc)
1078 cf(i,1)=0.5_r8
1079# else
1080 fc(i,0)=2.0_r8*t(i,j,1,nstp,itrc)
1081 cf(i,1)=1.0_r8
1082# endif
1083 END DO
1084 DO k=1,n(ng)-1
1085 DO i=istr,iend
1086 cff=1.0_r8/(2.0_r8*hz(i,j,k)+ &
1087 & hz(i,j,k+1)*(2.0_r8-cf(i,k)))
1088 cf(i,k+1)=cff*hz(i,j,k)
1089 fc(i,k)=cff*(3.0_r8*(hz(i,j,k )*t(i,j,k+1,nstp,itrc)+ &
1090 & hz(i,j,k+1)*t(i,j,k ,nstp,itrc))- &
1091 & hz(i,j,k+1)*fc(i,k-1))
1092 END DO
1093 END DO
1094 DO i=istr,iend
1095# ifdef NEUMANN
1096 fc(i,n(ng))=(3.0_r8*t(i,j,n(ng),nstp,itrc)- &
1097 & fc(i,n(ng)-1))/(2.0_r8-cf(i,n(ng)))
1098# else
1099 fc(i,n(ng))=(2.0_r8*t(i,j,n(ng),nstp,itrc)- &
1100 & fc(i,n(ng)-1))/(1.0_r8-cf(i,n(ng)))
1101# endif
1102 END DO
1103 DO k=n(ng)-1,0,-1
1104 DO i=istr,iend
1105 fc(i,k)=fc(i,k)-cf(i,k+1)*fc(i,k+1)
1106 fc(i,k+1)=w(i,j,k+1)*fc(i,k+1)
1107 END DO
1108 END DO
1109 DO i=istr,iend
1110 fc(i,n(ng))=0.0_r8
1111 fc(i,0)=0.0_r8
1112 END DO
1113!
1114 ELSE IF (ad_vadvection(itrc,ng)%AKIMA4) THEN
1115!
1116! Fourth-order, BASIC STATE Akima vertical advective flux.
1117!
1118 DO k=1,n(ng)-1
1119 DO i=istr,iend
1120 fc(i,k)=t(i,j,k+1,nstp,itrc)- &
1121 & t(i,j,k ,nstp,itrc)
1122 END DO
1123 END DO
1124 DO i=istr,iend
1125 fc(i,0)=fc(i,1)
1126 fc(i,n(ng))=fc(i,n(ng)-1)
1127 END DO
1128 DO k=1,n(ng)
1129 DO i=istr,iend
1130 cff=2.0_r8*fc(i,k)*fc(i,k-1)
1131 IF (cff.gt.eps) THEN
1132 cf(i,k)=cff/(fc(i,k)+fc(i,k-1))
1133 ELSE
1134 cf(i,k)=0.0_r8
1135 END IF
1136 END DO
1137 END DO
1138 cff1=1.0_r8/3.0_r8
1139 DO k=1,n(ng)-1
1140 DO i=istr,iend
1141 fc(i,k)=w(i,j,k)* &
1142 & 0.5_r8*(t(i,j,k ,nstp,itrc)+ &
1143 & t(i,j,k+1,nstp,itrc)- &
1144 & cff1*(cf(i,k+1)-cf(i,k)))
1145 END DO
1146 END DO
1147 DO i=istr,iend
1148 fc(i,0)=0.0_r8
1149 fc(i,n(ng))=0.0_r8
1150 END DO
1151!
1152 ELSE IF (ad_vadvection(itrc,ng)%CENTERED2) THEN
1153!
1154! Second-order, BASIC STATE central differences vertical advective
1155! flux.
1156!
1157 DO k=1,n(ng)-1
1158 DO i=istr,iend
1159 fc(i,k)=w(i,j,k)* &
1160 & 0.5_r8*(t(i,j,k ,nstp,itrc)+ &
1161 & t(i,j,k+1,nstp,itrc))
1162 END DO
1163 END DO
1164 DO i=istr,iend
1165 fc(i,0)=0.0_r8
1166 fc(i,n(ng))=0.0_r8
1167 END DO
1168!
1169 ELSE IF ((ad_vadvection(itrc,ng)%MPDATA).or. &
1170 & (ad_vadvection(itrc,ng)%HSIMT)) THEN
1171!
1172! First_order, BASIC STATE upstream differences vertical advective
1173! flux.
1174!
1175 DO i=istr,iend
1176 DO k=1,n(ng)-1
1177 cff1=max(w(i,j,k),0.0_r8)
1178 cff2=min(w(i,j,k),0.0_r8)
1179 fc(i,k)=cff1*t(i,j,k ,nstp,itrc)+ &
1180 & cff2*t(i,j,k+1,nstp,itrc)
1181 END DO
1182 fc(i,0)=0.0_r8
1183 fc(i,n(ng))=0.0_r8
1184 END DO
1185!
1186 ELSE IF ((ad_vadvection(itrc,ng)%CENTERED4).or. &
1187 & (ad_vadvection(itrc,ng)%SPLIT_U3)) THEN
1188!
1189! Fourth-order, BASIC STATE central differences vertical advective
1190! flux.
1191!
1192 cff1=0.5_r8
1193 cff2=7.0_r8/12.0_r8
1194 cff3=1.0_r8/12.0_r8
1195 DO k=2,n(ng)-2
1196 DO i=istr,iend
1197 fc(i,k)=w(i,j,k)* &
1198 & (cff2*(t(i,j,k ,nstp,itrc)+ &
1199 & t(i,j,k+1,nstp,itrc))- &
1200 & cff3*(t(i,j,k-1,nstp,itrc)+ &
1201 & t(i,j,k+2,nstp,itrc)))
1202 END DO
1203 END DO
1204 DO i=istr,iend
1205 fc(i,0)=0.0_r8
1206 fc(i,1)=w(i,j,1)* &
1207 & (cff1*t(i,j,1,nstp,itrc)+ &
1208 & cff2*t(i,j,2,nstp,itrc)- &
1209 & cff3*t(i,j,3,nstp,itrc))
1210 fc(i,n(ng)-1)=w(i,j,n(ng)-1)* &
1211 & (cff1*t(i,j,n(ng) ,nstp,itrc)+ &
1212 & cff2*t(i,j,n(ng)-1,nstp,itrc)- &
1213 & cff3*t(i,j,n(ng)-2,nstp,itrc))
1214 fc(i,n(ng))=0.0_r8
1215 END DO
1216 END IF vadv_flux1
1217!
1218! Compute BASIC STATE artificial continuity equation and load it
1219! into private array DC (1/m). It is needed to preserve tracer
1220! constancy. It is not the same for all the tracer advection
1221! schemes because of the Gamma value.
1222!
1223 IF ((vadvection(itrc,ng)%MPDATA).or. &
1224 & (vadvection(itrc,ng)%HSIMT)) THEN
1225 gamma=0.5_r8
1226 ELSE
1227 gamma=1.0_r8/6.0_r8
1228 END IF
1229 IF (iic(ng).eq.ntfirst(ng)) THEN
1230 cff=0.5_r8*dt(ng)
1231 ELSE
1232 cff=(1.0_r8-gamma)*dt(ng)
1233 END IF
1234 DO k=1,n(ng)
1235 DO i=istr,iend
1236 dc(i,k)=1.0_r8/(hz(i,j,k)- &
1237 & cff*pm(i,j)*pn(i,j)* &
1238 & (huon(i+1,j,k)-huon(i,j,k)+ &
1239 & hvom(i,j+1,k)-hvom(i,j,k)+ &
1240 & (w(i,j,k)-w(i,j,k-1))))
1241 END DO
1242 END DO
1243!
1244! Adjoint of Time-step vertical advection of tracers (Tunits).
1245!
1246! WARNING: t(:,:,:,3,itrc) at this point should be in units of
1247! ======= m Tunits. But, t(:,:,:,3,itrc) is read from a BASIC
1248! STATE file and is in Tunits, so we need to multiply
1249! by level thickness (Hz).
1250!
1251 DO k=1,n(ng)
1252 DO i=istr,iend
1253 cff1=cff*pm(i,j)*pn(i,j)
1254!^ tl_t(i,j,k,3,itrc)=tl_DC(i,k)* &
1255!^ & (t(i,j,k,3,itrc)*Hz(i,j,k)- &
1256!^ & cff1*(FC(i,k)-FC(i,k-1)))+ &
1257!^ & DC(i,k)* &
1258!^ & (tl_t(i,j,k,3,itrc)- &
1259!^ & cff1*(tl_FC(i,k)-tl_FC(i,k-1)))
1260!^
1261 adfac=dc(i,k)*ad_t(i,j,k,3,itrc)
1262 adfac1=adfac*cff1
1263 ad_dc(i,k)=ad_dc(i,k)+ &
1264 & (t(i,j,k,3,itrc)*hz(i,j,k)- &
1265 & cff1*(fc(i,k)-fc(i,k-1)))* &
1266 & ad_t(i,j,k,3,itrc)
1267 ad_fc(i,k-1)=ad_fc(i,k-1)+adfac1
1268 ad_fc(i,k )=ad_fc(i,k )-adfac1
1269 ad_t(i,j,k,3,itrc)=adfac
1270 END DO
1271 END DO
1272!
1273! Compute adjoint artificial continuity equation (same for all tracers)
1274! and load it into private array DC (1/m). Notice that "cff" has not
1275! changed from above.
1276!
1277 DO k=1,n(ng)
1278 DO i=istr,iend
1279!^ tl_DC(i,k)=-DC(i,k)*DC(i,k)* &
1280!^ & (tl_Hz(i,j,k)- &
1281!^ & cff*pm(i,j)*pn(i,j)* &
1282!^ & (tl_Huon(i+1,j,k)-tl_Huon(i,j,k)+ &
1283!^ & tl_Hvom(i,j+1,k)-tl_Hvom(i,j,k)+ &
1284!^ & (tl_W(i,j,k)-tl_W(i,j,k-1))))
1285!^
1286 adfac=-dc(i,k)*dc(i,k)*ad_dc(i,k)
1287 adfac1=adfac*cff*pm(i,j)*pn(i,j)
1288 ad_hz(i,j,k)=ad_hz(i,j,k)+adfac
1289 ad_huon(i ,j,k)=ad_huon(i ,j,k)+adfac1
1290 ad_huon(i+1,j,k)=ad_huon(i+1,j,k)-adfac1
1291 ad_hvom(i,j ,k)=ad_hvom(i,j ,k)+adfac1
1292 ad_hvom(i,j+1,k)=ad_hvom(i,j+1,k)-adfac1
1293 ad_w(i,j,k-1)=ad_w(i,j,k-1)+adfac1
1294 ad_w(i,j,k )=ad_w(i,j,k )-adfac1
1295 ad_dc(i,k)=0.0_r8
1296 END DO
1297 END DO
1298!
1299!-----------------------------------------------------------------------
1300! Compute adjoint time rate of change of intermediate tracer due to
1301! vertical advection.
1302!-----------------------------------------------------------------------
1303!
1304 vadv_flux2 : IF (ad_vadvection(itrc,ng)%SPLINES) THEN
1305!
1306! Construct BASIC STATE conservative parabolic splines for the vertical
1307! derivatives "CF" of the tracer.
1308!
1309 DO i=istr,iend
1310# ifdef NEUMANN
1311 fc(i,0)=1.5_r8*t(i,j,1,nstp,itrc)
1312 cf(i,1)=0.5_r8
1313# else
1314 fc(i,0)=2.0_r8*t(i,j,1,nstp,itrc)
1315 cf(i,1)=1.0_r8
1316# endif
1317 END DO
1318 DO k=1,n(ng)-1
1319 DO i=istr,iend
1320 cff=1.0_r8/(2.0_r8*hz(i,j,k)+ &
1321 & hz(i,j,k+1)*(2.0_r8-cf(i,k)))
1322 cf(i,k+1)=cff*hz(i,j,k)
1323 fc(i,k)=cff*(3.0_r8*(hz(i,j,k )*t(i,j,k+1,nstp,itrc)+ &
1324 & hz(i,j,k+1)*t(i,j,k ,nstp,itrc))- &
1325 & hz(i,j,k+1)*fc(i,k-1))
1326 END DO
1327 END DO
1328 DO i=istr,iend
1329# ifdef NEUMANN
1330 fc(i,n(ng))=(3.0_r8*t(i,j,n(ng),nstp,itrc)- &
1331 & fc(i,n(ng)-1))/(2.0_r8-cf(i,n(ng)))
1332# else
1333 fc(i,n(ng))=(2.0_r8*t(i,j,n(ng),nstp,itrc)- &
1334 & fc(i,n(ng)-1))/(1.0_r8-cf(i,n(ng)))
1335# endif
1336 END DO
1337 DO k=n(ng)-1,0,-1
1338 DO i=istr,iend
1339 fc(i,k)=fc(i,k)-cf(i,k+1)*fc(i,k+1)
1340 END DO
1341 END DO
1342!
1343! Now the adjoint spline code.
1344!
1345 DO i=istr,iend
1346!^ tl_FC(i,N(ng))=0.0_r8
1347!^
1348 ad_fc(i,n(ng))=0.0_r8
1349!^ tl_FC(i,0)=0.0_r8
1350!^
1351 ad_fc(i,0)=0.0_r8
1352 END DO
1353!
1354! Adjoint back substitution.
1355!
1356 DO k=0,n(ng)-1
1357 DO i=istr,iend
1358!^ tl_FC(i,k+1)=tl_W(i,j,k+1)*FC(i,k+1)+ &
1359!^ & W(i,j,k+1)*tl_FC(i,k+1)
1360!^
1361 ad_w(i,j,k+1)=ad_w(i,j,k+1)+fc(i,k+1)*ad_fc(i,k+1)
1362 ad_fc(i,k+1)=w(i,j,k+1)*ad_fc(i,k+1)
1363!^ tl_FC(i,k)=tl_FC(i,k)-CF(i,k+1)*tl_FC(i,k+1)
1364!^
1365 ad_fc(i,k+1)=ad_fc(i,k+1)-cf(i,k+1)*ad_fc(i,k)
1366 END DO
1367 END DO
1368!
1369 DO i=istr,iend
1370# ifdef NEUMANN
1371!^ tl_FC(i,N(ng))=(3.0_r8*tl_t(i,j,N(ng),nstp,itrc)- &
1372!^ & tl_FC(i,N(ng)-1))/ &
1373!^ & (2.0_r8-CF(i,N(ng)))
1374!^
1375 adfac=ad_fc(i,n(ng))/(2.0_r8-cf(i,n(ng)))
1376 ad_t(i,j,n(ng),nstp,itrc)=ad_t(i,j,n(ng),nstp,itrc)+ &
1377 & 3.0_r8*adfac
1378 ad_fc(i,n(ng)-1)=ad_fc(i,n(ng)-1)-adfac
1379 ad_fc(i,n(ng))=0.0_r8
1380# else
1381!^ tl_FC(i,N(ng))=(2.0_r8*tl_t(i,j,N(ng),nstp,itrc)- &
1382!^ & tl_FC(i,N(ng)-1))/ &
1383!^ & (1.0_r8-CF(i,N(ng)))
1384!^
1385 adfac=ad_fc(i,n(ng))/(1.0_r8-cf(i,n(ng)))
1386 ad_t(i,j,n(ng),nstp,itrc)=ad_t(i,j,n(ng),nstp,itrc)+ &
1387 & 2.0_r8*adfac
1388 ad_fc(i,n(ng)-1))=ad_fc(i,n(ng)-1))-adfac
1389 ad_fc(i,n(ng))=0.0_r8
1390# endif
1391 END DO
1392!
1393 DO k=n(ng)-1,1,-1
1394 DO i=istr,iend
1395 cff=1.0_r8/(2.0_r8*hz(i,j,k)+ &
1396 & hz(i,j,k+1)*(2.0_r8-cf(i,k)))
1397!^ tl_FC(i,k)=cff* &
1398!^ & (3.0_r8*(Hz(i,j,k )* &
1399!^ & tl_t(i,j,k+1,nstp,itrc)+ &
1400!^ & Hz(i,j,k+1)* &
1401!^ & tl_t(i,j,k ,nstp,itrc)+ &
1402!^ & tl_Hz(i,j,k )* &
1403!^ & t(i,j,k+1,nstp,itrc)+ &
1404!^ & tl_Hz(i,j,k+1)* &
1405!^ & t(i,j,k ,nstp,itrc))- &
1406!^ & (tl_Hz(i,j,k+1)*FC(i,k-1)+ &
1407!^ & 2.0_r8*(tl_Hz(i,j,k )+ &
1408!^ & tl_Hz(i,j,k+1))*FC(i,k)+ &
1409!^ & tl_Hz(i,j,k )*FC(i,k+1))- &
1410!^ & Hz(i,j,k+1)*tl_FC(i,k-1))
1411!^
1412 adfac=cff*ad_fc(i,k)
1413 adfac1=3.0_r8*adfac
1414 adfac2=2.0_r8*adfac
1415 ad_t(i,j,k ,nstp,itrc)=ad_t(i,j,k ,nstp,itrc)+ &
1416 & hz(i,j,k+1)*adfac1
1417 ad_t(i,j,k+1,nstp,itrc)=ad_t(i,j,k+1,nstp,itrc)+ &
1418 & hz(i,j,k )*adfac1
1419 ad_hz(i,j,k )=ad_hz(i,j,k )+ &
1420 & t(i,j,k+1,nstp,itrc)*adfac1- &
1421 & fc(i,k )*adfac2- &
1422 & fc(i,k+1)*adfac
1423 ad_hz(i,j,k+1)=ad_hz(i,j,k+1)+ &
1424 & t(i,j,k ,nstp,itrc)*adfac1- &
1425 & fc(i,k-1)*adfac- &
1426 & fc(i,k )*adfac2
1427 ad_fc(i,k-1)=ad_fc(i,k-1)-hz(i,j,k+1)*adfac
1428 ad_fc(i,k)=0.0_r8
1429 END DO
1430 END DO
1431!
1432 DO i=istr,iend
1433# ifdef NEUMANN
1434!^ tl_FC(i,0)=1.5_r8*tl_t(i,j,1,nstp,itrc)
1435!^
1436 ad_t(i,j,1,nstp,itrc)=ad_t(i,j,1,nstp,itrc)+ &
1437 & 1.5_r8*ad_fc(i,0)
1438 ad_fc(i,0)=0.0_r8
1439# else
1440!^ tl_FC(i,0)=2.0_r8*tl_t(i,j,1,nstp,itrc)
1441!^
1442 ad_t(i,j,1,nstp,itrc)=ad_t(i,j,1,nstp,itrc)+ &
1443 & 2.0_r8*ad_fc(i,0)
1444 ad_fc(i,0)=0.0_r8
1445# endif
1446 END DO
1447!
1448 ELSE IF (ad_vadvection(itrc,ng)%AKIMA4) THEN
1449!
1450! Fourth-order, adjoint Akima vertical advective flux.
1451!
1452 DO k=1,n(ng)-1
1453 DO i=istr,iend
1454 fc(i,k)=t(i,j,k+1,nstp,itrc)- &
1455 & t(i,j,k ,nstp,itrc)
1456 END DO
1457 END DO
1458 DO i=istr,iend
1459 fc(i,0)=fc(i,1)
1460 fc(i,n(ng))=fc(i,n(ng)-1)
1461 END DO
1462 DO k=1,n(ng)
1463 DO i=istr,iend
1464 cff=2.0_r8*fc(i,k)*fc(i,k-1)
1465 IF (cff.gt.eps) THEN
1466 cf(i,k)=cff/(fc(i,k)+fc(i,k-1))
1467 ELSE
1468 cf(i,k)=0.0_r8
1469 END IF
1470 END DO
1471 END DO
1472 DO i=istr,iend
1473!^ tl_FC(i,N(ng))=0.0_r8
1474!^
1475 ad_fc(i,n(ng))=0.0_r8
1476!^ tl_FC(i,0)=0.0_r8
1477!^
1478 ad_fc(i,0)=0.0_r8
1479 END DO
1480 cff1=1.0_r8/3.0_r8
1481 DO k=1,n(ng)-1
1482 DO i=istr,iend
1483!^ tl_FC(i,k)=0.5_r8* &
1484!^ & (tl_W(i,j,k)* &
1485!^ & (t(i,j,k ,nstp,itrc)+ &
1486!^ & t(i,j,k+1,nstp,itrc)- &
1487!^ & cff1*(CF(i,k+1)-CF(i,k)))+ &
1488!^ & W(i,j,k)* &
1489!^ & (tl_t(i,j,k ,nstp,itrc)+ &
1490!^ & tl_t(i,j,k+1,nstp,itrc)- &
1491!^ & cff1*(tl_CF(i,k+1)-tl_CF(i,k))))
1492!^
1493 adfac=0.5_r8*ad_fc(i,k)
1494 adfac1=adfac*w(i,j,k)
1495 adfac2=adfac1*cff1
1496 ad_cf(i,k )=ad_cf(i,k )+adfac2
1497 ad_cf(i,k+1)=ad_cf(i,k+1)-adfac2
1498 ad_t(i,j,k ,nstp,itrc)=ad_t(i,j,k ,nstp,itrc)+adfac1
1499 ad_t(i,j,k+1,nstp,itrc)=ad_t(i,j,k+1,nstp,itrc)+adfac1
1500 ad_w(i,j,k)=ad_w(i,j,k)+ &
1501 & (t(i,j,k ,nstp,itrc)+ &
1502 & t(i,j,k+1,nstp,itrc)- &
1503 & cff1*(cf(i,k+1)-cf(i,k)))*adfac
1504 ad_fc(i,k)=0.0_r8
1505 END DO
1506 END DO
1507 DO k=1,n(ng)
1508 DO i=istr,iend
1509 cff=2.0_r8*fc(i,k)*fc(i,k-1)
1510 IF (cff.gt.eps) THEN
1511!^ tl_CF(i,k)=((FC(i,k)+FC(i,k-1))*tl_cff- &
1512!^ & cff*(tl_FC(i,k)+tl_FC(i,k-1)))/ &
1513!^ & ((FC(i,k)+FC(i,k-1))*(FC(i,k)+FC(i,k-1)))
1514!^
1515 adfac=ad_cf(i,k)/ &
1516 & ((fc(i,k)+fc(i,k-1))*(fc(i,k)+fc(i,k-1)))
1517 adfac1=adfac*cff
1518 ad_fc(i,k-1)=ad_fc(i,k-1)-adfac1
1519 ad_fc(i,k )=ad_fc(i,k )-adfac1
1520 ad_cff=ad_cff+(fc(i,k)+fc(i,k-1))*adfac
1521 ad_cf(i,k)=0.0_r8
1522 ELSE
1523!^ tl_CF(i,k)=0.0_r8
1524!^
1525 ad_cf(i,k)=0.0_r8
1526 END IF
1527!^ tl_cff=2.0_r8*(tl_FC(i,k)*FC(i,k-1)+ &
1528!^ & FC(i,k)*tl_FC(i,k-1))
1529!^
1530 adfac=2.0_r8*ad_cff
1531 ad_fc(i,k-1)=ad_fc(i,k-1)+fc(i,k )*adfac
1532 ad_fc(i,k )=ad_fc(i,k )+fc(i,k-1)*adfac
1533 ad_cff=0.0_r8
1534 END DO
1535 END DO
1536 DO i=istr,iend
1537!^ tl_FC(i,N(ng))=tl_FC(i,N(ng)-1)
1538!^
1539 ad_fc(i,n(ng)-1)=ad_fc(i,n(ng)-1)+ad_fc(i,n(ng))
1540 ad_fc(i,n(ng))=0.0_r8
1541!^ tl_FC(i,0)=tl_FC(i,1)
1542!^
1543 ad_fc(i,1)=ad_fc(i,1)+ad_fc(i,0)
1544 ad_fc(i,0)=0.0_r8
1545 END DO
1546 DO k=1,n(ng)-1
1547 DO i=istr,iend
1548!^ tl_FC(i,k)=tl_t(i,j,k+1,nstp,itrc)- &
1549!^ & tl_t(i,j,k ,nstp,itrc)
1550!^
1551 ad_t(i,j,k ,nstp,itrc)=ad_t(i,j,k ,nstp,itrc)- &
1552 & ad_fc(i,k)
1553 ad_t(i,j,k+1,nstp,itrc)=ad_t(i,j,k+1,nstp,itrc)+ &
1554 & ad_fc(i,k)
1555 ad_fc(i,k)=0.0_r8
1556 END DO
1557 END DO
1558!
1559 ELSE IF (ad_vadvection(itrc,ng)%CENTERED2) THEN
1560!
1561! Second-order, central differences vertical advective flux.
1562!
1563 DO i=istr,iend
1564!^ tl_FC(i,N(ng))=0.0_r8
1565!^
1566 ad_fc(i,n(ng))=0.0_r8
1567!^ tl_FC(i,0)=0.0_r8
1568!^
1569 ad_fc(i,0)=0.0_r8
1570 END DO
1571 DO k=1,n(ng)-1
1572 DO i=istr,iend
1573!^ tl_FC(i,k)=0.5_r8* &
1574!^ & (tl_W(i,j,k)* &
1575!^ & (t(i,j,k ,nstp,itrc)+ &
1576!^ & t(i,j,k+1,nstp,itrc))+ &
1577!^ & W(i,j,k)* &
1578!^ & (tl_t(i,j,k ,nstp,itrc)+ &
1579!^ & tl_t(i,j,k+1,nstp,itrc)))
1580!^
1581 adfac=0.5_r8*ad_fc(i,k)
1582 adfac1=adfac*w(i,j,k)
1583 ad_w(i,j,k)=ad_w(i,j,k)+ &
1584 & (t(i,j,k ,nstp,itrc)+ &
1585 & t(i,j,k+1,nstp,itrc))*adfac
1586 ad_t(i,j,k ,nstp,itrc)=ad_t(i,j,k ,nstp,itrc)+adfac1
1587 ad_t(i,j,k+1,nstp,itrc)=ad_t(i,j,k+1,nstp,itrc)+adfac1
1588 ad_fc(i,k)=0.0_r8
1589 END DO
1590 END DO
1591!
1592 ELSE IF ((ad_vadvection(itrc,ng)%MPDATA).or. &
1593 & (ad_vadvection(itrc,ng)%HSIMT)) THEN
1594!
1595! First_order, upstream differences vertical advective flux.
1596!
1597 DO i=istr,iend
1598!^ tl_FC(i,N(ng))=0.0_r8
1599!^
1600 ad_fc(i,n(ng))=0.0_r8
1601!^ tl_FC(i,0)=0.0_r8
1602!^
1603 ad_fc(i,0)=0.0_r8
1604 END DO
1605 DO k=1,n(ng)-1
1606 DO i=istr,iend
1607 cff1=max(w(i,j,k),0.0_r8)
1608 cff2=min(w(i,j,k),0.0_r8)
1609!^ tl_FC(i,k)=tl_cff1*t(i,j,k ,nstp,itrc)+ &
1610!^ & cff1*tl_t(i,j,k ,nstp,itrc)+ &
1611!^ & tl_cff2*t(i,j,k+1,nstp,itrc)+ &
1612!^ & cff2*tl_t(i,j,k+1,nstp,itrc)
1613!^
1614 ad_t(i,j,k ,nstp,itrc)=ad_t(i,j,k ,nstp,itrc)+ &
1615 & cff1*ad_fc(i,k)
1616 ad_t(i,j,k+1,nstp,itrc)=ad_t(i,j,k+1,nstp,itrc)+ &
1617 & cff2*ad_fc(i,k)
1618 ad_cff1=ad_cff1+t(i,j,k ,nstp,itrc)*ad_fc(i,k)
1619 ad_cff2=ad_cff2+t(i,j,k+1,nstp,itrc)*ad_fc(i,k)
1620 ad_fc(i,k)=0.0_r8
1621!^ tl_cff1=(0.5_r8+SIGN(0.5_r8, W(i,j,k)))*tl_W(i,j,k)
1622!^ tl_cff2=(0.5_r8+SIGN(0.5_r8,-W(i,j,k)))*tl_W(i,j,k)
1623!^
1624 ad_w(i,j,k)=ad_w(i,j,k)+ &
1625 & (0.5_r8+sign(0.5_r8,-w(i,j,k)))*ad_cff2+ &
1626 & (0.5_r8+sign(0.5_r8, w(i,j,k)))*ad_cff1
1627 ad_cff2=0.0_r8
1628 ad_cff1=0.0_r8
1629 END DO
1630 END DO
1631!
1632 ELSE IF ((ad_vadvection(itrc,ng)%CENTERED4).or. &
1633 & (ad_vadvection(itrc,ng)%SPLIT_U3)) THEN
1634!
1635! Fourth-order, central differences vertical advective flux.
1636!
1637 cff1=0.5_r8
1638 cff2=7.0_r8/12.0_r8
1639 cff3=1.0_r8/12.0_r8
1640 DO i=istr,iend
1641!^ tl_FC(i,N(ng))=0.0_r8
1642!^
1643 ad_fc(i,n(ng))=0.0_r8
1644!^ tl_FC(i,N(ng)-1)=tl_W(i,j,N(ng)-1)* &
1645!^ & (cff1*t(i,j,N(ng) ,nstp,itrc)+ &
1646!^ & cff2*t(i,j,N(ng)-1,nstp,itrc)- &
1647!^ & cff3*t(i,j,N(ng)-2,nstp,itrc))+ &
1648!^ & W(i,j,N(ng)-1)* &
1649!^ & (cff1*tl_t(i,j,N(ng) ,nstp,itrc)+ &
1650!^ & cff2*tl_t(i,j,N(ng)-1,nstp,itrc)- &
1651!^ & cff3*tl_t(i,j,N(ng)-2,nstp,itrc))
1652!^
1653 adfac=w(i,j,n(ng)-1)*ad_fc(i,n(ng)-1)
1654 ad_w(i,j,n(ng)-1)=ad_w(i,j,n(ng)-1)+ &
1655 & (cff1*t(i,j,n(ng) ,nstp,itrc)+ &
1656 & cff2*t(i,j,n(ng)-1,nstp,itrc)- &
1657 & cff3*t(i,j,n(ng)-2,nstp,itrc))* &
1658 & ad_fc(i,n(ng)-1)
1659 ad_t(i,j,n(ng)-2,nstp,itrc)=ad_t(i,j,n(ng)-2,nstp,itrc)- &
1660 & cff3*adfac
1661 ad_t(i,j,n(ng)-1,nstp,itrc)=ad_t(i,j,n(ng)-1,nstp,itrc)+ &
1662 & cff2*adfac
1663 ad_t(i,j,n(ng) ,nstp,itrc)=ad_t(i,j,n(ng) ,nstp,itrc)+ &
1664 & cff1*adfac
1665 ad_fc(i,n(ng)-1)=0.0_r8
1666!^ tl_FC(i,1)=tl_W(i,j,1)* &
1667!^ & (cff1*t(i,j,1,nstp,itrc)+ &
1668!^ & cff2*t(i,j,2,nstp,itrc)- &
1669!^ & cff3*t(i,j,3,nstp,itrc))+ &
1670!^ & W(i,j,1)* &
1671!^ & (cff1*tl_t(i,j,1,nstp,itrc)+ &
1672!^ & cff2*tl_t(i,j,2,nstp,itrc)- &
1673!^ & cff3*tl_t(i,j,3,nstp,itrc))
1674!^
1675 adfac=w(i,j,1)*ad_fc(i,1)
1676 ad_w(i,j,1)=ad_w(i,j,1)+ &
1677 & (cff1*t(i,j,1,nstp,itrc)+ &
1678 & cff2*t(i,j,2,nstp,itrc)- &
1679 & cff3*t(i,j,3,nstp,itrc))*ad_fc(i,1)
1680 ad_t(i,j,1,nstp,itrc)=ad_t(i,j,1,nstp,itrc)+cff1*adfac
1681 ad_t(i,j,2,nstp,itrc)=ad_t(i,j,2,nstp,itrc)+cff2*adfac
1682 ad_t(i,j,3,nstp,itrc)=ad_t(i,j,3,nstp,itrc)-cff3*adfac
1683 ad_fc(i,1)=0.0_r8
1684!^ tl_FC(i,0)=0.0_r8
1685!^
1686 ad_fc(i,0)=0.0_r8
1687 END DO
1688 DO k=2,n(ng)-2
1689 DO i=istr,iend
1690!^ tl_FC(i,k)=tl_W(i,j,k)* &
1691!^ & (cff2*(t(i,j,k ,nstp,itrc)+ &
1692!^ & t(i,j,k+1,nstp,itrc))- &
1693!^ & cff3*(t(i,j,k-1,nstp,itrc)+ &
1694!^ & t(i,j,k+2,nstp,itrc)))+ &
1695!^ & W(i,j,k)* &
1696!^ & (cff2*(tl_t(i,j,k ,nstp,itrc)+ &
1697!^ & tl_t(i,j,k+1,nstp,itrc))- &
1698!^ & cff3*(tl_t(i,j,k-1,nstp,itrc)+ &
1699!^ & tl_t(i,j,k+2,nstp,itrc)))
1700!^
1701 adfac=w(i,j,k)*ad_fc(i,k)
1702 adfac1=adfac*cff2
1703 adfac2=adfac*cff3
1704 ad_w(i,j,k)=ad_w(i,j,k)+ &
1705 & (cff2*(t(i,j,k ,nstp,itrc)+ &
1706 & t(i,j,k+1,nstp,itrc))- &
1707 & cff3*(t(i,j,k-1,nstp,itrc)+ &
1708 & t(i,j,k+2,nstp,itrc)))*ad_fc(i,k)
1709 ad_t(i,j,k-1,nstp,itrc)=ad_t(i,j,k-1,nstp,itrc)-adfac2
1710 ad_t(i,j,k ,nstp,itrc)=ad_t(i,j,k ,nstp,itrc)+adfac1
1711 ad_t(i,j,k+1,nstp,itrc)=ad_t(i,j,k+1,nstp,itrc)+adfac1
1712 ad_t(i,j,k+2,nstp,itrc)=ad_t(i,j,k+2,nstp,itrc)-adfac2
1713 ad_fc(i,k)=0.0_r8
1714 END DO
1715 END DO
1716 END IF vadv_flux2
1717 END DO t_loop2
1718 END DO j_loop1
1719
1720# if defined AGE_MEAN && defined T_PASSIVE
1721!
1722! If inert passive tracer and Mean Age, compute age concentration (even
1723! inert index) forced by the right-hand-side term that is concentration
1724! of an associated conservative passive tracer (odd inert index).
1725! Implemented and tested by W.G. Zhang and J. Wilkin. (m Tunits)
1726!
1727 DO itrc=1,npt,2
1728 IF (.not.((ad_hadvection(inert(itrc),ng)%MPDATA).or. &
1729 & (ad_hadvection(inert(itrc),ng)%HSIMT))) THEN
1730 IF (iic(ng).eq.ntfirst(ng)) THEN
1731 cff=0.5_r8*dt(ng)
1732 ELSE
1733 gamma=1.0_r8/6.0_r8
1734 cff=(1.0_r8-gamma)*dt(ng)
1735 END IF
1736 iage=inert(itrc+1) ! even inert tracer index
1737 DO k=1,n(ng)
1738 DO j=jstr,jend
1739 DO i=istr,iend
1740!^ tl_t(i,j,k,3,iage)=tl_t(i,j,k,3,iage)+ &
1741!^ & cff* &
1742!^ & (Hz(i,j,k)* &
1743!^ & tl_t(i,j,k,nnew,inert(itrc))+ &
1744!^ & tl_Hz(i,j,k)* &
1745!^ & t(i,j,k,nnew,inert(itrc)))
1746!^
1747 adfac=cff*ad_t(i,j,k,3,iage)
1748 ad_t(i,j,k,nnew,inert(itrc))=ad_t(i,j,k,nnew, &
1749 & inert(itrc))+ &
1750 & hz(i,j,k)*adfac
1751 ad_hz(i,j,k)=ad_hz(i,j,k)+ &
1752 & t(i,j,k,nnew,inert(itrc))*adfac
1753 END DO
1754 END DO
1755 END DO
1756 END IF
1757 END DO
1758# endif
1759!
1760! Compute adjoint time rate of change of intermediate tracer due to
1761! horizontal advection.
1762!
1763 t_loop1 :DO itrc=1,nt(ng)
1764 k_loop: DO k=1,n(ng)
1765!
1766! Time-step horizontal advection (m Tunits).
1767!
1768 IF ((ad_hadvection(itrc,ng)%MPDATA).or. &
1769 & (ad_hadvection(itrc,ng)%HSIMT)) THEN
1770 gamma=0.5_r8
1771 ELSE
1772 gamma=1.0_r8/6.0_r8
1773 END IF
1774 IF (iic(ng).eq.ntfirst(ng)) THEN
1775 cff=0.5_r8*dt(ng)
1776 cff1=1.0_r8
1777 cff2=0.0_r8
1778 ELSE
1779 cff=(1.0_r8-gamma)*dt(ng)
1780 cff1=0.5_r8+gamma
1781 cff2=0.5_r8-gamma
1782 END IF
1783 DO j=jstr,jend
1784 DO i=istr,iend
1785!^ tl_t(i,j,k,3,itrc)=tl_Hz(i,j,k)* &
1786!^ & (cff1*t(i,j,k,nstp,itrc)+ &
1787!^ & cff2*t(i,j,k,nnew,itrc))+ &
1788!^ & Hz(i,j,k)* &
1789!^ & (cff1*tl_t(i,j,k,nstp,itrc)+ &
1790!^ & cff2*tl_t(i,j,k,nnew,itrc))- &
1791!^ & cff*pm(i,j)*pn(i,j)* &
1792!^ & (tl_FX(i+1,j)-tl_FX(i,j)+ &
1793!^ & tl_FE(i,j+1)-tl_FE(i,j))
1794!^
1795 adfac1=hz(i,j,k)*ad_t(i,j,k,3,itrc)
1796 adfac2=cff*pm(i,j)*pn(i,j)*ad_t(i,j,k,3,itrc)
1797 ad_hz(i,j,k)=ad_hz(i,j,k)+ &
1798 & (cff1*t(i,j,k,nstp,itrc)+ &
1799 & cff2*t(i,j,k,nnew,itrc))*ad_t(i,j,k,3,itrc)
1800 ad_t(i,j,k,nstp,itrc)=ad_t(i,j,k,nstp,itrc)+cff1*adfac1
1801 ad_t(i,j,k,nnew,itrc)=ad_t(i,j,k,nnew,itrc)+cff2*adfac1
1802 ad_fe(i,j )=ad_fe(i,j )+adfac2
1803 ad_fe(i,j+1)=ad_fe(i,j+1)-adfac2
1804 ad_fx(i ,j)=ad_fx(i ,j)+adfac2
1805 ad_fx(i+1,j)=ad_fx(i+1,j)-adfac2
1806 ad_t(i,j,k,3,itrc)=0.0_r8
1807 END DO
1808 END DO
1809!
1810! Apply tracers point sources to the horizontal advection terms,
1811! if any.
1812!
1813! Dsrc(is) = 0, flow across grid cell u-face (positive or negative)
1814! Dsrc(is) = 1, flow across grid cell v-face (positive or negative)
1815!
1816 IF (luvsrc(ng)) THEN
1817 DO is=1,nsrc(ng)
1818 isrc=sources(ng)%Isrc(is)
1819 jsrc=sources(ng)%Jsrc(is)
1820 IF (((istr.le.isrc).and.(isrc.le.iend+1)).and. &
1821 & ((jstr.le.jsrc).and.(jsrc.le.jend+1))) THEN
1822 IF (int(sources(ng)%Dsrc(is)).eq.0) THEN
1823 IF (ltracersrc(itrc,ng)) THEN
1824!^ tl_FX(Isrc,Jsrc)=tl_Huon(Isrc,Jsrc,k)* &
1825!^ & SOURCES(ng)%Tsrc(is,k,itrc)+ &
1826!^ & Huon(Isrc,Jsrc,k)* &
1827!^ & SOURCES(ng)%tl_Tsrc(is,k,itrc)
1828!^
1829 ad_huon(isrc,jsrc,k)=ad_huon(isrc,jsrc,k)+ &
1830 & sources(ng)%Tsrc(is,k,itrc)* &
1831 & ad_fx(isrc,jsrc)
1832 sources(ng)%ad_Tsrc(is,k,itrc)= &
1833 & sources(ng)%ad_Tsrc(is,k,itrc)+ &
1834 & huon(isrc,jsrc,k)* &
1835 & ad_fx(isrc,jsrc)
1836 ad_fx(isrc,jsrc)=0.0_r8
1837 ELSE
1838!^ tl_FX(Isrc,Jsrc)=0.0_r8
1839!^
1840 ad_fx(isrc,jsrc)=0.0_r8
1841 END IF
1842 ELSE IF (int(sources(ng)%Dsrc(is)).eq.1) THEN
1843 IF (ltracersrc(itrc,ng)) THEN
1844!^ tl_FE(Isrc,Jsrc)=tl_Hvom(Isrc,Jsrc,k)* &
1845!^ & SOURCES(ng)%Tsrc(is,k,itrc)+ &
1846!^ & Hvom(Isrc,Jsrc,k)* &
1847!^ & SOURCES(ng)%tl_Tsrc(is,k,itrc)
1848!^
1849 ad_hvom(isrc,jsrc,k)=ad_hvom(isrc,jsrc,k)+ &
1850 & sources(ng)%Tsrc(is,k,itrc)* &
1851 & ad_fe(isrc,jsrc)
1852 sources(ng)%ad_Tsrc(is,k,itrc)= &
1853 & sources(ng)%ad_Tsrc(is,k,itrc)+ &
1854 & hvom(isrc,jsrc,k)* &
1855 & ad_fe(isrc,jsrc)
1856 ad_fe(isrc,jsrc)=0.0_r8
1857 ELSE
1858!^ tl_FE(Isrc,Jsrc)=0.0_r8
1859!^
1860 ad_fe(isrc,jsrc)=0.0_r8
1861 END IF
1862 END IF
1863 END IF
1864 END DO
1865 END IF
1866!
1867! Compute adjoint of intermediate tracer horizontal advection fluxes.
1868!
1869 hadv_flux : IF (ad_hadvection(itrc,ng)%CENTERED2) THEN
1870!
1871! Second-order, centered differences horizontal advective fluxes.
1872!
1873 DO j=jstr,jend+1
1874 DO i=istr,iend
1875!^ tl_FE(i,j)=0.5_r8* &
1876!^ & (tl_Hvom(i,j,k)* &
1877!^ & (t(i,j-1,k,nstp,itrc)+ &
1878!^ & t(i,j ,k,nstp,itrc))+ &
1879!^ & Hvom(i,j,k)* &
1880!^ & (tl_t(i,j-1,k,nstp,itrc)+ &
1881!^ & tl_t(i,j ,k,nstp,itrc)))
1882!^
1883 adfac=0.5_r8*ad_fe(i,j)
1884 adfac1=adfac*hvom(i,j,k)
1885 ad_hvom(i,j,k)=ad_hvom(i,j,k)+ &
1886 & adfac*(t(i,j ,k,nstp,itrc)+ &
1887 & t(i,j-1,k,nstp,itrc))
1888 ad_t(i,j-1,k,nstp,itrc)=ad_t(i,j-1,k,nstp,itrc)+adfac1
1889 ad_t(i,j ,k,nstp,itrc)=ad_t(i,j ,k,nstp,itrc)+adfac1
1890 ad_fe(i,j)=0.0_r8
1891 END DO
1892 END DO
1893 DO j=jstr,jend
1894 DO i=istr,iend+1
1895!^ tl_FX(i,j)=0.5_r8* &
1896!^ & (tl_Huon(i,j,k)* &
1897!^ & (t(i-1,j,k,nstp,itrc)+ &
1898!^ & t(i ,j,k,nstp,itrc))+ &
1899!^ & Huon(i,j,k)* &
1900!^ & (tl_t(i-1,j,k,nstp,itrc)+ &
1901!^ & tl_t(i ,j,k,nstp,itrc)))
1902!^
1903 adfac=0.5_r8*ad_fx(i,j)
1904 adfac1=adfac*huon(i,j,k)
1905 ad_huon(i,j,k)=ad_huon(i,j,k)+ &
1906 & adfac*(t(i-1,j,k,nstp,itrc)+ &
1907 & t(i ,j,k,nstp,itrc))
1908 ad_t(i-1,j,k,nstp,itrc)=ad_t(i-1,j,k,nstp,itrc)+adfac1
1909 ad_t(i ,j,k,nstp,itrc)=ad_t(i ,j,k,nstp,itrc)+adfac1
1910 ad_fx(i,j)=0.0_r8
1911 END DO
1912 END DO
1913!
1914 ELSE IF ((ad_hadvection(itrc,ng)%MPDATA).or. &
1915 & (ad_hadvection(itrc,ng)%HSIMT)) THEN
1916!
1917! First-order, upstream differences horizontal advective fluxes.
1918!
1919 DO j=jstr,jend+1
1920 DO i=istr,iend
1921 cff1=max(hvom(i,j,k),0.0_r8)
1922 cff2=min(hvom(i,j,k),0.0_r8)
1923!^ tl_FE(i,j)=tl_cff1*t(i,j-1,k,nstp,itrc)+ &
1924!^ & cff1*tl_t(i,j-1,k,nstp,itrc)+ &
1925!^ & tl_cff2*t(i,j ,k,nstp,itrc)+ &
1926!^ & cff2*tl_t(i,j ,k,nstp,itrc)
1927!^
1928 ad_t(i,j-1,k,nstp,itrc)=ad_t(i,j-1,k,nstp,itrc)+ &
1929 & cff1*ad_fe(i,j)
1930 ad_t(i,j ,k,nstp,itrc)=ad_t(i,j ,k,nstp,itrc)+ &
1931 & cff2*ad_fe(i,j)
1932 ad_cff1=ad_cff1+t(i,j-1,k,nstp,itrc)*ad_fe(i,j)
1933 ad_cff2=ad_cff2+t(i,j ,k,nstp,itrc)*ad_fe(i,j)
1934 ad_fe(i,j)=0.0_r8
1935!^ tl_cff1=(0.5_r8+SIGN(0.5_r8, Hvom(i,j,k)))* &
1936!^ & tl_Hvom(i,j,k)
1937!^ tl_cff2=(0.5_r8+SIGN(0.5_r8,-Hvom(i,j,k)))* &
1938!^ & tl_Hvom(i,j,k)
1939!^
1940 ad_hvom(i,j,k)=ad_hvom(i,j,k)+ &
1941 & (0.5_r8+sign(0.5_r8,-hvom(i,j,k)))* &
1942 & ad_cff2+ &
1943 & (0.5_r8+sign(0.5_r8, hvom(i,j,k)))* &
1944 & ad_cff1
1945 ad_cff2=0.0_r8
1946 ad_cff1=0.0_r8
1947 END DO
1948 END DO
1949 DO j=jstr,jend
1950 DO i=istr,iend+1
1951 cff1=max(huon(i,j,k),0.0_r8)
1952 cff2=min(huon(i,j,k),0.0_r8)
1953!^ tl_FX(i,j)=tl_cff1*t(i-1,j,k,nstp,itrc)+ &
1954!^ & cff1*tl_t(i-1,j,k,nstp,itrc)+ &
1955!^ & tl_cff2*t(i ,j,k,nstp,itrc)+ &
1956!^ & cff2*tl_t(i ,j,k,nstp,itrc)
1957!^
1958 ad_t(i ,j,k,nstp,itrc)=ad_t(i ,j,k,nstp,itrc)+ &
1959 & cff2*ad_fx(i,j)
1960 ad_t(i-1,j,k,nstp,itrc)=ad_t(i-1,j,k,nstp,itrc)+ &
1961 & cff1*ad_fx(i,j)
1962 ad_cff1=ad_cff1+t(i-1,j,k,nstp,itrc)*ad_fx(i,j)
1963 ad_cff2=ad_cff2+t(i ,j,k,nstp,itrc)*ad_fx(i,j)
1964 ad_fx(i,j)=0.0_r8
1965!^ tl_cff1=(0.5_r8+SIGN(0.5_r8, Huon(i,j,k)))* &
1966!^ & tl_Huon(i,j,k)
1967!^ tl_cff2=(0.5_r8+SIGN(0.5_r8,-Huon(i,j,k)))* &
1968!^ & tl_Huon(i,j,k)
1969!^
1970 ad_huon(i,j,k)=ad_huon(i,j,k)+ &
1971 & (0.5_r8+sign(0.5_r8,-huon(i,j,k)))* &
1972 & ad_cff2+ &
1973 & (0.5_r8+sign(0.5_r8, huon(i,j,k)))* &
1974 & ad_cff1
1975 ad_cff2=0.0_r8
1976 ad_cff1=0.0_r8
1977 END DO
1978 END DO
1979!
1980 ELSE IF ((ad_hadvection(itrc,ng)%AKIMA4).or. &
1981 & (ad_hadvection(itrc,ng)%CENTERED4).or. &
1982 & (ad_hadvection(itrc,ng)%SPLIT_U3).or. &
1983 & (ad_hadvection(itrc,ng)%UPSTREAM3)) THEN
1984!
1985! Fourth-order Akima, fourth-order centered differences, or third-order
1986! upstream-biased horizontal advective fluxes.
1987! Compute BASIC STATE "curv" and "grad" scratch arrays.
1988!
1989 DO j=jstrm1,jendp2
1990 DO i=istr,iend
1991 fe(i,j)=t(i,j ,k,nstp,itrc)- &
1992 & t(i,j-1,k,nstp,itrc)
1993# ifdef MASKING
1994 fe(i,j)=fe(i,j)*vmask(i,j)
1995# endif
1996 END DO
1997 END DO
1998 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
1999 IF (domain(ng)%Southern_Edge(tile)) THEN
2000 DO i=istr,iend
2001 fe(i,jstr-1)=fe(i,jstr)
2002 END DO
2003 END IF
2004 END IF
2005 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
2006 IF (domain(ng)%Northern_Edge(tile)) THEN
2007 DO i=istr,iend
2008 fe(i,jend+2)=fe(i,jend+1)
2009 END DO
2010 END IF
2011 END IF
2012!
2013 DO j=jstr-1,jend+1
2014 DO i=istr,iend
2015 IF (ad_hadvection(itrc,ng)%UPSTREAM3) THEN
2016 curv(i,j)=fe(i,j+1)-fe(i,j)
2017 ELSE IF (ad_hadvection(itrc,ng)%AKIMA4) THEN
2018 cff=2.0_r8*fe(i,j+1)*fe(i,j)
2019 IF (cff.gt.eps) THEN
2020 grad(i,j)=cff/(fe(i,j+1)+fe(i,j))
2021 ELSE
2022 grad(i,j)=0.0_r8
2023 END IF
2024 ELSE IF ((ad_hadvection(itrc,ng)%CENTERED4).or. &
2025 & (ad_hadvection(itrc,ng)%SPLIT_U3)) THEN
2026 grad(i,j)=0.5_r8*(fe(i,j+1)+fe(i,j))
2027 END IF
2028 END DO
2029 END DO
2030!
2031 cff1=1.0_r8/6.0_r8
2032 cff2=1.0_r8/3.0_r8
2033 DO j=jstr,jend+1
2034 DO i=istr,iend
2035 IF (ad_hadvection(itrc,ng)%UPSTREAM3) THEN
2036!^ tl_FE(i,j)=0.5_r8* &
2037!^ & (tl_Hvom(i,j,k)* &
2038!^ & (t(i,j-1,k,nstp,itrc)+ &
2039!^ & t(i,j ,k,nstp,itrc))+ &
2040!^ & Hvom(i,j,k)* &
2041!^ & (tl_t(i,j-1,k,nstp,itrc)+ &
2042!^ & tl_t(i,j ,k,nstp,itrc)))- &
2043!^ & cff1* &
2044!^ & (tl_curv(i,j-1)*MAX(Hvom(i,j,k),0.0_r8)+ &
2045!^ & curv(i,j-1)* &
2046!^ & (0.5_r8+SIGN(0.5_r8, Hvom(i,j,k)))* &
2047!^ & tl_Hvom(i,j,k)+ &
2048!^ & tl_curv(i,j )*MIN(Hvom(i,j,k),0.0_r8)+ &
2049!^ & curv(i,j )* &
2050!^ & (0.5_r8+SIGN(0.5_r8,-Hvom(i,j,k)))* &
2051!^ & tl_Hvom(i,j,k))
2052!^
2053 adfac=0.5_r8*ad_fe(i,j)
2054 adfac1=adfac*hvom(i,j,k)
2055 adfac2=cff1*ad_fe(i,j)
2056 ad_hvom(i,j,k)=ad_hvom(i,j,k)+ &
2057 & (t(i,j-1,k,nstp,itrc)+ &
2058 & t(i,j ,k,nstp,itrc))*adfac- &
2059 & (curv(i,j-1)* &
2060 & (0.5_r8+sign(0.5_r8, hvom(i,j,k)))+ &
2061 & curv(i,j )* &
2062 & (0.5_r8+sign(0.5_r8,-hvom(i,j,k))))* &
2063 & adfac2
2064 ad_t(i,j-1,k,nstp,itrc)=ad_t(i,j-1,k,nstp,itrc)+adfac1
2065 ad_t(i,j ,k,nstp,itrc)=ad_t(i,j ,k,nstp,itrc)+adfac1
2066 ad_curv(i,j-1)=ad_curv(i,j-1)- &
2067 & max(hvom(i,j,k),0.0_r8)*adfac2
2068 ad_curv(i,j )=ad_curv(i,j )- &
2069 & min(hvom(i,j,k),0.0_r8)*adfac2
2070 ad_fe(i,j)=0.0_r8
2071 ELSE IF ((ad_hadvection(itrc,ng)%AKIMA4).or. &
2072 & (ad_hadvection(itrc,ng)%CENTERED4).or. &
2073 & (ad_hadvection(itrc,ng)%SPLIT_U3)) THEN
2074!^ tl_FE(i,j)=0.5_r8* &
2075!^ & (tl_Hvom(i,j,k)* &
2076!^ & (t(i,j-1,k,nstp,itrc)+ &
2077!^ & t(i,j ,k,nstp,itrc)- &
2078!^ & cff2*(grad(i,j )- &
2079!^ & grad(i,j-1)))+ &
2080!^ & Hvom(i,j,k)* &
2081!^ & (tl_t(i,j-1,k,nstp,itrc)+ &
2082!^ & tl_t(i,j ,k,nstp,itrc)- &
2083!^ & cff2*(tl_grad(i,j )- &
2084!^ & tl_grad(i,j-1))))
2085!^
2086 adfac=0.5_r8*ad_fe(i,j)
2087 adfac1=adfac*hvom(i,j,k)
2088 adfac2=adfac1*cff2
2089 ad_hvom(i,j,k)=ad_hvom(i,j,k)+ &
2090 & adfac*(t(i,j-1,k,nstp,itrc)+ &
2091 & t(i,j ,k,nstp,itrc)- &
2092 & cff2*(grad(i,j )- &
2093 & grad(i,j-1)))
2094 ad_t(i,j-1,k,nstp,itrc)=ad_t(i,j-1,k,nstp,itrc)+adfac1
2095 ad_t(i,j ,k,nstp,itrc)=ad_t(i,j ,k,nstp,itrc)+adfac1
2096 ad_grad(i,j-1)=ad_grad(i,j-1)+adfac2
2097 ad_grad(i,j )=ad_grad(i,j )-adfac2
2098 ad_fe(i,j)=0.0_r8
2099 END IF
2100 END DO
2101 END DO
2102!
2103 DO j=jstr-1,jend+1
2104 DO i=istr,iend
2105 IF (ad_hadvection(itrc,ng)%UPSTREAM3) THEN
2106!^ tl_curv(i,j)=tl_FE(i,j+1)-tl_FE(i,j)
2107!^
2108 ad_fe(i,j )=ad_fe(i,j )-ad_curv(i,j)
2109 ad_fe(i,j+1)=ad_fe(i,j+1)+ad_curv(i,j)
2110 ad_curv(i,j)=0.0_r8
2111 ELSE IF (ad_hadvection(itrc,ng)%AKIMA4) THEN
2112 cff=2.0_r8*fe(i,j+1)*fe(i,j)
2113 IF (cff.gt.eps) THEN
2114!^ tl_grad(i,j)=((FE(i,j+1)+FE(i,j))*tl_cff- &
2115!^ & cff*(tl_FE(i,j+1)+tl_FE(i,j)))/ &
2116!^ & ((FE(i,j+1)+FE(i,j))* &
2117!^ & (FE(i,j+1)+FE(i,j)))
2118!^
2119 adfac=ad_grad(i,j)/ &
2120 & ((fe(i,j+1)+fe(i,j))*(fe(i,j+1)+fe(i,j)))
2121 adfac1=adfac*cff
2122 ad_fe(i,j )=ad_fe(i,j )-adfac1
2123 ad_fe(i,j+1)=ad_fe(i,j+1)-adfac1
2124 ad_cff=ad_cff+(fe(i,j+1)+fe(i,j))*adfac
2125 ad_grad(i,j)=0.0_r8
2126 ELSE
2127!^ tl_grad(i,j)=0.0_r8
2128!^
2129 ad_grad(i,j)=0.0_r8
2130 END IF
2131!^ tl_cff=2.0_r8*(tl_FE(i,j+1)*FE(i,j)+ &
2132!^ & FE(i,j+1)*tl_FE(i,j))
2133!^
2134 adfac=2.0_r8*ad_cff
2135 ad_fe(i,j )=ad_fe(i,j )+fe(i,j+1)*adfac
2136 ad_fe(i,j+1)=ad_fe(i,j+1)+fe(i,j )*adfac
2137 ad_cff=0.0_r8
2138 ELSE IF ((ad_hadvection(itrc,ng)%CENTERED4).or. &
2139 & (ad_hadvection(itrc,ng)%SPLIT_U3)) THEN
2140!^ tl_grad(i,j)=0.5_r8*(tl_FE(i,j+1)+tl_FE(i,j))
2141!^
2142 adfac=0.5_r8*ad_grad(i,j)
2143 ad_fe(i,j )=ad_fe(i,j )+adfac
2144 ad_fe(i,j+1)=ad_fe(i,j+1)+adfac
2145 ad_grad(i,j)=0.0_r8
2146 END IF
2147 END DO
2148 END DO
2149 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
2150 IF (domain(ng)%Northern_Edge(tile)) THEN
2151 DO i=istr,iend
2152!^ tl_FE(i,Jend+2)=tl_FE(i,Jend+1)
2153!^
2154 ad_fe(i,jend+1)=ad_fe(i,jend+1)+ad_fe(i,jend+2)
2155 ad_fe(i,jend+2)=0.0_r8
2156 END DO
2157 END IF
2158 END IF
2159 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
2160 IF (domain(ng)%Southern_Edge(tile)) THEN
2161 DO i=istr,iend
2162!^ tl_FE(i,Jstr-1)=tl_FE(i,Jstr)
2163!^
2164 ad_fe(i,jstr)=ad_fe(i,jstr)+ad_fe(i,jstr-1)
2165 ad_fe(i,jstr-1)=0.0_r8
2166 END DO
2167 END IF
2168 END IF
2169!
2170 DO j=jstrm1,jendp2
2171 DO i=istr,iend
2172# ifdef MASKING
2173!^ tl_FE(i,j)=tl_FE(i,j)*vmask(i,j)
2174!^
2175 ad_fe(i,j)=ad_fe(i,j)*vmask(i,j)
2176# endif
2177!^ tl_FE(i,j)=tl_t(i,j ,k,nstp,itrc)- &
2178!^ & tl_t(i,j-1,k,nstp,itrc)
2179!^
2180 ad_t(i,j-1,k,nstp,itrc)=ad_t(i,j-1,k,nstp,itrc)- &
2181 & ad_fe(i,j)
2182 ad_t(i,j ,k,nstp,itrc)=ad_t(i,j ,k,nstp,itrc)+ &
2183 & ad_fe(i,j)
2184 ad_fe(i,j)=0.0_r8
2185 END DO
2186 END DO
2187!
2188 DO j=jstr,jend
2189 DO i=istrm1,iendp2
2190 fx(i,j)=t(i ,j,k,nstp,itrc)- &
2191 & t(i-1,j,k,nstp,itrc)
2192# ifdef MASKING
2193 fx(i,j)=fx(i,j)*umask(i,j)
2194# endif
2195 END DO
2196 END DO
2197 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
2198 IF (domain(ng)%Western_Edge(tile)) THEN
2199 DO j=jstr,jend
2200 fx(istr-1,j)=fx(istr,j)
2201 END DO
2202 END IF
2203 END IF
2204 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
2205 IF (domain(ng)%Eastern_Edge(tile)) THEN
2206 DO j=jstr,jend
2207 fx(iend+2,j)=fx(iend+1,j)
2208 END DO
2209 END IF
2210 END IF
2211!
2212 DO j=jstr,jend
2213 DO i=istr-1,iend+1
2214 IF (ad_hadvection(itrc,ng)%UPSTREAM3) THEN
2215 curv(i,j)=fx(i+1,j)-fx(i,j)
2216 ELSE IF (ad_hadvection(itrc,ng)%AKIMA4) THEN
2217 cff=2.0_r8*fx(i+1,j)*fx(i,j)
2218 IF (cff.gt.eps) THEN
2219 grad(i,j)=cff/(fx(i+1,j)+fx(i,j))
2220 ELSE
2221 grad(i,j)=0.0_r8
2222 END IF
2223 ELSE IF ((ad_hadvection(itrc,ng)%CENTERED4).or. &
2224 & (ad_hadvection(itrc,ng)%SPLIT_U3)) THEN
2225 grad(i,j)=0.5_r8*(fx(i+1,j)+fx(i,j))
2226 END IF
2227 END DO
2228 END DO
2229!
2230 cff1=1.0_r8/6.0_r8
2231 cff2=1.0_r8/3.0_r8
2232 DO j=jstr,jend
2233 DO i=istr,iend+1
2234 IF (ad_hadvection(itrc,ng)%UPSTREAM3) THEN
2235!^ tl_FX(i,j)=0.5_r8* &
2236!^ & (tl_Huon(i,j,k)* &
2237!^ & (t(i-1,j,k,nstp,itrc)+ &
2238!^ & t(i ,j,k,nstp,itrc))+ &
2239!^ & Huon(i,j,k)* &
2240!^ & (tl_t(i-1,j,k,nstp,itrc)+ &
2241!^ & tl_t(i ,j,k,nstp,itrc)))- &
2242!^ & cff1* &
2243!^ & (tl_curv(i-1,j)*MAX(Huon(i,j,k),0.0_r8)+ &
2244!^ & curv(i-1,j)* &
2245!^ & (0.5_r8+SIGN(0.5_r8, Huon(i,j,k)))* &
2246!^ & tl_Huon(i,j,k)+ &
2247!^ & tl_curv(i ,j)*MIN(Huon(i,j,k),0.0_r8)+ &
2248!^ & curv(i ,j)* &
2249!^ & (0.5_r8+SIGN(0.5_r8,-Huon(i,j,k)))* &
2250!^ & tl_Huon(i,j,k))
2251!^
2252 adfac=0.5_r8*ad_fx(i,j)
2253 adfac1=adfac*huon(i,j,k)
2254 adfac2=cff1*ad_fx(i,j)
2255 ad_huon(i,j,k)=ad_huon(i,j,k)+ &
2256 & (t(i-1,j,k,nstp,itrc)+ &
2257 & t(i ,j,k,nstp,itrc))*adfac- &
2258 & (curv(i-1,j)* &
2259 & (0.5_r8+sign(0.5_r8, huon(i,j,k)))+ &
2260 & curv(i ,j)* &
2261 & (0.5_r8+sign(0.5_r8,-huon(i,j,k))))* &
2262 & adfac2
2263 ad_t(i-1,j,k,nstp,itrc)=ad_t(i-1,j,k,nstp,itrc)+adfac1
2264 ad_t(i ,j,k,nstp,itrc)=ad_t(i ,j,k,nstp,itrc)+adfac1
2265 ad_curv(i-1,j)=ad_curv(i-1,j)- &
2266 & max(huon(i,j,k),0.0_r8)*adfac2
2267 ad_curv(i ,j)=ad_curv(i ,j)- &
2268 & min(huon(i,j,k),0.0_r8)*adfac2
2269 ad_fx(i,j)=0.0_r8
2270 ELSE IF ((ad_hadvection(itrc,ng)%AKIMA4).or. &
2271 & (ad_hadvection(itrc,ng)%CENTERED4).or. &
2272 & (ad_hadvection(itrc,ng)%SPLIT_U3)) THEN
2273!^ tl_FX(i,j)=0.5_r8* &
2274!^ & (tl_Huon(i,j,k)* &
2275!^ & (t(i-1,j,k,nstp,itrc)+ &
2276!^ & t(i ,j,k,nstp,itrc)- &
2277!^ & cff2*(grad(i ,j)- &
2278!^ & grad(i-1,j)))+ &
2279!^ & Huon(i,j,k)* &
2280!^ & (tl_t(i-1,j,k,nstp,itrc)+ &
2281!^ & tl_t(i ,j,k,nstp,itrc)- &
2282!^ & cff2*(tl_grad(i ,j)- &
2283!^ & tl_grad(i-1,j))))
2284!^
2285 adfac=0.5_r8*ad_fx(i,j)
2286 adfac1=adfac*huon(i,j,k)
2287 adfac2=adfac1*cff2
2288 ad_huon(i,j,k)=ad_huon(i,j,k)+ &
2289 & adfac*(t(i-1,j,k,nstp,itrc)+ &
2290 & t(i ,j,k,nstp,itrc)- &
2291 & cff2*(grad(i ,j)- &
2292 & grad(i-1,j)))
2293 ad_t(i-1,j,k,nstp,itrc)=ad_t(i-1,j,k,nstp,itrc)+adfac1
2294 ad_t(i ,j,k,nstp,itrc)=ad_t(i ,j,k,nstp,itrc)+adfac1
2295 ad_grad(i-1,j)=ad_grad(i-1,j)+adfac2
2296 ad_grad(i ,j)=ad_grad(i ,j)-adfac2
2297 ad_fx(i,j)=0.0_r8
2298 END IF
2299 END DO
2300 END DO
2301!
2302 DO j=jstr,jend
2303 DO i=istr-1,iend+1
2304 IF (ad_hadvection(itrc,ng)%UPSTREAM3) THEN
2305!^ tl_curv(i,j)=tl_FX(i+1,j)-tl_FX(i,j)
2306!^
2307 ad_fx(i ,j)=ad_fx(i ,j)-ad_curv(i,j)
2308 ad_fx(i+1,j)=ad_fx(i+1,j)+ad_curv(i,j)
2309 ad_curv(i,j)=0.0_r8
2310 ELSE IF (ad_hadvection(itrc,ng)%AKIMA4) THEN
2311 cff=2.0_r8*fx(i+1,j)*fx(i,j)
2312 IF (cff.gt.eps) THEN
2313!^ tl_grad(i,j)=((FX(i+1,j)+FX(i,j))*tl_cff- &
2314!^ & cff*(tl_FX(i+1,j)+tl_FX(i,j)))/ &
2315!^ & ((FX(i+1,j)+FX(i,j))* &
2316!^ & (FX(i+1,j)+FX(i,j)))
2317!^
2318 adfac=ad_grad(i,j)/ &
2319 & ((fx(i+1,j)+fx(i,j))*(fx(i+1,j)+fx(i,j)))
2320 adfac1=adfac*cff
2321 ad_fx(i ,j)=ad_fx(i ,j)-adfac1
2322 ad_fx(i+1,j)=ad_fx(i+1,j)-adfac1
2323 ad_cff=ad_cff+(fx(i+1,j)+fx(i,j))*adfac
2324 ad_grad(i,j)=0.0_r8
2325 ELSE
2326!^ tl_grad(i,j)=0.0_r8
2327!^
2328 ad_grad(i,j)=0.0_r8
2329 END IF
2330 ELSE IF ((ad_hadvection(itrc,ng)%CENTERED4).or. &
2331 & (ad_hadvection(itrc,ng)%SPLIT_U3)) THEN
2332!^ tl_grad(i,j)=0.5_r8*(tl_FX(i+1,j)+tl_FX(i,j))
2333!^
2334 adfac=0.5_r8*ad_grad(i,j)
2335 ad_fx(i ,j)=ad_fx(i ,j)+adfac
2336 ad_fx(i+1,j)=ad_fx(i+1,j)+adfac
2337 ad_grad(i,j)=0.0_r8
2338 END IF
2339!^ tl_cff=2.0_r8*(tl_FX(i+1,j)*FX(i,j)+ &
2340!^ & FX(i+1,j)*tl_FX(i,j))
2341!^
2342 adfac=2.0_r8*ad_cff
2343 ad_fx(i ,j)=ad_fx(i ,j)+fx(i+1,j)*adfac
2344 ad_fx(i+1,j)=ad_fx(i+1,j)+fx(i ,j)*adfac
2345 ad_cff=0.0_r8
2346 END DO
2347 END DO
2348 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
2349 IF (domain(ng)%Eastern_Edge(tile)) THEN
2350 DO j=jstr,jend
2351!^ tl_FX(Iend+2,j)=tl_FX(Iend+1,j)
2352!^
2353 ad_fx(iend+1,j)=ad_fx(iend+1,j)+ad_fx(iend+2,j)
2354 ad_fx(iend+2,j)=0.0_r8
2355 END DO
2356 END IF
2357 END IF
2358 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
2359 IF (domain(ng)%Western_Edge(tile)) THEN
2360 DO j=jstr,jend
2361!^ tl_FX(Istr-1,j)=tl_FX(Istr,j)
2362!^
2363 ad_fx(istr,j)=ad_fx(istr,j)+ad_fx(istr-1,j)
2364 ad_fx(istr-1,j)=0.0_r8
2365 END DO
2366 END IF
2367 END IF
2368!
2369 DO j=jstr,jend
2370 DO i=istrm1,iendp2
2371# ifdef MASKING
2372!^ tl_FX(i,j)=tl_FX(i,j)*umask(i,j)
2373!^
2374 ad_fx(i,j)=ad_fx(i,j)*umask(i,j)
2375# endif
2376!^ tl_FX(i,j)=tl_t(i ,j,k,nstp,itrc)- &
2377!^ & tl_t(i-1,j,k,nstp,itrc)
2378!^
2379 ad_t(i-1,j,k,nstp,itrc)=ad_t(i-1,j,k,nstp,itrc)- &
2380 & ad_fx(i,j)
2381 ad_t(i ,j,k,nstp,itrc)=ad_t(i ,j,k,nstp,itrc)+ &
2382 & ad_fx(i,j)
2383 ad_fx(i,j)=0.0_r8
2384 END DO
2385 END DO
2386 END IF hadv_flux
2387 END DO k_loop
2388 END DO t_loop1
2389
2390# ifdef SOLAR_SOURCE
2391!
2392! Compute fraction of the solar shortwave radiation, "swdk"
2393! (at vertical W-points) penetrating water column. First, compute
2394! BASIC STATE arrays.
2395!
2396 DO k=1,n(ng)-1
2397 DO j=jstr,jend
2398 DO i=istr,iend
2399 fx(i,j)=z_w(i,j,n(ng))-z_w(i,j,k)
2400 END DO
2401 END DO
2402 DO j=jstr,jend
2403 DO i=istr,iend
2404!^ tl_swdk(i,j,k)=tl_FE(i,j)
2405!^
2406 ad_fe(i,j)=ad_fe(i,j)+ad_swdk(i,j,k)
2407 ad_swdk(i,j,k)=0.0_r8
2408 END DO
2409 END DO
2410!^ CALL tl_lmd_swfrac_tile (ng, tile, &
2411!^ & LBi, UBi, LBj, UBj, &
2412!^ & IminS, ImaxS, JminS, JmaxS, &
2413!^ & -1.0_r8, FX, tl_FX, tl_FE)
2414!^
2415 CALL ad_lmd_swfrac_tile (ng, tile, &
2416 & lbi, ubi, lbj, ubj, &
2417 & imins, imaxs, jmins, jmaxs, &
2418 & -1.0_r8, fx, ad_fx, ad_fe)
2419 DO j=jstr,jend
2420 DO i=istr,iend
2421!^ tl_FX(i,j)=tl_z_w(i,j,N(ng))-tl_z_w(i,j,k)
2422!^
2423 ad_z_w(i,j,k )=ad_z_w(i,j,k )-ad_fx(i,j)
2424 ad_z_w(i,j,n(ng))=ad_z_w(i,j,n(ng))+ad_fx(i,j)
2425 ad_fx(i,j)=0.0_r8
2426 END DO
2427 END DO
2428 END DO
2429# endif
2430# endif /* !TS_FIXED */
2431!
2432 RETURN
2433 END SUBROUTINE ad_pre_step3d_tile
2434#endif
2435 END MODULE ad_pre_step3d_mod
ad_lmd_swfrac_tile
subroutine ad_lmd_swfrac_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, zscale, z, ad_z, ad_swdk)
Definition ad_lmd_swfrac.F:10
ad_exchange_3d_mod
Definition ad_exchange_3d.F:2
ad_exchange_3d_mod::ad_exchange_r3d_tile
subroutine ad_exchange_r3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, ad_a)
Definition ad_exchange_3d.F:461
ad_pre_step3d_mod
Definition ad_pre_step3d.F:2
ad_pre_step3d_mod::ad_pre_step3d_tile
subroutine ad_pre_step3d_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, kstp, knew, nrhs, nstp, nnew, rmask, umask, vmask, rmask_wet, om_v, on_u, pm, pn, hz, ad_hz, huon, ad_huon, hvom, ad_hvom, z_r, ad_z_r, z_w, ad_z_w, ad_btflx, ad_bustr, ad_bvstr, ad_stflx, ad_sustr, ad_svstr, srflx, akt, ad_akt, akv, ad_akv, ad_ubar, ad_vbar, w, ad_w, ad_ru, ad_rv, t, ad_t, u, ad_u, v, ad_v)
Definition ad_pre_step3d.F:196
ad_pre_step3d_mod::ad_pre_step3d
subroutine, public ad_pre_step3d(ng, tile)
Definition ad_pre_step3d.F:45
ad_t3dbc_mod
Definition ad_t3dbc_im.F:2
ad_t3dbc_mod::ad_t3dbc_tile
subroutine, public ad_t3dbc_tile(ng, tile, itrc, ic, lbi, ubi, lbj, ubj, ubk, ubt, imins, imaxs, jmins, jmaxs, nstp, nout, ad_t)
Definition ad_t3dbc_im.F:57
mod_forces
Definition mod_forces.F:2
mod_forces::forces
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
mod_grid
Definition mod_grid.F:2
mod_grid::grid
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
mod_mixing
Definition mod_mixing.F:2
mod_mixing::mixing
type(t_mixing), dimension(:), allocatable mixing
Definition mod_mixing.F:399
mod_ocean
Definition mod_ocean.F:2
mod_ocean::ocean
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
mod_param
Definition mod_param.F:2
mod_param::ad_hadvection
type(t_adv), dimension(:,:), allocatable ad_hadvection
Definition mod_param.F:407
mod_param::nghostpoints
integer nghostpoints
Definition mod_param.F:710
mod_param::iadm
integer, parameter iadm
Definition mod_param.F:665
mod_param::domain
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
mod_param::vadvection
type(t_adv), dimension(:,:), allocatable vadvection
Definition mod_param.F:404
mod_param::ad_vadvection
type(t_adv), dimension(:,:), allocatable ad_vadvection
Definition mod_param.F:408
mod_param::npt
integer npt
Definition mod_param.F:505
mod_scalars
Definition mod_scalars.F:2
mod_scalars::luvsrc
logical, dimension(:), allocatable luvsrc
Definition mod_scalars.F:542
mod_scalars::ltracersrc
logical, dimension(:,:), allocatable ltracersrc
Definition mod_scalars.F:544
mod_scalars::iic
integer, dimension(:), allocatable iic
Definition mod_scalars.F:249
mod_scalars::dt
real(dp), dimension(:), allocatable dt
Definition mod_scalars.F:269
mod_scalars::lambda
real(r8) lambda
Definition mod_scalars.F:750
mod_scalars::ewperiodic
logical, dimension(:), allocatable ewperiodic
Definition mod_scalars.F:510
mod_scalars::iwest
integer, parameter iwest
Definition mod_scalars.F:1432
mod_scalars::nsperiodic
logical, dimension(:), allocatable nsperiodic
Definition mod_scalars.F:511
mod_scalars::compositegrid
logical, dimension(:,:), allocatable compositegrid
Definition mod_scalars.F:493
mod_scalars::itemp
integer itemp
Definition mod_scalars.F:126
mod_scalars::isouth
integer, parameter isouth
Definition mod_scalars.F:1433
mod_scalars::inert
integer, dimension(:), pointer inert
Definition mod_scalars.F:129
mod_scalars::ntfirst
integer, dimension(:), allocatable ntfirst
Definition mod_scalars.F:377
mod_scalars::ieast
integer, parameter ieast
Definition mod_scalars.F:1434
mod_scalars::ltracerclm
logical, dimension(:,:), allocatable ltracerclm
Definition mod_scalars.F:528
mod_scalars::inorth
integer, parameter inorth
Definition mod_scalars.F:1435
mod_scalars::lnudgetclm
logical, dimension(:,:), allocatable lnudgetclm
Definition mod_scalars.F:535
mod_sources
Definition mod_sources.F:2
mod_sources::sources
type(t_sources), dimension(:), allocatable sources
Definition mod_sources.F:90
mod_sources::nsrc
integer, dimension(:), allocatable nsrc
Definition mod_sources.F:97
mod_stepping
Definition mod_stepping.F:2
mod_stepping::kstp
integer, dimension(:), allocatable kstp
Definition mod_stepping.F:66
mod_stepping::knew
integer, dimension(:), allocatable knew
Definition mod_stepping.F:64
mod_stepping::nrhs
integer, dimension(:), allocatable nrhs
Definition mod_stepping.F:70
mod_stepping::nnew
integer, dimension(:), allocatable nnew
Definition mod_stepping.F:69
mod_stepping::nstp
integer, dimension(:), allocatable nstp
Definition mod_stepping.F:71
mp_exchange_mod
Definition mp_exchange.F:2
mp_exchange_mod::ad_mp_exchange4d
subroutine ad_mp_exchange4d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, lbt, ubt, nghost, ew_periodic, ns_periodic, ad_a, ad_b, ad_c)
Definition mp_exchange.F:7058
wclock_off
recursive subroutine wclock_off(ng, model, region, line, routine)
Definition timers.F:148
wclock_on
recursive subroutine wclock_on(ng, model, region, line, routine)
Definition timers.F:3