ROMS
Loading...
Searching...
No Matches
ad_uv3dmix4_s.h
Go to the documentation of this file.
1 MODULE ad_uv3dmix4_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 biharmonic mixing of momentum, !
11! along constant S-surfaces, from the horizontal divergence of the !
12! stress tensor. A transverse isotropy is assumed so the stress !
13! tensor is split into vertical and horizontal subtensors. !
14! !
15! Reference: !
16! !
17! Wajsowicz, R.C, 1993: A consistent formulation of the !
18! anisotropic stress tensor for use in models of the !
19! large-scale ocean circulation, JCP, 105, 333-338. !
20! !
21! Sadourny, R. and K. Maynard, 1997: Formulations of !
22! lateral diffusion in geophysical fluid dynamics !
23! models, In Numerical Methods of Atmospheric and !
24! Oceanic Modelling. Lin, Laprise, and Ritchie, !
25! Eds., NRC Research Press, 547-556. !
26! !
27! Griffies, S.M. and R.W. Hallberg, 2000: Biharmonic !
28! friction with a Smagorinsky-like viscosity for !
29! use in large-scale eddy-permitting ocean models, !
30! Monthly Weather Rev., 128, 8, 2935-2946. !
31! !
32! Basic state variables required: visc4, u, v, Hz. !
33! !
34!=======================================================================
35!
36 implicit none
37!
38 PRIVATE
39 PUBLIC ad_uv3dmix4
40!
41 CONTAINS
42!
43!***********************************************************************
44 SUBROUTINE ad_uv3dmix4 (ng, tile)
45!***********************************************************************
46!
47 USE mod_param
48 USE mod_ncparam
49 USE mod_coupling
50!!#ifdef DIAGNOSTICS_UV
51!! USE mod_diags
52!!#endif
53 USE mod_grid
54 USE mod_mixing
55 USE mod_ocean
56 USE mod_stepping
57!
58! Imported variable declarations.
59!
60 integer, intent(in) :: ng, tile
61!
62! Local variable declarations.
63!
64 character (len=*), parameter :: MyFile = &
65 & __FILE__
66!
67#include "tile.h"
68!
69#ifdef PROFILE
70 CALL wclock_on (ng, iadm, 32, __line__, myfile)
71#endif
72 CALL ad_uv3dmix4_s_tile (ng, tile, &
73 & lbi, ubi, lbj, ubj, &
74 & imins, imaxs, jmins, jmaxs, &
75 & nrhs(ng), nnew(ng), &
76#ifdef MASKING
77 & grid(ng) % pmask, &
78#endif
79 & grid(ng) % Hz, &
80 & grid(ng) % ad_Hz, &
81 & grid(ng) % om_p, &
82 & grid(ng) % om_r, &
83 & grid(ng) % on_p, &
84 & grid(ng) % on_r, &
85 & grid(ng) % pm, &
86 & grid(ng) % pmon_p, &
87 & grid(ng) % pmon_r, &
88 & grid(ng) % pn, &
89 & grid(ng) % pnom_p, &
90 & grid(ng) % pnom_r, &
91 & mixing(ng) % visc4_p, &
92 & mixing(ng) % visc4_r, &
93!!#ifdef DIAGNOSTICS_UV
94!! & DIAGS(ng) % DiaRUfrc, &
95!! & DIAGS(ng) % DiaRVfrc, &
96!! & DIAGS(ng) % DiaU3wrk, &
97!! & DIAGS(ng) % DiaV3wrk, &
98!!#endif
99 & ocean(ng) % u, &
100 & ocean(ng) % v, &
101 & coupling(ng) % ad_rufrc, &
102 & coupling(ng) % ad_rvfrc, &
103 & ocean(ng) % ad_u, &
104 & ocean(ng) % ad_v)
105#ifdef PROFILE
106 CALL wclock_off (ng, iadm, 32, __line__, myfile)
107#endif
108!
109 RETURN
110 END SUBROUTINE ad_uv3dmix4
111
112!
113!***********************************************************************
114 SUBROUTINE ad_uv3dmix4_s_tile (ng, tile, &
115 & LBi, UBi, LBj, UBj, &
116 & IminS, ImaxS, JminS, JmaxS, &
117 & nrhs, nnew, &
118#ifdef MASKING
119 & pmask, &
120#endif
121 & Hz, ad_Hz, &
122 & om_p, om_r, on_p, on_r, &
123 & pm, pmon_p, pmon_r, &
124 & pn, pnom_p, pnom_r, &
125 & visc4_p, visc4_r, &
126!!#ifdef DIAGNOSTICS_UV
127!! & DiaRUfrc, DiaRVfrc, &
128!! & DiaU3wrk, DiaV3wrk, &
129!!#endif
130 & u, v, &
131 & ad_rufrc, ad_rvfrc, ad_u, ad_v)
132!***********************************************************************
133!
134 USE mod_param
135 USE mod_scalars
136!
137! Imported variable declarations.
138!
139 integer, intent(in) :: ng, tile
140 integer, intent(in) :: LBi, UBi, LBj, UBj
141 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
142 integer, intent(in) :: nrhs, nnew
143
144#ifdef ASSUMED_SHAPE
145# ifdef MASKING
146 real(r8), intent(in) :: pmask(LBi:,LBj:)
147# endif
148 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
149 real(r8), intent(in) :: om_p(LBi:,LBj:)
150 real(r8), intent(in) :: om_r(LBi:,LBj:)
151 real(r8), intent(in) :: on_p(LBi:,LBj:)
152 real(r8), intent(in) :: on_r(LBi:,LBj:)
153 real(r8), intent(in) :: pm(LBi:,LBj:)
154 real(r8), intent(in) :: pmon_p(LBi:,LBj:)
155 real(r8), intent(in) :: pmon_r(LBi:,LBj:)
156 real(r8), intent(in) :: pn(LBi:,LBj:)
157 real(r8), intent(in) :: pnom_p(LBi:,LBj:)
158 real(r8), intent(in) :: pnom_r(LBi:,LBj:)
159 real(r8), intent(in) :: visc4_p(LBi:,LBj:)
160 real(r8), intent(in) :: visc4_r(LBi:,LBj:)
161
162 real(r8), intent(in) :: u(LBi:,LBj:,:,:)
163 real(r8), intent(in) :: v(LBi:,LBj:,:,:)
164
165!!# ifdef DIAGNOSTICS_UV
166!! real(r8), intent(inout) :: DiaRUfrc(LBi:,LBj:,:,:)
167!! real(r8), intent(inout) :: DiaRVfrc(LBi:,LBj:,:,:)
168!! real(r8), intent(inout) :: DiaU3wrk(LBi:,LBj:,:,:)
169!! real(r8), intent(inout) :: DiaV3wrk(LBi:,LBj:,:,:)
170!!# endif
171
172 real(r8), intent(inout) :: ad_Hz(LBi:,LBj:,:)
173 real(r8), intent(inout) :: ad_rufrc(LBi:,LBj:)
174 real(r8), intent(inout) :: ad_rvfrc(LBi:,LBj:)
175 real(r8), intent(inout) :: ad_u(LBi:,LBj:,:,:)
176 real(r8), intent(inout) :: ad_v(LBi:,LBj:,:,:)
177
178#else
179
180# ifdef MASKING
181 real(r8), intent(in) :: pmask(LBi:UBi,LBj:UBj)
182# endif
183 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
184 real(r8), intent(in) :: om_p(LBi:UBi,LBj:UBj)
185 real(r8), intent(in) :: om_r(LBi:UBi,LBj:UBj)
186 real(r8), intent(in) :: on_p(LBi:UBi,LBj:UBj)
187 real(r8), intent(in) :: on_r(LBi:UBi,LBj:UBj)
188 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
189 real(r8), intent(in) :: pmon_p(LBi:UBi,LBj:UBj)
190 real(r8), intent(in) :: pmon_r(LBi:UBi,LBj:UBj)
191 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
192 real(r8), intent(in) :: pnom_p(LBi:UBi,LBj:UBj)
193 real(r8), intent(in) :: pnom_r(LBi:UBi,LBj:UBj)
194 real(r8), intent(in) :: visc4_p(LBi:UBi,LBj:UBj)
195 real(r8), intent(in) :: visc4_r(LBi:UBi,LBj:UBj)
196
197 real(r8), intent(inout) :: u(LBi:UBi,LBj:UBj,N(ng),2)
198 real(r8), intent(inout) :: v(LBi:UBi,LBj:UBj,N(ng),2)
199
200!!# ifdef DIAGNOSTICS_UV
201!! real(r8), intent(inout) :: DiaRUfrc(LBi:UBi,LBj:UBj,3,NDM2d-1)
202!! real(r8), intent(inout) :: DiaRVfrc(LBi:UBi,LBj:UBj,3,NDM2d-1)
203!! real(r8), intent(inout) :: DiaU3wrk(LBi:UBi,LBj:UBj,N(ng),NDM3d)
204!! real(r8), intent(inout) :: DiaV3wrk(LBi:UBi,LBj:UBj,N(ng),NDM3d)
205!!# endif
206
207 real(r8), intent(inout) :: ad_Hz(LBi:UBi,LBj:UBj,N(ng))
208 real(r8), intent(inout) :: ad_rufrc(LBi:UBi,LBj:UBj)
209 real(r8), intent(inout) :: ad_rvfrc(LBi:UBi,LBj:UBj)
210 real(r8), intent(inout) :: ad_u(LBi:UBi,LBj:UBj,N(ng),2)
211 real(r8), intent(inout) :: ad_v(LBi:UBi,LBj:UBj,N(ng),2)
212#endif
213!
214! Local variable declarations.
215!
216 integer :: IminU, IminV, ImaxU, ImaxV
217 integer :: JminU, JminV, JmaxU, JmaxV
218 integer :: i, j, k
219
220 real(r8) :: cff, cff1, cff2
221 real(r8) :: ad_cff, ad_cff1, ad_cff2
222 real(r8) :: adfac, adfac1, adfac2, adfac3, adfac4
223
224 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: LapU
225 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: LapV
226 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: UFe
227 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: VFe
228 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: UFx
229 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: VFx
230
231 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_LapU
232 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_LapV
233 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_UFe
234 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_VFe
235 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_UFx
236 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: ad_VFx
237
238#include "set_bounds.h"
239!
240!-----------------------------------------------------------------------
241! Initialize adjoint private variables.
242!-----------------------------------------------------------------------
243!
244 ad_cff=0.0_r8
245 ad_cff1=0.0_r8
246 ad_cff2=0.0_r8
247
248 ad_lapu(imins:imaxs,jmins:jmaxs)=0.0_r8
249 ad_lapv(imins:imaxs,jmins:jmaxs)=0.0_r8
250
251 ad_ufe(imins:imaxs,jmins:jmaxs)=0.0_r8
252 ad_vfe(imins:imaxs,jmins:jmaxs)=0.0_r8
253 ad_ufx(imins:imaxs,jmins:jmaxs)=0.0_r8
254 ad_vfx(imins:imaxs,jmins:jmaxs)=0.0_r8
255!
256!-----------------------------------------------------------------------
257! Compute adjoint horizontal biharmonic viscosity along constant
258! S-surfaces. The biharmonic operator is computed by applying the
259! harmonic operator twice.
260!-----------------------------------------------------------------------
261!
262! Set local I- and J-ranges.
263!
264 IF (ewperiodic(ng)) THEN
265 iminu=istr-1
266 imaxu=iend+1
267 iminv=istr-1
268 imaxv=iend+1
269 ELSE
270 iminu=max(2,istru-1)
271 imaxu=min(iend+1,lm(ng))
272 iminv=max(1,istr-1)
273 imaxv=min(iend+1,lm(ng))
274 END IF
275 IF (nsperiodic(ng)) THEN
276 jminu=jstr-1
277 jmaxu=jend+1
278 jminv=jstr-1
279 jmaxv=jend+1
280 ELSE
281 jminu=max(1,jstr-1)
282 jmaxu=min(jend+1,mm(ng))
283 jminv=max(2,jstrv-1)
284 jmaxv=min(jend+1,mm(ng))
285 END IF
286!
287! Computed required BASIC STATE flux-components of the horizontal
288! divergence of the stress tensor in XI- and ETA-directions.
289!
290 k_loop : DO k=1,n(ng)
291 DO j=-1+jminv,jmaxv
292 DO i=-1+iminu,imaxu
293 cff=visc4_r(i,j)*0.5_r8* &
294 & (pmon_r(i,j)* &
295 & ((pn(i ,j)+pn(i+1,j))*u(i+1,j,k,nrhs)- &
296 & (pn(i-1,j)+pn(i ,j))*u(i ,j,k,nrhs))- &
297 & pnom_r(i,j)* &
298 & ((pm(i,j )+pm(i,j+1))*v(i,j+1,k,nrhs)- &
299 & (pm(i,j-1)+pm(i,j ))*v(i,j ,k,nrhs)))
300 ufx(i,j)=on_r(i,j)*on_r(i,j)*cff
301 vfe(i,j)=om_r(i,j)*om_r(i,j)*cff
302 END DO
303 END DO
304 DO j=jminu,jmaxu+1
305 DO i=iminv,imaxv+1
306 cff=visc4_p(i,j)*0.5_r8* &
307 & (pmon_p(i,j)* &
308 & ((pn(i ,j-1)+pn(i ,j))*v(i ,j,k,nrhs)- &
309 & (pn(i-1,j-1)+pn(i-1,j))*v(i-1,j,k,nrhs))+ &
310 & pnom_p(i,j)* &
311 & ((pm(i-1,j )+pm(i,j ))*u(i,j ,k,nrhs)- &
312 & (pm(i-1,j-1)+pm(i,j-1))*u(i,j-1,k,nrhs)))
313#ifdef MASKING
314 cff=cff*pmask(i,j)
315#endif
316 ufe(i,j)=om_p(i,j)*om_p(i,j)*cff
317 vfx(i,j)=on_p(i,j)*on_p(i,j)*cff
318 END DO
319 END DO
320!
321! Compute BASIC STATE first harmonic operator (m s^-3/2).
322!
323 DO j=jminu,jmaxu
324 DO i=iminu,imaxu
325 lapu(i,j)=0.125_r8* &
326 & (pm(i-1,j)+pm(i,j))*(pn(i-1,j)+pn(i,j))* &
327 & ((pn(i-1,j)+pn(i,j))*(ufx(i,j )-ufx(i-1,j))+ &
328 & (pm(i-1,j)+pm(i,j))*(ufe(i,j+1)-ufe(i ,j)))
329 END DO
330 END DO
331 DO j=jminv,jmaxv
332 DO i=iminv,imaxv
333 lapv(i,j)=0.125_r8* &
334 & (pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))* &
335 & ((pn(i,j-1)+pn(i,j))*(vfx(i+1,j)-vfx(i,j ))- &
336 & (pm(i,j-1)+pm(i,j))*(vfe(i ,j)-vfe(i,j-1)))
337 END DO
338 END DO
339!
340! Apply boundary conditions (other than periodic) to the BASIC STATE
341! first harmonic operator. These are gradient or closed (free slip or
342! no slip) boundary conditions.
343!
344 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
345 IF (domain(ng)%Western_Edge(tile)) THEN
346 IF (ad_lbc(iwest,isuvel,ng)%closed) THEN
347 DO j=jminu,jmaxu
348 lapu(istr,j)=0.0_r8
349 END DO
350 ELSE
351 DO j=jminu,jmaxu
352 lapu(istr,j)=lapu(istr+1,j)
353 END DO
354 END IF
355 IF (ad_lbc(iwest,isvvel,ng)%closed) THEN
356 DO j=jminv,jmaxv
357 lapv(istr-1,j)=gamma2(ng)*lapv(istr,j)
358 END DO
359 ELSE
360 DO j=jminv,jmaxv
361 lapv(istr-1,j)=0.0_r8
362 END DO
363 END IF
364 END IF
365 END IF
366!
367 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
368 IF (domain(ng)%Eastern_Edge(tile)) THEN
369 IF (ad_lbc(ieast,isuvel,ng)%closed) THEN
370 DO j=jminu,jmaxu
371 lapu(iend+1,j)=0.0_r8
372 END DO
373 ELSE
374 DO j=jminu,jmaxu
375 lapu(iend+1,j)=lapu(iend,j)
376 END DO
377 END IF
378 IF (ad_lbc(ieast,isvvel,ng)%closed) THEN
379 DO j=jminv,jmaxv
380 lapv(iend+1,j)=gamma2(ng)*lapv(iend,j)
381 END DO
382 ELSE
383 DO j=jminv,jmaxv
384 lapv(iend+1,j)=0.0_r8
385 END DO
386 END IF
387 END IF
388 END IF
389!
390 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
391 IF (domain(ng)%Southern_Edge(tile)) THEN
392 IF (ad_lbc(isouth,isuvel,ng)%closed) THEN
393 DO i=iminu,imaxu
394 lapu(i,jstr-1)=gamma2(ng)*lapu(i,jstr)
395 END DO
396 ELSE
397 DO i=iminu,imaxu
398 lapu(i,jstr-1)=0.0_r8
399 END DO
400 END IF
401 IF (ad_lbc(isouth,isvvel,ng)%closed) THEN
402 DO i=iminv,imaxv
403 lapv(i,jstr)=0.0_r8
404 END DO
405 ELSE
406 DO i=iminv,imaxv
407 lapv(i,jstr)=lapv(i,jstr+1)
408 END DO
409 END IF
410 END IF
411 END IF
412!
413 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
414 IF (domain(ng)%Northern_Edge(tile)) THEN
415 IF (ad_lbc(inorth,isuvel,ng)%closed) THEN
416 DO i=iminu,imaxu
417 lapu(i,jend+1)=gamma2(ng)*lapu(i,jend)
418 END DO
419 ELSE
420 DO i=iminu,imaxu
421 lapu(i,jend+1)=0.0_r8
422 END DO
423 END IF
424 IF (ad_lbc(inorth,isvvel,ng)%closed) THEN
425 DO i=iminv,imaxv
426 lapv(i,jend+1)=0.0_r8
427 END DO
428 ELSE
429 DO i=iminv,imaxv
430 lapv(i,jend+1)=lapv(i,jend)
431 END DO
432 END IF
433 END IF
434 END IF
435!
436 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng).or. &
437 & compositegrid(iwest ,ng).or.ewperiodic(ng))) THEN
438 IF (domain(ng)%SouthWest_Corner(tile)) THEN
439 lapu(istr ,jstr-1)=0.5_r8* &
440 & (lapu(istr+1,jstr-1)+ &
441 & lapu(istr ,jstr ))
442 lapv(istr-1,jstr )=0.5_r8* &
443 & (lapv(istr-1,jstr+1)+ &
444 & lapv(istr ,jstr ))
445 END IF
446 END IF
447
448 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng).or. &
449 & compositegrid(ieast ,ng).or.ewperiodic(ng))) THEN
450 IF (domain(ng)%SouthEast_Corner(tile)) THEN
451 lapu(iend+1,jstr-1)=0.5_r8* &
452 & (lapu(iend ,jstr-1)+ &
453 & lapu(iend+1,jstr ))
454 lapv(iend+1,jstr )=0.5_r8* &
455 & (lapv(iend ,jstr )+ &
456 & lapv(iend+1,jstr+1))
457 END IF
458 END IF
459
460 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng).or. &
461 & compositegrid(iwest ,ng).or.ewperiodic(ng))) THEN
462 IF (domain(ng)%NorthWest_Corner(tile)) THEN
463 lapu(istr ,jend+1)=0.5_r8* &
464 & (lapu(istr+1,jend+1)+ &
465 & lapu(istr ,jend ))
466 lapv(istr-1,jend+1)=0.5_r8* &
467 & (lapv(istr ,jend+1)+ &
468 & lapv(istr-1,jend ))
469 END IF
470 END IF
471
472 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng).or. &
473 & compositegrid(ieast ,ng).or.ewperiodic(ng))) THEN
474 IF (domain(ng)%NorthEast_Corner(tile)) THEN
475 lapu(iend+1,jend+1)=0.5_r8* &
476 & (lapu(iend ,jend+1)+ &
477 & lapu(iend+1,jend ))
478 lapv(iend+1,jend+1)=0.5_r8* &
479 & (lapv(iend ,jend+1)+ &
480 & lapv(iend+1,jend ))
481 END IF
482 END IF
483!
484! Time-step adjoint biharmonic, S-surfaces viscosity term. Notice that
485! momentum at this stage is HzU and HzV and has units m2/s.
486!
487 DO j=jstrv,jend
488 DO i=istr,iend
489 cff=0.25_r8*(pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))
490!!#ifdef DIAGNOSTICS_UV
491!! DiaRVfrc(i,j,3,M2hvis)=DiaRVfrc(i,j,3,M2hvis)-cff1
492!! DiaV3wrk(i,j,k,M3hvis)=-cff2
493!!#endif
494!^ tl_v(i,j,k,nnew)=tl_v(i,j,k,nnew)-tl_cff2
495!^
496 ad_cff2=ad_cff2-ad_v(i,j,k,nnew)
497!^ tl_rvfrc(i,j)=tl_rvfrc(i,j)-tl_cff1
498!^
499 ad_cff1=ad_cff1-ad_rvfrc(i,j)
500!^ tl_cff2=dt(ng)*cff*tl_cff1
501!^
502 ad_cff1=ad_cff1+dt(ng)*cff*ad_cff2
503 ad_cff2=0.0_r8
504!^ tl_cff1=0.5_r8*((pn(i,j-1)+pn(i,j))* &
505!^ & (tl_VFx(i+1,j)-tl_VFx(i,j ))- &
506!^ & (pm(i,j-1)+pm(i,j))* &
507!^ & (tl_VFe(i ,j)-tl_VFe(i,j-1)))
508!^
509 adfac=0.5_r8*ad_cff1
510 adfac1=adfac*(pn(i,j-1)+pn(i,j))
511 adfac2=adfac*(pm(i,j-1)+pm(i,j))
512 ad_vfx(i ,j )=ad_vfx(i ,j )-adfac1
513 ad_vfx(i+1,j )=ad_vfx(i+1,j )+adfac1
514 ad_vfe(i ,j-1)=ad_vfe(i ,j-1)+adfac2
515 ad_vfe(i ,j )=ad_vfe(i ,j )-adfac2
516 ad_cff1=0.0_r8
517 END DO
518 END DO
519 DO j=jstr,jend
520 DO i=istru,iend
521 cff=0.25_r8*(pm(i-1,j)+pm(i,j))*(pn(i-1,j)+pn(i,j))
522!!#ifdef DIAGNOSTICS_UV
523!! DiaRUfrc(i,j,3,M2hvis)=DiaRUfrc(i,j,3,M2hvis)-cff1
524!! DiaU3wrk(i,j,k,M3hvis)=-cff2
525!!#endif
526!^ tl_u(i,j,k,nnew)=tl_u(i,j,k,nnew)-tl_cff2
527!^
528 ad_cff2=ad_cff2-ad_u(i,j,k,nnew)
529!^ tl_rufrc(i,j)=tl_rufrc(i,j)-tl_cff1
530!^
531 ad_cff1=ad_cff1-ad_rufrc(i,j)
532!^ tl_cff2=dt(ng)*cff*tl_cff1
533!^
534 ad_cff1=ad_cff1+dt(ng)*cff*ad_cff2
535 ad_cff2=0.0_r8
536!^ tl_cff1=0.5_r8*((pn(i-1,j)+pn(i,j))* &
537!^ & (tl_UFx(i,j )-tl_UFx(i-1,j))+ &
538!^ & (pm(i-1,j)+pm(i,j))* &
539!^ & (tl_UFe(i,j+1)-tl_UFe(i ,j)))
540!^
541 adfac=0.5_r8*ad_cff1
542 adfac1=adfac*(pn(i-1,j)+pn(i,j))
543 adfac2=adfac*(pm(i-1,j)+pm(i,j))
544 ad_ufx(i-1,j )=ad_ufx(i-1,j )-adfac1
545 ad_ufx(i ,j )=ad_ufx(i ,j )+adfac1
546 ad_ufe(i ,j )=ad_ufe(i ,j )-adfac2
547 ad_ufe(i ,j+1)=ad_ufe(i ,j+1)+adfac2
548 ad_cff1=0.0_r8
549 END DO
550 END DO
551!
552! Compute flux-components of the adjoint horizontal divergence of the
553! harmonic stress tensor (m4/s2) in XI- and ETA-directions.
554!
555 DO j=jstr,jend+1
556 DO i=istr,iend+1
557!^ tl_VFx(i,j)=on_p(i,j)*on_p(i,j)*tl_cff
558!^
559 ad_cff=ad_cff+on_p(i,j)*on_p(i,j)*ad_vfx(i,j)
560 ad_vfx(i,j)=0.0_r8
561!^ tl_UFe(i,j)=om_p(i,j)*om_p(i,j)*tl_cff
562!^
563 ad_cff=ad_cff+om_p(i,j)*om_p(i,j)*ad_ufe(i,j)
564 ad_ufe(i,j)=0.0_r8
565#ifdef MASKING
566!^ tl_cff=tl_cff*pmask(i,j)
567!^
568 ad_cff=ad_cff*pmask(i,j)
569#endif
570!^ tl_cff=visc4_p(i,j)*0.125_r8* &
571!^ & ((tl_Hz(i-1,j ,k)+tl_Hz(i,j ,k)+ &
572!^ & tl_Hz(i-1,j-1,k)+tl_Hz(i,j-1,k))* &
573!^ & (pmon_p(i,j)* &
574!^ & ((pn(i ,j-1)+pn(i ,j))*LapV(i ,j)- &
575!^ & (pn(i-1,j-1)+pn(i-1,j))*LapV(i-1,j))+ &
576!^ & pnom_p(i,j)* &
577!^ & ((pm(i-1,j )+pm(i,j ))*LapU(i,j )- &
578!^ & (pm(i-1,j-1)+pm(i,j-1))*LapU(i,j-1)))+ &
579!^ & (Hz(i-1,j ,k)+Hz(i,j ,k)+ &
580!^ & Hz(i-1,j-1,k)+Hz(i,j-1,k))* &
581!^ & (pmon_p(i,j)* &
582!^ & ((pn(i ,j-1)+pn(i ,j))*tl_LapV(i ,j)- &
583!^ & (pn(i-1,j-1)+pn(i-1,j))*tl_LapV(i-1,j))+ &
584!^ & pnom_p(i,j)* &
585!^ & ((pm(i-1,j )+pm(i,j ))*tl_LapU(i,j )- &
586!^ & (pm(i-1,j-1)+pm(i,j-1))*tl_LapU(i,j-1))))
587!^
588 adfac=visc4_p(i,j)*0.125_r8*ad_cff
589 adfac1=adfac*(pmon_p(i,j)* &
590 & ((pn(i ,j-1)+pn(i ,j))*lapv(i ,j)- &
591 & (pn(i-1,j-1)+pn(i-1,j))*lapv(i-1,j))+ &
592 & pnom_p(i,j)* &
593 & ((pm(i-1,j )+pm(i,j ))*lapu(i,j )- &
594 & (pm(i-1,j-1)+pm(i,j-1))*lapu(i,j-1)))
595 adfac2=adfac*(hz(i-1,j ,k)+hz(i,j ,k)+ &
596 & hz(i-1,j-1,k)+hz(i,j-1,k))
597 adfac3=adfac2*pmon_p(i,j)
598 adfac4=adfac2*pnom_p(i,j)
599 ad_hz(i-1,j-1,k)=ad_hz(i-1,j-1,k)+adfac1
600 ad_hz(i-1,j ,k)=ad_hz(i-1,j ,k)+adfac1
601 ad_hz(i ,j-1,k)=ad_hz(i ,j-1,k)+adfac1
602 ad_hz(i ,j ,k)=ad_hz(i ,j ,k)+adfac1
603 ad_lapv(i-1,j)=ad_lapv(i-1,j)- &
604 & (pn(i-1,j-1)+pn(i-1,j))*adfac3
605 ad_lapv(i ,j)=ad_lapv(i ,j)+ &
606 & (pn(i ,j-1)+pn(i ,j))*adfac3
607 ad_lapu(i,j-1)=ad_lapu(i,j-1)- &
608 & (pm(i-1,j-1)+pm(i,j-1))*adfac4
609 ad_lapu(i,j )=ad_lapu(i,j )+ &
610 & (pm(i-1,j )+pm(i,j ))*adfac4
611 ad_cff=0.0_r8
612 END DO
613 END DO
614!
615 DO j=jstrv-1,jend
616 DO i=istru-1,iend
617!^ tl_VFe(i,j)=om_r(i,j)*om_r(i,j)*tl_cff
618!^
619 ad_cff=ad_cff+om_r(i,j)*om_r(i,j)*ad_vfe(i,j)
620 ad_vfe(i,j)=0.0_r8
621!^ tl_UFx(i,j)=on_r(i,j)*on_r(i,j)*tl_cff
622!^
623 ad_cff=ad_cff+on_r(i,j)*on_r(i,j)*ad_ufx(i,j)
624 ad_ufx(i,j)=0.0_r8
625!^ tl_cff=visc4_r(i,j)*0.5_r8* &
626!^ & (tl_Hz(i,j,k)*
627!^ & (pmon_r(i,j)* &
628!^ & ((pn(i ,j)+pn(i+1,j))*LapU(i+1,j)- &
629!^ & (pn(i-1,j)+pn(i ,j))*LapU(i ,j))- &
630!^ & pnom_r(i,j)* &
631!^ & ((pm(i,j )+pm(i,j+1))*LapV(i,j+1)- &
632!^ & (pm(i,j-1)+pm(i,j ))*LapV(i,j )))+ &
633!^ & Hz(i,j,k)* &
634!^ & (pmon_r(i,j)* &
635!^ & ((pn(i ,j)+pn(i+1,j))*tl_LapU(i+1,j)- &
636!^ & (pn(i-1,j)+pn(i ,j))*tl_LapU(i ,j))- &
637!^ & pnom_r(i,j)* &
638!^ & ((pm(i,j )+pm(i,j+1))*tl_LapV(i,j+1)- &
639!^ & (pm(i,j-1)+pm(i,j ))*tl_LapV(i,j ))))
640!^
641 adfac=visc4_r(i,j)*0.5_r8*ad_cff
642 adfac1=adfac*hz(i,j,k)
643 adfac2=adfac1*pmon_r(i,j)
644 adfac3=adfac1*pnom_r(i,j)
645 ad_hz(i,j,k)=ad_hz(i,j,k)+ &
646 & (pmon_r(i,j)* &
647 & ((pn(i ,j)+pn(i+1,j))*lapu(i+1,j)- &
648 & (pn(i-1,j)+pn(i ,j))*lapu(i ,j))- &
649 & pnom_r(i,j)* &
650 & ((pm(i,j )+pm(i,j+1))*lapv(i,j+1)- &
651 & (pm(i,j-1)+pm(i,j ))*lapv(i,j )))*adfac
652 ad_lapu(i ,j)=ad_lapu(i ,j)- &
653 & (pn(i-1,j)+pn(i ,j))*adfac2
654 ad_lapu(i+1,j)=ad_lapu(i+1,j)+ &
655 & (pn(i ,j)+pn(i+1,j))*adfac2
656 ad_lapv(i,j )=ad_lapv(i,j )+ &
657 & (pm(i,j-1)+pm(i,j ))*adfac3
658 ad_lapv(i,j+1)=ad_lapv(i,j+1)- &
659 & (pm(i,j )+pm(i,j+1))*adfac3
660 ad_cff=0.0_r8
661 END DO
662 END DO
663!
664! Apply boundary conditions (other than periodic) to the first
665! adjoint harmonic operator. These are gradient or closed
666! (free slip or no slip) boundary conditions.
667!
668 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng).or. &
669 & compositegrid(ieast ,ng).or.ewperiodic(ng))) THEN
670 IF (domain(ng)%NorthEast_Corner(tile)) THEN
671!^ tl_LapV(Iend+1,Jend+1)=0.5_r8* &
672!^ & (tl_LapV(Iend ,Jend+1)+ &
673!^ & tl_LapV(Iend+1,Jend ))
674!^
675 adfac=0.5_r8*ad_lapv(iend+1,jend+1)
676 ad_lapv(iend ,jend+1)=ad_lapv(iend ,jend+1)+adfac
677 ad_lapv(iend+1,jend )=ad_lapv(iend+1,jend )+adfac
678 ad_lapv(iend+1,jend+1)=0.0_r8
679!^ tl_LapU(Iend+1,Jend+1)=0.5_r8* &
680!^ & (tl_LapU(Iend ,Jend+1)+ &
681!^ & tl_LapU(Iend+1,Jend ))
682!^
683 adfac=0.5_r8*ad_lapu(iend+1,jend+1)
684 ad_lapu(iend ,jend+1)=ad_lapu(iend ,jend+1)+adfac
685 ad_lapu(iend+1,jend )=ad_lapu(iend+1,jend )+adfac
686 ad_lapu(iend+1,jend+1)=0.0_r8
687 END IF
688 END IF
689
690 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng).or. &
691 & compositegrid(iwest ,ng).or.ewperiodic(ng))) THEN
692 IF (domain(ng)%NorthWest_Corner(tile)) THEN
693!^ tl_LapV(Istr-1,Jend+1)=0.5_r8* &
694!^ & (tl_LapV(Istr ,Jend+1)+ &
695!^ & tl_LapV(Istr-1,Jend ))
696!^
697 adfac=0.5_r8*ad_lapv(istr-1,jend+1)
698 ad_lapv(istr ,jend+1)=ad_lapv(istr ,jend+1)+adfac
699 ad_lapv(istr-1,jend )=ad_lapv(istr-1,jend )+adfac
700 ad_lapv(istr-1,jend+1)=0.0_r8
701!^ tl_LapU(Istr ,Jend+1)=0.5_r8* &
702!^ & (tl_LapU(Istr+1,Jend+1)+ &
703!^ & tl_LapU(Istr ,Jend ))
704!^
705 adfac=0.5_r8*ad_lapu(istr,jend+1)
706 ad_lapu(istr+1,jend+1)=ad_lapu(istr+1,jend+1)+adfac
707 ad_lapu(istr ,jend )=ad_lapu(istr ,jend )+adfac
708 ad_lapu(istr ,jend+1)=0.0_r8
709 END IF
710 END IF
711
712 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng).or. &
713 & compositegrid(ieast ,ng).or.ewperiodic(ng))) THEN
714 IF (domain(ng)%SouthEast_Corner(tile)) THEN
715!^ tl_LapV(Iend+1,Jstr )=0.5_r8* &
716!^ & (tl_LapV(Iend ,Jstr )+ &
717!^ & tl_LapV(Iend+1,Jstr+1))
718!^
719 adfac=0.5_r8*ad_lapv(iend+1,jstr)
720 ad_lapv(iend ,jstr )=ad_lapv(iend ,jstr )+adfac
721 ad_lapv(iend+1,jstr+1)=ad_lapv(iend+1,jstr+1)+adfac
722 ad_lapv(iend+1,jstr )=0.0_r8
723!^ tl_LapU(Iend+1,Jstr-1)=0.5_r8* &
724!^ & (tl_LapU(Iend ,Jstr-1)+ &
725!^ & tl_LapU(Iend+1,Jstr ))
726!^
727 adfac=0.5_r8*ad_lapu(iend+1,jstr-1)
728 ad_lapu(iend ,jstr-1)=ad_lapu(iend ,jstr-1)+adfac
729 ad_lapu(iend+1,jstr )=ad_lapu(iend+1,jstr )+adfac
730 ad_lapu(iend+1,jstr-1)=0.0_r8
731 END IF
732 END IF
733
734 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng).or. &
735 & compositegrid(iwest ,ng).or.ewperiodic(ng))) THEN
736 IF (domain(ng)%SouthWest_Corner(tile)) THEN
737!^ tl_LapV(Istr-1,Jstr )=0.5_r8* &
738!^ & (tl_LapV(Istr-1,Jstr+1)+ &
739!^ & tl_LapV(Istr ,Jstr ))
740!^
741 adfac=0.5_r8*ad_lapv(istr-1,jstr)
742 ad_lapv(istr-1,jstr+1)=ad_lapv(istr-1,jstr+1)+adfac
743 ad_lapv(istr ,jstr )=ad_lapv(istr ,jstr )+adfac
744 ad_lapv(istr-1,jstr )=0.0_r8
745!^ tl_LapU(Istr ,Jstr-1)=0.5_r8* &
746!^ & (tl_LapU(Istr+1,Jstr-1)+ &
747!^ & tl_LapU(Istr ,Jstr ))
748!^
749 adfac=0.5_r8*ad_lapu(istr,jstr-1)
750 ad_lapu(istr+1,jstr-1)=ad_lapu(istr+1,jstr-1)+adfac
751 ad_lapu(istr ,jstr )=ad_lapu(istr ,jstr )+adfac
752 ad_lapu(istr ,jstr-1)=0.0_r8
753 END IF
754 END IF
755!
756 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
757 IF (domain(ng)%Northern_Edge(tile)) THEN
758 IF (ad_lbc(inorth,isvvel,ng)%closed) THEN
759 DO i=iminv,imaxv
760!^ tl_LapV(i,Jend+1)=0.0_r8
761!^
762 ad_lapv(i,jend+1)=0.0_r8
763 END DO
764 ELSE
765 DO i=iminv,imaxv
766!^ tl_LapV(i,Jend+1)=tl_LapV(i,Jend)
767!^
768 ad_lapv(i,jend)=ad_lapv(i,jend)+ad_lapv(i,jend+1)
769 ad_lapv(i,jend+1)=0.0_r8
770 END DO
771 END IF
772 IF (ad_lbc(inorth,isuvel,ng)%closed) THEN
773 DO i=iminu,imaxu
774!^ tl_LapU(i,Jend+1)=gamma2(ng)*tl_LapU(i,Jend)
775!^
776 ad_lapu(i,jend)=ad_lapu(i,jend)+ &
777 & gamma2(ng)*ad_lapu(i,jend+1)
778 ad_lapu(i,jend+1)=0.0_r8
779 END DO
780 ELSE
781 DO i=iminu,imaxu
782!^ tl_LapU(i,Jend+1)=0.0_r8
783!^
784 ad_lapu(i,jend+1)=0.0_r8
785 END DO
786 END IF
787 END IF
788 END IF
789!
790 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
791 IF (domain(ng)%Southern_Edge(tile)) THEN
792 IF (ad_lbc(isouth,isvvel,ng)%closed) THEN
793 DO i=iminv,imaxv
794!^ tl_LapV(i,JstrV-1)=0.0_r8
795!^
796 ad_lapv(i,jstrv-1)=0.0_r8
797 END DO
798 ELSE
799 DO i=iminv,imaxv
800!^ tl_LapV(i,JstrV-1)=tl_LapV(i,JstrV)
801!^
802 ad_lapv(i,jstrv)=ad_lapv(i,jstrv)+ad_lapv(i,jstrv-1)
803 ad_lapv(i,jstrv-1)=0.0_r8
804 END DO
805 END IF
806 IF (ad_lbc(isouth,isuvel,ng)%closed) THEN
807 DO i=iminu,imaxu
808!^ tl_LapU(i,Jstr-1)=gamma2(ng)*tl_LapU(i,Jstr)
809!^
810 ad_lapu(i,jstr)=ad_lapu(i,jstr)+ &
811 & gamma2(ng)*ad_lapu(i,jstr-1)
812 ad_lapu(i,jstr-1)=0.0_r8
813 END DO
814 ELSE
815 DO i=iminu,imaxu
816!^ tl_LapU(i,Jstr-1)=0.0_r8
817!^
818 ad_lapu(i,jstr-1)=0.0_r8
819 END DO
820 END IF
821 END IF
822 END IF
823!
824 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
825 IF (domain(ng)%Eastern_Edge(tile)) THEN
826 IF (ad_lbc(ieast,isvvel,ng)%closed) THEN
827 DO j=jminv,jmaxv
828!^ tl_LapV(Iend+1,j)=gamma2(ng)*tl_LapV(Iend,j)
829!^
830 ad_lapv(iend,j)=ad_lapv(iend,j)+ &
831 & gamma2(ng)*ad_lapv(iend+1,j)
832 ad_lapv(iend+1,j)=0.0_r8
833 END DO
834 ELSE
835 DO j=jminv,jmaxv
836!^ tl_LapV(Iend+1,j)=0.0_r8
837!^
838 ad_lapv(iend+1,j)=0.0_r8
839 END DO
840 END IF
841 IF (ad_lbc(ieast,isuvel,ng)%closed) THEN
842 DO j=jminu,jmaxu
843!^ tl_LapU(Iend+1,j)=0.0_r8
844!^
845 ad_lapu(iend+1,j)=0.0_r8
846 END DO
847 ELSE
848 DO j=jminu,jmaxu
849!^ tl_LapU(Iend+1,j)=tl_LapU(Iend,j)
850!^
851 ad_lapu(iend,j)=ad_lapu(iend,j)+ad_lapu(iend+1,j)
852 ad_lapu(iend+1,j)=0.0_r8
853 END DO
854 END IF
855 END IF
856 END IF
857!
858 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
859 IF (domain(ng)%Western_Edge(tile)) THEN
860 IF (ad_lbc(iwest,isvvel,ng)%closed) THEN
861 DO j=jminv,jmaxv
862!^ tl_LapV(Istr-1,j)=gamma2(ng)*tl_LapV(Istr,j)
863!^
864 ad_lapv(istr,j)=ad_lapv(istr,j)+ &
865 & gamma2(ng)*ad_lapv(istr-1,j)
866 ad_lapv(istr-1,j)=0.0_r8
867 END DO
868 ELSE
869 DO j=jminv,jmaxv
870!^ tl_LapV(Istr-1,j)=0.0_r8
871!^
872 ad_lapv(istr-1,j)=0.0_r8
873 END DO
874 END IF
875 IF (ad_lbc(iwest,isuvel,ng)%closed) THEN
876 DO j=jminu,jmaxu
877!^ tl_LapU(IstrU-1,j)=0.0_r8
878!^
879 ad_lapu(istru-1,j)=0.0_r8
880 END DO
881 ELSE
882 DO j=jminu,jmaxu
883!^ tl_LapU(IstrU-1,j)=tl_LapU(IstrU,j)
884!^
885 ad_lapu(istru,j)=ad_lapu(istru,j)+ad_lapu(istru-1,j)
886 ad_lapu(istru-1,j)=0.0_r8
887 END DO
888 END IF
889 END IF
890 END IF
891!
892! Compute adjoint first harmonic operator (m s^-3/2).
893!
894 DO j=jminv,jmaxv
895 DO i=iminv,imaxv
896 cff=0.125_r8*(pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))
897!^ tl_LapV(i,j)=cff* &
898!^ & ((pn(i,j-1)+pn(i,j))* &
899!^ & (tl_VFx(i+1,j)-tl_VFx(i,j ))- &
900!^ & (pm(i,j-1)+pm(i,j))* &
901!^ & (tl_VFe(i ,j)-tl_VFe(i,j-1)))
902!^
903 adfac=cff*ad_lapv(i,j)
904 adfac1=adfac*(pn(i,j-1)+pn(i,j))
905 adfac2=adfac*(pm(i,j-1)+pm(i,j))
906 ad_vfx(i ,j)=ad_vfx(i ,j)-adfac1
907 ad_vfx(i+1,j)=ad_vfx(i+1,j)+adfac1
908 ad_vfe(i,j-1)=ad_vfe(i,j-1)+adfac2
909 ad_vfe(i, j)=ad_vfe(i, j)-adfac2
910 ad_lapv(i,j)=0.0_r8
911 END DO
912 END DO
913 DO j=jminu,jmaxu
914 DO i=iminu,imaxu
915 cff=0.125_r8*(pm(i-1,j)+pm(i,j))*(pn(i-1,j)+pn(i,j))
916!^ tl_LapU(i,j)=cff* &
917!^ & ((pn(i-1,j)+pn(i,j))* &
918!^ & (tl_UFx(i,j )-tl_UFx(i-1,j))+ &
919!^ & (pm(i-1,j)+pm(i,j))* &
920!^ & (tl_UFe(i,j+1)-tl_UFe(i ,j)))
921!^
922 adfac=cff*ad_lapu(i,j)
923 adfac1=adfac*(pn(i-1,j)+pn(i,j))
924 adfac2=adfac*(pm(i-1,j)+pm(i,j))
925 ad_ufx(i-1,j)=ad_ufx(i-1,j)-adfac1
926 ad_ufx(i ,j )=ad_ufx(i ,j)+adfac1
927 ad_ufe(i,j )=ad_ufe(i,j )-adfac2
928 ad_ufe(i,j+1)=ad_ufe(i,j+1)+adfac2
929 ad_lapu(i,j)=0.0_r8
930 END DO
931 END DO
932!
933! Compute flux-components of the horizontal divergence of the stress
934! tensor (m4 s^-3/2) in XI- and ETA-directions. It is assumed here
935! that "visc4_r" and "visc4_p" are the squared root of the biharmonic
936! viscosity coefficient. For momentum balance purposes, the
937! thickness "Hz" appears only when computing the second harmonic
938! operator.
939!
940 DO j=jminu,jmaxu+1
941 DO i=iminv,imaxv+1
942!^ tl_VFx(i,j)=on_p(i,j)*on_p(i,j)*tl_cff
943!^
944 ad_cff=ad_cff+on_p(i,j)*on_p(i,j)*ad_vfx(i,j)
945 ad_vfx(i,j)=0.0_r8
946!^ tl_UFe(i,j)=om_p(i,j)*om_p(i,j)*tl_cff
947!^
948 ad_cff=ad_cff+om_p(i,j)*om_p(i,j)*ad_ufe(i,j)
949 ad_ufe(i,j)=0.0_r8
950#ifdef MASKING
951!^ tl_cff=tl_cff*pmask(i,j)
952!^
953 ad_cff=ad_cff*pmask(i,j)
954#endif
955!^ tl_cff=visc4_p(i,j)*0.5_r8* &
956!^ & (pmon_p(i,j)* &
957!^ & ((pn(i ,j-1)+pn(i ,j))*tl_v(i ,j,k,nrhs)- &
958!^ & (pn(i-1,j-1)+pn(i-1,j))*tl_v(i-1,j,k,nrhs))+ &
959!^ & pnom_p(i,j)* &
960!^ & ((pm(i-1,j )+pm(i,j ))*tl_u(i,j ,k,nrhs)- &
961!^ & (pm(i-1,j-1)+pm(i,j-1))*tl_u(i,j-1,k,nrhs)))
962!^
963 adfac=visc4_p(i,j)*0.5_r8*ad_cff
964 adfac1=adfac*pmon_p(i,j)
965 adfac2=adfac*pnom_p(i,j)
966 ad_v(i-1,j,k,nrhs)=ad_v(i-1,j,k,nrhs)- &
967 & (pn(i-1,j-1)+pn(i-1,j))*adfac1
968 ad_v(i ,j,k,nrhs)=ad_v(i ,j,k,nrhs)+ &
969 & (pn(i ,j-1)+pn(i ,j))*adfac1
970 ad_u(i,j ,k,nrhs)=ad_u(i,j ,k,nrhs)+ &
971 & (pm(i-1,j )+pm(i,j ))*adfac2
972 ad_u(i,j-1,k,nrhs)=ad_u(i,j-1,k,nrhs)- &
973 & (pm(i-1,j-1)+pm(i,j-1))*adfac2
974 ad_cff=0.0_r8
975 END DO
976 END DO
977 DO j=-1+jminv,jmaxv
978 DO i=-1+iminu,imaxu
979!^ tl_VFe(i,j)=om_r(i,j)*om_r(i,j)*tl_cff
980!^
981 ad_cff=ad_cff+om_r(i,j)*om_r(i,j)*ad_vfe(i,j)
982 ad_vfe(i,j)=0.0_r8
983!^ tl_UFx(i,j)=on_r(i,j)*on_r(i,j)*tl_cff
984!^
985 ad_cff=ad_cff+on_r(i,j)*on_r(i,j)*ad_ufx(i,j)
986 ad_ufx(i,j)=0.0_r8
987!^ tl_cff=visc4_r(i,j)*0.5_r8* &
988!^ & (pmon_r(i,j)* &
989!^ & ((pn(i ,j)+pn(i+1,j))*tl_u(i+1,j,k,nrhs)- &
990!^ & (pn(i-1,j)+pn(i ,j))*tl_u(i ,j,k,nrhs))- &
991!^ & pnom_r(i,j)* &
992!^ & ((pm(i,j )+pm(i,j+1))*tl_v(i,j+1,k,nrhs)- &
993!^ & (pm(i,j-1)+pm(i,j ))*tl_v(i,j ,k,nrhs)))
994!^
995 adfac=visc4_r(i,j)*0.5_r8*ad_cff
996 adfac1=adfac*pmon_r(i,j)
997 adfac2=adfac*pnom_r(i,j)
998 ad_u(i ,j,k,nrhs)=ad_u(i ,j,k,nrhs)- &
999 & (pn(i-1,j)+pn(i ,j))*adfac1
1000 ad_u(i+1,j,k,nrhs)=ad_u(i+1,j,k,nrhs)+ &
1001 & (pn(i ,j)+pn(i+1,j))*adfac1
1002 ad_v(i,j ,k,nrhs)=ad_v(i,j ,k,nrhs)+ &
1003 & (pm(i,j-1)+pm(i,j ))*adfac2
1004 ad_v(i,j+1,k,nrhs)=ad_v(i,j+1,k,nrhs)- &
1005 & (pm(i,j )+pm(i,j+1))*adfac2
1006 ad_cff=0.0_r8
1007 END DO
1008 END DO
1009 END DO k_loop
1010!
1011 RETURN
1012 END SUBROUTINE ad_uv3dmix4_s_tile
1013
1014 END MODULE ad_uv3dmix4_mod
subroutine, public ad_uv3dmix4(ng, tile)
subroutine ad_uv3dmix4_s_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, nnew, pmask, hz, ad_hz, om_p, om_r, on_p, on_r, pm, pmon_p, pmon_r, pn, pnom_p, pnom_r, visc4_p, visc4_r, u, v, ad_rufrc, ad_rvfrc, ad_u, ad_v)
type(t_coupling), dimension(:), allocatable coupling
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
type(t_lbc), dimension(:,:,:), allocatable ad_lbc
Definition mod_param.F:378
integer, parameter iadm
Definition mod_param.F:665
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer, dimension(:), allocatable lm
Definition mod_param.F:455
integer, dimension(:), allocatable mm
Definition mod_param.F:456
real(dp), dimension(:), allocatable dt
logical, dimension(:), allocatable ewperiodic
integer, parameter iwest
logical, dimension(:), allocatable nsperiodic
real(r8), dimension(:), allocatable gamma2
logical, dimension(:,:), allocatable compositegrid
integer, parameter isouth
integer, parameter ieast
integer, parameter inorth
integer, dimension(:), allocatable nrhs
integer, dimension(:), allocatable nnew
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