ROMS
Loading...
Searching...
No Matches
pre_step3d_mod Module Reference

Functions/Subroutines

subroutine, public pre_step3d (ng, tile)
 
subroutine pre_step3d_tile (ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, nstp, nnew, rmask, umask, vmask, rmask_wet, pm, pn, hz, huon, hvom, z_r, z_w, btflx, bustr, bvstr, stflx, sustr, svstr, srflx, akt, akv, ghats, w, ru, rv, diatwrk, diau3wrk, diav3wrk, diaru, diarv, t, u, v)
 

Function/Subroutine Documentation

◆ pre_step3d()

subroutine, public pre_step3d_mod::pre_step3d ( integer, intent(in) ng,
integer, intent(in) tile )

Definition at line 39 of file pre_step3d.F.

40!***********************************************************************
41!
42 USE mod_param
43# ifdef DIAGNOSTICS
44 USE mod_diags
45# endif
46 USE mod_forces
47 USE mod_grid
48 USE mod_mixing
49 USE mod_ocean
50 USE mod_stepping
51!
52! Imported variable declarations.
53!
54 integer, intent(in) :: ng, tile
55!
56! Local variable declarations.
57!
58 character (len=*), parameter :: MyFile = &
59 & __FILE__
60!
61# include "tile.h"
62!
63# ifdef PROFILE
64 CALL wclock_on (ng, inlm, 22, __line__, myfile)
65# endif
66 CALL pre_step3d_tile (ng, tile, &
67 & lbi, ubi, lbj, ubj, &
68 & imins, imaxs, jmins, jmaxs, &
69 & nrhs(ng), nstp(ng), nnew(ng), &
70# ifdef MASKING
71 & grid(ng) % rmask, &
72 & grid(ng) % umask, &
73 & grid(ng) % vmask, &
74# if defined SOLAR_SOURCE && defined WET_DRY
75 & grid(ng) % rmask_wet, &
76# endif
77# endif
78
79 & grid(ng) % pm, &
80 & grid(ng) % pn, &
81 & grid(ng) % Hz, &
82 & grid(ng) % Huon, &
83 & grid(ng) % Hvom, &
84 & grid(ng) % z_r, &
85 & grid(ng) % z_w, &
86 & forces(ng) % btflx, &
87 & forces(ng) % bustr, &
88 & forces(ng) % bvstr, &
89 & forces(ng) % stflx, &
90 & forces(ng) % sustr, &
91 & forces(ng) % svstr, &
92# ifdef SOLAR_SOURCE
93 & forces(ng) % srflx, &
94# endif
95 & mixing(ng) % Akt, &
96 & mixing(ng) % Akv, &
97# ifdef LMD_NONLOCAL
98 & mixing(ng) % ghats, &
99# endif
100 & ocean(ng) % W, &
101# ifdef WEC_VF
102 & ocean(ng) % W_stokes, &
103# endif
104 & ocean(ng) % ru, &
105 & ocean(ng) % rv, &
106# ifdef DIAGNOSTICS_TS
107 & diags(ng) % DiaTwrk, &
108# endif
109# ifdef DIAGNOSTICS_UV
110 & diags(ng) % DiaU3wrk, &
111 & diags(ng) % DiaV3wrk, &
112 & diags(ng) % DiaRU, &
113 & diags(ng) % DiaRV, &
114# endif
115 & ocean(ng) % t, &
116 & ocean(ng) % u, &
117 & ocean(ng) % v)
118# ifdef PROFILE
119 CALL wclock_off (ng, inlm, 22, __line__, myfile)
120# endif
121!
122 RETURN
type(t_diags), dimension(:), allocatable diags
Definition mod_diags.F:100
type(t_forces), dimension(:), allocatable forces
Definition mod_forces.F:554
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
integer, dimension(:), allocatable nrhs
integer, dimension(:), allocatable nnew
integer, dimension(:), allocatable nstp
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

References mod_diags::diags, mod_forces::forces, mod_grid::grid, mod_param::inlm, mod_mixing::mixing, mod_stepping::nnew, mod_stepping::nrhs, mod_stepping::nstp, mod_ocean::ocean, pre_step3d_tile(), wclock_off(), and wclock_on().

Referenced by rhs3d_mod::rhs3d().

Here is the call graph for this function:
Here is the caller graph for this function:

◆ pre_step3d_tile()

subroutine pre_step3d_mod::pre_step3d_tile ( integer, intent(in) ng,
integer, intent(in) tile,
integer, intent(in) lbi,
integer, intent(in) ubi,
integer, intent(in) lbj,
integer, intent(in) ubj,
integer, intent(in) imins,
integer, intent(in) imaxs,
integer, intent(in) jmins,
integer, intent(in) jmaxs,
integer, intent(in) nrhs,
integer, intent(in) nstp,
integer, intent(in) nnew,
real(r8), dimension(lbi:,lbj:), intent(in) rmask,
real(r8), dimension(lbi:,lbj:), intent(in) umask,
real(r8), dimension(lbi:,lbj:), intent(in) vmask,
real(r8), dimension(lbi:,lbj:), intent(in) rmask_wet,
real(r8), dimension(lbi:,lbj:), intent(in) pm,
real(r8), dimension(lbi:,lbj:), intent(in) pn,
real(r8), dimension(lbi:,lbj:,:), intent(in) hz,
real(r8), dimension(lbi:,lbj:,:), intent(in) huon,
real(r8), dimension(lbi:,lbj:,:), intent(in) hvom,
real(r8), dimension(lbi:,lbj:,:), intent(in) z_r,
real(r8), dimension(lbi:,lbj:,0:), intent(in) z_w,
real(r8), dimension(lbi:,lbj:,:), intent(in) btflx,
real(r8), dimension(lbi:,lbj:), intent(in) bustr,
real(r8), dimension(lbi:,lbj:), intent(in) bvstr,
real(r8), dimension(lbi:,lbj:,:), intent(in) stflx,
real(r8), dimension(lbi:,lbj:), intent(in) sustr,
real(r8), dimension(lbi:,lbj:), intent(in) svstr,
real(r8), dimension(lbi:,lbj:), intent(in) srflx,
real(r8), dimension(lbi:,lbj:,0:,:), intent(in) akt,
real(r8), dimension(lbi:,lbj:,0:), intent(in) akv,
real(r8), dimension(lbi:,lbj:,0:,:), intent(in) ghats,
real(r8), dimension(lbi:,lbj:,0:), intent(in) w,
real(r8), dimension(lbi:,lbj:,0:,:), intent(in) ru,
real(r8), dimension(lbi:,lbj:,0:,:), intent(in) rv,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) diatwrk,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) diau3wrk,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) diav3wrk,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) diaru,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) diarv,
real(r8), dimension(lbi:,lbj:,:,:,:), intent(inout) t,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) u,
real(r8), dimension(lbi:,lbj:,:,:), intent(inout) v )
private

Definition at line 126 of file pre_step3d.F.

161!***********************************************************************
162!
163 USE mod_param
164 USE mod_scalars
165 USE mod_sources
166!
168# ifdef DISTRIBUTE
170# endif
171 USE t3dbc_mod, ONLY : t3dbc_tile
172!
173! Imported variable declarations.
174!
175 integer, intent(in) :: ng, tile
176 integer, intent(in) :: LBi, UBi, LBj, UBj
177 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
178 integer, intent(in) :: nrhs, nstp, nnew
179!
180# ifdef ASSUMED_SHAPE
181# ifdef MASKING
182 real(r8), intent(in) :: rmask(LBi:,LBj:)
183 real(r8), intent(in) :: umask(LBi:,LBj:)
184 real(r8), intent(in) :: vmask(LBi:,LBj:)
185# if defined SOLAR_SOURCE && defined WET_DRY
186 real(r8), intent(in) :: rmask_wet(LBi:,LBj:)
187# endif
188# endif
189 real(r8), intent(in) :: pm(LBi:,LBj:)
190 real(r8), intent(in) :: pn(LBi:,LBj:)
191 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
192 real(r8), intent(in) :: Huon(LBi:,LBj:,:)
193 real(r8), intent(in) :: Hvom(LBi:,LBj:,:)
194 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
195 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
196 real(r8), intent(in) :: btflx(LBi:,LBj:,:)
197 real(r8), intent(in) :: bustr(LBi:,LBj:)
198 real(r8), intent(in) :: bvstr(LBi:,LBj:)
199 real(r8), intent(in) :: stflx(LBi:,LBj:,:)
200 real(r8), intent(in) :: sustr(LBi:,LBj:)
201 real(r8), intent(in) :: svstr(LBi:,LBj:)
202# ifdef SOLAR_SOURCE
203 real(r8), intent(in) :: srflx(LBi:,LBj:)
204# endif
205# ifdef SUN
206 real(r8), intent(in) :: Akt(LBi:UBi,LBj:UBj,0:N(ng),NAT)
207# else
208 real(r8), intent(in) :: Akt(LBi:,LBj:,0:,:)
209# endif
210 real(r8), intent(in) :: Akv(LBi:,LBj:,0:)
211# ifdef LMD_NONLOCAL
212 real(r8), intent(in) :: ghats(LBi:,LBj:,0:,:)
213# endif
214 real(r8), intent(in) :: W(LBi:,LBj:,0:)
215# ifdef WEC_VF
216 real(r8), intent(in) :: W_stokes(LBi:,LBj:,0:)
217# endif
218 real(r8), intent(in) :: ru(LBi:,LBj:,0:,:)
219 real(r8), intent(in) :: rv(LBi:,LBj:,0:,:)
220# ifdef DIAGNOSTICS_TS
221 real(r8), intent(inout) :: DiaTwrk(LBi:,LBj:,:,:,:)
222# endif
223# ifdef DIAGNOSTICS_UV
224 real(r8), intent(inout) :: DiaU3wrk(LBi:,LBj:,:,:)
225 real(r8), intent(inout) :: DiaV3wrk(LBi:,LBj:,:,:)
226 real(r8), intent(inout) :: DiaRU(LBi:,LBj:,:,:,:)
227 real(r8), intent(inout) :: DiaRV(LBi:,LBj:,:,:,:)
228# endif
229# ifdef SUN
230 real(r8), intent(inout) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
231# else
232 real(r8), intent(inout) :: t(LBi:,LBj:,:,:,:)
233# endif
234 real(r8), intent(inout) :: u(LBi:,LBj:,:,:)
235 real(r8), intent(inout) :: v(LBi:,LBj:,:,:)
236
237# else
238
239# ifdef MASKING
240 real(r8), intent(in) :: rmask(LBi:UBi,LBj:UBj)
241 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
242 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
243# if defined SOLAR_SOURCE && defined WET_DRY
244 real(r8), intent(in) :: rmask_wet(LBi:UBi,LBj:UBj)
245# endif
246# endif
247 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
248 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
249 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
250 real(r8), intent(in) :: Huon(LBi:UBi,LBj:UBj,N(ng))
251 real(r8), intent(in) :: Hvom(LBi:UBi,LBj:UBj,N(ng))
252 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
253 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
254 real(r8), intent(in) :: btflx(LBi:UBi,LBj:UBj,NT(ng))
255 real(r8), intent(in) :: bustr(LBi:UBi,LBj:UBj)
256 real(r8), intent(in) :: bvstr(LBi:UBi,LBj:UBj)
257 real(r8), intent(in) :: stflx(LBi:UBi,LBj:UBj,NT(ng))
258 real(r8), intent(in) :: sustr(LBi:UBi,LBj:UBj)
259 real(r8), intent(in) :: svstr(LBi:UBi,LBj:UBj)
260# ifdef SOLAR_SOURCE
261 real(r8), intent(in) :: srflx(LBi:UBi,LBj:UBj)
262# endif
263 real(r8), intent(in) :: Akt(LBi:UBi,LBj:UBj,0:N(ng),NAT)
264 real(r8), intent(in) :: Akv(LBi:UBi,LBj:UBj,0:N(ng))
265# ifdef LMD_NONLOCAL
266 real(r8), intent(in) :: ghats(LBi:UBi,LBj:UBj,0:N(ng),NAT)
267# endif
268 real(r8), intent(in) :: W(LBi:UBi,LBj:UBj,0:N(ng))
269# ifdef WEC_VF
270 real(r8), intent(in) :: W_stokes(LBi:UBi,LBj:UBj,0:N(ng))
271# endif
272 real(r8), intent(in) :: ru(LBi:UBi,LBj:UBj,0:N(ng),2)
273 real(r8), intent(in) :: rv(LBi:UBi,LBj:UBj,0:N(ng),2)
274# ifdef DIAGNOSTICS_TS
275 real(r8), intent(inout) :: DiaTwrk(LBi:UBi,LBj:UBj,N(ng),NT(ng), &
276 & NDT)
277# endif
278# ifdef DIAGNOSTICS_UV
279 real(r8), intent(inout) :: DiaU3wrk(LBi:UBi,LBj:UBj,N(ng),NDM3d)
280 real(r8), intent(inout) :: DiaV3wrk(LBi:UBi,LBj:UBj,N(ng),NDM3d)
281 real(r8), intent(inout) :: DiaRU(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
282 real(r8), intent(inout) :: DiaRV(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
283# endif
284 real(r8), intent(inout) :: t(LBi:UBi,LBj:UBj,N(ng),3,NT(ng))
285 real(r8), intent(inout) :: u(LBi:UBi,LBj:UBj,N(ng),2)
286 real(r8), intent(inout) :: v(LBi:UBi,LBj:UBj,N(ng),2)
287# endif
288!
289! Local variable declarations.
290!
291 integer :: Isrc, Jsrc
292 integer :: i, ic, indx, is, itrc, j, k, ltrc
293# if defined AGE_MEAN && defined T_PASSIVE
294 integer :: iage
295# endif
296# if defined DIAGNOSTICS_TS || defined DIAGNOSTICS_UV
297 integer :: idiag
298# endif
299 real(r8), parameter :: eps = 1.0e-16_r8
300
301 real(r8) :: cff, cff1, cff2, cff3, cff4
302 real(r8) :: Gamma
303
304 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: CF
305 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: DC
306 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: FC
307
308# ifdef SOLAR_SOURCE
309 real(r8), dimension(IminS:ImaxS,JminS:JmaxS,0:N(ng)) :: swdk
310# endif
311 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FE
312 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: FX
313 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: curv
314 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: grad
315
316# include "set_bounds.h"
317
318# ifndef TS_FIXED
319!
320!=======================================================================
321! Tracer equation(s).
322!=======================================================================
323
324# ifdef SOLAR_SOURCE
325!
326! Compute fraction of the solar shortwave radiation, "swdk"
327! (at vertical W-points) penetrating water column.
328!
329 DO k=1,n(ng)-1
330 DO j=jstr,jend
331 DO i=istr,iend
332 fx(i,j)=z_w(i,j,n(ng))-z_w(i,j,k)
333 END DO
334 END DO
335 CALL lmd_swfrac_tile (ng, tile, &
336 & lbi, ubi, lbj, ubj, &
337 & imins, imaxs, jmins, jmaxs, &
338 & -1.0_r8, fx, fe)
339 DO j=jstr,jend
340 DO i=istr,iend
341 swdk(i,j,k)=fe(i,j)
342 END DO
343 END DO
344 END DO
345# endif
346!
347!-----------------------------------------------------------------------
348! Compute intermediate tracer at n+1/2 time-step, t(i,j,k,3,itrc).
349!-----------------------------------------------------------------------
350!
351! Compute time rate of change of intermediate tracer due to
352! horizontal advection.
353!
354 t_loop1 :DO itrc=1,nt(ng)
355 k_loop: DO k=1,n(ng)
356!
357 hadv_flux : IF (hadvection(itrc,ng)%CENTERED2) THEN
358!
359! Second-order, centered differences horizontal advective fluxes.
360!
361 DO j=jstr,jend
362 DO i=istr,iend+1
363 fx(i,j)=huon(i,j,k)* &
364 & 0.5_r8*(t(i-1,j,k,nstp,itrc)+ &
365 & t(i ,j,k,nstp,itrc))
366 END DO
367 END DO
368 DO j=jstr,jend+1
369 DO i=istr,iend
370 fe(i,j)=hvom(i,j,k)* &
371 & 0.5_r8*(t(i,j-1,k,nstp,itrc)+ &
372 & t(i,j ,k,nstp,itrc))
373 END DO
374 END DO
375!
376 ELSE IF ((hadvection(itrc,ng)%MPDATA).or. &
377 & (hadvection(itrc,ng)%HSIMT)) THEN
378!
379! First-order, upstream differences horizontal advective fluxes.
380!
381 DO j=jstr,jend
382 DO i=istr,iend+1
383 cff1=max(huon(i,j,k),0.0_r8)
384 cff2=min(huon(i,j,k),0.0_r8)
385 fx(i,j)=cff1*t(i-1,j,k,nstp,itrc)+ &
386 & cff2*t(i ,j,k,nstp,itrc)
387 END DO
388 END DO
389 DO j=jstr,jend+1
390 DO i=istr,iend
391 cff1=max(hvom(i,j,k),0.0_r8)
392 cff2=min(hvom(i,j,k),0.0_r8)
393 fe(i,j)=cff1*t(i,j-1,k,nstp,itrc)+ &
394 & cff2*t(i,j ,k,nstp,itrc)
395 END DO
396 END DO
397!
398 ELSE IF ((hadvection(itrc,ng)%AKIMA4).or. &
399 & (hadvection(itrc,ng)%CENTERED4).or. &
400 & (hadvection(itrc,ng)%SPLIT_U3).or. &
401 & (hadvection(itrc,ng)%UPSTREAM3)) THEN
402!
403! Fourth-order Akima, fourth-order centered differences, or third-order
404! upstream-biased horizontal advective fluxes.
405!
406 DO j=jstr,jend
407 DO i=istrm1,iendp2
408 fx(i,j)=t(i ,j,k,nstp,itrc)- &
409 & t(i-1,j,k,nstp,itrc)
410# ifdef MASKING
411 fx(i,j)=fx(i,j)*umask(i,j)
412# endif
413 END DO
414 END DO
415 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
416 IF (domain(ng)%Western_Edge(tile)) THEN
417 DO j=jstr,jend
418 fx(istr-1,j)=fx(istr,j)
419 END DO
420 END IF
421 END IF
422 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
423 IF (domain(ng)%Eastern_Edge(tile)) THEN
424 DO j=jstr,jend
425 fx(iend+2,j)=fx(iend+1,j)
426 END DO
427 END IF
428 END IF
429!
430 DO j=jstr,jend
431 DO i=istr-1,iend+1
432 IF (hadvection(itrc,ng)%UPSTREAM3) THEN
433 curv(i,j)=fx(i+1,j)-fx(i,j)
434 ELSE IF (hadvection(itrc,ng)%AKIMA4) THEN
435 cff=2.0_r8*fx(i+1,j)*fx(i,j)
436 IF (cff.gt.eps) THEN
437 grad(i,j)=cff/(fx(i+1,j)+fx(i,j))
438 ELSE
439 grad(i,j)=0.0_r8
440 END IF
441 ELSE IF ((hadvection(itrc,ng)%CENTERED4).or. &
442 & (hadvection(itrc,ng)%SPLIT_U3)) THEN
443 grad(i,j)=0.5_r8*(fx(i+1,j)+fx(i,j))
444 END IF
445 END DO
446 END DO
447!
448 cff1=1.0_r8/6.0_r8
449 cff2=1.0_r8/3.0_r8
450 DO j=jstr,jend
451 DO i=istr,iend+1
452 IF (hadvection(itrc,ng)%UPSTREAM3) THEN
453 fx(i,j)=huon(i,j,k)*0.5_r8* &
454 & (t(i-1,j,k,nstp,itrc)+ &
455 & t(i ,j,k,nstp,itrc))- &
456 & cff1*(curv(i-1,j)*max(huon(i,j,k),0.0_r8)+ &
457 & curv(i ,j)*min(huon(i,j,k),0.0_r8))
458 ELSE IF ((hadvection(itrc,ng)%AKIMA4).or. &
459 & (hadvection(itrc,ng)%CENTERED4).or. &
460 & (hadvection(itrc,ng)%SPLIT_U3)) THEN
461
462 fx(i,j)=huon(i,j,k)*0.5_r8* &
463 & (t(i-1,j,k,nstp,itrc)+ &
464 & t(i ,j,k,nstp,itrc)- &
465 & cff2*(grad(i ,j)- &
466 & grad(i-1,j)))
467 END IF
468 END DO
469 END DO
470!
471 DO j=jstrm1,jendp2
472 DO i=istr,iend
473 fe(i,j)=t(i,j ,k,nstp,itrc)- &
474 & t(i,j-1,k,nstp,itrc)
475# ifdef MASKING
476 fe(i,j)=fe(i,j)*vmask(i,j)
477# endif
478 END DO
479 END DO
480 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
481 IF (domain(ng)%Southern_Edge(tile)) THEN
482 DO i=istr,iend
483 fe(i,jstr-1)=fe(i,jstr)
484 END DO
485 END IF
486 END IF
487 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
488 IF (domain(ng)%Northern_Edge(tile)) THEN
489 DO i=istr,iend
490 fe(i,jend+2)=fe(i,jend+1)
491 END DO
492 END IF
493 END IF
494!
495 DO j=jstr-1,jend+1
496 DO i=istr,iend
497 IF (hadvection(itrc,ng)%UPSTREAM3) THEN
498 curv(i,j)=fe(i,j+1)-fe(i,j)
499 ELSE IF (hadvection(itrc,ng)%AKIMA4) THEN
500 cff=2.0_r8*fe(i,j+1)*fe(i,j)
501 IF (cff.gt.eps) THEN
502 grad(i,j)=cff/(fe(i,j+1)+fe(i,j))
503 ELSE
504 grad(i,j)=0.0_r8
505 END IF
506 ELSE IF ((hadvection(itrc,ng)%CENTERED4).or. &
507 & (hadvection(itrc,ng)%SPLIT_U3)) THEN
508 grad(i,j)=0.5_r8*(fe(i,j+1)+fe(i,j))
509 END IF
510 END DO
511 END DO
512!
513 cff1=1.0_r8/6.0_r8
514 cff2=1.0_r8/3.0_r8
515 DO j=jstr,jend+1
516 DO i=istr,iend
517 IF (hadvection(itrc,ng)%UPSTREAM3) THEN
518 fe(i,j)=hvom(i,j,k)*0.5_r8* &
519 & (t(i,j-1,k,nstp,itrc)+ &
520 & t(i,j ,k,nstp,itrc))- &
521 & cff1*(curv(i,j-1)*max(hvom(i,j,k),0.0_r8)+ &
522 & curv(i,j )*min(hvom(i,j,k),0.0_r8))
523 ELSE IF ((hadvection(itrc,ng)%AKIMA4).or. &
524 & (hadvection(itrc,ng)%CENTERED4).or. &
525 & (hadvection(itrc,ng)%SPLIT_U3)) THEN
526 fe(i,j)=hvom(i,j,k)*0.5_r8* &
527 & (t(i,j-1,k,nstp,itrc)+ &
528 & t(i,j ,k,nstp,itrc)- &
529 & cff2*(grad(i,j )- &
530 & grad(i,j-1)))
531 END IF
532 END DO
533 END DO
534 END IF hadv_flux
535!
536! Apply tracers point sources to the horizontal advection terms,
537! if any.
538!
539! Dsrc(is) = 0, flow across grid cell u-face (positive or negative)
540! Dsrc(is) = 1, flow across grid cell v-face (positive or negative)
541!
542 IF (luvsrc(ng)) THEN
543 DO is=1,nsrc(ng)
544 isrc=sources(ng)%Isrc(is)
545 jsrc=sources(ng)%Jsrc(is)
546 IF (((istr.le.isrc).and.(isrc.le.iend+1)).and. &
547 & ((jstr.le.jsrc).and.(jsrc.le.jend+1))) THEN
548 IF (int(sources(ng)%Dsrc(is)).eq.0) THEN
549 IF (ltracersrc(itrc,ng)) THEN
550 fx(isrc,jsrc)=huon(isrc,jsrc,k)* &
551 & sources(ng)%Tsrc(is,k,itrc)
552 ELSE
553 fx(isrc,jsrc)=0.0_r8
554 END IF
555 ELSE IF (int(sources(ng)%Dsrc(is)).eq.1) THEN
556 IF (ltracersrc(itrc,ng)) THEN
557 fe(isrc,jsrc)=hvom(isrc,jsrc,k)* &
558 & sources(ng)%Tsrc(is,k,itrc)
559 ELSE
560 fe(isrc,jsrc)=0.0_r8
561 END IF
562 END IF
563 END IF
564 END DO
565 END IF
566!
567! Time-step horizontal advection (m Tunits).
568!
569 IF ((hadvection(itrc,ng)%MPDATA).or. &
570 & (hadvection(itrc,ng)%HSIMT)) THEN
571 gamma=0.5_r8
572 ELSE
573 gamma=1.0_r8/6.0_r8
574 END IF
575 IF (iic(ng).eq.ntfirst(ng)) THEN
576 cff=0.5_r8*dt(ng)
577 cff1=1.0_r8
578 cff2=0.0_r8
579 ELSE
580 cff=(1.0_r8-gamma)*dt(ng)
581 cff1=0.5_r8+gamma
582 cff2=0.5_r8-gamma
583 END IF
584 DO j=jstr,jend
585 DO i=istr,iend
586 t(i,j,k,3,itrc)=hz(i,j,k)*(cff1*t(i,j,k,nstp,itrc)+ &
587 & cff2*t(i,j,k,nnew,itrc))- &
588 & cff*pm(i,j)*pn(i,j)* &
589 & (fx(i+1,j)-fx(i,j)+ &
590 & fe(i,j+1)-fe(i,j))
591 END DO
592 END DO
593 END DO k_loop
594 END DO t_loop1
595
596# if defined AGE_MEAN && defined T_PASSIVE
597!
598! If inert passive tracer and Mean Age, compute age concentration (even
599! inert index) forced by the right-hand-side term that is concentration
600! of an associated conservative passive tracer (odd inert index).
601! Implemented and tested by W.G. Zhang and J. Wilkin. (m Tunits)
602!
603 DO itrc=1,npt,2
604 IF (.not.((hadvection(inert(itrc),ng)%MPDATA).or. &
605 & (hadvection(inert(itrc),ng)%HSIMT))) THEN
606 IF (iic(ng).eq.ntfirst(ng)) THEN
607 cff=0.5_r8*dt(ng)
608 ELSE
609 gamma=1.0_r8/6.0_r8
610 cff=(1.0_r8-gamma)*dt(ng)
611 END IF
612 iage=inert(itrc+1) ! even inert tracer index
613 DO k=1,n(ng)
614 DO j=jstr,jend
615 DO i=istr,iend
616 t(i,j,k,3,iage)=t(i,j,k,3,iage)+ &
617 & cff*hz(i,j,k)* &
618 & t(i,j,k,nnew,inert(itrc))
619 END DO
620 END DO
621 END DO
622 END IF
623 END DO
624# endif
625!
626!-----------------------------------------------------------------------
627! Compute time rate of change of intermediate tracer due to vertical
628! advection. Impose artificial continuity equation.
629!-----------------------------------------------------------------------
630!
631 j_loop1 : DO j=jstr,jend
632 t_loop2 : DO itrc=1,nt(ng)
633!
634 vadv_flux : IF (vadvection(itrc,ng)%SPLINES) THEN
635!
636! Build conservative parabolic splines for the vertical derivatives
637! "FC" of the tracer. Then, the interfacial "FC" values are
638! converted to vertical advective flux.
639!
640 DO i=istr,iend
641# ifdef NEUMANN
642 fc(i,0)=1.5_r8*t(i,j,1,nstp,itrc)
643 cf(i,1)=0.5_r8
644# else
645 fc(i,0)=2.0_r8*t(i,j,1,nstp,itrc)
646 cf(i,1)=1.0_r8
647# endif
648 END DO
649 DO k=1,n(ng)-1
650 DO i=istr,iend
651 cff=1.0_r8/(2.0_r8*hz(i,j,k)+ &
652 & hz(i,j,k+1)*(2.0_r8-cf(i,k)))
653 cf(i,k+1)=cff*hz(i,j,k)
654 fc(i,k)=cff*(3.0_r8*(hz(i,j,k )*t(i,j,k+1,nstp,itrc)+ &
655 & hz(i,j,k+1)*t(i,j,k ,nstp,itrc))- &
656 & hz(i,j,k+1)*fc(i,k-1))
657 END DO
658 END DO
659 DO i=istr,iend
660# ifdef NEUMANN
661 fc(i,n(ng))=(3.0_r8*t(i,j,n(ng),nstp,itrc)- &
662 & fc(i,n(ng)-1))/(2.0_r8-cf(i,n(ng)))
663# else
664 fc(i,n(ng))=(2.0_r8*t(i,j,n(ng),nstp,itrc)- &
665 & fc(i,n(ng)-1))/(1.0_r8-cf(i,n(ng)))
666# endif
667 END DO
668 DO k=n(ng)-1,0,-1
669 DO i=istr,iend
670 fc(i,k)=fc(i,k)-cf(i,k+1)*fc(i,k+1)
671# ifdef WEC_VF
672 fc(i,k+1)=(w(i,j,k+1)+w_stokes(i,j,k+1))*fc(i,k+1)
673# else
674 fc(i,k+1)=w(i,j,k+1)*fc(i,k+1)
675# endif
676 END DO
677 END DO
678 DO i=istr,iend
679 fc(i,n(ng))=0.0_r8
680 fc(i,0)=0.0_r8
681 END DO
682!
683 ELSE IF (vadvection(itrc,ng)%AKIMA4) THEN
684!
685! Fourth-order, Akima vertical advective flux.
686!
687 DO k=1,n(ng)-1
688 DO i=istr,iend
689 fc(i,k)=t(i,j,k+1,nstp,itrc)- &
690 & t(i,j,k ,nstp,itrc)
691 END DO
692 END DO
693 DO i=istr,iend
694 fc(i,0)=fc(i,1)
695 fc(i,n(ng))=fc(i,n(ng)-1)
696 END DO
697 DO k=1,n(ng)
698 DO i=istr,iend
699 cff=2.0_r8*fc(i,k)*fc(i,k-1)
700 IF (cff.gt.eps) THEN
701 cf(i,k)=cff/(fc(i,k)+fc(i,k-1))
702 ELSE
703 cf(i,k)=0.0_r8
704 END IF
705 END DO
706 END DO
707 cff1=1.0_r8/3.0_r8
708 DO k=1,n(ng)-1
709 DO i=istr,iend
710# ifdef WEC_VF
711 fc(i,k)=(w(i,j,k)+w_stokes(i,j,k))* &
712# else
713 fc(i,k)=w(i,j,k)* &
714# endif
715 & 0.5_r8*(t(i,j,k ,nstp,itrc)+ &
716 & t(i,j,k+1,nstp,itrc)- &
717 & cff1*(cf(i,k+1)-cf(i,k)))
718 END DO
719 END DO
720 DO i=istr,iend
721 fc(i,0)=0.0_r8
722 fc(i,n(ng))=0.0_r8
723 END DO
724!
725 ELSE IF (vadvection(itrc,ng)%CENTERED2) THEN
726!
727! Second-order, central differences vertical advective flux.
728!
729 DO k=1,n(ng)-1
730 DO i=istr,iend
731# ifdef WEC_VF
732 fc(i,k)=(w(i,j,k)+w_stokes(i,j,k))* &
733# else
734 fc(i,k)=w(i,j,k)* &
735# endif
736 & 0.5_r8*(t(i,j,k ,nstp,itrc)+ &
737 & t(i,j,k+1,nstp,itrc))
738 END DO
739 END DO
740 DO i=istr,iend
741 fc(i,0)=0.0_r8
742 fc(i,n(ng))=0.0_r8
743 END DO
744!
745 ELSE IF ((vadvection(itrc,ng)%MPDATA).or. &
746 & (vadvection(itrc,ng)%HSIMT)) THEN
747!
748! First_order, upstream differences vertical advective flux.
749!
750 DO k=1,n(ng)-1
751 DO i=istr,iend
752# if defined WEC_VF
753 cff1=max(w(i,j,k)+w_stokes(i,j,k),0.0_r8)
754 cff2=min(w(i,j,k)+w_stokes(i,j,k),0.0_r8)
755# else
756 cff1=max(w(i,j,k),0.0_r8)
757 cff2=min(w(i,j,k),0.0_r8)
758# endif
759 fc(i,k)=cff1*t(i,j,k ,nstp,itrc)+ &
760 & cff2*t(i,j,k+1,nstp,itrc)
761 END DO
762 END DO
763 DO i=istr,iend
764 fc(i,0)=0.0_r8
765 fc(i,n(ng))=0.0_r8
766 END DO
767!
768 ELSE IF ((vadvection(itrc,ng)%CENTERED4).or. &
769 & (vadvection(itrc,ng)%SPLIT_U3)) THEN
770!
771! Fourth-order, central differences vertical advective flux.
772!
773 cff1=0.5_r8
774 cff2=7.0_r8/12.0_r8
775 cff3=1.0_r8/12.0_r8
776 DO k=2,n(ng)-2
777 DO i=istr,iend
778# ifdef WEC_VF
779 fc(i,k)=(w(i,j,k)+w_stokes(i,j,k))* &
780# else
781 fc(i,k)=w(i,j,k)* &
782# endif
783 & (cff2*(t(i,j,k ,nstp,itrc)+ &
784 & t(i,j,k+1,nstp,itrc))- &
785 & cff3*(t(i,j,k-1,nstp,itrc)+ &
786 & t(i,j,k+2,nstp,itrc)))
787 END DO
788 END DO
789 DO i=istr,iend
790 fc(i,0)=0.0_r8
791# ifdef WEC_VF
792 fc(i,1)=(w(i,j,1)+w_stokes(i,j,1))* &
793# else
794 fc(i,1)=w(i,j,1)* &
795# endif
796 & (cff1*t(i,j,1,nstp,itrc)+ &
797 & cff2*t(i,j,2,nstp,itrc)- &
798 & cff3*t(i,j,3,nstp,itrc))
799# ifdef WEC_VF
800 fc(i,n(ng)-1)=(w(i,j,n(ng)-1)+w_stokes(i,j,n(ng)-1))* &
801# else
802 fc(i,n(ng)-1)=w(i,j,n(ng)-1)* &
803# endif
804 & (cff1*t(i,j,n(ng) ,nstp,itrc)+ &
805 & cff2*t(i,j,n(ng)-1,nstp,itrc)- &
806 & cff3*t(i,j,n(ng)-2,nstp,itrc))
807 fc(i,n(ng))=0.0_r8
808 END DO
809 END IF vadv_flux
810!
811! Compute artificial continuity equation and load it into private
812! array DC (1/m). It is imposed to preserve tracer constancy. It is
813! not the same for all the tracer advection schemes because of the
814! Gamma value.
815!
816 IF ((vadvection(itrc,ng)%MPDATA).or. &
817 & (vadvection(itrc,ng)%HSIMT)) THEN
818 gamma=0.5_r8
819 ELSE
820 gamma=1.0_r8/6.0_r8
821 END IF
822 IF (iic(ng).eq.ntfirst(ng)) THEN
823 cff=0.5_r8*dt(ng)
824 ELSE
825 cff=(1.0_r8-gamma)*dt(ng)
826 END IF
827 DO k=1,n(ng)
828 DO i=istr,iend
829 dc(i,k)=1.0_r8/(hz(i,j,k)- &
830 & cff*pm(i,j)*pn(i,j)* &
831 & (huon(i+1,j,k)-huon(i,j,k)+ &
832 & hvom(i,j+1,k)-hvom(i,j,k)+ &
833# ifdef WEC_VF
834 & (w_stokes(i,j,k)-w_stokes(i,j,k-1))+ &
835# endif
836 & (w(i,j,k)-w(i,j,k-1))))
837 END DO
838 END DO
839!
840! Time-step vertical advection of tracers (Tunits). Impose artificial
841! continuity equation.
842!
843 DO k=1,n(ng)
844 DO i=istr,iend
845 cff1=cff*pm(i,j)*pn(i,j)
846 t(i,j,k,3,itrc)=dc(i,k)* &
847 & (t(i,j,k,3,itrc)- &
848 & cff1*(fc(i,k)-fc(i,k-1)))
849 END DO
850 END DO
851 END DO t_loop2
852 END DO j_loop1
853!
854!-----------------------------------------------------------------------
855! Start computation of tracers at n+1 time-step, t(i,j,k,nnew,itrc).
856!-----------------------------------------------------------------------
857!
858! Compute vertical diffusive fluxes "FC" of the tracer fields at
859! current time step n, and at horizontal RHO-points and vertical
860! W-points. Notice that the vertical diffusion coefficients for
861! passive tracers is the same as that for salinity (ltrc=NAT).
862!
863 DO j=jstr,jend
864 cff3=dt(ng)*(1.0_r8-lambda)
865 DO itrc=1,nt(ng)
866 ltrc=min(nat,itrc)
867 DO k=1,n(ng)-1
868 DO i=istr,iend
869 cff=1.0_r8/(z_r(i,j,k+1)-z_r(i,j,k))
870 fc(i,k)=cff3*cff*akt(i,j,k,ltrc)* &
871 & (t(i,j,k+1,nstp,itrc)- &
872 & t(i,j,k ,nstp,itrc))
873 END DO
874 END DO
875
876# ifdef LMD_NONLOCAL
877!
878! Add in the nonlocal transport flux for unstable (convective)
879! forcing conditions into matrix FC when using the Large et al.
880! KPP scheme. The nonlocal transport is only applied to active
881! tracers.
882!
883 IF (itrc.le.nat) THEN
884 DO k=1,n(ng)-1
885 DO i=istr,iend
886 fc(i,k)=fc(i,k)- &
887 & dt(ng)*akt(i,j,k,itrc)*ghats(i,j,k,itrc)
888 END DO
889 END DO
890 END IF
891# endif
892# ifdef SOLAR_SOURCE
893!
894! Add in incoming solar radiation at interior W-points using decay
895! decay penetration function based on Jerlov water type.
896!
897 IF (itrc.eq.itemp) THEN
898 DO k=1,n(ng)-1
899 DO i=istr,iend
900 fc(i,k)=fc(i,k)+ &
901 & dt(ng)*srflx(i,j)* &
902# ifdef WET_DRY
903 & rmask_wet(i,j)* &
904# endif
905 & swdk(i,j,k)
906 END DO
907 END DO
908 END IF
909# endif
910!
911! Apply bottom and surface tracer flux conditions.
912!
913 DO i=istr,iend
914 fc(i,0)=dt(ng)*btflx(i,j,itrc)
915 fc(i,n(ng))=dt(ng)*stflx(i,j,itrc)
916 END DO
917!
918! Compute new tracer field (m Tunits).
919!
920 DO k=1,n(ng)
921 DO i=istr,iend
922 cff1=hz(i,j,k)*t(i,j,k,nstp,itrc)
923 cff2=fc(i,k)-fc(i,k-1)
924 t(i,j,k,nnew,itrc)=cff1+cff2
925# ifdef DIAGNOSTICS_TS
926 diatwrk(i,j,k,itrc,itrate)=cff1
927 diatwrk(i,j,k,itrc,itvdif)=cff2
928# endif
929 END DO
930 END DO
931 END DO
932 END DO
933# endif /* !TS_FIXED */
934!
935!=======================================================================
936! 3D momentum equation in the XI-direction.
937!=======================================================================
938!
939! Compute U-component viscous vertical momentum fluxes "FC" at
940! current time-step n, and at horizontal U-points and vertical
941! W-points.
942!
943 j_loop2: DO j=jstr,jend
944 cff3=dt(ng)*(1.0_r8-lambda)
945 DO k=1,n(ng)-1
946 DO i=istru,iend
947 cff=1.0_r8/(z_r(i,j,k+1)+z_r(i-1,j,k+1)- &
948 & z_r(i,j,k )-z_r(i-1,j,k ))
949 fc(i,k)=cff3*cff*(u(i,j,k+1,nstp)-u(i,j,k,nstp))* &
950 & (akv(i,j,k)+akv(i-1,j,k))
951 END DO
952 END DO
953!
954! Apply bottom and surface stresses, if so is prescribed.
955!
956 DO i=istru,iend
957# ifdef BODYFORCE
958 fc(i,0)=0.0_r8
959 fc(i,n(ng))=0.0_r8
960# else
961 fc(i,0)=dt(ng)*bustr(i,j)
962 fc(i,n(ng))=dt(ng)*sustr(i,j)
963# endif
964 END DO
965!
966! Compute new U-momentum (m m/s).
967!
968 cff=dt(ng)*0.25_r8
969 DO i=istru,iend
970 dc(i,0)=cff*(pm(i,j)+pm(i-1,j))*(pn(i,j)+pn(i-1,j))
971 END DO
972 indx=3-nrhs
973 IF (iic(ng).eq.ntfirst(ng)) THEN
974 DO k=1,n(ng)
975 DO i=istru,iend
976 cff1=u(i,j,k,nstp)*0.5_r8*(hz(i,j,k)+hz(i-1,j,k))
977 cff2=fc(i,k)-fc(i,k-1)
978 u(i,j,k,nnew)=cff1+cff2
979# ifdef DIAGNOSTICS_UV
980 DO idiag=1,m3pgrd
981 diau3wrk(i,j,k,idiag)=0.0_r8
982 END DO
983 diau3wrk(i,j,k,m3vvis)=cff2
984 diau3wrk(i,j,k,m3rate)=cff1
985# endif
986 END DO
987 END DO
988 ELSE IF (iic(ng).eq.(ntfirst(ng)+1)) THEN
989 DO k=1,n(ng)
990 DO i=istru,iend
991 cff1=u(i,j,k,nstp)*0.5_r8*(hz(i,j,k)+hz(i-1,j,k))
992 cff2=fc(i,k)-fc(i,k-1)
993 cff3=0.5_r8*dc(i,0)
994 u(i,j,k,nnew)=cff1- &
995 & cff3*ru(i,j,k,indx)+ &
996 & cff2
997# ifdef DIAGNOSTICS_UV
998 DO idiag=1,m3pgrd
999 diau3wrk(i,j,k,idiag)=-cff3*diaru(i,j,k,indx,idiag)
1000 END DO
1001 diau3wrk(i,j,k,m3vvis)=cff2
1002# ifdef BODYFORCE
1003 diau3wrk(i,j,k,m3vvis)=diau3wrk(i,j,k,m3vvis)- &
1004 & cff3*diaru(i,j,k,indx,m3vvis)
1005# endif
1006 diau3wrk(i,j,k,m3rate)=cff1
1007# endif
1008 END DO
1009 END DO
1010 ELSE
1011 cff1= 5.0_r8/12.0_r8
1012 cff2=16.0_r8/12.0_r8
1013 DO k=1,n(ng)
1014 DO i=istru,iend
1015 cff3=u(i,j,k,nstp)*0.5_r8*(hz(i,j,k)+hz(i-1,j,k))
1016 cff4=fc(i,k)-fc(i,k-1)
1017 u(i,j,k,nnew)=cff3+ &
1018 & dc(i,0)*(cff1*ru(i,j,k,nrhs)- &
1019 & cff2*ru(i,j,k,indx))+ &
1020 & cff4
1021# ifdef DIAGNOSTICS_UV
1022 DO idiag=1,m3pgrd
1023 diau3wrk(i,j,k,idiag)=dc(i,0)* &
1024 & (cff1*diaru(i,j,k,nrhs,idiag)- &
1025 & cff2*diaru(i,j,k,indx,idiag))
1026 END DO
1027 diau3wrk(i,j,k,m3vvis)=cff4
1028# ifdef BODYFORCE
1029 diau3wrk(i,j,k,m3vvis)=diau3wrk(i,j,k,m3vvis)+ &
1030 & dc(i,0)* &
1031 & (cff1*diaru(i,j,k,nrhs,m3vvis)- &
1032 & cff2*diaru(i,j,k,indx,m3vvis))
1033# endif
1034 diau3wrk(i,j,k,m3rate)=cff3
1035# endif
1036 END DO
1037 END DO
1038 END IF
1039!
1040!=======================================================================
1041! 3D momentum equation in the ETA-direction.
1042!=======================================================================
1043!
1044! Compute V-component viscous vertical momentum fluxes "FC" at
1045! current time-step n, and at horizontal V-points and vertical
1046! W-points.
1047!
1048 IF (j.ge.jstrv) THEN
1049 cff3=dt(ng)*(1.0_r8-lambda)
1050 DO k=1,n(ng)-1
1051 DO i=istr,iend
1052 cff=1.0_r8/(z_r(i,j,k+1)+z_r(i,j-1,k+1)- &
1053 & z_r(i,j,k )-z_r(i,j-1,k ))
1054 fc(i,k)=cff3*cff*(v(i,j,k+1,nstp)-v(i,j,k,nstp))* &
1055 & (akv(i,j,k)+akv(i,j-1,k))
1056 END DO
1057 END DO
1058!
1059! Apply bottom and surface stresses, if so is prescribed.
1060!
1061 DO i=istr,iend
1062# ifdef BODYFORCE
1063 fc(i,0)=0.0_r8
1064 fc(i,n(ng))=0.0_r8
1065# else
1066 fc(i,0)=dt(ng)*bvstr(i,j)
1067 fc(i,n(ng))=dt(ng)*svstr(i,j)
1068# endif
1069 END DO
1070!
1071! Compute new V-momentum (m m/s).
1072!
1073 cff=dt(ng)*0.25_r8
1074 DO i=istr,iend
1075 dc(i,0)=cff*(pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))
1076 END DO
1077 IF (iic(ng).eq.ntfirst(ng)) THEN
1078 DO k=1,n(ng)
1079 DO i=istr,iend
1080 cff1=v(i,j,k,nstp)*0.5_r8*(hz(i,j,k)+hz(i,j-1,k))
1081 cff2=fc(i,k)-fc(i,k-1)
1082 v(i,j,k,nnew)=cff1+cff2
1083# ifdef DIAGNOSTICS_UV
1084 DO idiag=1,m3pgrd
1085 diav3wrk(i,j,k,idiag)=0.0_r8
1086 END DO
1087 diav3wrk(i,j,k,m3vvis)=cff2
1088 diav3wrk(i,j,k,m3rate)=cff1
1089# endif
1090 END DO
1091 END DO
1092 ELSE IF (iic(ng).eq.(ntfirst(ng)+1)) THEN
1093 DO k=1,n(ng)
1094 DO i=istr,iend
1095 cff1=v(i,j,k,nstp)*0.5_r8*(hz(i,j,k)+hz(i,j-1,k))
1096 cff2=fc(i,k)-fc(i,k-1)
1097 cff3=0.5_r8*dc(i,0)
1098 v(i,j,k,nnew)=cff1- &
1099 & cff3*rv(i,j,k,indx)+ &
1100 & cff2
1101# ifdef DIAGNOSTICS_UV
1102 DO idiag=1,m3pgrd
1103 diav3wrk(i,j,k,idiag)=-cff3*diarv(i,j,k,indx,idiag)
1104 END DO
1105 diav3wrk(i,j,k,m3vvis)=cff2
1106# ifdef BODYFORCE
1107 diav3wrk(i,j,k,m3vvis)=diav3wrk(i,j,k,m3vvis)- &
1108 & cff3*diarv(i,j,k,indx,m3vvis)
1109# endif
1110 diav3wrk(i,j,k,m3rate)=cff1
1111# endif
1112 END DO
1113 END DO
1114 ELSE
1115 cff1= 5.0_r8/12.0_r8
1116 cff2=16.0_r8/12.0_r8
1117 DO k=1,n(ng)
1118 DO i=istr,iend
1119 cff3=v(i,j,k,nstp)*0.5_r8*(hz(i,j,k)+hz(i,j-1,k))
1120 cff4=fc(i,k)-fc(i,k-1)
1121 v(i,j,k,nnew)=cff3+ &
1122 & dc(i,0)*(cff1*rv(i,j,k,nrhs)- &
1123 & cff2*rv(i,j,k,indx))+ &
1124 & cff4
1125# ifdef DIAGNOSTICS_UV
1126 DO idiag=1,m3pgrd
1127 diav3wrk(i,j,k,idiag)=dc(i,0)* &
1128 & (cff1*diarv(i,j,k,nrhs,idiag)- &
1129 & cff2*diarv(i,j,k,indx,idiag))
1130 END DO
1131 diav3wrk(i,j,k,m3vvis)=cff4
1132# ifdef BODYFORCE
1133 diav3wrk(i,j,k,m3vvis)=diav3wrk(i,j,k,m3vvis)+ &
1134 & dc(i,0)* &
1135 & (cff1*diarv(i,j,k,nrhs,m3vvis)- &
1136 & cff2*diarv(i,j,k,indx,m3vvis))
1137# endif
1138 diav3wrk(i,j,k,m3rate)=cff3
1139# endif
1140 END DO
1141 END DO
1142 END IF
1143 END IF
1144 END DO j_loop2
1145
1146# ifndef TS_FIXED
1147!
1148!=======================================================================
1149! Apply intermediate tracers lateral boundary conditions.
1150!=======================================================================
1151!
1152 ic=0
1153 DO itrc=1,nt(ng)
1154 IF (ltracerclm(itrc,ng).and.lnudgetclm(itrc,ng)) THEN
1155 ic=ic+1
1156 END IF
1157 CALL t3dbc_tile (ng, tile, itrc, ic, &
1158 & lbi, ubi, lbj, ubj, n(ng), nt(ng), &
1159 & imins, imaxs, jmins, jmaxs, &
1160 & nstp, 3, &
1161 & t)
1162
1163 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
1164 CALL exchange_r3d_tile (ng, tile, &
1165 & lbi, ubi, lbj, ubj, 1, n(ng), &
1166 & t(:,:,:,3,itrc))
1167 END IF
1168 END DO
1169
1170# ifdef DISTRIBUTE
1171 CALL mp_exchange4d (ng, tile, inlm, 1, &
1172 & lbi, ubi, lbj, ubj, 1, n(ng), 1, nt(ng), &
1173 & nghostpoints, &
1174 & ewperiodic(ng), nsperiodic(ng), &
1175 & t(:,:,:,3,:))
1176# endif
1177# endif
1178!
1179 RETURN
subroutine lmd_swfrac_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, zscale, z, swdk)
Definition lmd_swfrac.F:10
subroutine exchange_r3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
type(t_adv), dimension(:,:), allocatable hadvection
Definition mod_param.F:403
integer nat
Definition mod_param.F:499
integer, dimension(:), allocatable n
Definition mod_param.F:479
integer nghostpoints
Definition mod_param.F:710
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
type(t_adv), dimension(:,:), allocatable vadvection
Definition mod_param.F:404
integer, dimension(:), allocatable nt
Definition mod_param.F:489
integer npt
Definition mod_param.F:505
integer m3vvis
logical, dimension(:), allocatable luvsrc
logical, dimension(:,:), allocatable ltracersrc
integer, dimension(:), allocatable iic
real(dp), dimension(:), allocatable dt
real(r8) lambda
logical, dimension(:), allocatable ewperiodic
integer, parameter iwest
logical, dimension(:), allocatable nsperiodic
integer itrate
logical, dimension(:,:), allocatable compositegrid
integer itemp
integer m3rate
integer, parameter isouth
integer itvdif
integer, dimension(:), pointer inert
integer, dimension(:), allocatable ntfirst
integer, parameter ieast
logical, dimension(:,:), allocatable ltracerclm
integer, parameter inorth
integer m3pgrd
logical, dimension(:,:), allocatable lnudgetclm
type(t_sources), dimension(:), allocatable sources
Definition mod_sources.F:90
integer, dimension(:), allocatable nsrc
Definition mod_sources.F:97
subroutine mp_exchange4d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, lbt, ubt, nghost, ew_periodic, ns_periodic, a, b, c)
subroutine, public t3dbc_tile(ng, tile, itrc, ic, lbi, ubi, lbj, ubj, ubk, ubt, imins, imaxs, jmins, jmaxs, nstp, nout, t)
Definition t3dbc_im.F:55

References mod_scalars::compositegrid, mod_param::domain, mod_scalars::dt, mod_scalars::ewperiodic, exchange_3d_mod::exchange_r3d_tile(), mod_param::hadvection, mod_scalars::ieast, mod_scalars::iic, mod_scalars::inert, mod_param::inlm, mod_scalars::inorth, mod_scalars::isouth, mod_scalars::itemp, mod_scalars::itrate, mod_scalars::itvdif, mod_scalars::iwest, mod_scalars::lambda, lmd_swfrac_tile(), mod_scalars::lnudgetclm, mod_scalars::ltracerclm, mod_scalars::ltracersrc, mod_scalars::luvsrc, mod_scalars::m3pgrd, mod_scalars::m3rate, mod_scalars::m3vvis, mp_exchange_mod::mp_exchange4d(), mod_param::nghostpoints, mod_param::npt, mod_scalars::nsperiodic, mod_sources::nsrc, mod_scalars::ntfirst, mod_sources::sources, t3dbc_mod::t3dbc_tile(), and mod_param::vadvection.

Referenced by pre_step3d().

Here is the call graph for this function:
Here is the caller graph for this function: