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

Functions/Subroutines

subroutine, public t3dmix4 (ng, tile)
 
subroutine t3dmix4_geo_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, nstp, nnew, umask, vmask, umask_wet, vmask_wet, om_v, on_u, pm, pn, hz, z_r, diff3d_u, diff3d_v, diff3d_r, diff4, tclm, diatwrk, t)
 
subroutine t3dmix4_iso_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, nstp, nnew, umask, vmask, umask_wet, vmask_wet, om_v, on_u, pm, pn, hz, z_r, diff3d_u, diff3d_v, diff3d_r, diff4, pden, tclm, diatwrk, t)
 
subroutine t3dmix4_s_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, nstp, nnew, umask, vmask, umask_wet, vmask_wet, hz, pmon_u, pnom_v, pm, pn, diff3d_u, diff3d_v, diff3d_r, diff4, tclm, diatwrk, t)
 

Function/Subroutine Documentation

◆ t3dmix4()

subroutine public t3dmix4_mod::t3dmix4 ( integer, intent(in) ng,
integer, intent(in) tile )

Definition at line 23 of file t3dmix4_geo.h.

24!***********************************************************************
25!
26 USE mod_param
27#ifdef TS_MIX_CLIMA
28 USE mod_clima
29#endif
30#ifdef DIAGNOSTICS_TS
31 USE mod_diags
32#endif
33 USE mod_grid
34 USE mod_mixing
35 USE mod_ocean
36 USE mod_stepping
37!
38! Imported variable declarations.
39!
40 integer, intent(in) :: ng, tile
41!
42! Local variable declarations.
43!
44 character (len=*), parameter :: MyFile = &
45 & __FILE__
46!
47#include "tile.h"
48!
49#ifdef PROFILE
50 CALL wclock_on (ng, inlm, 28, __line__, myfile)
51#endif
52 CALL t3dmix4_geo_tile (ng, tile, &
53 & lbi, ubi, lbj, ubj, &
54 & imins, imaxs, jmins, jmaxs, &
55 & nrhs(ng), nstp(ng), nnew(ng), &
56#ifdef MASKING
57 & grid(ng) % umask, &
58 & grid(ng) % vmask, &
59#endif
60#ifdef WET_DRY
61 & grid(ng) % umask_wet, &
62 & grid(ng) % vmask_wet, &
63#endif
64 & grid(ng) % om_v, &
65 & grid(ng) % on_u, &
66 & grid(ng) % pm, &
67 & grid(ng) % pn, &
68 & grid(ng) % Hz, &
69 & grid(ng) % z_r, &
70#ifdef DIFF_3DCOEF
71# ifdef TS_U3ADV_SPLIT
72 & mixing(ng) % diff3d_u, &
73 & mixing(ng) % diff3d_v, &
74# else
75 & mixing(ng) % diff3d_r, &
76# endif
77#else
78 & mixing(ng) % diff4, &
79#endif
80#ifdef TS_MIX_CLIMA
81 & clima(ng) % tclm, &
82#endif
83#ifdef DIAGNOSTICS_TS
84 & diags(ng) % DiaTwrk, &
85#endif
86 & ocean(ng) % t)
87#ifdef PROFILE
88 CALL wclock_off (ng, inlm, 28, __line__, myfile)
89#endif
90!
91 RETURN
mod_clima
Definition mod_clima.F:2
mod_clima::clima
type(t_clima), dimension(:), allocatable clima
Definition mod_clima.F:153
mod_diags
Definition mod_diags.F:2
mod_diags::diags
type(t_diags), dimension(:), allocatable diags
Definition mod_diags.F:100
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::inlm
integer, parameter inlm
Definition mod_param.F:662
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
mod_stepping::nstp
integer, dimension(:), allocatable nstp
Definition mod_stepping.F:71
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_clima::clima, mod_diags::diags, mod_grid::grid, mod_param::inlm, mod_mixing::mixing, mod_stepping::nnew, mod_stepping::nrhs, mod_stepping::nstp, mod_ocean::ocean, t3dmix4_geo_tile(), wclock_off(), and wclock_on().

Referenced by rhs3d_mod::rhs3d().

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

◆ t3dmix4_geo_tile()

subroutine t3dmix4_mod::t3dmix4_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) nstp,
integer, intent(in) nnew,
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) umask_wet,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) vmask_wet,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) om_v,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) on_u,
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) z_r,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) diff3d_u,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) diff3d_v,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) diff3d_r,
real(r8), dimension(lbi:ubi,lbj:ubj,nt(ng)), intent(in) diff4,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),nt(ng)), intent(in) tclm,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),nt(ng), ndt), intent(inout) diatwrk,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),3,nt(ng)), intent(inout) t )
private

Definition at line 95 of file t3dmix4_geo.h.

123!***********************************************************************
124!
125 USE mod_param
126 USE mod_ncparam
127 USE mod_scalars
128!
129! Imported variable declarations.
130!
131 integer, intent(in) :: ng, tile
132 integer, intent(in) :: LBi, UBi, LBj, UBj
133 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
134 integer, intent(in) :: nrhs, nstp, nnew
135
136#ifdef ASSUMED_SHAPE
137# ifdef MASKING
138 real(r8), intent(in) :: umask(LBi:,LBj:)
139 real(r8), intent(in) :: vmask(LBi:,LBj:)
140# endif
141# ifdef WET_DRY
142 real(r8), intent(in) :: umask_wet(LBi:,LBj:)
143 real(r8), intent(in) :: vmask_wet(LBi:,LBj:)
144# endif
145# ifdef DIFF_3DCOEF
146# ifdef TS_U3ADV_SPLIT
147 real(r8), intent(in) :: diff3d_u(LBi:,LBj:,:)
148 real(r8), intent(in) :: diff3d_v(LBi:,LBj:,:)
149# else
150 real(r8), intent(in) :: diff3d_r(LBi:,LBj:,:)
151# endif
152# else
153 real(r8), intent(in) :: diff4(LBi:,LBj:,:)
154# endif
155 real(r8), intent(in) :: om_v(LBi:,LBj:)
156 real(r8), intent(in) :: on_u(LBi:,LBj:)
157 real(r8), intent(in) :: pm(LBi:,LBj:)
158 real(r8), intent(in) :: pn(LBi:,LBj:)
159 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
160 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
161# ifdef TS_MIX_CLIMA
162 real(r8), intent(in) :: tclm(LBi:,LBj:,:,:)
163# endif
164# ifdef DIAGNOSTICS_TS
165 real(r8), intent(inout) :: DiaTwrk(LBi:,LBj:,:,:,:)
166# endif
167 real(r8), intent(inout) :: t(LBi:,LBj:,:,:,:)
168#else
169# ifdef MASKING
170 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
171 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
172# endif
173# ifdef WET_DRY
174 real(r8), intent(in) :: umask_wet(LBi:UBi,LBj:UBj)
175 real(r8), intent(in) :: vmask_wet(LBi:UBi,LBj:UBj)
176# endif
177# ifdef DIFF_3DCOEF
178# ifdef TS_U3ADV_SPLIT
179 real(r8), intent(in) :: diff3d_u(LBi:UBi,LBj:UBj,N(ng))
180 real(r8), intent(in) :: diff3d_v(LBi:UBi,LBj:UBj,N(ng))
181# else
182 real(r8), intent(in) :: diff3d_r(LBi:UBi,LBj:UBj,N(ng))
183# endif
184# else
185 real(r8), intent(in) :: diff4(LBi:UBi,LBj:UBj,NT(ng))
186# endif
187 real(r8), intent(in) :: om_v(LBi:UBi,LBj:UBj)
188 real(r8), intent(in) :: on_u(LBi:UBi,LBj:UBj)
189 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
190 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
191 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
192 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
193# ifdef TS_MIX_CLIMA
194 real(r8), intent(in) :: tclm(LBi:UBi,LBj:UBj,N(ng),NT(ng))
195# endif
196# ifdef DIAGNOSTICS_TS
197 real(r8), intent(inout) :: DiaTwrk(LBi:UBi,LBj:UBj,N(ng),NT(ng), &
198 & NDT)
199# endif
200 real(r8), intent(inout) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
201#endif
202!
203! Local variable declarations.
204!
205 integer :: Imin, Imax, Jmin, Jmax
206 integer :: i, itrc, j, k, k1, k2
207
208 real(r8) :: cff, cff1, cff2, cff3, cff4, dife, difx
209
210 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,N(ng)) :: LapT
211
212 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FE
213 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FX
214
215 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: FS
216 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dTde
217 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dTdx
218 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dTdz
219 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dZde
220 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dZdx
221
222#include "set_bounds.h"
223!
224!-----------------------------------------------------------------------
225! Compute horizontal biharmonic diffusion along geopotential
226! surfaces. The biharmonic operator is computed by applying
227! the harmonic operator twice.
228!-----------------------------------------------------------------------
229!
230! Set local I- and J-ranges.
231!
232 IF (ewperiodic(ng)) THEN
233 imin=istr-1
234 imax=iend+1
235 ELSE
236 imin=max(istr-1,1)
237 imax=min(iend+1,lm(ng))
238 END IF
239 IF (nsperiodic(ng)) THEN
240 jmin=jstr-1
241 jmax=jend+1
242 ELSE
243 jmin=max(jstr-1,1)
244 jmax=min(jend+1,mm(ng))
245 END IF
246!
247! Compute horizontal and vertical gradients associated with the
248! first rotated harmonic operator. Notice the recursive blocking
249! sequence. The vertical placement of the gradients is:
250!
251! dTdx,dTde(:,:,k1) k rho-points
252! dTdx,dTde(:,:,k2) k+1 rho-points
253! FS,dTdz(:,:,k1) k-1/2 W-points
254! FS,dTdz(:,:,k2) k+1/2 W-points
255!
256#ifdef TS_MIX_STABILITY
257! In order to increase stability, the rotated biharmonic is applied
258! as: 3/4 t(:,:,:,nrhs,:) + 1/4 t(:,:,:,nstp,:).
259!
260#endif
261
262 t_loop : DO itrc=1,nt(ng)
263 k2=1
264 k_loop1 : DO k=0,n(ng)
265 k1=k2
266 k2=3-k1
267 IF (k.lt.n(ng)) THEN
268 DO j=jmin,jmax
269 DO i=imin,imax+1
270 cff=0.5_r8*(pm(i,j)+pm(i-1,j))
271#ifdef MASKING
272 cff=cff*umask(i,j)
273#endif
274#ifdef WET_DRY
275 cff=cff*umask_wet(i,j)
276#endif
277 dzdx(i,j,k2)=cff*(z_r(i ,j,k+1)- &
278 & z_r(i-1,j,k+1))
279#if defined TS_MIX_STABILITY
280 dtdx(i,j,k2)=cff*(0.75_r8*(t(i ,j,k+1,nrhs,itrc)- &
281 & t(i-1,j,k+1,nrhs,itrc))+ &
282 & 0.25_r8*(t(i ,j,k+1,nstp,itrc)- &
283 & t(i-1,j,k+1,nstp,itrc)))
284#elif defined TS_MIX_CLIMA
285 IF (ltracerclm(itrc,ng)) THEN
286 dtdx(i,j,k2)=cff*((t(i ,j,k+1,nrhs,itrc)- &
287 & tclm(i ,j,k+1,itrc))- &
288 & (t(i-1,j,k+1,nrhs,itrc)- &
289 & tclm(i-1,j,k+1,itrc)))
290 ELSE
291 dtdx(i,j,k2)=cff*(t(i ,j,k+1,nrhs,itrc)- &
292 & t(i-1,j,k+1,nrhs,itrc))
293 END IF
294#else
295 dtdx(i,j,k2)=cff*(t(i ,j,k+1,nrhs,itrc)- &
296 & t(i-1,j,k+1,nrhs,itrc))
297#endif
298 END DO
299 END DO
300 DO j=jmin,jmax+1
301 DO i=imin,imax
302 cff=0.5_r8*(pn(i,j)+pn(i,j-1))
303#ifdef MASKING
304 cff=cff*vmask(i,j)
305#endif
306#ifdef WET_DRY
307 cff=cff*vmask_wet(i,j)
308#endif
309 dzde(i,j,k2)=cff*(z_r(i,j ,k+1)- &
310 & z_r(i,j-1,k+1))
311#if defined TS_MIX_STABILITY
312 dtde(i,j,k2)=cff*(0.75_r8*(t(i,j ,k+1,nrhs,itrc)- &
313 & t(i,j-1,k+1,nrhs,itrc))+ &
314 & 0.25_r8*(t(i,j ,k+1,nstp,itrc)- &
315 & t(i,j-1,k+1,nstp,itrc)))
316#elif defined TS_MIX_CLIMA
317 IF (ltracerclm(itrc,ng)) THEN
318 dtde(i,j,k2)=cff*((t(i,j ,k+1,nrhs,itrc)- &
319 & tclm(i,j ,k+1,itrc))- &
320 & (t(i,j-1,k+1,nrhs,itrc)- &
321 & tclm(i,j-1,k+1,itrc)))
322 ELSE
323 dtde(i,j,k2)=cff*(t(i,j ,k+1,nrhs,itrc)- &
324 & t(i,j-1,k+1,nrhs,itrc))
325 END IF
326#else
327 dtde(i,j,k2)=cff*(t(i,j ,k+1,nrhs,itrc)- &
328 & t(i,j-1,k+1,nrhs,itrc))
329#endif
330 END DO
331 END DO
332 END IF
333 IF ((k.eq.0).or.(k.eq.n(ng))) THEN
334 DO j=jmin-1,jmax+1
335 DO i=imin-1,imax+1
336 dtdz(i,j,k2)=0.0_r8
337 fs(i,j,k2)=0.0_r8
338 END DO
339 END DO
340 ELSE
341 DO j=jmin-1,jmax+1
342 DO i=imin-1,imax+1
343 cff=1.0_r8/(z_r(i,j,k+1)- &
344 & z_r(i,j,k ))
345#if defined TS_MIX_STABILITY
346 dtdz(i,j,k2)=cff*(0.75_r8*(t(i,j,k+1,nrhs,itrc)- &
347 & t(i,j,k ,nrhs,itrc))+ &
348 & 0.25_r8*(t(i,j,k+1,nstp,itrc)- &
349 & t(i,j,k ,nstp,itrc)))
350#elif defined TS_MIX_CLIMA
351 IF (ltracerclm(itrc,ng)) THEN
352 dtdz(i,j,k2)=cff*((t(i,j,k+1,nrhs,itrc)- &
353 & tclm(i,j,k+1,itrc))- &
354 & (t(i,j,k ,nrhs,itrc)- &
355 & tclm(i,j,k ,itrc)))
356 ELSE
357 dtdz(i,j,k2)=cff*(t(i,j,k+1,nrhs,itrc)- &
358 & t(i,j,k ,nrhs,itrc))
359 END IF
360#else
361 dtdz(i,j,k2)=cff*(t(i,j,k+1,nrhs,itrc)- &
362 & t(i,j,k ,nrhs,itrc))
363#endif
364 END DO
365 END DO
366 END IF
367 IF (k.gt.0) THEN
368 DO j=jmin,jmax
369 DO i=imin,imax+1
370#ifdef DIFF_3DCOEF
371# ifdef TS_U3ADV_SPLIT
372 cff=0.5_r8*diff3d_u(i,j,k)*on_u(i,j)
373# else
374 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i-1,j,k))* &
375 & on_u(i,j)
376# endif
377#else
378 cff=0.25_r8*(diff4(i,j,itrc)+diff4(i-1,j,itrc))* &
379 & on_u(i,j)
380#endif
381 fx(i,j)=cff* &
382 & (hz(i,j,k)+hz(i-1,j,k))* &
383 & (dtdx(i,j,k1)- &
384 & 0.5_r8*(min(dzdx(i,j,k1),0.0_r8)* &
385 & (dtdz(i-1,j,k1)+ &
386 & dtdz(i ,j,k2))+ &
387 & max(dzdx(i,j,k1),0.0_r8)* &
388 & (dtdz(i-1,j,k2)+ &
389 & dtdz(i ,j,k1))))
390 END DO
391 END DO
392 DO j=jmin,jmax+1
393 DO i=imin,imax
394#ifdef DIFF_3DCOEF
395# ifdef TS_U3ADV_SPLIT
396 cff=0.5_r8*diff3d_v(i,j,k)*om_v(i,j)
397# else
398 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i,j-1,k))* &
399 & om_v(i,j)
400# endif
401#else
402 cff=0.25_r8*(diff4(i,j,itrc)+diff4(i,j-1,itrc))* &
403 & om_v(i,j)
404#endif
405 fe(i,j)=cff* &
406 & (hz(i,j,k)+hz(i,j-1,k))* &
407 & (dtde(i,j,k1)- &
408 & 0.5_r8*(min(dzde(i,j,k1),0.0_r8)* &
409 & (dtdz(i,j-1,k1)+ &
410 & dtdz(i,j ,k2))+ &
411 & max(dzde(i,j,k1),0.0_r8)* &
412 & (dtdz(i,j-1,k2)+ &
413 & dtdz(i,j ,k1))))
414 END DO
415 END DO
416 IF (k.lt.n(ng)) THEN
417 DO j=jmin,jmax
418 DO i=imin,imax
419#ifdef DIFF_3DCOEF
420# ifdef TS_U3ADV_SPLIT
421 difx=0.125_r8*(diff3d_u(i,j,k )+diff3d_u(i+1,j,k )+ &
422 & diff3d_u(i,j,k+1)+diff3d_u(i+1,j,k+1))
423 dife=0.125_r8*(diff3d_v(i,j,k )+diff3d_v(i,j+1,k )+ &
424 & diff3d_v(i,j,k+1)+diff3d_v(i,j+1,k+1))
425# else
426 difx=0.5_r8*diff3d_r(i,j,k)
427 dife=difx
428# endif
429#else
430 difx=0.5_r8*diff4(i,j,itrc)
431 dife=difx
432#endif
433 cff1=min(dzdx(i ,j,k1),0.0_r8)
434 cff2=min(dzdx(i+1,j,k2),0.0_r8)
435 cff3=max(dzdx(i ,j,k2),0.0_r8)
436 cff4=max(dzdx(i+1,j,k1),0.0_r8)
437 fs(i,j,k2)=difx* &
438 & (cff1*(cff1*dtdz(i,j,k2)- &
439 & dtdx(i ,j,k1))+ &
440 & cff2*(cff2*dtdz(i,j,k2)- &
441 & dtdx(i+1,j,k2))+ &
442 & cff3*(cff3*dtdz(i,j,k2)- &
443 & dtdx(i ,j,k2))+ &
444 & cff4*(cff4*dtdz(i,j,k2)- &
445 & dtdx(i+1,j,k1)))
446!
447 cff1=min(dzde(i,j ,k1),0.0_r8)
448 cff2=min(dzde(i,j+1,k2),0.0_r8)
449 cff3=max(dzde(i,j ,k2),0.0_r8)
450 cff4=max(dzde(i,j+1,k1),0.0_r8)
451 fs(i,j,k2)=fs(i,j,k2)+ &
452 & dife* &
453 & (cff1*(cff1*dtdz(i,j,k2)- &
454 & dtde(i,j ,k1))+ &
455 & cff2*(cff2*dtdz(i,j,k2)- &
456 & dtde(i,j+1,k2))+ &
457 & cff3*(cff3*dtdz(i,j,k2)- &
458 & dtde(i,j ,k2))+ &
459 & cff4*(cff4*dtdz(i,j,k2)- &
460 & dtde(i,j+1,k1)))
461 END DO
462 END DO
463 END IF
464!
465! Compute first harmonic operator, without mixing coefficient.
466! Multiply by the metrics of the second harmonic operator. Save
467! into work array "LapT".
468!
469 DO j=jmin,jmax
470 DO i=imin,imax
471 cff=pm(i,j)*pn(i,j)
472 cff1=1.0_r8/hz(i,j,k)
473 lapt(i,j,k)=cff1*(cff* &
474 & (fx(i+1,j)-fx(i,j)+ &
475 & fe(i,j+1)-fe(i,j))+ &
476 & (fs(i,j,k2)-fs(i,j,k1)))
477 END DO
478 END DO
479 END IF
480 END DO k_loop1
481!
482! Apply boundary conditions (except periodic; closed or gradient)
483! to the first harmonic operator.
484!
485 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
486 IF (domain(ng)%Western_Edge(tile)) THEN
487 IF (lbc(iwest,istvar(itrc),ng)%closed) THEN
488 DO k=1,n(ng)
489 DO j=jmin,jmax
490 lapt(istr-1,j,k)=0.0_r8
491 END DO
492 END DO
493 ELSE
494 DO k=1,n(ng)
495 DO j=jmin,jmax
496 lapt(istr-1,j,k)=lapt(istr,j,k)
497 END DO
498 END DO
499 END IF
500 END IF
501 END IF
502!
503 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
504 IF (domain(ng)%Eastern_Edge(tile)) THEN
505 IF (lbc(ieast,istvar(itrc),ng)%closed) THEN
506 DO k=1,n(ng)
507 DO j=jmin,jmax
508 lapt(iend+1,j,k)=0.0_r8
509 END DO
510 END DO
511 ELSE
512 DO k=1,n(ng)
513 DO j=jmin,jmax
514 lapt(iend+1,j,k)=lapt(iend,j,k)
515 END DO
516 END DO
517 END IF
518 END IF
519 END IF
520!
521 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
522 IF (domain(ng)%Southern_Edge(tile)) THEN
523 IF (lbc(isouth,istvar(itrc),ng)%closed) THEN
524 DO k=1,n(ng)
525 DO i=imin,imax
526 lapt(i,jstr-1,k)=0.0_r8
527 END DO
528 END DO
529 ELSE
530 DO k=1,n(ng)
531 DO i=imin,imax
532 lapt(i,jstr-1,k)=lapt(i,jstr,k)
533 END DO
534 END DO
535 END IF
536 END IF
537 END IF
538!
539 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
540 IF (domain(ng)%Northern_Edge(tile)) THEN
541 IF (lbc(inorth,istvar(itrc),ng)%closed) THEN
542 DO k=1,n(ng)
543 DO i=imin,imax
544 lapt(i,jend+1,k)=0.0_r8
545 END DO
546 END DO
547 ELSE
548 DO k=1,n(ng)
549 DO i=imin,imax
550 lapt(i,jend+1,k)=lapt(i,jend,k)
551 END DO
552 END DO
553 END IF
554 END IF
555 END IF
556!
557 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng).or. &
558 & compositegrid(iwest ,ng).or.ewperiodic(ng))) THEN
559 IF (domain(ng)%SouthWest_Corner(tile)) THEN
560 DO k=1,n(ng)
561 lapt(istr-1,jstr-1,k)=0.5_r8* &
562 & (lapt(istr ,jstr-1,k)+ &
563 & lapt(istr-1,jstr ,k))
564 END DO
565 END IF
566 END IF
567
568 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng).or. &
569 & compositegrid(ieast ,ng).or.ewperiodic(ng))) THEN
570 IF (domain(ng)%SouthEast_Corner(tile)) THEN
571 DO k=1,n(ng)
572 lapt(iend+1,jstr-1,k)=0.5_r8* &
573 & (lapt(iend ,jstr-1,k)+ &
574 & lapt(iend+1,jstr ,k))
575 END DO
576 END IF
577 END IF
578
579 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng).or. &
580 & compositegrid(iwest ,ng).or.ewperiodic(ng))) THEN
581 IF (domain(ng)%NorthWest_Corner(tile)) THEN
582 DO k=1,n(ng)
583 lapt(istr-1,jend+1,k)=0.5_r8* &
584 & (lapt(istr ,jend+1,k)+ &
585 & lapt(istr-1,jend ,k))
586 END DO
587 END IF
588 END IF
589
590 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng).or. &
591 & compositegrid(ieast ,ng).or.ewperiodic(ng))) THEN
592 IF (domain(ng)%NorthEast_Corner(tile)) THEN
593 DO k=1,n(ng)
594 lapt(iend+1,jend+1,k)=0.5_r8* &
595 & (lapt(iend ,jend+1,k)+ &
596 & lapt(iend+1,jend ,k))
597 END DO
598 END IF
599 END IF
600!
601! Compute horizontal and vertical gradients associated with the
602! second rotated harmonic operator.
603!
604 k2=1
605 k_loop2: DO k=0,n(ng)
606 k1=k2
607 k2=3-k1
608 IF (k.lt.n(ng)) THEN
609 DO j=jstr,jend
610 DO i=istr,iend+1
611 cff=0.5_r8*(pm(i,j)+pm(i-1,j))
612#ifdef MASKING
613 cff=cff*umask(i,j)
614#endif
615#ifdef WET_DRY
616 cff=cff*umask_wet(i,j)
617#endif
618 dzdx(i,j,k2)=cff*(z_r(i ,j,k+1)- &
619 & z_r(i-1,j,k+1))
620 dtdx(i,j,k2)=cff*(lapt(i ,j,k+1)- &
621 & lapt(i-1,j,k+1))
622 END DO
623 END DO
624 DO j=jstr,jend+1
625 DO i=istr,iend
626 cff=0.5_r8*(pn(i,j)+pn(i,j-1))
627#ifdef MASKING
628 cff=cff*vmask(i,j)
629#endif
630#ifdef WET_DRY
631 cff=cff*vmask_wet(i,j)
632#endif
633 dzde(i,j,k2)=cff*(z_r(i,j ,k+1)- &
634 & z_r(i,j-1,k+1))
635 dtde(i,j,k2)=cff*(lapt(i,j ,k+1)- &
636 & lapt(i,j-1,k+1))
637 END DO
638 END DO
639 END IF
640 IF ((k.eq.0).or.(k.eq.n(ng))) THEN
641 DO j=jstr-1,jend+1
642 DO i=istr-1,iend+1
643 dtdz(i,j,k2)=0.0_r8
644 fs(i,j,k2)=0.0_r8
645 END DO
646 END DO
647 ELSE
648 DO j=jstr-1,jend+1
649 DO i=istr-1,iend+1
650 cff=1.0_r8/(z_r(i,j,k+1)- &
651 & z_r(i,j,k ))
652 dtdz(i,j,k2)=cff*(lapt(i,j,k+1)- &
653 & lapt(i,j,k ))
654 END DO
655 END DO
656 END IF
657!
658! Compute components of the rotated tracer flux (T m4/s) along
659! geopotential surfaces.
660!
661 IF (k.gt.0) THEN
662 DO j=jstr,jend
663 DO i=istr,iend+1
664#ifdef DIFF_3DCOEF
665# ifdef TS_U3ADV_SPLIT
666 cff=0.5_r8*diff3d_u(i,j,k)*on_u(i,j)
667# else
668 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i-1,j,k))* &
669 & on_u(i,j)
670# endif
671#else
672 cff=0.25_r8*(diff4(i,j,itrc)+diff4(i-1,j,itrc))* &
673 & on_u(i,j)
674#endif
675 fx(i,j)=cff* &
676 & (hz(i,j,k)+hz(i-1,j,k))* &
677 & (dtdx(i ,j,k1)- &
678 & 0.5_r8*(min(dzdx(i,j,k1),0.0_r8)* &
679 & (dtdz(i-1,j,k1)+ &
680 & dtdz(i ,j,k2))+ &
681 & max(dzdx(i,j,k1),0.0_r8)* &
682 & (dtdz(i-1,j,k2)+ &
683 & dtdz(i ,j,k1))))
684 END DO
685 END DO
686 DO j=jstr,jend+1
687 DO i=istr,iend
688#ifdef DIFF_3DCOEF
689# ifdef TS_U3ADV_SPLIT
690 cff=0.5_r8*diff3d_v(i,j,k)*om_v(i,j)
691# else
692 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i,j-1,k))* &
693 & om_v(i,j)
694# endif
695#else
696 cff=0.25_r8*(diff4(i,j,itrc)+diff4(i,j-1,itrc))* &
697 & om_v(i,j)
698#endif
699 fe(i,j)=cff* &
700 & (hz(i,j,k)+hz(i,j-1,k))* &
701 & (dtde(i,j,k1)- &
702 & 0.5_r8*(min(dzde(i,j,k1),0.0_r8)* &
703 & (dtdz(i,j-1,k1)+ &
704 & dtdz(i,j ,k2))+ &
705 & max(dzde(i,j,k1),0.0_r8)* &
706 & (dtdz(i,j-1,k2)+ &
707 & dtdz(i,j ,k1))))
708 END DO
709 END DO
710 IF (k.lt.n(ng)) THEN
711 DO j=jstr,jend
712 DO i=istr,iend
713#ifdef DIFF_3DCOEF
714# ifdef TS_U3ADV_SPLIT
715 difx=0.125_r8*(diff3d_u(i,j,k )+diff3d_u(i+1,j,k )+ &
716 & diff3d_u(i,j,k+1)+diff3d_u(i+1,j,k+1))
717 dife=0.125_r8*(diff3d_v(i,j,k )+diff3d_v(i,j+1,k )+ &
718 & diff3d_v(i,j,k+1)+diff3d_v(i,j+1,k+1))
719# else
720 difx=0.5_r8*diff3d_r(i,j,k)
721 dife=difx
722# endif
723#else
724 difx=0.5_r8*diff4(i,j,itrc)
725 dife=difx
726#endif
727 cff1=min(dzdx(i ,j,k1),0.0_r8)
728 cff2=min(dzdx(i+1,j,k2),0.0_r8)
729 cff3=max(dzdx(i ,j,k2),0.0_r8)
730 cff4=max(dzdx(i+1,j,k1),0.0_r8)
731 fs(i,j,k2)=difx* &
732 & (cff1*(cff1*dtdz(i,j,k2)- &
733 & dtdx(i ,j,k1))+ &
734 & cff2*(cff2*dtdz(i,j,k2)- &
735 & dtdx(i+1,j,k2))+ &
736 & cff3*(cff3*dtdz(i,j,k2)- &
737 & dtdx(i ,j,k2))+ &
738 & cff4*(cff4*dtdz(i,j,k2)- &
739 & dtdx(i+1,j,k1)))
740!
741 cff1=min(dzde(i,j ,k1),0.0_r8)
742 cff2=min(dzde(i,j+1,k2),0.0_r8)
743 cff3=max(dzde(i,j ,k2),0.0_r8)
744 cff4=max(dzde(i,j+1,k1),0.0_r8)
745 fs(i,j,k2)=fs(i,j,k2)+ &
746 & dife* &
747 & (cff1*(cff1*dtdz(i,j,k2)- &
748 & dtde(i,j ,k1))+ &
749 & cff2*(cff2*dtdz(i,j,k2)- &
750 & dtde(i,j+1,k2))+ &
751 & cff3*(cff3*dtdz(i,j,k2)- &
752 & dtde(i,j ,k2))+ &
753 & cff4*(cff4*dtdz(i,j,k2)- &
754 & dtde(i,j+1,k1)))
755 END DO
756 END DO
757 END IF
758!
759! Time-step biharmonic, geopotential diffusion term (m Tunits).
760!
761 DO j=jstr,jend
762 DO i=istr,iend
763 cff=dt(ng)*pm(i,j)*pn(i,j)
764 cff1=cff*(fx(i+1,j )-fx(i,j))
765 cff2=cff*(fe(i ,j+1)-fe(i,j))
766 cff3=dt(ng)*(fs(i,j,k2)-fs(i,j,k1))
767 cff4=cff1+cff2+cff3
768 t(i,j,k,nnew,itrc)=t(i,j,k,nnew,itrc)-cff4
769#ifdef DIAGNOSTICS_TS
770 diatwrk(i,j,k,itrc,itxdif)=-cff1
771 diatwrk(i,j,k,itrc,itydif)=-cff2
772 diatwrk(i,j,k,itrc,itsdif)=-cff3
773 diatwrk(i,j,k,itrc,ithdif)=-cff4
774#endif
775 END DO
776 END DO
777 END IF
778 END DO k_loop2
779 END DO t_loop
780!
781 RETURN
mod_ncparam
Definition mod_ncparam.F:2
mod_ncparam::istvar
integer, dimension(:), allocatable istvar
Definition mod_ncparam.F:520
mod_param::n
integer, dimension(:), allocatable n
Definition mod_param.F:479
mod_param::lbc
type(t_lbc), dimension(:,:,:), allocatable lbc
Definition mod_param.F:375
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::nt
integer, dimension(:), allocatable nt
Definition mod_param.F:489
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::itxdif
integer itxdif
Definition mod_scalars.F:145
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::itsdif
integer itsdif
Definition mod_scalars.F:147
mod_scalars::ieast
integer, parameter ieast
Definition mod_scalars.F:1434
mod_scalars::itydif
integer itydif
Definition mod_scalars.F:146
mod_scalars::ltracerclm
logical, dimension(:,:), allocatable ltracerclm
Definition mod_scalars.F:528
mod_scalars::inorth
integer, parameter inorth
Definition mod_scalars.F:1435
mod_scalars::ithdif
integer ithdif
Definition mod_scalars.F:144

References mod_scalars::compositegrid, mod_param::domain, mod_scalars::dt, mod_scalars::ewperiodic, mod_scalars::ieast, mod_scalars::inorth, mod_scalars::isouth, mod_ncparam::istvar, mod_scalars::ithdif, mod_scalars::itsdif, mod_scalars::itxdif, mod_scalars::itydif, mod_scalars::iwest, mod_param::lbc, mod_param::lm, mod_scalars::ltracerclm, mod_param::mm, and mod_scalars::nsperiodic.

Referenced by t3dmix4().

Here is the caller graph for this function:

◆ t3dmix4_iso_tile()

subroutine t3dmix4_mod::t3dmix4_iso_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) nstp,
integer, intent(in) nnew,
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) umask_wet,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) vmask_wet,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) om_v,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) on_u,
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) z_r,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) diff3d_u,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) diff3d_v,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) diff3d_r,
real(r8), dimension(lbi:ubi,lbj:ubj,nt(ng)), intent(in) diff4,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) pden,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),nt(ng)), intent(in) tclm,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),nt(ng), ndt), intent(inout) diatwrk,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),3,nt(ng)), intent(inout) t )
private

Definition at line 96 of file t3dmix4_iso.h.

125!***********************************************************************
126!
127 USE mod_param
128 USE mod_ncparam
129 USE mod_scalars
130!
131! Imported variable declarations.
132!
133 integer, intent(in) :: ng, tile
134 integer, intent(in) :: LBi, UBi, LBj, UBj
135 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
136 integer, intent(in) :: nrhs, nstp, nnew
137
138#ifdef ASSUMED_SHAPE
139# ifdef MASKING
140 real(r8), intent(in) :: umask(LBi:,LBj:)
141 real(r8), intent(in) :: vmask(LBi:,LBj:)
142# endif
143# ifdef WET_DRY
144 real(r8), intent(in) :: umask_wet(LBi:,LBj:)
145 real(r8), intent(in) :: vmask_wet(LBi:,LBj:)
146# endif
147# ifdef DIFF_3DCOEF
148# ifdef TS_U3ADV_SPLIT
149 real(r8), intent(in) :: diff3d_u(LBi:,LBj:,:)
150 real(r8), intent(in) :: diff3d_v(LBi:,LBj:,:)
151# else
152 real(r8), intent(in) :: diff3d_r(LBi:,LBj:,:)
153# endif
154# else
155 real(r8), intent(in) :: diff4(LBi:,LBj:,:)
156# endif
157 real(r8), intent(in) :: om_v(LBi:,LBj:)
158 real(r8), intent(in) :: on_u(LBi:,LBj:)
159 real(r8), intent(in) :: pm(LBi:,LBj:)
160 real(r8), intent(in) :: pn(LBi:,LBj:)
161 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
162 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
163 real(r8), intent(in) :: pden(LBi:,LBj:,:)
164# ifdef TS_MIX_CLIMA
165 real(r8), intent(in) :: tclm(LBi:,LBj:,:,:)
166# endif
167# ifdef DIAGNOSTICS_TS
168 real(r8), intent(inout) :: DiaTwrk(LBi:,LBj:,:,:,:)
169# endif
170 real(r8), intent(inout) :: t(LBi:,LBj:,:,:,:)
171#else
172# ifdef MASKING
173 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
174 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
175# endif
176# ifdef WET_DRY
177 real(r8), intent(in) :: umask_wet(LBi:UBi,LBj:UBj)
178 real(r8), intent(in) :: vmask_wet(LBi:UBi,LBj:UBj)
179# endif
180# ifdef DIFF_3DCOEF
181# ifdef TS_U3ADV_SPLIT
182 real(r8), intent(in) :: diff3d_u(LBi:UBi,LBj:UBj,N(ng))
183 real(r8), intent(in) :: diff3d_v(LBi:UBi,LBj:UBj,N(ng))
184# else
185 real(r8), intent(in) :: diff3d_r(LBi:UBi,LBj:UBj,N(ng))
186# endif
187# else
188 real(r8), intent(in) :: diff4(LBi:UBi,LBj:UBj,NT(ng))
189# endif
190 real(r8), intent(in) :: om_v(LBi:UBi,LBj:UBj)
191 real(r8), intent(in) :: on_u(LBi:UBi,LBj:UBj)
192 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
193 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
194 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
195 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
196 real(r8), intent(in) :: pden(LBi:UBi,LBj:UBj,N(ng))
197# ifdef TS_MIX_CLIMA
198 real(r8), intent(in) :: tclm(LBi:UBi,LBj:UBj,N(ng),NT(ng))
199# endif
200# ifdef DIAGNOSTICS_TS
201 real(r8), intent(inout) :: DiaTwrk(LBi:UBi,LBj:UBj,N(ng),NT(ng), &
202 & NDT)
203# endif
204 real(r8), intent(inout) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
205#endif
206!
207! Local variable declarations.
208!
209 integer :: Imin, Imax, Jmin, Jmax
210 integer :: i, itrc, j, k, k1, k2
211
212 real(r8), parameter :: eps = 0.5_r8
213 real(r8), parameter :: small = 1.0e-14_r8
214 real(r8), parameter :: slope_max = 0.0001_r8
215 real(r8), parameter :: strat_min = 0.1_r8
216
217 real(r8) :: cff, cff1, cff2, cff3, cff4, dife, difx
218
219 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,N(ng)) :: LapT
220
221 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FE
222 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FX
223
224 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: FS
225 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dRde
226 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dRdx
227 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dTde
228 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dTdr
229 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dTdx
230
231#include "set_bounds.h"
232!
233!-----------------------------------------------------------------------
234! Compute horizontal biharmonic diffusion along isopycnic surfaces.
235! The biharmonic operator is computed by applying the harmonic
236! operator twice.
237!-----------------------------------------------------------------------
238!
239! Set local I- and J-ranges.
240!
241 IF (ewperiodic(ng)) THEN
242 imin=istr-1
243 imax=iend+1
244 ELSE
245 imin=max(istr-1,1)
246 imax=min(iend+1,lm(ng))
247 END IF
248 IF (nsperiodic(ng)) THEN
249 jmin=jstr-1
250 jmax=jend+1
251 ELSE
252 jmin=max(jstr-1,1)
253 jmax=min(jend+1,mm(ng))
254 END IF
255!
256! Compute horizontal and density gradients. Notice the recursive
257! blocking sequence. The vertical placement of the gradients is:
258!
259! dTdx,dTde(:,:,k1) k rho-points
260! dTdx,dTde(:,:,k2) k+1 rho-points
261! FS,dTdr(:,:,k1) k-1/2 W-points
262! FS,dTdr(:,:,k2) k+1/2 W-points
263!
264#ifdef TS_MIX_STABILITY
265! In order to increase stability, the biharmonic operator is applied
266! as: 3/4 t(:,:,:,nrhs,:) + 1/4 t(:,:,:,nstp,:).
267!
268#endif
269
270 t_loop : DO itrc=1,nt(ng)
271 k2=1
272 k_loop1 : DO k=0,n(ng)
273 k1=k2
274 k2=3-k1
275 IF (k.lt.n(ng)) THEN
276 DO j=jmin,jmax
277 DO i=imin,imax+1
278 cff=0.5_r8*(pm(i,j)+pm(i-1,j))
279#ifdef MASKING
280 cff=cff*umask(i,j)
281#endif
282#ifdef WET_DRY
283 cff=cff*umask_wet(i,j)
284#endif
285 drdx(i,j,k2)=cff*(pden(i ,j,k+1)- &
286 & pden(i-1,j,k+1))
287#if defined TS_MIX_STABILITY
288 dtdx(i,j,k2)=cff*(0.75_r8*(t(i ,j,k+1,nrhs,itrc)- &
289 & t(i-1,j,k+1,nrhs,itrc))+ &
290 & 0.25_r8*(t(i ,j,k+1,nstp,itrc)- &
291 & t(i-1,j,k+1,nstp,itrc)))
292#elif defined TS_MIX_CLIMA
293 IF (ltracerclm(itrc,ng)) THEN
294 dtdx(i,j,k2)=cff*((t(i ,j,k+1,nrhs,itrc)- &
295 & tclm(i ,j,k+1,itrc))- &
296 & (t(i-1,j,k+1,nrhs,itrc)- &
297 & tclm(i-1,j,k+1,itrc)))
298 ELSE
299 dtdx(i,j,k2)=cff*(t(i ,j,k+1,nrhs,itrc)- &
300 & t(i-1,j,k+1,nrhs,itrc))
301 END IF
302#else
303 dtdx(i,j,k2)=cff*(t(i ,j,k+1,nrhs,itrc)- &
304 & t(i-1,j,k+1,nrhs,itrc))
305#endif
306 END DO
307 END DO
308 DO j=jmin,jmax+1
309 DO i=imin,imax
310 cff=0.5_r8*(pn(i,j)+pn(i,j-1))
311#ifdef MASKING
312 cff=cff*vmask(i,j)
313#endif
314#ifdef WET_DRY
315 cff=cff*vmask_wet(i,j)
316#endif
317 drde(i,j,k2)=cff*(pden(i,j ,k+1)- &
318 & pden(i,j-1,k+1))
319#if defined TS_MIX_STABILITY
320 dtde(i,j,k2)=cff*(0.75_r8*(t(i,j ,k+1,nrhs,itrc)- &
321 & t(i,j-1,k+1,nrhs,itrc))+ &
322 & 0.25_r8*(t(i,j ,k+1,nstp,itrc)- &
323 & t(i,j-1,k+1,nstp,itrc)))
324#elif defined TS_MIX_CLIMA
325 IF (ltracerclm(itrc,ng)) THEN
326 dtde(i,j,k2)=cff*((t(i,j ,k+1,nrhs,itrc)- &
327 & tclm(i,j ,k+1,itrc))- &
328 & (t(i,j-1,k+1,nrhs,itrc)- &
329 & tclm(i,j-1,k+1,itrc)))
330 ELSE
331 dtde(i,j,k2)=cff*(t(i,j ,k+1,nrhs,itrc)- &
332 & t(i,j-1,k+1,nrhs,itrc))
333 END IF
334#else
335 dtde(i,j,k2)=cff*(t(i,j ,k+1,nrhs,itrc)- &
336 & t(i,j-1,k+1,nrhs,itrc))
337#endif
338 END DO
339 END DO
340 END IF
341 IF ((k.eq.0).or.(k.eq.n(ng))) THEN
342 DO j=jmin-1,jmax+1
343 DO i=imin-1,imax+1
344 dtdr(i,j,k2)=0.0_r8
345 fs(i,j,k2)=0.0_r8
346 END DO
347 END DO
348 ELSE
349 DO j=jmin-1,jmax+1
350 DO i=imin-1,imax+1
351#if defined TS_MIX_MAX_SLOPE
352 cff1=sqrt(drdx(i,j,k2)**2+drdx(i+1,j,k2)**2+ &
353 & drdx(i,j,k1)**2+drdx(i+1,j,k1)**2+ &
354 & drde(i,j,k2)**2+drde(i,j+1,k2)**2+ &
355 & drde(i,j,k1)**2+drde(i,j+1,k1)**2)
356 cff2=0.25_r8*slope_max* &
357 & (z_r(i,j,k+1)-z_r(i,j,k))*cff1
358 cff3=max(pden(i,j,k)-pden(i,j,k+1),small)
359 cff4=max(cff2,cff3)
360 cff=-1.0_r8/cff4
361#elif defined TS_MIX_MIN_STRAT
362 cff1=max(pden(i,j,k)-pden(i,j,k+1), &
363 & strat_min*(z_r(i,j,k+1)-z_r(i,j,k)))
364 cff=-1.0_r8/cff1
365#else
366 cff1=max(pden(i,j,k)-pden(i,j,k+1),eps)
367 cff=-1.0_r8/cff1
368#endif
369#if defined TS_MIX_STABILITY
370 dtdr(i,j,k2)=cff*(0.75_r8*(t(i,j,k+1,nrhs,itrc)- &
371 & t(i,j,k ,nrhs,itrc))+ &
372 & 0.25_r8*(t(i,j,k+1,nstp,itrc)- &
373 & t(i,j,k ,nstp,itrc)))
374#elif defined TS_MIX_CLIMA
375 IF (ltracerclm(itrc,ng)) THEN
376 dtdr(i,j,k2)=cff*((t(i,j,k+1,nrhs,itrc)- &
377 & tclm(i,j,k+1,itrc))- &
378 & (t(i,j,k ,nrhs,itrc)- &
379 & tclm(i,j,k ,itrc)))
380 ELSE
381 dtdr(i,j,k2)=cff*(t(i,j,k+1,nrhs,itrc)- &
382 & t(i,j,k ,nrhs,itrc))
383 END IF
384#else
385 dtdr(i,j,k2)=cff*(t(i,j,k+1,nrhs,itrc)- &
386 & t(i,j,k ,nrhs,itrc))
387#endif
388 fs(i,j,k2)=cff*(z_r(i,j,k+1)- &
389 & z_r(i,j,k ))
390 END DO
391 END DO
392 END IF
393 IF (k.gt.0) THEN
394 DO j=jmin,jmax
395 DO i=imin,imax+1
396#ifdef DIFF_3DCOEF
397# ifdef TS_U3ADV_SPLIT
398 cff=0.5_r8*diff3d_u(i,j,k)*on_u(i,j)
399# else
400 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i-1,j,k))* &
401 & on_u(i,j)
402# endif
403#else
404 cff=0.25_r8*(diff4(i,j,itrc)+diff4(i-1,j,itrc))* &
405 & on_u(i,j)
406#endif
407 fx(i,j)=cff* &
408 & (hz(i,j,k)+hz(i-1,j,k))* &
409 & (dtdx(i,j,k1)- &
410 & 0.5_r8*(max(drdx(i,j,k1),0.0_r8)* &
411 & (dtdr(i-1,j,k1)+ &
412 & dtdr(i ,j,k2))+ &
413 & min(drdx(i,j,k1),0.0_r8)* &
414 & (dtdr(i-1,j,k2)+ &
415 & dtdr(i ,j,k1))))
416 END DO
417 END DO
418 DO j=jmin,jmax+1
419 DO i=imin,imax
420#ifdef DIFF_3DCOEF
421# ifdef TS_U3ADV_SPLIT
422 cff=0.5_r8*diff3d_v(i,j,k)*om_v(i,j)
423# else
424 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i,j-1,k))* &
425 & om_v(i,j)
426# endif
427#else
428 cff=0.25_r8*(diff4(i,j,itrc)+diff4(i,j-1,itrc))* &
429 & om_v(i,j)
430#endif
431 fe(i,j)=cff* &
432 & (hz(i,j,k)+hz(i,j-1,k))* &
433 & (dtde(i,j,k1)- &
434 & 0.5_r8*(max(drde(i,j,k1),0.0_r8)* &
435 & (dtdr(i,j-1,k1)+ &
436 & dtdr(i,j ,k2))+ &
437 & min(drde(i,j,k1),0.0_r8)* &
438 & (dtdr(i,j-1,k2)+ &
439 & dtdr(i,j ,k1))))
440 END DO
441 END DO
442 IF (k.lt.n(ng)) THEN
443 DO j=jmin,jmax
444 DO i=imin,imax
445#ifdef DIFF_3DCOEF
446# ifdef TS_U3ADV_SPLIT
447 difx=0.125_r8*(diff3d_u(i,j,k )+diff3d_u(i+1,j,k )+ &
448 & diff3d_u(i,j,k+1)+diff3d_u(i+1,j,k+1))
449 dife=0.125_r8*(diff3d_v(i,j,k )+diff3d_v(i,j+1,k )+ &
450 & diff3d_v(i,j,k+1)+diff3d_v(i,j+1,k+1))
451# else
452 difx=0.5_r8*diff3d_r(i,j,k)
453 dife=difx
454# endif
455#else
456 difx=0.5_r8*diff4(i,j,itrc)
457 dife=difx
458#endif
459 cff1=max(drdx(i ,j,k1),0.0_r8)
460 cff2=max(drdx(i+1,j,k2),0.0_r8)
461 cff3=min(drdx(i ,j,k2),0.0_r8)
462 cff4=min(drdx(i+1,j,k1),0.0_r8)
463 cff=difx* &
464 & (cff1*(cff1*dtdr(i,j,k2)-dtdx(i ,j,k1))+ &
465 & cff2*(cff2*dtdr(i,j,k2)-dtdx(i+1,j,k2))+ &
466 & cff3*(cff3*dtdr(i,j,k2)-dtdx(i ,j,k2))+ &
467 & cff4*(cff4*dtdr(i,j,k2)-dtdx(i+1,j,k1)))
468!
469 cff1=max(drde(i,j ,k1),0.0_r8)
470 cff2=max(drde(i,j+1,k2),0.0_r8)
471 cff3=min(drde(i,j ,k2),0.0_r8)
472 cff4=min(drde(i,j+1,k1),0.0_r8)
473 cff=cff+ &
474 & dife* &
475 & (cff1*(cff1*dtdr(i,j,k2)-dtde(i,j ,k1))+ &
476 & cff2*(cff2*dtdr(i,j,k2)-dtde(i,j+1,k2))+ &
477 & cff3*(cff3*dtdr(i,j,k2)-dtde(i,j ,k2))+ &
478 & cff4*(cff4*dtdr(i,j,k2)-dtde(i,j+1,k1)))
479 fs(i,j,k2)=cff*fs(i,j,k2)
480 END DO
481 END DO
482 END IF
483!
484! Compute first harmonic operator, without mixing coefficient.
485! Multiply by the metrics of the second harmonic operator. Save
486! into work array "LapT".
487!
488 DO j=jmin,jmax
489 DO i=imin,imax
490 cff=pm(i,j)*pn(i,j)
491 cff1=1.0_r8/hz(i,j,k)
492 lapt(i,j,k)=cff1*(cff* &
493 & (fx(i+1,j)-fx(i,j)+ &
494 & fe(i,j+1)-fe(i,j))+ &
495 & (fs(i,j,k2)-fs(i,j,k1)))
496 END DO
497 END DO
498 END IF
499 END DO k_loop1
500!
501! Apply boundary conditions (except periodic; closed or gradient)
502! to the first harmonic operator.
503!
504 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
505 IF (domain(ng)%Western_Edge(tile)) THEN
506 IF (lbc(iwest,istvar(itrc),ng)%closed) THEN
507 DO k=1,n(ng)
508 DO j=jmin,jmax
509 lapt(istr-1,j,k)=0.0_r8
510 END DO
511 END DO
512 ELSE
513 DO k=1,n(ng)
514 DO j=jmin,jmax
515 lapt(istr-1,j,k)=lapt(istr,j,k)
516 END DO
517 END DO
518 END IF
519 END IF
520 END IF
521!
522 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
523 IF (domain(ng)%Eastern_Edge(tile)) THEN
524 IF (lbc(ieast,istvar(itrc),ng)%closed) THEN
525 DO k=1,n(ng)
526 DO j=jmin,jmax
527 lapt(iend+1,j,k)=0.0_r8
528 END DO
529 END DO
530 ELSE
531 DO k=1,n(ng)
532 DO j=jmin,jmax
533 lapt(iend+1,j,k)=lapt(iend,j,k)
534 END DO
535 END DO
536 END IF
537 END IF
538 END IF
539!
540 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
541 IF (domain(ng)%Southern_Edge(tile)) THEN
542 IF (lbc(isouth,istvar(itrc),ng)%closed) THEN
543 DO k=1,n(ng)
544 DO i=imin,imax
545 lapt(i,jstr-1,k)=0.0_r8
546 END DO
547 END DO
548 ELSE
549 DO k=1,n(ng)
550 DO i=imin,imax
551 lapt(i,jstr-1,k)=lapt(i,jstr,k)
552 END DO
553 END DO
554 END IF
555 END IF
556 END IF
557!
558 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
559 IF (domain(ng)%Northern_Edge(tile)) THEN
560 IF (lbc(inorth,istvar(itrc),ng)%closed) THEN
561 DO k=1,n(ng)
562 DO i=imin,imax
563 lapt(i,jend+1,k)=0.0_r8
564 END DO
565 END DO
566 ELSE
567 DO k=1,n(ng)
568 DO i=imin,imax
569 lapt(i,jend+1,k)=lapt(i,jend,k)
570 END DO
571 END DO
572 END IF
573 END IF
574 END IF
575!
576 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng).or. &
577 & compositegrid(iwest ,ng).or.ewperiodic(ng))) THEN
578 IF (domain(ng)%SouthWest_Corner(tile)) THEN
579 DO k=1,n(ng)
580 lapt(istr-1,jstr-1,k)=0.5_r8* &
581 & (lapt(istr ,jstr-1,k)+ &
582 & lapt(istr-1,jstr ,k))
583 END DO
584 END IF
585 END IF
586
587 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng).or. &
588 & compositegrid(ieast ,ng).or.ewperiodic(ng))) THEN
589 IF (domain(ng)%SouthEast_Corner(tile)) THEN
590 DO k=1,n(ng)
591 lapt(iend+1,jstr-1,k)=0.5_r8* &
592 & (lapt(iend ,jstr-1,k)+ &
593 & lapt(iend+1,jstr ,k))
594 END DO
595 END IF
596 END IF
597
598 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng).or. &
599 & compositegrid(iwest ,ng).or.ewperiodic(ng))) THEN
600 IF (domain(ng)%NorthWest_Corner(tile)) THEN
601 DO k=1,n(ng)
602 lapt(istr-1,jend+1,k)=0.5_r8* &
603 & (lapt(istr ,jend+1,k)+ &
604 & lapt(istr-1,jend ,k))
605 END DO
606 END IF
607 END IF
608
609 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng).or. &
610 & compositegrid(ieast ,ng).or.ewperiodic(ng))) THEN
611 IF (domain(ng)%NorthEast_Corner(tile)) THEN
612 DO k=1,n(ng)
613 lapt(iend+1,jend+1,k)=0.5_r8* &
614 & (lapt(iend ,jend+1,k)+ &
615 & lapt(iend+1,jend ,k))
616 END DO
617 END IF
618 END IF
619!
620! Compute horizontal and density gradients associated with the
621! second rotated harmonic operator.
622!
623 k2=1
624 k_loop2: DO k=0,n(ng)
625 k1=k2
626 k2=3-k1
627 IF (k.lt.n(ng)) THEN
628 DO j=jstr,jend
629 DO i=istr,iend+1
630 cff=0.5_r8*(pm(i,j)+pm(i-1,j))
631#ifdef MASKING
632 cff=cff*umask(i,j)
633#endif
634#ifdef WET_DRY
635 cff=cff*umask_wet(i,j)
636#endif
637 drdx(i,j,k2)=cff*(pden(i ,j,k+1)- &
638 & pden(i-1,j,k+1))
639 dtdx(i,j,k2)=cff*(lapt(i ,j,k+1)- &
640 & lapt(i-1,j,k+1))
641 END DO
642 END DO
643 DO j=jstr,jend+1
644 DO i=istr,iend
645 cff=0.5_r8*(pn(i,j)+pn(i,j-1))
646#ifdef MASKING
647 cff=cff*vmask(i,j)
648#endif
649#ifdef WET_DRY
650 cff=cff*vmask_wet(i,j)
651#endif
652 drde(i,j,k2)=cff*(pden(i,j ,k+1)- &
653 & pden(i,j-1,k+1))
654 dtde(i,j,k2)=cff*(lapt(i,j ,k+1)- &
655 & lapt(i,j-1,k+1))
656 END DO
657 END DO
658 END IF
659 IF ((k.eq.0).or.(k.eq.n(ng))) THEN
660 DO j=jstr-1,jend+1
661 DO i=istr-1,iend+1
662 dtdr(i,j,k2)=0.0_r8
663 fs(i,j,k2)=0.0_r8
664 END DO
665 END DO
666 ELSE
667 DO j=jstr-1,jend+1
668 DO i=istr-1,iend+1
669#if defined TS_MIX_MAX_SLOPE
670 cff1=sqrt(drdx(i,j,k2)**2+drdx(i+1,j,k2)**2+ &
671 & drdx(i,j,k1)**2+drdx(i+1,j,k1)**2+ &
672 & drde(i,j,k2)**2+drde(i,j+1,k2)**2+ &
673 & drde(i,j,k1)**2+drde(i,j+1,k1)**2)
674 cff2=0.25_r8*slope_max* &
675 & (z_r(i,j,k+1)-z_r(i,j,k))*cff1
676 cff3=max(pden(i,j,k)-pden(i,j,k+1),small)
677 cff4=max(cff2,cff3)
678 cff=-1.0_r8/cff4
679#elif defined TS_MIX_MIN_STRAT
680 cff1=max(pden(i,j,k)-pden(i,j,k+1), &
681 & strat_min*(z_r(i,j,k+1)-z_r(i,j,k)))
682 cff=-1.0_r8/cff1
683#else
684 cff1=max(pden(i,j,k)-pden(i,j,k+1),eps)
685 cff=-1.0_r8/cff1
686#endif
687 dtdr(i,j,k2)=cff*(lapt(i,j,k+1)- &
688 & lapt(i,j,k ))
689 fs(i,j,k2)=cff*(z_r(i,j,k+1)- &
690 & z_r(i,j,k ))
691 END DO
692 END DO
693 END IF
694!
695! Compute components of the rotated tracer flux (T m4/s) along
696! isopycnic surfaces.
697!
698 IF (k.gt.0) THEN
699 DO j=jstr,jend
700 DO i=istr,iend+1
701#ifdef DIFF_3DCOEF
702# ifdef TS_U3ADV_SPLIT
703 cff=0.5_r8*diff3d_u(i,j,k)*on_u(i,j)
704# else
705 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i-1,j,k))* &
706 & on_u(i,j)
707# endif
708#else
709 cff=0.25_r8*(diff4(i,j,itrc)+diff4(i-1,j,itrc))* &
710 & on_u(i,j)
711#endif
712 fx(i,j)=cff* &
713 & (hz(i,j,k)+hz(i-1,j,k))* &
714 & (dtdx(i,j,k1)- &
715 & 0.5_r8*(max(drdx(i,j,k1),0.0_r8)* &
716 & (dtdr(i-1,j,k1)+ &
717 & dtdr(i ,j,k2))+ &
718 & min(drdx(i,j,k1),0.0_r8)* &
719 & (dtdr(i-1,j,k2)+ &
720 & dtdr(i ,j,k1))))
721 END DO
722 END DO
723 DO j=jstr,jend+1
724 DO i=istr,iend
725#ifdef DIFF_3DCOEF
726# ifdef TS_U3ADV_SPLIT
727 cff=0.5_r8*diff3d_v(i,j,k)*om_v(i,j)
728# else
729 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i,j-1,k))* &
730 & om_v(i,j)
731# endif
732#else
733 cff=0.25_r8*(diff4(i,j,itrc)+diff4(i,j-1,itrc))* &
734 & om_v(i,j)
735#endif
736 fe(i,j)=cff* &
737 & (hz(i,j,k)+hz(i,j-1,k))* &
738 & (dtde(i,j,k1)- &
739 & 0.5_r8*(max(drde(i,j,k1),0.0_r8)* &
740 & (dtdr(i,j-1,k1)+ &
741 & dtdr(i,j ,k2))+ &
742 & min(drde(i,j,k1),0.0_r8)* &
743 & (dtdr(i,j-1,k2)+ &
744 & dtdr(i,j ,k1))))
745 END DO
746 END DO
747 IF (k.lt.n(ng)) THEN
748 DO j=jstr,jend
749 DO i=istr,iend
750#ifdef DIFF_3DCOEF
751# ifdef TS_U3ADV_SPLIT
752 difx=0.125_r8*(diff3d_u(i,j,k )+diff3d_u(i+1,j,k )+ &
753 & diff3d_u(i,j,k+1)+diff3d_u(i+1,j,k+1))
754 dife=0.125_r8*(diff3d_v(i,j,k )+diff3d_v(i,j+1,k )+ &
755 & diff3d_v(i,j,k+1)+diff3d_v(i,j+1,k+1))
756# else
757 difx=0.5_r8*diff3d_r(i,j,k)
758 dife=difx
759# endif
760#else
761 difx=0.5_r8*diff4(i,j,itrc)
762 dife=difx
763#endif
764 cff1=max(drdx(i ,j,k1),0.0_r8)
765 cff2=max(drdx(i+1,j,k2),0.0_r8)
766 cff3=min(drdx(i ,j,k2),0.0_r8)
767 cff4=min(drdx(i+1,j,k1),0.0_r8)
768 cff=difx* &
769 & (cff1*(cff1*dtdr(i,j,k2)-dtdx(i ,j,k1))+ &
770 & cff2*(cff2*dtdr(i,j,k2)-dtdx(i+1,j,k2))+ &
771 & cff3*(cff3*dtdr(i,j,k2)-dtdx(i ,j,k2))+ &
772 & cff4*(cff4*dtdr(i,j,k2)-dtdx(i+1,j,k1)))
773!
774 cff1=max(drde(i,j ,k1),0.0_r8)
775 cff2=max(drde(i,j+1,k2),0.0_r8)
776 cff3=min(drde(i,j ,k2),0.0_r8)
777 cff4=min(drde(i,j+1,k1),0.0_r8)
778 cff=cff+ &
779 & dife* &
780 & (cff1*(cff1*dtdr(i,j,k2)-dtde(i,j ,k1))+ &
781 & cff2*(cff2*dtdr(i,j,k2)-dtde(i,j+1,k2))+ &
782 & cff3*(cff3*dtdr(i,j,k2)-dtde(i,j ,k2))+ &
783 & cff4*(cff4*dtdr(i,j,k2)-dtde(i,j+1,k1)))
784 fs(i,j,k2)=cff*fs(i,j,k2)
785 END DO
786 END DO
787 END IF
788!
789! Time-step biharmonic, isopycnal diffusion term (m Tunits).
790!
791 DO j=jstr,jend
792 DO i=istr,iend
793 cff=dt(ng)*pm(i,j)*pn(i,j)
794 cff1=cff*(fx(i+1,j )-fx(i,j))
795 cff2=cff*(fe(i ,j+1)-fe(i,j))
796 cff3=dt(ng)*(fs(i,j,k2)-fs(i,j,k1))
797 cff4=cff1+cff2+cff3
798 t(i,j,k,nnew,itrc)=t(i,j,k,nnew,itrc)-cff4
799#ifdef DIAGNOSTICS_TS
800 diatwrk(i,j,k,itrc,itxdif)=-cff1
801 diatwrk(i,j,k,itrc,itydif)=-cff2
802 diatwrk(i,j,k,itrc,itsdif)=-cff3
803 diatwrk(i,j,k,itrc,ithdif)=-cff4
804#endif
805 END DO
806 END DO
807 END IF
808 END DO k_loop2
809 END DO t_loop
810!
811 RETURN

References mod_scalars::compositegrid, mod_param::domain, mod_scalars::dt, mod_scalars::ewperiodic, mod_scalars::ieast, mod_scalars::inorth, mod_scalars::isouth, mod_ncparam::istvar, mod_scalars::ithdif, mod_scalars::itsdif, mod_scalars::itxdif, mod_scalars::itydif, mod_scalars::iwest, mod_param::lbc, mod_param::lm, mod_scalars::ltracerclm, mod_param::mm, and mod_scalars::nsperiodic.

◆ t3dmix4_s_tile()

subroutine t3dmix4_mod::t3dmix4_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) nstp,
integer, intent(in) nnew,
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) umask_wet,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) vmask_wet,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) hz,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) pmon_u,
real(r8), dimension(lbi:ubi,lbj:ubj), intent(in) pnom_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) diff3d_u,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) diff3d_v,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng)), intent(in) diff3d_r,
real(r8), dimension(lbi:ubi,lbj:ubj,nt(ng)), intent(in) diff4,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),nt(ng)), intent(in) tclm,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),nt(ng), ndt), intent(inout) diatwrk,
real(r8), dimension(lbi:ubi,lbj:ubj,n(ng),3,nt(ng)), intent(inout) t )
private

Definition at line 94 of file t3dmix4_s.h.

121!***********************************************************************
122!
123 USE mod_param
124 USE mod_ncparam
125 USE mod_scalars
126!
127! Imported variable declarations.
128!
129 integer, intent(in) :: ng, tile
130 integer, intent(in) :: LBi, UBi, LBj, UBj
131 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
132 integer, intent(in) :: nrhs, nstp, nnew
133
134#ifdef ASSUMED_SHAPE
135# ifdef MASKING
136 real(r8), intent(in) :: umask(LBi:,LBj:)
137 real(r8), intent(in) :: vmask(LBi:,LBj:)
138# endif
139# ifdef WET_DRY
140 real(r8), intent(in) :: umask_wet(LBi:,LBj:)
141 real(r8), intent(in) :: vmask_wet(LBi:,LBj:)
142# endif
143# ifdef DIFF_3DCOEF
144# ifdef TS_U3ADV_SPLIT
145 real(r8), intent(in) :: diff3d_u(LBi:,LBj:,:)
146 real(r8), intent(in) :: diff3d_v(LBi:,LBj:,:)
147# else
148 real(r8), intent(in) :: diff3d_r(LBi:,LBj:,:)
149# endif
150# else
151 real(r8), intent(in) :: diff4(LBi:,LBj:,:)
152# endif
153 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
154 real(r8), intent(in) :: pmon_u(LBi:,LBj:)
155 real(r8), intent(in) :: pnom_v(LBi:,LBj:)
156 real(r8), intent(in) :: pm(LBi:,LBj:)
157 real(r8), intent(in) :: pn(LBi:,LBj:)
158# ifdef TS_MIX_CLIMA
159 real(r8), intent(in) :: tclm(LBi:,LBj:,:,:)
160# endif
161# ifdef DIAGNOSTICS_TS
162 real(r8), intent(inout) :: DiaTwrk(LBi:,LBj:,:,:,:)
163# endif
164 real(r8), intent(inout) :: t(LBi:,LBj:,:,:,:)
165#else
166# ifdef MASKING
167 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
168 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
169# endif
170# ifdef WET_DRY
171 real(r8), intent(in) :: umask_wet(LBi:UBi,LBj:UBj)
172 real(r8), intent(in) :: vmask_wet(LBi:UBi,LBj:UBj)
173# endif
174# ifdef DIFF_3DCOEF
175# ifdef TS_U3ADV_SPLIT
176 real(r8), intent(in) :: diff3d_u(LBi:UBi,LBj:UBj,N(ng))
177 real(r8), intent(in) :: diff3d_v(LBi:UBi,LBj:UBj,N(ng))
178# else
179 real(r8), intent(in) :: diff3d_r(LBi:UBi,LBj:UBj,N(ng))
180# endif
181# else
182 real(r8), intent(in) :: diff4(LBi:UBi,LBj:UBj,NT(ng))
183# endif
184 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
185 real(r8), intent(in) :: pmon_u(LBi:UBi,LBj:UBj)
186 real(r8), intent(in) :: pnom_v(LBi:UBi,LBj:UBj)
187 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
188 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
189# ifdef TS_MIX_CLIMA
190 real(r8), intent(in) :: tclm(LBi:UBi,LBj:UBj,N(ng),NT(ng))
191# endif
192# ifdef DIAGNOSTICS_TS
193 real(r8), intent(inout) :: DiaTwrk(LBi:UBi,LBj:UBj,N(ng),NT(ng), &
194 & NDT)
195# endif
196 real(r8), intent(inout) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
197#endif
198!
199! Local variable declarations.
200!
201 integer :: Imin, Imax, Jmin, Jmax
202 integer :: i, itrc, j, k
203
204 real(r8) :: cff, cff1, cff2, cff3
205
206 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FE
207 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FX
208 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: LapT
209
210#include "set_bounds.h"
211!
212!-----------------------------------------------------------------------
213! Compute horizontal biharmonic diffusion along constant S-surfaces.
214! The biharmonic operator is computed by applying the harmonic
215! operator twice.
216#ifdef TS_MIX_STABILITY
217! In order to increase stability, the biharmonic operator is applied
218! as: 3/4 t(:,:,:,nrhs,:) + 1/4 t(:,:,:,nstp,:).
219#endif
220!-----------------------------------------------------------------------
221!
222! Set local I- and J-ranges.
223!
224 IF (ewperiodic(ng)) THEN
225 imin=istr-1
226 imax=iend+1
227 ELSE
228 imin=max(istr-1,1)
229 imax=min(iend+1,lm(ng))
230 END IF
231 IF (nsperiodic(ng)) THEN
232 jmin=jstr-1
233 jmax=jend+1
234 ELSE
235 jmin=max(jstr-1,1)
236 jmax=min(jend+1,mm(ng))
237 END IF
238!
239! Compute horizontal tracer flux in the XI- and ETA-directions.
240!
241 DO itrc=1,nt(ng)
242 DO k=1,n(ng)
243 DO j=jmin,jmax
244 DO i=imin,imax+1
245#ifdef DIFF_3DCOEF
246# ifdef TS_U3ADV_SPLIT
247 cff=0.5_r8*diff3d_u(i,j,k)*pmon_u(i,j)
248# else
249 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i-1,j,k))* &
250 & pmon_u(i,j)
251# endif
252#else
253 cff=0.25_r8*(diff4(i,j,itrc)+diff4(i-1,j,itrc))* &
254 & pmon_u(i,j)
255#endif
256#ifdef MASKING
257 cff=cff*umask(i,j)
258#endif
259#ifdef WET_DRY
260 cff=cff*umask_wet(i,j)
261#endif
262#if defined TS_MIX_STABILITY
263 fx(i,j)=cff* &
264 & (hz(i,j,k)+hz(i-1,j,k))* &
265 & (0.75_r8*(t(i ,j,k,nrhs,itrc)- &
266 & t(i-1,j,k,nrhs,itrc))+ &
267 & 0.25_r8*(t(i ,j,k,nstp,itrc)- &
268 & t(i-1,j,k,nstp,itrc)))
269#elif defined TS_MIX_CLIMA
270 IF (ltracerclm(itrc,ng)) THEN
271 fx(i,j)=cff* &
272 & (hz(i,j,k)+hz(i-1,j,k))* &
273 & ((t(i ,j,k,nrhs,itrc)-tclm(i ,j,k,itrc))- &
274 & (t(i-1,j,k,nrhs,itrc)-tclm(i-1,j,k,itrc)))
275 ELSE
276 fx(i,j)=cff* &
277 & (hz(i,j,k)+hz(i-1,j,k))* &
278 & (t(i ,j,k,nrhs,itrc)- &
279 & t(i-1,j,k,nrhs,itrc))
280 END IF
281#else
282 fx(i,j)=cff* &
283 & (hz(i,j,k)+hz(i-1,j,k))* &
284 & (t(i ,j,k,nrhs,itrc)- &
285 & t(i-1,j,k,nrhs,itrc))
286#endif
287 END DO
288 END DO
289 DO j=jmin,jmax+1
290 DO i=imin,imax
291#ifdef DIFF_3DCOEF
292# ifdef TS_U3ADV_SPLIT
293 cff=0.5_r8*diff3d_v(i,j,k)*pnom_v(i,j)
294# else
295 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i,j-1,k))* &
296 & pnom_v(i,j)
297# endif
298#else
299 cff=0.25_r8*(diff4(i,j,itrc)+diff4(i,j-1,itrc))* &
300 & pnom_v(i,j)
301#endif
302#ifdef MASKING
303 cff=cff*vmask(i,j)
304#endif
305#ifdef WET_DRY
306 cff=cff*vmask_wet(i,j)
307#endif
308#if defined TS_MIX_STABILITY
309 fe(i,j)=cff* &
310 & (hz(i,j,k)+hz(i,j-1,k))* &
311 & (0.75_r8*(t(i,j ,k,nrhs,itrc)- &
312 & t(i,j-1,k,nrhs,itrc))+ &
313 & 0.25_r8*(t(i,j ,k,nstp,itrc)- &
314 & t(i,j-1,k,nstp,itrc)))
315#elif defined TS_MIX_CLIMA
316 IF (ltracerclm(itrc,ng)) THEN
317 fe(i,j)=cff* &
318 & (hz(i,j,k)+hz(i,j-1,k))* &
319 & ((t(i,j ,k,nrhs,itrc)-tclm(i,j ,k,itrc))- &
320 & (t(i,j-1,k,nrhs,itrc)-tclm(i,j-1,k,itrc)))
321 ELSE
322 fe(i,j)=cff* &
323 & (hz(i,j,k)+hz(i,j-1,k))* &
324 & (t(i,j ,k,nrhs,itrc)- &
325 & t(i,j-1,k,nrhs,itrc))
326 END IF
327#else
328 fe(i,j)=cff* &
329 & (hz(i,j,k)+hz(i,j-1,k))* &
330 & (t(i,j ,k,nrhs,itrc)- &
331 & t(i,j-1,k,nrhs,itrc))
332#endif
333 END DO
334 END DO
335!
336! Compute first harmonic operator and multiply by the metrics of the
337! second harmonic operator.
338!
339 DO j=jmin,jmax
340 DO i=imin,imax
341 cff=1.0_r8/hz(i,j,k)
342 lapt(i,j)=pm(i,j)*pn(i,j)*cff* &
343 & (fx(i+1,j)-fx(i,j)+ &
344 & fe(i,j+1)-fe(i,j))
345 END DO
346 END DO
347!
348! Apply boundary conditions (except periodic; closed or gradient)
349! to the first harmonic operator.
350!
351 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
352 IF (domain(ng)%Western_Edge(tile)) THEN
353 IF (lbc(iwest,istvar(itrc),ng)%closed) THEN
354 DO j=jmin,jmax
355 lapt(istr-1,j)=0.0_r8
356 END DO
357 ELSE
358 DO j=jmin,jmax
359 lapt(istr-1,j)=lapt(istr,j)
360 END DO
361 END IF
362 END IF
363 END IF
364!
365 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
366 IF (domain(ng)%Eastern_Edge(tile)) THEN
367 IF (lbc(ieast,istvar(itrc),ng)%closed) THEN
368 DO j=jmin,jmax
369 lapt(iend+1,j)=0.0_r8
370 END DO
371 ELSE
372 DO j=jmin,jmax
373 lapt(iend+1,j)=lapt(iend,j)
374 END DO
375 END IF
376 END IF
377 END IF
378!
379 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
380 IF (domain(ng)%Southern_Edge(tile)) THEN
381 IF (lbc(isouth,istvar(itrc),ng)%closed) THEN
382 DO i=imin,imax
383 lapt(i,jstr-1)=0.0_r8
384 END DO
385 ELSE
386 DO i=imin,imax
387 lapt(i,jstr-1)=lapt(i,jstr)
388 END DO
389 END IF
390 END IF
391 END IF
392!
393 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
394 IF (domain(ng)%Northern_Edge(tile)) THEN
395 IF (lbc(inorth,istvar(itrc),ng)%closed) THEN
396 DO i=imin,imax
397 lapt(i,jend+1)=0.0_r8
398 END DO
399 ELSE
400 DO i=imin,imax
401 lapt(i,jend+1)=lapt(i,jend)
402 END DO
403 END IF
404 END IF
405 END IF
406!
407! Compute FX=d(LapT)/d(xi) and FE=d(LapT)/d(eta) terms.
408!
409 DO j=jstr,jend
410 DO i=istr,iend+1
411#ifdef DIFF_3DCOEF
412# ifdef TS_U3ADV_SPLIT
413 cff=0.5_r8*diff3d_u(i,j,k)*pmon_u(i,j)
414# else
415 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i-1,j,k))* &
416 & pmon_u(i,j)
417# endif
418#else
419 cff=0.25_r8*(diff4(i,j,itrc)+diff4(i-1,j,itrc))* &
420 & pmon_u(i,j)
421#endif
422 fx(i,j)=cff* &
423 & (hz(i,j,k)+hz(i-1,j,k))* &
424 & (lapt(i,j)-lapt(i-1,j))
425#ifdef MASKING
426 fx(i,j)=fx(i,j)*umask(i,j)
427#endif
428#ifdef WET_DRY
429 fx(i,j)=fx(i,j)*umask_wet(i,j)
430#endif
431 END DO
432 END DO
433 DO j=jstr,jend+1
434 DO i=istr,iend
435#ifdef DIFF_3DCOEF
436# ifdef TS_U3ADV_SPLIT
437 cff=0.5_r8*diff3d_v(i,j,k)*pnom_v(i,j)
438# else
439 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i,j-1,k))* &
440 & pnom_v(i,j)
441# endif
442#else
443 cff=0.25_r8*(diff4(i,j,itrc)+diff4(i,j-1,itrc))* &
444 & pnom_v(i,j)
445#endif
446 fe(i,j)=cff* &
447 & (hz(i,j,k)+hz(i,j-1,k))* &
448 & (lapt(i,j)-lapt(i,j-1))
449#ifdef MASKING
450 fe(i,j)=fe(i,j)*vmask(i,j)
451#endif
452#ifdef WET_DRY
453 fe(i,j)=fe(i,j)*vmask_wet(i,j)
454#endif
455 END DO
456 END DO
457!
458! Time-step biharmonic, S-surfaces diffusion term (m Tunits).
459!
460 DO j=jstr,jend
461 DO i=istr,iend
462 cff=dt(ng)*pm(i,j)*pn(i,j)
463 cff1=cff*(fx(i+1,j )-fx(i,j))
464 cff2=cff*(fe(i ,j+1)-fe(i,j))
465 cff3=cff1+cff2
466 t(i,j,k,nnew,itrc)=t(i,j,k,nnew,itrc)-cff3
467#ifdef DIAGNOSTICS_TS
468 diatwrk(i,j,k,itrc,itxdif)=-cff1
469 diatwrk(i,j,k,itrc,itydif)=-cff2
470 diatwrk(i,j,k,itrc,ithdif)=-cff3
471#endif
472 END DO
473 END DO
474 END DO
475 END DO
476!
477 RETURN

References mod_scalars::compositegrid, mod_param::domain, mod_scalars::dt, mod_scalars::ewperiodic, mod_scalars::ieast, mod_scalars::inorth, mod_scalars::isouth, mod_ncparam::istvar, mod_scalars::ithdif, mod_scalars::itxdif, mod_scalars::itydif, mod_scalars::iwest, mod_param::lbc, mod_param::lm, mod_scalars::ltracerclm, mod_param::mm, and mod_scalars::nsperiodic.