ROMS
Loading...
Searching...
No Matches
tl_uv3dmix4_s.h
Go to the documentation of this file.
1 MODULE tl_uv3dmix4_mod
2!
3!git $Id$
4!================================================== Hernan G. Arango ===
5! Copyright (c) 2002-2025 The ROMS Group Andrew M. Moore !
6! Licensed under a MIT/X style license !
7! See License_ROMS.md !
8!=======================================================================
9! !
10! This subroutine computes tangent linear biharmonic mixing of !
11! momentum, along constant S-surfaces, from the horizontal !
12! divergence of the stress tensor. A transverse isotropy is !
13! assumed so the stress tensor is split into vertical and !
14! horizontal subtensors. !
15! !
16! Reference: !
17! !
18! Wajsowicz, R.C, 1993: A consistent formulation of the !
19! anisotropic stress tensor for use in models of the !
20! large-scale ocean circulation, JCP, 105, 333-338. !
21! !
22! Sadourny, R. and K. Maynard, 1997: Formulations of !
23! lateral diffusion in geophysical fluid dynamics !
24! models, In Numerical Methods of Atmospheric and !
25! Oceanic Modelling. Lin, Laprise, and Ritchie, !
26! Eds., NRC Research Press, 547-556. !
27! !
28! Griffies, S.M. and R.W. Hallberg, 2000: Biharmonic !
29! friction with a Smagorinsky-like viscosity for !
30! use in large-scale eddy-permitting ocean models, !
31! Monthly Weather Rev., 128, 8, 2935-2946. !
32! !
33! Basic state variables required: visc4, u, v, Hz. !
34! !
35!=======================================================================
36!
37 implicit none
38!
39 PRIVATE
40 PUBLIC tl_uv3dmix4
41!
42 CONTAINS
43!
44!***********************************************************************
45 SUBROUTINE tl_uv3dmix4 (ng, tile)
46!***********************************************************************
47!
48 USE mod_param
49 USE mod_coupling
50#ifdef DIAGNOSTICS_UV
51!! USE mod_diags
52#endif
53 USE mod_grid
54 USE mod_mixing
55 USE mod_ocean
56 USE mod_stepping
57!
58! Imported variable declarations.
59!
60 integer, intent(in) :: ng, tile
61!
62! Local variable declarations.
63!
64 character (len=*), parameter :: MyFile = &
65 & __FILE__
66!
67#include "tile.h"
68!
69#ifdef PROFILE
70 CALL wclock_on (ng, itlm, 32, __line__, myfile)
71#endif
72 CALL tl_uv3dmix4_s_tile (ng, tile, &
73 & lbi, ubi, lbj, ubj, &
74 & imins, imaxs, jmins, jmaxs, &
75 & nrhs(ng), nnew(ng), &
76#ifdef MASKING
77 & grid(ng) % pmask, &
78#endif
79 & grid(ng) % Hz, &
80 & grid(ng) % tl_Hz, &
81 & grid(ng) % om_p, &
82 & grid(ng) % om_r, &
83 & grid(ng) % on_p, &
84 & grid(ng) % on_r, &
85 & grid(ng) % pm, &
86 & grid(ng) % pmon_p, &
87 & grid(ng) % pmon_r, &
88 & grid(ng) % pn, &
89 & grid(ng) % pnom_p, &
90 & grid(ng) % pnom_r, &
91 & mixing(ng) % visc4_p, &
92 & mixing(ng) % visc4_r, &
93#ifdef DIAGNOSTICS_UV
94!! & DIAGS(ng) % DiaRUfrc, &
95!! & DIAGS(ng) % DiaRVfrc, &
96!! & DIAGS(ng) % DiaU3wrk, &
97!! & DIAGS(ng) % DiaV3wrk, &
98#endif
99 & ocean(ng) % u, &
100 & ocean(ng) % v, &
101 & coupling(ng) % tl_rufrc, &
102 & coupling(ng) % tl_rvfrc, &
103 & ocean(ng) % tl_u, &
104 & ocean(ng) % tl_v)
105#ifdef PROFILE
106 CALL wclock_off (ng, itlm, 32, __line__, myfile)
107#endif
108!
109 RETURN
110 END SUBROUTINE tl_uv3dmix4
111!
112!***********************************************************************
113 SUBROUTINE tl_uv3dmix4_s_tile (ng, tile, &
114 & LBi, UBi, LBj, UBj, &
115 & IminS, ImaxS, JminS, JmaxS, &
116 & nrhs, nnew, &
117#ifdef MASKING
118 & pmask, &
119#endif
120 & Hz, tl_Hz, &
121 & om_p, om_r, on_p, on_r, &
122 & pm, pmon_p, pmon_r, &
123 & pn, pnom_p, pnom_r, &
124 & visc4_p, visc4_r, &
125#ifdef DIAGNOSTICS_UV
126!! & DiaRUfrc, DiaRVfrc, &
127!! & DiaU3wrk, DiaV3wrk, &
128#endif
129 & u, v, &
130 & tl_rufrc, tl_rvfrc, tl_u, tl_v)
131!***********************************************************************
132!
133 USE mod_param
134 USE mod_ncparam
135 USE mod_scalars
136!
137! Imported variable declarations.
138!
139 integer, intent(in) :: ng, tile
140 integer, intent(in) :: LBi, UBi, LBj, UBj
141 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
142 integer, intent(in) :: nrhs, nnew
143
144#ifdef ASSUMED_SHAPE
145# ifdef MASKING
146 real(r8), intent(in) :: pmask(LBi:,LBj:)
147# endif
148 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
149 real(r8), intent(in) :: tl_Hz(LBi:,LBj:,:)
150 real(r8), intent(in) :: om_p(LBi:,LBj:)
151 real(r8), intent(in) :: om_r(LBi:,LBj:)
152 real(r8), intent(in) :: on_p(LBi:,LBj:)
153 real(r8), intent(in) :: on_r(LBi:,LBj:)
154 real(r8), intent(in) :: pm(LBi:,LBj:)
155 real(r8), intent(in) :: pmon_p(LBi:,LBj:)
156 real(r8), intent(in) :: pmon_r(LBi:,LBj:)
157 real(r8), intent(in) :: pn(LBi:,LBj:)
158 real(r8), intent(in) :: pnom_p(LBi:,LBj:)
159 real(r8), intent(in) :: pnom_r(LBi:,LBj:)
160 real(r8), intent(in) :: visc4_p(LBi:,LBj:)
161 real(r8), intent(in) :: visc4_r(LBi:,LBj:)
162
163 real(r8), intent(in) :: u(LBi:,LBj:,:,:)
164 real(r8), intent(in) :: v(LBi:,LBj:,:,:)
165
166# ifdef DIAGNOSTICS_UV
167!! real(r8), intent(inout) :: DiaRUfrc(LBi:,LBj:,:,:)
168!! real(r8), intent(inout) :: DiaRVfrc(LBi:,LBj:,:,:)
169!! real(r8), intent(inout) :: DiaU3wrk(LBi:,LBj:,:,:)
170!! real(r8), intent(inout) :: DiaV3wrk(LBi:,LBj:,:,:)
171# endif
172
173 real(r8), intent(inout) :: tl_rufrc(LBi:,LBj:)
174 real(r8), intent(inout) :: tl_rvfrc(LBi:,LBj:)
175 real(r8), intent(inout) :: tl_u(LBi:,LBj:,:,:)
176 real(r8), intent(inout) :: tl_v(LBi:,LBj:,:,:)
177
178#else
179
180# ifdef MASKING
181 real(r8), intent(in) :: pmask(LBi:UBi,LBj:UBj)
182# endif
183 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
184 real(r8), intent(in) :: tl_Hz(LBi:UBi,LBj:UBj,N(ng))
185 real(r8), intent(in) :: om_p(LBi:UBi,LBj:UBj)
186 real(r8), intent(in) :: om_r(LBi:UBi,LBj:UBj)
187 real(r8), intent(in) :: on_p(LBi:UBi,LBj:UBj)
188 real(r8), intent(in) :: on_r(LBi:UBi,LBj:UBj)
189 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
190 real(r8), intent(in) :: pmon_p(LBi:UBi,LBj:UBj)
191 real(r8), intent(in) :: pmon_r(LBi:UBi,LBj:UBj)
192 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
193 real(r8), intent(in) :: pnom_p(LBi:UBi,LBj:UBj)
194 real(r8), intent(in) :: pnom_r(LBi:UBi,LBj:UBj)
195 real(r8), intent(in) :: visc4_p(LBi:UBi,LBj:UBj)
196 real(r8), intent(in) :: visc4_r(LBi:UBi,LBj:UBj)
197
198 real(r8), intent(inout) :: u(LBi:UBi,LBj:UBj,N(ng),2)
199 real(r8), intent(inout) :: v(LBi:UBi,LBj:UBj,N(ng),2)
200
201# ifdef DIAGNOSTICS_UV
202!! real(r8), intent(inout) :: DiaRUfrc(LBi:UBi,LBj:UBj,3,NDM2d-1)
203!! real(r8), intent(inout) :: DiaRVfrc(LBi:UBi,LBj:UBj,3,NDM2d-1)
204!! real(r8), intent(inout) :: DiaU3wrk(LBi:UBi,LBj:UBj,N(ng),NDM3d)
205!! real(r8), intent(inout) :: DiaV3wrk(LBi:UBi,LBj:UBj,N(ng),NDM3d)
206# endif
207
208 real(r8), intent(inout) :: tl_rufrc(LBi:UBi,LBj:UBj)
209 real(r8), intent(inout) :: tl_rvfrc(LBi:UBi,LBj:UBj)
210 real(r8), intent(inout) :: tl_u(LBi:UBi,LBj:UBj,N(ng),2)
211 real(r8), intent(inout) :: tl_v(LBi:UBi,LBj:UBj,N(ng),2)
212#endif
213!
214! Local variable declarations.
215!
216 integer :: IminU, IminV, ImaxU, ImaxV
217 integer :: JminU, JminV, JmaxU, JmaxV
218 integer :: i, j, k
219
220 real(r8) :: cff, cff1, cff2
221 real(r8) :: tl_cff, tl_cff1, tl_cff2
222
223 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: LapU
224 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: LapV
225 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: UFe
226 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: VFe
227 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: UFx
228 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: VFx
229
230 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_LapU
231 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_LapV
232 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_UFe
233 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_VFe
234 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_UFx
235 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_VFx
236
237#include "set_bounds.h"
238!
239!-----------------------------------------------------------------------
240! Compute horizontal biharmonic viscosity along constant S-surfaces.
241! The biharmonic operator is computed by applying the harmonic
242! operator twice.
243!-----------------------------------------------------------------------
244!
245! Set local I- and J-ranges.
246!
247 IF (ewperiodic(ng)) THEN
248 iminu=istr-1
249 imaxu=iend+1
250 iminv=istr-1
251 imaxv=iend+1
252 ELSE
253 iminu=max(2,istru-1)
254 imaxu=min(iend+1,lm(ng))
255 iminv=max(1,istr-1)
256 imaxv=min(iend+1,lm(ng))
257 END IF
258 IF (nsperiodic(ng)) THEN
259 jminu=jstr-1
260 jmaxu=jend+1
261 jminv=jstr-1
262 jmaxv=jend+1
263 ELSE
264 jminu=max(1,jstr-1)
265 jmaxu=min(jend+1,mm(ng))
266 jminv=max(2,jstrv-1)
267 jmaxv=min(jend+1,mm(ng))
268 END IF
269!
270! Compute flux-components of the horizontal divergence of the stress
271! tensor (m4 s^-3/2) in XI- and ETA-directions. It is assumed here
272! that mixing coefficients are the squared root of the biharmonic
273! viscosity coefficient. For momentum balance purposes, the
274! thickness "Hz" appears only when computing the second harmonic
275! operator.
276!
277 k_loop : DO k=1,n(ng)
278 DO j=jminv-1,jmaxv
279 DO i=iminu-1,imaxu
280 cff=visc4_r(i,j)*0.5_r8* &
281 & (pmon_r(i,j)* &
282 & ((pn(i ,j)+pn(i+1,j))*u(i+1,j,k,nrhs)- &
283 & (pn(i-1,j)+pn(i ,j))*u(i ,j,k,nrhs))- &
284 & pnom_r(i,j)* &
285 & ((pm(i,j )+pm(i,j+1))*v(i,j+1,k,nrhs)- &
286 & (pm(i,j-1)+pm(i,j ))*v(i,j ,k,nrhs)))
287 tl_cff=visc4_r(i,j)*0.5_r8* &
288 & (pmon_r(i,j)* &
289 & ((pn(i ,j)+pn(i+1,j))*tl_u(i+1,j,k,nrhs)- &
290 & (pn(i-1,j)+pn(i ,j))*tl_u(i ,j,k,nrhs))- &
291 & pnom_r(i,j)* &
292 & ((pm(i,j )+pm(i,j+1))*tl_v(i,j+1,k,nrhs)- &
293 & (pm(i,j-1)+pm(i,j ))*tl_v(i,j ,k,nrhs)))
294 ufx(i,j)=on_r(i,j)*on_r(i,j)*cff
295 tl_ufx(i,j)=on_r(i,j)*on_r(i,j)*tl_cff
296 vfe(i,j)=om_r(i,j)*om_r(i,j)*cff
297 tl_vfe(i,j)=om_r(i,j)*om_r(i,j)*tl_cff
298 END DO
299 END DO
300 DO j=jminu,jmaxu+1
301 DO i=iminv,imaxv+1
302 cff=visc4_p(i,j)*0.5_r8* &
303 & (pmon_p(i,j)* &
304 & ((pn(i ,j-1)+pn(i ,j))*v(i ,j,k,nrhs)- &
305 & (pn(i-1,j-1)+pn(i-1,j))*v(i-1,j,k,nrhs))+ &
306 & pnom_p(i,j)* &
307 & ((pm(i-1,j )+pm(i,j ))*u(i,j ,k,nrhs)- &
308 & (pm(i-1,j-1)+pm(i,j-1))*u(i,j-1,k,nrhs)))
309 tl_cff=visc4_p(i,j)*0.5_r8* &
310 & (pmon_p(i,j)* &
311 & ((pn(i ,j-1)+pn(i ,j))*tl_v(i ,j,k,nrhs)- &
312 & (pn(i-1,j-1)+pn(i-1,j))*tl_v(i-1,j,k,nrhs))+ &
313 & pnom_p(i,j)* &
314 & ((pm(i-1,j )+pm(i,j ))*tl_u(i,j ,k,nrhs)- &
315 & (pm(i-1,j-1)+pm(i,j-1))*tl_u(i,j-1,k,nrhs)))
316#ifdef MASKING
317 cff=cff*pmask(i,j)
318 tl_cff=tl_cff*pmask(i,j)
319#endif
320 ufe(i,j)=om_p(i,j)*om_p(i,j)*cff
321 tl_ufe(i,j)=om_p(i,j)*om_p(i,j)*tl_cff
322 vfx(i,j)=on_p(i,j)*on_p(i,j)*cff
323 tl_vfx(i,j)=on_p(i,j)*on_p(i,j)*tl_cff
324 END DO
325 END DO
326!
327! Compute first harmonic operator (m s^-3/2).
328!
329 DO j=jminu,jmaxu
330 DO i=iminu,imaxu
331 lapu(i,j)=0.125_r8* &
332 & (pm(i-1,j)+pm(i,j))*(pn(i-1,j)+pn(i,j))* &
333 & ((pn(i-1,j)+pn(i,j))* &
334 & (ufx(i,j )-ufx(i-1,j))+ &
335 & (pm(i-1,j)+pm(i,j))* &
336 & (ufe(i,j+1)-ufe(i ,j)))
337 tl_lapu(i,j)=0.125_r8* &
338 & (pm(i-1,j)+pm(i,j))*(pn(i-1,j)+pn(i,j))* &
339 & ((pn(i-1,j)+pn(i,j))* &
340 & (tl_ufx(i,j )-tl_ufx(i-1,j))+ &
341 & (pm(i-1,j)+pm(i,j))* &
342 & (tl_ufe(i,j+1)-tl_ufe(i ,j)))
343 END DO
344 END DO
345 DO j=jminv,jmaxv
346 DO i=iminv,imaxv
347 lapv(i,j)=0.125_r8* &
348 & (pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))* &
349 & ((pn(i,j-1)+pn(i,j))* &
350 & (vfx(i+1,j)-vfx(i,j ))- &
351 & (pm(i,j-1)+pm(i,j))* &
352 & (vfe(i ,j)-vfe(i,j-1)))
353 tl_lapv(i,j)=0.125_r8* &
354 & (pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))* &
355 & ((pn(i,j-1)+pn(i,j))* &
356 & (tl_vfx(i+1,j)-tl_vfx(i,j ))- &
357 & (pm(i,j-1)+pm(i,j))* &
358 & (tl_vfe(i ,j)-tl_vfe(i,j-1)))
359 END DO
360 END DO
361!
362! Apply boundary conditions (other than periodic) to the first
363! harmonic operator. These are gradient or closed (free slip or
364! no slip) boundary conditions.
365!
366 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
367 IF (domain(ng)%Western_Edge(tile)) THEN
368 IF (tl_lbc(iwest,isuvel,ng)%closed) THEN
369 DO j=jminu,jmaxu
370 lapu(istr,j)=0.0_r8
371 tl_lapu(istr,j)=0.0_r8
372 END DO
373 ELSE
374 DO j=jminu,jmaxu
375 lapu(istr,j)=lapu(istr+1,j)
376 tl_lapu(istr,j)=tl_lapu(istr+1,j)
377 END DO
378 END IF
379 IF (tl_lbc(iwest,isvvel,ng)%closed) THEN
380 DO j=jminv,jmaxv
381 lapv(istr-1,j)=gamma2(ng)*lapv(istr,j)
382 tl_lapv(istr-1,j)=gamma2(ng)*tl_lapv(istr,j)
383 END DO
384 ELSE
385 DO j=jminv,jmaxv
386 lapv(istr-1,j)=0.0_r8
387 tl_lapv(istr-1,j)=0.0_r8
388 END DO
389 END IF
390 END IF
391 END IF
392!
393 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
394 IF (domain(ng)%Eastern_Edge(tile)) THEN
395 IF (tl_lbc(ieast,isuvel,ng)%closed) THEN
396 DO j=jminu,jmaxu
397 lapu(iend+1,j)=0.0_r8
398 tl_lapu(iend+1,j)=0.0_r8
399 END DO
400 ELSE
401 DO j=jminu,jmaxu
402 lapu(iend+1,j)=lapu(iend,j)
403 tl_lapu(iend+1,j)=tl_lapu(iend,j)
404 END DO
405 END IF
406 IF (tl_lbc(ieast,isvvel,ng)%closed) THEN
407 DO j=jminv,jmaxv
408 lapv(iend+1,j)=gamma2(ng)*lapv(iend,j)
409 tl_lapv(iend+1,j)=gamma2(ng)*tl_lapv(iend,j)
410 END DO
411 ELSE
412 DO j=jminv,jmaxv
413 lapv(iend+1,j)=0.0_r8
414 tl_lapv(iend+1,j)=0.0_r8
415 END DO
416 END IF
417 END IF
418 END IF
419!
420 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
421 IF (domain(ng)%Southern_Edge(tile)) THEN
422 IF (tl_lbc(isouth,isuvel,ng)%closed) THEN
423 DO i=iminu,imaxu
424 lapu(i,jstr-1)=gamma2(ng)*lapu(i,jstr)
425 tl_lapu(i,jstr-1)=gamma2(ng)*tl_lapu(i,jstr)
426 END DO
427 ELSE
428 DO i=iminu,imaxu
429 lapu(i,jstr-1)=0.0_r8
430 tl_lapu(i,jstr-1)=0.0_r8
431 END DO
432 END IF
433 IF (tl_lbc(isouth,isvvel,ng)%closed) THEN
434 DO i=iminv,imaxv
435 lapv(i,jstr)=0.0_r8
436 tl_lapv(i,jstr)=0.0_r8
437 END DO
438 ELSE
439 DO i=iminv,imaxv
440 lapv(i,jstr)=lapv(i,jstr+1)
441 tl_lapv(i,jstr)=tl_lapv(i,jstr+1)
442 END DO
443 END IF
444 END IF
445 END IF
446!
447 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
448 IF (domain(ng)%Northern_Edge(tile)) THEN
449 IF (tl_lbc(inorth,isuvel,ng)%closed) THEN
450 DO i=iminu,imaxu
451 lapu(i,jend+1)=gamma2(ng)*lapu(i,jend)
452 tl_lapu(i,jend+1)=gamma2(ng)*tl_lapu(i,jend)
453 END DO
454 ELSE
455 DO i=iminu,imaxu
456 lapu(i,jend+1)=0.0_r8
457 tl_lapu(i,jend+1)=0.0_r8
458 END DO
459 END IF
460 IF (tl_lbc(inorth,isvvel,ng)%closed) THEN
461 DO i=iminv,imaxv
462 lapv(i,jend+1)=0.0_r8
463 tl_lapv(i,jend+1)=0.0_r8
464 END DO
465 ELSE
466 DO i=iminv,imaxv
467 lapv(i,jend+1)=lapv(i,jend)
468 tl_lapv(i,jend+1)=tl_lapv(i,jend)
469 END DO
470 END IF
471 END IF
472 END IF
473!
474 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng).or. &
475 & compositegrid(iwest ,ng).or.ewperiodic(ng))) THEN
476 IF (domain(ng)%SouthWest_Corner(tile)) THEN
477 lapu(istr ,jstr-1)=0.5_r8* &
478 & (lapu(istr+1,jstr-1)+ &
479 & lapu(istr ,jstr ))
480 tl_lapu(istr ,jstr-1)=0.5_r8* &
481 & (tl_lapu(istr+1,jstr-1)+ &
482 & tl_lapu(istr ,jstr ))
483 lapv(istr-1,jstr )=0.5_r8* &
484 & (lapv(istr-1,jstr+1)+ &
485 & lapv(istr ,jstr ))
486 tl_lapv(istr-1,jstr )=0.5_r8* &
487 & (tl_lapv(istr-1,jstr+1)+ &
488 & tl_lapv(istr ,jstr ))
489 END IF
490 END IF
491
492 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng).or. &
493 & compositegrid(ieast ,ng).or.ewperiodic(ng))) THEN
494 IF (domain(ng)%SouthEast_Corner(tile)) THEN
495 lapu(iend+1,jstr-1)=0.5_r8* &
496 & (lapu(iend ,jstr-1)+ &
497 & lapu(iend+1,jstr ))
498 tl_lapu(iend+1,jstr-1)=0.5_r8* &
499 & (tl_lapu(iend ,jstr-1)+ &
500 & tl_lapu(iend+1,jstr ))
501 lapv(iend+1,jstr )=0.5_r8* &
502 & (lapv(iend ,jstr )+ &
503 & lapv(iend+1,jstr+1))
504 tl_lapv(iend+1,jstr )=0.5_r8* &
505 & (tl_lapv(iend ,jstr )+ &
506 & tl_lapv(iend+1,jstr+1))
507 END IF
508 END IF
509
510 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng).or. &
511 & compositegrid(iwest ,ng).or.ewperiodic(ng))) THEN
512 IF (domain(ng)%NorthWest_Corner(tile)) THEN
513 lapu(istr ,jend+1)=0.5_r8* &
514 & (lapu(istr+1,jend+1)+ &
515 & lapu(istr ,jend ))
516 tl_lapu(istr ,jend+1)=0.5_r8* &
517 & (tl_lapu(istr+1,jend+1)+ &
518 & tl_lapu(istr ,jend ))
519 lapv(istr-1,jend+1)=0.5_r8* &
520 & (lapv(istr ,jend+1)+ &
521 & lapv(istr-1,jend ))
522 tl_lapv(istr-1,jend+1)=0.5_r8* &
523 & (tl_lapv(istr ,jend+1)+ &
524 & tl_lapv(istr-1,jend ))
525 END IF
526 END IF
527
528 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng).or. &
529 & compositegrid(ieast ,ng).or.ewperiodic(ng))) THEN
530 IF (domain(ng)%NorthEast_Corner(tile)) THEN
531 lapu(iend+1,jend+1)=0.5_r8* &
532 & (lapu(iend ,jend+1)+ &
533 & lapu(iend+1,jend ))
534 tl_lapu(iend+1,jend+1)=0.5_r8* &
535 & (tl_lapu(iend ,jend+1)+ &
536 & tl_lapu(iend+1,jend ))
537 lapv(iend+1,jend+1)=0.5_r8* &
538 & (lapv(iend ,jend+1)+ &
539 & lapv(iend+1,jend ))
540 tl_lapv(iend+1,jend+1)=0.5_r8* &
541 & (tl_lapv(iend ,jend+1)+ &
542 & tl_lapv(iend+1,jend ))
543 END IF
544 END IF
545!
546! Compute flux-components of the horizontal divergence of the
547! harmonic stress tensor (m4/s2) in XI- and ETA-directions.
548!
549 DO j=jstrv-1,jend
550 DO i=istru-1,iend
551!^ cff=visc4_r(i,j)*Hz(i,j,k)*0.5_r8* &
552!^ & (pmon_r(i,j)* &
553!^ & ((pn(i ,j)+pn(i+1,j))*LapU(i+1,j)- &
554!^ & (pn(i-1,j)+pn(i ,j))*LapU(i ,j))- &
555!^ & pnom_r(i,j)* &
556!^ & ((pm(i,j )+pm(i,j+1))*LapV(i,j+1)- &
557!^ & (pm(i,j-1)+pm(i,j ))*LapV(i,j )))
558!^
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!^ UFx(i,j)=on_r(i,j)*on_r(i,j)*cff
575!^
576 tl_ufx(i,j)=on_r(i,j)*on_r(i,j)*tl_cff
577!^ VFe(i,j)=om_r(i,j)*om_r(i,j)*cff
578!^
579 tl_vfe(i,j)=om_r(i,j)*om_r(i,j)*tl_cff
580 END DO
581 END DO
582 DO j=jstr,jend+1
583 DO i=istr,iend+1
584!^ cff=visc4_p(i,j)*0.125_r8*(Hz(i-1,j ,k)+Hz(i,j ,k)+ &
585!^ & Hz(i-1,j-1,k)+Hz(i,j-1,k))* &
586!^ & (pmon_p(i,j)* &
587!^ & ((pn(i ,j-1)+pn(i ,j))*LapV(i ,j)- &
588!^ & (pn(i-1,j-1)+pn(i-1,j))*LapV(i-1,j))+ &
589!^ & pnom_p(i,j)* &
590!^ & ((pm(i-1,j )+pm(i,j ))*LapU(i,j )- &
591!^ & (pm(i-1,j-1)+pm(i,j-1))*LapU(i,j-1)))
592!^
593 tl_cff=visc4_p(i,j)*0.125_r8* &
594 & ((tl_hz(i-1,j ,k)+tl_hz(i,j ,k)+ &
595 & tl_hz(i-1,j-1,k)+tl_hz(i,j-1,k))* &
596 & (pmon_p(i,j)* &
597 & ((pn(i ,j-1)+pn(i ,j))*lapv(i ,j)- &
598 & (pn(i-1,j-1)+pn(i-1,j))*lapv(i-1,j))+ &
599 & pnom_p(i,j)* &
600 & ((pm(i-1,j )+pm(i,j ))*lapu(i,j )- &
601 & (pm(i-1,j-1)+pm(i,j-1))*lapu(i,j-1)))+ &
602 & (hz(i-1,j ,k)+hz(i,j ,k)+ &
603 & hz(i-1,j-1,k)+hz(i,j-1,k))* &
604 & (pmon_p(i,j)* &
605 & ((pn(i ,j-1)+pn(i ,j))*tl_lapv(i ,j)- &
606 & (pn(i-1,j-1)+pn(i-1,j))*tl_lapv(i-1,j))+ &
607 & pnom_p(i,j)* &
608 & ((pm(i-1,j )+pm(i,j ))*tl_lapu(i,j )- &
609 & (pm(i-1,j-1)+pm(i,j-1))*tl_lapu(i,j-1))))
610#ifdef MASKING
611!^ cff=cff*pmask(i,j)
612!^
613 tl_cff=tl_cff*pmask(i,j)
614#endif
615!^ UFe(i,j)=om_p(i,j)*om_p(i,j)*cff
616!^
617 tl_ufe(i,j)=om_p(i,j)*om_p(i,j)*tl_cff
618!^ VFx(i,j)=on_p(i,j)*on_p(i,j)*cff
619!^
620 tl_vfx(i,j)=on_p(i,j)*on_p(i,j)*tl_cff
621 END DO
622 END DO
623!
624! Time-step biharmonic, S-surfaces viscosity term. Notice that
625! momentum at this stage is HzU and HzV and has units m2/s. Add
626! contribution for barotropic forcing terms.
627!
628 DO j=jstr,jend
629 DO i=istru,iend
630 cff=0.25_r8*(pm(i-1,j)+pm(i,j))*(pn(i-1,j)+pn(i,j))
631!^ cff1=0.5_r8*((pn(i-1,j)+pn(i,j))* &
632!^ & (UFx(i,j )-UFx(i-1,j))+ &
633!^ & (pm(i-1,j)+pm(i,j))* &
634!^ & (UFe(i,j+1)-UFe(i ,j)))
635!^
636 tl_cff1=0.5_r8*((pn(i-1,j)+pn(i,j))* &
637 & (tl_ufx(i,j )-tl_ufx(i-1,j))+ &
638 & (pm(i-1,j)+pm(i,j))* &
639 & (tl_ufe(i,j+1)-tl_ufe(i ,j)))
640!^ cff2=dt(ng)*cff*cff1
641!^
642 tl_cff2=dt(ng)*cff*tl_cff1
643!^ rufrc(i,j)=rufrc(i,j)-cff1
644!^
645 tl_rufrc(i,j)=tl_rufrc(i,j)-tl_cff1
646!^ u(i,j,k,nnew)=u(i,j,k,nnew)-cff2
647!^
648 tl_u(i,j,k,nnew)=tl_u(i,j,k,nnew)-tl_cff2
649#ifdef DIAGNOSTICS_UV
650!! DiaRUfrc(i,j,3,M2hvis)=DiaRUfrc(i,j,3,M2hvis)-cff1
651!! DiaU3wrk(i,j,k,M3hvis)=-cff2
652#endif
653 END DO
654 END DO
655 DO j=jstrv,jend
656 DO i=istr,iend
657 cff=0.25_r8*(pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))
658!^ cff1=0.5_r8*((pn(i,j-1)+pn(i,j))* &
659!^ & (VFx(i+1,j)-VFx(i,j ))- &
660!^ & (pm(i,j-1)+pm(i,j))* &
661!^ & (VFe(i ,j)-VFe(i,j-1)))
662!^
663 tl_cff1=0.5_r8*((pn(i,j-1)+pn(i,j))* &
664 & (tl_vfx(i+1,j)-tl_vfx(i,j ))- &
665 & (pm(i,j-1)+pm(i,j))* &
666 & (tl_vfe(i ,j)-tl_vfe(i,j-1)))
667!^ cff2=dt(ng)*cff*cff1
668!^
669 tl_cff2=dt(ng)*cff*tl_cff1
670!^ rvfrc(i,j)=rvfrc(i,j)-cff1
671!^
672 tl_rvfrc(i,j)=tl_rvfrc(i,j)-tl_cff1
673!^ v(i,j,k,nnew)=v(i,j,k,nnew)-cff2
674!^
675 tl_v(i,j,k,nnew)=tl_v(i,j,k,nnew)-tl_cff2
676#ifdef DIAGNOSTICS_UV
677!! DiaRVfrc(i,j,3,M2hvis)=DiaRVfrc(i,j,3,M2hvis)-cff1
678!! DiaV3wrk(i,j,k,M3hvis)=-cff2
679#endif
680 END DO
681 END DO
682 END DO k_loop
683!
684 RETURN
685 END SUBROUTINE tl_uv3dmix4_s_tile
686
687 END MODULE tl_uv3dmix4_mod
mod_coupling
Definition mod_coupling.F:2
mod_coupling::coupling
type(t_coupling), dimension(:), allocatable coupling
Definition mod_coupling.F:102
mod_grid
Definition mod_grid.F:2
mod_grid::grid
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
mod_mixing
Definition mod_mixing.F:2
mod_mixing::mixing
type(t_mixing), dimension(:), allocatable mixing
Definition mod_mixing.F:399
mod_ncparam
Definition mod_ncparam.F:2
mod_ncparam::isvvel
integer isvvel
Definition mod_ncparam.F:507
mod_ncparam::isuvel
integer isuvel
Definition mod_ncparam.F:506
mod_ocean
Definition mod_ocean.F:2
mod_ocean::ocean
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
mod_param
Definition mod_param.F:2
mod_param::tl_lbc
type(t_lbc), dimension(:,:,:), allocatable tl_lbc
Definition mod_param.F:379
mod_param::domain
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
mod_param::lm
integer, dimension(:), allocatable lm
Definition mod_param.F:455
mod_param::itlm
integer, parameter itlm
Definition mod_param.F:663
mod_param::mm
integer, dimension(:), allocatable mm
Definition mod_param.F:456
mod_scalars
Definition mod_scalars.F:2
mod_scalars::dt
real(dp), dimension(:), allocatable dt
Definition mod_scalars.F:269
mod_scalars::ewperiodic
logical, dimension(:), allocatable ewperiodic
Definition mod_scalars.F:510
mod_scalars::iwest
integer, parameter iwest
Definition mod_scalars.F:1432
mod_scalars::nsperiodic
logical, dimension(:), allocatable nsperiodic
Definition mod_scalars.F:511
mod_scalars::gamma2
real(r8), dimension(:), allocatable gamma2
Definition mod_scalars.F:743
mod_scalars::compositegrid
logical, dimension(:,:), allocatable compositegrid
Definition mod_scalars.F:493
mod_scalars::isouth
integer, parameter isouth
Definition mod_scalars.F:1433
mod_scalars::ieast
integer, parameter ieast
Definition mod_scalars.F:1434
mod_scalars::inorth
integer, parameter inorth
Definition mod_scalars.F:1435
mod_stepping
Definition mod_stepping.F:2
mod_stepping::nrhs
integer, dimension(:), allocatable nrhs
Definition mod_stepping.F:70
mod_stepping::nnew
integer, dimension(:), allocatable nnew
Definition mod_stepping.F:69
tl_uv3dmix4_mod
Definition tl_uv3dmix4_geo.h:1
tl_uv3dmix4_mod::tl_uv3dmix4_s_tile
subroutine tl_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)
Definition tl_uv3dmix4_s.h:131
tl_uv3dmix4_mod::tl_uv3dmix4
subroutine, public tl_uv3dmix4(ng, tile)
Definition tl_uv3dmix4_geo.h:44
wclock_off
recursive subroutine wclock_off(ng, model, region, line, routine)
Definition timers.F:148
wclock_on
recursive subroutine wclock_on(ng, model, region, line, routine)
Definition timers.F:3