ROMS
Loading...
Searching...
No Matches
rp_t3dmix2_s.h
Go to the documentation of this file.
1 MODULE rp_t3dmix2_mod
2!
3!git $Id$
4!================================================== Hernan G. Arango ===
5! Copyright (c) 2002-2025 The ROMS Group Andrew M. Moore !
6! Licensed under a MIT/X style license !
7! See License_ROMS.md !
8!=======================================================================
9! !
10! This subroutine computes representers tangent linear horizontal !
11! harmonic mixing of tracers along S-coordinate levels surfaces. !
12! !
13! BASIC STATE variables needed: t, Hz !
14! !
15!=======================================================================
16!
17 implicit none
18!
19 PRIVATE
20 PUBLIC rp_t3dmix2
21!
22 CONTAINS
23!
24!***********************************************************************
25 SUBROUTINE rp_t3dmix2 (ng, tile)
26!***********************************************************************
27!
28 USE mod_param
29#ifdef TS_MIX_CLIMA
30 USE mod_clima
31#endif
32#ifdef DIAGNOSTICS_TS
33!! USE mod_diags
34#endif
35 USE mod_grid
36 USE mod_mixing
37 USE mod_ocean
38 USE mod_stepping
39!
40! Imported variable declarations.
41!
42 integer, intent(in) :: ng, tile
43!
44! Local variable declarations.
45!
46 character (len=*), parameter :: MyFile = &
47 & __FILE__
48!
49#include "tile.h"
50!
51#ifdef PROFILE
52 CALL wclock_on (ng, irpm, 24, __line__, myfile)
53#endif
54 CALL rp_t3dmix2_s_tile (ng, tile, &
55 & lbi, ubi, lbj, ubj, &
56 & imins, imaxs, jmins, jmaxs, &
57 & nrhs(ng), nstp(ng), nnew(ng), &
58#ifdef MASKING
59 & grid(ng) % umask, &
60 & grid(ng) % vmask, &
61#endif
62#ifdef WET_DRY_NOT_YET
63 & grid(ng) % umask_wet, &
64 & grid(ng) % vmask_wet, &
65#endif
66 & grid(ng) % Hz, &
67 & grid(ng) % tl_Hz, &
68 & grid(ng) % pmon_u, &
69 & grid(ng) % pnom_v, &
70 & grid(ng) % pm, &
71 & grid(ng) % pn, &
72#ifdef DIFF_3DCOEF
73 & mixing(ng) % diff3d_r, &
74#else
75 & mixing(ng) % diff2, &
76#endif
77#ifdef TS_MIX_CLIMA
78 & clima(ng) % tclm, &
79#endif
80#ifdef DIAGNOSTICS_TS
81!! & DIAGS(ng) % DiaTwrk, &
82#endif
83 & ocean(ng) % t, &
84 & ocean(ng) % tl_t)
85#ifdef PROFILE
86 CALL wclock_off (ng, irpm, 24, __line__, myfile)
87#endif
88!
89 RETURN
90 END SUBROUTINE rp_t3dmix2
91!
92!***********************************************************************
93 SUBROUTINE rp_t3dmix2_s_tile (ng, tile, &
94 & LBi, UBi, LBj, UBj, &
95 & IminS, ImaxS, JminS, JmaxS, &
96 & nrhs, nstp, nnew, &
97#ifdef MASKING
98 & umask, vmask, &
99#endif
100#ifdef WET_DRY_NOT_YET
101 & umask_wet, vmask_wet, &
102#endif
103 & Hz, tl_Hz, &
104 & pmon_u, pnom_v, pm, pn, &
105#ifdef DIFF_3DCOEF
106 & diff3d_r, &
107#else
108 & diff2, &
109#endif
110#ifdef TS_MIX_CLIMA
111 & tclm, &
112#endif
113#ifdef DIAGNOSTICS_TS
114!! & DiaTwrk, &
115#endif
116 & t, tl_t)
117!***********************************************************************
118!
119 USE mod_param
120 USE mod_scalars
121!
122! Imported variable declarations.
123!
124 integer, intent(in) :: ng, tile
125 integer, intent(in) :: LBi, UBi, LBj, UBj
126 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
127 integer, intent(in) :: nrhs, nstp, nnew
128
129#ifdef ASSUMED_SHAPE
130# ifdef MASKING
131 real(r8), intent(in) :: umask(LBi:,LBj:)
132 real(r8), intent(in) :: vmask(LBi:,LBj:)
133# endif
134# ifdef WET_DRY_NOT_YET
135 real(r8), intent(in) :: umask_wet(LBi:,LBj:)
136 real(r8), intent(in) :: vmask_wet(LBi:,LBj:)
137# endif
138# ifdef DIFF_3DCOEF
139 real(r8), intent(in) :: diff3d_r(LBi:,LBj:,:)
140# else
141 real(r8), intent(in) :: diff2(LBi:,LBj:,:)
142# endif
143 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
144 real(r8), intent(in) :: tl_Hz(LBi:,LBj:,:)
145 real(r8), intent(in) :: pmon_u(LBi:,LBj:)
146 real(r8), intent(in) :: pnom_v(LBi:,LBj:)
147 real(r8), intent(in) :: pm(LBi:,LBj:)
148 real(r8), intent(in) :: pn(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(in) :: t(LBi:,LBj:,:,:,:)
156
157 real(r8), intent(inout) :: tl_t(LBi:,LBj:,:,:,:)
158#else
159# ifdef MASKING
160 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
161 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
162# endif
163# ifdef WET_DRY_NOT_YET
164 real(r8), intent(in) :: umask_wet(LBi:UBi,LBj:UBj)
165 real(r8), intent(in) :: vmask_wet(LBi:UBi,LBj:UBj)
166# endif
167# ifdef DIFF_3DCOEF
168 real(r8), intent(in) :: diff3d_r(LBi:UBi,LBj:UBj,N(ng))
169# else
170 real(r8), intent(in) :: diff2(LBi:UBi,LBj:UBj,NT(ng))
171# endif
172 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
173 real(r8), intent(in) :: tl_Hz(LBi:UBi,LBj:UBj,N(ng))
174 real(r8), intent(in) :: pmon_u(LBi:UBi,LBj:UBj)
175 real(r8), intent(in) :: pnom_v(LBi:UBi,LBj:UBj)
176 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
177 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
178# ifdef TS_MIX_CLIMA
179 real(r8), intent(in) :: tclm(LBi:UBi,LBj:UBj,N(ng),NT(ng))
180# endif
181# ifdef DIAGNOSTICS_TS
182!! real(r8), intent(inout) :: DiaTwrk(LBi:UBi,LBj:UBj,N(ng),NT(ng), &
183!! & NDT)
184# endif
185 real(r8), intent(in) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
186
187 real(r8), intent(inout) :: tl_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
188#endif
189!
190! Local variable declarations.
191!
192 integer :: i, itrc, j, k
193
194 real(r8) :: cff, cff1, tl_cff, tl_cff1
195
196 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_FE
197 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_FX
198
199#include "set_bounds.h"
200!
201!-----------------------------------------------------------------------
202! Compute tangent linear horizontal harmonic diffusion along constant
203! S-surfaces.
204!-----------------------------------------------------------------------
205!
206 DO itrc=1,nt(ng)
207 DO k=1,n(ng)
208!
209! Compute XI- and ETA-components of diffusive tracer flux (T m3/s).
210!
211 DO j=jstr,jend
212 DO i=istr,iend+1
213#ifdef DIFF_3DCOEF
214 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i-1,j,k))* &
215 & pmon_u(i,j)
216#else
217 cff=0.25_r8*(diff2(i,j,itrc)+diff2(i-1,j,itrc))* &
218 & pmon_u(i,j)
219#endif
220#if defined TS_MIX_STABILITY
221!^ FX(i,j)=cff* &
222!^ & (Hz(i,j,k)+Hz(i-1,j,k))* &
223!^ & (0.75_r8*(t(i ,j,k,nrhs,itrc)- &
224!^ & t(i-1,j,k,nrhs,itrc))+ &
225!^ & 0.25_r8*(t(i ,j,k,nstp,itrc)- &
226!^ & t(i-1,j,k,nstp,itrc)))
227!^
228 tl_fx(i,j)=cff* &
229 & ((tl_hz(i,j,k)+tl_hz(i-1,j,k))* &
230 & (0.75_r8*(t(i ,j,k,nrhs,itrc)- &
231 & t(i-1,j,k,nrhs,itrc))+ &
232 & 0.25_r8*(t(i ,j,k,nstp,itrc)- &
233 & t(i-1,j,k,nstp,itrc)))+ &
234 & (hz(i,j,k)+hz(i-1,j,k))* &
235 & (0.75_r8*(tl_t(i ,j,k,nrhs,itrc)- &
236 & tl_t(i-1,j,k,nrhs,itrc))+ &
237 & 0.25_r8*(tl_t(i ,j,k,nstp,itrc)- &
238 & tl_t(i-1,j,k,nstp,itrc))))- &
239# ifdef TL_IOMS
240 & cff* &
241 & (hz(i,j,k)+hz(i-1,j,k))* &
242 & (0.75_r8*(t(i ,j,k,nrhs,itrc)- &
243 & t(i-1,j,k,nrhs,itrc))+ &
244 & 0.25_r8*(t(i ,j,k,nstp,itrc)- &
245 & t(i-1,j,k,nstp,itrc)))
246# endif
247#elif defined TS_MIX_CLIMA
248 IF (ltracerclm(itrc,ng)) THEN
249!^ FX(i,j)=cff* &
250!^ & (Hz(i,j,k)+Hz(i-1,j,k))* &
251!^ & ((t(i ,j,k,nrhs,itrc)-tclm(i ,j,k,itrc))- &
252!^ & (t(i-1,j,k,nrhs,itrc)-tclm(i-1,j,k,itrc)))
253!^
254 tl_fx(i,j)=cff* &
255 & ((tl_hz(i,j,k)+tl_hz(i-1,j,k))* &
256 & ((t(i ,j,k,nrhs,itrc)- &
257 & tclm(i ,j,k,itrc))- &
258 & (t(i-1,j,k,nrhs,itrc)- &
259 & tclm(i-1,j,k,itrc)))+ &
260 & (hz(i,j,k)+hz(i-1,j,k))* &
261 & (tl_t(i ,j,k,nrhs,itrc)- &
262 & tl_t(i-1,j,k,nrhs,itrc)))- &
263# ifdef TL_IOMS
264 & cff* &
265 & (hz(i,j,k)+hz(i-1,j,k))* &
266 & ((t(i ,j,k,nrhs,itrc)-tclm(i ,j,k,itrc))- &
267 & (t(i-1,j,k,nrhs,itrc)-tclm(i-1,j,k,itrc)))
268# endif
269 ELSE
270!^ FX(i,j)=cff* &
271!^ & (Hz(i,j,k)+Hz(i-1,j,k))* &
272!^ & (t(i,j,k,nrhs,itrc)-t(i-1,j,k,nrhs,itrc))
273!^
274 tl_fx(i,j)=cff* &
275 & ((tl_hz(i,j,k)+tl_hz(i-1,j,k))* &
276 & (t(i ,j,k,nrhs,itrc)- &
277 & t(i-1,j,k,nrhs,itrc))+ &
278 & (hz(i,j,k)+hz(i-1,j,k))* &
279 & (tl_t(i ,j,k,nrhs,itrc)- &
280 & tl_t(i-1,j,k,nrhs,itrc)))- &
281# ifdef TL_IOMS
282 & cff* &
283 & (hz(i,j,k)+hz(i-1,j,k))* &
284 & (t(i,j,k,nrhs,itrc)-t(i-1,j,k,nrhs,itrc))
285# endif
286 END IF
287#else
288!^ FX(i,j)=cff* &
289!^ & (Hz(i,j,k)+Hz(i-1,j,k))* &
290!^ & (t(i,j,k,nrhs,itrc)-t(i-1,j,k,nrhs,itrc))
291!^
292 tl_fx(i,j)=cff* &
293 & ((tl_hz(i,j,k)+tl_hz(i-1,j,k))* &
294 & (t(i ,j,k,nrhs,itrc)- &
295 & t(i-1,j,k,nrhs,itrc))+ &
296 & (hz(i,j,k)+hz(i-1,j,k))* &
297 & (tl_t(i ,j,k,nrhs,itrc)- &
298 & tl_t(i-1,j,k,nrhs,itrc)))- &
299# ifdef TL_IOMS
300 & cff* &
301 & (hz(i,j,k)+hz(i-1,j,k))* &
302 & (t(i,j,k,nrhs,itrc)-t(i-1,j,k,nrhs,itrc))
303# endif
304#endif
305#ifdef MASKING
306!^ FX(i,j)=FX(i,j)*umask(i,j)
307!^
308 tl_fx(i,j)=tl_fx(i,j)*umask(i,j)
309#endif
310#ifdef WET_DRY_NOT_YET
311 fx(i,j)=fx(i,j)*umask_wet(i,j)
312#endif
313 END DO
314 END DO
315 DO j=jstr,jend+1
316 DO i=istr,iend
317#ifdef DIFF_3DCOEF
318 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i,j-1,k))* &
319 & pnom_v(i,j)
320#else
321 cff=0.25_r8*(diff2(i,j,itrc)+diff2(i,j-1,itrc))* &
322 & pnom_v(i,j)
323#endif
324#if defined TS_MIX_STABILITY
325!^ FE(i,j)=cff* &
326!^ & (Hz(i,j,k)+Hz(i,j-1,k))* &
327!^ & (0.75_r8*(t(i,j ,k,nrhs,itrc)- &
328!^ & t(i,j-1,k,nrhs,itrc))+ &
329!^ & 0.25_r8*(t(i,j ,k,nstp,itrc)- &
330!^ & t(i,j-1,k,nstp,itrc)))
331!^
332 tl_fe(i,j)=cff* &
333 & ((tl_hz(i,j,k)+tl_hz(i,j-1,k))* &
334 & (0.75_r8*(t(i,j ,k,nrhs,itrc)- &
335 & t(i,j-1,k,nrhs,itrc))+ &
336 & 0.25_r8*(t(i,j ,k,nstp,itrc)- &
337 & t(i,j-1,k,nstp,itrc)))+ &
338 & (hz(i,j,k)+hz(i,j-1,k))* &
339 & (0.75_r8*(tl_t(i,j ,k,nrhs,itrc)- &
340 & tl_t(i,j-1,k,nrhs,itrc))+ &
341 & 0.25_r8*(tl_t(i,j ,k,nstp,itrc)- &
342 & tl_t(i,j-1,k,nstp,itrc))))- &
343# ifdef TL_IOMS
344 & cff* &
345 & (hz(i,j,k)+hz(i,j-1,k))* &
346 & (0.75_r8*(t(i,j ,k,nrhs,itrc)- &
347 & t(i,j-1,k,nrhs,itrc))+ &
348 & 0.25_r8*(t(i,j ,k,nstp,itrc)- &
349 & t(i,j-1,k,nstp,itrc)))
350# endif
351#elif defined TS_MIX_CLIMA
352 IF (ltracerclm(itrc,ng)) THEN
353!^ FE(i,j)=cff* &
354!^ & (Hz(i,j,k)+Hz(i,j-1,k))* &
355!^ & ((t(i,j ,k,nrhs,itrc)-tclm(i,j ,k,itrc))- &
356!^ & (t(i,j-1,k,nrhs,itrc)-tclm(i,j-1,k,itrc)))
357!^
358 tl_fe(i,j)=cff* &
359 & ((tl_hz(i,j,k)+tl_hz(i,j-1,k))* &
360 & ((t(i,j ,k,nrhs,itrc)- &
361 & tclm(i,j ,k,itrc))- &
362 & (t(i,j-1,k,nrhs,itrc)- &
363 & tclm(i,j-1,k,itrc)))+ &
364 & (hz(i,j,k)+hz(i,j-1,k))* &
365 & (tl_t(i,j ,k,nrhs,itrc)- &
366 & tl_t(i,j-1,k,nrhs,itrc)))- &
367# ifdef TL_IOMS
368 & cff* &
369 & (hz(i,j,k)+hz(i,j-1,k))* &
370 & ((t(i,j ,k,nrhs,itrc)-tclm(i,j ,k,itrc))- &
371 & (t(i,j-1,k,nrhs,itrc)-tclm(i,j-1,k,itrc)))
372# endif
373 ELSE
374!^ FE(i,j)=cff* &
375!^ & (Hz(i,j,k)+Hz(i,j-1,k))* &
376!^ & (t(i,j,k,nrhs,itrc)-t(i,j-1,k,nrhs,itrc))
377!^
378 tl_fe(i,j)=cff* &
379 & ((tl_hz(i,j,k)+tl_hz(i,j-1,k))* &
380 & (t(i,j ,k,nrhs,itrc)- &
381 & t(i,j-1,k,nrhs,itrc))+ &
382 & (hz(i,j,k)+hz(i,j-1,k))* &
383 & (tl_t(i,j ,k,nrhs,itrc)- &
384 & tl_t(i,j-1,k,nrhs,itrc)))- &
385# ifdef TL_IOMS
386 & cff* &
387 & (hz(i,j,k)+hz(i,j-1,k))* &
388 & (t(i,j,k,nrhs,itrc)-t(i,j-1,k,nrhs,itrc))
389# endif
390 END IF
391#else
392!^ FE(i,j)=cff* &
393!^ & (Hz(i,j,k)+Hz(i,j-1,k))* &
394!^ & (t(i,j,k,nrhs,itrc)-t(i,j-1,k,nrhs,itrc))
395!^
396 tl_fe(i,j)=cff* &
397 & ((tl_hz(i,j,k)+tl_hz(i,j-1,k))* &
398 & (t(i,j ,k,nrhs,itrc)- &
399 & t(i,j-1,k,nrhs,itrc))+ &
400 & (hz(i,j,k)+hz(i,j-1,k))* &
401 & (tl_t(i,j ,k,nrhs,itrc)- &
402 & tl_t(i,j-1,k,nrhs,itrc)))- &
403# ifdef TL_IOMS
404 & cff* &
405 & (hz(i,j,k)+hz(i,j-1,k))* &
406 & (t(i,j,k,nrhs,itrc)-t(i,j-1,k,nrhs,itrc))
407# endif
408#endif
409#ifdef MASKING
410!^ FE(i,j)=FE(i,j)*vmask(i,j)
411!^
412 tl_fe(i,j)=tl_fe(i,j)*vmask(i,j)
413#endif
414#ifdef WET_DRY_NOT_YET
415 fe(i,j)=fe(i,j)*vmask_wet(i,j)
416#endif
417 END DO
418 END DO
419!
420! Time-step harmonic, S-surfaces diffusion term (m Tunits).
421!
422 DO j=jstr,jend
423 DO i=istr,iend
424!^ cff=dt(ng)*pm(i,j)*pn(i,j)* &
425!^ & (FX(i+1,j)-FX(i,j)+ &
426!^ & FE(i,j+1)-FE(i,j))
427!^
428 tl_cff=dt(ng)*pm(i,j)*pn(i,j)* &
429 & (tl_fx(i+1,j)-tl_fx(i,j)+ &
430 & tl_fe(i,j+1)-tl_fe(i,j))
431!^ t(i,j,k,nnew,itrc)=t(i,j,k,nnew,itrc)+cff
432!^
433 tl_t(i,j,k,nnew,itrc)=tl_t(i,j,k,nnew,itrc)+tl_cff
434#ifdef DIAGNOSTICS_TS
435!! DiaTwrk(i,j,k,itrc,iThdif)=cff
436#endif
437 END DO
438 END DO
439 END DO
440 END DO
441!
442 RETURN
443 END SUBROUTINE rp_t3dmix2_s_tile
444
445 END MODULE rp_t3dmix2_mod
type(t_clima), dimension(:), allocatable clima
Definition mod_clima.F:153
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 irpm
Definition mod_param.F:664
real(dp), dimension(:), allocatable dt
logical, dimension(:,:), allocatable ltracerclm
integer, dimension(:), allocatable nrhs
integer, dimension(:), allocatable nnew
integer, dimension(:), allocatable nstp
subroutine rp_t3dmix2_s_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, nstp, nnew, umask, vmask, umask_wet, vmask_wet, hz, tl_hz, pmon_u, pnom_v, pm, pn, diff3d_r, diff2, tclm, t, tl_t)
subroutine, public rp_t3dmix2(ng, tile)
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