ROMS
Loading...
Searching...
No Matches
t3dmix4_s.h
Go to the documentation of this file.
1 MODULE t3dmix4_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 biharmonic mixing of tracers !
11! along S-coordinate levels surfaces. !
12! !
13!=======================================================================
14!
15 implicit none
16!
17 PRIVATE
18 PUBLIC t3dmix4
19!
20 CONTAINS
21!
22!***********************************************************************
23 SUBROUTINE t3dmix4 (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, 27, __line__, myfile)
51#endif
52 CALL t3dmix4_s_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) % Hz, &
65 & grid(ng) % pmon_u, &
66 & grid(ng) % pnom_v, &
67 & grid(ng) % pm, &
68 & grid(ng) % pn, &
69#ifdef DIFF_3DCOEF
70# ifdef TS_U3ADV_SPLIT
71 & mixing(ng) % diff3d_u, &
72 & mixing(ng) % diff3d_v, &
73# else
74 & mixing(ng) % diff3d_r, &
75# endif
76#else
77 & mixing(ng) % diff4, &
78#endif
79#ifdef TS_MIX_CLIMA
80 & clima(ng) % tclm, &
81#endif
82#ifdef DIAGNOSTICS_TS
83 & diags(ng) % DiaTwrk, &
84#endif
85 & ocean(ng) % t)
86#ifdef PROFILE
87 CALL wclock_off (ng, inlm, 27, __line__, myfile)
88#endif
89!
90 RETURN
91 END SUBROUTINE t3dmix4
92!
93!***********************************************************************
94 SUBROUTINE t3dmix4_s_tile (ng, tile, &
95 & LBi, UBi, LBj, UBj, &
96 & IminS, ImaxS, JminS, JmaxS, &
97 & nrhs, nstp, nnew, &
98#ifdef MASKING
99 & umask, vmask, &
100#endif
101#ifdef WET_DRY
102 & umask_wet, vmask_wet, &
103#endif
104 & Hz, pmon_u, pnom_v, pm, pn, &
105#ifdef DIFF_3DCOEF
106# ifdef TS_U3ADV_SPLIT
107 & diff3d_u, diff3d_v, &
108# else
109 & diff3d_r, &
110# endif
111#else
112 & diff4, &
113#endif
114#ifdef TS_MIX_CLIMA
115 & tclm, &
116#endif
117#ifdef DIAGNOSTICS_TS
118 & DiaTwrk, &
119#endif
120 & t)
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
478 END SUBROUTINE t3dmix4_s_tile
479
480 END MODULE t3dmix4_mod
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_ncparam
Definition mod_ncparam.F:2
mod_ncparam::istvar
integer, dimension(:), allocatable istvar
Definition mod_ncparam.F:520
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_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::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::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
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
t3dmix4_mod
Definition t3dmix4_geo.h:1
t3dmix4_mod::t3dmix4
subroutine, public t3dmix4(ng, tile)
Definition t3dmix4_geo.h:24
t3dmix4_mod::t3dmix4_s_tile
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)
Definition t3dmix4_s.h:121
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