ROMS
Loading...
Searching...
No Matches
Functions/Subroutines
rp_uv3dmix4_mod Module Reference

Functions/Subroutines

subroutine, public rp_uv3dmix4 (ng, tile)
 
subroutine rp_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)
 
subroutine rp_uv3dmix4_s_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, nnew, pmask, hz, tl_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, tl_rufrc, tl_rvfrc, tl_u, tl_v)
 

Function/Subroutine Documentation

◆ rp_uv3dmix4()

subroutine public rp_uv3dmix4_mod::rp_uv3dmix4 ( integer, intent(in) ng,
integer, intent(in) tile )

Definition at line 43 of file rp_uv3dmix4_geo.h.

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, irpm, 33, __line__, myfile)
69#endif
70 CALL rp_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, irpm, 33, __line__, myfile)
122#endif
123!
124 RETURN
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_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::irpm
integer, parameter irpm
Definition mod_param.F:664
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
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

References mod_coupling::coupling, mod_grid::grid, mod_param::irpm, mod_mixing::mixing, mod_stepping::nnew, mod_stepping::nrhs, mod_ocean::ocean, rp_uv3dmix4_geo_tile(), wclock_off(), and wclock_on().

Referenced by rp_rhs3d_mod::rp_rhs3d().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ rp_uv3dmix4_geo_tile()

subroutine rp_uv3dmix4_mod::rp_uv3dmix4_geo_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) imins,
integer, intent(in) imaxs,
integer, intent(in) jmins,
integer, intent(in) jmaxs,
integer, intent(in) nrhs,
integer, intent(in) nnew,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) pmask,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) rmask,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) umask,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) vmask,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) om_p,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) om_r,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) om_u,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) om_v,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) on_p,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) on_r,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) on_u,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) on_v,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) pm,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) pn,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) hz,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) tl_hz,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) z_r,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) tl_z_r,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) uvis3d_r,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) vvis3d_r,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(inout) tl_uvis3d_r,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(inout) tl_vvis3d_r,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) visc3d_r,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(inout) tl_visc3d_r,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) visc4_p,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) visc4_r,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),2), intent(in) u,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),2), intent(in) v,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),2), intent(inout) tl_u,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),2), intent(inout) tl_v,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(inout) tl_rufrc,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(inout) tl_rvfrc )
private

Definition at line 128 of file rp_uv3dmix4_geo.h.

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#ifdef TL_IOMS
538 & 2.0_r8*cff
539#endif
540 dudz(i,j,k2)=cff*(u(i,j,k+1,nrhs)- &
541 & u(i,j,k ,nrhs))
542 tl_dudz(i,j,k2)=tl_cff*(u(i,j,k+1,nrhs)- &
543 & u(i,j,k ,nrhs))+ &
544 & cff*(tl_u(i,j,k+1,nrhs)- &
545 & tl_u(i,j,k ,nrhs))- &
546#ifdef TL_IOMS
547 & dudz(i,j,k2)
548#endif
549 END DO
550 END DO
551
552 DO j=jstrvm2,jendp2
553 DO i=istrm2,iendp2
554 cff=1.0_r8/(0.5_r8*(z_r(i,j-1,k+1)- &
555 & z_r(i,j-1,k )+ &
556 & z_r(i,j ,k+1)- &
557 & z_r(i,j ,k )))
558 tl_cff=-cff*cff*(0.5_r8*(tl_z_r(i,j-1,k+1)- &
559 & tl_z_r(i,j-1,k )+ &
560 & tl_z_r(i,j ,k+1)- &
561 & tl_z_r(i,j ,k )))+ &
562#ifdef TL_IOMS
563 & 2.0_r8*cff
564#endif
565 dvdz(i,j,k2)=cff*(v(i,j,k+1,nrhs)- &
566 & v(i,j,k ,nrhs))
567 tl_dvdz(i,j,k2)=tl_cff*(v(i,j,k+1,nrhs)- &
568 & v(i,j,k ,nrhs))+ &
569 & cff*(tl_v(i,j,k+1,nrhs)- &
570 & tl_v(i,j,k ,nrhs))- &
571#ifdef TL_IOMS
572 & dvdz(i,j,k2)
573#endif
574 END DO
575 END DO
576 END IF
577!
578! Compute components of the rotated viscous flux (m^4 s-^3/2) along
579! geopotential surfaces in the XI- and ETA-directions.
580!
581 IF (k.gt.0) THEN
582 DO j=jstrvm2,jendp1
583 DO i=istrum2,iendp1
584 cff1=min(dzdx_r(i,j,k1),0.0_r8)
585 cff2=max(dzdx_r(i,j,k1),0.0_r8)
586 cff3=min(dzde_r(i,j,k1),0.0_r8)
587 cff4=max(dzde_r(i,j,k1),0.0_r8)
588 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_r(i,j,k1)))* &
589 & tl_dzdx_r(i,j,k1)
590 tl_cff2=(0.5_r8+sign(0.5_r8, dzdx_r(i,j,k1)))* &
591 & tl_dzdx_r(i,j,k1)
592 tl_cff3=(0.5_r8+sign(0.5_r8,-dzde_r(i,j,k1)))* &
593 & tl_dzde_r(i,j,k1)
594 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_r(i,j,k1)))* &
595 & tl_dzde_r(i,j,k1)
596 cff=on_r(i,j)*(dnudx(i,j,k1)- &
597 & 0.5_r8*pn(i,j)* &
598 & (cff1*(dudz(i ,j,k1)+ &
599 & dudz(i+1,j,k2))+ &
600 & cff2*(dudz(i ,j,k2)+ &
601 & dudz(i+1,j,k1))))- &
602 & om_r(i,j)*(dmvde(i,j,k1)- &
603 & 0.5_r8*pm(i,j)* &
604 & (cff3*(dvdz(i,j ,k1)+ &
605 & dvdz(i,j+1,k2))+ &
606 & cff4*(dvdz(i,j ,k2)+ &
607 & dvdz(i,j+1,k1))))
608 tl_cff=on_r(i,j)*(tl_dnudx(i,j,k1)- &
609 & 0.5_r8*pn(i,j)* &
610 & (tl_cff1*(dudz(i ,j,k1)+ &
611 & dudz(i+1,j,k2))+ &
612 & cff1*(tl_dudz(i ,j,k1)+ &
613 & tl_dudz(i+1,j,k2))+ &
614 & tl_cff2*(dudz(i ,j,k2)+ &
615 & dudz(i+1,j,k1))+ &
616 & cff2*(tl_dudz(i ,j,k2)+ &
617 & tl_dudz(i+1,j,k1))))- &
618 & om_r(i,j)*(tl_dmvde(i,j,k1)- &
619 & 0.5_r8*pm(i,j)* &
620 & (tl_cff3*(dvdz(i,j ,k1)+ &
621 & dvdz(i,j+1,k2))+ &
622 & cff3*(tl_dvdz(i,j ,k1)+ &
623 & tl_dvdz(i,j+1,k2))+ &
624 & tl_cff4*(dvdz(i,j ,k2)+ &
625 & dvdz(i,j+1,k1))+ &
626 & cff4*(tl_dvdz(i,j ,k2)+ &
627 & tl_dvdz(i,j+1,k1))))- &
628#ifdef TL_IOMS
629 & (-on_r(i,j)*0.5_r8*pn(i,j)* &
630 & (cff1*(dudz(i ,j,k1)+ &
631 & dudz(i+1,j,k2))+ &
632 & cff2*(dudz(i ,j,k2)+ &
633 & dudz(i+1,j,k1)))+ &
634 & om_r(i,j)*0.5_r8*pm(i,j)* &
635 & (cff3*(dvdz(i,j ,k1)+ &
636 & dvdz(i,j+1,k2))+ &
637 & cff4*(dvdz(i,j ,k2)+ &
638 & dvdz(i,j+1,k1))))
639#endif
640#ifdef MASKING
641 cff=cff*rmask(i,j)
642 tl_cff=tl_cff*rmask(i,j)
643#endif
644#ifdef VISC_3DCOEF
645# ifdef UV_U3ADV_SPLIT
646 ufx(i,j)=on_r(i,j)*on_r(i,j)*uvis3d_r(i,j,k)*cff
647 tl_ufx(i,j)=on_r(i,j)*on_r(i,j)* &
648 & (tl_uvis3d_r(i,j,k)*cff+ &
649 & uvis3d_r(i,j,k)*tl_cff)- &
650# ifdef TL_IOMS
651 & ufx(i,j)
652# endif
653 vfe(i,j)=om_r(i,j)*om_r(i,j)*vvis3d_r(i,j,k)*cff
654 tl_vfe(i,j)=om_r(i,j)*om_r(i,j)* &
655 & (tl_vvis3d_r(i,j,k)*cff+ &
656 & vvis3d_r(i,j,k)*tl_cff)- &
657# ifdef TL_IOMS
658 & vfe(i,j)
659# endif
660# else
661 ufx(i,j)=on_r(i,j)*on_r(i,j)*visc3d_r(i,j,k)*cff
662 tl_ufx(i,j)=on_r(i,j)*on_r(i,j)* &
663 & (tl_visc3d_r(i,j,k)*cff+ &
664 & visc3d_r(i,j,k)*tl_cff)- &
665# ifdef TL_IOMS
666 & ufx(i,j)
667# endif
668 vfe(i,j)=om_r(i,j)*om_r(i,j)*visc3d_r(i,j,k)*cff
669 tl_vfe(i,j)=om_r(i,j)*om_r(i,j)* &
670 & (tl_visc3d_r(i,j,k)*cff+ &
671 & visc3d_r(i,j,k)*tl_cff)- &
672# ifdef TL_IOMS
673 & vfe(i,j)
674# endif
675# endif
676#else
677 ufx(i,j)=on_r(i,j)*on_r(i,j)*visc4_r(i,j)*cff
678 tl_ufx(i,j)=on_r(i,j)*on_r(i,j)*visc4_r(i,j)*tl_cff
679 vfe(i,j)=om_r(i,j)*om_r(i,j)*visc4_r(i,j)*cff
680 tl_vfe(i,j)=om_r(i,j)*om_r(i,j)*visc4_r(i,j)*tl_cff
681#endif
682 END DO
683 END DO
684
685 DO j=jstrm1,jendp2
686 DO i=istrm1,iendp2
687 pm_p=0.25_r8*(pm(i-1,j-1)+pm(i-1,j)+ &
688 & pm(i ,j-1)+pm(i ,j))
689 pn_p=0.25_r8*(pn(i-1,j-1)+pn(i-1,j)+ &
690 & pn(i ,j-1)+pn(i ,j))
691 cff1=min(dzdx_p(i,j,k1),0.0_r8)
692 cff2=max(dzdx_p(i,j,k1),0.0_r8)
693 cff3=min(dzde_p(i,j,k1),0.0_r8)
694 cff4=max(dzde_p(i,j,k1),0.0_r8)
695 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_p(i,j,k1)))* &
696 & tl_dzdx_p(i,j,k1)
697 tl_cff2=(0.5_r8+sign(0.5_r8, dzdx_p(i,j,k1)))* &
698 & tl_dzdx_p(i,j,k1)
699 tl_cff3=(0.5_r8+sign(0.5_r8,-dzde_p(i,j,k1)))* &
700 & tl_dzde_p(i,j,k1)
701 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_p(i,j,k1)))* &
702 & tl_dzde_p(i,j,k1)
703 cff=on_p(i,j)*(dnvdx(i,j,k1)- &
704 & 0.5_r8*pn_p* &
705 & (cff1*(dvdz(i-1,j,k1)+ &
706 & dvdz(i ,j,k2))+ &
707 & cff2*(dvdz(i-1,j,k2)+ &
708 & dvdz(i ,j,k1))))+ &
709 & om_p(i,j)*(dmude(i,j,k1)- &
710 & 0.5_r8*pm_p* &
711 & (cff3*(dudz(i,j-1,k1)+ &
712 & dudz(i,j ,k2))+ &
713 & cff4*(dudz(i,j-1,k2)+ &
714 & dudz(i,j ,k1))))
715 tl_cff=on_p(i,j)*(tl_dnvdx(i,j,k1)- &
716 & 0.5_r8*pn_p* &
717 & (tl_cff1*(dvdz(i-1,j,k1)+ &
718 & dvdz(i ,j,k2))+ &
719 & cff1*(tl_dvdz(i-1,j,k1)+ &
720 & tl_dvdz(i ,j,k2))+ &
721 & tl_cff2*(dvdz(i-1,j,k2)+ &
722 & dvdz(i ,j,k1))+ &
723 & cff2*(tl_dvdz(i-1,j,k2)+ &
724 & tl_dvdz(i ,j,k1))))+ &
725 & om_p(i,j)*(tl_dmude(i,j,k1)- &
726 & 0.5_r8*pm_p* &
727 & (tl_cff3*(dudz(i,j-1,k1)+ &
728 & dudz(i,j ,k2))+ &
729 & cff3*(tl_dudz(i,j-1,k1)+ &
730 & tl_dudz(i,j ,k2))+ &
731 & tl_cff4*(dudz(i,j-1,k2)+ &
732 & dudz(i,j ,k1))+ &
733 & cff4*(tl_dudz(i,j-1,k2)+ &
734 & tl_dudz(i,j ,k1))))- &
735#ifdef TL_IOMS
736 & (-on_p(i,j)*0.5_r8*pn_p* &
737 & (cff1*(dvdz(i-1,j,k1)+ &
738 & dvdz(i ,j,k2))+ &
739 & cff2*(dvdz(i-1,j,k2)+ &
740 & dvdz(i ,j,k1)))- &
741 & om_p(i,j)*0.5_r8*pm_p* &
742 & (cff3*(dudz(i,j-1,k1)+ &
743 & dudz(i,j ,k2))+ &
744 & cff4*(dudz(i,j-1,k2)+ &
745 & dudz(i,j ,k1))))
746#endif
747#ifdef MASKING
748 cff=cff*pmask(i,j)
749 tl_cff=tl_cff*pmask(i,j)
750#endif
751#ifdef VISC_3DCOEF
752# ifdef UV_U3ADV_SPLIT
753 uvis_p=0.25_r8* &
754 & (uvis3d_r(i-1,j-1,k)+uvis3d_r(i-1,j,k)+ &
755 & uvis3d_r(i ,j-1,k)+uvis3d_r(i ,j,k))
756 tl_uvis_p=0.25_r8* &
757 & (tl_uvis3d_r(i-1,j-1,k)+tl_uvis3d_r(i-1,j,k)+ &
758 & tl_uvis3d_r(i ,j-1,k)+tl_uvis3d_r(i ,j,k))
759 vvis_p=0.25_r8* &
760 & (vvis3d_r(i-1,j-1,k)+vvis3d_r(i-1,j,k)+ &
761 & vvis3d_r(i ,j-1,k)+vvis3d_r(i ,j,k))
762 tl_vvis_p=0.25_r8* &
763 & (tl_vvis3d_r(i-1,j-1,k)+tl_vvis3d_r(i-1,j,k)+ &
764 & tl_vvis3d_r(i ,j-1,k)+tl_vvis3d_r(i ,j,k))
765 ufe(i,j)=om_p(i,j)*om_p(i,j)*uvis_p*cff
766 tl_ufe(i,j)=om_p(i,j)*om_p(i,j)* &
767 & (tl_uvis_p*cff+uvis_p*tl_cff)- &
768# ifdef TL_IOMS
769 & ufe(i,j)
770# endif
771 vfx(i,j)=on_p(i,j)*on_p(i,j)*vvis_p*cff
772 tl_vfx(i,j)=on_p(i,j)*on_p(i,j)* &
773 & (tl_vvis_p*cff+vvis_p*tl_cff)- &
774# ifdef TL_IOMS
775 & vfx(i,j)
776# endif
777# else
778 visc_p=0.25_r8* &
779 & (visc3d_r(i-1,j-1,k)+visc3d_r(i-1,j,k)+ &
780 & visc3d_r(i ,j-1,k)+visc3d_r(i ,j,k))
781 tl_visc_p=0.25_r8* &
782 & (tl_visc3d_r(i-1,j-1,k)+tl_visc3d_r(i-1,j,k)+ &
783 & tl_visc3d_r(i ,j-1,k)+tl_visc3d_r(i ,j,k))
784 ufe(i,j)=om_p(i,j)*om_p(i,j)*visc_p*cff
785 tl_ufe(i,j)=om_p(i,j)*om_p(i,j)* &
786 & (tl_visc_p*cff+visc_p*tl_cff)- &
787# ifdef TL_IOMS
788 & ufe(i,j)
789# endif
790 vfx(i,j)=on_p(i,j)*on_p(i,j)*visc_p*cff
791 tl_vfx(i,j)=on_p(i,j)*on_p(i,j)* &
792 & (tl_visc_p*cff+visc_p*tl_cff)- &
793# ifdef TL_IOMS
794 & vfx(i,j)
795# endif
796# endif
797#else
798 ufe(i,j)=om_p(i,j)*om_p(i,j)*visc4_p(i,j)*cff
799 tl_ufe(i,j)=om_p(i,j)*om_p(i,j)*visc4_p(i,j)*tl_cff
800 vfx(i,j)=on_p(i,j)*on_p(i,j)*visc4_p(i,j)*cff
801 tl_vfx(i,j)=on_p(i,j)*on_p(i,j)*visc4_p(i,j)*tl_cff
802#endif
803 END DO
804 END DO
805!
806! Compute vertical flux (m^2 s^-3/2) due to sloping terrain-following
807! surfaces.
808!
809 IF (k.lt.n(ng)) THEN
810 DO j=jstrm1,jendp1
811 DO i=istrum1,iendp1
812#ifdef VISC_3DCOEF
813# ifdef UV_U3ADV_SPLIT
814 cff=0.125_r8* &
815 & (uvis3d_r(i-1,j,k )+uvis3d_r(i,j,k )+ &
816 & uvis3d_r(i-1,j,k+1)+uvis3d_r(i,j,k+1))
817 tl_cff=0.125_r8* &
818 & (tl_uvis3d_r(i-1,j,k )+tl_uvis3d_r(i,j,k )+ &
819 & tl_uvis3d_r(i-1,j,k+1)+tl_uvis3d_r(i,j,k+1))
820# else
821 cff=0.125_r8* &
822 & (visc3d_r(i-1,j,k )+visc3d_r(i,j,k )+ &
823 & visc3d_r(i-1,j,k+1)+visc3d_r(i,j,k+1))
824 tl_cff=0.125_r8* &
825 & (tl_visc3d_r(i-1,j,k )+tl_visc3d_r(i,j,k )+ &
826 & tl_visc3d_r(i-1,j,k+1)+tl_visc3d_r(i,j,k+1))
827# endif
828 fac1=cff*on_u(i,j)
829 tl_fac1=tl_cff*on_u(i,j)
830 fac2=cff*om_u(i,j)
831 tl_fac2=tl_cff*om_u(i,j)
832#else
833 cff=0.25_r8*(visc4_r(i-1,j)+visc4_r(i,j))
834 fac1=cff*on_u(i,j)
835 fac2=cff*om_u(i,j)
836#endif
837 cff=0.5_r8*(pn(i-1,j)+pn(i,j))
838 dnudz=cff*dudz(i,j,k2)
839 tl_dnudz=cff*tl_dudz(i,j,k2)
840 dnvdz=cff*0.25_r8*(dvdz(i-1,j+1,k2)+ &
841 & dvdz(i ,j+1,k2)+ &
842 & dvdz(i-1,j ,k2)+ &
843 & dvdz(i ,j ,k2))
844 tl_dnvdz=cff*0.25_r8*(tl_dvdz(i-1,j+1,k2)+ &
845 & tl_dvdz(i ,j+1,k2)+ &
846 & tl_dvdz(i-1,j ,k2)+ &
847 & tl_dvdz(i ,j ,k2))
848 cff=0.5_r8*(pm(i-1,j)+pm(i,j))
849 dmudz=cff*dudz(i,j,k2)
850 tl_dmudz=cff*tl_dudz(i,j,k2)
851 dmvdz=cff*0.25_r8*(dvdz(i-1,j+1,k2)+ &
852 & dvdz(i ,j+1,k2)+ &
853 & dvdz(i-1,j ,k2)+ &
854 & dvdz(i ,j ,k2))
855 tl_dmvdz=cff*0.25_r8*(tl_dvdz(i-1,j+1,k2)+ &
856 & tl_dvdz(i ,j+1,k2)+ &
857 & tl_dvdz(i-1,j ,k2)+ &
858 & tl_dvdz(i ,j ,k2))
859
860 cff1=min(dzdx_r(i-1,j,k1),0.0_r8)
861 cff2=min(dzdx_r(i ,j,k2),0.0_r8)
862 cff3=max(dzdx_r(i-1,j,k2),0.0_r8)
863 cff4=max(dzdx_r(i ,j,k1),0.0_r8)
864 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_r(i-1,j,k1)))* &
865 & tl_dzdx_r(i-1,j,k1)
866 tl_cff2=(0.5_r8+sign(0.5_r8,-dzdx_r(i ,j,k2)))* &
867 & tl_dzdx_r(i ,j,k2)
868 tl_cff3=(0.5_r8+sign(0.5_r8, dzdx_r(i-1,j,k2)))* &
869 & tl_dzdx_r(i-1,j,k2)
870 tl_cff4=(0.5_r8+sign(0.5_r8, dzdx_r(i ,j,k1)))* &
871 & tl_dzdx_r(i ,j,k1)
872 ufsx(i,j,k2)=fac1* &
873 & (cff1*(cff1*dnudz-dnudx(i-1,j,k1))+ &
874 & cff2*(cff2*dnudz-dnudx(i ,j,k2))+ &
875 & cff3*(cff3*dnudz-dnudx(i-1,j,k2))+ &
876 & cff4*(cff4*dnudz-dnudx(i ,j,k1)))
877 tl_ufsx(i,j,k2)=fac1* &
878 & (tl_cff1*(cff1*dnudz-dnudx(i-1,j,k1))+ &
879 & tl_cff2*(cff2*dnudz-dnudx(i ,j,k2))+ &
880 & tl_cff3*(cff3*dnudz-dnudx(i-1,j,k2))+ &
881 & tl_cff4*(cff4*dnudz-dnudx(i ,j,k1))+ &
882 & cff1*(tl_cff1*dnudz+cff1*tl_dnudz- &
883 & tl_dnudx(i-1,j,k1))+ &
884 & cff2*(tl_cff2*dnudz+cff2*tl_dnudz- &
885 & tl_dnudx(i ,j,k2))+ &
886 & cff3*(tl_cff3*dnudz+cff3*tl_dnudz- &
887 & tl_dnudx(i-1,j,k2))+ &
888 & cff4*(tl_cff4*dnudz+cff4*tl_dnudz- &
889 & tl_dnudx(i ,j,k1)))- &
890#ifdef TL_IOMS
891 & fac1* &
892 & (cff1*(2.0_r8*cff1*dnudz- &
893 & dnudx(i-1,j,k1))+ &
894 & cff2*(2.0_r8*cff2*dnudz- &
895 & dnudx(i ,j,k2))+ &
896 & cff3*(2.0_r8*cff3*dnudz- &
897 & dnudx(i-1,j,k2))+ &
898 & cff4*(2.0_r8*cff4*dnudz- &
899 & dnudx(i ,j,k1)))
900#endif
901#ifdef VISC_3DCOEF
902 tl_ufsx(i,j,k2)=tl_ufsx(i,j,k2)+ &
903 & tl_fac1* &
904 & (cff1*(cff1*dnudz-dnudx(i-1,j,k1))+ &
905 & cff2*(cff2*dnudz-dnudx(i ,j,k2))+ &
906 & cff3*(cff3*dnudz-dnudx(i-1,j,k2))+ &
907 & cff4*(cff4*dnudz-dnudx(i ,j,k1)))
908#endif
909
910 cff1=min(dzde_p(i,j ,k1),0.0_r8)
911 cff2=min(dzde_p(i,j+1,k2),0.0_r8)
912 cff3=max(dzde_p(i,j ,k2),0.0_r8)
913 cff4=max(dzde_p(i,j+1,k1),0.0_r8)
914 tl_cff1=(0.5_r8+sign(0.5_r8,-dzde_p(i,j ,k1)))* &
915 & tl_dzde_p(i,j ,k1)
916 tl_cff2=(0.5_r8+sign(0.5_r8,-dzde_p(i,j+1,k2)))* &
917 & tl_dzde_p(i,j+1,k2)
918 tl_cff3=(0.5_r8+sign(0.5_r8, dzde_p(i,j ,k2)))* &
919 & tl_dzde_p(i,j ,k2)
920 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_p(i,j+1,k1)))* &
921 & tl_dzde_p(i,j+1,k1)
922 ufse(i,j,k2)=fac2* &
923 & (cff1*(cff1*dmudz-dmude(i,j ,k1))+ &
924 & cff2*(cff2*dmudz-dmude(i,j+1,k2))+ &
925 & cff3*(cff3*dmudz-dmude(i,j ,k2))+ &
926 & cff4*(cff4*dmudz-dmude(i,j+1,k1)))
927 tl_ufse(i,j,k2)=fac2* &
928 & (tl_cff1*(cff1*dmudz-dmude(i,j ,k1))+ &
929 & tl_cff2*(cff2*dmudz-dmude(i,j+1,k2))+ &
930 & tl_cff3*(cff3*dmudz-dmude(i,j ,k2))+ &
931 & tl_cff4*(cff4*dmudz-dmude(i,j+1,k1))+ &
932 & cff1*(tl_cff1*dmudz+cff1*tl_dmudz- &
933 & tl_dmude(i,j ,k1))+ &
934 & cff2*(tl_cff2*dmudz+cff2*tl_dmudz- &
935 & tl_dmude(i,j+1,k2))+ &
936 & cff3*(tl_cff3*dmudz+cff3*tl_dmudz- &
937 & tl_dmude(i,j ,k2))+ &
938 & cff4*(tl_cff4*dmudz+cff4*tl_dmudz- &
939 & tl_dmude(i,j+1,k1)))- &
940#ifdef TL_IOMS
941 & fac2* &
942 & (cff1*(2.0_r8*cff1*dmudz- &
943 & dmude(i,j ,k1))+ &
944 & cff2*(2.0_r8*cff2*dmudz- &
945 & dmude(i,j+1,k2))+ &
946 & cff3*(2.0_r8*cff3*dmudz- &
947 & dmude(i,j ,k2))+ &
948 & cff4*(2.0_r8*cff4*dmudz- &
949 & dmude(i,j+1,k1)))
950#endif
951#ifdef VISC_3DCOEF
952 tl_ufse(i,j,k2)=tl_ufse(i,j,k2)+ &
953 & tl_fac2* &
954 & (cff1*(cff1*dmudz-dmude(i,j ,k1))+ &
955 & cff2*(cff2*dmudz-dmude(i,j+1,k2))+ &
956 & cff3*(cff3*dmudz-dmude(i,j ,k2))+ &
957 & cff4*(cff4*dmudz-dmude(i,j+1,k1)))
958#endif
959
960 cff1=min(dzde_p(i,j ,k1),0.0_r8)
961 cff2=min(dzde_p(i,j+1,k2),0.0_r8)
962 cff3=max(dzde_p(i,j ,k2),0.0_r8)
963 cff4=max(dzde_p(i,j+1,k1),0.0_r8)
964 cff5=min(dzdx_p(i,j ,k1),0.0_r8)
965 cff6=min(dzdx_p(i,j+1,k2),0.0_r8)
966 cff7=max(dzdx_p(i,j ,k2),0.0_r8)
967 cff8=max(dzdx_p(i,j+1,k1),0.0_r8)
968 tl_cff1=(0.5_r8+sign(0.5_r8,-dzde_p(i,j ,k1)))* &
969 & tl_dzde_p(i,j ,k1)
970 tl_cff2=(0.5_r8+sign(0.5_r8,-dzde_p(i,j+1,k2)))* &
971 & tl_dzde_p(i,j+1,k2)
972 tl_cff3=(0.5_r8+sign(0.5_r8, dzde_p(i,j ,k2)))* &
973 & tl_dzde_p(i,j ,k2)
974 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_p(i,j+1,k1)))* &
975 & tl_dzde_p(i,j+1,k1)
976 tl_cff5=(0.5_r8+sign(0.5_r8,-dzdx_p(i,j ,k1)))* &
977 & tl_dzdx_p(i,j ,k1)
978 tl_cff6=(0.5_r8+sign(0.5_r8,-dzdx_p(i,j+1,k2)))* &
979 & tl_dzdx_p(i,j+1,k2)
980 tl_cff7=(0.5_r8+sign(0.5_r8, dzdx_p(i,j ,k2)))* &
981 & tl_dzdx_p(i,j ,k2)
982 tl_cff8=(0.5_r8+sign(0.5_r8, dzdx_p(i,j+1,k1)))* &
983 & tl_dzdx_p(i,j+1,k1)
984 ufsx(i,j,k2)=ufsx(i,j,k2)+ &
985 & fac1* &
986 & (cff1*(cff5*dnvdz-dnvdx(i,j ,k1))+ &
987 & cff2*(cff6*dnvdz-dnvdx(i,j+1,k2))+ &
988 & cff3*(cff7*dnvdz-dnvdx(i,j ,k2))+ &
989 & cff4*(cff8*dnvdz-dnvdx(i,j+1,k1)))
990 tl_ufsx(i,j,k2)=tl_ufsx(i,j,k2)+ &
991 & fac1* &
992 & (tl_cff1*(cff5*dnvdz-dnvdx(i,j ,k1))+ &
993 & tl_cff2*(cff6*dnvdz-dnvdx(i,j+1,k2))+ &
994 & tl_cff3*(cff7*dnvdz-dnvdx(i,j ,k2))+ &
995 & tl_cff4*(cff8*dnvdz-dnvdx(i,j+1,k1))+ &
996 & cff1*(tl_cff5*dnvdz+cff5*tl_dnvdz- &
997 & tl_dnvdx(i,j ,k1))+ &
998 & cff2*(tl_cff6*dnvdz+cff6*tl_dnvdz- &
999 & tl_dnvdx(i,j+1,k2))+ &
1000 & cff3*(tl_cff7*dnvdz+cff7*tl_dnvdz- &
1001 & tl_dnvdx(i,j ,k2))+ &
1002 & cff4*(tl_cff8*dnvdz+cff8*tl_dnvdz- &
1003 & tl_dnvdx(i,j+1,k1)))- &
1004#ifdef TL_IOMS
1005 & fac1* &
1006 & (cff1*(2.0_r8*cff5*dnvdz- &
1007 & dnvdx(i,j ,k1))+ &
1008 & cff2*(2.0_r8*cff6*dnvdz- &
1009 & dnvdx(i,j+1,k2))+ &
1010 & cff3*(2.0_r8*cff7*dnvdz- &
1011 & dnvdx(i,j ,k2))+ &
1012 & cff4*(2.0_r8*cff8*dnvdz- &
1013 & dnvdx(i,j+1,k1)))
1014#endif
1015#ifdef VISC_3DCOEF
1016 tl_ufsx(i,j,k2)=tl_ufsx(i,j,k2)+ &
1017 & tl_fac1* &
1018 & (cff1*(cff5*dnvdz-dnvdx(i,j ,k1))+ &
1019 & cff2*(cff6*dnvdz-dnvdx(i,j+1,k2))+ &
1020 & cff3*(cff7*dnvdz-dnvdx(i,j ,k2))+ &
1021 & cff4*(cff8*dnvdz-dnvdx(i,j+1,k1)))
1022#endif
1023
1024 cff1=min(dzdx_r(i-1,j,k1),0.0_r8)
1025 cff2=min(dzdx_r(i ,j,k2),0.0_r8)
1026 cff3=max(dzdx_r(i-1,j,k2),0.0_r8)
1027 cff4=max(dzdx_r(i ,j,k1),0.0_r8)
1028 cff5=min(dzde_r(i-1,j,k1),0.0_r8)
1029 cff6=min(dzde_r(i ,j,k2),0.0_r8)
1030 cff7=max(dzde_r(i-1,j,k2),0.0_r8)
1031 cff8=max(dzde_r(i ,j,k1),0.0_r8)
1032 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_r(i-1,j,k1)))* &
1033 & tl_dzdx_r(i-1,j,k1)
1034 tl_cff2=(0.5_r8+sign(0.5_r8,-dzdx_r(i ,j,k2)))* &
1035 & tl_dzdx_r(i ,j,k2)
1036 tl_cff3=(0.5_r8+sign(0.5_r8, dzdx_r(i-1,j,k2)))* &
1037 & tl_dzdx_r(i-1,j,k2)
1038 tl_cff4=(0.5_r8+sign(0.5_r8, dzdx_r(i ,j,k1)))* &
1039 & tl_dzdx_r(i ,j,k1)
1040 tl_cff5=(0.5_r8+sign(0.5_r8,-dzde_r(i-1,j,k1)))* &
1041 & tl_dzde_r(i-1,j,k1)
1042 tl_cff6=(0.5_r8+sign(0.5_r8,-dzde_r(i ,j,k2)))* &
1043 & tl_dzde_r(i ,j,k2)
1044 tl_cff7=(0.5_r8+sign(0.5_r8, dzde_r(i-1,j,k2)))* &
1045 & tl_dzde_r(i-1,j,k2)
1046 tl_cff8=(0.5_r8+sign(0.5_r8, dzde_r(i ,j,k1)))* &
1047 & tl_dzde_r(i ,j,k1)
1048 ufse(i,j,k2)=ufse(i,j,k2)- &
1049 & fac2* &
1050 & (cff1*(cff5*dmvdz-dmvde(i-1,j,k1))+ &
1051 & cff2*(cff6*dmvdz-dmvde(i ,j,k2))+ &
1052 & cff3*(cff7*dmvdz-dmvde(i-1,j,k2))+ &
1053 & cff4*(cff8*dmvdz-dmvde(i ,j,k1)))
1054 tl_ufse(i,j,k2)=tl_ufse(i,j,k2)- &
1055 & fac2* &
1056 & (tl_cff1*(cff5*dmvdz-dmvde(i-1,j,k1))+ &
1057 & tl_cff2*(cff6*dmvdz-dmvde(i ,j,k2))+ &
1058 & tl_cff3*(cff7*dmvdz-dmvde(i-1,j,k2))+ &
1059 & tl_cff4*(cff8*dmvdz-dmvde(i ,j,k1))+ &
1060 & cff1*(tl_cff5*dmvdz+cff5*tl_dmvdz- &
1061 & tl_dmvde(i-1,j,k1))+ &
1062 & cff2*(tl_cff6*dmvdz+cff6*tl_dmvdz- &
1063 & tl_dmvde(i ,j,k2))+ &
1064 & cff3*(tl_cff7*dmvdz+cff7*tl_dmvdz- &
1065 & tl_dmvde(i-1,j,k2))+ &
1066 & cff4*(tl_cff8*dmvdz+cff8*tl_dmvdz- &
1067 & tl_dmvde(i ,j,k1)))+ &
1068#ifdef TL_IOMS
1069 & fac2* &
1070 & (cff1*(2.0_r8*cff5*dmvdz- &
1071 & dmvde(i-1,j,k1))+ &
1072 & cff2*(2.0_r8*cff6*dmvdz- &
1073 & dmvde(i ,j,k2))+ &
1074 & cff3*(2.0_r8*cff7*dmvdz- &
1075 & dmvde(i-1,j,k2))+ &
1076 & cff4*(2.0_r8*cff8*dmvdz- &
1077 & dmvde(i ,j,k1)))
1078#endif
1079#ifdef VISC_3DCOEF
1080 tl_ufse(i,j,k2)=tl_ufse(i,j,k2)- &
1081 & tl_fac2* &
1082 & (cff1*(cff5*dmvdz-dmvde(i-1,j,k1))+ &
1083 & cff2*(cff6*dmvdz-dmvde(i ,j,k2))+ &
1084 & cff3*(cff7*dmvdz-dmvde(i-1,j,k2))+ &
1085 & cff4*(cff8*dmvdz-dmvde(i ,j,k1)))
1086#endif
1087 END DO
1088 END DO
1089!
1090 DO j=jstrvm1,jendp1
1091 DO i=istrm1,iendp1
1092#ifdef VISC_3DCOEF
1093# ifdef UV_U3ADV_SPLIT
1094 cff=0.125_r8* &
1095 & (vvis3d_r(i,j-1,k )+vvis3d_r(i,j,k )+ &
1096 & vvis3d_r(i,j-1,k+1)+vvis3d_r(i,j,k+1))
1097 tl_cff=0.125_r8* &
1098 & (tl_vvis3d_r(i,j-1,k )+tl_vvis3d_r(i,j,k )+ &
1099 & tl_vvis3d_r(i,j-1,k+1)+tl_vvis3d_r(i,j,k+1))
1100# else
1101 cff=0.125_r8* &
1102 & (visc3d_r(i,j-1,k )+visc3d_r(i,j,k )+ &
1103 & visc3d_r(i,j-1,k+1)+visc3d_r(i,j,k+1))
1104 tl_cff=0.125_r8* &
1105 & (tl_visc3d_r(i,j-1,k )+tl_visc3d_r(i,j,k )+ &
1106 & tl_visc3d_r(i,j-1,k+1)+tl_visc3d_r(i,j,k+1))
1107# endif
1108 fac1=cff*on_v(i,j)
1109 tl_fac1=tl_cff*on_v(i,j)
1110 fac2=cff*om_v(i,j)
1111 tl_fac2=tl_cff*om_v(i,j)
1112#else
1113 cff=0.25_r8*(visc4_r(i,j-1)+visc4_r(i,j))
1114 fac1=cff*on_v(i,j)
1115 fac2=cff*om_v(i,j)
1116#endif
1117 cff=0.5_r8*(pn(i,j-1)+pn(i,j))
1118 dnudz=cff*0.25_r8*(dudz(i ,j ,k2)+ &
1119 & dudz(i+1,j ,k2)+ &
1120 & dudz(i ,j-1,k2)+ &
1121 & dudz(i+1,j-1,k2))
1122 tl_dnudz=cff*0.25_r8*(tl_dudz(i ,j ,k2)+ &
1123 & tl_dudz(i+1,j ,k2)+ &
1124 & tl_dudz(i ,j-1,k2)+ &
1125 & tl_dudz(i+1,j-1,k2))
1126 dnvdz=cff*dvdz(i,j,k2)
1127 tl_dnvdz=cff*tl_dvdz(i,j,k2)
1128 cff=0.5_r8*(pm(i,j-1)+pm(i,j))
1129 dmudz=cff*0.25_r8*(dudz(i ,j ,k2)+ &
1130 & dudz(i+1,j ,k2)+ &
1131 & dudz(i ,j-1,k2)+ &
1132 & dudz(i+1,j-1,k2))
1133 tl_dmudz=cff*0.25_r8*(tl_dudz(i ,j ,k2)+ &
1134 & tl_dudz(i+1,j ,k2)+ &
1135 & tl_dudz(i ,j-1,k2)+ &
1136 & tl_dudz(i+1,j-1,k2))
1137 dmvdz=cff*dvdz(i,j,k2)
1138 tl_dmvdz=cff*tl_dvdz(i,j,k2)
1139
1140 cff1=min(dzdx_p(i ,j,k1),0.0_r8)
1141 cff2=min(dzdx_p(i+1,j,k2),0.0_r8)
1142 cff3=max(dzdx_p(i ,j,k2),0.0_r8)
1143 cff4=max(dzdx_p(i+1,j,k1),0.0_r8)
1144 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_p(i ,j,k1)))* &
1145 & tl_dzdx_p(i ,j,k1)
1146 tl_cff2=(0.5_r8+sign(0.5_r8,-dzdx_p(i+1,j,k2)))* &
1147 & tl_dzdx_p(i+1,j,k2)
1148 tl_cff3=(0.5_r8+sign(0.5_r8, dzdx_p(i ,j,k2)))* &
1149 & tl_dzdx_p(i ,j,k2)
1150 tl_cff4=(0.5_r8+sign(0.5_r8, dzdx_p(i+1,j,k1)))* &
1151 & tl_dzdx_p(i+1,j,k1)
1152 vfsx(i,j,k2)=fac1* &
1153 & (cff1*(cff1*dnvdz-dnvdx(i ,j,k1))+ &
1154 & cff2*(cff2*dnvdz-dnvdx(i+1,j,k2))+ &
1155 & cff3*(cff3*dnvdz-dnvdx(i ,j,k2))+ &
1156 & cff4*(cff4*dnvdz-dnvdx(i+1,j,k1)))
1157 tl_vfsx(i,j,k2)=fac1* &
1158 & (tl_cff1*(cff1*dnvdz-dnvdx(i ,j,k1))+ &
1159 & tl_cff2*(cff2*dnvdz-dnvdx(i+1,j,k2))+ &
1160 & tl_cff3*(cff3*dnvdz-dnvdx(i ,j,k2))+ &
1161 & tl_cff4*(cff4*dnvdz-dnvdx(i+1,j,k1))+ &
1162 & cff1*(tl_cff1*dnvdz+cff1*tl_dnvdz- &
1163 & tl_dnvdx(i ,j,k1))+ &
1164 & cff2*(tl_cff2*dnvdz+cff2*tl_dnvdz- &
1165 & tl_dnvdx(i+1,j,k2))+ &
1166 & cff3*(tl_cff3*dnvdz+cff3*tl_dnvdz- &
1167 & tl_dnvdx(i ,j,k2))+ &
1168 & cff4*(tl_cff4*dnvdz+cff4*tl_dnvdz- &
1169 & tl_dnvdx(i+1,j,k1)))- &
1170#ifdef TL_IOMS
1171 & fac1* &
1172 & (cff1*(2.0_r8*cff1*dnvdz- &
1173 & dnvdx(i ,j,k1))+ &
1174 & cff2*(2.0_r8*cff2*dnvdz- &
1175 & dnvdx(i+1,j,k2))+ &
1176 & cff3*(2.0_r8*cff3*dnvdz- &
1177 & dnvdx(i ,j,k2))+ &
1178 & cff4*(2.0_r8*cff4*dnvdz- &
1179 & dnvdx(i+1,j,k1)))
1180#endif
1181#ifdef VISC_3DCOEF
1182 tl_vfsx(i,j,k2)=tl_vfsx(i,j,k2)+ &
1183 & tl_fac1* &
1184 & (cff1*(cff1*dnvdz-dnvdx(i ,j,k1))+ &
1185 & cff2*(cff2*dnvdz-dnvdx(i+1,j,k2))+ &
1186 & cff3*(cff3*dnvdz-dnvdx(i ,j,k2))+ &
1187 & cff4*(cff4*dnvdz-dnvdx(i+1,j,k1)))
1188#endif
1189
1190 cff1=min(dzde_r(i,j-1,k1),0.0_r8)
1191 cff2=min(dzde_r(i,j ,k2),0.0_r8)
1192 cff3=max(dzde_r(i,j-1,k2),0.0_r8)
1193 cff4=max(dzde_r(i,j ,k1),0.0_r8)
1194 tl_cff1=(0.5_r8+sign(0.5_r8,-dzde_r(i,j-1,k1)))* &
1195 & tl_dzde_r(i,j-1,k1)
1196 tl_cff2=(0.5_r8+sign(0.5_r8,-dzde_r(i,j ,k2)))* &
1197 & tl_dzde_r(i,j ,k2)
1198 tl_cff3=(0.5_r8+sign(0.5_r8, dzde_r(i,j-1,k2)))* &
1199 & tl_dzde_r(i,j-1,k2)
1200 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_r(i,j ,k1)))* &
1201 & tl_dzde_r(i,j ,k1)
1202 vfse(i,j,k2)=fac2* &
1203 & (cff1*(cff1*dmvdz-dmvde(i,j-1,k1))+ &
1204 & cff2*(cff2*dmvdz-dmvde(i,j ,k2))+ &
1205 & cff3*(cff3*dmvdz-dmvde(i,j-1,k2))+ &
1206 & cff4*(cff4*dmvdz-dmvde(i,j ,k1)))
1207 tl_vfse(i,j,k2)=fac2* &
1208 & (tl_cff1*(cff1*dmvdz-dmvde(i,j-1,k1))+ &
1209 & tl_cff2*(cff2*dmvdz-dmvde(i,j ,k2))+ &
1210 & tl_cff3*(cff3*dmvdz-dmvde(i,j-1,k2))+ &
1211 & tl_cff4*(cff4*dmvdz-dmvde(i,j ,k1))+ &
1212 & cff1*(tl_cff1*dmvdz+cff1*tl_dmvdz- &
1213 & tl_dmvde(i,j-1,k1))+ &
1214 & cff2*(tl_cff2*dmvdz+cff2*tl_dmvdz- &
1215 & tl_dmvde(i,j ,k2))+ &
1216 & cff3*(tl_cff3*dmvdz+cff3*tl_dmvdz- &
1217 & tl_dmvde(i,j-1,k2))+ &
1218 & cff4*(tl_cff4*dmvdz+cff4*tl_dmvdz- &
1219 & tl_dmvde(i,j ,k1)))- &
1220#ifdef TL_IOMS
1221 & fac2* &
1222 & (cff1*(2.0_r8*cff1*dmvdz- &
1223 & dmvde(i,j-1,k1))+ &
1224 & cff2*(2.0_r8*cff2*dmvdz- &
1225 & dmvde(i,j ,k2))+ &
1226 & cff3*(2.0_r8*cff3*dmvdz- &
1227 & dmvde(i,j-1,k2))+ &
1228 & cff4*(2.0_r8*cff4*dmvdz- &
1229 & dmvde(i,j ,k1)))
1230#endif
1231#ifdef VISC_3DCOEF
1232 tl_vfse(i,j,k2)=tl_vfse(i,j,k2)+ &
1233 & tl_fac2* &
1234 & (cff1*(cff1*dmvdz-dmvde(i,j-1,k1))+ &
1235 & cff2*(cff2*dmvdz-dmvde(i,j ,k2))+ &
1236 & cff3*(cff3*dmvdz-dmvde(i,j-1,k2))+ &
1237 & cff4*(cff4*dmvdz-dmvde(i,j ,k1)))
1238#endif
1239
1240 cff1=min(dzde_r(i,j-1,k1),0.0_r8)
1241 cff2=min(dzde_r(i,j ,k2),0.0_r8)
1242 cff3=max(dzde_r(i,j-1,k2),0.0_r8)
1243 cff4=max(dzde_r(i,j ,k1),0.0_r8)
1244 cff5=min(dzdx_r(i,j-1,k1),0.0_r8)
1245 cff6=min(dzdx_r(i,j ,k2),0.0_r8)
1246 cff7=max(dzdx_r(i,j-1,k2),0.0_r8)
1247 cff8=max(dzdx_r(i,j ,k1),0.0_r8)
1248 tl_cff1=(0.5_r8+sign(0.5_r8,-dzde_r(i,j-1,k1)))* &
1249 & tl_dzde_r(i,j-1,k1)
1250 tl_cff2=(0.5_r8+sign(0.5_r8,-dzde_r(i,j ,k2)))* &
1251 & tl_dzde_r(i,j ,k2)
1252 tl_cff3=(0.5_r8+sign(0.5_r8, dzde_r(i,j-1,k2)))* &
1253 & tl_dzde_r(i,j-1,k2)
1254 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_r(i,j ,k1)))* &
1255 & tl_dzde_r(i,j ,k1)
1256 tl_cff5=(0.5_r8+sign(0.5_r8,-dzdx_r(i,j-1,k1)))* &
1257 & tl_dzdx_r(i,j-1,k1)
1258 tl_cff6=(0.5_r8+sign(0.5_r8,-dzdx_r(i,j ,k2)))* &
1259 & tl_dzdx_r(i,j ,k2)
1260 tl_cff7=(0.5_r8+sign(0.5_r8, dzdx_r(i,j-1,k2)))* &
1261 & tl_dzdx_r(i,j-1,k2)
1262 tl_cff8=(0.5_r8+sign(0.5_r8, dzdx_r(i,j ,k1)))* &
1263 & tl_dzdx_r(i,j ,k1)
1264 vfsx(i,j,k2)=vfsx(i,j,k2)- &
1265 & fac1* &
1266 & (cff1*(cff5*dnudz-dnudx(i,j-1,k1))+ &
1267 & cff2*(cff6*dnudz-dnudx(i,j ,k2))+ &
1268 & cff3*(cff7*dnudz-dnudx(i,j-1,k2))+ &
1269 & cff4*(cff8*dnudz-dnudx(i,j ,k1)))
1270 tl_vfsx(i,j,k2)=tl_vfsx(i,j,k2)- &
1271 & fac1* &
1272 & (tl_cff1*(cff5*dnudz-dnudx(i,j-1,k1))+ &
1273 & tl_cff2*(cff6*dnudz-dnudx(i,j ,k2))+ &
1274 & tl_cff3*(cff7*dnudz-dnudx(i,j-1,k2))+ &
1275 & tl_cff4*(cff8*dnudz-dnudx(i,j ,k1))+ &
1276 & cff1*(tl_cff5*dnudz+cff5*tl_dnudz- &
1277 & tl_dnudx(i,j-1,k1))+ &
1278 & cff2*(tl_cff6*dnudz+cff6*tl_dnudz- &
1279 & tl_dnudx(i,j ,k2))+ &
1280 & cff3*(tl_cff7*dnudz+cff7*tl_dnudz- &
1281 & tl_dnudx(i,j-1,k2))+ &
1282 & cff4*(tl_cff8*dnudz+cff8*tl_dnudz- &
1283 & tl_dnudx(i,j ,k1)))+ &
1284#ifdef TL_IOMS
1285 & fac1* &
1286 & (cff1*(2.0_r8*cff5*dnudz- &
1287 & dnudx(i,j-1,k1))+ &
1288 & cff2*(2.0_r8*cff6*dnudz- &
1289 & dnudx(i,j ,k2))+ &
1290 & cff3*(2.0_r8*cff7*dnudz- &
1291 & dnudx(i,j-1,k2))+ &
1292 & cff4*(2.0_r8*cff8*dnudz- &
1293 & dnudx(i,j ,k1)))
1294#endif
1295#ifdef VISC_3DCOEF
1296 tl_vfsx(i,j,k2)=tl_vfsx(i,j,k2)- &
1297 & tl_fac1* &
1298 & (cff1*(cff5*dnudz-dnudx(i,j-1,k1))+ &
1299 & cff2*(cff6*dnudz-dnudx(i,j ,k2))+ &
1300 & cff3*(cff7*dnudz-dnudx(i,j-1,k2))+ &
1301 & cff4*(cff8*dnudz-dnudx(i,j ,k1)))
1302#endif
1303
1304 cff1=min(dzdx_p(i ,j,k1),0.0_r8)
1305 cff2=min(dzdx_p(i+1,j,k2),0.0_r8)
1306 cff3=max(dzdx_p(i ,j,k2),0.0_r8)
1307 cff4=max(dzdx_p(i+1,j,k1),0.0_r8)
1308 cff5=min(dzde_p(i ,j,k1),0.0_r8)
1309 cff6=min(dzde_p(i+1,j,k2),0.0_r8)
1310 cff7=max(dzde_p(i ,j,k2),0.0_r8)
1311 cff8=max(dzde_p(i+1,j,k1),0.0_r8)
1312 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_p(i ,j,k1)))* &
1313 & tl_dzdx_p(i ,j,k1)
1314 tl_cff2=(0.5_r8+sign(0.5_r8,-dzdx_p(i+1,j,k2)))* &
1315 & tl_dzdx_p(i+1,j,k2)
1316 tl_cff3=(0.5_r8+sign(0.5_r8, dzdx_p(i ,j,k2)))* &
1317 & tl_dzdx_p(i ,j,k2)
1318 tl_cff4=(0.5_r8+sign(0.5_r8, dzdx_p(i+1,j,k1)))* &
1319 & tl_dzdx_p(i+1,j,k1)
1320 tl_cff5=(0.5_r8+sign(0.5_r8,-dzde_p(i ,j,k1)))* &
1321 & tl_dzde_p(i ,j,k1)
1322 tl_cff6=(0.5_r8+sign(0.5_r8,-dzde_p(i+1,j,k2)))* &
1323 & tl_dzde_p(i+1,j,k2)
1324 tl_cff7=(0.5_r8+sign(0.5_r8, dzde_p(i ,j,k2)))* &
1325 & tl_dzde_p(i ,j,k2)
1326 tl_cff8=(0.5_r8+sign(0.5_r8, dzde_p(i+1,j,k1)))* &
1327 & tl_dzde_p(i+1,j,k1)
1328 vfse(i,j,k2)=vfse(i,j,k2)+ &
1329 & fac2* &
1330 & (cff1*(cff5*dmudz-dmude(i ,j,k1))+ &
1331 & cff2*(cff6*dmudz-dmude(i+1,j,k2))+ &
1332 & cff3*(cff7*dmudz-dmude(i ,j,k2))+ &
1333 & cff4*(cff8*dmudz-dmude(i+1,j,k1)))
1334 tl_vfse(i,j,k2)=tl_vfse(i,j,k2)+ &
1335 & fac2* &
1336 & (tl_cff1*(cff5*dmudz-dmude(i ,j,k1))+ &
1337 & tl_cff2*(cff6*dmudz-dmude(i+1,j,k2))+ &
1338 & tl_cff3*(cff7*dmudz-dmude(i ,j,k2))+ &
1339 & tl_cff4*(cff8*dmudz-dmude(i+1,j,k1))+ &
1340 & cff1*(tl_cff5*dmudz+cff5*tl_dmudz- &
1341 & tl_dmude(i ,j,k1))+ &
1342 & cff2*(tl_cff6*dmudz+cff6*tl_dmudz- &
1343 & tl_dmude(i+1,j,k2))+ &
1344 & cff3*(tl_cff7*dmudz+cff7*tl_dmudz- &
1345 & tl_dmude(i ,j,k2))+ &
1346 & cff4*(tl_cff8*dmudz+cff8*tl_dmudz- &
1347 & tl_dmude(i+1,j,k1)))- &
1348#ifdef TL_IOMS
1349 & fac2* &
1350 & (cff1*(2.0_r8*cff5*dmudz- &
1351 & dmude(i ,j,k1))+ &
1352 & cff2*(2.0_r8*cff6*dmudz- &
1353 & dmude(i+1,j,k2))+ &
1354 & cff3*(2.0_r8*cff7*dmudz- &
1355 & dmude(i ,j,k2))+ &
1356 & cff4*(2.0_r8*cff8*dmudz- &
1357 & dmude(i+1,j,k1)))
1358#endif
1359#ifdef VISC_3DCOEF
1360 tl_vfse(i,j,k2)=tl_vfse(i,j,k2)+ &
1361 & tl_fac2* &
1362 & (cff1*(cff5*dmudz-dmude(i ,j,k1))+ &
1363 & cff2*(cff6*dmudz-dmude(i+1,j,k2))+ &
1364 & cff3*(cff7*dmudz-dmude(i ,j,k2))+ &
1365 & cff4*(cff8*dmudz-dmude(i+1,j,k1)))
1366#endif
1367 END DO
1368 END DO
1369 END IF
1370!
1371! Compute first harmonic operator (m s^-3/2).
1372!
1373 DO j=jstrm1,jendp1
1374 DO i=istrum1,iendp1
1375 cff=0.125_r8*(pm(i-1,j)+pm(i,j))* &
1376 & (pn(i-1,j)+pn(i,j))
1377 cff1=1.0_r8/(0.5_r8*(hz(i-1,j,k)+hz(i,j,k)))
1378 tl_cff1=-cff1*cff1* &
1379 & (0.5_r8*(tl_hz(i-1,j,k)+tl_hz(i,j,k)))+ &
1380#ifdef TL_IOMS
1381 & 2.0_r8*cff1
1382#endif
1383 lapu(i,j,k)=cff*((pn(i-1,j)+pn(i,j))* &
1384 (ufx(i,j)-ufx(i-1,j))+ &
1385 & (pm(i-1,j)+pm(i,j))* &
1386 & (ufe(i,j+1)-ufe(i,j)))+ &
1387 & cff1*((ufsx(i,j,k2)+ufse(i,j,k2))- &
1388 & (ufsx(i,j,k1)+ufse(i,j,k1)))
1389 tl_lapu(i,j,k)=cff*((pn(i-1,j)+pn(i,j))* &
1390 (tl_ufx(i,j)-tl_ufx(i-1,j))+ &
1391 & (pm(i-1,j)+pm(i,j))* &
1392 & (tl_ufe(i,j+1)-tl_ufe(i,j)))+ &
1393 & tl_cff1*((ufsx(i,j,k2)+ufse(i,j,k2))- &
1394 & (ufsx(i,j,k1)+ufse(i,j,k1)))+ &
1395 & cff1*((tl_ufsx(i,j,k2)+tl_ufse(i,j,k2))- &
1396 & (tl_ufsx(i,j,k1)+tl_ufse(i,j,k1)))- &
1397#ifdef TL_IOMS
1398 & cff1*((ufsx(i,j,k2)+ufse(i,j,k2))- &
1399 & (ufsx(i,j,k1)+ufse(i,j,k1)))
1400#endif
1401#ifdef MASKING
1402 lapu(i,j,k)=lapu(i,j,k)*umask(i,j)
1403 tl_lapu(i,j,k)=tl_lapu(i,j,k)*umask(i,j)
1404#endif
1405 END DO
1406 END DO
1407
1408 DO j=jstrvm1,jendp1
1409 DO i=istrm1,iendp1
1410 cff=0.125_r8*(pm(i,j)+pm(i,j-1))* &
1411 & (pn(i,j)+pn(i,j-1))
1412 cff1=1.0_r8/(0.5_r8*(hz(i,j-1,k)+hz(i,j,k)))
1413 tl_cff1=-cff1*cff1* &
1414 & (0.5_r8*(tl_hz(i,j-1,k)+tl_hz(i,j,k)))+ &
1415#ifdef TL_IOMS
1416 & 2.0_r8*cff1
1417#endif
1418 lapv(i,j,k)=cff*((pn(i,j-1)+pn(i,j))* &
1419 & (vfx(i+1,j)-vfx(i,j))- &
1420 & (pm(i,j-1)+pm(i,j))* &
1421 & (vfe(i,j)-vfe(i,j-1)))+ &
1422 & cff1*((vfsx(i,j,k2)+vfse(i,j,k2))- &
1423 & (vfsx(i,j,k1)+vfse(i,j,k1)))
1424 tl_lapv(i,j,k)=cff*((pn(i,j-1)+pn(i,j))* &
1425 & (tl_vfx(i+1,j)-tl_vfx(i,j))- &
1426 & (pm(i,j-1)+pm(i,j))* &
1427 & (tl_vfe(i,j)-tl_vfe(i,j-1)))+ &
1428 & tl_cff1*((vfsx(i,j,k2)+vfse(i,j,k2))- &
1429 & (vfsx(i,j,k1)+vfse(i,j,k1)))+ &
1430 & cff1*((tl_vfsx(i,j,k2)+tl_vfse(i,j,k2))- &
1431 & (tl_vfsx(i,j,k1)+tl_vfse(i,j,k1)))- &
1432#ifdef TL_IOMS
1433 & cff1*((vfsx(i,j,k2)+vfse(i,j,k2))- &
1434 & (vfsx(i,j,k1)+vfse(i,j,k1)))
1435#endif
1436#ifdef MASKING
1437 lapv(i,j,k)=lapv(i,j,k)*vmask(i,j)
1438 tl_lapv(i,j,k)=tl_lapv(i,j,k)*vmask(i,j)
1439#endif
1440 END DO
1441 END DO
1442 END IF
1443 END DO k_loop1
1444!
1445! Apply boundary conditions (closed or gradient; except periodic)
1446! to the first harmonic operator.
1447!
1448 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
1449 IF (domain(ng)%Western_Edge(tile)) THEN
1450 IF (tl_lbc(iwest,isuvel,ng)%closed) THEN
1451 DO k=1,n(ng)
1452 DO j=jstrm1,jendp1
1453 lapu(istru-1,j,k)=0.0_r8
1454 tl_lapu(istru-1,j,k)=0.0_r8
1455 END DO
1456 END DO
1457 ELSE
1458 DO k=1,n(ng)
1459 DO j=jstrm1,jendp1
1460 lapu(istru-1,j,k)=lapu(istru,j,k)
1461 tl_lapu(istru-1,j,k)=tl_lapu(istru,j,k)
1462 END DO
1463 END DO
1464 END IF
1465 IF (tl_lbc(iwest,isvvel,ng)%closed) THEN
1466 DO k=1,n(ng)
1467 DO j=jstrvm1,jendp1
1468 lapv(istr-1,j,k)=gamma2(ng)*lapv(istr,j,k)
1469 tl_lapv(istr-1,j,k)=gamma2(ng)*tl_lapv(istr,j,k)
1470 END DO
1471 END DO
1472 ELSE
1473 DO k=1,n(ng)
1474 DO j=jstrvm1,jendp1
1475 lapv(istr-1,j,k)=0.0_r8
1476 tl_lapv(istr-1,j,k)=0.0_r8
1477 END DO
1478 END DO
1479 END IF
1480 END IF
1481 END IF
1482!
1483 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
1484 IF (domain(ng)%Eastern_Edge(tile)) THEN
1485 IF (tl_lbc(ieast,isuvel,ng)%closed) THEN
1486 DO k=1,n(ng)
1487 DO j=jstrm1,jendp1
1488 lapu(iend+1,j,k)=0.0_r8
1489 tl_lapu(iend+1,j,k)=0.0_r8
1490 END DO
1491 END DO
1492 ELSE
1493 DO k=1,n(ng)
1494 DO j=jstrm1,jendp1
1495 lapu(iend+1,j,k)=lapu(iend,j,k)
1496 tl_lapu(iend+1,j,k)=tl_lapu(iend,j,k)
1497 END DO
1498 END DO
1499 END IF
1500 IF (tl_lbc(ieast,isvvel,ng)%closed) THEN
1501 DO k=1,n(ng)
1502 DO j=jstrvm1,jendp1
1503 lapv(iend+1,j,k)=gamma2(ng)*lapv(iend,j,k)
1504 tl_lapv(iend+1,j,k)=gamma2(ng)*tl_lapv(iend,j,k)
1505 END DO
1506 END DO
1507 ELSE
1508 DO k=1,n(ng)
1509 DO j=jstrvm1,jendp1
1510 lapv(iend+1,j,k)=0.0_r8
1511 tl_lapv(iend+1,j,k)=0.0_r8
1512 END DO
1513 END DO
1514 END IF
1515 END IF
1516 END IF
1517!
1518 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
1519 IF (domain(ng)%Southern_Edge(tile)) THEN
1520 IF (tl_lbc(isouth,isuvel,ng)%closed) THEN
1521 DO k=1,n(ng)
1522 DO i=istrum1,iendp1
1523 lapu(i,jstr-1,k)=gamma2(ng)*lapu(i,jstr,k)
1524 tl_lapu(i,jstr-1,k)=gamma2(ng)*tl_lapu(i,jstr,k)
1525 END DO
1526 END DO
1527 ELSE
1528 DO k=1,n(ng)
1529 DO i=istrum1,iendp1
1530 lapu(i,jstr-1,k)=0.0_r8
1531 tl_lapu(i,jstr-1,k)=0.0_r8
1532 END DO
1533 END DO
1534 END IF
1535 IF (tl_lbc(isouth,isvvel,ng)%closed) THEN
1536 DO k=1,n(ng)
1537 DO i=istrm1,iendp1
1538 lapv(i,jstrv-1,k)=0.0_r8
1539 tl_lapv(i,jstrv-1,k)=0.0_r8
1540 END DO
1541 END DO
1542 ELSE
1543 DO k=1,n(ng)
1544 DO i=istrm1,iendp1
1545 lapv(i,jstrv-1,k)=lapv(i,jstrv,k)
1546 tl_lapv(i,jstrv-1,k)=tl_lapv(i,jstrv,k)
1547 END DO
1548 END DO
1549 END IF
1550 END IF
1551 END IF
1552!
1553 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
1554 IF (domain(ng)%Northern_Edge(tile)) THEN
1555 IF (tl_lbc(inorth,isuvel,ng)%closed) THEN
1556 DO k=1,n(ng)
1557 DO i=istrum1,iendp1
1558 lapu(i,jend+1,k)=gamma2(ng)*lapu(i,jend,k)
1559 tl_lapu(i,jend+1,k)=gamma2(ng)*tl_lapu(i,jend,k)
1560 END DO
1561 END DO
1562 ELSE
1563 DO k=1,n(ng)
1564 DO i=istrum1,iendp1
1565 lapu(i,jend+1,k)=0.0_r8
1566 tl_lapu(i,jend+1,k)=0.0_r8
1567 END DO
1568 END DO
1569 END IF
1570 IF (tl_lbc(inorth,isvvel,ng)%closed) THEN
1571 DO k=1,n(ng)
1572 DO i=istrm1,iendp1
1573 lapv(i,jend+1,k)=0.0_r8
1574 tl_lapv(i,jend+1,k)=0.0_r8
1575 END DO
1576 END DO
1577 ELSE
1578 DO k=1,n(ng)
1579 DO i=istrm1,iendp1
1580 lapv(i,jend+1,k)=lapv(i,jend,k)
1581 tl_lapv(i,jend+1,k)=tl_lapv(i,jend,k)
1582 END DO
1583 END DO
1584 END IF
1585 END IF
1586 END IF
1587!
1588 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng).or. &
1589 & compositegrid(iwest ,ng).or.ewperiodic(ng))) THEN
1590 IF (domain(ng)%SouthWest_Corner(tile)) THEN
1591 DO k=1,n(ng)
1592 lapu(istr ,jstr-1,k)=0.5_r8* &
1593 & (lapu(istr+1,jstr-1,k)+ &
1594 & lapu(istr ,jstr ,k))
1595 tl_lapu(istr ,jstr-1,k)=0.5_r8* &
1596 & (tl_lapu(istr+1,jstr-1,k)+ &
1597 & tl_lapu(istr ,jstr ,k))
1598 lapv(istr-1,jstr ,k)=0.5_r8* &
1599 & (lapv(istr-1,jstr+1,k)+ &
1600 & lapv(istr ,jstr ,k))
1601 tl_lapv(istr-1,jstr ,k)=0.5_r8* &
1602 & (tl_lapv(istr-1,jstr+1,k)+ &
1603 & tl_lapv(istr ,jstr ,k))
1604 END DO
1605 END IF
1606 END IF
1607
1608 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng).or. &
1609 & compositegrid(ieast ,ng).or.ewperiodic(ng))) THEN
1610 IF (domain(ng)%SouthEast_Corner(tile)) THEN
1611 DO k=1,n(ng)
1612 lapu(iend+1,jstr-1,k)=0.5_r8* &
1613 & (lapu(iend ,jstr-1,k)+ &
1614 & lapu(iend+1,jstr ,k))
1615 tl_lapu(iend+1,jstr-1,k)=0.5_r8* &
1616 & (tl_lapu(iend ,jstr-1,k)+ &
1617 & tl_lapu(iend+1,jstr ,k))
1618 lapv(iend+1,jstr ,k)=0.5_r8* &
1619 & (lapv(iend ,jstr ,k)+ &
1620 & lapv(iend+1,jstr+1,k))
1621 tl_lapv(iend+1,jstr ,k)=0.5_r8* &
1622 & (tl_lapv(iend ,jstr ,k)+ &
1623 & tl_lapv(iend+1,jstr+1,k))
1624 END DO
1625 END IF
1626 END IF
1627
1628 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng).or. &
1629 & compositegrid(iwest ,ng).or.ewperiodic(ng))) THEN
1630 IF (domain(ng)%NorthWest_Corner(tile)) THEN
1631 DO k=1,n(ng)
1632 lapu(istr ,jend+1,k)=0.5_r8* &
1633 & (lapu(istr+1,jend+1,k)+ &
1634 & lapu(istr ,jend ,k))
1635 tl_lapu(istr ,jend+1,k)=0.5_r8* &
1636 & (tl_lapu(istr+1,jend+1,k)+ &
1637 & tl_lapu(istr ,jend ,k))
1638 lapv(istr-1,jend+1,k)=0.5_r8* &
1639 & (lapv(istr ,jend+1,k)+ &
1640 & lapv(istr-1,jend ,k))
1641 tl_lapv(istr-1,jend+1,k)=0.5_r8* &
1642 & (tl_lapv(istr ,jend+1,k)+ &
1643 & tl_lapv(istr-1,jend ,k))
1644 END DO
1645 END IF
1646 END IF
1647
1648 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng).or. &
1649 & compositegrid(ieast ,ng).or.ewperiodic(ng))) THEN
1650 IF (domain(ng)%NorthEast_Corner(tile)) THEN
1651 DO k=1,n(ng)
1652 lapu(iend+1,jend+1,k)=0.5_r8* &
1653 & (lapu(iend ,jend+1,k)+ &
1654 & lapu(iend+1,jend ,k))
1655 tl_lapu(iend+1,jend+1,k)=0.5_r8* &
1656 & (tl_lapu(iend ,jend+1,k)+ &
1657 & tl_lapu(iend+1,jend ,k))
1658 lapv(iend+1,jend+1,k)=0.5_r8* &
1659 & (lapv(iend ,jend+1,k)+ &
1660 & lapv(iend+1,jend ,k))
1661 tl_lapv(iend+1,jend+1,k)=0.5_r8* &
1662 & (tl_lapv(iend ,jend+1,k)+ &
1663 & tl_lapv(iend+1,jend ,k))
1664 END DO
1665 END IF
1666 END IF
1667!
1668! Compute horizontal and vertical gradients associated with the
1669! second rotated harmonic operator.
1670!
1671 k2=1
1672 k_loop2 : DO k=0,n(ng)
1673 k1=k2
1674 k2=3-k1
1675 IF (k.lt.n(ng)) THEN
1676!
1677! Compute slopes (nondimensional) at RHO- and PSI-points.
1678!
1679 DO j=jstr-1,jend+1
1680 DO i=istru-1,iend+1
1681 cff=0.5_r8*(pm(i-1,j)+pm(i,j))
1682#ifdef MASKING
1683 cff=cff*umask(i,j)
1684#endif
1685 ufx(i,j)=cff*(z_r(i ,j,k+1)- &
1686 & z_r(i-1,j,k+1))
1687 tl_ufx(i,j)=cff*(tl_z_r(i ,j,k+1)- &
1688 & tl_z_r(i-1,j,k+1))
1689 END DO
1690 END DO
1691 DO j=jstrv-1,jend+1
1692 DO i=istr-1,iend+1
1693 cff=0.5_r8*(pn(i,j-1)+pn(i,j))
1694#ifdef MASKING
1695 cff=cff*vmask(i,j)
1696#endif
1697 vfe(i,j)=cff*(z_r(i,j ,k+1)- &
1698 & z_r(i,j-1,k+1))
1699 tl_vfe(i,j)=cff*(tl_z_r(i,j ,k+1)- &
1700 & tl_z_r(i,j-1,k+1))
1701 END DO
1702 END DO
1703!
1704 DO j=jstr,jend+1
1705 DO i=istr,iend+1
1706 dzdx_p(i,j,k2)=0.5_r8*(ufx(i,j-1)+ &
1707 & ufx(i,j ))
1708 tl_dzdx_p(i,j,k2)=0.5_r8*(tl_ufx(i,j-1)+ &
1709 & tl_ufx(i,j ))
1710 dzde_p(i,j,k2)=0.5_r8*(vfe(i-1,j)+ &
1711 & vfe(i ,j))
1712 tl_dzde_p(i,j,k2)=0.5_r8*(tl_vfe(i-1,j)+ &
1713 & tl_vfe(i ,j))
1714 END DO
1715 END DO
1716 DO j=jstrv-1,jend
1717 DO i=istru-1,iend
1718 dzdx_r(i,j,k2)=0.5_r8*(ufx(i ,j)+ &
1719 & ufx(i+1,j))
1720 tl_dzdx_r(i,j,k2)=0.5_r8*(tl_ufx(i ,j)+ &
1721 & tl_ufx(i+1,j))
1722 dzde_r(i,j,k2)=0.5_r8*(vfe(i,j )+ &
1723 & vfe(i,j+1))
1724 tl_dzde_r(i,j,k2)=0.5_r8*(tl_vfe(i,j )+ &
1725 & tl_vfe(i,j+1))
1726 END DO
1727 END DO
1728!
1729! Compute momentum horizontal (m^-1 s^-3/2) and vertical (s^-3/2)
1730! gradients.
1731!
1732 DO j=jstrv-1,jend
1733 DO i=istru-1,iend
1734 cff=0.5_r8*pm(i,j)
1735#ifdef MASKING
1736 cff=cff*rmask(i,j)
1737#endif
1738 dnudx(i,j,k2)=cff*((pn(i ,j)+pn(i+1,j))* &
1739 & lapu(i+1,j,k+1)- &
1740 & (pn(i-1,j)+pn(i ,j))* &
1741 & lapu(i ,j,k+1))
1742 tl_dnudx(i,j,k2)=cff*((pn(i ,j)+pn(i+1,j))* &
1743 & tl_lapu(i+1,j,k+1)- &
1744 & (pn(i-1,j)+pn(i ,j))* &
1745 & tl_lapu(i ,j,k+1))
1746 END DO
1747 END DO
1748
1749 DO j=jstr,jend+1
1750 DO i=istr,iend+1
1751 cff=0.125_r8*(pn(i-1,j )+pn(i,j )+ &
1752 & pn(i-1,j-1)+pn(i,j-1))
1753#ifdef MASKING
1754 cff=cff*pmask(i,j)
1755#endif
1756 dmude(i,j,k2)=cff*((pm(i-1,j )+pm(i,j ))* &
1757 & lapu(i,j ,k+1)- &
1758 & (pm(i-1,j-1)+pm(i,j-1))* &
1759 & lapu(i,j-1,k+1))
1760 tl_dmude(i,j,k2)=cff*((pm(i-1,j )+pm(i,j ))* &
1761 & tl_lapu(i,j ,k+1)- &
1762 & (pm(i-1,j-1)+pm(i,j-1))* &
1763 & tl_lapu(i,j-1,k+1))
1764 END DO
1765 END DO
1766
1767 DO j=jstr,jend+1
1768 DO i=istr,iend+1
1769 cff=0.125_r8*(pm(i-1,j )+pm(i,j )+ &
1770 & pm(i-1,j-1)+pm(i,j-1))
1771#ifdef MASKING
1772 cff=cff*pmask(i,j)
1773#endif
1774 dnvdx(i,j,k2)=cff*((pn(i ,j-1)+pn(i ,j))* &
1775 & lapv(i ,j,k+1)- &
1776 & (pn(i-1,j-1)+pn(i-1,j))* &
1777 & lapv(i-1,j,k+1))
1778 tl_dnvdx(i,j,k2)=cff*((pn(i ,j-1)+pn(i ,j))* &
1779 & tl_lapv(i ,j,k+1)- &
1780 & (pn(i-1,j-1)+pn(i-1,j))* &
1781 & tl_lapv(i-1,j,k+1))
1782 END DO
1783 END DO
1784
1785 DO j=jstrv-1,jend
1786 DO i=istru-1,iend
1787 cff=0.5_r8*pn(i,j)
1788#ifdef MASKING
1789 cff=cff*rmask(i,j)
1790#endif
1791 dmvde(i,j,k2)=cff*((pm(i,j )+pm(i,j+1))* &
1792 & lapv(i,j+1,k+1)- &
1793 & (pm(i,j-1)+pm(i,j ))* &
1794 & lapv(i,j ,k+1))
1795 tl_dmvde(i,j,k2)=cff*((pm(i,j )+pm(i,j+1))* &
1796 & tl_lapv(i,j+1,k+1)- &
1797 & (pm(i,j-1)+pm(i,j ))* &
1798 & tl_lapv(i,j ,k+1))
1799 END DO
1800 END DO
1801 END IF
1802
1803 IF ((k.eq.0).or.(k.eq.n(ng))) THEN
1804 DO j=jstr-1,jend+1
1805 DO i=istru-1,iend+1
1806 dudz(i,j,k2)=0.0_r8
1807 tl_dudz(i,j,k2)=0.0_r8
1808 END DO
1809 END DO
1810 DO j=jstrv-1,jend+1
1811 DO i=istr-1,iend+1
1812 dvdz(i,j,k2)=0.0_r8
1813 tl_dvdz(i,j,k2)=0.0_r8
1814 END DO
1815 END DO
1816
1817 DO j=jstr,jend
1818 DO i=istru,iend
1819 ufsx(i,j,k2)=0.0_r8
1820 tl_ufsx(i,j,k2)=0.0_r8
1821 ufse(i,j,k2)=0.0_r8
1822 tl_ufse(i,j,k2)=0.0_r8
1823 END DO
1824 END DO
1825 DO j=jstrv,jend
1826 DO i=istr,iend
1827 vfsx(i,j,k2)=0.0_r8
1828 tl_vfsx(i,j,k2)=0.0_r8
1829 vfse(i,j,k2)=0.0_r8
1830 tl_vfse(i,j,k2)=0.0_r8
1831 END DO
1832 END DO
1833 ELSE
1834 DO j=jstr-1,jend+1
1835 DO i=istru-1,iend+1
1836 cff=1.0_r8/(0.5_r8*(z_r(i-1,j,k+1)- &
1837 & z_r(i-1,j,k )+ &
1838 & z_r(i ,j,k+1)- &
1839 & z_r(i ,j,k )))
1840 tl_cff=-cff*cff*(0.5_r8*(tl_z_r(i-1,j,k+1)- &
1841 & tl_z_r(i-1,j,k )+ &
1842 & tl_z_r(i ,j,k+1)- &
1843 & tl_z_r(i ,j,k )))+ &
1844#ifdef TL_IOMS
1845 & 2.0_r8*cff
1846#endif
1847 dudz(i,j,k2)=cff*(lapu(i,j,k+1)- &
1848 & lapu(i,j,k ))
1849 tl_dudz(i,j,k2)=tl_cff*(lapu(i,j,k+1)- &
1850 & lapu(i,j,k ))+ &
1851 & cff*(tl_lapu(i,j,k+1)- &
1852 & tl_lapu(i,j,k ))- &
1853#ifdef TL_IOMS
1854 & dudz(i,j,k2)
1855#endif
1856 END DO
1857 END DO
1858
1859 DO j=jstrv-1,jend+1
1860 DO i=istr-1,iend+1
1861 cff=1.0_r8/(0.5_r8*(z_r(i,j-1,k+1)- &
1862 & z_r(i,j-1,k )+ &
1863 & z_r(i,j ,k+1)- &
1864 & z_r(i,j ,k )))
1865 tl_cff=-cff*cff*(0.5_r8*(tl_z_r(i,j-1,k+1)- &
1866 & tl_z_r(i,j-1,k )+ &
1867 & tl_z_r(i,j ,k+1)- &
1868 & tl_z_r(i,j ,k )))+ &
1869#ifdef TL_IOMS
1870 & 2.0_r8*cff
1871#endif
1872 dvdz(i,j,k2)=cff*(lapv(i,j,k+1)- &
1873 & lapv(i,j,k ))
1874 tl_dvdz(i,j,k2)=tl_cff*(lapv(i,j,k+1)- &
1875 & lapv(i,j,k ))+ &
1876 & cff*(tl_lapv(i,j,k+1)- &
1877 & tl_lapv(i,j,k ))- &
1878#ifdef TL_IOMS
1879 & dvdz(i,j,k2)
1880#endif
1881 END DO
1882 END DO
1883 END IF
1884!
1885! Compute components of the rotated viscous flux (m5/s2) along
1886! geopotential surfaces in the XI- and ETA-directions.
1887!
1888 IF (k.gt.0) THEN
1889 DO j=jstrv-1,jend
1890 DO i=istru-1,iend
1891 cff1=min(dzdx_r(i,j,k1),0.0_r8)
1892 cff2=max(dzdx_r(i,j,k1),0.0_r8)
1893 cff3=min(dzde_r(i,j,k1),0.0_r8)
1894 cff4=max(dzde_r(i,j,k1),0.0_r8)
1895 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_r(i,j,k1)))* &
1896 & tl_dzdx_r(i,j,k1)
1897 tl_cff2=(0.5_r8+sign(0.5_r8, dzdx_r(i,j,k1)))* &
1898 & tl_dzdx_r(i,j,k1)
1899 tl_cff3=(0.5_r8+sign(0.5_r8,-dzde_r(i,j,k1)))* &
1900 & tl_dzde_r(i,j,k1)
1901 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_r(i,j,k1)))* &
1902 & tl_dzde_r(i,j,k1)
1903#ifdef VISC_3DCOEF
1904 cff=hz(i,j,k)* &
1905 & (on_r(i,j)*(dnudx(i,j,k1)- &
1906 & 0.5_r8*pn(i,j)* &
1907 & (cff1*(dudz(i ,j,k1)+ &
1908 & dudz(i+1,j,k2))+ &
1909 & cff2*(dudz(i ,j,k2)+ &
1910 & dudz(i+1,j,k1))))- &
1911 & om_r(i,j)*(dmvde(i,j,k1)- &
1912 & 0.5_r8*pm(i,j)* &
1913 & (cff3*(dvdz(i,j ,k1)+ &
1914 & dvdz(i,j+1,k2))+ &
1915 & cff4*(dvdz(i,j ,k2)+ &
1916 & dvdz(i,j+1,k1)))))
1917#else
1918!^ cff=Hz(i,j,k)* &
1919!^ & (on_r(i,j)*(dnUdx(i,j,k1)- &
1920!^ & 0.5_r8*pn(i,j)* &
1921!^ & (cff1*(dUdz(i ,j,k1)+ &
1922!^ & dUdz(i+1,j,k2))+ &
1923!^ & cff2*(dUdz(i ,j,k2)+ &
1924!^ & dUdz(i+1,j,k1))))- &
1925!^ & om_r(i,j)*(dmVde(i,j,k1)- &
1926!^ & 0.5_r8*pm(i,j)* &
1927!^ & (cff3*(dVdz(i,j ,k1)+ &
1928!^ & dVdz(i,j+1,k2))+ &
1929!^ & cff4*(dVdz(i,j ,k2)+ &
1930!^ & dVdz(i,j+1,k1)))))
1931!^
1932#endif
1933 tl_cff=tl_hz(i,j,k)* &
1934 & (on_r(i,j)*(dnudx(i,j,k1)- &
1935 & 0.5_r8*pn(i,j)* &
1936 & (cff1*(dudz(i ,j,k1)+ &
1937 & dudz(i+1,j,k2))+ &
1938 & cff2*(dudz(i ,j,k2)+ &
1939 & dudz(i+1,j,k1))))- &
1940 & om_r(i,j)*(dmvde(i,j,k1)- &
1941 & 0.5_r8*pm(i,j)* &
1942 & (cff3*(dvdz(i,j ,k1)+ &
1943 & dvdz(i,j+1,k2))+ &
1944 & cff4*(dvdz(i,j ,k2)+ &
1945 & dvdz(i,j+1,k1)))))+ &
1946 & hz(i,j,k)* &
1947 & (on_r(i,j)*(tl_dnudx(i,j,k1)- &
1948 & 0.5_r8*pn(i,j)* &
1949 & (tl_cff1*(dudz(i ,j,k1)+ &
1950 & dudz(i+1,j,k2))+ &
1951 & cff1*(tl_dudz(i ,j,k1)+ &
1952 & tl_dudz(i+1,j,k2))+ &
1953 & tl_cff2*(dudz(i ,j,k2)+ &
1954 & dudz(i+1,j,k1))+ &
1955 & cff2*(tl_dudz(i ,j,k2)+ &
1956 & tl_dudz(i+1,j,k1))))- &
1957 & om_r(i,j)*(tl_dmvde(i,j,k1)- &
1958 & 0.5_r8*pm(i,j)* &
1959 & (tl_cff3*(dvdz(i,j ,k1)+ &
1960 & dvdz(i,j+1,k2))+ &
1961 & cff3*(tl_dvdz(i,j ,k1)+ &
1962 & tl_dvdz(i,j+1,k2))+ &
1963 & tl_cff4*(dvdz(i,j ,k2)+ &
1964 & dvdz(i,j+1,k1))+ &
1965 & cff4*(tl_dvdz(i,j ,k2)+ &
1966 & tl_dvdz(i,j+1,k1)))))- &
1967#ifdef TL_IOMS
1968 & hz(i,j,k)* &
1969 & (on_r(i,j)*(dnudx(i,j,k1)- &
1970 & pn(i,j)* &
1971 & (cff1*(dudz(i ,j,k1)+ &
1972 & dudz(i+1,j,k2))+ &
1973 & cff2*(dudz(i ,j,k2)+ &
1974 & dudz(i+1,j,k1))))- &
1975 & om_r(i,j)*(dmvde(i,j,k1)- &
1976 & pm(i,j)* &
1977 & (cff3*(dvdz(i,j ,k1)+ &
1978 & dvdz(i,j+1,k2))+ &
1979 & cff4*(dvdz(i,j ,k2)+ &
1980 & dvdz(i,j+1,k1)))))
1981#endif
1982#ifdef MASKING
1983# ifdef VISC_3DCOEF
1984 cff=cff*rmask(i,j)
1985# else
1986!^ cff=cff*rmask(i,j)
1987!^
1988# endif
1989 tl_cff=tl_cff*rmask(i,j)
1990#endif
1991#ifdef VISC_3DCOEF
1992# ifdef UV_U3ADV_SPLIT
1993 ufx(i,j)=on_r(i,j)*on_r(i,j)*uvis3d_r(i,j,k)*cff
1994 tl_ufx(i,j)=on_r(i,j)*on_r(i,j)* &
1995 & (tl_uvis3d_r(i,j,k)*cff+ &
1996 & uvis3d_r(i,j,k)*tl_cff)- &
1997# ifdef TL_IOMS
1998 & ufx(i,j)
1999# endif
2000 vfe(i,j)=om_r(i,j)*om_r(i,j)*vvis3d_r(i,j,k)*cff
2001 tl_vfe(i,j)=om_r(i,j)*om_r(i,j)* &
2002 & (tl_vvis3d_r(i,j,k)*cff+ &
2003 & vvis3d_r(i,j,k)*tl_cff)- &
2004# ifdef TL_IOMS
2005 & vfe(i,j)
2006# endif
2007# else
2008 ufx(i,j)=on_r(i,j)*on_r(i,j)*visc3d_r(i,j,k)*cff
2009 tl_ufx(i,j)=on_r(i,j)*on_r(i,j)* &
2010 & (tl_visc3d_r(i,j,k)*cff+ &
2011 & visc3d_r(i,j,k)*tl_cff)- &
2012# ifdef TL_IOMS
2013 & ufx(i,j)
2014# endif
2015 vfe(i,j)=om_r(i,j)*om_r(i,j)*visc3d_r(i,j,k)*cff
2016 tl_vfe(i,j)=om_r(i,j)*om_r(i,j)* &
2017 & (tl_visc3d_r(i,j,k)*cff+ &
2018 & visc3d_r(i,j,k)*tl_cff)- &
2019# ifdef TL_IOMS
2020 & vfe(i,j)
2021# endif
2022# endif
2023#else
2024!^ UFx(i,j)=on_r(i,j)*on_r(i,j)*visc4_r(i,j)*cff
2025!^
2026 tl_ufx(i,j)=on_r(i,j)*on_r(i,j)*visc4_r(i,j)*tl_cff
2027!^ VFe(i,j)=om_r(i,j)*om_r(i,j)*visc4_r(i,j)*cff
2028!^
2029 tl_vfe(i,j)=om_r(i,j)*om_r(i,j)*visc4_r(i,j)*tl_cff
2030#endif
2031 END DO
2032 END DO
2033
2034 DO j=jstr,jend+1
2035 DO i=istr,iend+1
2036 pm_p=0.25_r8*(pm(i-1,j-1)+pm(i-1,j)+ &
2037 & pm(i ,j-1)+pm(i ,j))
2038 pn_p=0.25_r8*(pn(i-1,j-1)+pn(i-1,j)+ &
2039 & pn(i ,j-1)+pn(i ,j))
2040 cff1=min(dzdx_p(i,j,k1),0.0_r8)
2041 cff2=max(dzdx_p(i,j,k1),0.0_r8)
2042 cff3=min(dzde_p(i,j,k1),0.0_r8)
2043 cff4=max(dzde_p(i,j,k1),0.0_r8)
2044 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_p(i,j,k1)))* &
2045 & tl_dzdx_p(i,j,k1)
2046 tl_cff2=(0.5_r8+sign(0.5_r8, dzdx_p(i,j,k1)))* &
2047 & tl_dzdx_p(i,j,k1)
2048 tl_cff3=(0.5_r8+sign(0.5_r8,-dzde_p(i,j,k1)))* &
2049 & tl_dzde_p(i,j,k1)
2050 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_p(i,j,k1)))* &
2051 & tl_dzde_p(i,j,k1)
2052#ifdef VISC_3DCOEF
2053 cff=0.25_r8* &
2054 & (hz(i-1,j ,k)+hz(i,j ,k)+ &
2055 & hz(i-1,j-1,k)+hz(i,j-1,k))* &
2056 & (on_p(i,j)*(dnvdx(i,j,k1)- &
2057 & 0.5_r8*pn_p* &
2058 & (cff1*(dvdz(i-1,j,k1)+ &
2059 & dvdz(i ,j,k2))+ &
2060 & cff2*(dvdz(i-1,j,k2)+ &
2061 & dvdz(i ,j,k1))))+ &
2062 & om_p(i,j)*(dmude(i,j,k1)- &
2063 & 0.5_r8*pm_p* &
2064 & (cff3*(dudz(i,j-1,k1)+ &
2065 & dudz(i,j ,k2))+ &
2066 & cff4*(dudz(i,j-1,k2)+ &
2067 & dudz(i,j ,k1)))))
2068#else
2069!^ cff=0.25_r8* &
2070!^ & (Hz(i-1,j ,k)+Hz(i,j ,k)+ &
2071!^ & Hz(i-1,j-1,k)+Hz(i,j-1,k))* &
2072!^ & (on_p(i,j)*(dnVdx(i,j,k1)- &
2073!^ & 0.5_r8*pn_p* &
2074!^ & (cff1*(dVdz(i-1,j,k1)+ &
2075!^ & dVdz(i ,j,k2))+ &
2076!^ & cff2*(dVdz(i-1,j,k2)+ &
2077!^ & dVdz(i ,j,k1))))+ &
2078!^ & om_p(i,j)*(dmUde(i,j,k1)- &
2079!^ & 0.5_r8*pm_p* &
2080!^ & (cff3*(dUdz(i,j-1,k1)+ &
2081!^ & dUdz(i,j ,k2))+ &
2082!^ & cff4*(dUdz(i,j-1,k2)+ &
2083!^ & dUdz(i,j ,k1)))))
2084!^
2085#endif
2086 tl_cff=0.25_r8* &
2087 & ((tl_hz(i-1,j ,k)+tl_hz(i,j ,k)+ &
2088 & tl_hz(i-1,j-1,k)+tl_hz(i,j-1,k))* &
2089 & (on_p(i,j)*(dnvdx(i,j,k1)- &
2090 & 0.5_r8*pn_p* &
2091 & (cff1*(dvdz(i-1,j,k1)+ &
2092 & dvdz(i ,j,k2))+ &
2093 & cff2*(dvdz(i-1,j,k2)+ &
2094 & dvdz(i ,j,k1))))+ &
2095 & om_p(i,j)*(dmude(i,j,k1)- &
2096 & 0.5_r8*pm_p* &
2097 & (cff3*(dudz(i,j-1,k1)+ &
2098 & dudz(i,j ,k2))+ &
2099 & cff4*(dudz(i,j-1,k2)+ &
2100 & dudz(i,j ,k1)))))+ &
2101 & (hz(i-1,j ,k)+hz(i,j ,k)+ &
2102 & hz(i-1,j-1,k)+hz(i,j-1,k))* &
2103 & (on_p(i,j)*(tl_dnvdx(i,j,k1)- &
2104 & 0.5_r8*pn_p* &
2105 & (tl_cff1*(dvdz(i-1,j,k1)+ &
2106 & dvdz(i ,j,k2))+ &
2107 & cff1*(tl_dvdz(i-1,j,k1)+ &
2108 & tl_dvdz(i ,j,k2))+ &
2109 & tl_cff2*(dvdz(i-1,j,k2)+ &
2110 & dvdz(i ,j,k1))+ &
2111 & cff2*(tl_dvdz(i-1,j,k2)+ &
2112 & tl_dvdz(i ,j,k1))))+ &
2113 & om_p(i,j)*(tl_dmude(i,j,k1)- &
2114 & 0.5_r8*pm_p* &
2115 & (tl_cff3*(dudz(i,j-1,k1)+ &
2116 & dudz(i,j ,k2))+ &
2117 & cff3*(tl_dudz(i,j-1,k1)+ &
2118 & tl_dudz(i,j ,k2))+ &
2119 & tl_cff4*(dudz(i,j-1,k2)+ &
2120 & dudz(i,j ,k1))+ &
2121 & cff4*(tl_dudz(i,j-1,k2)+ &
2122 & tl_dudz(i,j ,k1))))))- &
2123#ifdef TL_IOMS
2124 & visc4_p(i,j)*0.25_r8* &
2125 & (hz(i-1,j ,k)+hz(i,j ,k)+ &
2126 & hz(i-1,j-1,k)+hz(i,j-1,k))* &
2127 & (on_p(i,j)*(dnvdx(i,j,k1)- &
2128 & pn_p* &
2129 & (cff1*(dvdz(i-1,j,k1)+ &
2130 & dvdz(i ,j,k2))+ &
2131 & cff2*(dvdz(i-1,j,k2)+ &
2132 & dvdz(i ,j,k1))))+ &
2133 & om_p(i,j)*(dmude(i,j,k1)- &
2134 & pm_p* &
2135 & (cff3*(dudz(i,j-1,k1)+ &
2136 & dudz(i,j ,k2))+ &
2137 & cff4*(dudz(i,j-1,k2)+ &
2138 & dudz(i,j ,k1)))))
2139#endif
2140#ifdef MASKING
2141# ifdef VISC_3DCOEF
2142 cff=cff*pmask(i,j)
2143# else
2144!^ cff=cff*pmask(i,j)
2145!^
2146# endif
2147 tl_cff=tl_cff*pmask(i,j)
2148#endif
2149#ifdef VISC_3DCOEF
2150# ifdef UV_U3ADV_SPLIT
2151 uvis_p=0.25_r8* &
2152 & (uvis3d_r(i-1,j-1,k)+uvis3d_r(i-1,j,k)+ &
2153 & uvis3d_r(i ,j-1,k)+uvis3d_r(i ,j,k))
2154 tl_uvis_p=0.25_r8* &
2155 & (tl_uvis3d_r(i-1,j-1,k)+tl_uvis3d_r(i-1,j,k)+ &
2156 & tl_uvis3d_r(i ,j-1,k)+tl_uvis3d_r(i ,j,k))
2157 vvis_p=0.25_r8* &
2158 & (vvis3d_r(i-1,j-1,k)+vvis3d_r(i-1,j,k)+ &
2159 & vvis3d_r(i ,j-1,k)+vvis3d_r(i ,j,k))
2160 tl_vvis_p=0.25_r8* &
2161 & (tl_vvis3d_r(i-1,j-1,k)+tl_vvis3d_r(i-1,j,k)+ &
2162 & tl_vvis3d_r(i ,j-1,k)+tl_vvis3d_r(i ,j,k))
2163 ufe(i,j)=om_p(i,j)*om_p(i,j)*uvis_p*cff
2164 tl_ufe(i,j)=om_p(i,j)*om_p(i,j)* &
2165 & (tl_uvis_p*cff+uvis_p*tl_cff)- &
2166# ifdef TL_IOMS
2167 & ufe(i,j)
2168# endif
2169 vfx(i,j)=on_p(i,j)*on_p(i,j)*vvis_p*cff
2170 tl_vfx(i,j)=on_p(i,j)*on_p(i,j)* &
2171 & (tl_vvis_p*cff+vvis_p*tl_cff)- &
2172# ifdef TL_IOMS
2173 & vfx(i,j)
2174# endif
2175# else
2176 visc_p=0.25_r8* &
2177 & (visc3d_r(i-1,j-1,k)+visc3d_r(i-1,j,k)+ &
2178 & visc3d_r(i ,j-1,k)+visc3d_r(i ,j,k))
2179 tl_visc_p=0.25_r8* &
2180 & (tl_visc3d_r(i-1,j-1,k)+tl_visc3d_r(i-1,j,k)+ &
2181 & tl_visc3d_r(i ,j-1,k)+tl_visc3d_r(i ,j,k))
2182 ufe(i,j)=om_p(i,j)*om_p(i,j)*visc_p*cff
2183 tl_ufe(i,j)=om_p(i,j)*om_p(i,j)* &
2184 & (tl_visc_p*cff+visc_p*tl_cff)- &
2185# ifdef TL_IOMS
2186 & ufe(i,j)
2187# endif
2188 vfx(i,j)=on_p(i,j)*on_p(i,j)*visc_p*cff
2189 tl_vfx(i,j)=on_p(i,j)*on_p(i,j)* &
2190 & (tl_visc_p*cff+visc_p*tl_cff)- &
2191# ifdef TL_IOMS
2192 & vfx(i,j)
2193# endif
2194# endif
2195#else
2196!^ UFe(i,j)=om_p(i,j)*om_p(i,j)*visc4_p(i,j)*cff
2197!^
2198 tl_ufe(i,j)=om_p(i,j)*om_p(i,j)*visc4_p(i,j)*tl_cff
2199!^ VFx(i,j)=on_p(i,j)*on_p(i,j)*visc4_p(i,j)*cff
2200!^
2201 tl_vfx(i,j)=on_p(i,j)*on_p(i,j)*visc4_p(i,j)*tl_cff
2202#endif
2203 END DO
2204 END DO
2205!
2206! Compute vertical flux (m2/s2) due to sloping terrain-following
2207! surfaces.
2208!
2209 IF (k.lt.n(ng)) THEN
2210 DO j=jstr,jend
2211 DO i=istru,iend
2212#ifdef VISC_3DCOEF
2213# ifdef UV_U3ADV_SPLIT
2214 cff=0.125_r8* &
2215 & (uvis3d_r(i-1,j,k )+uvis3d_r(i,j,k )+ &
2216 & uvis3d_r(i-1,j,k+1)+uvis3d_r(i,j,k+1))
2217 tl_cff=0.125_r8* &
2218 & (tl_uvis3d_r(i-1,j,k )+tl_uvis3d_r(i,j,k )+ &
2219 & tl_uvis3d_r(i-1,j,k+1)+tl_uvis3d_r(i,j,k+1))
2220# else
2221 cff=0.125_r8* &
2222 & (visc3d_r(i-1,j,k )+visc3d_r(i,j,k )+ &
2223 & visc3d_r(i-1,j,k+1)+visc3d_r(i,j,k+1))
2224 tl_cff=0.125_r8* &
2225 & (tl_visc3d_r(i-1,j,k )+tl_visc3d_r(i,j,k )+ &
2226 & tl_visc3d_r(i-1,j,k+1)+tl_visc3d_r(i,j,k+1))
2227# endif
2228 fac1=cff*on_u(i,j)
2229 tl_fac1=tl_cff*on_u(i,j)
2230 fac2=cff*om_u(i,j)
2231 tl_fac2=tl_cff*om_u(i,j)
2232#else
2233 cff=0.25_r8*(visc4_r(i-1,j)+visc4_r(i,j))
2234 fac1=cff*on_u(i,j)
2235 fac2=cff*om_u(i,j)
2236#endif
2237 cff=0.5_r8*(pn(i-1,j)+pn(i,j))
2238 dnudz=cff*dudz(i,j,k2)
2239 tl_dnudz=cff*tl_dudz(i,j,k2)
2240 dnvdz=cff*0.25_r8*(dvdz(i-1,j+1,k2)+ &
2241 & dvdz(i ,j+1,k2)+ &
2242 & dvdz(i-1,j ,k2)+ &
2243 & dvdz(i ,j ,k2))
2244 tl_dnvdz=cff*0.25_r8*(tl_dvdz(i-1,j+1,k2)+ &
2245 & tl_dvdz(i ,j+1,k2)+ &
2246 & tl_dvdz(i-1,j ,k2)+ &
2247 & tl_dvdz(i ,j ,k2))
2248 cff=0.5_r8*(pm(i-1,j)+pm(i,j))
2249 dmudz=cff*dudz(i,j,k2)
2250 tl_dmudz=cff*tl_dudz(i,j,k2)
2251 dmvdz=cff*0.25_r8*(dvdz(i-1,j+1,k2)+ &
2252 & dvdz(i ,j+1,k2)+ &
2253 & dvdz(i-1,j ,k2)+ &
2254 & dvdz(i ,j ,k2))
2255 tl_dmvdz=cff*0.25_r8*(tl_dvdz(i-1,j+1,k2)+ &
2256 & tl_dvdz(i ,j+1,k2)+ &
2257 & tl_dvdz(i-1,j ,k2)+ &
2258 & tl_dvdz(i ,j ,k2))
2259
2260 cff1=min(dzdx_r(i-1,j,k1),0.0_r8)
2261 cff2=min(dzdx_r(i ,j,k2),0.0_r8)
2262 cff3=max(dzdx_r(i-1,j,k2),0.0_r8)
2263 cff4=max(dzdx_r(i ,j,k1),0.0_r8)
2264 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_r(i-1,j,k1)))* &
2265 & tl_dzdx_r(i-1,j,k1)
2266 tl_cff2=(0.5_r8+sign(0.5_r8,-dzdx_r(i ,j,k2)))* &
2267 & tl_dzdx_r(i ,j,k2)
2268 tl_cff3=(0.5_r8+sign(0.5_r8, dzdx_r(i-1,j,k2)))* &
2269 & tl_dzdx_r(i-1,j,k2)
2270 tl_cff4=(0.5_r8+sign(0.5_r8, dzdx_r(i ,j,k1)))* &
2271 tl_dzdx_r(i ,j,k1)
2272!^ UFsx(i,j,k2)=fac1* &
2273!^ & (cff1*(cff1*dnUdz-dnUdx(i-1,j,k1))+ &
2274!^ & cff2*(cff2*dnUdz-dnUdx(i ,j,k2))+ &
2275!^ & cff3*(cff3*dnUdz-dnUdx(i-1,j,k2))+ &
2276!^ & cff4*(cff4*dnUdz-dnUdx(i ,j,k1)))
2277!^
2278 tl_ufsx(i,j,k2)=fac1* &
2279 & (tl_cff1*(cff1*dnudz-dnudx(i-1,j,k1))+ &
2280 & tl_cff2*(cff2*dnudz-dnudx(i ,j,k2))+ &
2281 & tl_cff3*(cff3*dnudz-dnudx(i-1,j,k2))+ &
2282 & tl_cff4*(cff4*dnudz-dnudx(i ,j,k1))+ &
2283 & cff1*(tl_cff1*dnudz+cff1*tl_dnudz- &
2284 & tl_dnudx(i-1,j,k1))+ &
2285 & cff2*(tl_cff2*dnudz+cff2*tl_dnudz- &
2286 & tl_dnudx(i ,j,k2))+ &
2287 & cff3*(tl_cff3*dnudz+cff3*tl_dnudz- &
2288 & tl_dnudx(i-1,j,k2))+ &
2289 & cff4*(tl_cff4*dnudz+cff4*tl_dnudz- &
2290 & tl_dnudx(i ,j,k1)))- &
2291#ifdef TL_IOMS
2292 & fac1* &
2293 & (cff1*(2.0_r8*cff1*dnudz- &
2294 & dnudx(i-1,j,k1))+ &
2295 & cff2*(2.0_r8*cff2*dnudz- &
2296 & dnudx(i ,j,k2))+ &
2297 & cff3*(2.0_r8*cff3*dnudz- &
2298 & dnudx(i-1,j,k2))+ &
2299 & cff4*(2.0_r8*cff4*dnudz- &
2300 & dnudx(i ,j,k1)))
2301#endif
2302#ifdef VISC_3DCOEF
2303 tl_ufsx(i,j,k2)=tl_ufsx(i,j,k2)+ &
2304 & tl_fac1* &
2305 & (cff1*(cff1*dnudz-dnudx(i-1,j,k1))+ &
2306 & cff2*(cff2*dnudz-dnudx(i ,j,k2))+ &
2307 & cff3*(cff3*dnudz-dnudx(i-1,j,k2))+ &
2308 & cff4*(cff4*dnudz-dnudx(i ,j,k1)))
2309#endif
2310
2311 cff1=min(dzde_p(i,j ,k1),0.0_r8)
2312 cff2=min(dzde_p(i,j+1,k2),0.0_r8)
2313 cff3=max(dzde_p(i,j ,k2),0.0_r8)
2314 cff4=max(dzde_p(i,j+1,k1),0.0_r8)
2315 tl_cff1=(0.5_r8+sign(0.5_r8,-dzde_p(i,j ,k1)))* &
2316 & tl_dzde_p(i,j ,k1)
2317 tl_cff2=(0.5_r8+sign(0.5_r8,-dzde_p(i,j+1,k2)))* &
2318 & tl_dzde_p(i,j+1,k2)
2319 tl_cff3=(0.5_r8+sign(0.5_r8, dzde_p(i,j ,k2)))* &
2320 & tl_dzde_p(i,j ,k2)
2321 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_p(i,j+1,k1)))* &
2322 tl_dzde_p(i,j+1,k1)
2323!^ UFse(i,j,k2)=fac2* &
2324!^ & (cff1*(cff1*dmUdz-dmUde(i,j ,k1))+ &
2325!^ & cff2*(cff2*dmUdz-dmUde(i,j+1,k2))+ &
2326!^ & cff3*(cff3*dmUdz-dmUde(i,j ,k2))+ &
2327!^ & cff4*(cff4*dmUdz-dmUde(i,j+1,k1)))
2328!^
2329 tl_ufse(i,j,k2)=fac2* &
2330 & (tl_cff1*(cff1*dmudz-dmude(i,j ,k1))+ &
2331 & tl_cff2*(cff2*dmudz-dmude(i,j+1,k2))+ &
2332 & tl_cff3*(cff3*dmudz-dmude(i,j ,k2))+ &
2333 & tl_cff4*(cff4*dmudz-dmude(i,j+1,k1))+ &
2334 & cff1*(tl_cff1*dmudz+cff1*tl_dmudz- &
2335 & tl_dmude(i,j ,k1))+ &
2336 & cff2*(tl_cff2*dmudz+cff2*tl_dmudz- &
2337 & tl_dmude(i,j+1,k2))+ &
2338 & cff3*(tl_cff3*dmudz+cff3*tl_dmudz- &
2339 & tl_dmude(i,j ,k2))+ &
2340 & cff4*(tl_cff4*dmudz+cff4*tl_dmudz- &
2341 & tl_dmude(i,j+1,k1)))- &
2342#ifdef TL_IOMS
2343 & fac2* &
2344 & (cff1*(2.0_r8*cff1*dmudz- &
2345 & dmude(i,j ,k1))+ &
2346 & cff2*(2.0_r8*cff2*dmudz- &
2347 & dmude(i,j+1,k2))+ &
2348 & cff3*(2.0_r8*cff3*dmudz- &
2349 & dmude(i,j ,k2))+ &
2350 & cff4*(2.0_r8*cff4*dmudz- &
2351 & dmude(i,j+1,k1)))
2352#endif
2353#ifdef VISC_3DCOEF
2354 tl_ufse(i,j,k2)=tl_ufse(i,j,k2)+ &
2355 & tl_fac2* &
2356 & (cff1*(cff1*dmudz-dmude(i,j ,k1))+ &
2357 & cff2*(cff2*dmudz-dmude(i,j+1,k2))+ &
2358 & cff3*(cff3*dmudz-dmude(i,j ,k2))+ &
2359 & cff4*(cff4*dmudz-dmude(i,j+1,k1)))
2360#endif
2361
2362 cff1=min(dzde_p(i,j ,k1),0.0_r8)
2363 cff2=min(dzde_p(i,j+1,k2),0.0_r8)
2364 cff3=max(dzde_p(i,j ,k2),0.0_r8)
2365 cff4=max(dzde_p(i,j+1,k1),0.0_r8)
2366 cff5=min(dzdx_p(i,j ,k1),0.0_r8)
2367 cff6=min(dzdx_p(i,j+1,k2),0.0_r8)
2368 cff7=max(dzdx_p(i,j ,k2),0.0_r8)
2369 cff8=max(dzdx_p(i,j+1,k1),0.0_r8)
2370 tl_cff1=(0.5_r8+sign(0.5_r8,-dzde_p(i,j ,k1)))* &
2371 & tl_dzde_p(i,j ,k1)
2372 tl_cff2=(0.5_r8+sign(0.5_r8,-dzde_p(i,j+1,k2)))* &
2373 & tl_dzde_p(i,j+1,k2)
2374 tl_cff3=(0.5_r8+sign(0.5_r8, dzde_p(i,j ,k2)))* &
2375 & tl_dzde_p(i,j ,k2)
2376 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_p(i,j+1,k1)))* &
2377 & tl_dzde_p(i,j+1,k1)
2378 tl_cff5=(0.5_r8+sign(0.5_r8,-dzdx_p(i,j ,k1)))* &
2379 & tl_dzdx_p(i,j ,k1)
2380 tl_cff6=(0.5_r8+sign(0.5_r8,-dzdx_p(i,j+1,k2)))* &
2381 & tl_dzdx_p(i,j+1,k2)
2382 tl_cff7=(0.5_r8+sign(0.5_r8, dzdx_p(i,j ,k2)))* &
2383 & tl_dzdx_p(i,j ,k2)
2384 tl_cff8=(0.5_r8+sign(0.5_r8, dzdx_p(i,j+1,k1)))* &
2385 & tl_dzdx_p(i,j+1,k1)
2386!^ UFsx(i,j,k2)=UFsx(i,j,k2)+ &
2387!^ & fac1* &
2388!^ & (cff1*(cff5*dnVdz-dnVdx(i,j ,k1))+ &
2389!^ & cff2*(cff6*dnVdz-dnVdx(i,j+1,k2))+ &
2390!^ & cff3*(cff7*dnVdz-dnVdx(i,j ,k2))+ &
2391!^ & cff4*(cff8*dnVdz-dnVdx(i,j+1,k1)))
2392!^
2393 tl_ufsx(i,j,k2)=tl_ufsx(i,j,k2)+ &
2394 & fac1* &
2395 & (tl_cff1*(cff5*dnvdz-dnvdx(i,j ,k1))+ &
2396 & tl_cff2*(cff6*dnvdz-dnvdx(i,j+1,k2))+ &
2397 & tl_cff3*(cff7*dnvdz-dnvdx(i,j ,k2))+ &
2398 & tl_cff4*(cff8*dnvdz-dnvdx(i,j+1,k1))+ &
2399 & cff1*(tl_cff5*dnvdz+cff5*tl_dnvdz- &
2400 & tl_dnvdx(i,j ,k1))+ &
2401 & cff2*(tl_cff6*dnvdz+cff6*tl_dnvdz- &
2402 & tl_dnvdx(i,j+1,k2))+ &
2403 & cff3*(tl_cff7*dnvdz+cff7*tl_dnvdz- &
2404 & tl_dnvdx(i,j ,k2))+ &
2405 & cff4*(tl_cff8*dnvdz+cff8*tl_dnvdz- &
2406 & tl_dnvdx(i,j+1,k1)))- &
2407#ifdef TL_IOMS
2408 & fac1* &
2409 & (cff1*(2.0_r8*cff5*dnvdz- &
2410 & dnvdx(i,j ,k1))+ &
2411 & cff2*(2.0_r8*cff6*dnvdz- &
2412 & dnvdx(i,j+1,k2))+ &
2413 & cff3*(2.0_r8*cff7*dnvdz- &
2414 & dnvdx(i,j ,k2))+ &
2415 & cff4*(2.0_r8*cff8*dnvdz- &
2416 & dnvdx(i,j+1,k1)))
2417#endif
2418#ifdef VISC_3DCOEF
2419 tl_ufsx(i,j,k2)=tl_ufsx(i,j,k2)+ &
2420 & tl_fac1* &
2421 & (cff1*(cff5*dnvdz-dnvdx(i,j ,k1))+ &
2422 & cff2*(cff6*dnvdz-dnvdx(i,j+1,k2))+ &
2423 & cff3*(cff7*dnvdz-dnvdx(i,j ,k2))+ &
2424 & cff4*(cff8*dnvdz-dnvdx(i,j+1,k1)))
2425#endif
2426
2427 cff1=min(dzdx_r(i-1,j,k1),0.0_r8)
2428 cff2=min(dzdx_r(i ,j,k2),0.0_r8)
2429 cff3=max(dzdx_r(i-1,j,k2),0.0_r8)
2430 cff4=max(dzdx_r(i ,j,k1),0.0_r8)
2431 cff5=min(dzde_r(i-1,j,k1),0.0_r8)
2432 cff6=min(dzde_r(i ,j,k2),0.0_r8)
2433 cff7=max(dzde_r(i-1,j,k2),0.0_r8)
2434 cff8=max(dzde_r(i ,j,k1),0.0_r8)
2435 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_r(i-1,j,k1)))* &
2436 & tl_dzdx_r(i-1,j,k1)
2437 tl_cff2=(0.5_r8+sign(0.5_r8,-dzdx_r(i ,j,k2)))* &
2438 & tl_dzdx_r(i ,j,k2)
2439 tl_cff3=(0.5_r8+sign(0.5_r8, dzdx_r(i-1,j,k2)))* &
2440 & tl_dzdx_r(i-1,j,k2)
2441 tl_cff4=(0.5_r8+sign(0.5_r8, dzdx_r(i ,j,k1)))* &
2442 & tl_dzdx_r(i ,j,k1)
2443 tl_cff5=(0.5_r8+sign(0.5_r8,-dzde_r(i-1,j,k1)))* &
2444 & tl_dzde_r(i-1,j,k1)
2445 tl_cff6=(0.5_r8+sign(0.5_r8,-dzde_r(i ,j,k2)))* &
2446 & tl_dzde_r(i ,j,k2)
2447 tl_cff7=(0.5_r8+sign(0.5_r8, dzde_r(i-1,j,k2)))* &
2448 & tl_dzde_r(i-1,j,k2)
2449 tl_cff8=(0.5_r8+sign(0.5_r8, dzde_r(i ,j,k1)))* &
2450 & tl_dzde_r(i ,j,k1)
2451!^ UFse(i,j,k2)=UFse(i,j,k2)- &
2452!^ & fac2* &
2453!^ & (cff1*(cff5*dmVdz-dmVde(i-1,j,k1))+ &
2454!^ & cff2*(cff6*dmVdz-dmVde(i ,j,k2))+ &
2455!^ & cff3*(cff7*dmVdz-dmVde(i-1,j,k2))+ &
2456!^ & cff4*(cff8*dmVdz-dmVde(i ,j,k1)))
2457!^
2458 tl_ufse(i,j,k2)=tl_ufse(i,j,k2)- &
2459 & fac2* &
2460 & (tl_cff1*(cff5*dmvdz-dmvde(i-1,j,k1))+ &
2461 & tl_cff2*(cff6*dmvdz-dmvde(i ,j,k2))+ &
2462 & tl_cff3*(cff7*dmvdz-dmvde(i-1,j,k2))+ &
2463 & tl_cff4*(cff8*dmvdz-dmvde(i ,j,k1))+ &
2464 & cff1*(tl_cff5*dmvdz+cff5*tl_dmvdz- &
2465 & tl_dmvde(i-1,j,k1))+ &
2466 & cff2*(tl_cff6*dmvdz+cff6*tl_dmvdz- &
2467 & tl_dmvde(i ,j,k2))+ &
2468 & cff3*(tl_cff7*dmvdz+cff7*tl_dmvdz- &
2469 & tl_dmvde(i-1,j,k2))+ &
2470 & cff4*(tl_cff8*dmvdz+cff8*tl_dmvdz- &
2471 & tl_dmvde(i ,j,k1)))+ &
2472#ifdef TL_IOMS
2473 & fac2* &
2474 & (cff1*(2.0_r8*cff5*dmvdz- &
2475 & dmvde(i-1,j,k1))+ &
2476 & cff2*(2.0_r8*cff6*dmvdz- &
2477 & dmvde(i ,j,k2))+ &
2478 & cff3*(2.0_r8*cff7*dmvdz- &
2479 & dmvde(i-1,j,k2))+ &
2480 & cff4*(2.0_r8*cff8*dmvdz- &
2481 & dmvde(i ,j,k1)))
2482#endif
2483#ifdef VISC_3DCOEF
2484 tl_ufse(i,j,k2)=tl_ufse(i,j,k2)- &
2485 & tl_fac2* &
2486 & (cff1*(cff5*dmvdz-dmvde(i-1,j,k1))+ &
2487 & cff2*(cff6*dmvdz-dmvde(i ,j,k2))+ &
2488 & cff3*(cff7*dmvdz-dmvde(i-1,j,k2))+ &
2489 & cff4*(cff8*dmvdz-dmvde(i ,j,k1)))
2490#endif
2491 END DO
2492 END DO
2493!
2494 DO j=jstrv,jend
2495 DO i=istr,iend
2496#ifdef VISC_3DCOEF
2497# ifdef UV_U3ADV_SPLIT
2498 cff=0.125_r8* &
2499 & (vvis3d_r(i,j-1,k )+vvis3d_r(i,j,k )+ &
2500 & vvis3d_r(i,j-1,k+1)+vvis3d_r(i,j,k+1))
2501 tl_cff=0.125_r8* &
2502 & (tl_vvis3d_r(i,j-1,k )+tl_vvis3d_r(i,j,k )+ &
2503 & tl_vvis3d_r(i,j-1,k+1)+tl_vvis3d_r(i,j,k+1))
2504# else
2505 cff=0.125_r8* &
2506 & (visc3d_r(i,j-1,k )+visc3d_r(i,j,k )+ &
2507 & visc3d_r(i,j-1,k+1)+visc3d_r(i,j,k+1))
2508 tl_cff=0.125_r8* &
2509 & (tl_visc3d_r(i,j-1,k )+tl_visc3d_r(i,j,k )+ &
2510 & tl_visc3d_r(i,j-1,k+1)+tl_visc3d_r(i,j,k+1))
2511# endif
2512 fac1=cff*on_v(i,j)
2513 tl_fac1=tl_cff*on_v(i,j)
2514 fac2=cff*om_v(i,j)
2515 tl_fac2=tl_cff*om_v(i,j)
2516#else
2517 cff=0.25_r8*(visc4_r(i,j-1)+visc4_r(i,j))
2518 fac1=cff*on_v(i,j)
2519 fac2=cff*om_v(i,j)
2520#endif
2521 cff=0.5_r8*(pn(i,j-1)+pn(i,j))
2522 dnudz=cff*0.25_r8*(dudz(i ,j ,k2)+ &
2523 & dudz(i+1,j ,k2)+ &
2524 & dudz(i ,j-1,k2)+ &
2525 & dudz(i+1,j-1,k2))
2526 tl_dnudz=cff*0.25_r8*(tl_dudz(i ,j ,k2)+ &
2527 & tl_dudz(i+1,j ,k2)+ &
2528 & tl_dudz(i ,j-1,k2)+ &
2529 & tl_dudz(i+1,j-1,k2))
2530 dnvdz=cff*dvdz(i,j,k2)
2531 tl_dnvdz=cff*tl_dvdz(i,j,k2)
2532 cff=0.5_r8*(pm(i,j-1)+pm(i,j))
2533 dmudz=cff*0.25_r8*(dudz(i ,j ,k2)+ &
2534 & dudz(i+1,j ,k2)+ &
2535 & dudz(i ,j-1,k2)+ &
2536 & dudz(i+1,j-1,k2))
2537 tl_dmudz=cff*0.25_r8*(tl_dudz(i ,j ,k2)+ &
2538 & tl_dudz(i+1,j ,k2)+ &
2539 & tl_dudz(i ,j-1,k2)+ &
2540 & tl_dudz(i+1,j-1,k2))
2541 dmvdz=cff*dvdz(i,j,k2)
2542 tl_dmvdz=cff*tl_dvdz(i,j,k2)
2543
2544 cff1=min(dzdx_p(i ,j,k1),0.0_r8)
2545 cff2=min(dzdx_p(i+1,j,k2),0.0_r8)
2546 cff3=max(dzdx_p(i ,j,k2),0.0_r8)
2547 cff4=max(dzdx_p(i+1,j,k1),0.0_r8)
2548 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_p(i ,j,k1)))* &
2549 & tl_dzdx_p(i ,j,k1)
2550 tl_cff2=(0.5_r8+sign(0.5_r8,-dzdx_p(i+1,j,k2)))* &
2551 & tl_dzdx_p(i+1,j,k2)
2552 tl_cff3=(0.5_r8+sign(0.5_r8, dzdx_p(i ,j,k2)))* &
2553 & tl_dzdx_p(i ,j,k2)
2554 tl_cff4=(0.5_r8+sign(0.5_r8, dzdx_p(i+1,j,k1)))* &
2555 & tl_dzdx_p(i+1,j,k1)
2556!^ VFsx(i,j,k2)=fac1* &
2557!^ & (cff1*(cff1*dnVdz-dnVdx(i ,j,k1))+ &
2558!^ & cff2*(cff2*dnVdz-dnVdx(i+1,j,k2))+ &
2559!^ & cff3*(cff3*dnVdz-dnVdx(i ,j,k2))+ &
2560!^ & cff4*(cff4*dnVdz-dnVdx(i+1,j,k1)))
2561!^
2562 tl_vfsx(i,j,k2)=fac1* &
2563 & (tl_cff1*(cff1*dnvdz-dnvdx(i ,j,k1))+ &
2564 & tl_cff2*(cff2*dnvdz-dnvdx(i+1,j,k2))+ &
2565 & tl_cff3*(cff3*dnvdz-dnvdx(i ,j,k2))+ &
2566 & tl_cff4*(cff4*dnvdz-dnvdx(i+1,j,k1))+ &
2567 & cff1*(tl_cff1*dnvdz+cff1*tl_dnvdz- &
2568 & tl_dnvdx(i ,j,k1))+ &
2569 & cff2*(tl_cff2*dnvdz+cff2*tl_dnvdz- &
2570 & tl_dnvdx(i+1,j,k2))+ &
2571 & cff3*(tl_cff3*dnvdz+cff3*tl_dnvdz- &
2572 & tl_dnvdx(i ,j,k2))+ &
2573 & cff4*(tl_cff4*dnvdz+cff4*tl_dnvdz- &
2574 & tl_dnvdx(i+1,j,k1)))- &
2575#ifdef TL_IOMS
2576 & fac1* &
2577 & (cff1*(2.0_r8*cff1*dnvdz- &
2578 & dnvdx(i ,j,k1))+ &
2579 & cff2*(2.0_r8*cff2*dnvdz- &
2580 & dnvdx(i+1,j,k2))+ &
2581 & cff3*(2.0_r8*cff3*dnvdz- &
2582 & dnvdx(i ,j,k2))+ &
2583 & cff4*(2.0_r8*cff4*dnvdz- &
2584 & dnvdx(i+1,j,k1)))
2585#endif
2586#ifdef VISC_3DCOEF
2587 tl_vfsx(i,j,k2)=tl_vfsx(i,j,k2)+ &
2588 & tl_fac1* &
2589 & (cff1*(cff1*dnvdz-dnvdx(i ,j,k1))+ &
2590 & cff2*(cff2*dnvdz-dnvdx(i+1,j,k2))+ &
2591 & cff3*(cff3*dnvdz-dnvdx(i ,j,k2))+ &
2592 & cff4*(cff4*dnvdz-dnvdx(i+1,j,k1)))
2593#endif
2594
2595 cff1=min(dzde_r(i,j-1,k1),0.0_r8)
2596 cff2=min(dzde_r(i,j ,k2),0.0_r8)
2597 cff3=max(dzde_r(i,j-1,k2),0.0_r8)
2598 cff4=max(dzde_r(i,j ,k1),0.0_r8)
2599 tl_cff1=(0.5_r8+sign(0.5_r8,-dzde_r(i,j-1,k1)))* &
2600 & tl_dzde_r(i,j-1,k1)
2601 tl_cff2=(0.5_r8+sign(0.5_r8,-dzde_r(i,j ,k2)))* &
2602 & tl_dzde_r(i,j ,k2)
2603 tl_cff3=(0.5_r8+sign(0.5_r8, dzde_r(i,j-1,k2)))* &
2604 & tl_dzde_r(i,j-1,k2)
2605 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_r(i,j ,k1)))* &
2606 & tl_dzde_r(i,j ,k1)
2607!^ VFse(i,j,k2)=fac2* &
2608!^ & (cff1*(cff1*dmVdz-dmVde(i,j-1,k1))+ &
2609!^ & cff2*(cff2*dmVdz-dmVde(i,j ,k2))+ &
2610!^ & cff3*(cff3*dmVdz-dmVde(i,j-1,k2))+ &
2611!^ & cff4*(cff4*dmVdz-dmVde(i,j ,k1)))
2612!^
2613 tl_vfse(i,j,k2)=fac2* &
2614 & (tl_cff1*(cff1*dmvdz-dmvde(i,j-1,k1))+ &
2615 & tl_cff2*(cff2*dmvdz-dmvde(i,j ,k2))+ &
2616 & tl_cff3*(cff3*dmvdz-dmvde(i,j-1,k2))+ &
2617 & tl_cff4*(cff4*dmvdz-dmvde(i,j ,k1))+ &
2618 & cff1*(tl_cff1*dmvdz+cff1*tl_dmvdz- &
2619 & tl_dmvde(i,j-1,k1))+ &
2620 & cff2*(tl_cff2*dmvdz+cff2*tl_dmvdz- &
2621 & tl_dmvde(i,j ,k2))+ &
2622 & cff3*(tl_cff3*dmvdz+cff3*tl_dmvdz- &
2623 & tl_dmvde(i,j-1,k2))+ &
2624 & cff4*(tl_cff4*dmvdz+cff4*tl_dmvdz- &
2625 & tl_dmvde(i,j ,k1)))- &
2626#ifdef TL_IOMS
2627 & fac2* &
2628 & (cff1*(2.0_r8*cff1*dmvdz- &
2629 & dmvde(i,j-1,k1))+ &
2630 & cff2*(2.0_r8*cff2*dmvdz- &
2631 & dmvde(i,j ,k2))+ &
2632 & cff3*(2.0_r8*cff3*dmvdz- &
2633 & dmvde(i,j-1,k2))+ &
2634 & cff4*(2.0_r8*cff4*dmvdz- &
2635 & dmvde(i,j ,k1)))
2636#endif
2637#ifdef VISC_3DCOEF
2638 tl_vfse(i,j,k2)=tl_vfse(i,j,k2)+ &
2639 & tl_fac2* &
2640 & (cff1*(cff1*dmvdz-dmvde(i,j-1,k1))+ &
2641 & cff2*(cff2*dmvdz-dmvde(i,j ,k2))+ &
2642 & cff3*(cff3*dmvdz-dmvde(i,j-1,k2))+ &
2643 & cff4*(cff4*dmvdz-dmvde(i,j ,k1)))
2644#endif
2645
2646 cff1=min(dzde_r(i,j-1,k1),0.0_r8)
2647 cff2=min(dzde_r(i,j ,k2),0.0_r8)
2648 cff3=max(dzde_r(i,j-1,k2),0.0_r8)
2649 cff4=max(dzde_r(i,j ,k1),0.0_r8)
2650 cff5=min(dzdx_r(i,j-1,k1),0.0_r8)
2651 cff6=min(dzdx_r(i,j ,k2),0.0_r8)
2652 cff7=max(dzdx_r(i,j-1,k2),0.0_r8)
2653 cff8=max(dzdx_r(i,j ,k1),0.0_r8)
2654 tl_cff1=(0.5_r8+sign(0.5_r8,-dzde_r(i,j-1,k1)))* &
2655 & tl_dzde_r(i,j-1,k1)
2656 tl_cff2=(0.5_r8+sign(0.5_r8,-dzde_r(i,j ,k2)))* &
2657 & tl_dzde_r(i,j ,k2)
2658 tl_cff3=(0.5_r8+sign(0.5_r8, dzde_r(i,j-1,k2)))* &
2659 & tl_dzde_r(i,j-1,k2)
2660 tl_cff4=(0.5_r8+sign(0.5_r8, dzde_r(i,j ,k1)))* &
2661 & tl_dzde_r(i,j ,k1)
2662 tl_cff5=(0.5_r8+sign(0.5_r8,-dzdx_r(i,j-1,k1)))* &
2663 & tl_dzdx_r(i,j-1,k1)
2664 tl_cff6=(0.5_r8+sign(0.5_r8,-dzdx_r(i,j ,k2)))* &
2665 & tl_dzdx_r(i,j ,k2)
2666 tl_cff7=(0.5_r8+sign(0.5_r8, dzdx_r(i,j-1,k2)))* &
2667 & tl_dzdx_r(i,j-1,k2)
2668 tl_cff8=(0.5_r8+sign(0.5_r8, dzdx_r(i,j ,k1)))* &
2669 & tl_dzdx_r(i,j ,k1)
2670!^ VFsx(i,j,k2)=VFsx(i,j,k2)- &
2671!^ & fac1* &
2672!^ & (cff1*(cff5*dnUdz-dnUdx(i,j-1,k1))+ &
2673!^ & cff2*(cff6*dnUdz-dnUdx(i,j ,k2))+ &
2674!^ & cff3*(cff7*dnUdz-dnUdx(i,j-1,k2))+ &
2675!^ & cff4*(cff8*dnUdz-dnUdx(i,j ,k1)))
2676!^
2677 tl_vfsx(i,j,k2)=tl_vfsx(i,j,k2)- &
2678 & fac1* &
2679 & (tl_cff1*(cff5*dnudz-dnudx(i,j-1,k1))+ &
2680 & tl_cff2*(cff6*dnudz-dnudx(i,j ,k2))+ &
2681 & tl_cff3*(cff7*dnudz-dnudx(i,j-1,k2))+ &
2682 & tl_cff4*(cff8*dnudz-dnudx(i,j ,k1))+ &
2683 & cff1*(tl_cff5*dnudz+cff5*tl_dnudz- &
2684 & tl_dnudx(i,j-1,k1))+ &
2685 & cff2*(tl_cff6*dnudz+cff6*tl_dnudz- &
2686 & tl_dnudx(i,j ,k2))+ &
2687 & cff3*(tl_cff7*dnudz+cff7*tl_dnudz- &
2688 & tl_dnudx(i,j-1,k2))+ &
2689 & cff4*(tl_cff8*dnudz+cff8*tl_dnudz- &
2690 & tl_dnudx(i,j ,k1)))+ &
2691#ifdef TL_IOMS
2692 & fac1* &
2693 & (cff1*(2.0_r8*cff5*dnudz- &
2694 & dnudx(i,j-1,k1))+ &
2695 & cff2*(2.0_r8*cff6*dnudz- &
2696 & dnudx(i,j ,k2))+ &
2697 & cff3*(2.0_r8*cff7*dnudz- &
2698 & dnudx(i,j-1,k2))+ &
2699 & cff4*(2.0_r8*cff8*dnudz- &
2700 & dnudx(i,j ,k1)))
2701#endif
2702#ifdef VISC_3DCOEF
2703 tl_vfsx(i,j,k2)=tl_vfsx(i,j,k2)- &
2704 & tl_fac1* &
2705 & (cff1*(cff5*dnudz-dnudx(i,j-1,k1))+ &
2706 & cff2*(cff6*dnudz-dnudx(i,j ,k2))+ &
2707 & cff3*(cff7*dnudz-dnudx(i,j-1,k2))+ &
2708 & cff4*(cff8*dnudz-dnudx(i,j ,k1)))
2709#endif
2710
2711 cff1=min(dzdx_p(i ,j,k1),0.0_r8)
2712 cff2=min(dzdx_p(i+1,j,k2),0.0_r8)
2713 cff3=max(dzdx_p(i ,j,k2),0.0_r8)
2714 cff4=max(dzdx_p(i+1,j,k1),0.0_r8)
2715 cff5=min(dzde_p(i ,j,k1),0.0_r8)
2716 cff6=min(dzde_p(i+1,j,k2),0.0_r8)
2717 cff7=max(dzde_p(i ,j,k2),0.0_r8)
2718 cff8=max(dzde_p(i+1,j,k1),0.0_r8)
2719 tl_cff1=(0.5_r8+sign(0.5_r8,-dzdx_p(i ,j,k1)))* &
2720 & tl_dzdx_p(i ,j,k1)
2721 tl_cff2=(0.5_r8+sign(0.5_r8,-dzdx_p(i+1,j,k2)))* &
2722 & tl_dzdx_p(i+1,j,k2)
2723 tl_cff3=(0.5_r8+sign(0.5_r8, dzdx_p(i ,j,k2)))* &
2724 & tl_dzdx_p(i ,j,k2)
2725 tl_cff4=(0.5_r8+sign(0.5_r8, dzdx_p(i+1,j,k1)))* &
2726 & tl_dzdx_p(i+1,j,k1)
2727 tl_cff5=(0.5_r8+sign(0.5_r8,-dzde_p(i ,j,k1)))* &
2728 & tl_dzde_p(i ,j,k1)
2729 tl_cff6=(0.5_r8+sign(0.5_r8,-dzde_p(i+1,j,k2)))* &
2730 & tl_dzde_p(i+1,j,k2)
2731 tl_cff7=(0.5_r8+sign(0.5_r8, dzde_p(i ,j,k2)))* &
2732 & tl_dzde_p(i ,j,k2)
2733 tl_cff8=(0.5_r8+sign(0.5_r8, dzde_p(i+1,j,k1)))* &
2734 & tl_dzde_p(i+1,j,k1)
2735!^ VFse(i,j,k2)=VFse(i,j,k2)+ &
2736!^ & fac2* &
2737!^ & (cff1*(cff5*dmUdz-dmUde(i ,j,k1))+ &
2738!^ & cff2*(cff6*dmUdz-dmUde(i+1,j,k2))+ &
2739!^ & cff3*(cff7*dmUdz-dmUde(i ,j,k2))+ &
2740!^ & cff4*(cff8*dmUdz-dmUde(i+1,j,k1)))
2741!^
2742 tl_vfse(i,j,k2)=tl_vfse(i,j,k2)+ &
2743 & fac2* &
2744 & (tl_cff1*(cff5*dmudz-dmude(i ,j,k1))+ &
2745 & tl_cff2*(cff6*dmudz-dmude(i+1,j,k2))+ &
2746 & tl_cff3*(cff7*dmudz-dmude(i ,j,k2))+ &
2747 & tl_cff4*(cff8*dmudz-dmude(i+1,j,k1))+ &
2748 & cff1*(tl_cff5*dmudz+cff5*tl_dmudz- &
2749 & tl_dmude(i ,j,k1))+ &
2750 & cff2*(tl_cff6*dmudz+cff6*tl_dmudz- &
2751 & tl_dmude(i+1,j,k2))+ &
2752 & cff3*(tl_cff7*dmudz+cff7*tl_dmudz- &
2753 & tl_dmude(i ,j,k2))+ &
2754 & cff4*(tl_cff8*dmudz+cff8*tl_dmudz- &
2755 & tl_dmude(i+1,j,k1)))- &
2756#ifdef TL_IOMS
2757 & fac2* &
2758 & (cff1*(2.0_r8*cff5*dmudz- &
2759 & dmude(i ,j,k1))+ &
2760 & cff2*(2.0_r8*cff6*dmudz- &
2761 & dmude(i+1,j,k2))+ &
2762 & cff3*(2.0_r8*cff7*dmudz- &
2763 & dmude(i ,j,k2))+ &
2764 & cff4*(2.0_r8*cff8*dmudz- &
2765 & dmude(i+1,j,k1)))
2766#endif
2767#ifdef VISC_3DCOEF
2768 tl_vfse(i,j,k2)=tl_vfse(i,j,k2)+ &
2769 & tl_fac2* &
2770 & (cff1*(cff5*dmudz-dmude(i ,j,k1))+ &
2771 & cff2*(cff6*dmudz-dmude(i+1,j,k2))+ &
2772 & cff3*(cff7*dmudz-dmude(i ,j,k2))+ &
2773 & cff4*(cff8*dmudz-dmude(i+1,j,k1)))
2774#endif
2775 END DO
2776 END DO
2777 END IF
2778!
2779! Time-step biharmonic, geopotential viscosity term. Notice that
2780! momentum at this stage is HzU and HzV and has m2/s units. Add
2781! contribution for barotropic forcing terms.
2782#ifdef DIAGNOSTICS_UV
2783!! The rotated vertical term cannot be split from the horizontal
2784!! terms because of the 2D/3D momentum coupling.
2785#endif
2786!
2787 DO j=jstr,jend
2788 DO i=istru,iend
2789 cff=dt(ng)*0.25_r8*(pm(i-1,j)+pm(i,j))*(pn(i-1,j)+pn(i,j))
2790!^ cff1=0.5_r8*(pn(i-1,j)+pn(i,j))*(UFx(i,j )-UFx(i-1,j))
2791!^
2792 tl_cff1=0.5_r8*(pn(i-1,j)+pn(i,j))* &
2793 & (tl_ufx(i,j )-tl_ufx(i-1,j))
2794!^ cff2=0.5_r8*(pm(i-1,j)+pm(i,j))*(UFe(i,j+1)-UFe(i ,j))
2795!^
2796 tl_cff2=0.5_r8*(pm(i-1,j)+pm(i,j))* &
2797 & (tl_ufe(i,j+1)-tl_ufe(i ,j))
2798!^ cff3=UFsx(i,j,k2)-UFsx(i,j,k1)
2799!^
2800 tl_cff3=tl_ufsx(i,j,k2)-tl_ufsx(i,j,k1)
2801!^ cff4=UFse(i,j,k2)-UFse(i,j,k1)
2802!^
2803 tl_cff4=tl_ufse(i,j,k2)-tl_ufse(i,j,k1)
2804!^ cff5=cff*(cff1+cff2)
2805!^
2806 tl_cff5=cff*(tl_cff1+tl_cff2)
2807!^ cff6=dt(ng)*(cff3+cff4)
2808!^
2809 tl_cff6=dt(ng)*(tl_cff3+tl_cff4)
2810!^ rufrc(i,j)=rufrc(i,j)-cff1-cff2-cff3-cff4
2811!^
2812 tl_rufrc(i,j)=tl_rufrc(i,j)- &
2813 & tl_cff1-tl_cff2-tl_cff3-tl_cff4
2814!^ u(i,j,k,nnew)=u(i,j,k,nnew)-cff5-cff6
2815!^
2816 tl_u(i,j,k,nnew)=tl_u(i,j,k,nnew)-tl_cff5-tl_cff6
2817#ifdef DIAGNOSTICS_UV
2818!! DiaRUfrc(i,j,3,M2hvis)=DiaRUfrc(i,j,3,M2hvis)-cff1-cff2- &
2819!! & cff3-cff4
2820!! DiaRUfrc(i,j,3,M2xvis)=DiaRUfrc(i,j,3,M2xvis)-cff1-cff3
2821!! DiaRUfrc(i,j,3,M2yvis)=DiaRUfrc(i,j,3,M2yvis)-cff2-cff4
2822!! DiaU3wrk(i,j,k,M3hvis)=-cff5-cff6
2823!! DiaU3wrk(i,j,k,M3xvis)=-cff*cff1-dt(ng)*cff3
2824!! DiaU3wrk(i,j,k,M3yvis)=-cff*cff2-dt(ng)*cff4
2825#endif
2826 END DO
2827 END DO
2828
2829 DO j=jstrv,jend
2830 DO i=istr,iend
2831 cff=dt(ng)*0.25_r8*(pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))
2832!^ cff1=0.5_r8*(pn(i,j-1)+pn(i,j))*(VFx(i+1,j)-VFx(i,j ))
2833!^
2834 tl_cff1=0.5_r8*(pn(i,j-1)+pn(i,j))* &
2835 & (tl_vfx(i+1,j)-tl_vfx(i,j ))
2836!^ cff2=0.5_r8*(pm(i,j-1)+pm(i,j))*(VFe(i ,j)-VFe(i,j-1))
2837!^
2838 tl_cff2=0.5_r8*(pm(i,j-1)+pm(i,j))* &
2839 & (tl_vfe(i ,j)-tl_vfe(i,j-1))
2840!^ cff3=VFsx(i,j,k2)-VFsx(i,j,k1)
2841!^
2842 tl_cff3=tl_vfsx(i,j,k2)-tl_vfsx(i,j,k1)
2843!^ cff4=VFse(i,j,k2)-VFse(i,j,k1)
2844!^
2845 tl_cff4=tl_vfse(i,j,k2)-tl_vfse(i,j,k1)
2846!^ cff5=cff*(cff1-cff2)
2847!^
2848 tl_cff5=cff*(tl_cff1-tl_cff2)
2849!^ cff6=dt(ng)*(cff3+cff4)
2850!^
2851 tl_cff6=dt(ng)*(tl_cff3+tl_cff4)
2852!^ rvfrc(i,j)=rvfrc(i,j)-cff1+cff2-cff3-cff4
2853!^
2854 tl_rvfrc(i,j)=tl_rvfrc(i,j)- &
2855 & tl_cff1+tl_cff2-tl_cff3-tl_cff4
2856!^ v(i,j,k,nnew)=v(i,j,k,nnew)-cff5-cff6
2857!^
2858 tl_v(i,j,k,nnew)=tl_v(i,j,k,nnew)-tl_cff5-tl_cff6
2859#ifdef DIAGNOSTICS_UV
2860!! DiaRVfrc(i,j,3,M2hvis)=DiaRVfrc(i,j,3,M2hvis)-cff1+cff2- &
2861!! & cff3-cff4
2862!! DiaRVfrc(i,j,3,M2xvis)=DiaRVfrc(i,j,3,M2xvis)-cff1-cff3
2863!! DiaRVfrc(i,j,3,M2yvis)=DiaRVfrc(i,j,3,M2yvis)+cff2-cff4
2864!! DiaV3wrk(i,j,k,M3hvis)=-cff5-cff6
2865!! DiaV3wrk(i,j,k,M3xvis)=-cff*cff1-dt(ng)*cff3
2866!! DiaV3wrk(i,j,k,M3yvis)= cff*cff2-dt(ng)*cff4
2867#endif
2868 END DO
2869 END DO
2870 END IF
2871 END DO k_loop2
2872!
2873 RETURN
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_param::n
integer, dimension(:), allocatable n
Definition mod_param.F:479
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_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

References mod_scalars::compositegrid, mod_param::domain, mod_scalars::dt, mod_scalars::ewperiodic, mod_scalars::gamma2, mod_scalars::ieast, mod_scalars::inorth, mod_scalars::isouth, mod_ncparam::isuvel, mod_ncparam::isvvel, mod_scalars::iwest, mod_scalars::nsperiodic, and mod_param::tl_lbc.

Referenced by rp_uv3dmix4().

Here is the caller graph for this function:

◆ rp_uv3dmix4_s_tile()

subroutine rp_uv3dmix4_mod::rp_uv3dmix4_s_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) imins,
integer, intent(in) imaxs,
integer, intent(in) jmins,
integer, intent(in) jmaxs,
integer, intent(in) nrhs,
integer, intent(in) nnew,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) pmask,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) hz,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) tl_hz,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) om_p,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) om_r,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) on_p,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) on_r,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) pm,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) pmon_p,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) pmon_r,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) pn,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) pnom_p,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) pnom_r,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) visc4_p,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) visc4_r,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),2), intent(inout) u,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),2), intent(inout) v,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(inout) tl_rufrc,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(inout) tl_rvfrc,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),2), intent(inout) tl_u,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),2), intent(inout) tl_v )
private

Definition at line 114 of file rp_uv3dmix4_s.h.

132!***********************************************************************
133!
134 USE mod_param
135 USE mod_ncparam
136 USE mod_scalars
137!
138! Imported variable declarations.
139!
140 integer, intent(in) :: ng, tile
141 integer, intent(in) :: LBi, UBi, LBj, UBj
142 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
143 integer, intent(in) :: nrhs, nnew
144
145#ifdef ASSUMED_SHAPE
146# ifdef MASKING
147 real(r8), intent(in) :: pmask(LBi:,LBj:)
148# endif
149 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
150 real(r8), intent(in) :: tl_Hz(LBi:,LBj:,:)
151 real(r8), intent(in) :: om_p(LBi:,LBj:)
152 real(r8), intent(in) :: om_r(LBi:,LBj:)
153 real(r8), intent(in) :: on_p(LBi:,LBj:)
154 real(r8), intent(in) :: on_r(LBi:,LBj:)
155 real(r8), intent(in) :: pm(LBi:,LBj:)
156 real(r8), intent(in) :: pmon_p(LBi:,LBj:)
157 real(r8), intent(in) :: pmon_r(LBi:,LBj:)
158 real(r8), intent(in) :: pn(LBi:,LBj:)
159 real(r8), intent(in) :: pnom_p(LBi:,LBj:)
160 real(r8), intent(in) :: pnom_r(LBi:,LBj:)
161 real(r8), intent(in) :: visc4_p(LBi:,LBj:)
162 real(r8), intent(in) :: visc4_r(LBi:,LBj:)
163
164 real(r8), intent(in) :: u(LBi:,LBj:,:,:)
165 real(r8), intent(in) :: v(LBi:,LBj:,:,:)
166
167# ifdef DIAGNOSTICS_UV
168!! real(r8), intent(inout) :: DiaRUfrc(LBi:,LBj:,:,:)
169!! real(r8), intent(inout) :: DiaRVfrc(LBi:,LBj:,:,:)
170!! real(r8), intent(inout) :: DiaU3wrk(LBi:,LBj:,:,:)
171!! real(r8), intent(inout) :: DiaV3wrk(LBi:,LBj:,:,:)
172# endif
173
174 real(r8), intent(inout) :: tl_rufrc(LBi:,LBj:)
175 real(r8), intent(inout) :: tl_rvfrc(LBi:,LBj:)
176 real(r8), intent(inout) :: tl_u(LBi:,LBj:,:,:)
177 real(r8), intent(inout) :: tl_v(LBi:,LBj:,:,:)
178
179#else
180
181# ifdef MASKING
182 real(r8), intent(in) :: pmask(LBi:UBi,LBj:UBj)
183# endif
184 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
185 real(r8), intent(in) :: tl_Hz(LBi:UBi,LBj:UBj,N(ng))
186 real(r8), intent(in) :: om_p(LBi:UBi,LBj:UBj)
187 real(r8), intent(in) :: om_r(LBi:UBi,LBj:UBj)
188 real(r8), intent(in) :: on_p(LBi:UBi,LBj:UBj)
189 real(r8), intent(in) :: on_r(LBi:UBi,LBj:UBj)
190 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
191 real(r8), intent(in) :: pmon_p(LBi:UBi,LBj:UBj)
192 real(r8), intent(in) :: pmon_r(LBi:UBi,LBj:UBj)
193 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
194 real(r8), intent(in) :: pnom_p(LBi:UBi,LBj:UBj)
195 real(r8), intent(in) :: pnom_r(LBi:UBi,LBj:UBj)
196 real(r8), intent(in) :: visc4_p(LBi:UBi,LBj:UBj)
197 real(r8), intent(in) :: visc4_r(LBi:UBi,LBj:UBj)
198
199 real(r8), intent(inout) :: u(LBi:UBi,LBj:UBj,N(ng),2)
200 real(r8), intent(inout) :: v(LBi:UBi,LBj:UBj,N(ng),2)
201
202# ifdef DIAGNOSTICS_UV
203!! real(r8), intent(inout) :: DiaRUfrc(LBi:UBi,LBj:UBj,3,NDM2d-1)
204!! real(r8), intent(inout) :: DiaRVfrc(LBi:UBi,LBj:UBj,3,NDM2d-1)
205!! real(r8), intent(inout) :: DiaU3wrk(LBi:UBi,LBj:UBj,N(ng),NDM3d)
206!! real(r8), intent(inout) :: DiaV3wrk(LBi:UBi,LBj:UBj,N(ng),NDM3d)
207# endif
208
209 real(r8), intent(inout) :: tl_rufrc(LBi:UBi,LBj:UBj)
210 real(r8), intent(inout) :: tl_rvfrc(LBi:UBi,LBj:UBj)
211 real(r8), intent(inout) :: tl_u(LBi:UBi,LBj:UBj,N(ng),2)
212 real(r8), intent(inout) :: tl_v(LBi:UBi,LBj:UBj,N(ng),2)
213#endif
214!
215! Local variable declarations.
216!
217 integer :: IminU, IminV, ImaxU, ImaxV
218 integer :: JminU, JminV, JmaxU, JmaxV
219 integer :: i, j, k
220
221 real(r8) :: cff, cff1, cff2
222 real(r8) :: tl_cff, tl_cff1, tl_cff2
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) :: tl_LapU
232 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_LapV
233 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_UFe
234 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_VFe
235 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_UFx
236 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_VFx
237
238#include "set_bounds.h"
239!
240!-----------------------------------------------------------------------
241! Compute horizontal biharmonic viscosity along constant S-surfaces.
242! The biharmonic operator is computed by applying the harmonic
243! operator twice.
244!-----------------------------------------------------------------------
245!
246! Set local I- and J-ranges.
247!
248 IF (ewperiodic(ng)) THEN
249 iminu=istr-1
250 imaxu=iend+1
251 iminv=istr-1
252 imaxv=iend+1
253 ELSE
254 iminu=max(2,istru-1)
255 imaxu=min(iend+1,lm(ng))
256 iminv=max(1,istr-1)
257 imaxv=min(iend+1,lm(ng))
258 END IF
259 IF (nsperiodic(ng)) THEN
260 jminu=jstr-1
261 jmaxu=jend+1
262 jminv=jstr-1
263 jmaxv=jend+1
264 ELSE
265 jminu=max(1,jstr-1)
266 jmaxu=min(jend+1,mm(ng))
267 jminv=max(2,jstrv-1)
268 jmaxv=min(jend+1,mm(ng))
269 END IF
270!
271! Compute flux-components of the horizontal divergence of the stress
272! tensor (m4 s^-3/2) in XI- and ETA-directions. It is assumed here
273! that mixing coefficients are the squared root of the biharmonic
274! viscosity coefficient. For momentum balance purposes, the
275! thickness "Hz" appears only when computing the second harmonic
276! operator.
277!
278 k_loop : DO k=1,n(ng)
279 DO j=jminv-1,jmaxv
280 DO i=iminu-1,imaxu
281 cff=visc4_r(i,j)*0.5_r8* &
282 & (pmon_r(i,j)* &
283 & ((pn(i ,j)+pn(i+1,j))*u(i+1,j,k,nrhs)- &
284 & (pn(i-1,j)+pn(i ,j))*u(i ,j,k,nrhs))- &
285 & pnom_r(i,j)* &
286 & ((pm(i,j )+pm(i,j+1))*v(i,j+1,k,nrhs)- &
287 & (pm(i,j-1)+pm(i,j ))*v(i,j ,k,nrhs)))
288 tl_cff=visc4_r(i,j)*0.5_r8* &
289 & (pmon_r(i,j)* &
290 & ((pn(i ,j)+pn(i+1,j))*tl_u(i+1,j,k,nrhs)- &
291 & (pn(i-1,j)+pn(i ,j))*tl_u(i ,j,k,nrhs))- &
292 & pnom_r(i,j)* &
293 & ((pm(i,j )+pm(i,j+1))*tl_v(i,j+1,k,nrhs)- &
294 & (pm(i,j-1)+pm(i,j ))*tl_v(i,j ,k,nrhs)))
295 ufx(i,j)=on_r(i,j)*on_r(i,j)*cff
296 tl_ufx(i,j)=on_r(i,j)*on_r(i,j)*tl_cff
297 vfe(i,j)=om_r(i,j)*om_r(i,j)*cff
298 tl_vfe(i,j)=om_r(i,j)*om_r(i,j)*tl_cff
299 END DO
300 END DO
301 DO j=jminu,jmaxu+1
302 DO i=iminv,imaxv+1
303 cff=visc4_p(i,j)*0.5_r8* &
304 & (pmon_p(i,j)* &
305 & ((pn(i ,j-1)+pn(i ,j))*v(i ,j,k,nrhs)- &
306 & (pn(i-1,j-1)+pn(i-1,j))*v(i-1,j,k,nrhs))+ &
307 & pnom_p(i,j)* &
308 & ((pm(i-1,j )+pm(i,j ))*u(i,j ,k,nrhs)- &
309 & (pm(i-1,j-1)+pm(i,j-1))*u(i,j-1,k,nrhs)))
310 tl_cff=visc4_p(i,j)*0.5_r8* &
311 & (pmon_p(i,j)* &
312 & ((pn(i ,j-1)+pn(i ,j))*tl_v(i ,j,k,nrhs)- &
313 & (pn(i-1,j-1)+pn(i-1,j))*tl_v(i-1,j,k,nrhs))+ &
314 & pnom_p(i,j)* &
315 & ((pm(i-1,j )+pm(i,j ))*tl_u(i,j ,k,nrhs)- &
316 & (pm(i-1,j-1)+pm(i,j-1))*tl_u(i,j-1,k,nrhs)))
317#ifdef MASKING
318 cff=cff*pmask(i,j)
319 tl_cff=tl_cff*pmask(i,j)
320#endif
321 ufe(i,j)=om_p(i,j)*om_p(i,j)*cff
322 tl_ufe(i,j)=om_p(i,j)*om_p(i,j)*tl_cff
323 vfx(i,j)=on_p(i,j)*on_p(i,j)*cff
324 tl_vfx(i,j)=on_p(i,j)*on_p(i,j)*tl_cff
325 END DO
326 END DO
327!
328! Compute first harmonic operator (m s^-3/2).
329!
330 DO j=jminu,jmaxu
331 DO i=iminu,imaxu
332 lapu(i,j)=0.125_r8* &
333 & (pm(i-1,j)+pm(i,j))*(pn(i-1,j)+pn(i,j))* &
334 & ((pn(i-1,j)+pn(i,j))* &
335 & (ufx(i,j )-ufx(i-1,j))+ &
336 & (pm(i-1,j)+pm(i,j))* &
337 & (ufe(i,j+1)-ufe(i ,j)))
338 tl_lapu(i,j)=0.125_r8* &
339 & (pm(i-1,j)+pm(i,j))*(pn(i-1,j)+pn(i,j))* &
340 & ((pn(i-1,j)+pn(i,j))* &
341 & (tl_ufx(i,j )-tl_ufx(i-1,j))+ &
342 & (pm(i-1,j)+pm(i,j))* &
343 & (tl_ufe(i,j+1)-tl_ufe(i ,j)))
344 END DO
345 END DO
346 DO j=jminv,jmaxv
347 DO i=iminv,imaxv
348 lapv(i,j)=0.125_r8* &
349 & (pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))* &
350 & ((pn(i,j-1)+pn(i,j))* &
351 & (vfx(i+1,j)-vfx(i,j ))- &
352 & (pm(i,j-1)+pm(i,j))* &
353 & (vfe(i ,j)-vfe(i,j-1)))
354 tl_lapv(i,j)=0.125_r8* &
355 & (pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))* &
356 & ((pn(i,j-1)+pn(i,j))* &
357 & (tl_vfx(i+1,j)-tl_vfx(i,j ))- &
358 & (pm(i,j-1)+pm(i,j))* &
359 & (tl_vfe(i ,j)-tl_vfe(i,j-1)))
360 END DO
361 END DO
362!
363! Apply boundary conditions (other than periodic) to the first
364! harmonic operator. These are gradient or closed (free slip or
365! no slip) boundary conditions.
366!
367 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
368 IF (domain(ng)%Western_Edge(tile)) THEN
369 IF (tl_lbc(iwest,isuvel,ng)%closed) THEN
370 DO j=jminu,jmaxu
371 lapu(istr,j)=0.0_r8
372 tl_lapu(istr,j)=0.0_r8
373 END DO
374 ELSE
375 DO j=jminu,jmaxu
376 lapu(istr,j)=lapu(istr+1,j)
377 tl_lapu(istr,j)=tl_lapu(istr+1,j)
378 END DO
379 END IF
380 IF (tl_lbc(iwest,isvvel,ng)%closed) THEN
381 DO j=jminv,jmaxv
382 lapv(istr-1,j)=gamma2(ng)*lapv(istr,j)
383 tl_lapv(istr-1,j)=gamma2(ng)*tl_lapv(istr,j)
384 END DO
385 ELSE
386 DO j=jminv,jmaxv
387 lapv(istr-1,j)=0.0_r8
388 tl_lapv(istr-1,j)=0.0_r8
389 END DO
390 END IF
391 END IF
392 END IF
393!
394 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
395 IF (domain(ng)%Eastern_Edge(tile)) THEN
396 IF (tl_lbc(ieast,isuvel,ng)%closed) THEN
397 DO j=jminu,jmaxu
398 lapu(iend+1,j)=0.0_r8
399 tl_lapu(iend+1,j)=0.0_r8
400 END DO
401 ELSE
402 DO j=jminu,jmaxu
403 lapu(iend+1,j)=lapu(iend,j)
404 tl_lapu(iend+1,j)=tl_lapu(iend,j)
405 END DO
406 END IF
407 IF (tl_lbc(ieast,isvvel,ng)%closed) THEN
408 DO j=jminv,jmaxv
409 lapv(iend+1,j)=gamma2(ng)*lapv(iend,j)
410 tl_lapv(iend+1,j)=gamma2(ng)*tl_lapv(iend,j)
411 END DO
412 ELSE
413 DO j=jminv,jmaxv
414 lapv(iend+1,j)=0.0_r8
415 tl_lapv(iend+1,j)=0.0_r8
416 END DO
417 END IF
418 END IF
419 END IF
420!
421 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
422 IF (domain(ng)%Southern_Edge(tile)) THEN
423 IF (tl_lbc(isouth,isuvel,ng)%closed) THEN
424 DO i=iminu,imaxu
425 lapu(i,jstr-1)=gamma2(ng)*lapu(i,jstr)
426 tl_lapu(i,jstr-1)=gamma2(ng)*tl_lapu(i,jstr)
427 END DO
428 ELSE
429 DO i=iminu,imaxu
430 lapu(i,jstr-1)=0.0_r8
431 tl_lapu(i,jstr-1)=0.0_r8
432 END DO
433 END IF
434 IF (tl_lbc(isouth,isvvel,ng)%closed) THEN
435 DO i=iminv,imaxv
436 lapv(i,jstr)=0.0_r8
437 tl_lapv(i,jstr)=0.0_r8
438 END DO
439 ELSE
440 DO i=iminv,imaxv
441 lapv(i,jstr)=lapv(i,jstr+1)
442 tl_lapv(i,jstr)=tl_lapv(i,jstr+1)
443 END DO
444 END IF
445 END IF
446 END IF
447!
448 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
449 IF (domain(ng)%Northern_Edge(tile)) THEN
450 IF (tl_lbc(inorth,isuvel,ng)%closed) THEN
451 DO i=iminu,imaxu
452 lapu(i,jend+1)=gamma2(ng)*lapu(i,jend)
453 tl_lapu(i,jend+1)=gamma2(ng)*tl_lapu(i,jend)
454 END DO
455 ELSE
456 DO i=iminu,imaxu
457 lapu(i,jend+1)=0.0_r8
458 tl_lapu(i,jend+1)=0.0_r8
459 END DO
460 END IF
461 IF (tl_lbc(inorth,isvvel,ng)%closed) THEN
462 DO i=iminv,imaxv
463 lapv(i,jend+1)=0.0_r8
464 tl_lapv(i,jend+1)=0.0_r8
465 END DO
466 ELSE
467 DO i=iminv,imaxv
468 lapv(i,jend+1)=lapv(i,jend)
469 tl_lapv(i,jend+1)=tl_lapv(i,jend)
470 END DO
471 END IF
472 END IF
473 END IF
474!
475 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng).or. &
476 & compositegrid(iwest ,ng).or.ewperiodic(ng))) THEN
477 IF (domain(ng)%SouthWest_Corner(tile)) THEN
478 lapu(istr ,jstr-1)=0.5_r8* &
479 & (lapu(istr+1,jstr-1)+ &
480 & lapu(istr ,jstr ))
481 tl_lapu(istr ,jstr-1)=0.5_r8* &
482 & (tl_lapu(istr+1,jstr-1)+ &
483 & tl_lapu(istr ,jstr ))
484 lapv(istr-1,jstr )=0.5_r8* &
485 & (lapv(istr-1,jstr+1)+ &
486 & lapv(istr ,jstr ))
487 tl_lapv(istr-1,jstr )=0.5_r8* &
488 & (tl_lapv(istr-1,jstr+1)+ &
489 & tl_lapv(istr ,jstr ))
490 END IF
491 END IF
492
493 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng).or. &
494 & compositegrid(ieast ,ng).or.ewperiodic(ng))) THEN
495 IF (domain(ng)%SouthEast_Corner(tile)) THEN
496 lapu(iend+1,jstr-1)=0.5_r8* &
497 & (lapu(iend ,jstr-1)+ &
498 & lapu(iend+1,jstr ))
499 tl_lapu(iend+1,jstr-1)=0.5_r8* &
500 & (tl_lapu(iend ,jstr-1)+ &
501 & tl_lapu(iend+1,jstr ))
502 lapv(iend+1,jstr )=0.5_r8* &
503 & (lapv(iend ,jstr )+ &
504 & lapv(iend+1,jstr+1))
505 tl_lapv(iend+1,jstr )=0.5_r8* &
506 & (tl_lapv(iend ,jstr )+ &
507 & tl_lapv(iend+1,jstr+1))
508 END IF
509 END IF
510
511 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng).or. &
512 & compositegrid(iwest ,ng).or.ewperiodic(ng))) THEN
513 IF (domain(ng)%NorthWest_Corner(tile)) THEN
514 lapu(istr ,jend+1)=0.5_r8* &
515 & (lapu(istr+1,jend+1)+ &
516 & lapu(istr ,jend ))
517 tl_lapu(istr ,jend+1)=0.5_r8* &
518 & (tl_lapu(istr+1,jend+1)+ &
519 & tl_lapu(istr ,jend ))
520 lapv(istr-1,jend+1)=0.5_r8* &
521 & (lapv(istr ,jend+1)+ &
522 & lapv(istr-1,jend ))
523 tl_lapv(istr-1,jend+1)=0.5_r8* &
524 & (tl_lapv(istr ,jend+1)+ &
525 & tl_lapv(istr-1,jend ))
526 END IF
527 END IF
528
529 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng).or. &
530 & compositegrid(ieast ,ng).or.ewperiodic(ng))) THEN
531 IF (domain(ng)%NorthEast_Corner(tile)) THEN
532 lapu(iend+1,jend+1)=0.5_r8* &
533 & (lapu(iend ,jend+1)+ &
534 & lapu(iend+1,jend ))
535 tl_lapu(iend+1,jend+1)=0.5_r8* &
536 & (tl_lapu(iend ,jend+1)+ &
537 & tl_lapu(iend+1,jend ))
538 lapv(iend+1,jend+1)=0.5_r8* &
539 & (lapv(iend ,jend+1)+ &
540 & lapv(iend+1,jend ))
541 tl_lapv(iend+1,jend+1)=0.5_r8* &
542 & (tl_lapv(iend ,jend+1)+ &
543 & tl_lapv(iend+1,jend ))
544 END IF
545 END IF
546!
547! Compute flux-components of the horizontal divergence of the
548! harmonic stress tensor (m4/s2) in XI- and ETA-directions.
549!
550 DO j=jstrv-1,jend
551 DO i=istru-1,iend
552 cff=visc4_r(i,j)*hz(i,j,k)*0.5_r8* &
553 & (pmon_r(i,j)* &
554 & ((pn(i ,j)+pn(i+1,j))*lapu(i+1,j)- &
555 & (pn(i-1,j)+pn(i ,j))*lapu(i ,j))- &
556 & pnom_r(i,j)* &
557 & ((pm(i,j )+pm(i,j+1))*lapv(i,j+1)- &
558 & (pm(i,j-1)+pm(i,j ))*lapv(i,j )))
559 tl_cff=visc4_r(i,j)*0.5_r8* &
560 & (tl_hz(i,j,k)* &
561 & (pmon_r(i,j)* &
562 & ((pn(i ,j)+pn(i+1,j))*lapu(i+1,j)- &
563 & (pn(i-1,j)+pn(i ,j))*lapu(i ,j))- &
564 & pnom_r(i,j)* &
565 & ((pm(i,j )+pm(i,j+1))*lapv(i,j+1)- &
566 & (pm(i,j-1)+pm(i,j ))*lapv(i,j )))+ &
567 & hz(i,j,k)* &
568 & (pmon_r(i,j)* &
569 & ((pn(i ,j)+pn(i+1,j))*tl_lapu(i+1,j)- &
570 & (pn(i-1,j)+pn(i ,j))*tl_lapu(i ,j))- &
571 & pnom_r(i,j)* &
572 & ((pm(i,j )+pm(i,j+1))*tl_lapv(i,j+1)- &
573 & (pm(i,j-1)+pm(i,j ))*tl_lapv(i,j ))))- &
574#ifdef TL_IOMS
575 & cff
576#endif
577!^ UFx(i,j)=on_r(i,j)*on_r(i,j)*cff
578!^
579 tl_ufx(i,j)=on_r(i,j)*on_r(i,j)*tl_cff
580!^ VFe(i,j)=om_r(i,j)*om_r(i,j)*cff
581!^
582 tl_vfe(i,j)=om_r(i,j)*om_r(i,j)*tl_cff
583 END DO
584 END DO
585 DO j=jstr,jend+1
586 DO i=istr,iend+1
587 cff=visc4_p(i,j)*0.125_r8*(hz(i-1,j ,k)+hz(i,j ,k)+ &
588 & hz(i-1,j-1,k)+hz(i,j-1,k))* &
589 & (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 tl_cff=visc4_p(i,j)*0.125_r8* &
596 & ((tl_hz(i-1,j ,k)+tl_hz(i,j ,k)+ &
597 & tl_hz(i-1,j-1,k)+tl_hz(i,j-1,k))* &
598 & (pmon_p(i,j)* &
599 & ((pn(i ,j-1)+pn(i ,j))*lapv(i ,j)- &
600 & (pn(i-1,j-1)+pn(i-1,j))*lapv(i-1,j))+ &
601 & pnom_p(i,j)* &
602 & ((pm(i-1,j )+pm(i,j ))*lapu(i,j )- &
603 & (pm(i-1,j-1)+pm(i,j-1))*lapu(i,j-1)))+ &
604 & (hz(i-1,j ,k)+hz(i,j ,k)+ &
605 & hz(i-1,j-1,k)+hz(i,j-1,k))* &
606 & (pmon_p(i,j)* &
607 & ((pn(i ,j-1)+pn(i ,j))*tl_lapv(i ,j)- &
608 & (pn(i-1,j-1)+pn(i-1,j))*tl_lapv(i-1,j))+ &
609 & pnom_p(i,j)* &
610 & ((pm(i-1,j )+pm(i,j ))*tl_lapu(i,j )- &
611 & (pm(i-1,j-1)+pm(i,j-1))*tl_lapu(i,j-1))))- &
612#ifdef TL_IOMS
613 & cff
614#endif
615#ifdef MASKING
616!^ cff=cff*pmask(i,j)
617!^
618 tl_cff=tl_cff*pmask(i,j)
619#endif
620!^ UFe(i,j)=om_p(i,j)*om_p(i,j)*cff
621!^
622 tl_ufe(i,j)=om_p(i,j)*om_p(i,j)*tl_cff
623!^ VFx(i,j)=on_p(i,j)*on_p(i,j)*cff
624!^
625 tl_vfx(i,j)=on_p(i,j)*on_p(i,j)*tl_cff
626 END DO
627 END DO
628!
629! Time-step biharmonic, S-surfaces viscosity term. Notice that
630! momentum at this stage is HzU and HzV and has units m2/s. Add
631! contribution for barotropic forcing terms.
632!
633 DO j=jstr,jend
634 DO i=istru,iend
635 cff=0.25_r8*(pm(i-1,j)+pm(i,j))*(pn(i-1,j)+pn(i,j))
636!^ cff1=0.5_r8*((pn(i-1,j)+pn(i,j))* &
637!^ & (UFx(i,j )-UFx(i-1,j))+ &
638!^ & (pm(i-1,j)+pm(i,j))* &
639!^ & (UFe(i,j+1)-UFe(i ,j)))
640!^
641 tl_cff1=0.5_r8*((pn(i-1,j)+pn(i,j))* &
642 & (tl_ufx(i,j )-tl_ufx(i-1,j))+ &
643 & (pm(i-1,j)+pm(i,j))* &
644 & (tl_ufe(i,j+1)-tl_ufe(i ,j)))
645!^ cff2=dt(ng)*cff*cff1
646!^
647 tl_cff2=dt(ng)*cff*tl_cff1
648!^ rufrc(i,j)=rufrc(i,j)-cff1
649!^
650 tl_rufrc(i,j)=tl_rufrc(i,j)-tl_cff1
651!^ u(i,j,k,nnew)=u(i,j,k,nnew)-cff2
652!^
653 tl_u(i,j,k,nnew)=tl_u(i,j,k,nnew)-tl_cff2
654#ifdef DIAGNOSTICS_UV
655!! DiaRUfrc(i,j,3,M2hvis)=DiaRUfrc(i,j,3,M2hvis)-cff1
656!! DiaU3wrk(i,j,k,M3hvis)=-cff2
657#endif
658 END DO
659 END DO
660 DO j=jstrv,jend
661 DO i=istr,iend
662 cff=0.25_r8*(pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))
663!^ cff1=0.5_r8*((pn(i,j-1)+pn(i,j))* &
664!^ & (VFx(i+1,j)-VFx(i,j ))- &
665!^ & (pm(i,j-1)+pm(i,j))* &
666!^ & (VFe(i ,j)-VFe(i,j-1)))
667!^
668 tl_cff1=0.5_r8*((pn(i,j-1)+pn(i,j))* &
669 & (tl_vfx(i+1,j)-tl_vfx(i,j ))- &
670 & (pm(i,j-1)+pm(i,j))* &
671 & (tl_vfe(i ,j)-tl_vfe(i,j-1)))
672!^ cff2=dt(ng)*cff*cff1
673!^
674 tl_cff2=dt(ng)*cff*tl_cff1
675!^ rvfrc(i,j)=rvfrc(i,j)-cff1
676!^
677 tl_rvfrc(i,j)=tl_rvfrc(i,j)-tl_cff1
678!^ v(i,j,k,nnew)=v(i,j,k,nnew)-cff2
679!^
680 tl_v(i,j,k,nnew)=tl_v(i,j,k,nnew)-tl_cff2
681#ifdef DIAGNOSTICS_UV
682!! DiaRVfrc(i,j,3,M2hvis)=DiaRVfrc(i,j,3,M2hvis)-cff1
683!! DiaV3wrk(i,j,k,M3hvis)=-cff2
684#endif
685 END DO
686 END DO
687 END DO k_loop
688!
689 RETURN
mod_param::lm
integer, dimension(:), allocatable lm
Definition mod_param.F:455
mod_param::mm
integer, dimension(:), allocatable mm
Definition mod_param.F:456

References mod_scalars::compositegrid, mod_param::domain, mod_scalars::dt, mod_scalars::ewperiodic, mod_scalars::gamma2, mod_scalars::ieast, mod_scalars::inorth, mod_scalars::isouth, mod_ncparam::isuvel, mod_ncparam::isvvel, mod_scalars::iwest, mod_param::lm, mod_param::mm, mod_scalars::nsperiodic, and mod_param::tl_lbc.