ROMS
Loading...
Searching...
No Matches
t3dmix2_iso.h
Go to the documentation of this file.
1 MODULE t3dmix2_mod
2!
3!git $Id$
4!=======================================================================
5! Copyright (c) 2002-2025 The ROMS Group !
6! Licensed under a MIT/X style license !
7! See License_ROMS.md Hernan G. Arango !
8!========================================== Alexander F. Shchepetkin ===
9! !
10! This subroutine computes horizontal harmonic mixing of tracers !
11! along isopycnic surfaces. !
12! !
13!=======================================================================
14!
15 implicit none
16!
17 PRIVATE
18 PUBLIC t3dmix2
19!
20 CONTAINS
21!
22!***********************************************************************
23 SUBROUTINE t3dmix2 (ng, tile)
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, 26, __line__, myfile)
51#endif
52 CALL t3dmix2_iso_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 & mixing(ng) % diff3d_r, &
72#else
73 & mixing(ng) % diff2, &
74#endif
75 & ocean(ng) % pden, &
76#ifdef TS_MIX_CLIMA
77 & clima(ng) % tclm, &
78#endif
79#ifdef DIAGNOSTICS_TS
80 & diags(ng) % DiaTwrk, &
81#endif
82 & ocean(ng) % t)
83#ifdef PROFILE
84 CALL wclock_off (ng, inlm, 26, __line__, myfile)
85#endif
86!
87 RETURN
88 END SUBROUTINE t3dmix2
89!
90!***********************************************************************
91 SUBROUTINE t3dmix2_iso_tile (ng, tile, &
92 & LBi, UBi, LBj, UBj, &
93 & IminS, ImaxS, JminS, JmaxS, &
94 & nrhs, nstp, nnew, &
95#ifdef MASKING
96 & umask, vmask, &
97#endif
98#ifdef WET_DRY
99 & umask_wet, vmask_wet, &
100#endif
101 & om_v, on_u, pm, pn, &
102 & Hz, z_r, &
103#ifdef DIFF_3DCOEF
104 & diff3d_r, &
105#else
106 & diff2, &
107#endif
108 & pden, &
109#ifdef TS_MIX_CLIMA
110 & tclm, &
111#endif
112#ifdef DIAGNOSTICS_TS
113 & DiaTwrk, &
114#endif
115 & t)
116!***********************************************************************
117!
118 USE mod_param
119 USE mod_scalars
120!
121! Imported variable declarations.
122!
123 integer, intent(in) :: ng, tile
124 integer, intent(in) :: LBi, UBi, LBj, UBj
125 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
126 integer, intent(in) :: nrhs, nstp, nnew
127
128#ifdef ASSUMED_SHAPE
129# ifdef MASKING
130 real(r8), intent(in) :: umask(LBi:,LBj:)
131 real(r8), intent(in) :: vmask(LBi:,LBj:)
132# endif
133# ifdef WET_DRY
134 real(r8), intent(in) :: umask_wet(LBi:,LBj:)
135 real(r8), intent(in) :: vmask_wet(LBi:,LBj:)
136# endif
137# ifdef DIFF_3DCOEF
138 real(r8), intent(in) :: diff3d_r(LBi:,LBj:,:)
139# else
140 real(r8), intent(in) :: diff2(LBi:,LBj:,:)
141# endif
142 real(r8), intent(in) :: om_v(LBi:,LBj:)
143 real(r8), intent(in) :: on_u(LBi:,LBj:)
144 real(r8), intent(in) :: pm(LBi:,LBj:)
145 real(r8), intent(in) :: pn(LBi:,LBj:)
146 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
147 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
148 real(r8), intent(in) :: pden(LBi:,LBj:,:)
149# ifdef TS_MIX_CLIMA
150 real(r8), intent(in) :: tclm(LBi:,LBj:,:,:)
151# endif
152# ifdef DIAGNOSTICS_TS
153 real(r8), intent(inout) :: DiaTwrk(LBi:,LBj:,:,:,:)
154# endif
155 real(r8), intent(inout) :: t(LBi:,LBj:,:,:,:)
156#else
157# ifdef MASKING
158 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
159 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
160# endif
161# ifdef WET_DRY
162 real(r8), intent(in) :: umask_wet(LBi:UBi,LBj:UBj)
163 real(r8), intent(in) :: vmask_wet(LBi:UBi,LBj:UBj)
164# endif
165# ifdef DIFF_3DCOEF
166 real(r8), intent(in) :: diff3d_r(LBi:UBi,LBj:UBj,N(ng))
167# else
168 real(r8), intent(in) :: diff2(LBi:UBi,LBj:UBj,NT(ng))
169# endif
170 real(r8), intent(in) :: om_v(LBi:UBi,LBj:UBj)
171 real(r8), intent(in) :: on_u(LBi:UBi,LBj:UBj)
172 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
173 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
174 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
175 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
176 real(r8), intent(in) :: pden(LBi:UBi,LBj:UBj,N(ng))
177# ifdef TS_MIX_CLIMA
178 real(r8), intent(in) :: tclm(LBi:UBi,LBj:UBj,N(ng),NT(ng))
179# endif
180# ifdef DIAGNOSTICS_TS
181 real(r8), intent(inout) :: DiaTwrk(LBi:UBi,LBj:UBj,N(ng),NT(ng), &
182 & NDT)
183# endif
184 real(r8), intent(inout) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
185#endif
186!
187! Local variable declarations.
188!
189 integer :: i, itrc, j, k, k1, k2
190
191 real(r8), parameter :: eps = 0.5_r8
192 real(r8), parameter :: small = 1.0e-14_r8
193 real(r8), parameter :: slope_max = 0.0001_r8
194 real(r8), parameter :: strat_min = 0.1_r8
195
196 real(r8) :: cff, cff1, cff2, cff3, cff4
197
198 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FE
199 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FX
200
201 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: FS
202 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dRde
203 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dRdx
204 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dTde
205 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dTdr
206 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,2) :: dTdx
207
208#include "set_bounds.h"
209!
210!-----------------------------------------------------------------------
211! Compute horizontal harmonic diffusion along isopycnic surfaces.
212!-----------------------------------------------------------------------
213!
214! Compute horizontal and density gradients. Notice the recursive
215! blocking sequence. The vertical placement of the gradients is:
216!
217! dTdx,dTde(:,:,k1) k rho-points
218! dTdx,dTde(:,:,k2) k+1 rho-points
219! FS,dTdr(:,:,k1) k-1/2 W-points
220! FS,dTdr(:,:,k2) k+1/2 W-points
221!
222 t_loop : DO itrc=1,nt(ng)
223 k2=1
224 k_loop : DO k=0,n(ng)
225 k1=k2
226 k2=3-k1
227 IF (k.lt.n(ng)) THEN
228 DO j=jstr,jend
229 DO i=istr,iend+1
230 cff=0.5_r8*(pm(i,j)+pm(i-1,j))
231#ifdef MASKING
232 cff=cff*umask(i,j)
233#endif
234#ifdef WET_DRY
235 cff=cff*umask_wet(i,j)
236#endif
237 drdx(i,j,k2)=cff*(pden(i ,j,k+1)- &
238 & pden(i-1,j,k+1))
239#if defined TS_MIX_STABILITY
240 dtdx(i,j,k2)=cff*(0.75_r8*(t(i ,j,k+1,nrhs,itrc)- &
241 & t(i-1,j,k+1,nrhs,itrc))+ &
242 & 0.25_r8*(t(i ,j,k+1,nstp,itrc)- &
243 & t(i-1,j,k+1,nstp,itrc)))
244#elif defined TS_MIX_CLIMA
245 IF (ltracerclm(itrc,ng)) THEN
246 dtdx(i,j,k2)=cff*((t(i ,j,k+1,nrhs,itrc)- &
247 & tclm(i ,j,k+1,itrc))- &
248 & (t(i-1,j,k+1,nrhs,itrc)- &
249 & tclm(i-1,j,k+1,itrc)))
250 ELSE
251 dtdx(i,j,k2)=cff*(t(i ,j,k+1,nrhs,itrc)- &
252 & t(i-1,j,k+1,nrhs,itrc))
253 END IF
254#else
255 dtdx(i,j,k2)=cff*(t(i ,j,k+1,nrhs,itrc)- &
256 & t(i-1,j,k+1,nrhs,itrc))
257#endif
258 END DO
259 END DO
260 DO j=jstr,jend+1
261 DO i=istr,iend
262 cff=0.5_r8*(pn(i,j)+pn(i,j-1))
263#ifdef MASKING
264 cff=cff*vmask(i,j)
265#endif
266#ifdef WET_DRY
267 cff=cff*vmask_wet(i,j)
268#endif
269 drde(i,j,k2)=cff*(pden(i,j ,k+1)- &
270 & pden(i,j-1,k+1))
271#if defined TS_MIX_STABILITY
272 dtde(i,j,k2)=cff*(0.75_r8*(t(i,j ,k+1,nrhs,itrc)- &
273 & t(i,j-1,k+1,nrhs,itrc))+ &
274 & 0.25_r8*(t(i,j ,k+1,nstp,itrc)- &
275 & t(i,j-1,k+1,nstp,itrc)))
276#elif defined TS_MIX_CLIMA
277 IF (ltracerclm(itrc,ng)) THEN
278 dtde(i,j,k2)=cff*((t(i,j ,k+1,nrhs,itrc)- &
279 & tclm(i,j ,k+1,itrc))- &
280 & (t(i,j-1,k+1,nrhs,itrc)- &
281 & tclm(i,j-1,k+1,itrc)))
282 ELSE
283 dtde(i,j,k2)=cff*(t(i,j ,k+1,nrhs,itrc)- &
284 & t(i,j-1,k+1,nrhs,itrc))
285 END IF
286#else
287 dtde(i,j,k2)=cff*(t(i,j ,k+1,nrhs,itrc)- &
288 & t(i,j-1,k+1,nrhs,itrc))
289#endif
290 END DO
291 END DO
292 END IF
293 IF ((k.eq.0).or.(k.eq.n(ng))) THEN
294 DO j=jstr-1,jend+1
295 DO i=istr-1,iend+1
296 dtdr(i,j,k2)=0.0_r8
297 fs(i,j,k2)=0.0_r8
298 END DO
299 END DO
300 ELSE
301 DO j=jstr-1,jend+1
302 DO i=istr-1,iend+1
303#if defined TS_MIX_MAX_SLOPE
304 cff1=sqrt(drdx(i,j,k2)**2+drdx(i+1,j,k2)**2+ &
305 & drdx(i,j,k1)**2+drdx(i+1,j,k1)**2+ &
306 & drde(i,j,k2)**2+drde(i,j+1,k2)**2+ &
307 & drde(i,j,k1)**2+drde(i,j+1,k1)**2)
308 cff2=0.25_r8*slope_max* &
309 & (z_r(i,j,k+1)-z_r(i,j,k))*cff1
310 cff3=max(pden(i,j,k)-pden(i,j,k+1),small)
311 cff4=max(cff2,cff3)
312 cff=-1.0_r8/cff4
313#elif defined TS_MIX_MIN_STRAT
314 cff1=max(pden(i,j,k)-pden(i,j,k+1), &
315 & strat_min*(z_r(i,j,k+1)-z_r(i,j,k)))
316 cff=-1.0_r8/cff1
317#else
318 cff1=max(pden(i,j,k)-pden(i,j,k+1),eps)
319 cff=-1.0_r8/cff1
320#endif
321#if defined TS_MIX_STABILITY
322 dtdr(i,j,k2)=cff*(0.75_r8*(t(i,j,k+1,nrhs,itrc)- &
323 & t(i,j,k ,nrhs,itrc))+ &
324 & 0.25_r8*(t(i,j,k+1,nstp,itrc)- &
325 & t(i,j,k ,nstp,itrc)))
326#elif defined TS_MIX_CLIMA
327 IF (ltracerclm(itrc,ng)) THEN
328 dtdr(i,j,k2)=cff*((t(i,j,k+1,nrhs,itrc)- &
329 & tclm(i,j,k+1,itrc))- &
330 & (t(i,j,k ,nrhs,itrc)- &
331 & tclm(i,j,k ,itrc)))
332 ELSE
333 dtdr(i,j,k2)=cff*(t(i,j,k+1,nrhs,itrc)- &
334 & t(i,j,k ,nrhs,itrc))
335 END IF
336#else
337 dtdr(i,j,k2)=cff*(t(i,j,k+1,nrhs,itrc)- &
338 & t(i,j,k ,nrhs,itrc))
339#endif
340 fs(i,j,k2)=cff*(z_r(i,j,k+1)-z_r(i,j,k))
341 END DO
342 END DO
343 END IF
344!
345! Compute components of the rotated tracer flux (T m4/s) along
346! isopycnic surfaces.
347!
348 IF (k.gt.0) THEN
349 DO j=jstr,jend
350 DO i=istr,iend+1
351#ifdef DIFF_3DCOEF
352 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i-1,j,k))* &
353 & on_u(i,j)
354#else
355 cff=0.25_r8*(diff2(i,j,itrc)+diff2(i-1,j,itrc))* &
356 & on_u(i,j)
357#endif
358 fx(i,j)=cff* &
359 & (hz(i,j,k)+hz(i-1,j,k))* &
360 & (dtdx(i,j,k1)- &
361 & 0.5_r8*(max(drdx(i,j,k1),0.0_r8)* &
362 & (dtdr(i-1,j,k1)+ &
363 & dtdr(i ,j,k2))+ &
364 & min(drdx(i,j,k1),0.0_r8)* &
365 & (dtdr(i-1,j,k2)+ &
366 & dtdr(i ,j,k1))))
367 END DO
368 END DO
369 DO j=jstr,jend+1
370 DO i=istr,iend
371#ifdef DIFF_3DCOEF
372 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i,j-1,k))* &
373 & om_v(i,j)
374#else
375 cff=0.25_r8*(diff2(i,j,itrc)+diff2(i,j-1,itrc))* &
376 & om_v(i,j)
377#endif
378 fe(i,j)=cff* &
379 & (hz(i,j,k)+hz(i,j-1,k))* &
380 & (dtde(i,j,k1)- &
381 & 0.5_r8*(max(drde(i,j,k1),0.0_r8)* &
382 & (dtdr(i,j-1,k1)+ &
383 & dtdr(i,j ,k2))+ &
384 & min(drde(i,j,k1),0.0_r8)* &
385 & (dtdr(i,j-1,k2)+ &
386 & dtdr(i,j ,k1))))
387 END DO
388 END DO
389 IF (k.lt.n(ng)) THEN
390 DO j=jstr,jend
391 DO i=istr,iend
392 cff1=max(drdx(i ,j,k1),0.0_r8)
393 cff2=max(drdx(i+1,j,k2),0.0_r8)
394 cff3=min(drdx(i ,j,k2),0.0_r8)
395 cff4=min(drdx(i+1,j,k1),0.0_r8)
396 cff=cff1*(cff1*dtdr(i,j,k2)-dtdx(i ,j,k1))+ &
397 & cff2*(cff2*dtdr(i,j,k2)-dtdx(i+1,j,k2))+ &
398 & cff3*(cff3*dtdr(i,j,k2)-dtdx(i ,j,k2))+ &
399 & cff4*(cff4*dtdr(i,j,k2)-dtdx(i+1,j,k1))
400 cff1=max(drde(i,j ,k1),0.0_r8)
401 cff2=max(drde(i,j+1,k2),0.0_r8)
402 cff3=min(drde(i,j ,k2),0.0_r8)
403 cff4=min(drde(i,j+1,k1),0.0_r8)
404 cff=cff+ &
405 & cff1*(cff1*dtdr(i,j,k2)-dtde(i,j ,k1))+ &
406 & cff2*(cff2*dtdr(i,j,k2)-dtde(i,j+1,k2))+ &
407 & cff3*(cff3*dtdr(i,j,k2)-dtde(i,j ,k2))+ &
408 & cff4*(cff4*dtdr(i,j,k2)-dtde(i,j+1,k1))
409#ifdef DIFF_3DCOEF
410 fs(i,j,k2)=0.5_r8*cff*diff3d_r(i,j,k)*fs(i,j,k2)
411#else
412 fs(i,j,k2)=0.5_r8*cff*diff2(i,j,itrc)*fs(i,j,k2)
413#endif
414 END DO
415 END DO
416 END IF
417!
418! Time-step harmonic, isopycnic diffusion term (m Tunits).
419!
420 DO j=jstr,jend
421 DO i=istr,iend
422 cff=dt(ng)*pm(i,j)*pn(i,j)
423 cff1=cff*(fx(i+1,j )-fx(i,j))
424 cff2=cff*(fe(i ,j+1)-fe(i,j))
425 cff3=dt(ng)*(fs(i,j,k2)-fs(i,j,k1))
426 cff4=cff1+cff2+cff3
427 t(i,j,k,nnew,itrc)=t(i,j,k,nnew,itrc)+cff4
428#ifdef DIAGNOSTICS_TS
429 diatwrk(i,j,k,itrc,itxdif)=cff1
430 diatwrk(i,j,k,itrc,itydif)=cff2
431 diatwrk(i,j,k,itrc,itsdif)=cff3
432 diatwrk(i,j,k,itrc,ithdif)=cff4
433#endif
434 END DO
435 END DO
436 END IF
437 END DO k_loop
438 END DO t_loop
439!
440 RETURN
441 END SUBROUTINE t3dmix2_iso_tile
442
443 END MODULE t3dmix2_mod
type(t_clima), dimension(:), allocatable clima
Definition mod_clima.F:153
type(t_diags), dimension(:), allocatable diags
Definition mod_diags.F:100
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_mixing), dimension(:), allocatable mixing
Definition mod_mixing.F:399
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
integer, parameter inlm
Definition mod_param.F:662
real(dp), dimension(:), allocatable dt
integer itxdif
integer itsdif
integer itydif
logical, dimension(:,:), allocatable ltracerclm
integer ithdif
integer, dimension(:), allocatable nrhs
integer, dimension(:), allocatable nnew
integer, dimension(:), allocatable nstp
subroutine t3dmix2_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_r, diff2, pden, tclm, diatwrk, t)
subroutine, public t3dmix2(ng, tile)
Definition t3dmix2_geo.h:24
recursive subroutine wclock_off(ng, model, region, line, routine)
Definition timers.F:148
recursive subroutine wclock_on(ng, model, region, line, routine)
Definition timers.F:3