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