ROMS
Loading...
Searching...
No Matches
rp_rhs3d.F
Go to the documentation of this file.
1#include "cppdefs.h"
3#if defined TL_IOMS && defined SOLVE3D
4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group Andrew M. Moore !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This subroutine evaluates representers tangent linear right-hand- !
13! side terms for 3D momentum and tracers equations !
14! !
15! BASIC STATE variables needed: Hz, Huon, HVom, u, v, W, uclm, vclm, !
16! sustr, svstr, bustr, bvstr !
17! !
18!=======================================================================
19!
20 implicit none
21!
22 PRIVATE
23 PUBLIC :: rp_rhs3d
24!
25 CONTAINS
26!
27!***********************************************************************
28 SUBROUTINE rp_rhs3d (ng, tile)
29!***********************************************************************
30!
31 USE mod_param
32 USE mod_coupling
33# ifdef DIAGNOSTICS_UV
34 USE mod_diags
35# endif
36 USE mod_forces
37 USE mod_grid
38# ifdef WEC_MELLOR
39 USE mod_mixing
40# endif
41 USE mod_ocean
42 USE mod_stepping
43!
45 USE rp_prsgrd_mod, ONLY : rp_prsgrd
46# ifndef TS_FIXED
47# ifdef TS_DIF2
48 USE rp_t3dmix2_mod, ONLY : rp_t3dmix2
49# endif
50# ifdef TS_DIF4
51 USE rp_t3dmix4_mod, ONLY : rp_t3dmix4
52# endif
53# endif
54# ifdef RPM_RELAXATION
56# endif
57# ifdef UV_VIS2
59# endif
60# ifdef UV_VIS4
62# endif
63# ifdef RPM_RELAXATION
65# endif
66!
67! Imported variable declarations.
68!
69 integer, intent(in) :: ng, tile
70!
71! Local variable declarations.
72!
73 character (len=*), parameter :: myfile = &
74 & __FILE__
75!
76# include "tile.h"
77!
78!-----------------------------------------------------------------------
79! Initialize computations for new time step of the 3D primitive
80! variables.
81!-----------------------------------------------------------------------
82!
83 CALL rp_pre_step3d (ng, tile)
84!
85!-----------------------------------------------------------------------
86! Compute baroclinic pressure gradient.
87!-----------------------------------------------------------------------
88!
89 CALL rp_prsgrd (ng, tile)
90# ifndef TS_FIXED
91# ifdef TS_DIF2
92!
93!-----------------------------------------------------------------------
94! Compute horizontal harmonic mixing of tracer type variables.
95!-----------------------------------------------------------------------
96!
97 CALL rp_t3dmix2 (ng, tile)
98# endif
99# ifdef TS_DIF4
100!
101!-----------------------------------------------------------------------
102! Compute horizontal biharmonic mixing of tracer type variables.
103!-----------------------------------------------------------------------
104!
105 CALL rp_t3dmix4 (ng, tile)
106# endif
107# endif
108# ifdef RPM_RELAXATION
109!
110!-----------------------------------------------------------------------
111! Improve stability and convergence of the tangent linear representer
112! model tracer type variables by a "diffusive relaxation" to previous
113! Picard iteration solution.
114!-----------------------------------------------------------------------
115!
116 CALL rp_t3drelax (ng, tile)
117# endif
118!
119!-----------------------------------------------------------------------
120! Compute right-hand-side terms for the 3D momentum equations.
121!-----------------------------------------------------------------------
122!
123# ifdef PROFILE
124 CALL wclock_on (ng, irpm, 21, __line__, myfile)
125# endif
126 CALL rp_rhs3d_tile (ng, tile, &
127 & lbi, ubi, lbj, ubj, &
128 & imins, imaxs, jmins, jmaxs, &
129 & nrhs(ng), &
130 & grid(ng) % Hz, &
131 & grid(ng) % tl_Hz, &
132 & grid(ng) % Huon, &
133 & grid(ng) % tl_Huon, &
134 & grid(ng) % Hvom, &
135 & grid(ng) % tl_Hvom, &
136# if defined CURVGRID && defined UV_ADV
137 & grid(ng) % dmde, &
138 & grid(ng) % dndx, &
139# endif
140 & grid(ng) % fomn, &
141 & grid(ng) % om_u, &
142 & grid(ng) % om_v, &
143 & grid(ng) % on_u, &
144 & grid(ng) % on_v, &
145 & grid(ng) % pm, &
146 & grid(ng) % pn, &
147# ifdef WET_DRY_NOT_YET
148 & grid(ng)%umask_wet, &
149 & grid(ng)%vmask_wet, &
150# endif
151 & forces(ng) % bustr, &
152 & forces(ng) % tl_bustr, &
153 & forces(ng) % bvstr, &
154 & forces(ng) % tl_bvstr, &
155 & forces(ng) % sustr, &
156 & forces(ng) % tl_sustr, &
157 & forces(ng) % svstr, &
158 & forces(ng) % tl_svstr, &
159 & ocean(ng) % u, &
160 & ocean(ng) % tl_u, &
161 & ocean(ng) % v, &
162 & ocean(ng) % tl_v, &
163 & ocean(ng) % W, &
164 & ocean(ng) % tl_W, &
165# ifdef WEC_MELLOR
166 & ocean(ng) % u_stokes, &
167 & ocean(ng) % tl_u_stokes, &
168 & ocean(ng) % v_stokes, &
169 & ocean(ng) % tl_v_stokes, &
170 & ocean(ng) % tl_rulag3d, &
171 & ocean(ng) % tl_rvlag3d, &
172 & mixing(ng) % tl_rustr3d, &
173 & mixing(ng) % tl_rvstr3d, &
174# endif
175 & coupling(ng) % tl_rufrc, &
176 & coupling(ng) % tl_rvfrc, &
177# ifdef DIAGNOSTICS_UV
178!! & DIAGS(ng) % DiaRUfrc, &
179!! & DIAGS(ng) % DiaRVfrc, &
180!! & DIAGS(ng) % DiaRU, &
181!! & DIAGS(ng) % DiaRV, &
182# endif
183 & ocean(ng) % tl_ru, &
184 & ocean(ng) % tl_rv)
185# ifdef PROFILE
186 CALL wclock_off (ng, irpm, 21, __line__, myfile)
187# endif
188# ifdef UV_VIS2
189!
190!-----------------------------------------------------------------------
191! Compute horizontal, harmonic mixing of momentum.
192!-----------------------------------------------------------------------
193!
194 CALL rp_uv3dmix2 (ng, tile)
195# endif
196# ifdef UV_VIS4
197!
198!-----------------------------------------------------------------------
199! Compute horizontal, biharmonic mixing of momentum.
200!-----------------------------------------------------------------------
201!
202 CALL rp_uv3dmix4 (ng, tile)
203# endif
204# ifdef RPM_RELAXATION
205!
206!-----------------------------------------------------------------------
207! Improve stability and convergence of the tangent linear representer
208! model 3D momentum by a "diffusive relaxation" to previous Picard
209! iteration solution.
210!-----------------------------------------------------------------------
211!
212 CALL rp_uv3drelax (ng, tile)
213# endif
214
215 RETURN
216 END SUBROUTINE rp_rhs3d
217!
218!***********************************************************************
219 SUBROUTINE rp_rhs3d_tile (ng, tile, &
220 & LBi, UBi, LBj, UBj, &
221 & IminS, ImaxS, JminS, JmaxS, &
222 & nrhs, &
223 & Hz, tl_Hz, &
224 & Huon, tl_Huon, &
225 & Hvom, tl_Hvom, &
226# if defined CURVGRID && defined UV_ADV
227 & dmde, dndx, &
228# endif
229 & fomn, &
230 & om_u, om_v, on_u, on_v, pm, pn, &
231# ifdef WET_DRY_NOT_YET
232 & umask_wet, vmask_wet, &
233# endif
234 & bustr, tl_bustr, &
235 & bvstr, tl_bvstr, &
236 & sustr, tl_sustr, &
237 & svstr, tl_svstr, &
238 & u, tl_u, &
239 & v, tl_v, &
240 & W, tl_W, &
241# ifdef WEC_MELLOR
242 & u_stokes, tl_u_stokes, &
243 & v_stokes, tl_v_stokes, &
244 & tl_rulag3d, tl_rvlag3d, &
245 & tl_rustr3d, tl_rvstr3d, &
246# endif
247 & tl_rufrc, &
248 & tl_rvfrc, &
249# ifdef DIAGNOSTICS_UV
250!! & DiaRUfrc, DiaRVfrc, &
251!! & DiaRU, DiaRV, &
252# endif
253 & tl_ru, tl_rv)
254!***********************************************************************
255!
256 USE mod_param
257 USE mod_clima
258 USE mod_scalars
259!
260! Imported variable declarations.
261!
262 integer, intent(in) :: ng, tile
263 integer, intent(in) :: LBi, UBi, LBj, UBj
264 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
265 integer, intent(in) :: nrhs
266!
267# ifdef ASSUMED_SHAPE
268 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
269 real(r8), intent(in) :: Huon(LBi:,LBj:,:)
270 real(r8), intent(in) :: Hvom(LBi:,LBj:,:)
271# if defined CURVGRID && defined UV_ADV
272 real(r8), intent(in) :: dmde(LBi:,LBj:)
273 real(r8), intent(in) :: dndx(LBi:,LBj:)
274# endif
275 real(r8), intent(in) :: fomn(LBi:,LBj:)
276 real(r8), intent(in) :: om_u(LBi:,LBj:)
277 real(r8), intent(in) :: om_v(LBi:,LBj:)
278 real(r8), intent(in) :: on_u(LBi:,LBj:)
279 real(r8), intent(in) :: on_v(LBi:,LBj:)
280 real(r8), intent(in) :: pm(LBi:,LBj:)
281 real(r8), intent(in) :: pn(LBi:,LBj:)
282# ifdef WET_DRY_NOT_YET
283 real(r8), intent(in) :: umask_wet(LBi:,LBj:)
284 real(r8), intent(in) :: vmask_wet(LBi:,LBj:)
285# endif
286 real(r8), intent(in) :: bustr(LBi:,LBj:)
287 real(r8), intent(in) :: bvstr(LBi:,LBj:)
288 real(r8), intent(in) :: sustr(LBi:,LBj:)
289 real(r8), intent(in) :: svstr(LBi:,LBj:)
290 real(r8), intent(in) :: u(LBi:,LBj:,:,:)
291 real(r8), intent(in) :: v(LBi:,LBj:,:,:)
292 real(r8), intent(in) :: W(LBi:,LBj:,0:)
293
294 real(r8), intent(in) :: tl_Hz(LBi:,LBj:,:)
295 real(r8), intent(in) :: tl_Huon(LBi:,LBj:,:)
296 real(r8), intent(in) :: tl_Hvom(LBi:,LBj:,:)
297 real(r8), intent(in) :: tl_bustr(LBi:,LBj:)
298 real(r8), intent(in) :: tl_bvstr(LBi:,LBj:)
299 real(r8), intent(in) :: tl_sustr(LBi:,LBj:)
300 real(r8), intent(in) :: tl_svstr(LBi:,LBj:)
301 real(r8), intent(in) :: tl_u(LBi:,LBj:,:,:)
302 real(r8), intent(in) :: tl_v(LBi:,LBj:,:,:)
303 real(r8), intent(in) :: tl_W(LBi:,LBj:,0:)
304# ifdef WEC_MELLOR
305 real(r8), intent(in) :: u_stokes(LBi:,LBj:,:)
306 real(r8), intent(in) :: v_stokes(LBi:,LBj:,:)
307 real(r8), intent(in) :: tl_u_stokes(LBi:,LBj:,:)
308 real(r8), intent(in) :: tl_v_stokes(LBi:,LBj:,:)
309 real(r8), intent(in) :: tl_rulag3d(LBi:,LBj:,:)
310 real(r8), intent(in) :: tl_rvlag3d(LBi:,LBj:,:)
311 real(r8), intent(in) :: tl_rustr3d(LBi:,LBj:,:)
312 real(r8), intent(in) :: tl_rvstr3d(LBi:,LBj:,:)
313# endif
314# ifdef DIAGNOSTICS_UV
315!! real(r8), intent(inout) :: DiaRUfrc(LBi:,LBj:,:,:)
316!! real(r8), intent(inout) :: DiaRVfrc(LBi:,LBj:,:,:)
317!! real(r8), intent(inout) :: DiaRU(LBi:,LBj:,:,:,:)
318!! real(r8), intent(inout) :: DiaRV(LBi:,LBj:,:,:,:)
319# endif
320 real(r8), intent(inout) :: tl_ru(LBi:,LBj:,0:,:)
321 real(r8), intent(inout) :: tl_rv(LBi:,LBj:,0:,:)
322
323 real(r8), intent(out) :: tl_rufrc(LBi:,LBj:)
324 real(r8), intent(out) :: tl_rvfrc(LBi:,LBj:)
325# else
326 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
327 real(r8), intent(in) :: Huon(LBi:UBi,LBj:UBj,N(ng))
328 real(r8), intent(in) :: Hvom(LBi:UBi,LBj:UBj,N(ng))
329# if defined CURVGRID && defined UV_ADV
330 real(r8), intent(in) :: dmde(LBi:UBi,LBj:UBj)
331 real(r8), intent(in) :: dndx(LBi:UBi,LBj:UBj)
332# endif
333 real(r8), intent(in) :: fomn(LBi:UBi,LBj:UBj)
334 real(r8), intent(in) :: om_u(LBi:UBi,LBj:UBj)
335 real(r8), intent(in) :: om_v(LBi:UBi,LBj:UBj)
336 real(r8), intent(in) :: on_u(LBi:UBi,LBj:UBj)
337 real(r8), intent(in) :: on_v(LBi:UBi,LBj:UBj)
338 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
339 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
340# ifdef WET_DRY_NOT_YET
341 real(r8), intent(in) :: umask_wet(LBi:UBi,LBj:UBj)
342 real(r8), intent(in) :: vmask_wet(LBi:UBi,LBj:UBj)
343# endif
344 real(r8), intent(in) :: bustr(LBi:UBi,LBj:UBj)
345 real(r8), intent(in) :: bvstr(LBi:UBi,LBj:UBj)
346 real(r8), intent(in) :: sustr(LBi:UBi,LBj:UBj)
347 real(r8), intent(in) :: svstr(LBi:UBi,LBj:UBj)
348 real(r8), intent(in) :: u(LBi:UBi,LBj:UBj,N(ng),2)
349 real(r8), intent(in) :: v(LBi:UBi,LBj:UBj,N(ng),2)
350 real(r8), intent(in) :: W(LBi:UBi,LBj:UBj,0:N(ng))
351
352 real(r8), intent(in) :: tl_Hz(LBi:UBi,LBj:UBj,N(ng))
353 real(r8), intent(in) :: tl_Huon(LBi:UBi,LBj:UBj,N(ng))
354 real(r8), intent(in) :: tl_Hvom(LBi:UBi,LBj:UBj,N(ng))
355 real(r8), intent(in) :: tl_bustr(LBi:UBi,LBj:UBj)
356 real(r8), intent(in) :: tl_bvstr(LBi:UBi,LBj:UBj)
357 real(r8), intent(in) :: tl_sustr(LBi:UBi,LBj:UBj)
358 real(r8), intent(in) :: tl_svstr(LBi:UBi,LBj:UBj)
359 real(r8), intent(in) :: tl_u(LBi:UBi,LBj:UBj,N(ng),2)
360 real(r8), intent(in) :: tl_v(LBi:UBi,LBj:UBj,N(ng),2)
361 real(r8), intent(in) :: tl_W(LBi:UBi,LBj:UBj,0:N(ng))
362# ifdef WEC_MELLOR
363 real(r8), intent(in) :: u_stokes(LBi:UBi,LBj:UBj,N(ng))
364 real(r8), intent(in) :: v_stokes(LBi:UBi,LBj:UBj,N(ng))
365 real(r8), intent(in) :: tl_u_stokes(LBi:UBi,LBj:UBj,N(ng))
366 real(r8), intent(in) :: tl_v_stokes(LBi:UBi,LBj:UBj,N(ng))
367 real(r8), intent(in) :: tl_rulag3d(LBi:UBi,LBj:UBj,N(ng))
368 real(r8), intent(in) :: tl_rvlag3d(LBi:UBi,LBj:UBj,N(ng))
369 real(r8), intent(in) :: tl_rustr3d(LBi:UBi,LBj:UBj,N(ng))
370 real(r8), intent(in) :: tl_rvstr3d(LBi:UBi,LBj:UBj,N(ng))
371# endif
372# ifdef DIAGNOSTICS_UV
373!! real(r8), intent(inout) :: DiaRUfrc(LBi:UBi,LBj:UBj,3,NDM2d-1)
374!! real(r8), intent(inout) :: DiaRVfrc(LBi:UBi,LBj:UBj,3,NDM2d-1)
375!! real(r8), intent(inout) :: DiaRU(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
376!! real(r8), intent(inout) :: DiaRV(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
377# endif
378 real(r8), intent(inout) :: tl_ru(LBi:UBi,LBj:UBj,0:N(ng),2)
379 real(r8), intent(inout) :: tl_rv(LBi:UBi,LBj:UBj,0:N(ng),2)
380
381 real(r8), intent(out) :: tl_rufrc(LBi:UBi,LBj:UBj)
382 real(r8), intent(out) :: tl_rvfrc(LBi:UBi,LBj:UBj)
383# endif
384!
385! Local variable declarations.
386!
387 integer :: i, j, k
388
389 real(r8), parameter :: Gadv = -0.25_r8
390
391 real(r8) :: cff, cff1, cff2, cff3, cff4
392 real(r8) :: tl_cff, tl_cff1, tl_cff2, tl_cff3, tl_cff4
393
394 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: CF
395 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: DC
396 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: FC
397
398 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: tl_CF
399 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: tl_DC
400 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: tl_FC
401
402 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Huee
403 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Huxx
404 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Hvee
405 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Hvxx
406 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: UFx
407 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: UFe
408 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Uwrk
409 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: VFx
410 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: VFe
411 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: Vwrk
412 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: uee
413 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: uxx
414 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: vee
415 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: vxx
416 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: wrk
417
418 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_Huee
419 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_Huxx
420 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_Hvee
421 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_Hvxx
422 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_UFx
423 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_UFe
424 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_Uwrk
425 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_VFx
426 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_VFe
427 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_Vwrk
428 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_uee
429 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_uxx
430 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_vee
431 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_vxx
432 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: tl_wrk
433
434# include "set_bounds.h"
435
436# ifdef BODYFORCE
437!
438!-----------------------------------------------------------------------
439! Apply surface stress as a bodyforce: determine the thickness (m)
440! of the surface layer; then add in surface stress as a bodyfoce.
441!-----------------------------------------------------------------------
442!
443# ifdef DIAGNOSTICS_UV
444!! DO k=1,N(ng)
445!! DO j=Jstr,Jend
446!! DO i=Istr,Iend
447!! DiaRU(i,j,k,nrhs,M3vvis)=0.0_r8
448!! DiaRV(i,j,k,nrhs,M3vvis)=0.0_r8
449!! END DO
450!! END DO
451!! END DO
452!! DO j=Jstr,Jend
453!! DO i=IstrU,Iend
454!! DiaRUfrc(i,j,3,M2sstr)=0.0_r8
455!! DiaRUfrc(i,j,3,M2bstr)=0.0_r8
456!! END DO
457!! END DO
458!! DO j=JstrV,Jend
459!! DO i=Istr,Iend
460!! DiaRVfrc(i,j,3,M2sstr)=0.0_r8
461!! DiaRVfrc(i,j,3,M2bstr)=0.0_r8
462!! END DO
463!! END DO
464# endif
465 DO j=jstrv-1,jend
466 DO i=istru-1,iend
467 wrk(i,j)=0.0_r8
468 tl_wrk(i,j)=0.0_r8
469 END DO
470 END DO
471 DO k=n(ng),levsfrc(ng),-1
472 DO j=jstrv-1,jend
473 DO i=istru-1,iend
474 wrk(i,j)=wrk(i,j)+hz(i,j,k)
475 tl_wrk(i,j)=tl_wrk(i,j)+tl_hz(i,j,k)
476 END DO
477 END DO
478 END DO
479 DO j=jstr,jend
480 DO i=istru,iend
481 cff=0.25_r8*(pm(i-1,j)+pm(i,j))* &
482 & (pn(i-1,j)+pn(i,j))
483 cff1=1.0_r8/(cff*(wrk(i-1,j)+wrk(i,j)))
484 tl_cff1=-cff1*cff1*cff*(tl_wrk(i-1,j)+tl_wrk(i,j))+ &
485# ifdef TL_IOMS
486 & 2.0_r8*cff1
487# endif
488 uwrk(i,j)=sustr(i,j)*cff1
489 tl_uwrk(i,j)=tl_sustr(i,j)*cff1+ &
490 & sustr(i,j)*tl_cff1- &
491# ifdef TL_IOMS
492 & uwrk(i,j)
493# endif
494 END DO
495 END DO
496 DO j=jstrv,jend
497 DO i=istr,iend
498 cff=0.25*(pm(i,j-1)+pm(i,j))* &
499 & (pn(i,j-1)+pn(i,j))
500 cff1=1.0_r8/(cff*(wrk(i,j-1)+wrk(i,j)))
501 tl_cff1=-cff1*cff1*cff*(tl_wrk(i,j-1)+tl_wrk(i,j))+ &
502# ifdef TL_IOMS
503 & 2.0_r8*cff1
504# endif
505 vwrk(i,j)=svstr(i,j)*cff1
506 tl_vwrk(i,j)=tl_svstr(i,j)*cff1+ &
507 & svstr(i,j)*tl_cff1- &
508# ifdef TL_IOMS
509 & vwrk(i,j)
510# endif
511 END DO
512 END DO
513 DO k=levsfrc(ng),n(ng)
514 DO j=jstr,jend
515 DO i=istru,iend
516 cff=uwrk(i,j)*(hz(i ,j,k)+ &
517 & hz(i-1,j,k))
518 tl_cff=tl_uwrk(i,j)*(hz(i ,j,k)+ &
519 & hz(i-1,j,k))+ &
520 & uwrk(i,j)*(tl_hz(i ,j,k)+ &
521 & tl_hz(i-1,j,k))- &
522# ifdef TL_IOMS
523 & cff
524# endif
525
526!^ ru(i,j,k,nrhs)=ru(i,j,k,nrhs)+cff
527!^
528 tl_ru(i,j,k,nrhs)=tl_ru(i,j,k,nrhs)+tl_cff
529# ifdef DIAGNOSTICS_UV
530!! DiaRU(i,j,k,nrhs,M3vvis)=DiaRU(i,j,k,nrhs,M3vvis)+cff
531!! DiaRUfrc(i,j,3,M2sstr)=DiaRUfrc(i,j,3,M2sstr)+cff
532# endif
533 END DO
534 END DO
535 DO j=jstrv,jend
536 DO i=istr,iend
537 cff=vwrk(i,j)*(hz(i,j ,k)+ &
538 & hz(i,j-1,k))
539 tl_cff=tl_vwrk(i,j)*(hz(i,j ,k)+ &
540 & hz(i,j-1,k))+ &
541 & vwrk(i,j)*(tl_hz(i,j ,k)+ &
542 & tl_hz(i,j-1,k))- &
543# ifdef TL_IOMS
544 & cff
545# endif
546!^ rv(i,j,k,nrhs)=rv(i,j,k,nrhs)+cff
547!^
548 tl_rv(i,j,k,nrhs)=tl_rv(i,j,k,nrhs)+tl_cff
549# ifdef DIAGNOSTICS_UV
550!! DiaRV(i,j,k,nrhs,M3vvis)=DiaRV(i,j,k,nrhs,M3vvis)+cff
551!! DiaRVfrc(i,j,3,M2sstr)=DiaRVfrc(i,j,3,M2sstr)+cff
552# endif
553 END DO
554 END DO
555 END DO
556!
557! Apply bottom stress as a bodyforce: determine the thickness (m)
558! of the bottom layer; then add in bottom stress as a bodyfoce.
559!
560 DO j=jstrv-1,jend
561 DO i=istru-1,iend
562 wrk(i,j)=0.0_r8
563 tl_wrk(i,j)=0.0_r8
564 END DO
565 END DO
566 DO k=1,levbfrc(ng)
567 DO j=jstrv-1,jend
568 DO i=istru-1,iend
569 wrk(i,j)=wrk(i,j)+hz(i,j,k)
570 tl_wrk(i,j)=tl_wrk(i,j)+tl_hz(i,j,k)
571 END DO
572 END DO
573 END DO
574 DO j=jstr,jend
575 DO i=istru,iend
576 cff=0.25_r8*(pm(i-1,j)+pm(i,j))* &
577 & (pn(i-1,j)+pn(i,j))
578 cff1=1.0_r8/(cff*(wrk(i-1,j)+wrk(i,j)))
579 tl_cff1=-cff1*cff1*cff*(tl_wrk(i-1,j)+tl_wrk(i,j))+ &
580# ifdef TL_IOMS
581 & 2.0_r8*cff1
582# endif
583 uwrk(i,j)=bustr(i,j)*cff1
584 tl_uwrk(i,j)=tl_bustr(i,j)*cff1+ &
585 & bustr(i,j)*tl_cff1- &
586# ifdef TL_IOMS
587 & uwrk(i,j)
588# endif
589 END DO
590 END DO
591 DO j=jstrv,jend
592 DO i=istr,iend
593 cff=0.25_r8*(pm(i,j-1)+pm(i,j))* &
594 & (pn(i,j-1)+pn(i,j))
595 cff1=1.0_r8/(cff*(wrk(i,j-1)+wrk(i,j)))
596 tl_cff1=-cff1*cff1*cff*(tl_wrk(i,j-1)+tl_wrk(i,j))+ &
597# ifdef TL_IOMS
598 & 2.0_r8*cff1
599# endif
600 vwrk(i,j)=bvstr(i,j)*cff1
601 tl_vwrk(i,j)=tl_bvstr(i,j)*cff1+ &
602 & bvstr(i,j)*tl_cff1- &
603# ifdef TL_IOMS
604 & vwrk(i,j)
605# endif
606 END DO
607 END DO
608 DO k=1,levbfrc(ng)
609 DO j=jstr,jend
610 DO i=istru,iend
611 cff=uwrk(i,j)*(hz(i ,j,k)+ &
612 & hz(i-1,j,k))
613 tl_cff=tl_uwrk(i,j)*(hz(i ,j,k)+ &
614 & hz(i-1,j,k))+ &
615 & uwrk(i,j)*(tl_hz(i ,j,k)+ &
616 & tl_hz(i-1,j,k))- &
617# ifdef TL_IOMS
618 & cff
619# endif
620!^ ru(i,j,k,nrhs)=ru(i,j,k,nrhs)-cff
621!^
622 tl_ru(i,j,k,nrhs)=tl_ru(i,j,k,nrhs)-tl_cff
623# ifdef DIAGNOSTICS_UV
624!! DiaRU(i,j,k,nrhs,M3vvis)=DiaRU(i,j,k,nrhs,M3vvis)-cff
625!! DiaRUfrc(i,j,3,M2bstr)=DiaRUfrc(i,j,3,M2bstr)-cff
626# endif
627 END DO
628 END DO
629 DO j=jstrv,jend
630 DO i=istr,iend
631 cff=vwrk(i,j)*(hz(i,j ,k)+ &
632 & hz(i,j-1,k))
633 tl_cff=tl_vwrk(i,j)*(hz(i,j ,k)+ &
634 & hz(i,j-1,k))+ &
635 & vwrk(i,j)*(tl_hz(i,j ,k)+ &
636 & tl_hz(i,j-1,k))- &
637# ifdef TL_IOMS
638 & cff
639# endif
640!^ rv(i,j,k,nrhs)=rv(i,j,k,nrhs)-cff
641!^
642 tl_rv(i,j,k,nrhs)=tl_rv(i,j,k,nrhs)-tl_cff
643# ifdef DIAGNOSTICS_UV
644!! DiaRV(i,j,k,nrhs,M3vvis)=DiaRV(i,j,k,nrhs,M3vvis)-cff
645!! DiaRVfrc(i,j,3,M2bstr)=DiaRVfrc(i,j,3,M2bstr)-cff
646# endif
647 END DO
648 END DO
649 END DO
650# endif
651!
652 k_loop : DO k=1,n(ng)
653
654# ifdef UV_COR
655!
656!-----------------------------------------------------------------------
657! Add in Coriolis terms.
658!-----------------------------------------------------------------------
659!
660 DO j=jstrv-1,jend
661 DO i=istru-1,iend
662 cff=0.5_r8*hz(i,j,k)*fomn(i,j)
663 tl_cff=0.5_r8*tl_hz(i,j,k)*fomn(i,j)
664 ufx(i,j)=cff*(v(i,j ,k,nrhs)+ &
665# ifdef WEC_MELLOR
666 & v_stokes(i,j ,k)+ &
667 & v_stokes(i,j+1,k)+ &
668# endif
669 & v(i,j+1,k,nrhs))
670 tl_ufx(i,j)=tl_cff*(v(i,j ,k,nrhs)+ &
671# ifdef WEC_MELLOR
672 & v_stokes(i,j ,k)+ &
673 & v_stokes(i,j+1,k)+ &
674# endif
675 & v(i,j+1,k,nrhs))+ &
676 & cff*(tl_v(i,j ,k,nrhs)+ &
677# ifdef WEC_MELLOR
678 & tl_v_stokes(i,j ,k)+ &
679 & tl_v_stokes(i,j+1,k)+ &
680# endif
681 & tl_v(i,j+1,k,nrhs))- &
682# ifdef TL_IOMS
683 & ufx(i,j)
684# endif
685 vfe(i,j)=cff*(u(i ,j,k,nrhs)+ &
686# ifdef WEC_MELLOR
687 & u_stokes(i ,j,k)+ &
688 & u_stokes(i+1,j,k)+ &
689# endif
690 & u(i+1,j,k,nrhs))
691 tl_vfe(i,j)=tl_cff*(u(i ,j,k,nrhs)+ &
692# ifdef WEC_MELLOR
693 & u_stokes(i ,j,k)+ &
694 & u_stokes(i+1,j,k)+ &
695# endif
696 & u(i+1,j,k,nrhs))+ &
697 & cff*(tl_u(i ,j,k,nrhs)+ &
698# ifdef WEC_MELLOR
699 & tl_u_stokes(i ,j,k)+ &
700 & tl_u_stokes(i+1,j,k)+ &
701# endif
702 & tl_u(i+1,j,k,nrhs))- &
703# ifdef TL_IOMS
704 & vfe(i,j)
705# endif
706 END DO
707 END DO
708 DO j=jstr,jend
709 DO i=istru,iend
710!^ cff1=0.5_r8*(UFx(i,j)+UFx(i-1,j))
711!^
712 tl_cff1=0.5_r8*(tl_ufx(i,j)+tl_ufx(i-1,j))
713!^ ru(i,j,k,nrhs)=ru(i,j,k,nrhs)+cff1
714!^
715 tl_ru(i,j,k,nrhs)=tl_ru(i,j,k,nrhs)+tl_cff1
716# ifdef DIAGNOSTICS_UV
717!! DiaRU(i,j,k,nrhs,M3fcor)=cff1
718# endif
719 END DO
720 END DO
721 DO j=jstrv,jend
722 DO i=istr,iend
723!^ cff1=0.5_r8*(VFe(i,j)+VFe(i,j-1))
724!^
725 tl_cff1=0.5_r8*(tl_vfe(i,j)+tl_vfe(i,j-1))
726!^ rv(i,j,k,nrhs)=rv(i,j,k,nrhs)-cff1
727!^
728 tl_rv(i,j,k,nrhs)=tl_rv(i,j,k,nrhs)-tl_cff1
729# ifdef DIAGNOSTICS_UV
730!! DiaRV(i,j,k,nrhs,M3fcor)=-cff1
731# endif
732 END DO
733 END DO
734# endif
735# if defined CURVGRID && defined UV_ADV
736!
737!-----------------------------------------------------------------------
738! Add in curvilinear transformation terms.
739!-----------------------------------------------------------------------
740!
741 DO j=jstrv-1,jend
742 DO i=istru-1,iend
743 cff1=0.5_r8*(v(i,j ,k,nrhs)+ &
744# ifdef WEC_MELLOR
745 & v_stokes(i,j ,k)+ &
746 & v_stokes(i,j+1,k)+ &
747# endif
748 & v(i,j+1,k,nrhs))
749 tl_cff1=0.5_r8*(tl_v(i,j ,k,nrhs)+ &
750# ifdef WEC_MELLOR
751 & tl_v_stokes(i,j ,k)+ &
752 & tl_v_stokes(i,j+1,k)+ &
753# endif
754 & tl_v(i,j+1,k,nrhs))
755 cff2=0.5_r8*(u(i ,j,k,nrhs)+ &
756# ifdef WEC_MELLOR
757 & u_stokes(i ,j,k)+ &
758 & u_stokes(i+1,j,k)+ &
759# endif
760 & u(i+1,j,k,nrhs))
761 tl_cff2=0.5_r8*(tl_u(i ,j,k,nrhs)+ &
762# ifdef WEC_MELLOR
763 & tl_u_stokes(i ,j,k)+ &
764 & tl_u_stokes(i+1,j,k)+ &
765# endif
766 & tl_u(i+1,j,k,nrhs))
767 cff3=cff1*dndx(i,j)
768 tl_cff3=tl_cff1*dndx(i,j)
769 cff4=cff2*dmde(i,j)
770 tl_cff4=tl_cff2*dmde(i,j)
771 cff=hz(i,j,k)*(cff3-cff4)
772 tl_cff=tl_hz(i,j,k)*(cff3-cff4)+ &
773 & hz(i,j,k)*(tl_cff3-tl_cff4)- &
774# ifdef TL_IOMS
775 & cff
776# endif
777!^ UFx(i,j)=cff*cff1
778!^
779 tl_ufx(i,j)=tl_cff*cff1+cff*tl_cff1- &
780# ifdef TL_IOMS
781 & cff*cff1
782# endif
783!^ VFe(i,j)=cff*cff2
784!^
785 tl_vfe(i,j)=tl_cff*cff2+cff*tl_cff2- &
786# ifdef TL_IOMS
787 & cff*cff2
788# endif
789# if defined DIAGNOSTICS_UV
790!! cff=Hz(i,j,k)*cff4
791!! Uwrk(i,j)=-cff*cff1 ! u equation, ETA-term
792!! Vwrk(i,j)=-cff*cff2 ! v equation, ETA-term
793# endif
794 END DO
795 END DO
796 DO j=jstr,jend
797 DO i=istru,iend
798!^ cff1=0.5_r8*(UFx(i,j)+UFx(i-1,j))
799!^
800 tl_cff1=0.5_r8*(tl_ufx(i,j)+tl_ufx(i-1,j))
801!^ ru(i,j,k,nrhs)=ru(i,j,k,nrhs)+cff1
802!^
803 tl_ru(i,j,k,nrhs)=tl_ru(i,j,k,nrhs)+tl_cff1
804# ifdef DIAGNOSTICS_UV
805!! cff2=0.5_r8*(Uwrk(i,j)+Uwrk(i-1,j))
806!! DiaRU(i,j,k,nrhs,M3xadv)=cff1-cff2
807!! DiaRU(i,j,k,nrhs,M3yadv)=cff2
808!! DiaRU(i,j,k,nrhs,M3hadv)=cff1
809# endif
810 END DO
811 END DO
812 DO j=jstrv,jend
813 DO i=istr,iend
814!^ cff1=0.5_r8*(VFe(i,j)+VFe(i,j-1))
815!^
816 tl_cff1=0.5_r8*(tl_vfe(i,j)+tl_vfe(i,j-1))
817!^ rv(i,j,k,nrhs)=rv(i,j,k,nrhs)-cff1
818!^
819 tl_rv(i,j,k,nrhs)=tl_rv(i,j,k,nrhs)-tl_cff1
820# ifdef DIAGNOSTICS_UV
821!! cff2=0.5_r8*(Vwrk(i,j)+Vwrk(i,j-1))
822!! DiaRV(i,j,k,nrhs,M3xadv)=-cff1+cff2
823!! DiaRV(i,j,k,nrhs,M3yadv)=-cff2
824!! DiaRV(i,j,k,nrhs,M3hadv)=-cff1
825# endif
826 END DO
827 END DO
828# endif
829!
830!-----------------------------------------------------------------------
831! Add in nudging of 3D momentum climatology.
832!-----------------------------------------------------------------------
833!
834 IF (lnudgem3clm(ng)) THEN
835 DO j=jstr,jend
836 DO i=istru,iend
837 cff=0.25_r8*(clima(ng)%M3nudgcof(i-1,j,k)+ &
838 & clima(ng)%M3nudgcof(i ,j,k))* &
839 & om_u(i,j)*on_u(i,j)
840!^ ru(i,j,k,nrhs)=ru(i,j,k,nrhs)+ &
841!^ & cff*(Hz(i-1,j,k)+Hz(i,j,k))* &
842!^ & (CLIMA(ng)%uclm(i,j,k)- &
843!^ & u(i,j,k,nrhs))
844!^
845 tl_ru(i,j,k,nrhs)=tl_ru(i,j,k,nrhs)+ &
846 & cff*((hz(i-1,j,k)+hz(i,j,k))* &
847 & (-tl_u(i,j,k,nrhs))+ &
848 & (tl_hz(i-1,j,k)+tl_hz(i,j,k))* &
849 & (clima(ng)%uclm(i,j,k)- &
850 & u(i,j,k,nrhs)))+ &
851# ifdef TL_IOMS
852 & cff*(hz(i-1,j,k)+hz(i,j,k))* &
853 & u(i,j,k,nrhs)
854# endif
855 END DO
856 END DO
857 DO j=jstrv,jend
858 DO i=istr,iend
859 cff=0.25_r8*(clima(ng)%M3nudgcof(i,j-1,k)+ &
860 & clima(ng)%M3nudgcof(i,j ,k))* &
861 & om_v(i,j)*on_v(i,j)
862!^ rv(i,j,k,nrhs)=rv(i,j,k,nrhs)+ &
863!^ & cff*(Hz(i,j-1,k)+Hz(i,j,k))* &
864!^ & (CLIMA(ng)%vclm(i,j,k)- &
865!^ & v(i,j,k,nrhs))
866!^
867 tl_rv(i,j,k,nrhs)=tl_rv(i,j,k,nrhs)+ &
868 & cff*((hz(i,j-1,k)+hz(i,j,k))* &
869 & (-tl_v(i,j,k,nrhs))+ &
870 & (tl_hz(i,j-1,k)+tl_hz(i,j,k))* &
871 & (clima(ng)%vclm(i,j,k)- &
872 & v(i,j,k,nrhs)))+ &
873# ifdef TL_IOMS
874 & cff*(hz(i,j-1,k)+hz(i,j,k))* &
875 & v(i,j,k,nrhs)
876# endif
877 END DO
878 END DO
879 END IF
880
881# ifdef UV_ADV
882!
883!-----------------------------------------------------------------------
884! Add in horizontal advection of momentum.
885!-----------------------------------------------------------------------
886!
887! Compute diagonal [UFx,VFe] and off-diagonal [UFe,VFx] components
888! of tensor of momentum flux due to horizontal advection.
889!
890# ifdef UV_C2ADVECTION
891!
892! Second-order, centered differences advection.
893!
894 DO j=jstr,jend
895 DO i=istru-1,iend
896 ufx(i,j)=0.25_r8*(u(i ,j,k,nrhs)+ &
897# ifdef WEC_MELLOR
898 & u_stokes(i ,j,k)+ &
899 & u_stokes(i+1,j,k)+ &
900# endif
901 & u(i+1,j,k,nrhs))* &
902 & (huon(i ,j,k)+ &
903 & huon(i+1,j,k))
904 tl_ufx(i,j)=0.25_r8* &
905 & ((tl_u(i ,j,k,nrhs)+ &
906# ifdef WEC_MELLOR
907 & tl_u_stokes(i ,j,k)+ &
908 & tl_u_stokes(i+1,j,k)+ &
909# endif
910 & tl_u(i+1,j,k,nrhs))* &
911 & (huon(i ,j,k)+ &
912 & huon(i+1,j,k))+ &
913 & (u(i ,j,k,nrhs)+ &
914# ifdef WEC_MELLOR
915 & u_stokes(i ,j,k)+ &
916 & u_stokes(i+1,j,k)+ &
917# endif
918 & u(i+1,j,k,nrhs))* &
919 & (tl_huon(i ,j,k)+ &
920 & tl_huon(i+1,j,k)))- &
921# ifdef TL_IOMS
922 & ufx(i,j)
923# endif
924 END DO
925 END DO
926 DO j=jstr,jend+1
927 DO i=istru,iend
928 ufe(i,j)=0.25_r8*(u(i,j-1,k,nrhs)+ &
929# ifdef WEC_MELLOR
930 & u_stokes(i,j-1,k)+ &
931 & u_stokes(i,j ,k)+ &
932# endif
933 & u(i,j ,k,nrhs))* &
934 & (hvom(i-1,j,k)+ &
935 & hvom(i ,j,k))
936 tl_ufe(i,j)=0.25_r8* &
937 & ((tl_u(i,j-1,k,nrhs)+ &
938# ifdef WEC_MELLOR
939 & tl_u_stokes(i,j-1,k)+ &
940 & tl_u_stokes(i,j ,k)+ &
941# endif
942 & tl_u(i,j ,k,nrhs))* &
943 & (hvom(i-1,j,k)+ &
944 & hvom(i ,j,k))+ &
945 & (u(i,j-1,k,nrhs)+ &
946# ifdef WEC_MELLOR
947 & u_stokes(i,j-1,k)+ &
948 & u_stokes(i,j ,k)+ &
949# endif
950 & u(i,j ,k,nrhs))* &
951 & (tl_hvom(i-1,j,k)+ &
952 & tl_hvom(i ,j,k)))- &
953# ifdef TL_IOMS
954 & ufe(i,j)
955# endif
956 END DO
957 END DO
958 DO j=jstrv,jend
959 DO i=istr,iend+1
960 vfx(i,j)=0.25_r8*(v(i-1,j,k,nrhs)+ &
961# ifdef WEC_MELLOR
962 & v_stokes(i-1,j,k)+ &
963 & v_stokes(i ,j,k)+ &
964# endif
965 & v(i ,j,k,nrhs))* &
966 & (huon(i,j-1,k)+ &
967 huon(i,j ,k))
968 tl_vfx(i,j)=0.25_r8* &
969 & ((tl_v(i-1,j,k,nrhs)+ &
970# ifdef WEC_MELLOR
971 & tl_v_stokes(i-1,j,k)+ &
972 & tl_v_stokes(i ,j,k)+ &
973# endif
974 & tl_v(i ,j,k,nrhs))* &
975 & (huon(i,j-1,k)+ &
976 & huon(i,j ,k))+ &
977 & (v(i-1,j,k,nrhs)+ &
978# ifdef WEC_MELLOR
979 & v_stokes(i-1,j,k)+ &
980 & v_stokes(i ,j,k)+ &
981# endif
982 & v(i ,j,k,nrhs))* &
983 & (tl_huon(i,j-1,k)+ &
984 & tl_huon(i,j ,k)))- &
985# ifdef TL_IOMS
986 & vfx(i,j)
987# endif
988 END DO
989 END DO
990 DO j=jstrv-1,jend
991 DO i=istr,iend
992 vfe(i,j)=0.25_r8*(v(i,j ,k,nrhs)+ &
993# ifdef WEC_MELLOR
994 & v_stokes(i,j ,k)+ &
995 & v_stokes(i,j+1,k)+ &
996# endif
997 & v(i,j+1,k,nrhs))* &
998 & (hvom(i,j ,k)+ &
999 & hvom(i,j+1,k))
1000 tl_vfe(i,j)=0.25_r8* &
1001 & ((tl_v(i,j ,k,nrhs)+ &
1002# ifdef WEC_MELLOR
1003 & tl_v_stokes(i,j ,k)+ &
1004 & tl_v_stokes(i,j+1,k)+ &
1005# endif
1006 & tl_v(i,j+1,k,nrhs))* &
1007 & (hvom(i,j ,k)+ &
1008 & hvom(i,j+1,k))+ &
1009 & (v(i,j ,k,nrhs)+ &
1010# ifdef WEC_MELLOR
1011 & v_stokes(i,j ,k)+ &
1012 & v_stokes(i,j+1,k)+ &
1013# endif
1014 & v(i,j+1,k,nrhs))* &
1015 & (tl_hvom(i,j ,k)+ &
1016 & tl_hvom(i,j+1,k)))- &
1017# ifdef TL_IOMS
1018 & vfe(i,j)
1019# endif
1020 END DO
1021 END DO
1022# else
1023 DO j=jstr,jend
1024 DO i=istrum1,iendp1
1025 uxx(i,j)=u(i-1,j,k,nrhs)-2.0_r8*u(i,j,k,nrhs)+ &
1026# ifdef WEC_MELLOR
1027 & u_stokes(i-1,j,k)-2.0_r8*u_stokes(i,j,k)+ &
1028 & u_stokes(i+1,j,k)+ &
1029# endif
1030 & u(i+1,j,k,nrhs)
1031 tl_uxx(i,j)=tl_u(i-1,j,k,nrhs)-2.0_r8*tl_u(i,j,k,nrhs)+ &
1032# ifdef WEC_MELLOR
1033 & tl_u_stokes(i-1,j,k)-2.0_r8*tl_u_stokes(i,j,k)+ &
1034 & tl_u_stokes(i+1,j,k)+ &
1035# endif
1036 & tl_u(i+1,j,k,nrhs)
1037 huxx(i,j)=huon(i-1,j,k)-2.0_r8*huon(i,j,k)+huon(i+1,j,k)
1038 tl_huxx(i,j)=tl_huon(i-1,j,k)-2.0_r8*tl_huon(i,j,k)+ &
1039 & tl_huon(i+1,j,k)
1040 END DO
1041 END DO
1042 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
1043 IF (domain(ng)%Western_Edge(tile)) THEN
1044 DO j=jstr,jend
1045 uxx(istr,j)=uxx(istr+1,j)
1046 tl_uxx(istr,j)=tl_uxx(istr+1,j)
1047 huxx(istr,j)=huxx(istr+1,j)
1048 tl_huxx(istr,j)=tl_huxx(istr+1,j)
1049 END DO
1050 END IF
1051 END IF
1052 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
1053 IF (domain(ng)%Eastern_Edge(tile)) THEN
1054 DO j=jstr,jend
1055 uxx(iend+1,j)=uxx(iend,j)
1056 tl_uxx(iend+1,j)=tl_uxx(iend,j)
1057 huxx(iend+1,j)=huxx(iend,j)
1058 tl_huxx(iend+1,j)=tl_huxx(iend,j)
1059 END DO
1060 END IF
1061 END IF
1062# ifdef UV_C4ADVECTION
1063!
1064! Fourth-order, centered differences u-momentum horizontal advection.
1065!
1066 cff=1.0_r8/6.0_r8
1067 DO j=jstr,jend
1068 DO i=istru-1,iend
1069 ufx(i,j)=0.25_r8*(u(i ,j,k,nrhs)+ &
1070# ifdef WEC_MELLOR
1071 & u_stokes(i ,j,k)+ &
1072 & u_stokes(i+1,j,k)+ &
1073# endif
1074 & u(i+1,j,k,nrhs)- &
1075 & cff*(uxx(i ,j)+ &
1076 & uxx(i+1,j)))* &
1077 & (huon(i ,j,k)+ &
1078 & huon(i+1,j,k)- &
1079 & cff*(huxx(i ,j)+ &
1080 & huxx(i+1,j)))
1081 tl_ufx(i,j)=0.25_r8*((tl_u(i ,j,k,nrhs)+ &
1082# ifdef WEC_MELLOR
1083 & tl_u_stokes(i ,j,k)+ &
1084 & tl_u_stokes(i+1,j,k)+ &
1085# endif
1086 & tl_u(i+1,j,k,nrhs)- &
1087 & cff*(tl_uxx(i ,j)+ &
1088 & tl_uxx(i+1,j)))* &
1089 & (huon(i ,j,k)+ &
1090 & huon(i+1,j,k)- &
1091 & cff*(huxx(i ,j)+ &
1092 & huxx(i+1,j)))+ &
1093 & (u(i ,j,k,nrhs)+ &
1094# ifdef WEC_MELLOR
1095 & u_stokes(i ,j,k)+ &
1096 & u_stokes(i+1,j,k)+ &
1097# endif
1098 & u(i+1,j,k,nrhs)- &
1099 & cff*(uxx(i ,j)+ &
1100 & uxx(i+1,j)))* &
1101 & (tl_huon(i ,j,k)+ &
1102 & tl_huon(i+1,j,k)- &
1103 & cff*(tl_huxx(i ,j)+ &
1104 & tl_huxx(i+1,j))))- &
1105# ifdef TL_IOMS
1106 & ufx(i,j)
1107# endif
1108 END DO
1109 END DO
1110# else
1111!
1112! Third-order, upstream bias u-momentum advection with velocity
1113! dependent hyperdiffusion.
1114!
1115 DO j=jstr,jend
1116 DO i=istru-1,iend
1117 cff1=u(i ,j,k,nrhs)+ &
1118# ifdef WEC_MELLOR
1119 & u_stokes(i ,j,k)+ &
1120 & u_stokes(i+1,j,k)+ &
1121# endif
1122 & u(i+1,j,k,nrhs)
1123 tl_cff1=tl_u(i ,j,k,nrhs)+ &
1124# ifdef WEC_MELLOR
1125 & tl_u_stokes(i ,j,k)+ &
1126 & tl_u_stokes(i+1,j,k)+ &
1127# endif
1128 & tl_u(i+1,j,k,nrhs)
1129 IF (cff1.gt.0.0_r8) THEN
1130 cff=uxx(i,j)
1131 tl_cff=tl_uxx(i,j)
1132 ELSE
1133 cff=uxx(i+1,j)
1134 tl_cff=tl_uxx(i+1,j)
1135 END IF
1136 ufx(i,j)=0.25_r8*(cff1+gadv*cff)* &
1137 & (huon(i ,j,k)+ &
1138 & huon(i+1,j,k)+ &
1139 & gadv*0.5_r8*(huxx(i ,j)+ &
1140 & huxx(i+1,j)))
1141 tl_ufx(i,j)=0.25_r8* &
1142 & ((tl_cff1+gadv*tl_cff)* &
1143 & (huon(i ,j,k)+ &
1144 & huon(i+1,j,k)+ &
1145 & gadv*0.5_r8*(huxx(i ,j)+ &
1146 & huxx(i+1,j)))+ &
1147 & (cff1+gadv*cff)* &
1148 & (tl_huon(i ,j,k)+ &
1149 & tl_huon(i+1,j,k)+ &
1150 & gadv*0.5_r8*(tl_huxx(i ,j)+ &
1151 & tl_huxx(i+1,j))))- &
1152# ifdef TL_IOMS
1153 & ufx(i,j)
1154# endif
1155 END DO
1156 END DO
1157# endif
1158 DO j=jstrm1,jendp1
1159 DO i=istru,iend
1160 uee(i,j)=u(i,j-1,k,nrhs)-2.0_r8*u(i,j,k,nrhs)+ &
1161# ifdef WEC_MELLOR
1162 & u_stokes(i,j-1,k)-2.0_r8*u_stokes(i,j,k)+ &
1163 & u_stokes(i,j+1,k)+ &
1164# endif
1165 & u(i,j+1,k,nrhs)
1166 tl_uee(i,j)=tl_u(i,j-1,k,nrhs)-2.0_r8*tl_u(i,j,k,nrhs)+ &
1167# ifdef WEC_MELLOR
1168 & tl_u_stokes(i,j-1,k)-2.0_r8*tl_u_stokes(i,j,k)+ &
1169 & tl_u_stokes(i,j+1,k)+ &
1170# endif
1171 & tl_u(i,j+1,k,nrhs)
1172 END DO
1173 END DO
1174 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
1175 IF (domain(ng)%Southern_Edge(tile)) THEN
1176 DO i=istru,iend
1177 uee(i,jstr-1)=uee(i,jstr)
1178 tl_uee(i,jstr-1)=tl_uee(i,jstr)
1179 END DO
1180 END IF
1181 END IF
1182 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
1183 IF (domain(ng)%Northern_Edge(tile)) THEN
1184 DO i=istru,iend
1185 uee(i,jend+1)=uee(i,jend)
1186 tl_uee(i,jend+1)=tl_uee(i,jend)
1187 END DO
1188 END IF
1189 END IF
1190 DO j=jstr,jend+1
1191 DO i=istru-1,iend
1192 hvxx(i,j)=hvom(i-1,j,k)-2.0_r8*hvom(i,j,k)+hvom(i+1,j,k)
1193 tl_hvxx(i,j)=tl_hvom(i-1,j,k)-2.0_r8*tl_hvom(i,j,k)+ &
1194 & tl_hvom(i+1,j,k)
1195 END DO
1196 END DO
1197# ifdef UV_C4ADVECTION
1198 cff=1.0_r8/6.0_r8
1199 DO j=jstr,jend+1
1200 DO i=istru,iend
1201 ufe(i,j)=0.25_r8*(u(i,j ,k,nrhs)+ &
1202# ifdef WEC_MELLOR
1203 & u_stokes(i,j ,k)+ &
1204 & u_stokes(i,j-1,k)+ &
1205# endif
1206 & u(i,j-1,k,nrhs)- &
1207 & cff*(uee(i,j )+ &
1208 & uee(i,j-1)))* &
1209 & (hvom(i ,j,k)+ &
1210 & hvom(i-1,j,k)- &
1211 & cff*(hvxx(i ,j)+ &
1212 & hvxx(i-1,j)))
1213 tl_ufe(i,j)=0.25_r8*((tl_u(i,j ,k,nrhs)+ &
1214# ifdef WEC_MELLOR
1215 & tl_u_stokes(i,j ,k)+ &
1216 & tl_u_stokes(i,j-1,k)+ &
1217# endif
1218 & tl_u(i,j-1,k,nrhs)- &
1219 & cff*(tl_uee(i,j )+ &
1220 & tl_uee(i,j-1)))* &
1221 & (hvom(i ,j,k)+ &
1222 & hvom(i-1,j,k)- &
1223 & cff*(hvxx(i ,j)+ &
1224 & hvxx(i-1,j)))+ &
1225 & (u(i,j ,k,nrhs)+ &
1226# ifdef WEC_MELLOR
1227 & u_stokes(i,j ,k)+ &
1228 & u_stokes(i,j-1,k)+ &
1229# endif
1230 & u(i,j-1,k,nrhs)- &
1231 & cff*(uee(i,j )+ &
1232 & uee(i,j-1)))* &
1233 & (tl_hvom(i ,j,k)+ &
1234 & tl_hvom(i-1,j,k)- &
1235 & cff*(tl_hvxx(i ,j)+ &
1236 & tl_hvxx(i-1,j))))- &
1237# ifdef TL_IOMS
1238 & ufe(i,j)
1239# endif
1240 END DO
1241 END DO
1242# else
1243 DO j=jstr,jend+1
1244 DO i=istru,iend
1245 cff1=u(i,j ,k,nrhs)+ &
1246# ifdef WEC_MELLOR
1247 & u_stokes(i,j ,k)+ &
1248 & u_stokes(i,j-1,k)+ &
1249# endif
1250 & u(i,j-1,k,nrhs)
1251 tl_cff1=tl_u(i,j,k,nrhs)+ &
1252# ifdef WEC_MELLOR
1253 & tl_u_stokes(i,j ,k)+ &
1254 & tl_u_stokes(i,j-1,k)+ &
1255# endif
1256 & tl_u(i,j-1,k,nrhs)
1257 cff2=hvom(i,j,k)+hvom(i-1,j,k)
1258 tl_cff2=tl_hvom(i,j,k)+tl_hvom(i-1,j,k)
1259 IF (cff2.gt.0.0_r8) THEN
1260 cff=uee(i,j-1)
1261 tl_cff=tl_uee(i,j-1)
1262 ELSE
1263 cff=uee(i,j)
1264 tl_cff=tl_uee(i,j)
1265 END IF
1266 ufe(i,j)=0.25_r8*(cff1+gadv*cff)* &
1267 & (cff2+gadv*0.5_r8*(hvxx(i ,j)+ &
1268 & hvxx(i-1,j)))
1269 tl_ufe(i,j)=0.25_r8* &
1270 & ((tl_cff1+gadv*tl_cff)* &
1271 & (cff2+gadv*0.5_r8*(hvxx(i ,j)+ &
1272 & hvxx(i-1,j)))+ &
1273 & (cff1+gadv*cff)* &
1274 & (tl_cff2+gadv*0.5_r8*(tl_hvxx(i ,j)+ &
1275 & tl_hvxx(i-1,j))))- &
1276# ifdef TL_IOMS
1277 & ufe(i,j)
1278# endif
1279 END DO
1280 END DO
1281# endif
1282 DO j=jstrv,jend
1283 DO i=istrm1,iendp1
1284 vxx(i,j)=v(i-1,j,k,nrhs)-2.0_r8*v(i,j,k,nrhs)+ &
1285# ifdef WEC_MELLOR
1286 & v_stokes(i-1,j,k)-2.0_r8*v_stokes(i,j,k)+ &
1287 & v_stokes(i+1,j,k)+ &
1288# endif
1289 & v(i+1,j,k,nrhs)
1290 tl_vxx(i,j)=tl_v(i-1,j,k,nrhs)-2.0_r8*tl_v(i,j,k,nrhs)+ &
1291# ifdef WEC_MELLOR
1292 & tl_v_stokes(i-1,j,k)-2.0_r8*tl_v_stokes(i,j,k)+ &
1293 & tl_v_stokes(i+1,j,k)+ &
1294# endif
1295 & tl_v(i+1,j,k,nrhs)
1296 END DO
1297 END DO
1298 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
1299 IF (domain(ng)%Western_Edge(tile)) THEN
1300 DO j=jstrv,jend
1301 vxx(istr-1,j)=vxx(istr,j)
1302 tl_vxx(istr-1,j)=tl_vxx(istr,j)
1303 END DO
1304 END IF
1305 END IF
1306 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
1307 IF (domain(ng)%Eastern_Edge(tile)) THEN
1308 DO j=jstrv,jend
1309 vxx(iend+1,j)=vxx(iend,j)
1310 tl_vxx(iend+1,j)=tl_vxx(iend,j)
1311 END DO
1312 END IF
1313 END IF
1314 DO j=jstrv-1,jend
1315 DO i=istr,iend+1
1316 huee(i,j)=huon(i,j-1,k)-2.0_r8*huon(i,j,k)+huon(i,j+1,k)
1317 tl_huee(i,j)=tl_huon(i,j-1,k)-2.0_r8*tl_huon(i,j,k)+ &
1318 & tl_huon(i,j+1,k)
1319 END DO
1320 END DO
1321# ifdef UV_C4ADVECTION
1322!
1323! Fourth-order, centered differences v-momentum horizontal advection.
1324!
1325 cff=1.0_r8/6.0_r8
1326 DO j=jstrv,jend
1327 DO i=istr,iend+1
1328 vfx(i,j)=0.25_r8*(v(i ,j,k,nrhs)+ &
1329# ifdef WEC_MELLOR
1330 & v_stokes(i ,j,k)+ &
1331 & v_stokes(i-1,j,k)+ &
1332# endif
1333 & v(i-1,j,k,nrhs)- &
1334 & cff*(vxx(i ,j)+ &
1335 & vxx(i-1,j)))* &
1336 & (huon(i,j ,k)+ &
1337 & huon(i,j-1,k)- &
1338 & cff*(huee(i,j )+ &
1339 & huee(i,j-1)))
1340 tl_vfx(i,j)=0.25_r8*((tl_v(i ,j,k,nrhs)+ &
1341# ifdef WEC_MELLOR
1342 & tl_v_stokes(i ,j,k)+ &
1343 & tl_v_stokes(i-1,j,k)+ &
1344# endif
1345 & tl_v(i-1,j,k,nrhs)- &
1346 & cff*(tl_vxx(i ,j)+ &
1347 & tl_vxx(i-1,j)))* &
1348 & (huon(i,j ,k)+ &
1349 & huon(i,j-1,k)- &
1350 & cff*(huee(i,j )+ &
1351 & huee(i,j-1)))+ &
1352 & (v(i ,j,k,nrhs)+ &
1353# ifdef WEC_MELLOR
1354 & v_stokes(i ,j,k)+ &
1355 & v_stokes(i-1,j,k)+ &
1356# endif
1357 & v(i-1,j,k,nrhs)- &
1358 & cff*(vxx(i ,j)+ &
1359 & vxx(i-1,j)))* &
1360 & (tl_huon(i,j ,k)+ &
1361 & tl_huon(i,j-1,k)- &
1362 & cff*(tl_huee(i,j )+ &
1363 & tl_huee(i,j-1))))- &
1364# ifdef TL_IOMS
1365 & vfx(i,j)
1366# endif
1367 END DO
1368 END DO
1369# else
1370!
1371! Third-order, upstream bias v-momentum advection with velocity
1372! dependent hyperdiffusion.
1373!
1374 DO j=jstrv,jend
1375 DO i=istr,iend+1
1376 cff1=v(i ,j,k,nrhs)+ &
1377# ifdef WEC_MELLOR
1378 & v_stokes(i ,j,k)+ &
1379 & v_stokes(i-1,j,k)+ &
1380# endif
1381 & v(i-1,j,k,nrhs)
1382 tl_cff1=tl_v(i ,j,k,nrhs)+ &
1383# ifdef WEC_MELLOR
1384 & tl_v_stokes(i ,j,k)+ &
1385 & tl_v_stokes(i-1,j,k)+ &
1386# endif
1387 & tl_v(i-1,j,k,nrhs)
1388 cff2=huon(i,j,k)+huon(i,j-1,k)
1389 tl_cff2=tl_huon(i,j,k)+tl_huon(i,j-1,k)
1390 IF (cff2.gt.0.0_r8) THEN
1391 cff=vxx(i-1,j)
1392 tl_cff=tl_vxx(i-1,j)
1393 ELSE
1394 cff=vxx(i,j)
1395 tl_cff=tl_vxx(i,j)
1396 END IF
1397 vfx(i,j)=0.25_r8*(cff1+gadv*cff)* &
1398 & (cff2+gadv*0.5_r8*(huee(i,j )+ &
1399 & huee(i,j-1)))
1400 tl_vfx(i,j)=0.25_r8* &
1401 & ((tl_cff1+gadv*tl_cff)* &
1402 & (cff2+gadv*0.5_r8*(huee(i,j )+ &
1403 & huee(i,j-1)))+ &
1404 & (cff1+gadv*cff)* &
1405 & (tl_cff2+gadv*0.5_r8*(tl_huee(i,j )+ &
1406 & tl_huee(i,j-1))))- &
1407# ifdef TL_IOMS
1408 & vfx(i,j)
1409# endif
1410 END DO
1411 END DO
1412# endif
1413 DO j=jstrvm1,jendp1
1414 DO i=istr,iend
1415 vee(i,j)=v(i,j-1,k,nrhs)-2.0_r8*v(i,j,k,nrhs)+ &
1416# ifdef WEC_MELLOR
1417 & v_stokes(i,j-1,k)-2.0_r8*v_stokes(i,j,k)+ &
1418 & v_stokes(i,j+1,k)+ &
1419# endif
1420 & v(i,j+1,k,nrhs)
1421 tl_vee(i,j)=tl_v(i,j-1,k,nrhs)-2.0_r8*tl_v(i,j,k,nrhs)+ &
1422# ifdef WEC_MELLOR
1423 & tl_v_stokes(i,j-1,k)-2.0_r8*tl_v_stokes(i,j,k)+ &
1424 & tl_v_stokes(i,j+1,k)+ &
1425# endif
1426 & tl_v(i,j+1,k,nrhs)
1427 hvee(i,j)=hvom(i,j-1,k)-2.0_r8*hvom(i,j,k)+hvom(i,j+1,k)
1428 tl_hvee(i,j)=tl_hvom(i,j-1,k)-2.0_r8*tl_hvom(i,j,k)+ &
1429 & tl_hvom(i,j+1,k)
1430 END DO
1431 END DO
1432 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
1433 IF (domain(ng)%Southern_Edge(tile)) THEN
1434 DO i=istr,iend
1435 vee(i,jstr)=vee(i,jstr+1)
1436 tl_vee(i,jstr)=tl_vee(i,jstr+1)
1437 hvee(i,jstr)=hvee(i,jstr+1)
1438 tl_hvee(i,jstr)=tl_hvee(i,jstr+1)
1439 END DO
1440 END IF
1441 END IF
1442 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
1443 IF (domain(ng)%Northern_Edge(tile)) THEN
1444 DO i=istr,iend
1445 vee(i,jend+1)=vee(i,jend)
1446 tl_vee(i,jend+1)=tl_vee(i,jend)
1447 hvee(i,jend+1)=hvee(i,jend)
1448 tl_hvee(i,jend+1)=tl_hvee(i,jend)
1449 END DO
1450 END IF
1451 END IF
1452# ifdef UV_C4ADVECTION
1453 cff=1.0_r8/6.0_r8
1454 DO j=jstrv-1,jend
1455 DO i=istr,iend
1456 vfe(i,j)=0.25_r8*(v(i,j,k,nrhs)+ &
1457# ifdef WEC_MELLOR
1458 & v_stokes(i,j ,k)+ &
1459 & v_stokes(i,j+1,k)+ &
1460# endif
1461 & v(i,j+1,k,nrhs)- &
1462 & cff*(vee(i,j )+ &
1463 & vee(i,j+1)))* &
1464 & (hvom(i,j ,k)+ &
1465 & hvom(i,j+1,k)- &
1466 & cff*(hvee(i,j )+ &
1467 & hvee(i,j+1)))
1468 tl_vfe(i,j)=0.25_r8*((tl_v(i,j ,k,nrhs)+ &
1469# ifdef WEC_MELLOR
1470 & tl_v_stokes(i,j ,k)+ &
1471 & tl_v_stokes(i,j+1,k)+ &
1472# endif
1473 & tl_v(i,j+1,k,nrhs)- &
1474 & cff*(tl_vee(i,j )+ &
1475 & tl_vee(i,j+1)))* &
1476 & (hvom(i,j ,k)+ &
1477 & hvom(i,j+1,k)- &
1478 & cff*(hvee(i,j )+ &
1479 & hvee(i,j+1)))+ &
1480 & (v(i,j ,k,nrhs)+ &
1481# ifdef WEC_MELLOR
1482 & v_stokes(i,j ,k)+ &
1483 & v_stokes(i,j+1,k)+ &
1484# endif
1485 & v(i,j+1,k,nrhs)- &
1486 & cff*(vee(i,j )+ &
1487 & vee(i,j+1)))* &
1488 & (tl_hvom(i,j ,k)+ &
1489 & tl_hvom(i,j+1,k)- &
1490 & cff*(tl_hvee(i,j )+ &
1491 & tl_hvee(i,j+1))))- &
1492# ifdef TL_IOMS
1493 & vfe(i,j)
1494# endif
1495 END DO
1496 END DO
1497# else
1498 DO j=jstrv-1,jend
1499 DO i=istr,iend
1500 cff1=v(i,j ,k,nrhs)+ &
1501# ifdef WEC_MELLOR
1502 & v_stokes(i,j ,k)+ &
1503 & v_stokes(i,j+1,k)+ &
1504# endif
1505 & v(i,j+1,k,nrhs)
1506 tl_cff1=tl_v(i,j ,k,nrhs)+ &
1507# ifdef WEC_MELLOR
1508 & tl_v_stokes(i,j ,k)+ &
1509 & tl_v_stokes(i,j+1,k)+ &
1510# endif
1511 & tl_v(i,j+1,k,nrhs)
1512 IF (cff1.gt.0.0_r8) THEN
1513 cff=vee(i,j)
1514 tl_cff=tl_vee(i,j)
1515 ELSE
1516 cff=vee(i,j+1)
1517 tl_cff=tl_vee(i,j+1)
1518 END IF
1519 vfe(i,j)=0.25_r8*(cff1+gadv*cff)* &
1520 & (hvom(i,j ,k)+ &
1521 & hvom(i,j+1,k)+ &
1522 & gadv*0.5_r8*(hvee(i,j )+ &
1523 & hvee(i,j+1)))
1524 tl_vfe(i,j)=0.25_r8* &
1525 & ((tl_cff1+gadv*tl_cff)* &
1526 & (hvom(i,j ,k)+ &
1527 & hvom(i,j+1,k)+ &
1528 & gadv*0.5_r8*(hvee(i,j )+ &
1529 & hvee(i,j+1)))+ &
1530 & (cff1+gadv*cff)* &
1531 & (tl_hvom(i,j ,k)+ &
1532 & tl_hvom(i,j+1,k)+ &
1533 & gadv*0.5_r8*(tl_hvee(i,j )+ &
1534 & tl_hvee(i,j+1))))- &
1535# ifdef TL_IOMS
1536 & vfe(i,j)
1537# endif
1538 END DO
1539 END DO
1540# endif
1541# endif
1542!
1543! Add in horizontal advection.
1544!
1545 DO j=jstr,jend
1546 DO i=istru,iend
1547!^ cff1=UFx(i,j)-UFx(i-1,j)
1548!^
1549 tl_cff1=tl_ufx(i,j)-tl_ufx(i-1,j)
1550!^ cff2=UFe(i,j+1)-UFe(i,j)
1551!^
1552 tl_cff2=tl_ufe(i,j+1)-tl_ufe(i,j)
1553!^ cff=cff1+cff2
1554!^
1555 tl_cff=tl_cff1+tl_cff2
1556!^ ru(i,j,k,nrhs)=ru(i,j,k,nrhs)-cff
1557!^
1558 tl_ru(i,j,k,nrhs)=tl_ru(i,j,k,nrhs)-tl_cff
1559# ifdef DIAGNOSTICS_UV
1560# ifdef CURVGRID
1561!! DiaRU(i,j,k,nrhs,M3xadv)=DiaRU(i,j,k,nrhs,M3xadv)-cff1
1562!! DiaRU(i,j,k,nrhs,M3yadv)=DiaRU(i,j,k,nrhs,M3yadv)-cff2
1563!! DiaRU(i,j,k,nrhs,M3hadv)=DiaRU(i,j,k,nrhs,M3hadv)-cff
1564# else
1565!! DiaRU(i,j,k,nrhs,M3xadv)=-cff1
1566!! DiaRU(i,j,k,nrhs,M3yadv)=-cff2
1567!! DiaRU(i,j,k,nrhs,M3hadv)=-cff
1568# endif
1569# endif
1570 END DO
1571 END DO
1572 DO j=jstrv,jend
1573 DO i=istr,iend
1574!^ cff1=VFx(i+1,j)-VFx(i,j)
1575!^
1576 tl_cff1=tl_vfx(i+1,j)-tl_vfx(i,j)
1577!^ cff2=VFe(i,j)-VFe(i,j-1)
1578!^
1579 tl_cff2=tl_vfe(i,j)-tl_vfe(i,j-1)
1580!^ cff=cff1+cff2
1581!^
1582 tl_cff=tl_cff1+tl_cff2
1583!^ rv(i,j,k,nrhs)=rv(i,j,k,nrhs)-cff
1584!^
1585 tl_rv(i,j,k,nrhs)=tl_rv(i,j,k,nrhs)-tl_cff
1586# ifdef DIAGNOSTICS_UV
1587# ifdef CURVGRID
1588!! DiaRV(i,j,k,nrhs,M3xadv)=DiaRV(i,j,k,nrhs,M3xadv)-cff1
1589!! DiaRV(i,j,k,nrhs,M3yadv)=DiaRV(i,j,k,nrhs,M3yadv)-cff2
1590!! DiaRV(i,j,k,nrhs,M3hadv)=DiaRV(i,j,k,nrhs,M3hadv)-cff
1591# else
1592!! DiaRV(i,j,k,nrhs,M3xadv)=-cff1
1593!! DiaRV(i,j,k,nrhs,M3yadv)=-cff2
1594!! DiaRV(i,j,k,nrhs,M3hadv)=-cff
1595# endif
1596# endif
1597 END DO
1598 END DO
1599# endif
1600# ifdef WEC_MELLOR
1601!
1602!-----------------------------------------------------------------------
1603! Add in radiation stress terms. Convert stresses to m4/s2.
1604!-----------------------------------------------------------------------
1605!
1606 DO j=jstr,jend
1607 DO i=istru,iend
1608!^ ru(i,j,k,nrhs)=ru(i,j,k,nrhs)- &
1609!^ & rustr3d(i,j,k)*om_u(i,j)*on_u(i,j)- &
1610!^ & rulag3d(i,j,k)
1611!^
1612 tl_ru(i,j,k,nrhs)=tl_ru(i,j,k,nrhs)- &
1613 & tl_rustr3d(i,j,k)*om_u(i,j)*on_u(i,j)- &
1614 & tl_rulag3d(i,j,k)
1615 END DO
1616 END DO
1617 DO j=jstrv,jend
1618 DO i=istr,iend
1619!^ rv(i,j,k,nrhs)=rv(i,j,k,nrhs)- &
1620!^ & rvstr3d(i,j,k)*om_v(i,j)*on_v(i,j)- &
1621!^ & rvlag3d(i,j,k)
1622!^
1623 tl_rv(i,j,k,nrhs)=tl_rv(i,j,k,nrhs)- &
1624 & tl_rvstr3d(i,j,k)*om_v(i,j)*on_v(i,j)- &
1625 & tl_rvlag3d(i,j,k)
1626 END DO
1627 END DO
1628# endif
1629
1630 END DO k_loop
1631!
1632 j_loop : DO j=jstr,jend
1633# ifdef UV_ADV
1634!
1635!-----------------------------------------------------------------------
1636! Add in vertical advection.
1637!-----------------------------------------------------------------------
1638!
1639# ifdef UV_SADVECTION
1640!
1641! Apply spline code to BASIC STATE u-momentum which should be in
1642! units of m/s. CF will be used by the tangent linear spline code.
1643!
1644 cff1=9.0_r8/16.0_r8
1645 cff2=1.0_r8/16.0_r8
1646 DO k=1,n(ng)
1647 DO i=istru,iend
1648 dc(i,k)=cff1*(hz(i ,j,k)+ &
1649 & hz(i-1,j,k))- &
1650 & cff2*(hz(i+1,j,k)+ &
1651 & hz(i-2,j,k))
1652 END DO
1653 END DO
1654 DO i=istru,iend
1655 fc(i,0)=0.0_r8
1656 cf(i,0)=0.0_r8
1657 END DO
1658 DO k=1,n(ng)-1
1659 DO i=istru,iend
1660 cff=1.0_r8/(2.0_r8*dc(i,k+1)+dc(i,k)*(2.0_r8-fc(i,k-1)))
1661 fc(i,k)=cff*dc(i,k+1)
1662 cf(i,k)=cff*(6.0_r8*(u(i,j,k+1,nrhs)- &
1663# ifdef WEC_MELLOR
1664 & u_stokes(i,j,k )+ &
1665 & u_stokes(i,j,k+1)- &
1666# endif
1667 & u(i,j,k ,nrhs))- &
1668 & dc(i,k)*cf(i,k-1))
1669 END DO
1670 END DO
1671 DO i=istru,iend
1672 cf(i,n(ng))=0.0_r8
1673 END DO
1674 DO k=n(ng)-1,1,-1
1675 DO i=istru,iend
1676 cf(i,k)=cf(i,k)-fc(i,k)*cf(i,k+1)
1677 END DO
1678 END DO
1679!
1680! Construct tangent linear conservative parabolic splines for the
1681! vertical derivatives "tl_CF" of u-momentum.
1682!
1683 cff1=9.0_r8/16.0_r8
1684 cff2=1.0_r8/16.0_r8
1685 DO k=1,n(ng)
1686 DO i=istru,iend
1687 dc(i,k)=cff1*(hz(i ,j,k)+ &
1688 & hz(i-1,j,k))- &
1689 & cff2*(hz(i+1,j,k)+ &
1690 & hz(i-2,j,k))
1691 tl_dc(i,k)=cff1*(tl_hz(i ,j,k)+ &
1692 & tl_hz(i-1,j,k))- &
1693 & cff2*(tl_hz(i+1,j,k)+ &
1694 & tl_hz(i-2,j,k))
1695 END DO
1696 END DO
1697 DO i=istru,iend
1698 fc(i,0)=0.0_r8
1699 tl_cf(i,0)=0.0_r8
1700 END DO
1701 DO k=1,n(ng)-1
1702 DO i=istru,iend
1703 cff=1.0_r8/(2.0_r8*dc(i,k+1)+dc(i,k)*(2.0_r8-fc(i,k-1)))
1704 fc(i,k)=cff*dc(i,k+1)
1705# ifdef TL_IOMS
1706 tl_cf(i,k)=cff*(6.0_r8*(tl_u(i,j,k+1,nrhs)- &
1707# ifdef WEC_MELLOR
1708 & tl_u_stokes(i,j,k )+ &
1709 & tl_u_stokes(i,j,k+1)- &
1710# endif
1711 & tl_u(i,j,k ,nrhs))- &
1712 & ((tl_dc(i,k)-dc(i,k))*cf(i,k-1)+ &
1713 & 2.0_r8*(tl_dc(i,k )-dc(i,k )+ &
1714 & tl_dc(i,k+1)-dc(i,k+1))*cf(i,k)+ &
1715 & (tl_dc(i,k+1)-dc(i,k+1))*cf(i,k+1))- &
1716 & dc(i,k)*tl_cf(i,k-1))
1717# else
1718 tl_cf(i,k)=cff*(6.0_r8*(tl_u(i,j,k+1,nrhs)- &
1719# ifdef WEC_MELLOR
1720 & tl_u_stokes(i,j,k )+ &
1721 & tl_u_stokes(i,j,k+1)- &
1722# endif
1723 & tl_u(i,j,k ,nrhs))- &
1724 & (tl_dc(i,k)*cf(i,k-1)+ &
1725 & 2.0_r8*(tl_dc(i,k)+tl_dc(i,k+1))*cf(i,k)+ &
1726 & tl_dc(i,k+1)*cf(i,k+1))- &
1727 & dc(i,k)*tl_cf(i,k-1))
1728# endif
1729 END DO
1730 END DO
1731 DO i=istru,iend
1732 tl_cf(i,n(ng))=0.0_r8
1733 END DO
1734 DO k=n(ng)-1,1,-1
1735 DO i=istru,iend
1736 tl_cf(i,k)=tl_cf(i,k)-fc(i,k)*tl_cf(i,k+1)
1737 END DO
1738 END DO
1739!
1740! Compute spline-interpolated, vertical advective u-momentum flux.
1741!
1742 cff3=1.0_r8/3.0_r8
1743 cff4=1.0_r8/6.0_r8
1744 DO k=1,n(ng)-1
1745 DO i=istru,iend
1746 fc(i,k)=(cff1*(w(i ,j,k)+ &
1747 & w(i-1,j,k))- &
1748 & cff2*(w(i+1,j,k)+ &
1749 & w(i-2,j,k)))* &
1750 & (u(i,j,k,nrhs)+ &
1751# ifdef WEC_MELLOR
1752 & u_stokes(i,j,k)+ &
1753# endif
1754 & dc(i,k)*(cff3*cf(i,k )+ &
1755 & cff4*cf(i,k-1)))
1756 tl_fc(i,k)=(cff1*(tl_w(i ,j,k)+ &
1757 & tl_w(i-1,j,k))- &
1758 & cff2*(tl_w(i+1,j,k)+ &
1759 & tl_w(i-2,j,k)))* &
1760 & (u(i,j,k,nrhs)+ &
1761# ifdef WEC_MELLOR
1762 & u_stokes(i,j,k)+ &
1763# endif
1764 & dc(i,k)*(cff3*cf(i,k )+ &
1765 & cff4*cf(i,k-1)))+ &
1766 & (cff1*(w(i ,j,k)+ &
1767 & w(i-1,j,k))- &
1768 & cff2*(w(i+1,j,k)+ &
1769 & w(i-2,j,k)))* &
1770 & (tl_u(i,j,k,nrhs)+ &
1771# ifdef WEC_MELLOR
1772 & tl_u_stokes(i,j,k)+ &
1773# endif
1774 & dc(i,k)*(cff3*tl_cf(i,k )+ &
1775 & cff4*tl_cf(i,k-1))+ &
1776 & tl_dc(i,k)*(cff3*cf(i,k )+ &
1777 & cff4*cf(i,k-1)))- &
1778# ifdef TL_IOMS
1779 & 2.0_r8*fc(i,k)+(cff1*(w(i ,j,k)+ &
1780 & w(i-1,j,k))- &
1781 & cff2*(w(i+1,j,k)+ &
1782 & w(i-2,j,k)))* &
1783# ifdef WEC_MELLOR
1784 & (u(i,j,k,nrhs)+u_stokes(i,j,k))
1785# else
1786 & u(i,j,k,nrhs)
1787# endif
1788# endif
1789 END DO
1790 END DO
1791 DO i=istru,iend
1792!^ FC(i,N(ng))=0.0_r8
1793!^
1794 tl_fc(i,n(ng))=0.0_r8
1795!^ FC(i,0)=0.0_r8
1796!^
1797 tl_fc(i,0)=0.0_r8
1798 END DO
1799# elif defined UV_C2ADVECTION
1800 DO k=1,n(ng)-1
1801 DO i=istru,iend
1802 fc(i,k)=0.25_r8*(u(i,j,k ,nrhs)+ &
1803# ifdef WEC_MELLOR
1804 & u_stokes(i,j,k )+ &
1805 & u_stokes(i,j,k+1)+ &
1806# endif
1807 & u(i,j,k+1,nrhs))* &
1808 & (w(i ,j,k)+ &
1809 & w(i-1,j,k))
1810 tl_fc(i,k)=0.25_r8*((tl_u(i,j,k ,nrhs)+ &
1811# ifdef WEC_MELLOR
1812 & tl_u_stokes(i,j,k )+ &
1813 & tl_u_stokes(i,j,k+1)+ &
1814# endif
1815 & tl_u(i,j,k+1,nrhs))* &
1816 & (w(i ,j,k)+ &
1817 & w(i-1,j,k))+ &
1818 & (u(i,j,k ,nrhs)+ &
1819# ifdef WEC_MELLOR
1820 & u_stokes(i,j,k )+ &
1821 & u_stokes(i,j,k+1)+ &
1822# endif
1823 & u(i,j,k+1,nrhs))* &
1824 & (tl_w(i ,j,k)+ &
1825 & tl_w(i-1,j,k)))- &
1826# ifdef TL_IOMS
1827 & fc(i,k)
1828# endif
1829 END DO
1830 END DO
1831 DO i=istru,iend
1832!^ FC(i,0)=0.0_r8
1833!^
1834 tl_fc(i,0)=0.0_r8
1835!^ FC(i,N(ng))=0.0_r8
1836!^
1837 tl_fc(i,n(ng))=0.0_r8
1838 END DO
1839# elif defined UV_C4ADVECTION
1840 cff1=9.0_r8/32.0_r8
1841 cff2=1.0_r8/32.0_r8
1842 DO k=2,n(ng)-2
1843 DO i=istru,iend
1844 fc(i,k)=(cff1*(u(i,j,k ,nrhs)+ &
1845# ifdef WEC_MELLOR
1846 & u_stokes(i,j,k )+ &
1847 & u_stokes(i,j,k+1)+ &
1848# endif
1849 & u(i,j,k+1,nrhs))- &
1850 & cff2*(u(i,j,k-1,nrhs)+ &
1851# ifdef WEC_MELLOR
1852 & u_stokes(i,j,k-1)+ &
1853 & u_stokes(i,j,k+2)+ &
1854# endif
1855 & u(i,j,k+2,nrhs)))* &
1856 & (w(i ,j,k)+ &
1857 & w(i-1,j,k))
1858 tl_fc(i,k)=(cff1*(tl_u(i,j,k ,nrhs)+ &
1859# ifdef WEC_MELLOR
1860 & tl_u_stokes(i,j,k )+ &
1861 & tl_u_stokes(i,j,k+1)+ &
1862# endif
1863 & tl_u(i,j,k+1,nrhs))- &
1864 & cff2*(tl_u(i,j,k-1,nrhs)+ &
1865# ifdef WEC_MELLOR
1866 & tl_u_stokes(i,j,k-1)+ &
1867 & tl_u_stokes(i,j,k+2)+ &
1868# endif
1869 & tl_u(i,j,k+2,nrhs)))* &
1870 & (w(i ,j,k)+ &
1871 & w(i-1,j,k))+ &
1872 & (cff1*(u(i,j,k ,nrhs)+ &
1873# ifdef WEC_MELLOR
1874 & u_stokes(i,j,k )+ &
1875 & u_stokes(i,j,k+1)+ &
1876# endif
1877 & u(i,j,k+1,nrhs))- &
1878 & cff2*(u(i,j,k-1,nrhs)+ &
1879# ifdef WEC_MELLOR
1880 & u_stokes(i,j,k-1)+ &
1881 & u_stokes(i,j,k+2)+ &
1882# endif
1883 & u(i,j,k+2,nrhs)))* &
1884 & (tl_w(i ,j,k)+ &
1885 & tl_w(i-1,j,k))- &
1886# ifdef TL_IOMS
1887 & fc(i,k)
1888# endif
1889 END DO
1890 END DO
1891 DO i=istru,iend
1892!^ FC(i,N(ng))=0.0_r8
1893!^
1894 tl_fc(i,n(ng))=0.0_r8
1895 fc(i,n(ng)-1)=(cff1*(u(i,j,n(ng)-1,nrhs)+ &
1896# ifdef WEC_MELLOR
1897 & u_stokes(i,j,n(ng)-1)+ &
1898 & u_stokes(i,j,n(ng) )+ &
1899# endif
1900 & u(i,j,n(ng) ,nrhs))- &
1901 & cff2*(u(i,j,n(ng)-2,nrhs)+ &
1902# ifdef WEC_MELLOR
1903 & u_stokes(i,j,n(ng)-2)+ &
1904 & u_stokes(i,j,n(ng) )+ &
1905# endif
1906 & u(i,j,n(ng) ,nrhs)))* &
1907 & (w(i ,j,n(ng)-1)+ &
1908 & w(i-1,j,n(ng)-1))
1909 tl_fc(i,n(ng)-1)=(cff1*(tl_u(i,j,n(ng)-1,nrhs)+ &
1910# ifdef WEC_MELLOR
1911 & tl_u_stokes(i,j,n(ng)-1)+ &
1912 & tl_u_stokes(i,j,n(ng) )+ &
1913# endif
1914 & tl_u(i,j,n(ng) ,nrhs))- &
1915 & cff2*(tl_u(i,j,n(ng)-2,nrhs)+ &
1916# ifdef WEC_MELLOR
1917 & tl_u_stokes(i,j,n(ng)-2)+ &
1918 & tl_u_stokes(i,j,n(ng) )+ &
1919# endif
1920 & tl_u(i,j,n(ng) ,nrhs)))* &
1921 & (w(i ,j,n(ng)-1)+ &
1922 & w(i-1,j,n(ng)-1))+ &
1923 & (cff1*(u(i,j,n(ng)-1,nrhs)+ &
1924# ifdef WEC_MELLOR
1925 & u_stokes(i,j,n(ng)-1)+ &
1926 & u_stokes(i,j,n(ng) )+ &
1927# endif
1928 & u(i,j,n(ng) ,nrhs))- &
1929 & cff2*(u(i,j,n(ng)-2,nrhs)+ &
1930# ifdef WEC_MELLOR
1931 & u_stokes(i,j,n(ng)-2)+ &
1932 & u_stokes(i,j,n(ng) )+ &
1933# endif
1934 & u(i,j,n(ng) ,nrhs)))* &
1935 & (tl_w(i ,j,n(ng)-1)+ &
1936 & tl_w(i-1,j,n(ng)-1))- &
1937# ifdef TL_IOMS
1938 & fc(i,n(ng)-1)
1939# endif
1940 fc(i,1)=(cff1*(u(i,j,1,nrhs)+ &
1941# ifdef WEC_MELLOR
1942 & u_stokes(i,j,1)+ &
1943 & u_stokes(i,j,2)+ &
1944# endif
1945 & u(i,j,2,nrhs))- &
1946 & cff2*(u(i,j,1,nrhs)+ &
1947# ifdef WEC_MELLOR
1948 & u_stokes(i,j,1)+ &
1949 & u_stokes(i,j,3)+ &
1950# endif
1951 & u(i,j,3,nrhs)))* &
1952 & (w(i ,j,1)+ &
1953 & w(i-1,j,1))
1954 tl_fc(i,1)=(cff1*(tl_u(i,j,1,nrhs)+ &
1955# ifdef WEC_MELLOR
1956 & tl_u_stokes(i,j,1)+ &
1957 & tl_u_stokes(i,j,2)+ &
1958# endif
1959 & tl_u(i,j,2,nrhs))- &
1960 & cff2*(tl_u(i,j,1,nrhs)+ &
1961# ifdef WEC_MELLOR
1962 & tl_u_stokes(i,j,1)+ &
1963 & tl_u_stokes(i,j,3)+ &
1964# endif
1965 & tl_u(i,j,3,nrhs)))* &
1966 & (w(i ,j,1)+ &
1967 & w(i-1,j,1))+ &
1968 & (cff1*(u(i,j,1,nrhs)+ &
1969# ifdef WEC_MELLOR
1970 & u_stokes(i,j,1)+ &
1971 & u_stokes(i,j,2)+ &
1972# endif
1973 & u(i,j,2,nrhs))- &
1974 & cff2*(u(i,j,1,nrhs)+ &
1975# ifdef WEC_MELLOR
1976 & u_stokes(i,j,1)+ &
1977 & u_stokes(i,j,3)+ &
1978# endif
1979 & u(i,j,3,nrhs)))* &
1980 & (tl_w(i ,j,1)+ &
1981 & tl_w(i-1,j,1))- &
1982# ifdef TL_IOMS
1983 & fc(i,1)
1984# endif
1985!^ FC(i,0)=0.0_r8
1986!^
1987 tl_fc(i,0)=0.0_r8
1988 END DO
1989# else
1990 cff1=9.0_r8/16.0_r8
1991 cff2=1.0_r8/16.0_r8
1992 DO k=2,n(ng)-2
1993 DO i=istru,iend
1994 fc(i,k)=(cff1*(u(i,j,k ,nrhs)+ &
1995# ifdef WEC_MELLOR
1996 & u_stokes(i,j,k )+ &
1997 & u_stokes(i,j,k+1)+ &
1998# endif
1999 & u(i,j,k+1,nrhs))- &
2000 & cff2*(u(i,j,k-1,nrhs)+ &
2001# ifdef WEC_MELLOR
2002 & u_stokes(i,j,k-1)+ &
2003 & u_stokes(i,j,k+2)+ &
2004# endif
2005 & u(i,j,k+2,nrhs)))* &
2006 & (cff1*(w(i ,j,k)+ &
2007 & w(i-1,j,k))- &
2008 & cff2*(w(i+1,j,k)+ &
2009 & w(i-2,j,k)))
2010 tl_fc(i,k)=(cff1*(tl_u(i,j,k ,nrhs)+ &
2011# ifdef WEC_MELLOR
2012 & tl_u_stokes(i,j,k )+ &
2013 & tl_u_stokes(i,j,k+1)+ &
2014# endif
2015 & tl_u(i,j,k+1,nrhs))- &
2016 & cff2*(tl_u(i,j,k-1,nrhs)+ &
2017# ifdef WEC_MELLOR
2018 & tl_u_stokes(i,j,k-1)+ &
2019 & tl_u_stokes(i,j,k+2)+ &
2020# endif
2021 & tl_u(i,j,k+2,nrhs)))* &
2022 & (cff1*(w(i ,j,k)+ &
2023 & w(i-1,j,k))- &
2024 & cff2*(w(i+1,j,k)+ &
2025 & w(i-2,j,k)))+ &
2026 & (cff1*(u(i,j,k ,nrhs)+ &
2027# ifdef WEC_MELLOR
2028 & u_stokes(i,j,k )+ &
2029 & u_stokes(i,j,k+1)+ &
2030# endif
2031 & u(i,j,k+1,nrhs))- &
2032 & cff2*(u(i,j,k-1,nrhs)+ &
2033# ifdef WEC_MELLOR
2034 & u_stokes(i,j,k-1)+ &
2035 & u_stokes(i,j,k+2)+ &
2036# endif
2037 & u(i,j,k+2,nrhs)))* &
2038 & (cff1*(tl_w(i ,j,k)+ &
2039 & tl_w(i-1,j,k))- &
2040 & cff2*(tl_w(i+1,j,k)+ &
2041 & tl_w(i-2,j,k)))- &
2042# ifdef TL_IOMS
2043 & fc(i,k)
2044# endif
2045 END DO
2046 END DO
2047 DO i=istru,iend
2048!^ FC(i,N(ng))=0.0_r8
2049!^
2050 tl_fc(i,n(ng))=0.0_r8
2051 fc(i,n(ng)-1)=(cff1*(u(i,j,n(ng)-1,nrhs)+ &
2052# ifdef WEC_MELLOR
2053 & u_stokes(i,j,n(ng)-1)+ &
2054 & u_stokes(i,j,n(ng) )+ &
2055# endif
2056 & u(i,j,n(ng) ,nrhs))- &
2057 & cff2*(u(i,j,n(ng)-2,nrhs)+ &
2058# ifdef WEC_MELLOR
2059 & u_stokes(i,j,n(ng)-2)+ &
2060 & u_stokes(i,j,n(ng) )+ &
2061# endif
2062 & u(i,j,n(ng) ,nrhs)))* &
2063 & (cff1*(w(i ,j,n(ng)-1)+ &
2064 & w(i-1,j,n(ng)-1))- &
2065 & cff2*(w(i+1,j,n(ng)-1)+ &
2066 & w(i-2,j,n(ng)-1)))
2067 tl_fc(i,n(ng)-1)=(cff1*(tl_u(i,j,n(ng)-1,nrhs)+ &
2068# ifdef WEC_MELLOR
2069 & tl_u_stokes(i,j,n(ng)-1)+ &
2070 & tl_u_stokes(i,j,n(ng) )+ &
2071# endif
2072 & tl_u(i,j,n(ng) ,nrhs))- &
2073 & cff2*(tl_u(i,j,n(ng)-2,nrhs)+ &
2074# ifdef WEC_MELLOR
2075 & tl_u_stokes(i,j,n(ng)-2)+ &
2076 & tl_u_stokes(i,j,n(ng) )+ &
2077# endif
2078 & tl_u(i,j,n(ng) ,nrhs)))* &
2079 & (cff1*(w(i ,j,n(ng)-1)+ &
2080 & w(i-1,j,n(ng)-1))- &
2081 & cff2*(w(i+1,j,n(ng)-1)+ &
2082 & w(i-2,j,n(ng)-1)))+ &
2083 & (cff1*(u(i,j,n(ng)-1,nrhs)+ &
2084# ifdef WEC_MELLOR
2085 & u_stokes(i,j,n(ng)-1)+ &
2086 & u_stokes(i,j,n(ng) )+ &
2087# endif
2088 & u(i,j,n(ng) ,nrhs))- &
2089 & cff2*(u(i,j,n(ng)-2,nrhs)+ &
2090# ifdef WEC_MELLOR
2091 & u_stokes(i,j,n(ng)-2)+ &
2092 & u_stokes(i,j,n(ng) )+ &
2093# endif
2094 & u(i,j,n(ng) ,nrhs)))* &
2095 & (cff1*(tl_w(i ,j,n(ng)-1)+ &
2096 & tl_w(i-1,j,n(ng)-1))- &
2097 & cff2*(tl_w(i+1,j,n(ng)-1)+ &
2098 & tl_w(i-2,j,n(ng)-1)))- &
2099# ifdef TL_IOMS
2100 & fc(i,n(ng)-1)
2101# endif
2102 fc(i,1)=(cff1*(u(i,j,1,nrhs)+ &
2103# ifdef WEC_MELLOR
2104 & u_stokes(i,j,1)+ &
2105 & u_stokes(i,j,2)+ &
2106# endif
2107 & u(i,j,2,nrhs))- &
2108 & cff2*(u(i,j,1,nrhs)+ &
2109# ifdef WEC_MELLOR
2110 & u_stokes(i,j,1)+ &
2111 & u_stokes(i,j,3)+ &
2112# endif
2113 & u(i,j,3,nrhs)))* &
2114 & (cff1*(w(i ,j,1)+ &
2115 & w(i-1,j,1))- &
2116 & cff2*(w(i+1,j,1)+ &
2117 & w(i-2,j,1)))
2118 tl_fc(i,1)=(cff1*(tl_u(i,j,1,nrhs)+ &
2119# ifdef WEC_MELLOR
2120 & tl_u_stokes(i,j,1)+ &
2121 & tl_u_stokes(i,j,2)+ &
2122# endif
2123 & tl_u(i,j,2,nrhs))- &
2124 & cff2*(tl_u(i,j,1,nrhs)+ &
2125# ifdef WEC_MELLOR
2126 & tl_u_stokes(i,j,1)+ &
2127 & tl_u_stokes(i,j,3)+ &
2128# endif
2129 & tl_u(i,j,3,nrhs)))* &
2130 & (cff1*(w(i ,j,1)+ &
2131 & w(i-1,j,1))- &
2132 & cff2*(w(i+1,j,1)+ &
2133 & w(i-2,j,1)))+ &
2134 & (cff1*(u(i,j,1,nrhs)+ &
2135# ifdef WEC_MELLOR
2136 & u_stokes(i,j,1)+ &
2137 & u_stokes(i,j,2)+ &
2138# endif
2139 & u(i,j,2,nrhs))- &
2140 & cff2*(u(i,j,1,nrhs)+ &
2141# ifdef WEC_MELLOR
2142 & u_stokes(i,j,1)+ &
2143 & u_stokes(i,j,3)+ &
2144# endif
2145 & u(i,j,3,nrhs)))* &
2146 & (cff1*(tl_w(i ,j,1)+ &
2147 & tl_w(i-1,j,1))- &
2148 & cff2*(tl_w(i+1,j,1)+ &
2149 & tl_w(i-2,j,1)))- &
2150# ifdef TL_IOMS
2151 & fc(i,1)
2152# endif
2153!^ FC(i,0)=0.0_r8
2154!^
2155 tl_fc(i,0)=0.0_r8
2156 END DO
2157# endif
2158 DO k=1,n(ng)
2159 DO i=istru,iend
2160!^ cff=FC(i,k)-FC(i,k-1)
2161!^
2162 tl_cff=tl_fc(i,k)-tl_fc(i,k-1)
2163!^ ru(i,j,k,nrhs)=ru(i,j,k,nrhs)-cff
2164!^
2165 tl_ru(i,j,k,nrhs)=tl_ru(i,j,k,nrhs)-tl_cff
2166# ifdef DIAGNOSTICS_UV
2167!! DiaRU(i,j,k,nrhs,M3vadv)=-cff
2168# endif
2169 END DO
2170 END DO
2171 IF (j.ge.jstrv) THEN
2172# ifdef UV_SADVECTION
2173!
2174! Apply spline code to BASIC STATE v-momentum which should be in
2175! units of m/s. CF will be used by the tangent linear spline code.
2176!
2177 cff1=9.0_r8/16.0_r8
2178 cff2=1.0_r8/16.0_r8
2179 DO k=1,n(ng)
2180 DO i=istr,iend
2181 dc(i,k)=(cff1*(hz(i,j ,k)+ &
2182 & hz(i,j-1,k))- &
2183 & cff2*(hz(i,j+1,k)+ &
2184 & hz(i,j-2,k)))
2185 END DO
2186 END DO
2187 DO i=istr,iend
2188 fc(i,0)=0.0_r8
2189 cf(i,0)=0.0_r8
2190 END DO
2191 DO k=1,n(ng)-1
2192 DO i=istr,iend
2193 cff=1.0_r8/(2.0_r8*dc(i,k+1)+dc(i,k)*(2.0_r8-fc(i,k-1)))
2194 fc(i,k)=cff*dc(i,k+1)
2195 cf(i,k)=cff*(6.0_r8*(v(i,j,k+1,nrhs)- &
2196# ifdef WEC_MELLOR
2197 & v_stokes(i,j,k )+ &
2198 & v_stokes(i,j,k+1)- &
2199# endif
2200 & v(i,j,k ,nrhs))- &
2201 & dc(i,k)*cf(i,k-1))
2202 END DO
2203 END DO
2204 DO i=istr,iend
2205 cf(i,n(ng))=0.0_r8
2206 END DO
2207 DO k=n(ng)-1,1,-1
2208 DO i=istr,iend
2209 cf(i,k)=cf(i,k)-fc(i,k)*cf(i,k+1)
2210 END DO
2211 END DO
2212!
2213! Construct tangent linear conservative parabolic splines for the
2214! vertical derivatives "tl_CF" of v-momentum.
2215!
2216 cff1=9.0_r8/16.0_r8
2217 cff2=1.0_r8/16.0_r8
2218 DO k=1,n(ng)
2219 DO i=istr,iend
2220 dc(i,k)=(cff1*(hz(i,j ,k)+ &
2221 & hz(i,j-1,k))- &
2222 & cff2*(hz(i,j+1,k)+ &
2223 & hz(i,j-2,k)))
2224 tl_dc(i,k)=(cff1*(tl_hz(i,j ,k)+ &
2225 & tl_hz(i,j-1,k))- &
2226 & cff2*(tl_hz(i,j+1,k)+ &
2227 & tl_hz(i,j-2,k)))
2228 END DO
2229 END DO
2230 DO i=istr,iend
2231 fc(i,0)=0.0_r8
2232 tl_cf(i,0)=0.0_r8
2233 END DO
2234 DO k=1,n(ng)-1
2235 DO i=istr,iend
2236 cff=1.0_r8/(2.0_r8*dc(i,k+1)+dc(i,k)*(2.0_r8-fc(i,k-1)))
2237 fc(i,k)=cff*dc(i,k+1)
2238# ifdef TL_IOMS
2239 tl_cf(i,k)=cff*(6.0_r8*(tl_v(i,j,k+1,nrhs)- &
2240# ifdef WEC_MELLOR
2241 & tl_v_stokes(i,j,k )+ &
2242 & tl_v_stokes(i,j,k+1)- &
2243# endif
2244 & tl_v(i,j,k ,nrhs))- &
2245 & ((tl_dc(i,k)-dc(i,k))*cf(i,k-1)+ &
2246 & 2.0_r8*(tl_dc(i,k )-dc(i,k )+ &
2247 & tl_dc(i,k+1)-dc(i,k+1))*cf(i,k)+ &
2248 & (tl_dc(i,k+1)-dc(i,k+1))*cf(i,k+1))- &
2249 & dc(i,k)*tl_cf(i,k-1))
2250# else
2251 tl_cf(i,k)=cff*(6.0_r8*(tl_v(i,j,k+1,nrhs)- &
2252# ifdef WEC_MELLOR
2253 & tl_v_stokes(i,j,k )+ &
2254 & tl_v_stokes(i,j,k+1)- &
2255# endif
2256 & tl_v(i,j,k ,nrhs))- &
2257 & (tl_dc(i,k)*cf(i,k-1)+ &
2258 & 2.0_r8*(tl_dc(i,k )+ &
2259 & tl_dc(i,k+1))*cf(i,k)+ &
2260 & tl_dc(i,k+1)*cf(i,k+1))- &
2261 & dc(i,k)*tl_cf(i,k-1))
2262# endif
2263 END DO
2264 END DO
2265 DO i=istr,iend
2266 tl_cf(i,n(ng))=0.0_r8
2267 END DO
2268 DO k=n(ng)-1,1,-1
2269 DO i=istr,iend
2270 tl_cf(i,k)=tl_cf(i,k)-fc(i,k)*tl_cf(i,k+1)
2271 END DO
2272 END DO
2273!
2274! Compute spline-interpolated, vertical advective v-momentum flux.
2275!
2276 cff3=1.0_r8/3.0_r8
2277 cff4=1.0_r8/6.0_r8
2278 DO k=1,n(ng)-1
2279 DO i=istr,iend
2280 fc(i,k)=(cff1*(w(i,j ,k)+ &
2281 & w(i,j-1,k))- &
2282 & cff2*(w(i,j+1,k)+ &
2283 & w(i,j-2,k)))* &
2284 & (v(i,j,k,nrhs)+ &
2285# ifdef WEC_MELLOR
2286 & v_stokes(i,j,k)+ &
2287# endif
2288 & dc(i,k)*(cff3*cf(i,k )+ &
2289 & cff4*cf(i,k-1)))
2290 tl_fc(i,k)=(cff1*(tl_w(i,j ,k)+ &
2291 & tl_w(i,j-1,k))- &
2292 & cff2*(tl_w(i,j+1,k)+ &
2293 & tl_w(i,j-2,k)))* &
2294 & (v(i,j,k,nrhs)+ &
2295# ifdef WEC_MELLOR
2296 & v_stokes(i,j,k)+ &
2297# endif
2298 & dc(i,k)*(cff3*cf(i,k )+ &
2299 & cff4*cf(i,k-1)))+ &
2300 & (cff1*(w(i,j ,k)+ &
2301 & w(i,j-1,k))- &
2302 & cff2*(w(i,j+1,k)+ &
2303 & w(i,j-2,k)))* &
2304 & (tl_v(i,j,k,nrhs)+ &
2305# ifdef WEC_MELLOR
2306 & tl_v_stokes(i,j,k)+ &
2307# endif
2308 & dc(i,k)*(cff3*tl_cf(i,k )+ &
2309 & cff4*tl_cf(i,k-1))+ &
2310 & tl_dc(i,k)*(cff3*cf(i,k )+ &
2311 & cff4*cf(i,k-1)))- &
2312# ifdef TL_IOMS
2313 & 2.0_r8*fc(i,k)+(cff1*(w(i,j ,k)+ &
2314 & w(i,j-1,k))- &
2315 & cff2*(w(i,j+1,k)+ &
2316 & w(i,j-2,k)))* &
2317# ifdef WEC_MELLOR
2318 & (v(i,j,k,nrhs)+v_stokes(i,j,k))
2319# else
2320 & v(i,j,k,nrhs)
2321# endif
2322# endif
2323 END DO
2324 END DO
2325 DO i=istr,iend
2326!^ FC(i,N(ng))=0.0_r8
2327!^
2328 tl_fc(i,n(ng))=0.0_r8
2329!^ FC(i,0)=0.0_r8
2330!^
2331 tl_fc(i,0)=0.0_r8
2332 END DO
2333# elif defined UV_C2ADVECTION
2334!
2335! Second-order, centered differences vertical advection.
2336!
2337 DO k=1,n(ng)-1
2338 DO i=istr,iend
2339 fc(i,k)=0.25_r8*(v(i,j,k ,nrhs)+ &
2340# ifdef WEC_MELLOR
2341 & v_stokes(i,j,k )+ &
2342 & v_stokes(i,j,k+1)+ &
2343# endif
2344 & v(i,j,k+1,nrhs))* &
2345 & (w(i,j ,k)+ &
2346 & w(i,j-1,k))
2347 tl_fc(i,k)=0.25_r8*((tl_v(i,j,k ,nrhs)+ &
2348# ifdef WEC_MELLOR
2349 & tl_v_stokes(i,j,k )+ &
2350 & tl_v_stokes(i,j,k+1)+ &
2351# endif
2352 & tl_v(i,j,k+1,nrhs))* &
2353 & (w(i,j ,k)+ &
2354 & w(i,j-1,k))+ &
2355 & (v(i,j,k ,nrhs)+ &
2356# ifdef WEC_MELLOR
2357 & v_stokes(i,j,k )+ &
2358 & v_stokes(i,j,k+1)+ &
2359# endif
2360 & v(i,j,k+1,nrhs))* &
2361 & (tl_w(i,j ,k)+ &
2362 & tl_w(i,j-1,k)))- &
2363# ifdef TL_IOMS
2364 & fc(i,k)
2365# endif
2366 END DO
2367 END DO
2368 DO i=istr,iend
2369!^ FC(i,0)=0.0_r8
2370!^
2371 tl_fc(i,0)=0.0_r8
2372!^ FC(i,N(ng))=0.0_r8
2373!^
2374 tl_fc(i,n(ng))=0.0_r8
2375 END DO
2376# elif defined UV_C4ADVECTION
2377!
2378! Forth-order, centered differences vertical advection.
2379!
2380 cff1=9.0_r8/32.0_r8
2381 cff2=1.0_r8/32.0_r8
2382 DO k=2,n(ng)-2
2383 DO i=istr,iend
2384 fc(i,k)=(cff1*(v(i,j,k ,nrhs)+ &
2385# ifdef WEC_MELLOR
2386 & v_stokes(i,j,k )+ &
2387 & v_stokes(i,j,k+1)+ &
2388# endif
2389 & v(i,j,k+1,nrhs))- &
2390 & cff2*(v(i,j,k-1,nrhs)+ &
2391# ifdef WEC_MELLOR
2392 & v_stokes(i,j,k-1)+ &
2393 & v_stokes(i,j,k+2)+ &
2394# endif
2395 & v(i,j,k+2,nrhs)))* &
2396 & (w(i,j ,k)+ &
2397 & w(i,j-1,k))
2398 tl_fc(i,k)=(cff1*(tl_v(i,j,k ,nrhs)+ &
2399# ifdef WEC_MELLOR
2400 & tl_v_stokes(i,j,k )+ &
2401 & tl_v_stokes(i,j,k+1)+ &
2402# endif
2403 & tl_v(i,j,k+1,nrhs))- &
2404 & cff2*(tl_v(i,j,k-1,nrhs)+ &
2405# ifdef WEC_MELLOR
2406 & tl_v_stokes(i,j,k-1)+ &
2407 & tl_v_stokes(i,j,k+2)+ &
2408# endif
2409 & tl_v(i,j,k+2,nrhs)))* &
2410 & (w(i,j ,k)+ &
2411 & w(i,j-1,k))+ &
2412 & (cff1*(v(i,j,k ,nrhs)+ &
2413# ifdef WEC_MELLOR
2414 & v_stokes(i,j,k )+ &
2415 & v_stokes(i,j,k+1)+ &
2416# endif
2417 & v(i,j,k+1,nrhs))- &
2418 & cff2*(v(i,j,k-1,nrhs)+ &
2419# ifdef WEC_MELLOR
2420 & v_stokes(i,j,k-1)+ &
2421 & v_stokes(i,j,k+2)+ &
2422# endif
2423 & v(i,j,k+2,nrhs)))* &
2424 & (tl_w(i,j ,k)+ &
2425 & tl_w(i,j-1,k))- &
2426# ifdef TL_IOMS
2427 & fc(i,k)
2428# endif
2429 END DO
2430 END DO
2431 DO i=istr,iend
2432!^ FC(i,N(ng))=0.0_r8
2433!^
2434 tl_fc(i,n(ng))=0.0_r8
2435 fc(i,n(ng)-1)=(cff1*(v(i,j,n(ng)-1,nrhs)+ &
2436# ifdef WEC_MELLOR
2437 & v_stokes(i,j,n(ng)-1)+ &
2438 & v_stokes(i,j,n(ng) )+ &
2439# endif
2440 & v(i,j,n(ng) ,nrhs))- &
2441 & cff2*(v(i,j,n(ng)-2,nrhs)+ &
2442# ifdef WEC_MELLOR
2443 & v_stokes(i,j,n(ng)-2)+ &
2444 & v_stokes(i,j,n(ng) )+ &
2445# endif
2446 & v(i,j,n(ng) ,nrhs)))* &
2447 & (w(i,j ,n(ng)-1)+ &
2448 & w(i,j-1,n(ng)-1))
2449 tl_fc(i,n(ng)-1)=(cff1*(tl_v(i,j,n(ng)-1,nrhs)+ &
2450# ifdef WEC_MELLOR
2451 & tl_v_stokes(i,j,n(ng)-1)+ &
2452 & tl_v_stokes(i,j,n(ng) )+ &
2453# endif
2454 & tl_v(i,j,n(ng) ,nrhs))- &
2455 & cff2*(tl_v(i,j,n(ng)-2,nrhs)+ &
2456# ifdef WEC_MELLOR
2457 & tl_v_stokes(i,j,n(ng)-2)+ &
2458 & tl_v_stokes(i,j,n(ng) )+ &
2459# endif
2460 & tl_v(i,j,n(ng) ,nrhs)))* &
2461 & (w(i,j ,n(ng)-1)+ &
2462 & w(i,j-1,n(ng)-1))+ &
2463 & (cff1*(v(i,j,n(ng)-1,nrhs)+ &
2464# ifdef WEC_MELLOR
2465 & v_stokes(i,j,n(ng)-1)+ &
2466 & v_stokes(i,j,n(ng) )+ &
2467# endif
2468 & v(i,j,n(ng) ,nrhs))- &
2469 & cff2*(v(i,j,n(ng)-2,nrhs)+ &
2470# ifdef WEC_MELLOR
2471 & v_stokes(i,j,n(ng)-2)+ &
2472 & v_stokes(i,j,n(ng) )+ &
2473# endif
2474 & v(i,j,n(ng) ,nrhs)))* &
2475 & (tl_w(i,j ,n(ng)-1)+ &
2476 & tl_w(i,j-1,n(ng)-1))- &
2477# ifdef TL_IOMS
2478 & fc(i,n(ng)-1)
2479# endif
2480 fc(i,1)=(cff1*(v(i,j,1,nrhs)+ &
2481# ifdef WEC_MELLOR
2482 & v_stokes(i,j,1)+ &
2483 & v_stokes(i,j,2)+ &
2484# endif
2485 & v(i,j,2,nrhs))- &
2486 & cff2*(v(i,j,1,nrhs)+ &
2487# ifdef WEC_MELLOR
2488 & v_stokes(i,j,1)+ &
2489 & v_stokes(i,j,3)+ &
2490# endif
2491 & v(i,j,3,nrhs)))* &
2492 & (w(i,j ,1)+ &
2493 & w(i,j-1,1))
2494 tl_fc(i,1)=(cff1*(tl_v(i,j,1,nrhs)+ &
2495# ifdef WEC_MELLOR
2496 & tl_v_stokes(i,j,1)+ &
2497 & tl_v_stokes(i,j,2)+ &
2498# endif
2499 & tl_v(i,j,2,nrhs))- &
2500 & cff2*(tl_v(i,j,1,nrhs)+ &
2501# ifdef WEC_MELLOR
2502 & tl_v_stokes(i,j,1)+ &
2503 & tl_v_stokes(i,j,3)+ &
2504# endif
2505 & tl_v(i,j,3,nrhs)))* &
2506 & (w(i,j ,1)+ &
2507 & w(i,j-1,1))+ &
2508 & (cff1*(v(i,j,1,nrhs)+ &
2509# ifdef WEC_MELLOR
2510 & v_stokes(i,j,1)+ &
2511 & v_stokes(i,j,2)+ &
2512# endif
2513 & v(i,j,2,nrhs))- &
2514 & cff2*(v(i,j,1,nrhs)+ &
2515# ifdef WEC_MELLOR
2516 & v_stokes(i,j,1)+ &
2517 & v_stokes(i,j,3)+ &
2518# endif
2519 & v(i,j,3,nrhs)))* &
2520 & (tl_w(i,j ,1)+ &
2521 & tl_w(i,j-1,1))- &
2522# ifdef TL_IOMS
2523 & fc(i,1)
2524# endif
2525!^ FC(i,0)=0.0_r8
2526!^
2527 tl_fc(i,0)=0.0_r8
2528 END DO
2529# else
2530 cff1=9.0_r8/16.0_r8
2531 cff2=1.0_r8/16.0_r8
2532 DO k=2,n(ng)-2
2533 DO i=istr,iend
2534 fc(i,k)=(cff1*(v(i,j,k ,nrhs)+ &
2535# ifdef WEC_MELLOR
2536 & v_stokes(i,j,k )+ &
2537 & v_stokes(i,j,k+1)+ &
2538# endif
2539 & v(i,j,k+1,nrhs))- &
2540 & cff2*(v(i,j,k-1,nrhs)+ &
2541# ifdef WEC_MELLOR
2542 & v_stokes(i,j,k-1)+ &
2543 & v_stokes(i,j,k+2)+ &
2544# endif
2545 & v(i,j,k+2,nrhs)))* &
2546 & (cff1*(w(i,j ,k)+ &
2547 & w(i,j-1,k))- &
2548 & cff2*(w(i,j+1,k)+ &
2549 & w(i,j-2,k)))
2550 tl_fc(i,k)=(cff1*(tl_v(i,j,k ,nrhs)+ &
2551# ifdef WEC_MELLOR
2552 & tl_v_stokes(i,j,k )+ &
2553 & tl_v_stokes(i,j,k+1)+ &
2554# endif
2555 & tl_v(i,j,k+1,nrhs))- &
2556 & cff2*(tl_v(i,j,k-1,nrhs)+ &
2557# ifdef WEC_MELLOR
2558 & tl_v_stokes(i,j,k-1)+ &
2559 & tl_v_stokes(i,j,k+2)+ &
2560# endif
2561 & tl_v(i,j,k+2,nrhs)))* &
2562 & (cff1*(w(i,j ,k)+ &
2563 & w(i,j-1,k))- &
2564 & cff2*(w(i,j+1,k)+ &
2565 & w(i,j-2,k)))+ &
2566 & (cff1*(v(i,j,k ,nrhs)+ &
2567# ifdef WEC_MELLOR
2568 & v_stokes(i,j,k )+ &
2569 & v_stokes(i,j,k+1)+ &
2570# endif
2571 & v(i,j,k+1,nrhs))- &
2572 & cff2*(v(i,j,k-1,nrhs)+ &
2573# ifdef WEC_MELLOR
2574 & v_stokes(i,j,k-1)+ &
2575 & v_stokes(i,j,k+2)+ &
2576# endif
2577 & v(i,j,k+2,nrhs)))* &
2578 & (cff1*(tl_w(i,j ,k)+ &
2579 & tl_w(i,j-1,k))- &
2580 & cff2*(tl_w(i,j+1,k)+ &
2581 & tl_w(i,j-2,k)))- &
2582# ifdef TL_IOMS
2583 & fc(i,k)
2584# endif
2585 END DO
2586 END DO
2587 DO i=istr,iend
2588!^ FC(i,N(ng))=0.0_r8
2589!^
2590 tl_fc(i,n(ng))=0.0_r8
2591 fc(i,n(ng)-1)=(cff1*(v(i,j,n(ng)-1,nrhs)+ &
2592# ifdef WEC_MELLOR
2593 & v_stokes(i,j,n(ng)-1)+ &
2594 & v_stokes(i,j,n(ng) )+ &
2595# endif
2596 & v(i,j,n(ng) ,nrhs))- &
2597 & cff2*(v(i,j,n(ng)-2,nrhs)+ &
2598# ifdef WEC_MELLOR
2599 & v_stokes(i,j,n(ng)-2)+ &
2600 & v_stokes(i,j,n(ng) )+ &
2601# endif
2602 & v(i,j,n(ng) ,nrhs)))* &
2603 & (cff1*(w(i,j ,n(ng)-1)+ &
2604 & w(i,j-1,n(ng)-1))- &
2605 & cff2*(w(i,j+1,n(ng)-1)+ &
2606 & w(i,j-2,n(ng)-1)))
2607 tl_fc(i,n(ng)-1)=(cff1*(tl_v(i,j,n(ng)-1,nrhs)+ &
2608# ifdef WEC_MELLOR
2609 & tl_v_stokes(i,j,n(ng)-1)+ &
2610 & tl_v_stokes(i,j,n(ng) )+ &
2611# endif
2612 & tl_v(i,j,n(ng) ,nrhs))- &
2613 & cff2*(tl_v(i,j,n(ng)-2,nrhs)+ &
2614# ifdef WEC_MELLOR
2615 & tl_v_stokes(i,j,n(ng)-2)+ &
2616 & tl_v_stokes(i,j,n(ng) )+ &
2617# endif
2618 & tl_v(i,j,n(ng) ,nrhs)))* &
2619 & (cff1*(w(i,j ,n(ng)-1)+ &
2620 & w(i,j-1,n(ng)-1))- &
2621 & cff2*(w(i,j+1,n(ng)-1)+ &
2622 & w(i,j-2,n(ng)-1)))+ &
2623 & (cff1*(v(i,j,n(ng)-1,nrhs)+ &
2624# ifdef WEC_MELLOR
2625 & v_stokes(i,j,n(ng)-1)+ &
2626 & v_stokes(i,j,n(ng) )+ &
2627# endif
2628 & v(i,j,n(ng) ,nrhs))- &
2629 & cff2*(v(i,j,n(ng)-2,nrhs)+ &
2630# ifdef WEC_MELLOR
2631 & v_stokes(i,j,n(ng)-2)+ &
2632 & v_stokes(i,j,n(ng) )+ &
2633# endif
2634 & v(i,j,n(ng) ,nrhs)))* &
2635 & (cff1*(tl_w(i,j ,n(ng)-1)+ &
2636 & tl_w(i,j-1,n(ng)-1))- &
2637 & cff2*(tl_w(i,j+1,n(ng)-1)+ &
2638 & tl_w(i,j-2,n(ng)-1)))- &
2639# ifdef TL_IOMS
2640 & fc(i,n(ng)-1)
2641# endif
2642 fc(i,1)=(cff1*(v(i,j,1,nrhs)+ &
2643# ifdef WEC_MELLOR
2644 & v_stokes(i,j,1)+ &
2645 & v_stokes(i,j,2)+ &
2646# endif
2647 & v(i,j,2,nrhs))- &
2648 & cff2*(v(i,j,1,nrhs)+ &
2649# ifdef WEC_MELLOR
2650 & v_stokes(i,j,1)+ &
2651 & v_stokes(i,j,3)+ &
2652# endif
2653 & v(i,j,3,nrhs)))* &
2654 & (cff1*(w(i,j ,1)+ &
2655 & w(i,j-1,1))- &
2656 & cff2*(w(i,j+1,1)+ &
2657 & w(i,j-2,1)))
2658 tl_fc(i,1)=(cff1*(tl_v(i,j,1,nrhs)+ &
2659# ifdef WEC_MELLOR
2660 & tl_v_stokes(i,j,1)+ &
2661 & tl_v_stokes(i,j,2)+ &
2662# endif
2663 & tl_v(i,j,2,nrhs))- &
2664 & cff2*(tl_v(i,j,1,nrhs)+ &
2665# ifdef WEC_MELLOR
2666 & tl_v_stokes(i,j,1)+ &
2667 & tl_v_stokes(i,j,3)+ &
2668# endif
2669 & tl_v(i,j,3,nrhs)))* &
2670 & (cff1*(w(i,j ,1)+ &
2671 & w(i,j-1,1))- &
2672 & cff2*(w(i,j+1,1)+ &
2673 & w(i,j-2,1)))+ &
2674 & (cff1*(v(i,j,1,nrhs)+ &
2675# ifdef WEC_MELLOR
2676 & v_stokes(i,j,1)+ &
2677 & v_stokes(i,j,2)+ &
2678# endif
2679 & v(i,j,2,nrhs))- &
2680 & cff2*(v(i,j,1,nrhs)+ &
2681# ifdef WEC_MELLOR
2682 & v_stokes(i,j,1)+ &
2683 & v_stokes(i,j,3)+ &
2684# endif
2685 & v(i,j,3,nrhs)))* &
2686 & (cff1*(tl_w(i,j ,1)+ &
2687 & tl_w(i,j-1,1))- &
2688 & cff2*(tl_w(i,j+1,1)+ &
2689 & tl_w(i,j-2,1)))- &
2690# ifdef TL_IOMS
2691 & fc(i,1)
2692# endif
2693!^ FC(i,0)=0.0_r8
2694!^
2695 tl_fc(i,0)=0.0_r8
2696 END DO
2697# endif
2698 DO k=1,n(ng)
2699 DO i=istr,iend
2700!^ cff=FC(i,k)-FC(i,k-1)
2701!^
2702 tl_cff=tl_fc(i,k)-tl_fc(i,k-1)
2703!^ rv(i,j,k,nrhs)=rv(i,j,k,nrhs)-cff
2704!^
2705 tl_rv(i,j,k,nrhs)=tl_rv(i,j,k,nrhs)-tl_cff
2706# ifdef DIAGNOSTICS_UV
2707!! DiaRV(i,j,k,nrhs,M3vadv)=-cff
2708# endif
2709 END DO
2710 END DO
2711 END IF
2712# endif
2713!
2714!-----------------------------------------------------------------------
2715! Compute forcing term for the 2D momentum equations.
2716!-----------------------------------------------------------------------
2717!
2718! Vertically integrate baroclinic right-hand-side terms. If not
2719! body force stresses, add in the difference between surface and
2720! bottom stresses.
2721!
2722 DO i=istru,iend
2723# ifdef WET_DRY_NOT_YET
2724!^ ru(i,j,1,nrhs)=ru(i,j,1,nrhs)*umask_wet(i,j)
2725!^
2726 tl_ru(i,j,1,nrhs)=tl_ru(i,j,1,nrhs)*umask_wet(i,j)
2727# endif
2728!^ rufrc(i,j)=ru(i,j,1,nrhs)
2729!^
2730 tl_rufrc(i,j)=tl_ru(i,j,1,nrhs)
2731# ifdef DIAGNOSTICS_UV
2732!! DiaRUfrc(i,j,3,M2pgrd)=DiaRU(i,j,1,nrhs,M3pgrd)
2733# ifdef UV_COR
2734!! DiaRUfrc(i,j,3,M2fcor)=DiaRU(i,j,1,nrhs,M3fcor)
2735# endif
2736# ifdef UV_ADV
2737!! DiaRUfrc(i,j,3,M2xadv)=DiaRU(i,j,1,nrhs,M3xadv)
2738!! DiaRUfrc(i,j,3,M2yadv)=DiaRU(i,j,1,nrhs,M3yadv)
2739!! DiaRUfrc(i,j,3,M2hadv)=DiaRU(i,j,1,nrhs,M3hadv)
2740# endif
2741# ifdef WEC_MELLOR
2742!! DiaRUfrc(i,j,3,M2hrad)=DiaRU(i,j,1,nrhs,M3hrad)
2743# endif
2744# if defined UV_VIS2 || defined UV_VIS4
2745!! DiaRUfrc(i,j,3,M2xvis)=0.0_r8
2746!! DiaRUfrc(i,j,3,M2yvis)=0.0_r8
2747!! DiaRUfrc(i,j,3,M2hvis)=0.0_r8
2748# endif
2749# ifdef BODYFORCE
2750!! DiaRUfrc(i,j,3,M2strs)=DiaRU(i,j,1,nrhs,M3vvis)
2751# endif
2752# endif
2753 END DO
2754 DO k=2,n(ng)
2755 DO i=istru,iend
2756# ifdef WET_DRY_NOT_YET
2757!^ ru(i,j,k,nrhs)=ru(i,j,k,nrhs)*umask_wet(i,j)
2758!^
2759 tl_ru(i,j,k,nrhs)=tl_ru(i,j,k,nrhs)*umask_wet(i,j)
2760# endif
2761!^ rufrc(i,j)=rufrc(i,j)+ru(i,j,k,nrhs)
2762!^
2763 tl_rufrc(i,j)=tl_rufrc(i,j)+tl_ru(i,j,k,nrhs)
2764# ifdef DIAGNOSTICS_UV
2765!! DiaRUfrc(i,j,3,M2pgrd)=DiaRUfrc(i,j,3,M2pgrd)+ &
2766!! & DiaRU(i,j,k,nrhs,M3pgrd)
2767# ifdef UV_COR
2768!! DiaRUfrc(i,j,3,M2fcor)=DiaRUfrc(i,j,3,M2fcor)+ &
2769!! & DiaRU(i,j,k,nrhs,M3fcor)
2770# endif
2771# ifdef UV_ADV
2772!! DiaRUfrc(i,j,3,M2xadv)=DiaRUfrc(i,j,3,M2xadv)+ &
2773!! & DiaRU(i,j,k,nrhs,M3xadv)
2774!! DiaRUfrc(i,j,3,M2yadv)=DiaRUfrc(i,j,3,M2yadv)+ &
2775!! & DiaRU(i,j,k,nrhs,M3yadv)
2776!! DiaRUfrc(i,j,3,M2hadv)=DiaRUfrc(i,j,3,M2hadv)+ &
2777!! & DiaRU(i,j,k,nrhs,M3hadv)
2778# endif
2779# ifdef WEC_MELLOR
2780!! DiaRUfrc(i,j,3,M2hrad)=DiaRUfrc(i,j,3,M2hrad)+ &
2781!! & DiaRU(i,j,k,nrhs,M3hrad)
2782# endif
2783# ifdef BODYFORCE
2784!! DiaRUfrc(i,j,3,M2strs)=DiaRUfrc(i,j,3,M2strs)+ &
2785!! & DiaRU(i,j,k,nrhs,M3vvis)
2786# endif
2787# endif
2788 END DO
2789 END DO
2790# ifndef BODYFORCE
2791 DO i=istru,iend
2792 cff=om_u(i,j)*on_u(i,j)
2793!^ cff1= sustr(i,j)*cff
2794!^
2795 tl_cff1= tl_sustr(i,j)*cff
2796!^ cff2=-bustr(i,j)*cff
2797!^
2798 tl_cff2=-tl_bustr(i,j)*cff
2799!^ rufrc(i,j)=rufrc(i,j)+cff1+cff2
2800!^
2801 tl_rufrc(i,j)=tl_rufrc(i,j)+tl_cff1+tl_cff2
2802# ifdef WET_DRY_NOT_YET
2803!^ rufrc(i,j)=rufrc(i,j)*umask_wet(i,j)
2804!^
2805 tl_rufrc(i,j)=tl_rufrc(i,j)*umask_wet(i,j)
2806# endif
2807# ifdef DIAGNOSTICS_UV
2808!! DiaRUfrc(i,j,3,M2sstr)=cff1
2809!! DiaRUfrc(i,j,3,M2bstr)=cff2
2810# endif
2811 END DO
2812# endif
2813 IF (j.ge.jstrv) THEN
2814 DO i=istr,iend
2815# ifdef WET_DRY_NOT_YET
2816!^ rv(i,j,1,nrhs)=rv(i,j,1,nrhs)*vmask_wet(i,j)
2817!^
2818 tl_rv(i,j,1,nrhs)=tl_rv(i,j,1,nrhs)*vmask_wet(i,j)
2819# endif
2820!^ rvfrc(i,j)=rv(i,j,1,nrhs)
2821!^
2822 tl_rvfrc(i,j)=tl_rv(i,j,1,nrhs)
2823# ifdef DIAGNOSTICS_UV
2824!! DiaRVfrc(i,j,3,M2pgrd)=DiaRV(i,j,1,nrhs,M3pgrd)
2825# ifdef UV_COR
2826!! DiaRVfrc(i,j,3,M2fcor)=DiaRV(i,j,1,nrhs,M3fcor)
2827# endif
2828# ifdef UV_ADV
2829!! DiaRVfrc(i,j,3,M2xadv)=DiaRV(i,j,1,nrhs,M3xadv)
2830!! DiaRVfrc(i,j,3,M2yadv)=DiaRV(i,j,1,nrhs,M3yadv)
2831!! DiaRVfrc(i,j,3,M2hadv)=DiaRV(i,j,1,nrhs,M3hadv)
2832# endif
2833# ifdef WEC_MELLOR
2834!! DiaRVfrc(i,j,3,M2hrad)=DiaRV(i,j,1,nrhs,M3hrad)
2835# endif
2836# if defined UV_VIS2 || defined UV_VIS4
2837!! DiaRVfrc(i,j,3,M2hvis)=0.0_r8
2838!! DiaRVfrc(i,j,3,M2xvis)=0.0_r8
2839!! DiaRVfrc(i,j,3,M2yvis)=0.0_r8
2840# endif
2841# ifdef BODYFORCE
2842!! DiaRVfrc(i,j,3,M2strs)=DiaRV(i,j,1,nrhs,M3vvis)
2843# endif
2844# endif
2845 END DO
2846 DO k=2,n(ng)
2847 DO i=istr,iend
2848# ifdef WET_DRY_NOT_YET
2849!^ rv(i,j,k,nrhs)=rv(i,j,k,nrhs)*vmask_wet(i,j)
2850!^
2851 tl_rv(i,j,k,nrhs)=tl_rv(i,j,k,nrhs)*vmask_wet(i,j)
2852# endif
2853!^ rvfrc(i,j)=rvfrc(i,j)+rv(i,j,k,nrhs)
2854!^
2855 tl_rvfrc(i,j)=tl_rvfrc(i,j)+tl_rv(i,j,k,nrhs)
2856# ifdef DIAGNOSTICS_UV
2857!! DiaRVfrc(i,j,3,M2pgrd)=DiaRVfrc(i,j,3,M2pgrd)+ &
2858!! & DiaRV(i,j,k,nrhs,M3pgrd)
2859# ifdef UV_COR
2860!! DiaRVfrc(i,j,3,M2fcor)=DiaRVfrc(i,j,3,M2fcor)+ &
2861!! & DiaRV(i,j,k,nrhs,M3fcor)
2862# endif
2863# ifdef UV_ADV
2864!! DiaRVfrc(i,j,3,M2xadv)=DiaRVfrc(i,j,3,M2xadv)+ &
2865!! & DiaRV(i,j,k,nrhs,M3xadv)
2866!! DiaRVfrc(i,j,3,M2yadv)=DiaRVfrc(i,j,3,M2yadv)+ &
2867!! & DiaRV(i,j,k,nrhs,M3yadv)
2868!! DiaRVfrc(i,j,3,M2hadv)=DiaRVfrc(i,j,3,M2hadv)+ &
2869!! & DiaRV(i,j,k,nrhs,M3hadv)
2870# endif
2871# ifdef WEC_MELLOR
2872!! DiaRVfrc(i,j,3,M2hrad)=DiaRVfrc(i,j,3,M2hrad)+ &
2873!! & DiaRV(i,j,k,nrhs,M3hrad)
2874# endif
2875# ifdef BODYFORCE
2876!! DiaRVfrc(i,j,3,M2strs)=DiaRVfrc(i,j,3,M2strs)+ &
2877!! & DiaRV(i,j,k,nrhs,M3vvis)
2878# endif
2879# endif
2880 END DO
2881 END DO
2882# ifndef BODYFORCE
2883 DO i=istr,iend
2884 cff=om_v(i,j)*on_v(i,j)
2885!^ cff1= svstr(i,j)*cff
2886!^
2887 tl_cff1= tl_svstr(i,j)*cff
2888!^ cff2=-bvstr(i,j)*cff
2889!^
2890 tl_cff2=-tl_bvstr(i,j)*cff
2891!^ rvfrc(i,j)=rvfrc(i,j)+cff1+cff2
2892!^
2893 tl_rvfrc(i,j)=tl_rvfrc(i,j)+tl_cff1+tl_cff2
2894# ifdef WET_DRY_NOT_YET
2895!^ rvfrc(i,j)=rvfrc(i,j)*vmask_wet(i,j)
2896!^
2897 tl_rvfrc(i,j)=tl_rvfrc(i,j)*vmask_wet(i,j)
2898# endif
2899# ifdef DIAGNOSTICS_UV
2900!! DiaRVfrc(i,j,3,M2sstr)=cff1
2901!! DiaRVfrc(i,j,3,M2bstr)=cff2
2902# endif
2903 END DO
2904# endif
2905 END IF
2906 END DO j_loop
2907
2908 RETURN
2909 END SUBROUTINE rp_rhs3d_tile
2910#endif
2911 END MODULE rp_rhs3d_mod
type(t_clima), dimension(:), allocatable clima
Definition mod_clima.F:153
type(t_coupling), dimension(:), allocatable coupling
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 irpm
Definition mod_param.F:664
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer, dimension(:), allocatable levbfrc
logical, dimension(:), allocatable ewperiodic
integer, parameter iwest
logical, dimension(:), allocatable nsperiodic
logical, dimension(:), allocatable lnudgem3clm
logical, dimension(:,:), allocatable compositegrid
integer, parameter isouth
integer, parameter ieast
integer, parameter inorth
integer, dimension(:), allocatable levsfrc
integer, dimension(:), allocatable nrhs
subroutine, public rp_pre_step3d(ng, tile)
subroutine, public rp_prsgrd(ng, tile)
Definition rp_prsgrd31.h:36
subroutine, public rp_rhs3d(ng, tile)
Definition rp_rhs3d.F:29
subroutine rp_rhs3d_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, hz, tl_hz, huon, tl_huon, hvom, tl_hvom, dmde, dndx, fomn, om_u, om_v, on_u, on_v, pm, pn, umask_wet, vmask_wet, bustr, tl_bustr, bvstr, tl_bvstr, sustr, tl_sustr, svstr, tl_svstr, u, tl_u, v, tl_v, w, tl_w, tl_rufrc, tl_rvfrc, tl_ru, tl_rv)
Definition rp_rhs3d.F:254
subroutine, public rp_t3dmix2(ng, tile)
subroutine, public rp_t3dmix4(ng, tile)
subroutine, public rp_t3drelax(ng, tile)
Definition rp_t3drelax.F:29
subroutine, public rp_uv3dmix2(ng, tile)
subroutine, public rp_uv3dmix4(ng, tile)
subroutine, public rp_uv3drelax(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