ROMS
Loading...
Searching...
No Matches
t3dmix2_s.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 S-coordinate levels 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, 24, __line__, myfile)
51#endif
52 CALL t3dmix2_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 & mixing(ng) % diff3d_r, &
71#else
72 & mixing(ng) % diff2, &
73#endif
74#ifdef TS_MIX_CLIMA
75 & clima(ng) % tclm, &
76#endif
77#ifdef DIAGNOSTICS_TS
78 & diags(ng) % DiaTwrk, &
79#endif
80 & ocean(ng) % t)
81#ifdef PROFILE
82 CALL wclock_off (ng, inlm, 24, __line__, myfile)
83#endif
84!
85 RETURN
86 END SUBROUTINE t3dmix2
87!
88!***********************************************************************
89 SUBROUTINE t3dmix2_s_tile (ng, tile, &
90 & LBi, UBi, LBj, UBj, &
91 & IminS, ImaxS, JminS, JmaxS, &
92 & nrhs, nstp, nnew, &
93#ifdef MASKING
94 & umask, vmask, &
95#endif
96#ifdef WET_DRY
97 & umask_wet, vmask_wet, &
98#endif
99 & Hz, pmon_u, pnom_v, pm, pn, &
100#ifdef DIFF_3DCOEF
101 & diff3d_r, &
102#else
103 & diff2, &
104#endif
105#ifdef TS_MIX_CLIMA
106 & tclm, &
107#endif
108#ifdef DIAGNOSTICS_TS
109 & DiaTwrk, &
110#endif
111 & t)
112!***********************************************************************
113!
114 USE mod_param
115 USE mod_scalars
116!
117! Imported variable declarations.
118!
119 integer, intent(in) :: ng, tile
120 integer, intent(in) :: LBi, UBi, LBj, UBj
121 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
122 integer, intent(in) :: nrhs, nstp, nnew
123
124#ifdef ASSUMED_SHAPE
125# ifdef MASKING
126 real(r8), intent(in) :: umask(LBi:,LBj:)
127 real(r8), intent(in) :: vmask(LBi:,LBj:)
128# endif
129# ifdef WET_DRY
130 real(r8), intent(in) :: umask_wet(LBi:,LBj:)
131 real(r8), intent(in) :: vmask_wet(LBi:,LBj:)
132# endif
133# ifdef DIFF_3DCOEF
134 real(r8), intent(in) :: diff3d_r(LBi:,LBj:,:)
135# else
136 real(r8), intent(in) :: diff2(LBi:,LBj:,:)
137# endif
138 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
139 real(r8), intent(in) :: pmon_u(LBi:,LBj:)
140 real(r8), intent(in) :: pnom_v(LBi:,LBj:)
141 real(r8), intent(in) :: pm(LBi:,LBj:)
142 real(r8), intent(in) :: pn(LBi:,LBj:)
143# ifdef TS_MIX_CLIMA
144 real(r8), intent(in) :: tclm(LBi:,LBj:,:,:)
145# endif
146# ifdef DIAGNOSTICS_TS
147 real(r8), intent(inout) :: DiaTwrk(LBi:,LBj:,:,:,:)
148# endif
149 real(r8), intent(inout) :: t(LBi:,LBj:,:,:,:)
150#else
151# ifdef MASKING
152 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
153 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
154# endif
155# ifdef WET_DRY
156 real(r8), intent(in) :: umask_wet(LBi:UBi,LBj:UBj)
157 real(r8), intent(in) :: vmask_wet(LBi:UBi,LBj:UBj)
158# endif
159# ifdef DIFF_3DCOEF
160 real(r8), intent(in) :: diff3d_r(LBi:UBi,LBj:UBj,N(ng))
161# else
162 real(r8), intent(in) :: diff2(LBi:UBi,LBj:UBj,NT(ng))
163# endif
164 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
165 real(r8), intent(in) :: pmon_u(LBi:UBi,LBj:UBj)
166 real(r8), intent(in) :: pnom_v(LBi:UBi,LBj:UBj)
167 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
168 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
169# ifdef TS_MIX_CLIMA
170 real(r8), intent(in) :: tclm(LBi:UBi,LBj:UBj,N(ng),NT(ng))
171# endif
172# ifdef DIAGNOSTICS_TS
173 real(r8), intent(inout) :: DiaTwrk(LBi:UBi,LBj:UBj,N(ng),NT(ng), &
174 & NDT)
175# endif
176 real(r8), intent(inout) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
177#endif
178!
179! Local variable declarations.
180!
181 integer :: i, itrc, j, k
182
183 real(r8) :: cff, cff1, cff2, cff3
184
185 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FE
186 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FX
187
188#include "set_bounds.h"
189!
190!-----------------------------------------------------------------------
191! Compute horizontal harmonic diffusion along constant S-surfaces.
192#ifdef TS_MIX_STABILITY
193! In order to increase stability, the harmonic operator is applied
194! as: 3/4 t(:,:,:,nrhs,:) + 1/4 t(:,:,:,nstp,:).
195#endif
196!-----------------------------------------------------------------------
197!
198 DO itrc=1,nt(ng)
199 DO k=1,n(ng)
200!
201! Compute XI- and ETA-components of diffusive tracer flux (T m3/s).
202!
203 DO j=jstr,jend
204 DO i=istr,iend+1
205#ifdef DIFF_3DCOEF
206 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i-1,j,k))* &
207 & pmon_u(i,j)
208#else
209 cff=0.25_r8*(diff2(i,j,itrc)+diff2(i-1,j,itrc))* &
210 & pmon_u(i,j)
211#endif
212#if defined TS_MIX_STABILITY
213 fx(i,j)=cff* &
214 & (hz(i,j,k)+hz(i-1,j,k))* &
215 & (0.75_r8*(t(i ,j,k,nrhs,itrc)- &
216 & t(i-1,j,k,nrhs,itrc))+ &
217 & 0.25_r8*(t(i ,j,k,nstp,itrc)- &
218 & t(i-1,j,k,nstp,itrc)))
219#elif defined TS_MIX_CLIMA
220 IF (ltracerclm(itrc,ng)) THEN
221 fx(i,j)=cff* &
222 & (hz(i,j,k)+hz(i-1,j,k))* &
223 & ((t(i ,j,k,nrhs,itrc)-tclm(i ,j,k,itrc))- &
224 & (t(i-1,j,k,nrhs,itrc)-tclm(i-1,j,k,itrc)))
225 ELSE
226 fx(i,j)=cff* &
227 & (hz(i,j,k)+hz(i-1,j,k))* &
228 & (t(i,j,k,nrhs,itrc)-t(i-1,j,k,nrhs,itrc))
229 END IF
230#else
231 fx(i,j)=cff* &
232 & (hz(i,j,k)+hz(i-1,j,k))* &
233 & (t(i,j,k,nrhs,itrc)-t(i-1,j,k,nrhs,itrc))
234#endif
235#ifdef MASKING
236 fx(i,j)=fx(i,j)*umask(i,j)
237#endif
238#ifdef WET_DRY
239 fx(i,j)=fx(i,j)*umask_wet(i,j)
240#endif
241 END DO
242 END DO
243 DO j=jstr,jend+1
244 DO i=istr,iend
245#ifdef DIFF_3DCOEF
246 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i,j-1,k))* &
247 & pnom_v(i,j)
248#else
249 cff=0.25_r8*(diff2(i,j,itrc)+diff2(i,j-1,itrc))* &
250 & pnom_v(i,j)
251#endif
252#if defined TS_MIX_STABILITY
253 fe(i,j)=cff* &
254 & (hz(i,j,k)+hz(i,j-1,k))* &
255 & (0.75_r8*(t(i,j ,k,nrhs,itrc)- &
256 & t(i,j-1,k,nrhs,itrc))+ &
257 & 0.25_r8*(t(i,j ,k,nstp,itrc)- &
258 & t(i,j-1,k,nstp,itrc)))
259#elif defined TS_MIX_CLIMA
260 IF (ltracerclm(itrc,ng)) THEN
261 fe(i,j)=cff* &
262 & (hz(i,j,k)+hz(i,j-1,k))* &
263 & ((t(i,j ,k,nrhs,itrc)-tclm(i,j ,k,itrc))- &
264 & (t(i,j-1,k,nrhs,itrc)-tclm(i,j-1,k,itrc)))
265 ELSE
266 fe(i,j)=cff* &
267 & (hz(i,j,k)+hz(i,j-1,k))* &
268 & (t(i,j,k,nrhs,itrc)-t(i,j-1,k,nrhs,itrc))
269 END IF
270#else
271 fe(i,j)=cff* &
272 & (hz(i,j,k)+hz(i,j-1,k))* &
273 & (t(i,j,k,nrhs,itrc)-t(i,j-1,k,nrhs,itrc))
274#endif
275#ifdef MASKING
276 fe(i,j)=fe(i,j)*vmask(i,j)
277#endif
278#ifdef WET_DRY
279 fe(i,j)=fe(i,j)*vmask_wet(i,j)
280#endif
281 END DO
282 END DO
283!
284! Time-step harmonic, S-surfaces diffusion term (m Tunits).
285!
286 DO j=jstr,jend
287 DO i=istr,iend
288 cff=dt(ng)*pm(i,j)*pn(i,j)
289 cff1=cff*(fx(i+1,j )-fx(i,j))
290 cff2=cff*(fe(i ,j+1)-fe(i,j))
291 cff3=cff1+cff2
292 t(i,j,k,nnew,itrc)=t(i,j,k,nnew,itrc)+cff3
293#ifdef DIAGNOSTICS_TS
294 diatwrk(i,j,k,itrc,itxdif)=cff1
295 diatwrk(i,j,k,itrc,itydif)=cff2
296 diatwrk(i,j,k,itrc,ithdif)=cff3
297#endif
298 END DO
299 END DO
300 END DO
301 END DO
302!
303 RETURN
304 END SUBROUTINE t3dmix2_s_tile
305
306 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 itydif
logical, dimension(:,:), allocatable ltracerclm
integer ithdif
integer, dimension(:), allocatable nrhs
integer, dimension(:), allocatable nnew
integer, dimension(:), allocatable nstp
subroutine t3dmix2_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_r, diff2, tclm, diatwrk, t)
Definition t3dmix2_s.h:112
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