ROMS
Loading...
Searching...
No Matches
tl_uv3dmix4_geo.h
Go to the documentation of this file.
1 MODULE tl_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 tangent linear biharmonic mixing of !
11! momentum, rotated along geopotentials, from the horizontal !
12! divergence of the stress tensor. A transverse isotropy is !
13! assumed so the stress tensor is split into vertical and !
14! horizontal subtensors. !
15! !
16! Reference: !
17! !
18! Wajsowicz, R.C, 1993: A consistent formulation of the !
19! anisotropic stress tensor for use in models of the !
20! large-scale ocean circulation, JCP, 105, 333-338. !
21! !
22! Sadourny, R. and K. Maynard, 1997: Formulations of !
23! lateral diffusion in geophysical fluid dynamics !
24! models, In Numerical Methods of Atmospheric and !
25! Oceanic Modelling. Lin, Laprise, and Ritchie, !
26! Eds., NRC Research Press, 547-556. !
27! !
28! Griffies, S.M. and R.W. Hallberg, 2000: Biharmonic !
29! friction with a Smagorinsky-like viscosity for !
30! use in large-scale eddy-permitting ocean models, !
31! Monthly Weather Rev., 128, 8, 2935-2946. !
32! !
33!=======================================================================
34!
35 implicit none
36!
37 PRIVATE
38 PUBLIC tl_uv3dmix4
39!
40 CONTAINS
41!
42!***********************************************************************
43 SUBROUTINE tl_uv3dmix4 (ng, tile)
44!***********************************************************************
45!
46 USE mod_param
47 USE mod_coupling
48#ifdef DIAGNOSTICS_UV
49!! USE mod_diags
50#endif
51 USE mod_grid
52 USE mod_mixing
53 USE mod_ocean
54 USE mod_stepping
55!
56! Imported variable declarations.
57!
58 integer, intent(in) :: ng, tile
59!
60! Local variable declarations.
61!
62 character (len=*), parameter :: myfile = &
63 & __FILE__
64!
65#include "tile.h"
66!
67#ifdef PROFILE
68 CALL wclock_on (ng, itlm, 33, __line__, myfile)
69#endif
70 CALL tl_uv3dmix4_geo_tile (ng, tile, &
71 & lbi, ubi, lbj, ubj, &
72 & imins, imaxs, jmins, jmaxs, &
73 & nrhs(ng), nnew(ng), &
74#ifdef MASKING
75 & grid(ng) % pmask, &
76 & grid(ng) % rmask, &
77 & grid(ng) % umask, &
78 & grid(ng) % vmask, &
79#endif
80 & grid(ng) % om_p, &
81 & grid(ng) % om_r, &
82 & grid(ng) % om_u, &
83 & grid(ng) % om_v, &
84 & grid(ng) % on_p, &
85 & grid(ng) % on_r, &
86 & grid(ng) % on_u, &
87 & grid(ng) % on_v, &
88 & grid(ng) % pm, &
89 & grid(ng) % pn, &
90 & grid(ng) % Hz, &
91 & grid(ng) % tl_Hz, &
92 & grid(ng) % z_r, &
93 & grid(ng) % tl_z_r, &
94#ifdef VISC_3DCOEF
95# ifdef UV_U3ADV_SPLIT
96 & mixing(ng) % Uvis3d_r, &
97 & mixing(ng) % Vvis3d_r, &
98 & mixing(ng) % tl_Uvis3d_r, &
99 & mixing(ng) % tl_Vvis3d_r, &
100# else
101 & mixing(ng) % visc3d_r, &
102 & mixing(ng) % tl_visc3d_r, &
103# endif
104#else
105 & mixing(ng) % visc4_p, &
106 & mixing(ng) % visc4_r, &
107#endif
108#ifdef DIAGNOSTICS_UV
109!! & DIAGS(ng) % DiaRUfrc, &
110!! & DIAGS(ng) % DiaRVfrc, &
111!! & DIAGS(ng) % DiaU3wrk, &
112!! & DIAGS(ng) % DiaV3wrk, &
113#endif
114 & ocean(ng) % u, &
115 & ocean(ng) % v, &
116 & ocean(ng) % tl_u, &
117 & ocean(ng) % tl_v, &
118 & coupling(ng) % tl_rufrc, &
119 & coupling(ng) % tl_rvfrc)
120#ifdef PROFILE
121 CALL wclock_off (ng, itlm, 33, __line__, myfile)
122#endif
123!
124 RETURN
125 END SUBROUTINE tl_uv3dmix4
126!
127!***********************************************************************
128 SUBROUTINE tl_uv3dmix4_geo_tile (ng, tile, &
129 & LBi, UBi, LBj, UBj, &
130 & IminS, ImaxS, JminS, JmaxS, &
131 & nrhs, nnew, &
132#ifdef MASKING
133 & pmask, rmask, umask, vmask, &
134#endif
135 & om_p, om_r, om_u, om_v, &
136 & on_p, on_r, on_u, on_v, &
137 & pm, pn, &
138 & Hz, tl_Hz, &
139 & z_r, tl_z_r, &
140#ifdef VISC_3DCOEF
141# ifdef UV_U3ADV_SPLIT
142 & Uvis3d_r, Vvis3d_r, &
143 & tl_Uvis3d_r, tl_Vvis3d_r, &
144# else
145 & visc3d_r, tl_visc3d_r, &
146# endif
147#else
148 & visc4_p, visc4_r, &
149#endif
150#ifdef DIAGNOSTICS_UV
151!! & DiaRUfrc, DiaRVfrc, &
152!! & DiaU3wrk, DiaV3wrk, &
153#endif
154 & u, v, &
155 & tl_u, tl_v, &
156 & tl_rufrc, tl_rvfrc)
157!***********************************************************************
158!
159 USE mod_param
160 USE mod_ncparam
161 USE mod_scalars
162!
163! Imported variable declarations.
164!
165 integer, intent(in) :: ng, tile
166 integer, intent(in) :: LBi, UBi, LBj, UBj
167 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
168 integer, intent(in) :: nrhs, nnew
169
170#ifdef ASSUMED_SHAPE
171# ifdef MASKING
172 real(r8), intent(in) :: pmask(LBi:,LBj:)
173 real(r8), intent(in) :: rmask(LBi:,LBj:)
174 real(r8), intent(in) :: umask(LBi:,LBj:)
175 real(r8), intent(in) :: vmask(LBi:,LBj:)
176# endif
177 real(r8), intent(in) :: om_p(LBi:,LBj:)
178 real(r8), intent(in) :: om_r(LBi:,LBj:)
179 real(r8), intent(in) :: om_u(LBi:,LBj:)
180 real(r8), intent(in) :: om_v(LBi:,LBj:)
181 real(r8), intent(in) :: on_p(LBi:,LBj:)
182 real(r8), intent(in) :: on_r(LBi:,LBj:)
183 real(r8), intent(in) :: on_u(LBi:,LBj:)
184 real(r8), intent(in) :: on_v(LBi:,LBj:)
185 real(r8), intent(in) :: pm(LBi:,LBj:)
186 real(r8), intent(in) :: pn(LBi:,LBj:)
187 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
188 real(r8), intent(in) :: tl_Hz(LBi:,LBj:,:)
189 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
190 real(r8), intent(in) :: tl_z_r(LBi:,LBj:,:)
191# ifdef VISC_3DCOEF
192# ifdef UV_U3ADV_SPLIT
193 real(r8), intent(in) :: Uvis3d_r(LBi:,LBj:,:)
194 real(r8), intent(in) :: Vvis3d_r(LBi:,LBj:,:)
195# else
196 real(r8), intent(in) :: visc3d_r(LBi:,LBj:,:)
197# endif
198# else
199 real(r8), intent(in) :: visc4_p(LBi:,LBj:)
200 real(r8), intent(in) :: visc4_r(LBi:,LBj:)
201# endif
202 real(r8), intent(in) :: u(LBi:,LBj:,:,:)
203 real(r8), intent(in) :: v(LBi:,LBj:,:,:)
204
205# ifdef DIAGNOSTICS_UV
206!! real(r8), intent(inout) :: DiaRUfrc(LBi:,LBj:,:,:)
207!! real(r8), intent(inout) :: DiaRVfrc(LBi:,LBj:,:,:)
208!! real(r8), intent(inout) :: DiaU3wrk(LBi:,LBj:,:,:)
209!! real(r8), intent(inout) :: DiaV3wrk(LBi:,LBj:,:,:)
210# endif
211# ifdef VISC_3DCOEF
212# ifdef UV_U3ADV_SPLIT
213 real(r8), intent(inout) :: tl_Uvis3d_r(LBi:,LBj:,:)
214 real(r8), intent(inout) :: tl_Vvis3d_r(LBi:,LBj:,:)
215# else
216 real(r8), intent(inout) :: tl_visc3d_r(LBi:,LBj:,:)
217# endif
218# endif
219 real(r8), intent(inout) :: tl_rufrc(LBi:,LBj:)
220 real(r8), intent(inout) :: tl_rvfrc(LBi:,LBj:)
221 real(r8), intent(inout) :: tl_u(LBi:,LBj:,:,:)
222 real(r8), intent(inout) :: tl_v(LBi:,LBj:,:,:)
223#else
224# ifdef MASKING
225 real(r8), intent(in) :: pmask(LBi:UBi,LBj:UBj)
226 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
227 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
228 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
229# endif
230 real(r8), intent(in) :: om_p(LBi:UBi,LBj:UBj)
231 real(r8), intent(in) :: om_r(LBi:UBi,LBj:UBj)
232 real(r8), intent(in) :: om_u(LBi:UBi,LBj:UBj)
233 real(r8), intent(in) :: om_v(LBi:UBi,LBj:UBj)
234 real(r8), intent(in) :: on_p(LBi:UBi,LBj:UBj)
235 real(r8), intent(in) :: on_r(LBi:UBi,LBj:UBj)
236 real(r8), intent(in) :: on_u(LBi:UBi,LBj:UBj)
237 real(r8), intent(in) :: on_v(LBi:UBi,LBj:UBj)
238 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
239 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
240 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
241 real(r8), intent(in) :: tl_Hz(LBi:UBi,LBj:UBj,N(ng))
242 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
243 real(r8), intent(in) :: tl_z_r(LBi:UBi,LBj:UBj,N(ng))
244# ifdef VISC_3DCOEF
245# ifdef UV_U3ADV_SPLIT
246 real(r8), intent(in) :: Uvis3d_r(LBi:UBi,LBj:UBj,N(ng))
247 real(r8), intent(in) :: Vvis3d_r(LBi:UBi,LBj:UBj,N(ng))
248# else
249 real(r8), intent(in) :: visc3d_r(LBi:UBi,LBj:UBj,N(ng))
250# endif
251# else
252 real(r8), intent(in) :: visc4_p(LBi:UBi,LBj:UBj)
253 real(r8), intent(in) :: visc4_r(LBi:UBi,LBj:UBj)
254# endif
255 real(r8), intent(in) :: u(LBi:UBi,LBj:UBj,N(ng),2)
256 real(r8), intent(in) :: v(LBi:UBi,LBj:UBj,N(ng),2)
257
258# ifdef DIAGNOSTICS_UV
259!! real(r8), intent(inout) :: DiaRUfrc(LBi:UBi,LBj:UBj,3,NDM2d-1)
260!! real(r8), intent(inout) :: DiaRVfrc(LBi:UBi,LBj:UBj,3,NDM2d-1)
261!! real(r8), intent(inout) :: DiaU3wrk(LBi:UBi,LBj:UBj,N(ng),NDM3d)
262!! real(r8), intent(inout) :: DiaV3wrk(LBi:UBi,LBj:UBj,N(ng),NDM3d)
263# endif
264# ifdef VISC_3DCOEF
265# ifdef UV_U3ADV_SPLIT
266 real(r8), intent(inout) :: tl_Uvis3d_r(LBi:UBi,LBj:UBj,N(ng))
267 real(r8), intent(inout) :: tl_Vvis3d_r(LBi:UBi,LBj:UBj,N(ng))
268# else
269 real(r8), intent(inout) :: tl_visc3d_r(LBi:UBi,LBj:UBj,N(ng))
270# endif
271# endif
272 real(r8), intent(inout) :: tl_rufrc(LBi:UBi,LBj:UBj)
273 real(r8), intent(inout) :: tl_rvfrc(LBi:UBi,LBj:UBj)
274 real(r8), intent(inout) :: tl_u(LBi:UBi,LBj:UBj,N(ng),2)
275 real(r8), intent(inout) :: tl_v(LBi:UBi,LBj:UBj,N(ng),2)
276#endif
277!
278! Local variable declarations.
279!
280 integer :: i, j, k, k1, k2
281
282 real(r8) :: cff, fac1, fac2, pm_p, pn_p
283 real(r8) :: cff1, cff2, cff3, cff4
284 real(r8) :: cff5, cff6, cff7, cff8
285 real(r8) :: dmUdz, dnUdz, dmVdz, dnVdz
286#ifdef VISC_3DCOEF
287 real(r8) :: Uvis_p, Vvis_p, visc_p
288 real(r8) :: tl_fac1, tl_fac2, tl_Uvis_p, tl_Vvis_p, tl_visc_p
289#endif
290 real(r8) :: tl_cff
291 real(r8) :: tl_cff1, tl_cff2, tl_cff3, tl_cff4
292 real(r8) :: tl_cff5, tl_cff6, tl_cff7, tl_cff8
293 real(r8) :: tl_dmUdz, tl_dnUdz, tl_dmVdz, tl_dnVdz
294
295 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,N(ng)) :: LapU
296 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,N(ng)) :: LapV
297
298 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,N(ng)) :: tl_LapU
299 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,N(ng)) :: tl_LapV
300
301 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: UFe
302 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: UFx
303 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: VFe
304 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: VFx
305
306 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_UFe
307 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_UFx
308 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_VFe
309 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_VFx
310
311 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: UFse
312 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: UFsx
313 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: VFse
314 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: VFsx
315 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dmUde
316 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dmVde
317 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dnUdx
318 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dnVdx
319 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dUdz
320 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dVdz
321 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dZde_p
322 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dZde_r
323 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dZdx_p
324 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dZdx_r
325
326 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: tl_UFse
327 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: tl_UFsx
328 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: tl_VFse
329 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: tl_VFsx
330 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: tl_dmUde
331 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: tl_dmVde
332 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: tl_dnUdx
333 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: tl_dnVdx
334 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: tl_dUdz
335 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: tl_dVdz
336 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: tl_dZde_p
337 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: tl_dZde_r
338 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: tl_dZdx_p
339 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: tl_dZdx_r
340
341#include "set_bounds.h"
342!
343!-----------------------------------------------------------------------
344! Compute horizontal biharmonic viscosity along geopotential
345! surfaces. The biharmonic operator is computed by applying
346! the harmonic operator twice.
347!-----------------------------------------------------------------------
348!
349! Compute horizontal and vertical gradients. Notice the recursive
350! blocking sequence. It is assumed here that the mixing coefficients
351! are the squared root of the biharmonic viscosity coefficient. For
352! momentum balance purposes, the thickness "Hz" appears only when
353! computing the second harmonic operator. The vertical placement of
354! the gradients is:
355!
356! dZdx_r, dZde_r, dnUdx, dmVde(:,:,k1) k rho-points
357! dZdx_r, dZde_r, dnUdx, dmVde(:,:,k2) k+1 rho-points
358! dZdx_p, dZde_p, dnVdx, dmUde(:,:,k1) k psi-points
359! dZdx_p, dZde_p, dnVdx, dmUde(:,:,k2) k+1 psi-points
360! UFse, UFsx, dUdz(:,:,k1) k-1/2 WU-points
361! UFse, UFsx, dUdz(:,:,k2) k+1/2 WU-points
362! VFse, VFsx, dVdz(:,:,k1) k-1/2 WV-points
363! VFse, VFsx, dVdz(:,:,k2) k+1/2 WV-points
364!
365 k2=1
366 k_loop1 : DO k=0,n(ng)
367 k1=k2
368 k2=3-k1
369 IF (k.lt.n(ng)) THEN
370!
371! Compute slopes (nondimensional) at RHO- and PSI-points.
372!
373 DO j=jstrm2,jendp2
374 DO i=istrum2,iendp2
375 cff=0.5_r8*(pm(i-1,j)+pm(i,j))
376#ifdef MASKING
377 cff=cff*umask(i,j)
378#endif
379 ufx(i,j)=cff*(z_r(i ,j,k+1)- &
380 & z_r(i-1,j,k+1))
381 tl_ufx(i,j)=cff*(tl_z_r(i ,j,k+1)- &
382 & tl_z_r(i-1,j,k+1))
383 END DO
384 END DO
385 DO j=jstrvm2,jendp2
386 DO i=istrm2,iendp2
387 cff=0.5_r8*(pn(i,j-1)+pn(i,j))
388#ifdef MASKING
389 cff=cff*vmask(i,j)
390#endif
391 vfe(i,j)=cff*(z_r(i,j ,k+1)- &
392 & z_r(i,j-1,k+1))
393 tl_vfe(i,j)=cff*(tl_z_r(i,j ,k+1)- &
394 & tl_z_r(i,j-1,k+1))
395 END DO
396 END DO
397!
398 DO j=jstrm1,jendp2
399 DO i=istrm1,iendp2
400 dzdx_p(i,j,k2)=0.5_r8*(ufx(i,j-1)+ &
401 & ufx(i,j ))
402 tl_dzdx_p(i,j,k2)=0.5_r8*(tl_ufx(i,j-1)+ &
403 & tl_ufx(i,j ))
404 dzde_p(i,j,k2)=0.5_r8*(vfe(i-1,j)+ &
405 & vfe(i ,j))
406 tl_dzde_p(i,j,k2)=0.5_r8*(tl_vfe(i-1,j)+ &
407 & tl_vfe(i ,j))
408 END DO
409 END DO
410 DO j=jstrvm2,jendp1
411 DO i=istrum2,iendp1
412 dzdx_r(i,j,k2)=0.5_r8*(ufx(i ,j)+ &
413 & ufx(i+1,j))
414 tl_dzdx_r(i,j,k2)=0.5_r8*(tl_ufx(i ,j)+ &
415 & tl_ufx(i+1,j))
416 dzde_r(i,j,k2)=0.5_r8*(vfe(i,j )+ &
417 & vfe(i,j+1))
418 tl_dzde_r(i,j,k2)=0.5_r8*(tl_vfe(i,j )+ &
419 & tl_vfe(i,j+1))
420 END DO
421 END DO
422!
423! Compute momentum horizontal (1/m/s) and vertical (1/s) gradients.
424!
425 DO j=jstrvm2,jendp1
426 DO i=istrum2,iendp1
427 cff=0.5_r8*pm(i,j)
428#ifdef MASKING
429 cff=cff*rmask(i,j)
430#endif
431 dnudx(i,j,k2)=cff*((pn(i ,j)+pn(i+1,j))* &
432 & u(i+1,j,k+1,nrhs)- &
433 & (pn(i-1,j)+pn(i ,j))* &
434 & u(i ,j,k+1,nrhs))
435 tl_dnudx(i,j,k2)=cff*((pn(i ,j)+pn(i+1,j))* &
436 & tl_u(i+1,j,k+1,nrhs)- &
437 & (pn(i-1,j)+pn(i ,j))* &
438 & tl_u(i ,j,k+1,nrhs))
439 END DO
440 END DO
441
442 DO j=jstrm1,jendp2
443 DO i=istrm1,iendp2
444 cff=0.125_r8*(pn(i-1,j )+pn(i,j )+ &
445 & pn(i-1,j-1)+pn(i,j-1))
446#ifdef MASKING
447 cff=cff*pmask(i,j)
448#endif
449 dmude(i,j,k2)=cff*((pm(i-1,j )+pm(i,j ))* &
450 & u(i,j ,k+1,nrhs)- &
451 & (pm(i-1,j-1)+pm(i,j-1))* &
452 & u(i,j-1,k+1,nrhs))
453 tl_dmude(i,j,k2)=cff*((pm(i-1,j )+pm(i,j ))* &
454 & tl_u(i,j ,k+1,nrhs)- &
455 & (pm(i-1,j-1)+pm(i,j-1))* &
456 & tl_u(i,j-1,k+1,nrhs))
457 END DO
458 END DO
459
460 DO j=jstrm1,jendp2
461 DO i=istrm1,iendp2
462 cff=0.125_r8*(pm(i-1,j )+pm(i,j )+ &
463 & pm(i-1,j-1)+pm(i,j-1))
464#ifdef MASKING
465 cff=cff*pmask(i,j)
466#endif
467 dnvdx(i,j,k2)=cff*((pn(i ,j-1)+pn(i ,j))* &
468 & v(i ,j,k+1,nrhs)- &
469 & (pn(i-1,j-1)+pn(i-1,j))* &
470 & v(i-1,j,k+1,nrhs))
471 tl_dnvdx(i,j,k2)=cff*((pn(i ,j-1)+pn(i ,j))* &
472 & tl_v(i ,j,k+1,nrhs)- &
473 & (pn(i-1,j-1)+pn(i-1,j))* &
474 & tl_v(i-1,j,k+1,nrhs))
475 END DO
476 END DO
477
478 DO j=jstrvm2,jendp1
479 DO i=istrum2,iendp1
480 cff=0.5_r8*pn(i,j)
481#ifdef MASKING
482 cff=cff*rmask(i,j)
483#endif
484 dmvde(i,j,k2)=cff*((pm(i,j )+pm(i,j+1))* &
485 & v(i,j+1,k+1,nrhs)- &
486 & (pm(i,j-1)+pm(i,j ))* &
487 & v(i,j ,k+1,nrhs))
488 tl_dmvde(i,j,k2)=cff*((pm(i,j )+pm(i,j+1))* &
489 & tl_v(i,j+1,k+1,nrhs)- &
490 & (pm(i,j-1)+pm(i,j ))* &
491 & tl_v(i,j ,k+1,nrhs))
492 END DO
493 END DO
494 END IF
495
496 IF ((k.eq.0).or.(k.eq.n(ng))) THEN
497 DO j=jstrm2,jendp2
498 DO i=istrum2,iendp2
499 dudz(i,j,k2)=0.0_r8
500 tl_dudz(i,j,k2)=0.0_r8
501 END DO
502 END DO
503 DO j=jstrvm2,jendp2
504 DO i=istrm2,iendp2
505 dvdz(i,j,k2)=0.0_r8
506 tl_dvdz(i,j,k2)=0.0_r8
507 END DO
508 END DO
509
510 DO j=jstrm1,jendp1
511 DO i=istrum1,iendp1
512 ufsx(i,j,k2)=0.0_r8
513 tl_ufsx(i,j,k2)=0.0_r8
514 ufse(i,j,k2)=0.0_r8
515 tl_ufse(i,j,k2)=0.0_r8
516 END DO
517 END DO
518 DO j=jstrvm1,jendp1
519 DO i=istrm1,iendp1
520 vfsx(i,j,k2)=0.0_r8
521 tl_vfsx(i,j,k2)=0.0_r8
522 vfse(i,j,k2)=0.0_r8
523 tl_vfse(i,j,k2)=0.0_r8
524 END DO
525 END DO
526 ELSE
527 DO j=jstrm2,jendp2
528 DO i=istrum2,iendp2
529 cff=1.0_r8/(0.5_r8*(z_r(i-1,j,k+1)- &
530 & z_r(i-1,j,k )+ &
531 & z_r(i ,j,k+1)- &
532 & z_r(i ,j,k )))
533 tl_cff=-cff*cff*(0.5_r8*(tl_z_r(i-1,j,k+1)- &
534 & tl_z_r(i-1,j,k )+ &
535 & tl_z_r(i ,j,k+1)- &
536 & tl_z_r(i ,j,k )))
537 dudz(i,j,k2)=cff*(u(i,j,k+1,nrhs)- &
538 & u(i,j,k ,nrhs))
539 tl_dudz(i,j,k2)=tl_cff*(u(i,j,k+1,nrhs)- &
540 & u(i,j,k ,nrhs))+ &
541 & cff*(tl_u(i,j,k+1,nrhs)- &
542 & tl_u(i,j,k ,nrhs))
543 END DO
544 END DO
545
546 DO j=jstrvm2,jendp2
547 DO i=istrm2,iendp2
548 cff=1.0_r8/(0.5_r8*(z_r(i,j-1,k+1)- &
549 & z_r(i,j-1,k )+ &
550 & z_r(i,j ,k+1)- &
551 & z_r(i,j ,k )))
552 tl_cff=-cff*cff*(0.5_r8*(tl_z_r(i,j-1,k+1)- &
553 & tl_z_r(i,j-1,k )+ &
554 & tl_z_r(i,j ,k+1)- &
555 & tl_z_r(i,j ,k )))
556 dvdz(i,j,k2)=cff*(v(i,j,k+1,nrhs)- &
557 & v(i,j,k ,nrhs))
558 tl_dvdz(i,j,k2)=tl_cff*(v(i,j,k+1,nrhs)- &
559 & v(i,j,k ,nrhs))+ &
560 & cff*(tl_v(i,j,k+1,nrhs)- &
561 & tl_v(i,j,k ,nrhs))
562 END DO
563 END DO
564 END IF
565!
566! Compute components of the rotated viscous flux (m^4 s-^3/2) along
567! geopotential surfaces in the XI- and ETA-directions.
568!
569 IF (k.gt.0) THEN
570 DO j=jstrvm2,jendp1
571 DO i=istrum2,iendp1
572 cff1=min(dzdx_r(i,j,k1),0.0_r8)
573 cff2=max(dzdx_r(i,j,k1),0.0_r8)
574 cff3=min(dzde_r(i,j,k1),0.0_r8)
575 cff4=max(dzde_r(i,j,k1),0.0_r8)
576 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_r(i,j,k1)))* &
577 & tl_dzdx_r(i,j,k1)
578 tl_cff2=(0.5_r8+sign(0.5_r8, dzdx_r(i,j,k1)))* &
579 & tl_dzdx_r(i,j,k1)
580 tl_cff3=(0.5_r8+sign(0.5_r8,-dzde_r(i,j,k1)))* &
581 & tl_dzde_r(i,j,k1)
582 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_r(i,j,k1)))* &
583 & tl_dzde_r(i,j,k1)
584 cff=on_r(i,j)*(dnudx(i,j,k1)- &
585 & 0.5_r8*pn(i,j)* &
586 & (cff1*(dudz(i ,j,k1)+ &
587 & dudz(i+1,j,k2))+ &
588 & cff2*(dudz(i ,j,k2)+ &
589 & dudz(i+1,j,k1))))- &
590 & om_r(i,j)*(dmvde(i,j,k1)- &
591 & 0.5_r8*pm(i,j)* &
592 & (cff3*(dvdz(i,j ,k1)+ &
593 & dvdz(i,j+1,k2))+ &
594 & cff4*(dvdz(i,j ,k2)+ &
595 & dvdz(i,j+1,k1))))
596 tl_cff=on_r(i,j)*(tl_dnudx(i,j,k1)- &
597 & 0.5_r8*pn(i,j)* &
598 & (tl_cff1*(dudz(i ,j,k1)+ &
599 & dudz(i+1,j,k2))+ &
600 & cff1*(tl_dudz(i ,j,k1)+ &
601 & tl_dudz(i+1,j,k2))+ &
602 & tl_cff2*(dudz(i ,j,k2)+ &
603 & dudz(i+1,j,k1))+ &
604 & cff2*(tl_dudz(i ,j,k2)+ &
605 & tl_dudz(i+1,j,k1))))- &
606 & om_r(i,j)*(tl_dmvde(i,j,k1)- &
607 & 0.5_r8*pm(i,j)* &
608 & (tl_cff3*(dvdz(i,j ,k1)+ &
609 & dvdz(i,j+1,k2))+ &
610 & cff3*(tl_dvdz(i,j ,k1)+ &
611 & tl_dvdz(i,j+1,k2))+ &
612 & tl_cff4*(dvdz(i,j ,k2)+ &
613 & dvdz(i,j+1,k1))+ &
614 & cff4*(tl_dvdz(i,j ,k2)+ &
615 & tl_dvdz(i,j+1,k1))))
616#ifdef MASKING
617 cff=cff*rmask(i,j)
618 tl_cff=tl_cff*rmask(i,j)
619#endif
620#ifdef VISC_3DCOEF
621# ifdef UV_U3ADV_SPLIT
622 ufx(i,j)=on_r(i,j)*on_r(i,j)*uvis3d_r(i,j,k)*cff
623 tl_ufx(i,j)=on_r(i,j)*on_r(i,j)* &
624 & (tl_uvis3d_r(i,j,k)*cff+ &
625 & uvis3d_r(i,j,k)*tl_cff)
626 vfe(i,j)=om_r(i,j)*om_r(i,j)*vvis3d_r(i,j,k)*cff
627 tl_vfe(i,j)=om_r(i,j)*om_r(i,j)* &
628 & (tl_vvis3d_r(i,j,k)*cff+ &
629 & vvis3d_r(i,j,k)*tl_cff)
630# else
631 ufx(i,j)=on_r(i,j)*on_r(i,j)*visc3d_r(i,j,k)*cff
632 tl_ufx(i,j)=on_r(i,j)*on_r(i,j)* &
633 & (tl_visc3d_r(i,j,k)*cff+ &
634 & visc3d_r(i,j,k)*tl_cff)
635 vfe(i,j)=om_r(i,j)*om_r(i,j)*visc3d_r(i,j,k)*cff
636 tl_vfe(i,j)=om_r(i,j)*om_r(i,j)* &
637 & (tl_visc3d_r(i,j,k)*cff+ &
638 & visc3d_r(i,j,k)*tl_cff)
639# endif
640#else
641 ufx(i,j)=on_r(i,j)*on_r(i,j)*visc4_r(i,j)*cff
642 tl_ufx(i,j)=on_r(i,j)*on_r(i,j)*visc4_r(i,j)*tl_cff
643 vfe(i,j)=om_r(i,j)*om_r(i,j)*visc4_r(i,j)*cff
644 tl_vfe(i,j)=om_r(i,j)*om_r(i,j)*visc4_r(i,j)*tl_cff
645#endif
646 END DO
647 END DO
648
649 DO j=jstrm1,jendp2
650 DO i=istrm1,iendp2
651 pm_p=0.25_r8*(pm(i-1,j-1)+pm(i-1,j)+ &
652 & pm(i ,j-1)+pm(i ,j))
653 pn_p=0.25_r8*(pn(i-1,j-1)+pn(i-1,j)+ &
654 & pn(i ,j-1)+pn(i ,j))
655 cff1=min(dzdx_p(i,j,k1),0.0_r8)
656 cff2=max(dzdx_p(i,j,k1),0.0_r8)
657 cff3=min(dzde_p(i,j,k1),0.0_r8)
658 cff4=max(dzde_p(i,j,k1),0.0_r8)
659 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_p(i,j,k1)))* &
660 & tl_dzdx_p(i,j,k1)
661 tl_cff2=(0.5_r8+sign(0.5_r8, dzdx_p(i,j,k1)))* &
662 & tl_dzdx_p(i,j,k1)
663 tl_cff3=(0.5_r8+sign(0.5_r8,-dzde_p(i,j,k1)))* &
664 & tl_dzde_p(i,j,k1)
665 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_p(i,j,k1)))* &
666 & tl_dzde_p(i,j,k1)
667 cff=on_p(i,j)*(dnvdx(i,j,k1)- &
668 & 0.5_r8*pn_p* &
669 & (cff1*(dvdz(i-1,j,k1)+ &
670 & dvdz(i ,j,k2))+ &
671 & cff2*(dvdz(i-1,j,k2)+ &
672 & dvdz(i ,j,k1))))+ &
673 & om_p(i,j)*(dmude(i,j,k1)- &
674 & 0.5_r8*pm_p* &
675 & (cff3*(dudz(i,j-1,k1)+ &
676 & dudz(i,j ,k2))+ &
677 & cff4*(dudz(i,j-1,k2)+ &
678 & dudz(i,j ,k1))))
679 tl_cff=on_p(i,j)*(tl_dnvdx(i,j,k1)- &
680 & 0.5_r8*pn_p* &
681 & (tl_cff1*(dvdz(i-1,j,k1)+ &
682 & dvdz(i ,j,k2))+ &
683 & cff1*(tl_dvdz(i-1,j,k1)+ &
684 & tl_dvdz(i ,j,k2))+ &
685 & tl_cff2*(dvdz(i-1,j,k2)+ &
686 & dvdz(i ,j,k1))+ &
687 & cff2*(tl_dvdz(i-1,j,k2)+ &
688 & tl_dvdz(i ,j,k1))))+ &
689 & om_p(i,j)*(tl_dmude(i,j,k1)- &
690 & 0.5_r8*pm_p* &
691 & (tl_cff3*(dudz(i,j-1,k1)+ &
692 & dudz(i,j ,k2))+ &
693 & cff3*(tl_dudz(i,j-1,k1)+ &
694 & tl_dudz(i,j ,k2))+ &
695 & tl_cff4*(dudz(i,j-1,k2)+ &
696 & dudz(i,j ,k1))+ &
697 & cff4*(tl_dudz(i,j-1,k2)+ &
698 & tl_dudz(i,j ,k1))))
699#ifdef MASKING
700 cff=cff*pmask(i,j)
701 tl_cff=tl_cff*pmask(i,j)
702#endif
703#ifdef VISC_3DCOEF
704# ifdef UV_U3ADV_SPLIT
705 uvis_p=0.25_r8* &
706 & (uvis3d_r(i-1,j-1,k)+uvis3d_r(i-1,j,k)+ &
707 & uvis3d_r(i ,j-1,k)+uvis3d_r(i ,j,k))
708 tl_uvis_p=0.25_r8* &
709 & (tl_uvis3d_r(i-1,j-1,k)+tl_uvis3d_r(i-1,j,k)+ &
710 & tl_uvis3d_r(i ,j-1,k)+tl_uvis3d_r(i ,j,k))
711 vvis_p=0.25_r8* &
712 & (vvis3d_r(i-1,j-1,k)+vvis3d_r(i-1,j,k)+ &
713 & vvis3d_r(i ,j-1,k)+vvis3d_r(i ,j,k))
714 tl_vvis_p=0.25_r8* &
715 & (tl_vvis3d_r(i-1,j-1,k)+tl_vvis3d_r(i-1,j,k)+ &
716 & tl_vvis3d_r(i ,j-1,k)+tl_vvis3d_r(i ,j,k))
717 ufe(i,j)=om_p(i,j)*om_p(i,j)*uvis_p*cff
718 tl_ufe(i,j)=om_p(i,j)*om_p(i,j)* &
719 & (tl_uvis_p*cff+uvis_p*tl_cff)
720 vfx(i,j)=on_p(i,j)*on_p(i,j)*vvis_p*cff
721 tl_vfx(i,j)=on_p(i,j)*on_p(i,j)* &
722 & (tl_vvis_p*cff+vvis_p*tl_cff)
723# else
724 visc_p=0.25_r8* &
725 & (visc3d_r(i-1,j-1,k)+visc3d_r(i-1,j,k)+ &
726 & visc3d_r(i ,j-1,k)+visc3d_r(i ,j,k))
727 tl_visc_p=0.25_r8* &
728 & (tl_visc3d_r(i-1,j-1,k)+tl_visc3d_r(i-1,j,k)+ &
729 & tl_visc3d_r(i ,j-1,k)+tl_visc3d_r(i ,j,k))
730 ufe(i,j)=om_p(i,j)*om_p(i,j)*visc_p*cff
731 tl_ufe(i,j)=om_p(i,j)*om_p(i,j)* &
732 & (tl_visc_p*cff+visc_p*tl_cff)
733 vfx(i,j)=on_p(i,j)*on_p(i,j)*visc_p*cff
734 tl_vfx(i,j)=on_p(i,j)*on_p(i,j)* &
735 & (tl_visc_p*cff+visc_p*tl_cff)
736# endif
737#else
738 ufe(i,j)=om_p(i,j)*om_p(i,j)*visc4_p(i,j)*cff
739 tl_ufe(i,j)=om_p(i,j)*om_p(i,j)*visc4_p(i,j)*tl_cff
740 vfx(i,j)=on_p(i,j)*on_p(i,j)*visc4_p(i,j)*cff
741 tl_vfx(i,j)=on_p(i,j)*on_p(i,j)*visc4_p(i,j)*tl_cff
742#endif
743 END DO
744 END DO
745!
746! Compute vertical flux (m^2 s^-3/2) due to sloping terrain-following
747! surfaces.
748!
749 IF (k.lt.n(ng)) THEN
750 DO j=jstrm1,jendp1
751 DO i=istrum1,iendp1
752#ifdef VISC_3DCOEF
753# ifdef UV_U3ADV_SPLIT
754 cff=0.125_r8* &
755 & (uvis3d_r(i-1,j,k )+uvis3d_r(i,j,k )+ &
756 & uvis3d_r(i-1,j,k+1)+uvis3d_r(i,j,k+1))
757 tl_cff=0.125_r8* &
758 & (tl_uvis3d_r(i-1,j,k )+tl_uvis3d_r(i,j,k )+ &
759 & tl_uvis3d_r(i-1,j,k+1)+tl_uvis3d_r(i,j,k+1))
760# else
761 cff=0.125_r8* &
762 & (visc3d_r(i-1,j,k )+visc3d_r(i,j,k )+ &
763 & visc3d_r(i-1,j,k+1)+visc3d_r(i,j,k+1))
764 tl_cff=0.125_r8* &
765 & (tl_visc3d_r(i-1,j,k )+tl_visc3d_r(i,j,k )+ &
766 & tl_visc3d_r(i-1,j,k+1)+tl_visc3d_r(i,j,k+1))
767# endif
768 fac1=cff*on_u(i,j)
769 tl_fac1=tl_cff*on_u(i,j)
770 fac2=cff*om_u(i,j)
771 tl_fac2=tl_cff*om_u(i,j)
772#else
773 cff=0.25_r8*(visc4_r(i-1,j)+visc4_r(i,j))
774 fac1=cff*on_u(i,j)
775 fac2=cff*om_u(i,j)
776#endif
777 cff=0.5_r8*(pn(i-1,j)+pn(i,j))
778 dnudz=cff*dudz(i,j,k2)
779 tl_dnudz=cff*tl_dudz(i,j,k2)
780 dnvdz=cff*0.25_r8*(dvdz(i-1,j+1,k2)+ &
781 & dvdz(i ,j+1,k2)+ &
782 & dvdz(i-1,j ,k2)+ &
783 & dvdz(i ,j ,k2))
784 tl_dnvdz=cff*0.25_r8*(tl_dvdz(i-1,j+1,k2)+ &
785 & tl_dvdz(i ,j+1,k2)+ &
786 & tl_dvdz(i-1,j ,k2)+ &
787 & tl_dvdz(i ,j ,k2))
788 cff=0.5_r8*(pm(i-1,j)+pm(i,j))
789 dmudz=cff*dudz(i,j,k2)
790 tl_dmudz=cff*tl_dudz(i,j,k2)
791 dmvdz=cff*0.25_r8*(dvdz(i-1,j+1,k2)+ &
792 & dvdz(i ,j+1,k2)+ &
793 & dvdz(i-1,j ,k2)+ &
794 & dvdz(i ,j ,k2))
795 tl_dmvdz=cff*0.25_r8*(tl_dvdz(i-1,j+1,k2)+ &
796 & tl_dvdz(i ,j+1,k2)+ &
797 & tl_dvdz(i-1,j ,k2)+ &
798 & tl_dvdz(i ,j ,k2))
799
800 cff1=min(dzdx_r(i-1,j,k1),0.0_r8)
801 cff2=min(dzdx_r(i ,j,k2),0.0_r8)
802 cff3=max(dzdx_r(i-1,j,k2),0.0_r8)
803 cff4=max(dzdx_r(i ,j,k1),0.0_r8)
804 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_r(i-1,j,k1)))* &
805 & tl_dzdx_r(i-1,j,k1)
806 tl_cff2=(0.5_r8+sign(0.5_r8,-dzdx_r(i ,j,k2)))* &
807 & tl_dzdx_r(i ,j,k2)
808 tl_cff3=(0.5_r8+sign(0.5_r8, dzdx_r(i-1,j,k2)))* &
809 & tl_dzdx_r(i-1,j,k2)
810 tl_cff4=(0.5_r8+sign(0.5_r8, dzdx_r(i ,j,k1)))* &
811 & tl_dzdx_r(i ,j,k1)
812 ufsx(i,j,k2)=fac1* &
813 & (cff1*(cff1*dnudz-dnudx(i-1,j,k1))+ &
814 & cff2*(cff2*dnudz-dnudx(i ,j,k2))+ &
815 & cff3*(cff3*dnudz-dnudx(i-1,j,k2))+ &
816 & cff4*(cff4*dnudz-dnudx(i ,j,k1)))
817 tl_ufsx(i,j,k2)=fac1* &
818 & (tl_cff1*(cff1*dnudz-dnudx(i-1,j,k1))+ &
819 & tl_cff2*(cff2*dnudz-dnudx(i ,j,k2))+ &
820 & tl_cff3*(cff3*dnudz-dnudx(i-1,j,k2))+ &
821 & tl_cff4*(cff4*dnudz-dnudx(i ,j,k1))+ &
822 & cff1*(tl_cff1*dnudz+cff1*tl_dnudz- &
823 & tl_dnudx(i-1,j,k1))+ &
824 & cff2*(tl_cff2*dnudz+cff2*tl_dnudz- &
825 & tl_dnudx(i ,j,k2))+ &
826 & cff3*(tl_cff3*dnudz+cff3*tl_dnudz- &
827 & tl_dnudx(i-1,j,k2))+ &
828 & cff4*(tl_cff4*dnudz+cff4*tl_dnudz- &
829 & tl_dnudx(i ,j,k1)))
830#ifdef VISC_3DCOEF
831 tl_ufsx(i,j,k2)=tl_ufsx(i,j,k2)+ &
832 & tl_fac1* &
833 & (cff1*(cff1*dnudz-dnudx(i-1,j,k1))+ &
834 & cff2*(cff2*dnudz-dnudx(i ,j,k2))+ &
835 & cff3*(cff3*dnudz-dnudx(i-1,j,k2))+ &
836 & cff4*(cff4*dnudz-dnudx(i ,j,k1)))
837#endif
838
839 cff1=min(dzde_p(i,j ,k1),0.0_r8)
840 cff2=min(dzde_p(i,j+1,k2),0.0_r8)
841 cff3=max(dzde_p(i,j ,k2),0.0_r8)
842 cff4=max(dzde_p(i,j+1,k1),0.0_r8)
843 tl_cff1=(0.5_r8+sign(0.5_r8,-dzde_p(i,j ,k1)))* &
844 & tl_dzde_p(i,j ,k1)
845 tl_cff2=(0.5_r8+sign(0.5_r8,-dzde_p(i,j+1,k2)))* &
846 & tl_dzde_p(i,j+1,k2)
847 tl_cff3=(0.5_r8+sign(0.5_r8, dzde_p(i,j ,k2)))* &
848 & tl_dzde_p(i,j ,k2)
849 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_p(i,j+1,k1)))* &
850 & tl_dzde_p(i,j+1,k1)
851 ufse(i,j,k2)=fac2* &
852 & (cff1*(cff1*dmudz-dmude(i,j ,k1))+ &
853 & cff2*(cff2*dmudz-dmude(i,j+1,k2))+ &
854 & cff3*(cff3*dmudz-dmude(i,j ,k2))+ &
855 & cff4*(cff4*dmudz-dmude(i,j+1,k1)))
856 tl_ufse(i,j,k2)=fac2* &
857 & (tl_cff1*(cff1*dmudz-dmude(i,j ,k1))+ &
858 & tl_cff2*(cff2*dmudz-dmude(i,j+1,k2))+ &
859 & tl_cff3*(cff3*dmudz-dmude(i,j ,k2))+ &
860 & tl_cff4*(cff4*dmudz-dmude(i,j+1,k1))+ &
861 & cff1*(tl_cff1*dmudz+cff1*tl_dmudz- &
862 & tl_dmude(i,j ,k1))+ &
863 & cff2*(tl_cff2*dmudz+cff2*tl_dmudz- &
864 & tl_dmude(i,j+1,k2))+ &
865 & cff3*(tl_cff3*dmudz+cff3*tl_dmudz- &
866 & tl_dmude(i,j ,k2))+ &
867 & cff4*(tl_cff4*dmudz+cff4*tl_dmudz- &
868 & tl_dmude(i,j+1,k1)))
869#ifdef VISC_3DCOEF
870 tl_ufse(i,j,k2)=tl_ufse(i,j,k2)+ &
871 & tl_fac2* &
872 & (cff1*(cff1*dmudz-dmude(i,j ,k1))+ &
873 & cff2*(cff2*dmudz-dmude(i,j+1,k2))+ &
874 & cff3*(cff3*dmudz-dmude(i,j ,k2))+ &
875 & cff4*(cff4*dmudz-dmude(i,j+1,k1)))
876#endif
877
878 cff1=min(dzde_p(i,j ,k1),0.0_r8)
879 cff2=min(dzde_p(i,j+1,k2),0.0_r8)
880 cff3=max(dzde_p(i,j ,k2),0.0_r8)
881 cff4=max(dzde_p(i,j+1,k1),0.0_r8)
882 cff5=min(dzdx_p(i,j ,k1),0.0_r8)
883 cff6=min(dzdx_p(i,j+1,k2),0.0_r8)
884 cff7=max(dzdx_p(i,j ,k2),0.0_r8)
885 cff8=max(dzdx_p(i,j+1,k1),0.0_r8)
886 tl_cff1=(0.5_r8+sign(0.5_r8,-dzde_p(i,j ,k1)))* &
887 & tl_dzde_p(i,j ,k1)
888 tl_cff2=(0.5_r8+sign(0.5_r8,-dzde_p(i,j+1,k2)))* &
889 & tl_dzde_p(i,j+1,k2)
890 tl_cff3=(0.5_r8+sign(0.5_r8, dzde_p(i,j ,k2)))* &
891 & tl_dzde_p(i,j ,k2)
892 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_p(i,j+1,k1)))* &
893 & tl_dzde_p(i,j+1,k1)
894 tl_cff5=(0.5_r8+sign(0.5_r8,-dzdx_p(i,j ,k1)))* &
895 & tl_dzdx_p(i,j ,k1)
896 tl_cff6=(0.5_r8+sign(0.5_r8,-dzdx_p(i,j+1,k2)))* &
897 & tl_dzdx_p(i,j+1,k2)
898 tl_cff7=(0.5_r8+sign(0.5_r8, dzdx_p(i,j ,k2)))* &
899 & tl_dzdx_p(i,j ,k2)
900 tl_cff8=(0.5_r8+sign(0.5_r8, dzdx_p(i,j+1,k1)))* &
901 & tl_dzdx_p(i,j+1,k1)
902 ufsx(i,j,k2)=ufsx(i,j,k2)+ &
903 & fac1* &
904 & (cff1*(cff5*dnvdz-dnvdx(i,j ,k1))+ &
905 & cff2*(cff6*dnvdz-dnvdx(i,j+1,k2))+ &
906 & cff3*(cff7*dnvdz-dnvdx(i,j ,k2))+ &
907 & cff4*(cff8*dnvdz-dnvdx(i,j+1,k1)))
908 tl_ufsx(i,j,k2)=tl_ufsx(i,j,k2)+ &
909 & fac1* &
910 & (tl_cff1*(cff5*dnvdz-dnvdx(i,j ,k1))+ &
911 & tl_cff2*(cff6*dnvdz-dnvdx(i,j+1,k2))+ &
912 & tl_cff3*(cff7*dnvdz-dnvdx(i,j ,k2))+ &
913 & tl_cff4*(cff8*dnvdz-dnvdx(i,j+1,k1))+ &
914 & cff1*(tl_cff5*dnvdz+cff5*tl_dnvdz- &
915 & tl_dnvdx(i,j ,k1))+ &
916 & cff2*(tl_cff6*dnvdz+cff6*tl_dnvdz- &
917 & tl_dnvdx(i,j+1,k2))+ &
918 & cff3*(tl_cff7*dnvdz+cff7*tl_dnvdz- &
919 & tl_dnvdx(i,j ,k2))+ &
920 & cff4*(tl_cff8*dnvdz+cff8*tl_dnvdz- &
921 & tl_dnvdx(i,j+1,k1)))
922#ifdef VISC_3DCOEF
923 tl_ufsx(i,j,k2)=tl_ufsx(i,j,k2)+ &
924 & tl_fac1* &
925 & (cff1*(cff5*dnvdz-dnvdx(i,j ,k1))+ &
926 & cff2*(cff6*dnvdz-dnvdx(i,j+1,k2))+ &
927 & cff3*(cff7*dnvdz-dnvdx(i,j ,k2))+ &
928 & cff4*(cff8*dnvdz-dnvdx(i,j+1,k1)))
929#endif
930
931 cff1=min(dzdx_r(i-1,j,k1),0.0_r8)
932 cff2=min(dzdx_r(i ,j,k2),0.0_r8)
933 cff3=max(dzdx_r(i-1,j,k2),0.0_r8)
934 cff4=max(dzdx_r(i ,j,k1),0.0_r8)
935 cff5=min(dzde_r(i-1,j,k1),0.0_r8)
936 cff6=min(dzde_r(i ,j,k2),0.0_r8)
937 cff7=max(dzde_r(i-1,j,k2),0.0_r8)
938 cff8=max(dzde_r(i ,j,k1),0.0_r8)
939 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_r(i-1,j,k1)))* &
940 & tl_dzdx_r(i-1,j,k1)
941 tl_cff2=(0.5_r8+sign(0.5_r8,-dzdx_r(i ,j,k2)))* &
942 & tl_dzdx_r(i ,j,k2)
943 tl_cff3=(0.5_r8+sign(0.5_r8, dzdx_r(i-1,j,k2)))* &
944 & tl_dzdx_r(i-1,j,k2)
945 tl_cff4=(0.5_r8+sign(0.5_r8, dzdx_r(i ,j,k1)))* &
946 & tl_dzdx_r(i ,j,k1)
947 tl_cff5=(0.5_r8+sign(0.5_r8,-dzde_r(i-1,j,k1)))* &
948 & tl_dzde_r(i-1,j,k1)
949 tl_cff6=(0.5_r8+sign(0.5_r8,-dzde_r(i ,j,k2)))* &
950 & tl_dzde_r(i ,j,k2)
951 tl_cff7=(0.5_r8+sign(0.5_r8, dzde_r(i-1,j,k2)))* &
952 & tl_dzde_r(i-1,j,k2)
953 tl_cff8=(0.5_r8+sign(0.5_r8, dzde_r(i ,j,k1)))* &
954 & tl_dzde_r(i ,j,k1)
955 ufse(i,j,k2)=ufse(i,j,k2)- &
956 & fac2* &
957 & (cff1*(cff5*dmvdz-dmvde(i-1,j,k1))+ &
958 & cff2*(cff6*dmvdz-dmvde(i ,j,k2))+ &
959 & cff3*(cff7*dmvdz-dmvde(i-1,j,k2))+ &
960 & cff4*(cff8*dmvdz-dmvde(i ,j,k1)))
961 tl_ufse(i,j,k2)=tl_ufse(i,j,k2)- &
962 & fac2* &
963 & (tl_cff1*(cff5*dmvdz-dmvde(i-1,j,k1))+ &
964 & tl_cff2*(cff6*dmvdz-dmvde(i ,j,k2))+ &
965 & tl_cff3*(cff7*dmvdz-dmvde(i-1,j,k2))+ &
966 & tl_cff4*(cff8*dmvdz-dmvde(i ,j,k1))+ &
967 & cff1*(tl_cff5*dmvdz+cff5*tl_dmvdz- &
968 & tl_dmvde(i-1,j,k1))+ &
969 & cff2*(tl_cff6*dmvdz+cff6*tl_dmvdz- &
970 & tl_dmvde(i ,j,k2))+ &
971 & cff3*(tl_cff7*dmvdz+cff7*tl_dmvdz- &
972 & tl_dmvde(i-1,j,k2))+ &
973 & cff4*(tl_cff8*dmvdz+cff8*tl_dmvdz- &
974 & tl_dmvde(i ,j,k1)))
975#ifdef VISC_3DCOEF
976 tl_ufse(i,j,k2)=tl_ufse(i,j,k2)- &
977 & tl_fac2* &
978 & (cff1*(cff5*dmvdz-dmvde(i-1,j,k1))+ &
979 & cff2*(cff6*dmvdz-dmvde(i ,j,k2))+ &
980 & cff3*(cff7*dmvdz-dmvde(i-1,j,k2))+ &
981 & cff4*(cff8*dmvdz-dmvde(i ,j,k1)))
982#endif
983 END DO
984 END DO
985!
986 DO j=jstrvm1,jendp1
987 DO i=istrm1,iendp1
988#ifdef VISC_3DCOEF
989# ifdef UV_U3ADV_SPLIT
990 cff=0.125_r8* &
991 & (vvis3d_r(i,j-1,k )+vvis3d_r(i,j,k )+ &
992 & vvis3d_r(i,j-1,k+1)+vvis3d_r(i,j,k+1))
993 tl_cff=0.125_r8* &
994 & (tl_vvis3d_r(i,j-1,k )+tl_vvis3d_r(i,j,k )+ &
995 & tl_vvis3d_r(i,j-1,k+1)+tl_vvis3d_r(i,j,k+1))
996# else
997 cff=0.125_r8* &
998 & (visc3d_r(i,j-1,k )+visc3d_r(i,j,k )+ &
999 & visc3d_r(i,j-1,k+1)+visc3d_r(i,j,k+1))
1000 tl_cff=0.125_r8* &
1001 & (tl_visc3d_r(i,j-1,k )+tl_visc3d_r(i,j,k )+ &
1002 & tl_visc3d_r(i,j-1,k+1)+tl_visc3d_r(i,j,k+1))
1003# endif
1004 fac1=cff*on_v(i,j)
1005 tl_fac1=tl_cff*on_v(i,j)
1006 fac2=cff*om_v(i,j)
1007 tl_fac2=tl_cff*om_v(i,j)
1008#else
1009 cff=0.25_r8*(visc4_r(i,j-1)+visc4_r(i,j))
1010 fac1=cff*on_v(i,j)
1011 fac2=cff*om_v(i,j)
1012#endif
1013 cff=0.5_r8*(pn(i,j-1)+pn(i,j))
1014 dnudz=cff*0.25_r8*(dudz(i ,j ,k2)+ &
1015 & dudz(i+1,j ,k2)+ &
1016 & dudz(i ,j-1,k2)+ &
1017 & dudz(i+1,j-1,k2))
1018 tl_dnudz=cff*0.25_r8*(tl_dudz(i ,j ,k2)+ &
1019 & tl_dudz(i+1,j ,k2)+ &
1020 & tl_dudz(i ,j-1,k2)+ &
1021 & tl_dudz(i+1,j-1,k2))
1022 dnvdz=cff*dvdz(i,j,k2)
1023 tl_dnvdz=cff*tl_dvdz(i,j,k2)
1024 cff=0.5_r8*(pm(i,j-1)+pm(i,j))
1025 dmudz=cff*0.25_r8*(dudz(i ,j ,k2)+ &
1026 & dudz(i+1,j ,k2)+ &
1027 & dudz(i ,j-1,k2)+ &
1028 & dudz(i+1,j-1,k2))
1029 tl_dmudz=cff*0.25_r8*(tl_dudz(i ,j ,k2)+ &
1030 & tl_dudz(i+1,j ,k2)+ &
1031 & tl_dudz(i ,j-1,k2)+ &
1032 & tl_dudz(i+1,j-1,k2))
1033 dmvdz=cff*dvdz(i,j,k2)
1034 tl_dmvdz=cff*tl_dvdz(i,j,k2)
1035
1036 cff1=min(dzdx_p(i ,j,k1),0.0_r8)
1037 cff2=min(dzdx_p(i+1,j,k2),0.0_r8)
1038 cff3=max(dzdx_p(i ,j,k2),0.0_r8)
1039 cff4=max(dzdx_p(i+1,j,k1),0.0_r8)
1040 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_p(i ,j,k1)))* &
1041 & tl_dzdx_p(i ,j,k1)
1042 tl_cff2=(0.5_r8+sign(0.5_r8,-dzdx_p(i+1,j,k2)))* &
1043 & tl_dzdx_p(i+1,j,k2)
1044 tl_cff3=(0.5_r8+sign(0.5_r8, dzdx_p(i ,j,k2)))* &
1045 & tl_dzdx_p(i ,j,k2)
1046 tl_cff4=(0.5_r8+sign(0.5_r8, dzdx_p(i+1,j,k1)))* &
1047 & tl_dzdx_p(i+1,j,k1)
1048 vfsx(i,j,k2)=fac1* &
1049 & (cff1*(cff1*dnvdz-dnvdx(i ,j,k1))+ &
1050 & cff2*(cff2*dnvdz-dnvdx(i+1,j,k2))+ &
1051 & cff3*(cff3*dnvdz-dnvdx(i ,j,k2))+ &
1052 & cff4*(cff4*dnvdz-dnvdx(i+1,j,k1)))
1053 tl_vfsx(i,j,k2)=fac1* &
1054 & (tl_cff1*(cff1*dnvdz-dnvdx(i ,j,k1))+ &
1055 & tl_cff2*(cff2*dnvdz-dnvdx(i+1,j,k2))+ &
1056 & tl_cff3*(cff3*dnvdz-dnvdx(i ,j,k2))+ &
1057 & tl_cff4*(cff4*dnvdz-dnvdx(i+1,j,k1))+ &
1058 & cff1*(tl_cff1*dnvdz+cff1*tl_dnvdz- &
1059 & tl_dnvdx(i ,j,k1))+ &
1060 & cff2*(tl_cff2*dnvdz+cff2*tl_dnvdz- &
1061 & tl_dnvdx(i+1,j,k2))+ &
1062 & cff3*(tl_cff3*dnvdz+cff3*tl_dnvdz- &
1063 & tl_dnvdx(i ,j,k2))+ &
1064 & cff4*(tl_cff4*dnvdz+cff4*tl_dnvdz- &
1065 & tl_dnvdx(i+1,j,k1)))
1066#ifdef VISC_3DCOEF
1067 tl_vfsx(i,j,k2)=tl_vfsx(i,j,k2)+ &
1068 & tl_fac1* &
1069 & (cff1*(cff1*dnvdz-dnvdx(i ,j,k1))+ &
1070 & cff2*(cff2*dnvdz-dnvdx(i+1,j,k2))+ &
1071 & cff3*(cff3*dnvdz-dnvdx(i ,j,k2))+ &
1072 & cff4*(cff4*dnvdz-dnvdx(i+1,j,k1)))
1073#endif
1074
1075 cff1=min(dzde_r(i,j-1,k1),0.0_r8)
1076 cff2=min(dzde_r(i,j ,k2),0.0_r8)
1077 cff3=max(dzde_r(i,j-1,k2),0.0_r8)
1078 cff4=max(dzde_r(i,j ,k1),0.0_r8)
1079 tl_cff1=(0.5_r8+sign(0.5_r8,-dzde_r(i,j-1,k1)))* &
1080 & tl_dzde_r(i,j-1,k1)
1081 tl_cff2=(0.5_r8+sign(0.5_r8,-dzde_r(i,j ,k2)))* &
1082 & tl_dzde_r(i,j ,k2)
1083 tl_cff3=(0.5_r8+sign(0.5_r8, dzde_r(i,j-1,k2)))* &
1084 & tl_dzde_r(i,j-1,k2)
1085 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_r(i,j ,k1)))* &
1086 & tl_dzde_r(i,j ,k1)
1087 vfse(i,j,k2)=fac2* &
1088 & (cff1*(cff1*dmvdz-dmvde(i,j-1,k1))+ &
1089 & cff2*(cff2*dmvdz-dmvde(i,j ,k2))+ &
1090 & cff3*(cff3*dmvdz-dmvde(i,j-1,k2))+ &
1091 & cff4*(cff4*dmvdz-dmvde(i,j ,k1)))
1092 tl_vfse(i,j,k2)=fac2* &
1093 & (tl_cff1*(cff1*dmvdz-dmvde(i,j-1,k1))+ &
1094 & tl_cff2*(cff2*dmvdz-dmvde(i,j ,k2))+ &
1095 & tl_cff3*(cff3*dmvdz-dmvde(i,j-1,k2))+ &
1096 & tl_cff4*(cff4*dmvdz-dmvde(i,j ,k1))+ &
1097 & cff1*(tl_cff1*dmvdz+cff1*tl_dmvdz- &
1098 & tl_dmvde(i,j-1,k1))+ &
1099 & cff2*(tl_cff2*dmvdz+cff2*tl_dmvdz- &
1100 & tl_dmvde(i,j ,k2))+ &
1101 & cff3*(tl_cff3*dmvdz+cff3*tl_dmvdz- &
1102 & tl_dmvde(i,j-1,k2))+ &
1103 & cff4*(tl_cff4*dmvdz+cff4*tl_dmvdz- &
1104 & tl_dmvde(i,j ,k1)))
1105#ifdef VISC_3DCOEF
1106 tl_vfse(i,j,k2)=tl_vfse(i,j,k2)+ &
1107 & tl_fac2* &
1108 & (cff1*(cff1*dmvdz-dmvde(i,j-1,k1))+ &
1109 & cff2*(cff2*dmvdz-dmvde(i,j ,k2))+ &
1110 & cff3*(cff3*dmvdz-dmvde(i,j-1,k2))+ &
1111 & cff4*(cff4*dmvdz-dmvde(i,j ,k1)))
1112#endif
1113
1114 cff1=min(dzde_r(i,j-1,k1),0.0_r8)
1115 cff2=min(dzde_r(i,j ,k2),0.0_r8)
1116 cff3=max(dzde_r(i,j-1,k2),0.0_r8)
1117 cff4=max(dzde_r(i,j ,k1),0.0_r8)
1118 cff5=min(dzdx_r(i,j-1,k1),0.0_r8)
1119 cff6=min(dzdx_r(i,j ,k2),0.0_r8)
1120 cff7=max(dzdx_r(i,j-1,k2),0.0_r8)
1121 cff8=max(dzdx_r(i,j ,k1),0.0_r8)
1122 tl_cff1=(0.5_r8+sign(0.5_r8,-dzde_r(i,j-1,k1)))* &
1123 & tl_dzde_r(i,j-1,k1)
1124 tl_cff2=(0.5_r8+sign(0.5_r8,-dzde_r(i,j ,k2)))* &
1125 & tl_dzde_r(i,j ,k2)
1126 tl_cff3=(0.5_r8+sign(0.5_r8, dzde_r(i,j-1,k2)))* &
1127 & tl_dzde_r(i,j-1,k2)
1128 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_r(i,j ,k1)))* &
1129 & tl_dzde_r(i,j ,k1)
1130 tl_cff5=(0.5_r8+sign(0.5_r8,-dzdx_r(i,j-1,k1)))* &
1131 & tl_dzdx_r(i,j-1,k1)
1132 tl_cff6=(0.5_r8+sign(0.5_r8,-dzdx_r(i,j ,k2)))* &
1133 & tl_dzdx_r(i,j ,k2)
1134 tl_cff7=(0.5_r8+sign(0.5_r8, dzdx_r(i,j-1,k2)))* &
1135 & tl_dzdx_r(i,j-1,k2)
1136 tl_cff8=(0.5_r8+sign(0.5_r8, dzdx_r(i,j ,k1)))* &
1137 & tl_dzdx_r(i,j ,k1)
1138 vfsx(i,j,k2)=vfsx(i,j,k2)- &
1139 & fac1* &
1140 & (cff1*(cff5*dnudz-dnudx(i,j-1,k1))+ &
1141 & cff2*(cff6*dnudz-dnudx(i,j ,k2))+ &
1142 & cff3*(cff7*dnudz-dnudx(i,j-1,k2))+ &
1143 & cff4*(cff8*dnudz-dnudx(i,j ,k1)))
1144 tl_vfsx(i,j,k2)=tl_vfsx(i,j,k2)- &
1145 & fac1* &
1146 & (tl_cff1*(cff5*dnudz-dnudx(i,j-1,k1))+ &
1147 & tl_cff2*(cff6*dnudz-dnudx(i,j ,k2))+ &
1148 & tl_cff3*(cff7*dnudz-dnudx(i,j-1,k2))+ &
1149 & tl_cff4*(cff8*dnudz-dnudx(i,j ,k1))+ &
1150 & cff1*(tl_cff5*dnudz+cff5*tl_dnudz- &
1151 & tl_dnudx(i,j-1,k1))+ &
1152 & cff2*(tl_cff6*dnudz+cff6*tl_dnudz- &
1153 & tl_dnudx(i,j ,k2))+ &
1154 & cff3*(tl_cff7*dnudz+cff7*tl_dnudz- &
1155 & tl_dnudx(i,j-1,k2))+ &
1156 & cff4*(tl_cff8*dnudz+cff8*tl_dnudz- &
1157 & tl_dnudx(i,j ,k1)))
1158#ifdef VISC_3DCOEF
1159 tl_vfsx(i,j,k2)=tl_vfsx(i,j,k2)- &
1160 & tl_fac1* &
1161 & (cff1*(cff5*dnudz-dnudx(i,j-1,k1))+ &
1162 & cff2*(cff6*dnudz-dnudx(i,j ,k2))+ &
1163 & cff3*(cff7*dnudz-dnudx(i,j-1,k2))+ &
1164 & cff4*(cff8*dnudz-dnudx(i,j ,k1)))
1165#endif
1166
1167 cff1=min(dzdx_p(i ,j,k1),0.0_r8)
1168 cff2=min(dzdx_p(i+1,j,k2),0.0_r8)
1169 cff3=max(dzdx_p(i ,j,k2),0.0_r8)
1170 cff4=max(dzdx_p(i+1,j,k1),0.0_r8)
1171 cff5=min(dzde_p(i ,j,k1),0.0_r8)
1172 cff6=min(dzde_p(i+1,j,k2),0.0_r8)
1173 cff7=max(dzde_p(i ,j,k2),0.0_r8)
1174 cff8=max(dzde_p(i+1,j,k1),0.0_r8)
1175 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_p(i ,j,k1)))* &
1176 & tl_dzdx_p(i ,j,k1)
1177 tl_cff2=(0.5_r8+sign(0.5_r8,-dzdx_p(i+1,j,k2)))* &
1178 & tl_dzdx_p(i+1,j,k2)
1179 tl_cff3=(0.5_r8+sign(0.5_r8, dzdx_p(i ,j,k2)))* &
1180 & tl_dzdx_p(i ,j,k2)
1181 tl_cff4=(0.5_r8+sign(0.5_r8, dzdx_p(i+1,j,k1)))* &
1182 & tl_dzdx_p(i+1,j,k1)
1183 tl_cff5=(0.5_r8+sign(0.5_r8,-dzde_p(i ,j,k1)))* &
1184 & tl_dzde_p(i ,j,k1)
1185 tl_cff6=(0.5_r8+sign(0.5_r8,-dzde_p(i+1,j,k2)))* &
1186 & tl_dzde_p(i+1,j,k2)
1187 tl_cff7=(0.5_r8+sign(0.5_r8, dzde_p(i ,j,k2)))* &
1188 & tl_dzde_p(i ,j,k2)
1189 tl_cff8=(0.5_r8+sign(0.5_r8, dzde_p(i+1,j,k1)))* &
1190 & tl_dzde_p(i+1,j,k1)
1191 vfse(i,j,k2)=vfse(i,j,k2)+ &
1192 & fac2* &
1193 & (cff1*(cff5*dmudz-dmude(i ,j,k1))+ &
1194 & cff2*(cff6*dmudz-dmude(i+1,j,k2))+ &
1195 & cff3*(cff7*dmudz-dmude(i ,j,k2))+ &
1196 & cff4*(cff8*dmudz-dmude(i+1,j,k1)))
1197 tl_vfse(i,j,k2)=tl_vfse(i,j,k2)+ &
1198 & fac2* &
1199 & (tl_cff1*(cff5*dmudz-dmude(i ,j,k1))+ &
1200 & tl_cff2*(cff6*dmudz-dmude(i+1,j,k2))+ &
1201 & tl_cff3*(cff7*dmudz-dmude(i ,j,k2))+ &
1202 & tl_cff4*(cff8*dmudz-dmude(i+1,j,k1))+ &
1203 & cff1*(tl_cff5*dmudz+cff5*tl_dmudz- &
1204 & tl_dmude(i ,j,k1))+ &
1205 & cff2*(tl_cff6*dmudz+cff6*tl_dmudz- &
1206 & tl_dmude(i+1,j,k2))+ &
1207 & cff3*(tl_cff7*dmudz+cff7*tl_dmudz- &
1208 & tl_dmude(i ,j,k2))+ &
1209 & cff4*(tl_cff8*dmudz+cff8*tl_dmudz- &
1210 & tl_dmude(i+1,j,k1)))
1211#ifdef VISC_3DCOEF
1212 tl_vfse(i,j,k2)=tl_vfse(i,j,k2)+ &
1213 & tl_fac2* &
1214 & (cff1*(cff5*dmudz-dmude(i ,j,k1))+ &
1215 & cff2*(cff6*dmudz-dmude(i+1,j,k2))+ &
1216 & cff3*(cff7*dmudz-dmude(i ,j,k2))+ &
1217 & cff4*(cff8*dmudz-dmude(i+1,j,k1)))
1218#endif
1219 END DO
1220 END DO
1221 END IF
1222!
1223! Compute first harmonic operator (m s^-3/2).
1224!
1225 DO j=jstrm1,jendp1
1226 DO i=istrum1,iendp1
1227 cff=0.125_r8*(pm(i-1,j)+pm(i,j))* &
1228 & (pn(i-1,j)+pn(i,j))
1229 cff1=1.0_r8/(0.5_r8*(hz(i-1,j,k)+hz(i,j,k)))
1230 tl_cff1=-cff1*cff1* &
1231 & (0.5_r8*(tl_hz(i-1,j,k)+tl_hz(i,j,k)))
1232 lapu(i,j,k)=cff*((pn(i-1,j)+pn(i,j))* &
1233 (ufx(i,j)-ufx(i-1,j))+ &
1234 & (pm(i-1,j)+pm(i,j))* &
1235 & (ufe(i,j+1)-ufe(i,j)))+ &
1236 & cff1*((ufsx(i,j,k2)+ufse(i,j,k2))- &
1237 & (ufsx(i,j,k1)+ufse(i,j,k1)))
1238 tl_lapu(i,j,k)=cff*((pn(i-1,j)+pn(i,j))* &
1239 (tl_ufx(i,j)-tl_ufx(i-1,j))+ &
1240 & (pm(i-1,j)+pm(i,j))* &
1241 & (tl_ufe(i,j+1)-tl_ufe(i,j)))+ &
1242 & tl_cff1*((ufsx(i,j,k2)+ufse(i,j,k2))- &
1243 & (ufsx(i,j,k1)+ufse(i,j,k1)))+ &
1244 & cff1*((tl_ufsx(i,j,k2)+tl_ufse(i,j,k2))- &
1245 & (tl_ufsx(i,j,k1)+tl_ufse(i,j,k1)))
1246#ifdef MASKING
1247 lapu(i,j,k)=lapu(i,j,k)*umask(i,j)
1248 tl_lapu(i,j,k)=tl_lapu(i,j,k)*umask(i,j)
1249#endif
1250 END DO
1251 END DO
1252
1253 DO j=jstrvm1,jendp1
1254 DO i=istrm1,iendp1
1255 cff=0.125_r8*(pm(i,j)+pm(i,j-1))* &
1256 & (pn(i,j)+pn(i,j-1))
1257 cff1=1.0_r8/(0.5_r8*(hz(i,j-1,k)+hz(i,j,k)))
1258 tl_cff1=-cff1*cff1* &
1259 & (0.5_r8*(tl_hz(i,j-1,k)+tl_hz(i,j,k)))
1260 lapv(i,j,k)=cff*((pn(i,j-1)+pn(i,j))* &
1261 & (vfx(i+1,j)-vfx(i,j))- &
1262 & (pm(i,j-1)+pm(i,j))* &
1263 & (vfe(i,j)-vfe(i,j-1)))+ &
1264 & cff1*((vfsx(i,j,k2)+vfse(i,j,k2))- &
1265 & (vfsx(i,j,k1)+vfse(i,j,k1)))
1266 tl_lapv(i,j,k)=cff*((pn(i,j-1)+pn(i,j))* &
1267 & (tl_vfx(i+1,j)-tl_vfx(i,j))- &
1268 & (pm(i,j-1)+pm(i,j))* &
1269 & (tl_vfe(i,j)-tl_vfe(i,j-1)))+ &
1270 & tl_cff1*((vfsx(i,j,k2)+vfse(i,j,k2))- &
1271 & (vfsx(i,j,k1)+vfse(i,j,k1)))+ &
1272 & cff1*((tl_vfsx(i,j,k2)+tl_vfse(i,j,k2))- &
1273 & (tl_vfsx(i,j,k1)+tl_vfse(i,j,k1)))
1274#ifdef MASKING
1275 lapv(i,j,k)=lapv(i,j,k)*vmask(i,j)
1276 tl_lapv(i,j,k)=tl_lapv(i,j,k)*vmask(i,j)
1277#endif
1278 END DO
1279 END DO
1280 END IF
1281 END DO k_loop1
1282!
1283! Apply boundary conditions (closed or gradient; except periodic)
1284! to the first harmonic operator.
1285!
1286 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
1287 IF (domain(ng)%Western_Edge(tile)) THEN
1288 IF (tl_lbc(iwest,isuvel,ng)%closed) THEN
1289 DO k=1,n(ng)
1290 DO j=jstrm1,jendp1
1291 lapu(istru-1,j,k)=0.0_r8
1292 tl_lapu(istru-1,j,k)=0.0_r8
1293 END DO
1294 END DO
1295 ELSE
1296 DO k=1,n(ng)
1297 DO j=jstrm1,jendp1
1298 lapu(istru-1,j,k)=lapu(istru,j,k)
1299 tl_lapu(istru-1,j,k)=tl_lapu(istru,j,k)
1300 END DO
1301 END DO
1302 END IF
1303 IF (tl_lbc(iwest,isvvel,ng)%closed) THEN
1304 DO k=1,n(ng)
1305 DO j=jstrvm1,jendp1
1306 lapv(istr-1,j,k)=gamma2(ng)*lapv(istr,j,k)
1307 tl_lapv(istr-1,j,k)=gamma2(ng)*tl_lapv(istr,j,k)
1308 END DO
1309 END DO
1310 ELSE
1311 DO k=1,n(ng)
1312 DO j=jstrvm1,jendp1
1313 lapv(istr-1,j,k)=0.0_r8
1314 tl_lapv(istr-1,j,k)=0.0_r8
1315 END DO
1316 END DO
1317 END IF
1318 END IF
1319 END IF
1320!
1321 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
1322 IF (domain(ng)%Eastern_Edge(tile)) THEN
1323 IF (tl_lbc(ieast,isuvel,ng)%closed) THEN
1324 DO k=1,n(ng)
1325 DO j=jstrm1,jendp1
1326 lapu(iend+1,j,k)=0.0_r8
1327 tl_lapu(iend+1,j,k)=0.0_r8
1328 END DO
1329 END DO
1330 ELSE
1331 DO k=1,n(ng)
1332 DO j=jstrm1,jendp1
1333 lapu(iend+1,j,k)=lapu(iend,j,k)
1334 tl_lapu(iend+1,j,k)=tl_lapu(iend,j,k)
1335 END DO
1336 END DO
1337 END IF
1338 IF (tl_lbc(ieast,isvvel,ng)%closed) THEN
1339 DO k=1,n(ng)
1340 DO j=jstrvm1,jendp1
1341 lapv(iend+1,j,k)=gamma2(ng)*lapv(iend,j,k)
1342 tl_lapv(iend+1,j,k)=gamma2(ng)*tl_lapv(iend,j,k)
1343 END DO
1344 END DO
1345 ELSE
1346 DO k=1,n(ng)
1347 DO j=jstrvm1,jendp1
1348 lapv(iend+1,j,k)=0.0_r8
1349 tl_lapv(iend+1,j,k)=0.0_r8
1350 END DO
1351 END DO
1352 END IF
1353 END IF
1354 END IF
1355!
1356 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
1357 IF (domain(ng)%Southern_Edge(tile)) THEN
1358 IF (tl_lbc(isouth,isuvel,ng)%closed) THEN
1359 DO k=1,n(ng)
1360 DO i=istrum1,iendp1
1361 lapu(i,jstr-1,k)=gamma2(ng)*lapu(i,jstr,k)
1362 tl_lapu(i,jstr-1,k)=gamma2(ng)*tl_lapu(i,jstr,k)
1363 END DO
1364 END DO
1365 ELSE
1366 DO k=1,n(ng)
1367 DO i=istrum1,iendp1
1368 lapu(i,jstr-1,k)=0.0_r8
1369 tl_lapu(i,jstr-1,k)=0.0_r8
1370 END DO
1371 END DO
1372 END IF
1373 IF (tl_lbc(isouth,isvvel,ng)%closed) THEN
1374 DO k=1,n(ng)
1375 DO i=istrm1,iendp1
1376 lapv(i,jstrv-1,k)=0.0_r8
1377 tl_lapv(i,jstrv-1,k)=0.0_r8
1378 END DO
1379 END DO
1380 ELSE
1381 DO k=1,n(ng)
1382 DO i=istrm1,iendp1
1383 lapv(i,jstrv-1,k)=lapv(i,jstrv,k)
1384 tl_lapv(i,jstrv-1,k)=tl_lapv(i,jstrv,k)
1385 END DO
1386 END DO
1387 END IF
1388 END IF
1389 END IF
1390!
1391 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
1392 IF (domain(ng)%Northern_Edge(tile)) THEN
1393 IF (tl_lbc(inorth,isuvel,ng)%closed) THEN
1394 DO k=1,n(ng)
1395 DO i=istrum1,iendp1
1396 lapu(i,jend+1,k)=gamma2(ng)*lapu(i,jend,k)
1397 tl_lapu(i,jend+1,k)=gamma2(ng)*tl_lapu(i,jend,k)
1398 END DO
1399 END DO
1400 ELSE
1401 DO k=1,n(ng)
1402 DO i=istrum1,iendp1
1403 lapu(i,jend+1,k)=0.0_r8
1404 tl_lapu(i,jend+1,k)=0.0_r8
1405 END DO
1406 END DO
1407 END IF
1408 IF (tl_lbc(inorth,isvvel,ng)%closed) THEN
1409 DO k=1,n(ng)
1410 DO i=istrm1,iendp1
1411 lapv(i,jend+1,k)=0.0_r8
1412 tl_lapv(i,jend+1,k)=0.0_r8
1413 END DO
1414 END DO
1415 ELSE
1416 DO k=1,n(ng)
1417 DO i=istrm1,iendp1
1418 lapv(i,jend+1,k)=lapv(i,jend,k)
1419 tl_lapv(i,jend+1,k)=tl_lapv(i,jend,k)
1420 END DO
1421 END DO
1422 END IF
1423 END IF
1424 END IF
1425!
1426 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng).or. &
1427 & compositegrid(iwest ,ng).or.ewperiodic(ng))) THEN
1428 IF (domain(ng)%SouthWest_Corner(tile)) THEN
1429 DO k=1,n(ng)
1430 lapu(istr ,jstr-1,k)=0.5_r8* &
1431 & (lapu(istr+1,jstr-1,k)+ &
1432 & lapu(istr ,jstr ,k))
1433 tl_lapu(istr ,jstr-1,k)=0.5_r8* &
1434 & (tl_lapu(istr+1,jstr-1,k)+ &
1435 & tl_lapu(istr ,jstr ,k))
1436 lapv(istr-1,jstr ,k)=0.5_r8* &
1437 & (lapv(istr-1,jstr+1,k)+ &
1438 & lapv(istr ,jstr ,k))
1439 tl_lapv(istr-1,jstr ,k)=0.5_r8* &
1440 & (tl_lapv(istr-1,jstr+1,k)+ &
1441 & tl_lapv(istr ,jstr ,k))
1442 END DO
1443 END IF
1444 END IF
1445
1446 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng).or. &
1447 & compositegrid(ieast ,ng).or.ewperiodic(ng))) THEN
1448 IF (domain(ng)%SouthEast_Corner(tile)) THEN
1449 DO k=1,n(ng)
1450 lapu(iend+1,jstr-1,k)=0.5_r8* &
1451 & (lapu(iend ,jstr-1,k)+ &
1452 & lapu(iend+1,jstr ,k))
1453 tl_lapu(iend+1,jstr-1,k)=0.5_r8* &
1454 & (tl_lapu(iend ,jstr-1,k)+ &
1455 & tl_lapu(iend+1,jstr ,k))
1456 lapv(iend+1,jstr ,k)=0.5_r8* &
1457 & (lapv(iend ,jstr ,k)+ &
1458 & lapv(iend+1,jstr+1,k))
1459 tl_lapv(iend+1,jstr ,k)=0.5_r8* &
1460 & (tl_lapv(iend ,jstr ,k)+ &
1461 & tl_lapv(iend+1,jstr+1,k))
1462 END DO
1463 END IF
1464 END IF
1465
1466 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng).or. &
1467 & compositegrid(iwest ,ng).or.ewperiodic(ng))) THEN
1468 IF (domain(ng)%NorthWest_Corner(tile)) THEN
1469 DO k=1,n(ng)
1470 lapu(istr ,jend+1,k)=0.5_r8* &
1471 & (lapu(istr+1,jend+1,k)+ &
1472 & lapu(istr ,jend ,k))
1473 tl_lapu(istr ,jend+1,k)=0.5_r8* &
1474 & (tl_lapu(istr+1,jend+1,k)+ &
1475 & tl_lapu(istr ,jend ,k))
1476 lapv(istr-1,jend+1,k)=0.5_r8* &
1477 & (lapv(istr ,jend+1,k)+ &
1478 & lapv(istr-1,jend ,k))
1479 tl_lapv(istr-1,jend+1,k)=0.5_r8* &
1480 & (tl_lapv(istr ,jend+1,k)+ &
1481 & tl_lapv(istr-1,jend ,k))
1482 END DO
1483 END IF
1484 END IF
1485
1486 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng).or. &
1487 & compositegrid(ieast ,ng).or.ewperiodic(ng))) THEN
1488 IF (domain(ng)%NorthEast_Corner(tile)) THEN
1489 DO k=1,n(ng)
1490 lapu(iend+1,jend+1,k)=0.5_r8* &
1491 & (lapu(iend ,jend+1,k)+ &
1492 & lapu(iend+1,jend ,k))
1493 tl_lapu(iend+1,jend+1,k)=0.5_r8* &
1494 & (tl_lapu(iend ,jend+1,k)+ &
1495 & tl_lapu(iend+1,jend ,k))
1496 lapv(iend+1,jend+1,k)=0.5_r8* &
1497 & (lapv(iend ,jend+1,k)+ &
1498 & lapv(iend+1,jend ,k))
1499 tl_lapv(iend+1,jend+1,k)=0.5_r8* &
1500 & (tl_lapv(iend ,jend+1,k)+ &
1501 & tl_lapv(iend+1,jend ,k))
1502 END DO
1503 END IF
1504 END IF
1505!
1506! Compute horizontal and vertical gradients associated with the
1507! second rotated harmonic operator.
1508!
1509 k2=1
1510 k_loop2 : DO k=0,n(ng)
1511 k1=k2
1512 k2=3-k1
1513 IF (k.lt.n(ng)) THEN
1514!
1515! Compute slopes (nondimensional) at RHO- and PSI-points.
1516!
1517 DO j=jstr-1,jend+1
1518 DO i=istru-1,iend+1
1519 cff=0.5_r8*(pm(i-1,j)+pm(i,j))
1520#ifdef MASKING
1521 cff=cff*umask(i,j)
1522#endif
1523 ufx(i,j)=cff*(z_r(i ,j,k+1)- &
1524 & z_r(i-1,j,k+1))
1525 tl_ufx(i,j)=cff*(tl_z_r(i ,j,k+1)- &
1526 & tl_z_r(i-1,j,k+1))
1527 END DO
1528 END DO
1529 DO j=jstrv-1,jend+1
1530 DO i=istr-1,iend+1
1531 cff=0.5_r8*(pn(i,j-1)+pn(i,j))
1532#ifdef MASKING
1533 cff=cff*vmask(i,j)
1534#endif
1535 vfe(i,j)=cff*(z_r(i,j ,k+1)- &
1536 & z_r(i,j-1,k+1))
1537 tl_vfe(i,j)=cff*(tl_z_r(i,j ,k+1)- &
1538 & tl_z_r(i,j-1,k+1))
1539 END DO
1540 END DO
1541!
1542 DO j=jstr,jend+1
1543 DO i=istr,iend+1
1544 dzdx_p(i,j,k2)=0.5_r8*(ufx(i,j-1)+ &
1545 & ufx(i,j ))
1546 tl_dzdx_p(i,j,k2)=0.5_r8*(tl_ufx(i,j-1)+ &
1547 & tl_ufx(i,j ))
1548 dzde_p(i,j,k2)=0.5_r8*(vfe(i-1,j)+ &
1549 & vfe(i ,j))
1550 tl_dzde_p(i,j,k2)=0.5_r8*(tl_vfe(i-1,j)+ &
1551 & tl_vfe(i ,j))
1552 END DO
1553 END DO
1554 DO j=jstrv-1,jend
1555 DO i=istru-1,iend
1556 dzdx_r(i,j,k2)=0.5_r8*(ufx(i ,j)+ &
1557 & ufx(i+1,j))
1558 tl_dzdx_r(i,j,k2)=0.5_r8*(tl_ufx(i ,j)+ &
1559 & tl_ufx(i+1,j))
1560 dzde_r(i,j,k2)=0.5_r8*(vfe(i,j )+ &
1561 & vfe(i,j+1))
1562 tl_dzde_r(i,j,k2)=0.5_r8*(tl_vfe(i,j )+ &
1563 & tl_vfe(i,j+1))
1564 END DO
1565 END DO
1566!
1567! Compute momentum horizontal (m^-1 s^-3/2) and vertical (s^-3/2)
1568! gradients.
1569!
1570 DO j=jstrv-1,jend
1571 DO i=istru-1,iend
1572 cff=0.5_r8*pm(i,j)
1573#ifdef MASKING
1574 cff=cff*rmask(i,j)
1575#endif
1576 dnudx(i,j,k2)=cff*((pn(i ,j)+pn(i+1,j))* &
1577 & lapu(i+1,j,k+1)- &
1578 & (pn(i-1,j)+pn(i ,j))* &
1579 & lapu(i ,j,k+1))
1580 tl_dnudx(i,j,k2)=cff*((pn(i ,j)+pn(i+1,j))* &
1581 & tl_lapu(i+1,j,k+1)- &
1582 & (pn(i-1,j)+pn(i ,j))* &
1583 & tl_lapu(i ,j,k+1))
1584 END DO
1585 END DO
1586
1587 DO j=jstr,jend+1
1588 DO i=istr,iend+1
1589 cff=0.125_r8*(pn(i-1,j )+pn(i,j )+ &
1590 & pn(i-1,j-1)+pn(i,j-1))
1591#ifdef MASKING
1592 cff=cff*pmask(i,j)
1593#endif
1594 dmude(i,j,k2)=cff*((pm(i-1,j )+pm(i,j ))* &
1595 & lapu(i,j ,k+1)- &
1596 & (pm(i-1,j-1)+pm(i,j-1))* &
1597 & lapu(i,j-1,k+1))
1598 tl_dmude(i,j,k2)=cff*((pm(i-1,j )+pm(i,j ))* &
1599 & tl_lapu(i,j ,k+1)- &
1600 & (pm(i-1,j-1)+pm(i,j-1))* &
1601 & tl_lapu(i,j-1,k+1))
1602 END DO
1603 END DO
1604
1605 DO j=jstr,jend+1
1606 DO i=istr,iend+1
1607 cff=0.125_r8*(pm(i-1,j )+pm(i,j )+ &
1608 & pm(i-1,j-1)+pm(i,j-1))
1609#ifdef MASKING
1610 cff=cff*pmask(i,j)
1611#endif
1612 dnvdx(i,j,k2)=cff*((pn(i ,j-1)+pn(i ,j))* &
1613 & lapv(i ,j,k+1)- &
1614 & (pn(i-1,j-1)+pn(i-1,j))* &
1615 & lapv(i-1,j,k+1))
1616 tl_dnvdx(i,j,k2)=cff*((pn(i ,j-1)+pn(i ,j))* &
1617 & tl_lapv(i ,j,k+1)- &
1618 & (pn(i-1,j-1)+pn(i-1,j))* &
1619 & tl_lapv(i-1,j,k+1))
1620 END DO
1621 END DO
1622
1623 DO j=jstrv-1,jend
1624 DO i=istru-1,iend
1625 cff=0.5_r8*pn(i,j)
1626#ifdef MASKING
1627 cff=cff*rmask(i,j)
1628#endif
1629 dmvde(i,j,k2)=cff*((pm(i,j )+pm(i,j+1))* &
1630 & lapv(i,j+1,k+1)- &
1631 & (pm(i,j-1)+pm(i,j ))* &
1632 & lapv(i,j ,k+1))
1633 tl_dmvde(i,j,k2)=cff*((pm(i,j )+pm(i,j+1))* &
1634 & tl_lapv(i,j+1,k+1)- &
1635 & (pm(i,j-1)+pm(i,j ))* &
1636 & tl_lapv(i,j ,k+1))
1637 END DO
1638 END DO
1639 END IF
1640
1641 IF ((k.eq.0).or.(k.eq.n(ng))) THEN
1642 DO j=jstr-1,jend+1
1643 DO i=istru-1,iend+1
1644 dudz(i,j,k2)=0.0_r8
1645 tl_dudz(i,j,k2)=0.0_r8
1646 END DO
1647 END DO
1648 DO j=jstrv-1,jend+1
1649 DO i=istr-1,iend+1
1650 dvdz(i,j,k2)=0.0_r8
1651 tl_dvdz(i,j,k2)=0.0_r8
1652 END DO
1653 END DO
1654
1655 DO j=jstr,jend
1656 DO i=istru,iend
1657 ufsx(i,j,k2)=0.0_r8
1658 tl_ufsx(i,j,k2)=0.0_r8
1659 ufse(i,j,k2)=0.0_r8
1660 tl_ufse(i,j,k2)=0.0_r8
1661 END DO
1662 END DO
1663 DO j=jstrv,jend
1664 DO i=istr,iend
1665 vfsx(i,j,k2)=0.0_r8
1666 tl_vfsx(i,j,k2)=0.0_r8
1667 vfse(i,j,k2)=0.0_r8
1668 tl_vfse(i,j,k2)=0.0_r8
1669 END DO
1670 END DO
1671 ELSE
1672 DO j=jstr-1,jend+1
1673 DO i=istru-1,iend+1
1674 cff=1.0_r8/(0.5_r8*(z_r(i-1,j,k+1)- &
1675 & z_r(i-1,j,k )+ &
1676 & z_r(i ,j,k+1)- &
1677 & z_r(i ,j,k )))
1678 tl_cff=-cff*cff*(0.5_r8*(tl_z_r(i-1,j,k+1)- &
1679 & tl_z_r(i-1,j,k )+ &
1680 & tl_z_r(i ,j,k+1)- &
1681 & tl_z_r(i ,j,k )))
1682 dudz(i,j,k2)=cff*(lapu(i,j,k+1)- &
1683 & lapu(i,j,k ))
1684 tl_dudz(i,j,k2)=tl_cff*(lapu(i,j,k+1)- &
1685 & lapu(i,j,k ))+ &
1686 & cff*(tl_lapu(i,j,k+1)- &
1687 & tl_lapu(i,j,k ))
1688 END DO
1689 END DO
1690 DO j=jstrv-1,jend+1
1691 DO i=istr-1,iend+1
1692 cff=1.0_r8/(0.5_r8*(z_r(i,j-1,k+1)- &
1693 & z_r(i,j-1,k )+ &
1694 & z_r(i,j ,k+1)- &
1695 & z_r(i,j ,k )))
1696 tl_cff=-cff*cff*(0.5_r8*(tl_z_r(i,j-1,k+1)- &
1697 & tl_z_r(i,j-1,k )+ &
1698 & tl_z_r(i,j ,k+1)- &
1699 & tl_z_r(i,j ,k )))
1700 dvdz(i,j,k2)=cff*(lapv(i,j,k+1)- &
1701 & lapv(i,j,k ))
1702 tl_dvdz(i,j,k2)=tl_cff*(lapv(i,j,k+1)- &
1703 & lapv(i,j,k ))+ &
1704 & cff*(tl_lapv(i,j,k+1)- &
1705 & tl_lapv(i,j,k ))
1706 END DO
1707 END DO
1708 END IF
1709!
1710! Compute components of the rotated viscous flux (m5/s2) along
1711! geopotential surfaces in the XI- and ETA-directions.
1712!
1713 IF (k.gt.0) THEN
1714 DO j=jstrv-1,jend
1715 DO i=istru-1,iend
1716 cff1=min(dzdx_r(i,j,k1),0.0_r8)
1717 cff2=max(dzdx_r(i,j,k1),0.0_r8)
1718 cff3=min(dzde_r(i,j,k1),0.0_r8)
1719 cff4=max(dzde_r(i,j,k1),0.0_r8)
1720 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_r(i,j,k1)))* &
1721 & tl_dzdx_r(i,j,k1)
1722 tl_cff2=(0.5_r8+sign(0.5_r8, dzdx_r(i,j,k1)))* &
1723 & tl_dzdx_r(i,j,k1)
1724 tl_cff3=(0.5_r8+sign(0.5_r8,-dzde_r(i,j,k1)))* &
1725 & tl_dzde_r(i,j,k1)
1726 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_r(i,j,k1)))* &
1727 & tl_dzde_r(i,j,k1)
1728#ifdef VISC_3DCOEF
1729 cff=hz(i,j,k)* &
1730 & (on_r(i,j)*(dnudx(i,j,k1)- &
1731 & 0.5_r8*pn(i,j)* &
1732 & (cff1*(dudz(i ,j,k1)+ &
1733 & dudz(i+1,j,k2))+ &
1734 & cff2*(dudz(i ,j,k2)+ &
1735 & dudz(i+1,j,k1))))- &
1736 & om_r(i,j)*(dmvde(i,j,k1)- &
1737 & 0.5_r8*pm(i,j)* &
1738 & (cff3*(dvdz(i,j ,k1)+ &
1739 & dvdz(i,j+1,k2))+ &
1740 & cff4*(dvdz(i,j ,k2)+ &
1741 & dvdz(i,j+1,k1)))))
1742#else
1743!^ cff=Hz(i,j,k)* &
1744!^ & (on_r(i,j)*(dnUdx(i,j,k1)- &
1745!^ & 0.5_r8*pn(i,j)* &
1746!^ & (cff1*(dUdz(i ,j,k1)+ &
1747!^ & dUdz(i+1,j,k2))+ &
1748!^ & cff2*(dUdz(i ,j,k2)+ &
1749!^ & dUdz(i+1,j,k1))))- &
1750!^ & om_r(i,j)*(dmVde(i,j,k1)- &
1751!^ & 0.5_r8*pm(i,j)* &
1752!^ & (cff3*(dVdz(i,j ,k1)+ &
1753!^ & dVdz(i,j+1,k2))+ &
1754!^ & cff4*(dVdz(i,j ,k2)+ &
1755!^ & dVdz(i,j+1,k1)))))
1756!^
1757#endif
1758 tl_cff=tl_hz(i,j,k)* &
1759 & (on_r(i,j)*(dnudx(i,j,k1)- &
1760 & 0.5_r8*pn(i,j)* &
1761 & (cff1*(dudz(i ,j,k1)+ &
1762 & dudz(i+1,j,k2))+ &
1763 & cff2*(dudz(i ,j,k2)+ &
1764 & dudz(i+1,j,k1))))- &
1765 & om_r(i,j)*(dmvde(i,j,k1)- &
1766 & 0.5_r8*pm(i,j)* &
1767 & (cff3*(dvdz(i,j ,k1)+ &
1768 & dvdz(i,j+1,k2))+ &
1769 & cff4*(dvdz(i,j ,k2)+ &
1770 & dvdz(i,j+1,k1)))))+ &
1771 & hz(i,j,k)* &
1772 & (on_r(i,j)*(tl_dnudx(i,j,k1)- &
1773 & 0.5_r8*pn(i,j)* &
1774 & (tl_cff1*(dudz(i ,j,k1)+ &
1775 & dudz(i+1,j,k2))+ &
1776 & cff1*(tl_dudz(i ,j,k1)+ &
1777 & tl_dudz(i+1,j,k2))+ &
1778 & tl_cff2*(dudz(i ,j,k2)+ &
1779 & dudz(i+1,j,k1))+ &
1780 & cff2*(tl_dudz(i ,j,k2)+ &
1781 & tl_dudz(i+1,j,k1))))- &
1782 & om_r(i,j)*(tl_dmvde(i,j,k1)- &
1783 & 0.5_r8*pm(i,j)* &
1784 & (tl_cff3*(dvdz(i,j ,k1)+ &
1785 & dvdz(i,j+1,k2))+ &
1786 & cff3*(tl_dvdz(i,j ,k1)+ &
1787 & tl_dvdz(i,j+1,k2))+ &
1788 & tl_cff4*(dvdz(i,j ,k2)+ &
1789 & dvdz(i,j+1,k1))+ &
1790 & cff4*(tl_dvdz(i,j ,k2)+ &
1791 & tl_dvdz(i,j+1,k1)))))
1792#ifdef MASKING
1793# ifdef VISC_3DCOEF
1794 cff=cff*rmask(i,j)
1795# else
1796!^ cff=cff*rmask(i,j)
1797!^
1798# endif
1799 tl_cff=tl_cff*rmask(i,j)
1800#endif
1801#ifdef VISC_3DCOEF
1802# ifdef UV_U3ADV_SPLIT
1803!^ UFx(i,j)=on_r(i,j)*on_r(i,j)*Uvis3d_r(i,j,k)*cff
1804!^
1805 tl_ufx(i,j)=on_r(i,j)*on_r(i,j)* &
1806 & (tl_uvis3d_r(i,j,k)*cff+ &
1807 & uvis3d_r(i,j,k)*tl_cff)
1808!^ VFe(i,j)=om_r(i,j)*om_r(i,j)*Vvis3d_r(i,j,k)*cff
1809!^
1810 tl_vfe(i,j)=om_r(i,j)*om_r(i,j)* &
1811 & (tl_vvis3d_r(i,j,k)*cff+ &
1812 & vvis3d_r(i,j,k)*tl_cff)
1813# else
1814!^ UFx(i,j)=on_r(i,j)*on_r(i,j)*visc3d_r(i,j,k)*cff
1815!^
1816 tl_ufx(i,j)=on_r(i,j)*on_r(i,j)* &
1817 & (tl_visc3d_r(i,j,k)*cff+ &
1818 & visc3d_r(i,j,k)*tl_cff)
1819!^ VFe(i,j)=om_r(i,j)*om_r(i,j)*visc3d_r(i,j,k)*cff
1820!^
1821 tl_vfe(i,j)=om_r(i,j)*om_r(i,j)* &
1822 & (tl_visc3d_r(i,j,k)*cff+ &
1823 & visc3d_r(i,j,k)*tl_cff)
1824# endif
1825#else
1826!^ UFx(i,j)=on_r(i,j)*on_r(i,j)*visc4_r(i,j)*cff
1827!^
1828 tl_ufx(i,j)=on_r(i,j)*on_r(i,j)*visc4_r(i,j)*tl_cff
1829!^ VFe(i,j)=om_r(i,j)*om_r(i,j)*visc4_r(i,j)*cff
1830!^
1831 tl_vfe(i,j)=om_r(i,j)*om_r(i,j)*visc4_r(i,j)*tl_cff
1832#endif
1833 END DO
1834 END DO
1835
1836 DO j=jstr,jend+1
1837 DO i=istr,iend+1
1838 pm_p=0.25_r8*(pm(i-1,j-1)+pm(i-1,j)+ &
1839 & pm(i ,j-1)+pm(i ,j))
1840 pn_p=0.25_r8*(pn(i-1,j-1)+pn(i-1,j)+ &
1841 & pn(i ,j-1)+pn(i ,j))
1842 cff1=min(dzdx_p(i,j,k1),0.0_r8)
1843 cff2=max(dzdx_p(i,j,k1),0.0_r8)
1844 cff3=min(dzde_p(i,j,k1),0.0_r8)
1845 cff4=max(dzde_p(i,j,k1),0.0_r8)
1846 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_p(i,j,k1)))* &
1847 & tl_dzdx_p(i,j,k1)
1848 tl_cff2=(0.5_r8+sign(0.5_r8, dzdx_p(i,j,k1)))* &
1849 & tl_dzdx_p(i,j,k1)
1850 tl_cff3=(0.5_r8+sign(0.5_r8,-dzde_p(i,j,k1)))* &
1851 & tl_dzde_p(i,j,k1)
1852 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_p(i,j,k1)))* &
1853 & tl_dzde_p(i,j,k1)
1854#ifdef VISC_3DCOEF
1855 cff=0.25_r8* &
1856 & (hz(i-1,j ,k)+hz(i,j ,k)+ &
1857 & hz(i-1,j-1,k)+hz(i,j-1,k))* &
1858 & (on_p(i,j)*(dnvdx(i,j,k1)- &
1859 & 0.5_r8*pn_p* &
1860 & (cff1*(dvdz(i-1,j,k1)+ &
1861 & dvdz(i ,j,k2))+ &
1862 & cff2*(dvdz(i-1,j,k2)+ &
1863 & dvdz(i ,j,k1))))+ &
1864 & om_p(i,j)*(dmude(i,j,k1)- &
1865 & 0.5_r8*pm_p* &
1866 & (cff3*(dudz(i,j-1,k1)+ &
1867 & dudz(i,j ,k2))+ &
1868 & cff4*(dudz(i,j-1,k2)+ &
1869 & dudz(i,j ,k1)))))
1870#else
1871!^ cff=0.25_r8* &
1872!^ & (Hz(i-1,j ,k)+Hz(i,j ,k)+ &
1873!^ & Hz(i-1,j-1,k)+Hz(i,j-1,k))* &
1874!^ & (on_p(i,j)*(dnVdx(i,j,k1)- &
1875!^ & 0.5_r8*pn_p* &
1876!^ & (cff1*(dVdz(i-1,j,k1)+ &
1877!^ & dVdz(i ,j,k2))+ &
1878!^ & cff2*(dVdz(i-1,j,k2)+ &
1879!^ & dVdz(i ,j,k1))))+ &
1880!^ & om_p(i,j)*(dmUde(i,j,k1)- &
1881!^ & 0.5_r8*pm_p* &
1882!^ & (cff3*(dUdz(i,j-1,k1)+ &
1883!^ & dUdz(i,j ,k2))+ &
1884!^ & cff4*(dUdz(i,j-1,k2)+ &
1885!^ & dUdz(i,j ,k1)))))
1886!^
1887#endif
1888 tl_cff=0.25_r8* &
1889 & ((tl_hz(i-1,j ,k)+tl_hz(i,j ,k)+ &
1890 & tl_hz(i-1,j-1,k)+tl_hz(i,j-1,k))* &
1891 & (on_p(i,j)*(dnvdx(i,j,k1)- &
1892 & 0.5_r8*pn_p* &
1893 & (cff1*(dvdz(i-1,j,k1)+ &
1894 & dvdz(i ,j,k2))+ &
1895 & cff2*(dvdz(i-1,j,k2)+ &
1896 & dvdz(i ,j,k1))))+ &
1897 & om_p(i,j)*(dmude(i,j,k1)- &
1898 & 0.5_r8*pm_p* &
1899 & (cff3*(dudz(i,j-1,k1)+ &
1900 & dudz(i,j ,k2))+ &
1901 & cff4*(dudz(i,j-1,k2)+ &
1902 & dudz(i,j ,k1)))))+ &
1903 & (hz(i-1,j ,k)+hz(i,j ,k)+ &
1904 & hz(i-1,j-1,k)+hz(i,j-1,k))* &
1905 & (on_p(i,j)*(tl_dnvdx(i,j,k1)- &
1906 & 0.5_r8*pn_p* &
1907 & (tl_cff1*(dvdz(i-1,j,k1)+ &
1908 & dvdz(i ,j,k2))+ &
1909 & cff1*(tl_dvdz(i-1,j,k1)+ &
1910 & tl_dvdz(i ,j,k2))+ &
1911 & tl_cff2*(dvdz(i-1,j,k2)+ &
1912 & dvdz(i ,j,k1))+ &
1913 & cff2*(tl_dvdz(i-1,j,k2)+ &
1914 & tl_dvdz(i ,j,k1))))+ &
1915 & om_p(i,j)*(tl_dmude(i,j,k1)- &
1916 & 0.5_r8*pm_p* &
1917 & (tl_cff3*(dudz(i,j-1,k1)+ &
1918 & dudz(i,j ,k2))+ &
1919 & cff3*(tl_dudz(i,j-1,k1)+ &
1920 & tl_dudz(i,j ,k2))+ &
1921 & tl_cff4*(dudz(i,j-1,k2)+ &
1922 & dudz(i,j ,k1))+ &
1923 & cff4*(tl_dudz(i,j-1,k2)+ &
1924 & tl_dudz(i,j ,k1))))))
1925#ifdef MASKING
1926# ifdef VISC_3DCOEF
1927 cff=cff*pmask(i,j)
1928# else
1929!^ cff=cff*pmask(i,j)
1930!^
1931# endif
1932 tl_cff=tl_cff*pmask(i,j)
1933#endif
1934#ifdef VISC_3DCOEF
1935# ifdef UV_U3ADV_SPLIT
1936 uvis_p=0.25_r8* &
1937 & (uvis3d_r(i-1,j-1,k)+uvis3d_r(i-1,j,k)+ &
1938 & uvis3d_r(i ,j-1,k)+uvis3d_r(i ,j,k))
1939 tl_uvis_p=0.25_r8* &
1940 & (tl_uvis3d_r(i-1,j-1,k)+tl_uvis3d_r(i-1,j,k)+ &
1941 & tl_uvis3d_r(i ,j-1,k)+tl_uvis3d_r(i ,j,k))
1942 vvis_p=0.25_r8* &
1943 & (vvis3d_r(i-1,j-1,k)+vvis3d_r(i-1,j,k)+ &
1944 & vvis3d_r(i ,j-1,k)+vvis3d_r(i ,j,k))
1945 tl_vvis_p=0.25_r8* &
1946 & (tl_vvis3d_r(i-1,j-1,k)+tl_vvis3d_r(i-1,j,k)+ &
1947 & tl_vvis3d_r(i ,j-1,k)+tl_vvis3d_r(i ,j,k))
1948!^ UFe(i,j)=om_p(i,j)*om_p(i,j)*Uvis_p*cff
1949!^
1950 tl_ufe(i,j)=om_p(i,j)*om_p(i,j)* &
1951 & (tl_uvis_p*cff+uvis_p*tl_cff)
1952!^ VFx(i,j)=on_p(i,j)*on_p(i,j)*Vvis_p*cff
1953!^
1954 tl_vfx(i,j)=on_p(i,j)*on_p(i,j)* &
1955 & (tl_vvis_p*cff+vvis_p*tl_cff)
1956# else
1957 visc_p=0.25_r8* &
1958 & (visc3d_r(i-1,j-1,k)+visc3d_r(i-1,j,k)+ &
1959 & visc3d_r(i ,j-1,k)+visc3d_r(i ,j,k))
1960 tl_visc_p=0.25_r8* &
1961 & (tl_visc3d_r(i-1,j-1,k)+tl_visc3d_r(i-1,j,k)+ &
1962 & tl_visc3d_r(i ,j-1,k)+tl_visc3d_r(i ,j,k))
1963!^ UFe(i,j)=om_p(i,j)*om_p(i,j)*visc_p*cff
1964!^
1965 tl_ufe(i,j)=om_p(i,j)*om_p(i,j)* &
1966 & (tl_visc_p*cff+visc_p*tl_cff)
1967!^ VFx(i,j)=on_p(i,j)*on_p(i,j)*visc_p*cff
1968!^
1969 tl_vfx(i,j)=on_p(i,j)*on_p(i,j)* &
1970 & (tl_visc_p*cff+visc_p*tl_cff)
1971# endif
1972#else
1973!^ UFe(i,j)=om_p(i,j)*om_p(i,j)*visc4_p(i,j)*cff
1974!^
1975 tl_ufe(i,j)=om_p(i,j)*om_p(i,j)*visc4_p(i,j)*tl_cff
1976!^ VFx(i,j)=on_p(i,j)*on_p(i,j)*visc4_p(i,j)*cff
1977!^
1978 tl_vfx(i,j)=on_p(i,j)*on_p(i,j)*visc4_p(i,j)*tl_cff
1979#endif
1980 END DO
1981 END DO
1982!
1983! Compute vertical flux (m2/s2) due to sloping terrain-following
1984! surfaces.
1985!
1986 IF (k.lt.n(ng)) THEN
1987 DO j=jstr,jend
1988 DO i=istru,iend
1989#ifdef VISC_3DCOEF
1990# ifdef UV_U3ADV_SPLIT
1991 cff=0.125_r8* &
1992 & (uvis3d_r(i-1,j,k )+uvis3d_r(i,j,k )+ &
1993 & uvis3d_r(i-1,j,k+1)+uvis3d_r(i,j,k+1))
1994 tl_cff=0.125_r8* &
1995 & (tl_uvis3d_r(i-1,j,k )+tl_uvis3d_r(i,j,k )+ &
1996 & tl_uvis3d_r(i-1,j,k+1)+tl_uvis3d_r(i,j,k+1))
1997# else
1998 cff=0.125_r8* &
1999 & (visc3d_r(i-1,j,k )+visc3d_r(i,j,k )+ &
2000 & visc3d_r(i-1,j,k+1)+visc3d_r(i,j,k+1))
2001 tl_cff=0.125_r8* &
2002 & (tl_visc3d_r(i-1,j,k )+tl_visc3d_r(i,j,k )+ &
2003 & tl_visc3d_r(i-1,j,k+1)+tl_visc3d_r(i,j,k+1))
2004# endif
2005 fac1=cff*on_u(i,j)
2006 tl_fac1=tl_cff*on_u(i,j)
2007 fac2=cff*om_u(i,j)
2008 tl_fac2=tl_cff*om_u(i,j)
2009#else
2010 cff=0.25_r8*(visc4_r(i-1,j)+visc4_r(i,j))
2011 fac1=cff*on_u(i,j)
2012 fac2=cff*om_u(i,j)
2013#endif
2014 cff=0.5_r8*(pn(i-1,j)+pn(i,j))
2015 dnudz=cff*dudz(i,j,k2)
2016 tl_dnudz=cff*tl_dudz(i,j,k2)
2017 dnvdz=cff*0.25_r8*(dvdz(i-1,j+1,k2)+ &
2018 & dvdz(i ,j+1,k2)+ &
2019 & dvdz(i-1,j ,k2)+ &
2020 & dvdz(i ,j ,k2))
2021 tl_dnvdz=cff*0.25_r8*(tl_dvdz(i-1,j+1,k2)+ &
2022 & tl_dvdz(i ,j+1,k2)+ &
2023 & tl_dvdz(i-1,j ,k2)+ &
2024 & tl_dvdz(i ,j ,k2))
2025 cff=0.5_r8*(pm(i-1,j)+pm(i,j))
2026 dmudz=cff*dudz(i,j,k2)
2027 tl_dmudz=cff*tl_dudz(i,j,k2)
2028 dmvdz=cff*0.25_r8*(dvdz(i-1,j+1,k2)+ &
2029 & dvdz(i ,j+1,k2)+ &
2030 & dvdz(i-1,j ,k2)+ &
2031 & dvdz(i ,j ,k2))
2032 tl_dmvdz=cff*0.25_r8*(tl_dvdz(i-1,j+1,k2)+ &
2033 & tl_dvdz(i ,j+1,k2)+ &
2034 & tl_dvdz(i-1,j ,k2)+ &
2035 & tl_dvdz(i ,j ,k2))
2036
2037 cff1=min(dzdx_r(i-1,j,k1),0.0_r8)
2038 cff2=min(dzdx_r(i ,j,k2),0.0_r8)
2039 cff3=max(dzdx_r(i-1,j,k2),0.0_r8)
2040 cff4=max(dzdx_r(i ,j,k1),0.0_r8)
2041 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_r(i-1,j,k1)))* &
2042 & tl_dzdx_r(i-1,j,k1)
2043 tl_cff2=(0.5_r8+sign(0.5_r8,-dzdx_r(i ,j,k2)))* &
2044 & tl_dzdx_r(i ,j,k2)
2045 tl_cff3=(0.5_r8+sign(0.5_r8, dzdx_r(i-1,j,k2)))* &
2046 & tl_dzdx_r(i-1,j,k2)
2047 tl_cff4=(0.5_r8+sign(0.5_r8, dzdx_r(i ,j,k1)))* &
2048 tl_dzdx_r(i ,j,k1)
2049!^ UFsx(i,j,k2)=fac1* &
2050!^ & (cff1*(cff1*dnUdz-dnUdx(i-1,j,k1))+ &
2051!^ & cff2*(cff2*dnUdz-dnUdx(i ,j,k2))+ &
2052!^ & cff3*(cff3*dnUdz-dnUdx(i-1,j,k2))+ &
2053!^ & cff4*(cff4*dnUdz-dnUdx(i ,j,k1)))
2054!^
2055 tl_ufsx(i,j,k2)=fac1* &
2056 & (tl_cff1*(cff1*dnudz-dnudx(i-1,j,k1))+ &
2057 & tl_cff2*(cff2*dnudz-dnudx(i ,j,k2))+ &
2058 & tl_cff3*(cff3*dnudz-dnudx(i-1,j,k2))+ &
2059 & tl_cff4*(cff4*dnudz-dnudx(i ,j,k1))+ &
2060 & cff1*(tl_cff1*dnudz+cff1*tl_dnudz- &
2061 & tl_dnudx(i-1,j,k1))+ &
2062 & cff2*(tl_cff2*dnudz+cff2*tl_dnudz- &
2063 & tl_dnudx(i ,j,k2))+ &
2064 & cff3*(tl_cff3*dnudz+cff3*tl_dnudz- &
2065 & tl_dnudx(i-1,j,k2))+ &
2066 & cff4*(tl_cff4*dnudz+cff4*tl_dnudz- &
2067 & tl_dnudx(i ,j,k1)))
2068#ifdef VISC_3DCOEF
2069 tl_ufsx(i,j,k2)=tl_ufsx(i,j,k2)+ &
2070 & tl_fac1* &
2071 & (cff1*(cff1*dnudz-dnudx(i-1,j,k1))+ &
2072 & cff2*(cff2*dnudz-dnudx(i ,j,k2))+ &
2073 & cff3*(cff3*dnudz-dnudx(i-1,j,k2))+ &
2074 & cff4*(cff4*dnudz-dnudx(i ,j,k1)))
2075#endif
2076
2077 cff1=min(dzde_p(i,j ,k1),0.0_r8)
2078 cff2=min(dzde_p(i,j+1,k2),0.0_r8)
2079 cff3=max(dzde_p(i,j ,k2),0.0_r8)
2080 cff4=max(dzde_p(i,j+1,k1),0.0_r8)
2081 tl_cff1=(0.5_r8+sign(0.5_r8,-dzde_p(i,j ,k1)))* &
2082 & tl_dzde_p(i,j ,k1)
2083 tl_cff2=(0.5_r8+sign(0.5_r8,-dzde_p(i,j+1,k2)))* &
2084 & tl_dzde_p(i,j+1,k2)
2085 tl_cff3=(0.5_r8+sign(0.5_r8, dzde_p(i,j ,k2)))* &
2086 & tl_dzde_p(i,j ,k2)
2087 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_p(i,j+1,k1)))* &
2088 tl_dzde_p(i,j+1,k1)
2089!^ UFse(i,j,k2)=fac2* &
2090!^ & (cff1*(cff1*dmUdz-dmUde(i,j ,k1))+ &
2091!^ & cff2*(cff2*dmUdz-dmUde(i,j+1,k2))+ &
2092!^ & cff3*(cff3*dmUdz-dmUde(i,j ,k2))+ &
2093!^ & cff4*(cff4*dmUdz-dmUde(i,j+1,k1)))
2094!^
2095 tl_ufse(i,j,k2)=fac2* &
2096 & (tl_cff1*(cff1*dmudz-dmude(i,j ,k1))+ &
2097 & tl_cff2*(cff2*dmudz-dmude(i,j+1,k2))+ &
2098 & tl_cff3*(cff3*dmudz-dmude(i,j ,k2))+ &
2099 & tl_cff4*(cff4*dmudz-dmude(i,j+1,k1))+ &
2100 & cff1*(tl_cff1*dmudz+cff1*tl_dmudz- &
2101 & tl_dmude(i,j ,k1))+ &
2102 & cff2*(tl_cff2*dmudz+cff2*tl_dmudz- &
2103 & tl_dmude(i,j+1,k2))+ &
2104 & cff3*(tl_cff3*dmudz+cff3*tl_dmudz- &
2105 & tl_dmude(i,j ,k2))+ &
2106 & cff4*(tl_cff4*dmudz+cff4*tl_dmudz- &
2107 & tl_dmude(i,j+1,k1)))
2108#ifdef VISC_3DCOEF
2109 tl_ufse(i,j,k2)=tl_ufse(i,j,k2)+ &
2110 & tl_fac2* &
2111 & (cff1*(cff1*dmudz-dmude(i,j ,k1))+ &
2112 & cff2*(cff2*dmudz-dmude(i,j+1,k2))+ &
2113 & cff3*(cff3*dmudz-dmude(i,j ,k2))+ &
2114 & cff4*(cff4*dmudz-dmude(i,j+1,k1)))
2115#endif
2116
2117 cff1=min(dzde_p(i,j ,k1),0.0_r8)
2118 cff2=min(dzde_p(i,j+1,k2),0.0_r8)
2119 cff3=max(dzde_p(i,j ,k2),0.0_r8)
2120 cff4=max(dzde_p(i,j+1,k1),0.0_r8)
2121 cff5=min(dzdx_p(i,j ,k1),0.0_r8)
2122 cff6=min(dzdx_p(i,j+1,k2),0.0_r8)
2123 cff7=max(dzdx_p(i,j ,k2),0.0_r8)
2124 cff8=max(dzdx_p(i,j+1,k1),0.0_r8)
2125 tl_cff1=(0.5_r8+sign(0.5_r8,-dzde_p(i,j ,k1)))* &
2126 & tl_dzde_p(i,j ,k1)
2127 tl_cff2=(0.5_r8+sign(0.5_r8,-dzde_p(i,j+1,k2)))* &
2128 & tl_dzde_p(i,j+1,k2)
2129 tl_cff3=(0.5_r8+sign(0.5_r8, dzde_p(i,j ,k2)))* &
2130 & tl_dzde_p(i,j ,k2)
2131 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_p(i,j+1,k1)))* &
2132 & tl_dzde_p(i,j+1,k1)
2133 tl_cff5=(0.5_r8+sign(0.5_r8,-dzdx_p(i,j ,k1)))* &
2134 & tl_dzdx_p(i,j ,k1)
2135 tl_cff6=(0.5_r8+sign(0.5_r8,-dzdx_p(i,j+1,k2)))* &
2136 & tl_dzdx_p(i,j+1,k2)
2137 tl_cff7=(0.5_r8+sign(0.5_r8, dzdx_p(i,j ,k2)))* &
2138 & tl_dzdx_p(i,j ,k2)
2139 tl_cff8=(0.5_r8+sign(0.5_r8, dzdx_p(i,j+1,k1)))* &
2140 & tl_dzdx_p(i,j+1,k1)
2141!^ UFsx(i,j,k2)=UFsx(i,j,k2)+ &
2142!^ & fac1* &
2143!^ & (cff1*(cff5*dnVdz-dnVdx(i,j ,k1))+ &
2144!^ & cff2*(cff6*dnVdz-dnVdx(i,j+1,k2))+ &
2145!^ & cff3*(cff7*dnVdz-dnVdx(i,j ,k2))+ &
2146!^ & cff4*(cff8*dnVdz-dnVdx(i,j+1,k1)))
2147!^
2148 tl_ufsx(i,j,k2)=tl_ufsx(i,j,k2)+ &
2149 & fac1* &
2150 & (tl_cff1*(cff5*dnvdz-dnvdx(i,j ,k1))+ &
2151 & tl_cff2*(cff6*dnvdz-dnvdx(i,j+1,k2))+ &
2152 & tl_cff3*(cff7*dnvdz-dnvdx(i,j ,k2))+ &
2153 & tl_cff4*(cff8*dnvdz-dnvdx(i,j+1,k1))+ &
2154 & cff1*(tl_cff5*dnvdz+cff5*tl_dnvdz- &
2155 & tl_dnvdx(i,j ,k1))+ &
2156 & cff2*(tl_cff6*dnvdz+cff6*tl_dnvdz- &
2157 & tl_dnvdx(i,j+1,k2))+ &
2158 & cff3*(tl_cff7*dnvdz+cff7*tl_dnvdz- &
2159 & tl_dnvdx(i,j ,k2))+ &
2160 & cff4*(tl_cff8*dnvdz+cff8*tl_dnvdz- &
2161 & tl_dnvdx(i,j+1,k1)))
2162#ifdef VISC_3DCOEF
2163 tl_ufsx(i,j,k2)=tl_ufsx(i,j,k2)+ &
2164 & tl_fac1* &
2165 & (cff1*(cff5*dnvdz-dnvdx(i,j ,k1))+ &
2166 & cff2*(cff6*dnvdz-dnvdx(i,j+1,k2))+ &
2167 & cff3*(cff7*dnvdz-dnvdx(i,j ,k2))+ &
2168 & cff4*(cff8*dnvdz-dnvdx(i,j+1,k1)))
2169#endif
2170
2171 cff1=min(dzdx_r(i-1,j,k1),0.0_r8)
2172 cff2=min(dzdx_r(i ,j,k2),0.0_r8)
2173 cff3=max(dzdx_r(i-1,j,k2),0.0_r8)
2174 cff4=max(dzdx_r(i ,j,k1),0.0_r8)
2175 cff5=min(dzde_r(i-1,j,k1),0.0_r8)
2176 cff6=min(dzde_r(i ,j,k2),0.0_r8)
2177 cff7=max(dzde_r(i-1,j,k2),0.0_r8)
2178 cff8=max(dzde_r(i ,j,k1),0.0_r8)
2179 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_r(i-1,j,k1)))* &
2180 & tl_dzdx_r(i-1,j,k1)
2181 tl_cff2=(0.5_r8+sign(0.5_r8,-dzdx_r(i ,j,k2)))* &
2182 & tl_dzdx_r(i ,j,k2)
2183 tl_cff3=(0.5_r8+sign(0.5_r8, dzdx_r(i-1,j,k2)))* &
2184 & tl_dzdx_r(i-1,j,k2)
2185 tl_cff4=(0.5_r8+sign(0.5_r8, dzdx_r(i ,j,k1)))* &
2186 & tl_dzdx_r(i ,j,k1)
2187 tl_cff5=(0.5_r8+sign(0.5_r8,-dzde_r(i-1,j,k1)))* &
2188 & tl_dzde_r(i-1,j,k1)
2189 tl_cff6=(0.5_r8+sign(0.5_r8,-dzde_r(i ,j,k2)))* &
2190 & tl_dzde_r(i ,j,k2)
2191 tl_cff7=(0.5_r8+sign(0.5_r8, dzde_r(i-1,j,k2)))* &
2192 & tl_dzde_r(i-1,j,k2)
2193 tl_cff8=(0.5_r8+sign(0.5_r8, dzde_r(i ,j,k1)))* &
2194 & tl_dzde_r(i ,j,k1)
2195!^ UFse(i,j,k2)=UFse(i,j,k2)- &
2196!^ & fac2* &
2197!^ & (cff1*(cff5*dmVdz-dmVde(i-1,j,k1))+ &
2198!^ & cff2*(cff6*dmVdz-dmVde(i ,j,k2))+ &
2199!^ & cff3*(cff7*dmVdz-dmVde(i-1,j,k2))+ &
2200!^ & cff4*(cff8*dmVdz-dmVde(i ,j,k1)))
2201!^
2202 tl_ufse(i,j,k2)=tl_ufse(i,j,k2)- &
2203 & fac2* &
2204 & (tl_cff1*(cff5*dmvdz-dmvde(i-1,j,k1))+ &
2205 & tl_cff2*(cff6*dmvdz-dmvde(i ,j,k2))+ &
2206 & tl_cff3*(cff7*dmvdz-dmvde(i-1,j,k2))+ &
2207 & tl_cff4*(cff8*dmvdz-dmvde(i ,j,k1))+ &
2208 & cff1*(tl_cff5*dmvdz+cff5*tl_dmvdz- &
2209 & tl_dmvde(i-1,j,k1))+ &
2210 & cff2*(tl_cff6*dmvdz+cff6*tl_dmvdz- &
2211 & tl_dmvde(i ,j,k2))+ &
2212 & cff3*(tl_cff7*dmvdz+cff7*tl_dmvdz- &
2213 & tl_dmvde(i-1,j,k2))+ &
2214 & cff4*(tl_cff8*dmvdz+cff8*tl_dmvdz- &
2215 & tl_dmvde(i ,j,k1)))
2216#ifdef VISC_3DCOEF
2217 tl_ufse(i,j,k2)=tl_ufse(i,j,k2)- &
2218 & tl_fac2* &
2219 & (cff1*(cff5*dmvdz-dmvde(i-1,j,k1))+ &
2220 & cff2*(cff6*dmvdz-dmvde(i ,j,k2))+ &
2221 & cff3*(cff7*dmvdz-dmvde(i-1,j,k2))+ &
2222 & cff4*(cff8*dmvdz-dmvde(i ,j,k1)))
2223#endif
2224 END DO
2225 END DO
2226!
2227 DO j=jstrv,jend
2228 DO i=istr,iend
2229#ifdef VISC_3DCOEF
2230# ifdef UV_U3ADV_SPLIT
2231 cff=0.125_r8* &
2232 & (vvis3d_r(i,j-1,k )+vvis3d_r(i,j,k )+ &
2233 & vvis3d_r(i,j-1,k+1)+vvis3d_r(i,j,k+1))
2234 tl_cff=0.125_r8* &
2235 & (tl_vvis3d_r(i,j-1,k )+tl_vvis3d_r(i,j,k )+ &
2236 & tl_vvis3d_r(i,j-1,k+1)+tl_vvis3d_r(i,j,k+1))
2237# else
2238 cff=0.125_r8* &
2239 & (visc3d_r(i,j-1,k )+visc3d_r(i,j,k )+ &
2240 & visc3d_r(i,j-1,k+1)+visc3d_r(i,j,k+1))
2241 tl_cff=0.125_r8* &
2242 & (tl_visc3d_r(i,j-1,k )+tl_visc3d_r(i,j,k )+ &
2243 & tl_visc3d_r(i,j-1,k+1)+tl_visc3d_r(i,j,k+1))
2244# endif
2245 fac1=cff*on_v(i,j)
2246 tl_fac1=tl_cff*on_v(i,j)
2247 fac2=cff*om_v(i,j)
2248 tl_fac2=tl_cff*om_v(i,j)
2249#else
2250 cff=0.25_r8*(visc4_r(i,j-1)+visc4_r(i,j))
2251 fac1=cff*on_v(i,j)
2252 fac2=cff*om_v(i,j)
2253#endif
2254 cff=0.5_r8*(pn(i,j-1)+pn(i,j))
2255 dnudz=cff*0.25_r8*(dudz(i ,j ,k2)+ &
2256 & dudz(i+1,j ,k2)+ &
2257 & dudz(i ,j-1,k2)+ &
2258 & dudz(i+1,j-1,k2))
2259 tl_dnudz=cff*0.25_r8*(tl_dudz(i ,j ,k2)+ &
2260 & tl_dudz(i+1,j ,k2)+ &
2261 & tl_dudz(i ,j-1,k2)+ &
2262 & tl_dudz(i+1,j-1,k2))
2263 dnvdz=cff*dvdz(i,j,k2)
2264 tl_dnvdz=cff*tl_dvdz(i,j,k2)
2265 cff=0.5_r8*(pm(i,j-1)+pm(i,j))
2266 dmudz=cff*0.25_r8*(dudz(i ,j ,k2)+ &
2267 & dudz(i+1,j ,k2)+ &
2268 & dudz(i ,j-1,k2)+ &
2269 & dudz(i+1,j-1,k2))
2270 tl_dmudz=cff*0.25_r8*(tl_dudz(i ,j ,k2)+ &
2271 & tl_dudz(i+1,j ,k2)+ &
2272 & tl_dudz(i ,j-1,k2)+ &
2273 & tl_dudz(i+1,j-1,k2))
2274 dmvdz=cff*dvdz(i,j,k2)
2275 tl_dmvdz=cff*tl_dvdz(i,j,k2)
2276
2277 cff1=min(dzdx_p(i ,j,k1),0.0_r8)
2278 cff2=min(dzdx_p(i+1,j,k2),0.0_r8)
2279 cff3=max(dzdx_p(i ,j,k2),0.0_r8)
2280 cff4=max(dzdx_p(i+1,j,k1),0.0_r8)
2281 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_p(i ,j,k1)))* &
2282 & tl_dzdx_p(i ,j,k1)
2283 tl_cff2=(0.5_r8+sign(0.5_r8,-dzdx_p(i+1,j,k2)))* &
2284 & tl_dzdx_p(i+1,j,k2)
2285 tl_cff3=(0.5_r8+sign(0.5_r8, dzdx_p(i ,j,k2)))* &
2286 & tl_dzdx_p(i ,j,k2)
2287 tl_cff4=(0.5_r8+sign(0.5_r8, dzdx_p(i+1,j,k1)))* &
2288 & tl_dzdx_p(i+1,j,k1)
2289!^ VFsx(i,j,k2)=fac1* &
2290!^ & (cff1*(cff1*dnVdz-dnVdx(i ,j,k1))+ &
2291!^ & cff2*(cff2*dnVdz-dnVdx(i+1,j,k2))+ &
2292!^ & cff3*(cff3*dnVdz-dnVdx(i ,j,k2))+ &
2293!^ & cff4*(cff4*dnVdz-dnVdx(i+1,j,k1)))
2294!^
2295 tl_vfsx(i,j,k2)=fac1* &
2296 & (tl_cff1*(cff1*dnvdz-dnvdx(i ,j,k1))+ &
2297 & tl_cff2*(cff2*dnvdz-dnvdx(i+1,j,k2))+ &
2298 & tl_cff3*(cff3*dnvdz-dnvdx(i ,j,k2))+ &
2299 & tl_cff4*(cff4*dnvdz-dnvdx(i+1,j,k1))+ &
2300 & cff1*(tl_cff1*dnvdz+cff1*tl_dnvdz- &
2301 & tl_dnvdx(i ,j,k1))+ &
2302 & cff2*(tl_cff2*dnvdz+cff2*tl_dnvdz- &
2303 & tl_dnvdx(i+1,j,k2))+ &
2304 & cff3*(tl_cff3*dnvdz+cff3*tl_dnvdz- &
2305 & tl_dnvdx(i ,j,k2))+ &
2306 & cff4*(tl_cff4*dnvdz+cff4*tl_dnvdz- &
2307 & tl_dnvdx(i+1,j,k1)))
2308#ifdef VISC_3DCOEF
2309 tl_vfsx(i,j,k2)=tl_vfsx(i,j,k2)+ &
2310 & tl_fac1* &
2311 & (cff1*(cff1*dnvdz-dnvdx(i ,j,k1))+ &
2312 & cff2*(cff2*dnvdz-dnvdx(i+1,j,k2))+ &
2313 & cff3*(cff3*dnvdz-dnvdx(i ,j,k2))+ &
2314 & cff4*(cff4*dnvdz-dnvdx(i+1,j,k1)))
2315#endif
2316
2317 cff1=min(dzde_r(i,j-1,k1),0.0_r8)
2318 cff2=min(dzde_r(i,j ,k2),0.0_r8)
2319 cff3=max(dzde_r(i,j-1,k2),0.0_r8)
2320 cff4=max(dzde_r(i,j ,k1),0.0_r8)
2321 tl_cff1=(0.5_r8+sign(0.5_r8,-dzde_r(i,j-1,k1)))* &
2322 & tl_dzde_r(i,j-1,k1)
2323 tl_cff2=(0.5_r8+sign(0.5_r8,-dzde_r(i,j ,k2)))* &
2324 & tl_dzde_r(i,j ,k2)
2325 tl_cff3=(0.5_r8+sign(0.5_r8, dzde_r(i,j-1,k2)))* &
2326 & tl_dzde_r(i,j-1,k2)
2327 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_r(i,j ,k1)))* &
2328 & tl_dzde_r(i,j ,k1)
2329!^ VFse(i,j,k2)=fac2* &
2330!^ & (cff1*(cff1*dmVdz-dmVde(i,j-1,k1))+ &
2331!^ & cff2*(cff2*dmVdz-dmVde(i,j ,k2))+ &
2332!^ & cff3*(cff3*dmVdz-dmVde(i,j-1,k2))+ &
2333!^ & cff4*(cff4*dmVdz-dmVde(i,j ,k1)))
2334!^
2335 tl_vfse(i,j,k2)=fac2* &
2336 & (tl_cff1*(cff1*dmvdz-dmvde(i,j-1,k1))+ &
2337 & tl_cff2*(cff2*dmvdz-dmvde(i,j ,k2))+ &
2338 & tl_cff3*(cff3*dmvdz-dmvde(i,j-1,k2))+ &
2339 & tl_cff4*(cff4*dmvdz-dmvde(i,j ,k1))+ &
2340 & cff1*(tl_cff1*dmvdz+cff1*tl_dmvdz- &
2341 & tl_dmvde(i,j-1,k1))+ &
2342 & cff2*(tl_cff2*dmvdz+cff2*tl_dmvdz- &
2343 & tl_dmvde(i,j ,k2))+ &
2344 & cff3*(tl_cff3*dmvdz+cff3*tl_dmvdz- &
2345 & tl_dmvde(i,j-1,k2))+ &
2346 & cff4*(tl_cff4*dmvdz+cff4*tl_dmvdz- &
2347 & tl_dmvde(i,j ,k1)))
2348#ifdef VISC_3DCOEF
2349 tl_vfse(i,j,k2)=tl_vfse(i,j,k2)+ &
2350 & tl_fac2* &
2351 & (cff1*(cff1*dmvdz-dmvde(i,j-1,k1))+ &
2352 & cff2*(cff2*dmvdz-dmvde(i,j ,k2))+ &
2353 & cff3*(cff3*dmvdz-dmvde(i,j-1,k2))+ &
2354 & cff4*(cff4*dmvdz-dmvde(i,j ,k1)))
2355#endif
2356
2357 cff1=min(dzde_r(i,j-1,k1),0.0_r8)
2358 cff2=min(dzde_r(i,j ,k2),0.0_r8)
2359 cff3=max(dzde_r(i,j-1,k2),0.0_r8)
2360 cff4=max(dzde_r(i,j ,k1),0.0_r8)
2361 cff5=min(dzdx_r(i,j-1,k1),0.0_r8)
2362 cff6=min(dzdx_r(i,j ,k2),0.0_r8)
2363 cff7=max(dzdx_r(i,j-1,k2),0.0_r8)
2364 cff8=max(dzdx_r(i,j ,k1),0.0_r8)
2365 tl_cff1=(0.5_r8+sign(0.5_r8,-dzde_r(i,j-1,k1)))* &
2366 & tl_dzde_r(i,j-1,k1)
2367 tl_cff2=(0.5_r8+sign(0.5_r8,-dzde_r(i,j ,k2)))* &
2368 & tl_dzde_r(i,j ,k2)
2369 tl_cff3=(0.5_r8+sign(0.5_r8, dzde_r(i,j-1,k2)))* &
2370 & tl_dzde_r(i,j-1,k2)
2371 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_r(i,j ,k1)))* &
2372 & tl_dzde_r(i,j ,k1)
2373 tl_cff5=(0.5_r8+sign(0.5_r8,-dzdx_r(i,j-1,k1)))* &
2374 & tl_dzdx_r(i,j-1,k1)
2375 tl_cff6=(0.5_r8+sign(0.5_r8,-dzdx_r(i,j ,k2)))* &
2376 & tl_dzdx_r(i,j ,k2)
2377 tl_cff7=(0.5_r8+sign(0.5_r8, dzdx_r(i,j-1,k2)))* &
2378 & tl_dzdx_r(i,j-1,k2)
2379 tl_cff8=(0.5_r8+sign(0.5_r8, dzdx_r(i,j ,k1)))* &
2380 & tl_dzdx_r(i,j ,k1)
2381!^ VFsx(i,j,k2)=VFsx(i,j,k2)- &
2382!^ & fac1* &
2383!^ & (cff1*(cff5*dnUdz-dnUdx(i,j-1,k1))+ &
2384!^ & cff2*(cff6*dnUdz-dnUdx(i,j ,k2))+ &
2385!^ & cff3*(cff7*dnUdz-dnUdx(i,j-1,k2))+ &
2386!^ & cff4*(cff8*dnUdz-dnUdx(i,j ,k1)))
2387!^
2388 tl_vfsx(i,j,k2)=tl_vfsx(i,j,k2)- &
2389 & fac1* &
2390 & (tl_cff1*(cff5*dnudz-dnudx(i,j-1,k1))+ &
2391 & tl_cff2*(cff6*dnudz-dnudx(i,j ,k2))+ &
2392 & tl_cff3*(cff7*dnudz-dnudx(i,j-1,k2))+ &
2393 & tl_cff4*(cff8*dnudz-dnudx(i,j ,k1))+ &
2394 & cff1*(tl_cff5*dnudz+cff5*tl_dnudz- &
2395 & tl_dnudx(i,j-1,k1))+ &
2396 & cff2*(tl_cff6*dnudz+cff6*tl_dnudz- &
2397 & tl_dnudx(i,j ,k2))+ &
2398 & cff3*(tl_cff7*dnudz+cff7*tl_dnudz- &
2399 & tl_dnudx(i,j-1,k2))+ &
2400 & cff4*(tl_cff8*dnudz+cff8*tl_dnudz- &
2401 & tl_dnudx(i,j ,k1)))
2402#ifdef VISC_3DCOEF
2403 tl_vfsx(i,j,k2)=tl_vfsx(i,j,k2)- &
2404 & tl_fac1* &
2405 & (cff1*(cff5*dnudz-dnudx(i,j-1,k1))+ &
2406 & cff2*(cff6*dnudz-dnudx(i,j ,k2))+ &
2407 & cff3*(cff7*dnudz-dnudx(i,j-1,k2))+ &
2408 & cff4*(cff8*dnudz-dnudx(i,j ,k1)))
2409#endif
2410
2411 cff1=min(dzdx_p(i ,j,k1),0.0_r8)
2412 cff2=min(dzdx_p(i+1,j,k2),0.0_r8)
2413 cff3=max(dzdx_p(i ,j,k2),0.0_r8)
2414 cff4=max(dzdx_p(i+1,j,k1),0.0_r8)
2415 cff5=min(dzde_p(i ,j,k1),0.0_r8)
2416 cff6=min(dzde_p(i+1,j,k2),0.0_r8)
2417 cff7=max(dzde_p(i ,j,k2),0.0_r8)
2418 cff8=max(dzde_p(i+1,j,k1),0.0_r8)
2419 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_p(i ,j,k1)))* &
2420 & tl_dzdx_p(i ,j,k1)
2421 tl_cff2=(0.5_r8+sign(0.5_r8,-dzdx_p(i+1,j,k2)))* &
2422 & tl_dzdx_p(i+1,j,k2)
2423 tl_cff3=(0.5_r8+sign(0.5_r8, dzdx_p(i ,j,k2)))* &
2424 & tl_dzdx_p(i ,j,k2)
2425 tl_cff4=(0.5_r8+sign(0.5_r8, dzdx_p(i+1,j,k1)))* &
2426 & tl_dzdx_p(i+1,j,k1)
2427 tl_cff5=(0.5_r8+sign(0.5_r8,-dzde_p(i ,j,k1)))* &
2428 & tl_dzde_p(i ,j,k1)
2429 tl_cff6=(0.5_r8+sign(0.5_r8,-dzde_p(i+1,j,k2)))* &
2430 & tl_dzde_p(i+1,j,k2)
2431 tl_cff7=(0.5_r8+sign(0.5_r8, dzde_p(i ,j,k2)))* &
2432 & tl_dzde_p(i ,j,k2)
2433 tl_cff8=(0.5_r8+sign(0.5_r8, dzde_p(i+1,j,k1)))* &
2434 & tl_dzde_p(i+1,j,k1)
2435!^ VFse(i,j,k2)=VFse(i,j,k2)+ &
2436!^ & fac2* &
2437!^ & (cff1*(cff5*dmUdz-dmUde(i ,j,k1))+ &
2438!^ & cff2*(cff6*dmUdz-dmUde(i+1,j,k2))+ &
2439!^ & cff3*(cff7*dmUdz-dmUde(i ,j,k2))+ &
2440!^ & cff4*(cff8*dmUdz-dmUde(i+1,j,k1)))
2441!^
2442 tl_vfse(i,j,k2)=tl_vfse(i,j,k2)+ &
2443 & fac2* &
2444 & (tl_cff1*(cff5*dmudz-dmude(i ,j,k1))+ &
2445 & tl_cff2*(cff6*dmudz-dmude(i+1,j,k2))+ &
2446 & tl_cff3*(cff7*dmudz-dmude(i ,j,k2))+ &
2447 & tl_cff4*(cff8*dmudz-dmude(i+1,j,k1))+ &
2448 & cff1*(tl_cff5*dmudz+cff5*tl_dmudz- &
2449 & tl_dmude(i ,j,k1))+ &
2450 & cff2*(tl_cff6*dmudz+cff6*tl_dmudz- &
2451 & tl_dmude(i+1,j,k2))+ &
2452 & cff3*(tl_cff7*dmudz+cff7*tl_dmudz- &
2453 & tl_dmude(i ,j,k2))+ &
2454 & cff4*(tl_cff8*dmudz+cff8*tl_dmudz- &
2455 & tl_dmude(i+1,j,k1)))
2456#ifdef VISC_3DCOEF
2457 tl_vfse(i,j,k2)=tl_vfse(i,j,k2)+ &
2458 & tl_fac2* &
2459 & (cff1*(cff5*dmudz-dmude(i ,j,k1))+ &
2460 & cff2*(cff6*dmudz-dmude(i+1,j,k2))+ &
2461 & cff3*(cff7*dmudz-dmude(i ,j,k2))+ &
2462 & cff4*(cff8*dmudz-dmude(i+1,j,k1)))
2463#endif
2464 END DO
2465 END DO
2466 END IF
2467!
2468! Time-step biharmonic, geopotential viscosity term. Notice that
2469! momentum at this stage is HzU and HzV and has m2/s units. Add
2470! contribution for barotropic forcing terms.
2471#ifdef DIAGNOSTICS_UV
2472!! The rotated vertical term cannot be split from the horizontal
2473!! terms because of the 2D/3D momentum coupling.
2474#endif
2475!
2476 DO j=jstr,jend
2477 DO i=istru,iend
2478 cff=dt(ng)*0.25_r8*(pm(i-1,j)+pm(i,j))*(pn(i-1,j)+pn(i,j))
2479!^ cff1=0.5_r8*(pn(i-1,j)+pn(i,j))*(UFx(i,j )-UFx(i-1,j))
2480!^
2481 tl_cff1=0.5_r8*(pn(i-1,j)+pn(i,j))* &
2482 & (tl_ufx(i,j )-tl_ufx(i-1,j))
2483!^ cff2=0.5_r8*(pm(i-1,j)+pm(i,j))*(UFe(i,j+1)-UFe(i ,j))
2484!^
2485 tl_cff2=0.5_r8*(pm(i-1,j)+pm(i,j))* &
2486 & (tl_ufe(i,j+1)-tl_ufe(i ,j))
2487!^ cff3=UFsx(i,j,k2)-UFsx(i,j,k1)
2488!^
2489 tl_cff3=tl_ufsx(i,j,k2)-tl_ufsx(i,j,k1)
2490!^ cff4=UFse(i,j,k2)-UFse(i,j,k1)
2491!^
2492 tl_cff4=tl_ufse(i,j,k2)-tl_ufse(i,j,k1)
2493!^ cff5=cff*(cff1+cff2)
2494!^
2495 tl_cff5=cff*(tl_cff1+tl_cff2)
2496!^ cff6=dt(ng)*(cff3+cff4)
2497!^
2498 tl_cff6=dt(ng)*(tl_cff3+tl_cff4)
2499!^ rufrc(i,j)=rufrc(i,j)-cff1-cff2-cff3-cff4
2500!^
2501 tl_rufrc(i,j)=tl_rufrc(i,j)- &
2502 & tl_cff1-tl_cff2-tl_cff3-tl_cff4
2503!^ u(i,j,k,nnew)=u(i,j,k,nnew)-cff5-cff6
2504!^
2505 tl_u(i,j,k,nnew)=tl_u(i,j,k,nnew)-tl_cff5-tl_cff6
2506#ifdef DIAGNOSTICS_UV
2507!! DiaRUfrc(i,j,3,M2hvis)=DiaRUfrc(i,j,3,M2hvis)-cff1-cff2- &
2508!! & cff3-cff4
2509!! DiaRUfrc(i,j,3,M2xvis)=DiaRUfrc(i,j,3,M2xvis)-cff1-cff3
2510!! DiaRUfrc(i,j,3,M2yvis)=DiaRUfrc(i,j,3,M2yvis)-cff2-cff4
2511!! DiaU3wrk(i,j,k,M3hvis)=-cff5-cff6
2512!! DiaU3wrk(i,j,k,M3xvis)=-cff*cff1-dt(ng)*cff3
2513!! DiaU3wrk(i,j,k,M3yvis)=-cff*cff2-dt(ng)*cff4
2514#endif
2515 END DO
2516 END DO
2517
2518 DO j=jstrv,jend
2519 DO i=istr,iend
2520 cff=dt(ng)*0.25_r8*(pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))
2521!^ cff1=0.5_r8*(pn(i,j-1)+pn(i,j))*(VFx(i+1,j)-VFx(i,j ))
2522!^
2523 tl_cff1=0.5_r8*(pn(i,j-1)+pn(i,j))* &
2524 & (tl_vfx(i+1,j)-tl_vfx(i,j ))
2525!^ cff2=0.5_r8*(pm(i,j-1)+pm(i,j))*(VFe(i ,j)-VFe(i,j-1))
2526!^
2527 tl_cff2=0.5_r8*(pm(i,j-1)+pm(i,j))* &
2528 & (tl_vfe(i ,j)-tl_vfe(i,j-1))
2529!^ cff3=VFsx(i,j,k2)-VFsx(i,j,k1)
2530!^
2531 tl_cff3=tl_vfsx(i,j,k2)-tl_vfsx(i,j,k1)
2532!^ cff4=VFse(i,j,k2)-VFse(i,j,k1)
2533!^
2534 tl_cff4=tl_vfse(i,j,k2)-tl_vfse(i,j,k1)
2535!^ cff5=cff*(cff1-cff2)
2536!^
2537 tl_cff5=cff*(tl_cff1-tl_cff2)
2538!^ cff6=dt(ng)*(cff3+cff4)
2539!^
2540 tl_cff6=dt(ng)*(tl_cff3+tl_cff4)
2541!^ rvfrc(i,j)=rvfrc(i,j)-cff1+cff2-cff3-cff4
2542!^
2543 tl_rvfrc(i,j)=tl_rvfrc(i,j)- &
2544 & tl_cff1+tl_cff2-tl_cff3-tl_cff4
2545!^ v(i,j,k,nnew)=v(i,j,k,nnew)-cff5-cff6
2546!^
2547 tl_v(i,j,k,nnew)=tl_v(i,j,k,nnew)-tl_cff5-tl_cff6
2548#ifdef DIAGNOSTICS_UV
2549!! DiaRVfrc(i,j,3,M2hvis)=DiaRVfrc(i,j,3,M2hvis)-cff1+cff2- &
2550!! & cff3-cff4
2551!! DiaRVfrc(i,j,3,M2xvis)=DiaRVfrc(i,j,3,M2xvis)-cff1-cff3
2552!! DiaRVfrc(i,j,3,M2yvis)=DiaRVfrc(i,j,3,M2yvis)+cff2-cff4
2553!! DiaV3wrk(i,j,k,M3hvis)=-cff5-cff6
2554!! DiaV3wrk(i,j,k,M3xvis)=-cff*cff1-dt(ng)*cff3
2555!! DiaV3wrk(i,j,k,M3yvis)= cff*cff2-dt(ng)*cff4
2556#endif
2557 END DO
2558 END DO
2559 END IF
2560 END DO k_loop2
2561!
2562 RETURN
2563 END SUBROUTINE tl_uv3dmix4_geo_tile
2564
2565 END MODULE tl_uv3dmix4_mod
mod_coupling
Definition mod_coupling.F:2
mod_coupling::coupling
type(t_coupling), dimension(:), allocatable coupling
Definition mod_coupling.F:102
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_ncparam
Definition mod_ncparam.F:2
mod_ncparam::isvvel
integer isvvel
Definition mod_ncparam.F:507
mod_ncparam::isuvel
integer isuvel
Definition mod_ncparam.F:506
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::tl_lbc
type(t_lbc), dimension(:,:,:), allocatable tl_lbc
Definition mod_param.F:379
mod_param::domain
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
mod_param::itlm
integer, parameter itlm
Definition mod_param.F:663
mod_scalars
Definition mod_scalars.F:2
mod_scalars::dt
real(dp), dimension(:), allocatable dt
Definition mod_scalars.F:269
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::gamma2
real(r8), dimension(:), allocatable gamma2
Definition mod_scalars.F:743
mod_scalars::compositegrid
logical, dimension(:,:), allocatable compositegrid
Definition mod_scalars.F:493
mod_scalars::isouth
integer, parameter isouth
Definition mod_scalars.F:1433
mod_scalars::ieast
integer, parameter ieast
Definition mod_scalars.F:1434
mod_scalars::inorth
integer, parameter inorth
Definition mod_scalars.F:1435
mod_stepping
Definition mod_stepping.F:2
mod_stepping::nrhs
integer, dimension(:), allocatable nrhs
Definition mod_stepping.F:70
mod_stepping::nnew
integer, dimension(:), allocatable nnew
Definition mod_stepping.F:69
tl_uv3dmix4_mod
Definition tl_uv3dmix4_geo.h:1
tl_uv3dmix4_mod::tl_uv3dmix4_geo_tile
subroutine tl_uv3dmix4_geo_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, nnew, pmask, rmask, umask, vmask, om_p, om_r, om_u, om_v, on_p, on_r, on_u, on_v, pm, pn, hz, tl_hz, z_r, tl_z_r, uvis3d_r, vvis3d_r, tl_uvis3d_r, tl_vvis3d_r, visc3d_r, tl_visc3d_r, visc4_p, visc4_r, u, v, tl_u, tl_v, tl_rufrc, tl_rvfrc)
Definition tl_uv3dmix4_geo.h:157
tl_uv3dmix4_mod::tl_uv3dmix4
subroutine, public tl_uv3dmix4(ng, tile)
Definition tl_uv3dmix4_geo.h:44
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