ROMS
Loading...
Searching...
No Matches
tl_t3dmix4_s.h
Go to the documentation of this file.
1 MODULE tl_t3dmix4_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 tangent linear horizontal biharmonic !
11! mixing of tracers along S-coordinate levels surfaces. !
12! !
13!=======================================================================
14!
15 implicit none
16!
17 PRIVATE
18 PUBLIC tl_t3dmix4
19!
20 CONTAINS
21!
22!***********************************************************************
23 SUBROUTINE tl_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, itlm, 27, __line__, myfile)
51#endif
52 CALL tl_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_NOT_YET
61 & grid(ng) % umask_wet, &
62 & grid(ng) % vmask_wet, &
63#endif
64 & grid(ng) % Hz, &
65 & grid(ng) % tl_Hz, &
66 & grid(ng) % pmon_u, &
67 & grid(ng) % pnom_v, &
68 & grid(ng) % pm, &
69 & grid(ng) % pn, &
70#ifdef DIFF_3DCOEF
71# ifdef TS_U3ADV_SPLIT_NOT_YET
72 & mixing(ng) % diff3d_u, &
73 & mixing(ng) % diff3d_v, &
74# else
75 & mixing(ng) % diff3d_r, &
76# endif
77#else
78 & mixing(ng) % diff4, &
79#endif
80#ifdef TS_MIX_CLIMA
81 & clima(ng) % tclm, &
82#endif
83#ifdef DIAGNOSTICS_TS
84!! & DIAGS(ng) % DiaTwrk, &
85#endif
86 & ocean(ng) % t, &
87 & ocean(ng) % tl_t)
88#ifdef PROFILE
89 CALL wclock_off (ng, itlm, 27, __line__, myfile)
90#endif
91!
92 RETURN
93 END SUBROUTINE tl_t3dmix4
94!
95!***********************************************************************
96 SUBROUTINE tl_t3dmix4_s_tile (ng, tile, &
97 & LBi, UBi, LBj, UBj, &
98 & IminS, ImaxS, JminS, JmaxS, &
99 & nrhs, nstp, nnew, &
100#ifdef MASKING
101 & umask, vmask, &
102#endif
103#ifdef WET_DRY_NOT_YET
104 & umask_wet, vmask_wet, &
105#endif
106 & Hz, tl_Hz, &
107 & pmon_u, pnom_v, pm, pn, &
108#ifdef DIFF_3DCOEF
109# ifdef TS_U3ADV_SPLIT_NOT_YET
110 & diff3d_u, diff3d_v, &
111# else
112 & diff3d_r, &
113# endif
114#else
115 & diff4, &
116#endif
117#ifdef TS_MIX_CLIMA
118 & tclm, &
119#endif
120#ifdef DIAGNOSTICS_TS
121!! & DiaTwrk, &
122#endif
123 & t, tl_t)
124!***********************************************************************
125!
126 USE mod_param
127 USE mod_ncparam
128 USE mod_scalars
129!
130! Imported variable declarations.
131!
132 integer, intent(in) :: ng, tile
133 integer, intent(in) :: LBi, UBi, LBj, UBj
134 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
135 integer, intent(in) :: nrhs, nstp, nnew
136
137#ifdef ASSUMED_SHAPE
138# ifdef MASKING
139 real(r8), intent(in) :: umask(LBi:,LBj:)
140 real(r8), intent(in) :: vmask(LBi:,LBj:)
141# endif
142# ifdef WET_DRY_NOT_YET
143 real(r8), intent(in) :: umask_wet(LBi:,LBj:)
144 real(r8), intent(in) :: vmask_wet(LBi:,LBj:)
145# endif
146# ifdef DIFF_3DCOEF
147# ifdef TS_U3ADV_SPLIT_NOT_YET
148 real(r8), intent(in) :: diff3d_u(LBi:,LBj:,:)
149 real(r8), intent(in) :: diff3d_v(LBi:,LBj:,:)
150# else
151 real(r8), intent(in) :: diff3d_r(LBi:,LBj:,:)
152# endif
153# else
154 real(r8), intent(in) :: diff4(LBi:,LBj:,:)
155# endif
156 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
157 real(r8), intent(in) :: tl_Hz(LBi:,LBj:,:)
158 real(r8), intent(in) :: pmon_u(LBi:,LBj:)
159 real(r8), intent(in) :: pnom_v(LBi:,LBj:)
160 real(r8), intent(in) :: pm(LBi:,LBj:)
161 real(r8), intent(in) :: pn(LBi:,LBj:)
162 real(r8), intent(in) :: t(LBi:,LBj:,:,:,:)
163# ifdef TS_MIX_CLIMA
164 real(r8), intent(in) :: tclm(LBi:,LBj:,:,:)
165# endif
166# ifdef DIAGNOSTICS_TS
167!! real(r8), intent(inout) :: DiaTwrk(LBi:,LBj:,:,:,:)
168# endif
169 real(r8), intent(inout) :: tl_t(LBi:,LBj:,:,:,:)
170#else
171# ifdef MASKING
172 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
173 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
174# endif
175# ifdef WET_DRY_NOT_YET
176 real(r8), intent(in) :: umask_wet(LBi:UBi,LBj:UBj)
177 real(r8), intent(in) :: vmask_wet(LBi:UBi,LBj:UBj)
178# endif
179# ifdef DIFF_3DCOEF
180# ifdef TS_U3ADV_SPLIT_NOT_YET
181 real(r8), intent(in) :: diff3d_u(LBi:UBi,LBj:UBj,N(ng))
182 real(r8), intent(in) :: diff3d_v(LBi:UBi,LBj:UBj,N(ng))
183# else
184 real(r8), intent(in) :: diff3d_r(LBi:UBi,LBj:UBj,N(ng))
185# endif
186# else
187 real(r8), intent(in) :: diff4(LBi:UBi,LBj:UBj,NT(ng))
188# endif
189 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
190 real(r8), intent(in) :: tl_Hz(LBi:UBi,LBj:UBj,N(ng))
191 real(r8), intent(in) :: pmon_u(LBi:UBi,LBj:UBj)
192 real(r8), intent(in) :: pnom_v(LBi:UBi,LBj:UBj)
193 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
194 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
195 real(r8), intent(in) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
196# ifdef TS_MIX_CLIMA
197 real(r8), intent(in) :: tclm(LBi:UBi,LBj:UBj,N(ng),NT(ng))
198# endif
199# ifdef DIAGNOSTICS_TS
200!! real(r8), intent(inout) :: DiaTwrk(LBi:UBi,LBj:UBj,N(ng),NT(ng), &
201!! & NDT)
202# endif
203 real(r8), intent(inout) :: tl_t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
204#endif
205!
206! Local variable declarations.
207!
208 integer :: Imin, Imax, Jmin, Jmax
209 integer :: i, itrc, j, k
210
211 real(r8) :: cff, cff1, tl_cff, tl_cff1
212
213 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FE
214 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FX
215 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: LapT
216
217 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_FE
218 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_FX
219 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_LapT
220
221#include "set_bounds.h"
222!
223!-----------------------------------------------------------------------
224! Compute horizontal biharmonic diffusion along constant S-surfaces.
225! The biharmonic operator is computed by applying the harmonic
226! operator twice.
227#ifdef TS_MIX_STABILITY
228! In order to increase stability, the biharmonic operator is applied
229! as: 3/4 t(:,:,:,nrhs,:) + 1/4 t(:,:,:,nstp,:).
230#endif
231!-----------------------------------------------------------------------
232!
233! Set local I- and J-ranges.
234!
235 IF (ewperiodic(ng)) THEN
236 imin=istr-1
237 imax=iend+1
238 ELSE
239 imin=max(istr-1,1)
240 imax=min(iend+1,lm(ng))
241 END IF
242 IF (nsperiodic(ng)) THEN
243 jmin=jstr-1
244 jmax=jend+1
245 ELSE
246 jmin=max(jstr-1,1)
247 jmax=min(jend+1,mm(ng))
248 END IF
249!
250! Compute horizontal tracer flux in the XI- and ETA-directions.
251!
252 DO itrc=1,nt(ng)
253 DO k=1,n(ng)
254 DO j=jmin,jmax
255 DO i=imin,imax+1
256#ifdef DIFF_3DCOEF
257# ifdef TS_U3ADV_SPLIT_NOT_YET
258 cff=0.5_r8*diff3d_u(i,j,k)*pmon_u(i,j)
259# else
260 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i-1,j,k))* &
261 & pmon_u(i,j)
262# endif
263#else
264 cff=0.25_r8*(diff4(i,j,itrc)+diff4(i-1,j,itrc))* &
265 & pmon_u(i,j)
266#endif
267#ifdef MASKING
268 cff=cff*umask(i,j)
269#endif
270#ifdef WET_DRY_NOT_YET
271 cff=cff*umask_wet(i,j)
272#endif
273#if defined TS_MIX_STABILITY
274 fx(i,j)=cff* &
275 & (hz(i,j,k)+hz(i-1,j,k))* &
276 & (0.75_r8*(t(i ,j,k,nrhs,itrc)- &
277 & t(i-1,j,k,nrhs,itrc))+ &
278 & 0.25_r8*(t(i ,j,k,nstp,itrc)- &
279 & t(i-1,j,k,nstp,itrc)))
280 tl_fx(i,j)=cff* &
281 & ((tl_hz(i,j,k)+tl_hz(i-1,j,k))* &
282 & (0.75_r8*(t(i ,j,k,nrhs,itrc)- &
283 & t(i-1,j,k,nrhs,itrc))+ &
284 & 0.25_r8*(t(i ,j,k,nstp,itrc)- &
285 & t(i-1,j,k,nstp,itrc)))+ &
286 & (hz(i,j,k)+hz(i-1,j,k))* &
287 & (0.75_r8*(tl_t(i ,j,k,nrhs,itrc)- &
288 & tl_t(i-1,j,k,nrhs,itrc))+ &
289 & 0.25_r8*(tl_t(i ,j,k,nstp,itrc)- &
290 & tl_t(i-1,j,k,nstp,itrc))))
291#elif defined TS_MIX_CLIMA
292 IF (ltracerclm(itrc,ng)) THEN
293 fx(i,j)=cff* &
294 & (hz(i,j,k)+hz(i-1,j,k))* &
295 & ((t(i ,j,k,nrhs,itrc)-tclm(i ,j,k,itrc))- &
296 & (t(i-1,j,k,nrhs,itrc)-tclm(i-1,j,k,itrc)))
297 tl_fx(i,j)=cff* &
298 & ((tl_hz(i,j,k)+tl_hz(i-1,j,k))* &
299 & ((t(i ,j,k,nrhs,itrc)- &
300 & tclm(i ,j,k,itrc))- &
301 & (t(i-1,j,k,nrhs,itrc)- &
302 & tclm(i-1,j,k,itrc)))+ &
303 & (hz(i,j,k)+hz(i-1,j,k))* &
304 & (tl_t(i ,j,k,nrhs,itrc)- &
305 & tl_t(i-1,j,k,nrhs,itrc)))
306 ELSE
307 fx(i,j)=cff* &
308 & (hz(i,j,k)+hz(i-1,j,k))* &
309 & (t(i ,j,k,nrhs,itrc)- &
310 & t(i-1,j,k,nrhs,itrc))
311 tl_fx(i,j)=cff* &
312 & ((tl_hz(i,j,k)+tl_hz(i-1,j,k))* &
313 & (t(i ,j,k,nrhs,itrc)- &
314 & t(i-1,j,k,nrhs,itrc))+ &
315 & (hz(i,j,k)+hz(i-1,j,k))* &
316 & (tl_t(i ,j,k,nrhs,itrc)- &
317 & tl_t(i-1,j,k,nrhs,itrc)))
318 END IF
319#else
320 fx(i,j)=cff* &
321 & (hz(i,j,k)+hz(i-1,j,k))* &
322 & (t(i ,j,k,nrhs,itrc)- &
323 & t(i-1,j,k,nrhs,itrc))
324 tl_fx(i,j)=cff* &
325 & ((tl_hz(i,j,k)+tl_hz(i-1,j,k))* &
326 & (t(i ,j,k,nrhs,itrc)- &
327 & t(i-1,j,k,nrhs,itrc))+ &
328 & (hz(i,j,k)+hz(i-1,j,k))* &
329 & (tl_t(i ,j,k,nrhs,itrc)- &
330 & tl_t(i-1,j,k,nrhs,itrc)))
331#endif
332 END DO
333 END DO
334 DO j=jmin,jmax+1
335 DO i=imin,imax
336#ifdef DIFF_3DCOEF
337# ifdef TS_U3ADV_SPLIT_NOT_YET
338 cff=0.5_r8*diff3d_v(i,j,k)*pnom_v(i,j)
339# else
340 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i,j-1,k))* &
341 & pnom_v(i,j)
342# endif
343#else
344 cff=0.25_r8*(diff4(i,j,itrc)+diff4(i,j-1,itrc))* &
345 & pnom_v(i,j)
346#endif
347#ifdef MASKING
348 cff=cff*vmask(i,j)
349#endif
350#ifdef WET_DRY_NOT_YET
351 cff=cff*vmask_wet(i,j)
352#endif
353#if defined TS_MIX_STABILITY
354 fe(i,j)=cff* &
355 & (hz(i,j,k)+hz(i,j-1,k))* &
356 & (0.75_r8*(t(i,j ,k,nrhs,itrc)- &
357 & t(i,j-1,k,nrhs,itrc))+ &
358 & 0.25_r8*(t(i,j ,k,nstp,itrc)- &
359 & t(i,j-1,k,nstp,itrc)))
360 tl_fe(i,j)=cff* &
361 & ((tl_hz(i,j,k)+tl_hz(i,j-1,k))* &
362 & (0.75_r8*(t(i,j ,k,nrhs,itrc)- &
363 & t(i,j-1,k,nrhs,itrc))+ &
364 & 0.25_r8*(t(i,j ,k,nstp,itrc)- &
365 & t(i,j-1,k,nstp,itrc)))+ &
366 & (hz(i,j,k)+hz(i,j-1,k))* &
367 & (0.75_r8*(tl_t(i,j ,k,nrhs,itrc)- &
368 & tl_t(i,j-1,k,nrhs,itrc))+ &
369 & 0.25_r8*(tl_t(i,j ,k,nstp,itrc)- &
370 & tl_t(i,j-1,k,nstp,itrc))))
371#elif defined TS_MIX_CLIMA
372 IF (ltracerclm(itrc,ng)) THEN
373 fe(i,j)=cff* &
374 & (hz(i,j,k)+hz(i,j-1,k))* &
375 & ((t(i,j ,k,nrhs,itrc)-tclm(i,j ,k,itrc))- &
376 & (t(i,j-1,k,nrhs,itrc)-tclm(i,j-1,k,itrc)))
377 tl_fe(i,j)=cff* &
378 & ((tl_hz(i,j,k)+tl_hz(i,j-1,k))* &
379 & ((t(i,j ,k,nrhs,itrc)- &
380 & tclm(i,j ,k,itrc))- &
381 & (t(i,j-1,k,nrhs,itrc)- &
382 & tclm(i,j-1,k,itrc)))+ &
383 & (hz(i,j,k)+hz(i,j-1,k))* &
384 & (tl_t(i,j ,k,nrhs,itrc)- &
385 & tl_t(i,j-1,k,nrhs,itrc)))
386 ELSE
387 fe(i,j)=cff* &
388 & (hz(i,j,k)+hz(i,j-1,k))* &
389 & (t(i,j ,k,nrhs,itrc)- &
390 & t(i,j-1,k,nrhs,itrc))
391 tl_fe(i,j)=cff* &
392 & ((tl_hz(i,j,k)+tl_hz(i,j-1,k))* &
393 & (t(i,j ,k,nrhs,itrc)- &
394 & t(i,j-1,k,nrhs,itrc))+ &
395 & (hz(i,j,k)+hz(i,j-1,k))* &
396 & (tl_t(i,j ,k,nrhs,itrc)- &
397 & tl_t(i,j-1,k,nrhs,itrc)))
398 END IF
399#else
400 fe(i,j)=cff* &
401 & (hz(i,j,k)+hz(i,j-1,k))* &
402 & (t(i,j ,k,nrhs,itrc)- &
403 & t(i,j-1,k,nrhs,itrc))
404 tl_fe(i,j)=cff* &
405 & ((tl_hz(i,j,k)+tl_hz(i,j-1,k))* &
406 & (t(i,j ,k,nrhs,itrc)- &
407 & t(i,j-1,k,nrhs,itrc))+ &
408 & (hz(i,j,k)+hz(i,j-1,k))* &
409 & (tl_t(i,j ,k,nrhs,itrc)- &
410 & tl_t(i,j-1,k,nrhs,itrc)))
411#endif
412 END DO
413 END DO
414!
415! Compute first harmonic operator and multiply by the metrics of the
416! second harmonic operator.
417!
418 DO j=jmin,jmax
419 DO i=imin,imax
420 cff=1.0_r8/hz(i,j,k)
421 tl_cff=-cff*cff*tl_hz(i,j,k)
422 lapt(i,j)=pm(i,j)*pn(i,j)*cff* &
423 & (fx(i+1,j)-fx(i,j)+ &
424 & fe(i,j+1)-fe(i,j))
425 tl_lapt(i,j)=pm(i,j)*pn(i,j)* &
426 & (tl_cff* &
427 & (fx(i+1,j)-fx(i,j)+ &
428 & fe(i,j+1)-fe(i,j))+ &
429 & cff* &
430 & (tl_fx(i+1,j)-tl_fx(i,j)+ &
431 & tl_fe(i,j+1)-tl_fe(i,j)))
432 END DO
433 END DO
434!
435! Apply boundary conditions (except periodic; closed or gradient)
436! to the first harmonic operator.
437!
438 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
439 IF (domain(ng)%Western_Edge(tile)) THEN
440 IF (tl_lbc(iwest,istvar(itrc),ng)%closed) THEN
441 DO j=jmin,jmax
442 lapt(istr-1,j)=0.0_r8
443 tl_lapt(istr-1,j)=0.0_r8
444 END DO
445 ELSE
446 DO j=jmin,jmax
447 lapt(istr-1,j)=lapt(istr,j)
448 tl_lapt(istr-1,j)=tl_lapt(istr,j)
449 END DO
450 END IF
451 END IF
452 END IF
453!
454 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
455 IF (domain(ng)%Eastern_Edge(tile)) THEN
456 IF (tl_lbc(ieast,istvar(itrc),ng)%closed) THEN
457 DO j=jmin,jmax
458 lapt(iend+1,j)=0.0_r8
459 tl_lapt(iend+1,j)=0.0_r8
460 END DO
461 ELSE
462 DO j=jmin,jmax
463 lapt(iend+1,j)=lapt(iend,j)
464 tl_lapt(iend+1,j)=tl_lapt(iend,j)
465 END DO
466 END IF
467 END IF
468 END IF
469!
470 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
471 IF (domain(ng)%Southern_Edge(tile)) THEN
472 IF (tl_lbc(isouth,istvar(itrc),ng)%closed) THEN
473 DO i=imin,imax
474 lapt(i,jstr-1)=0.0_r8
475 tl_lapt(i,jstr-1)=0.0_r8
476 END DO
477 ELSE
478 DO i=imin,imax
479 lapt(i,jstr-1)=lapt(i,jstr)
480 tl_lapt(i,jstr-1)=tl_lapt(i,jstr)
481 END DO
482 END IF
483 END IF
484 END IF
485!
486 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
487 IF (domain(ng)%Northern_Edge(tile)) THEN
488 IF (tl_lbc(inorth,istvar(itrc),ng)%closed) THEN
489 DO i=imin,imax
490 lapt(i,jend+1)=0.0_r8
491 tl_lapt(i,jend+1)=0.0_r8
492 END DO
493 ELSE
494 DO i=imin,imax
495 lapt(i,jend+1)=lapt(i,jend)
496 tl_lapt(i,jend+1)=tl_lapt(i,jend)
497 END DO
498 END IF
499 END IF
500 END IF
501!
502! Compute FX=d(LapT)/d(xi) and FE=d(LapT)/d(eta) terms.
503!
504 DO j=jstr,jend
505 DO i=istr,iend+1
506#ifdef DIFF_3DCOEF
507# ifdef TS_U3ADV_SPLIT_NOT_YET
508 cff=0.5_r8*diff3d_u(i,j,k)*pmon_u(i,j)
509# else
510 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i-1,j,k))* &
511 & pmon_u(i,j)
512# endif
513#else
514 cff=0.25_r8*(diff4(i,j,itrc)+diff4(i-1,j,itrc))* &
515 & pmon_u(i,j)
516#endif
517!^ FX(i,j)=cff* &
518!^ & (Hz(i,j,k)+Hz(i-1,j,k))* &
519!^ & (LapT(i,j)-LapT(i-1,j))
520!^
521 tl_fx(i,j)=cff* &
522 & ((tl_hz(i,j,k)+tl_hz(i-1,j,k))* &
523 & (lapt(i,j)-lapt(i-1,j))+ &
524 & (hz(i,j,k)+hz(i-1,j,k))* &
525 & (tl_lapt(i,j)-tl_lapt(i-1,j)))
526#ifdef MASKING
527!^ FX(i,j)=FX(i,j)*umask(i,j)
528!^
529 tl_fx(i,j)=tl_fx(i,j)*umask(i,j)
530#endif
531#ifdef WET_DRY_NOT_YET
532 fx(i,j)=fx(i,j)*umask_wet(i,j)
533#endif
534 END DO
535 END DO
536 DO j=jstr,jend+1
537 DO i=istr,iend
538#ifdef DIFF_3DCOEF
539# ifdef TS_U3ADV_SPLIT_NOT_YET
540 cff=0.5_r8*diff3d_v(i,j,k)*pnom_v(i,j)
541# else
542 cff=0.25_r8*(diff3d_r(i,j,k)+diff3d_r(i,j-1,k))* &
543 & pnom_v(i,j)
544# endif
545#else
546 cff=0.25_r8*(diff4(i,j,itrc)+diff4(i,j-1,itrc))* &
547 & pnom_v(i,j)
548#endif
549!^ FE(i,j)=cff* &
550!^ & (Hz(i,j,k)+Hz(i,j-1,k))* &
551!^ & (LapT(i,j)-LapT(i,j-1))
552!^
553 tl_fe(i,j)=cff* &
554 & ((tl_hz(i,j,k)+tl_hz(i,j-1,k))* &
555 & (lapt(i,j)-lapt(i,j-1))+ &
556 & (hz(i,j,k)+hz(i,j-1,k))* &
557 & (tl_lapt(i,j)-tl_lapt(i,j-1)))
558#ifdef MASKING
559!^ FE(i,j)=FE(i,j)*vmask(i,j)
560!^
561 tl_fe(i,j)=tl_fe(i,j)*vmask(i,j)
562#endif
563#ifdef WET_DRY_NOT_YET
564 fe(i,j)=fe(i,j)*vmask_wet(i,j)
565#endif
566 END DO
567 END DO
568!
569! Time-step biharmonic, S-surfaces diffusion term (m Tunits).
570!
571 DO j=jstr,jend
572 DO i=istr,iend
573!^ cff=dt(ng)*pm(i,j)*pn(i,j)* &
574!^ & (FX(i+1,j)-FX(i,j)+ &
575!^ & FE(i,j+1)-FE(i,j))
576!^
577 tl_cff=dt(ng)*pm(i,j)*pn(i,j)* &
578 & (tl_fx(i+1,j)-tl_fx(i,j)+ &
579 & tl_fe(i,j+1)-tl_fe(i,j))
580!^ t(i,j,k,nnew,itrc)=t(i,j,k,nnew,itrc)-cff
581!^
582 tl_t(i,j,k,nnew,itrc)=tl_t(i,j,k,nnew,itrc)-tl_cff
583#ifdef DIAGNOSTICS_TS
584!! DiaTwrk(i,j,k,itrc,iThdif)=-cff
585#endif
586 END DO
587 END DO
588 END DO
589 END DO
590!
591 RETURN
592 END SUBROUTINE tl_t3dmix4_s_tile
593
594 END MODULE tl_t3dmix4_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
integer, dimension(:), allocatable istvar
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
type(t_lbc), dimension(:,:,:), allocatable tl_lbc
Definition mod_param.F:379
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer, dimension(:), allocatable lm
Definition mod_param.F:455
integer, parameter itlm
Definition mod_param.F:663
integer, dimension(:), allocatable mm
Definition mod_param.F:456
real(dp), dimension(:), allocatable dt
logical, dimension(:), allocatable ewperiodic
integer, parameter iwest
logical, dimension(:), allocatable nsperiodic
logical, dimension(:,:), allocatable compositegrid
integer, parameter isouth
integer, parameter ieast
logical, dimension(:,:), allocatable ltracerclm
integer, parameter inorth
integer, dimension(:), allocatable nrhs
integer, dimension(:), allocatable nnew
integer, dimension(:), allocatable nstp
subroutine tl_t3dmix4_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_u, diff3d_v, diff3d_r, diff4, tclm, t, tl_t)
subroutine, public tl_t3dmix4(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