ROMS
Loading...
Searching...
No Matches
rp_step3d_uv.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 routine time-steps the representers tangent linear horizontal !
13! momentum equations. The vertical viscosity terms are time-stepped !
14! using an implicit algorithm. !
15! !
16! BASIC STATE variables needed: u, v, ru, rv, Akv, !
17! Hz, Huon, Hvom, z_r, z_w, !
18! DU_avg1, DU_avg2, DV_avg1, DV_avg2, !
19! Qsrc !
20! !
21!=======================================================================
22!
23 USE mod_param
24 USE mod_coupling
25# ifdef DIAGNOSTICS_UV
26!! USE mod_diags
27# endif
28 USE mod_forces
29 USE mod_grid
30 USE mod_mixing
31 USE mod_ncparam
32 USE mod_ocean
33 USE mod_scalars
34 USE mod_sources
35!
38# ifdef DISTRIBUTE
40# endif
41 USE rp_u3dbc_mod, ONLY : rp_u3dbc_tile
42 USE rp_v3dbc_mod, ONLY : rp_v3dbc_tile
43# ifdef UV_DESTAGGERED
45# endif
46!
47 implicit none
48!
49 PRIVATE
50 PUBLIC :: rp_step3d_uv
51!
52 CONTAINS
53!
54!***********************************************************************
55 SUBROUTINE rp_step3d_uv (ng, tile)
56!***********************************************************************
57!
58 USE mod_stepping
59!
60! Imported variable declarations.
61!
62 integer, intent(in) :: ng, tile
63!
64! Local variable declarations.
65!
66 character (len=*), parameter :: myfile = &
67 & __FILE__
68!
69# include "tile.h"
70!
71# ifdef PROFILE
72 CALL wclock_on (ng, irpm, 34, __line__, myfile)
73# endif
74 CALL rp_step3d_uv_tile (ng, tile, &
75 & lbi, ubi, lbj, ubj, &
76 & imins, imaxs, jmins, jmaxs, &
77 & nrhs(ng), nstp(ng), nnew(ng), &
78# ifdef MASKING
79 & grid(ng) % umask, &
80 & grid(ng) % vmask, &
81# endif
82# ifdef WET_DRY_NOT_YET
83 & grid(ng) % umask_wet, &
84 & grid(ng) % vmask_wet, &
85# endif
86 & grid(ng) % om_v, &
87 & grid(ng) % on_u, &
88 & grid(ng) % pm, &
89 & grid(ng) % pn, &
90 & grid(ng) % Hz, &
91 & grid(ng) % tl_Hz, &
92 & grid(ng) % z_r, &
93 & grid(ng) % tl_z_r, &
94 & grid(ng) % z_w, &
95 & grid(ng) % tl_z_w, &
96 & mixing(ng) % Akv, &
97 & mixing(ng) % tl_Akv, &
98 & coupling(ng) % DU_avg1, &
99 & coupling(ng) % tl_DU_avg1, &
100 & coupling(ng) % DV_avg1, &
101 & coupling(ng) % tl_DV_avg1, &
102 & coupling(ng) % DU_avg2, &
103 & coupling(ng) % tl_DU_avg2, &
104 & coupling(ng) % DV_avg2, &
105 & coupling(ng) % tl_DV_avg2, &
106 & ocean(ng) % tl_ru, &
107 & ocean(ng) % tl_rv, &
108# ifdef DIAGNOSTICS_UV
109!! & DIAGS(ng) % DiaU2wrk, &
110!! & DIAGS(ng) % DiaV2wrk, &
111!! & DIAGS(ng) % DiaU2int, &
112!! & DIAGS(ng) % DiaV2int, &
113!! & DIAGS(ng) % DiaU3wrk, &
114!! & DIAGS(ng) % DiaV3wrk, &
115!! & DIAGS(ng) % DiaRU, &
116!! & DIAGS(ng) % DiaRV, &
117# endif
118 & ocean(ng) % u, &
119 & ocean(ng) % tl_u, &
120 & ocean(ng) % v, &
121 & ocean(ng) % tl_v, &
122 & ocean(ng) % tl_ubar, &
123 & ocean(ng) % tl_vbar, &
124# ifdef NEARSHORE_MELLOR
125 & ocean(ng) % ubar_stokes, &
126 & ocean(ng) % tl_ubar_stokes, &
127 & ocean(ng) % vbar_stokes, &
128 & ocean(ng) % tl_vbar_stokes, &
129 & ocean(ng) % u_stokes, &
130 & ocean(ng) % tl_u_stokes, &
131 & ocean(ng) % v_stokes, &
132 & ocean(ng) % tl_v_stokes, &
133# endif
134 & grid(ng) % Huon, &
135 & grid(ng) % tl_Huon, &
136 & grid(ng) % Hvom, &
137 & grid(ng) % tl_Hvom)
138# ifdef PROFILE
139 CALL wclock_off (ng, irpm, 34, __line__, myfile)
140# endif
141!
142 RETURN
143 END SUBROUTINE rp_step3d_uv
144!
145!***********************************************************************
146 SUBROUTINE rp_step3d_uv_tile (ng, tile, &
147 & LBi, UBi, LBj, UBj, &
148 & IminS, ImaxS, JminS, JmaxS, &
149 & nrhs, nstp, nnew, &
150# ifdef MASKING
151 & umask, vmask, &
152# endif
153# ifdef WET_DRY_NOT_YET
154 & umask_wet, vmask_wet, &
155# endif
156 & om_v, on_u, pm, pn, &
157 & Hz, tl_Hz, &
158 & z_r, tl_z_r, &
159 & z_w, tl_z_w, &
160 & Akv, tl_Akv, &
161 & DU_avg1, tl_DU_avg1, &
162 & DV_avg1, tl_DV_avg1, &
163 & DU_avg2, tl_DU_avg2, &
164 & DV_avg2, tl_DV_avg2, &
165 & tl_ru, tl_rv, &
166# ifdef DIAGNOSTICS_UV
167!! & DiaU2wrk, DiaV2wrk, &
168!! & DiaU2int, DiaV2int, &
169!! & DiaU3wrk, DiaV3wrk, &
170!! & DiaRU, DiaRV, &
171# endif
172 & u, tl_u, &
173 & v, tl_v, &
174 & tl_ubar, tl_vbar, &
175# ifdef NEARSHORE_MELLOR
176 & ubar_stokes, tl_ubar_stokes, &
177 & vbar_stokes, tl_vbar_stokes, &
178 & u_stokes, tl_u_stokes, &
179 & v_stokes, tl_v_stokes, &
180# endif
181 & Huon, tl_Huon, &
182 & Hvom, tl_Hvom)
183!***********************************************************************
184!
185! Imported variable declarations.
186!
187 integer, intent(in) :: ng, tile
188 integer, intent(in) :: LBi, UBi, LBj, UBj
189 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
190 integer, intent(in) :: nrhs, nstp, nnew
191!
192# ifdef ASSUMED_SHAPE
193# ifdef MASKING
194 real(r8), intent(in) :: umask(LBi:,LBj:)
195 real(r8), intent(in) :: vmask(LBi:,LBj:)
196# endif
197# ifdef WET_DRY_NOT_YET
198 real(r8), intent(in) :: umask_wet(LBi:,LBj:)
199 real(r8), intent(in) :: vmask_wet(LBi:,LBj:)
200# endif
201 real(r8), intent(in) :: om_v(LBi:,LBj:)
202 real(r8), intent(in) :: on_u(LBi:,LBj:)
203 real(r8), intent(in) :: pm(LBi:,LBj:)
204 real(r8), intent(in) :: pn(LBi:,LBj:)
205 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
206 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
207 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
208 real(r8), intent(in) :: Akv(LBi:,LBj:,0:)
209 real(r8), intent(in) :: DU_avg1(LBi:,LBj:)
210 real(r8), intent(in) :: DV_avg1(LBi:,LBj:)
211 real(r8), intent(in) :: DU_avg2(LBi:,LBj:)
212 real(r8), intent(in) :: DV_avg2(LBi:,LBj:)
213 real(r8), intent(in) :: u(LBi:,LBj:,:,:)
214 real(r8), intent(in) :: v(LBi:,LBj:,:,:)
215# ifdef NEARSHORE_MELLOR
216 real(r8), intent(in) :: ubar_stokes(LBi:,LBj:)
217 real(r8), intent(in) :: vbar_stokes(LBi:,LBj:)
218 real(r8), intent(in) :: tl_ubar_stokes(LBi:,LBj:)
219 real(r8), intent(in) :: tl_vbar_stokes(LBi:,LBj:)
220# endif
221 real(r8), intent(in) :: tl_Hz(LBi:,LBj:,:)
222 real(r8), intent(in) :: tl_z_r(LBi:,LBj:,:)
223 real(r8), intent(in) :: tl_z_w(LBi:,LBj:,0:)
224 real(r8), intent(in) :: tl_Akv(LBi:,LBj:,0:)
225 real(r8), intent(in) :: tl_DU_avg1(LBi:,LBj:)
226 real(r8), intent(in) :: tl_DV_avg1(LBi:,LBj:)
227 real(r8), intent(in) :: tl_DU_avg2(LBi:,LBj:)
228 real(r8), intent(in) :: tl_DV_avg2(LBi:,LBj:)
229 real(r8), intent(in) :: tl_ru(LBi:,LBj:,0:,:)
230 real(r8), intent(in) :: tl_rv(LBi:,LBj:,0:,:)
231
232# ifdef DIAGNOSTICS_UV
233!! real(r8), intent(inout) :: DiaU2wrk(LBi:,LBj:,:)
234!! real(r8), intent(inout) :: DiaV2wrk(LBi:,LBj:,:)
235!! real(r8), intent(inout) :: DiaU2int(LBi:,LBj:,:)
236!! real(r8), intent(inout) :: DiaV2int(LBi:,LBj:,:)
237!! real(r8), intent(inout) :: DiaU3wrk(LBi:,LBj:,:,:)
238!! real(r8), intent(inout) :: DiaV3wrk(LBi:,LBj:,:,:)
239!! real(r8), intent(inout) :: DiaRU(LBi:,LBj:,:,:,:)
240!! real(r8), intent(inout) :: DiaRV(LBi:,LBj:,:,:,:)
241# endif
242 real(r8), intent(inout) :: Huon(LBi:,LBj:,:)
243 real(r8), intent(inout) :: Hvom(LBi:,LBj:,:)
244 real(r8), intent(inout) :: tl_u(LBi:,LBj:,:,:)
245 real(r8), intent(inout) :: tl_v(LBi:,LBj:,:,:)
246# ifdef NEARSHORE_MELLOR
247 real(r8), intent(inout) :: u_stokes(LBi:,LBj:,:)
248 real(r8), intent(inout) :: v_stokes(LBi:,LBj:,:)
249 real(r8), intent(inout) :: tl_u_stokes(LBi:,LBj:,:)
250 real(r8), intent(inout) :: tl_v_stokes(LBi:,LBj:,:)
251# endif
252 real(r8), intent(out) :: tl_ubar(LBi:,LBj:,:)
253 real(r8), intent(out) :: tl_vbar(LBi:,LBj:,:)
254 real(r8), intent(out) :: tl_Huon(LBi:,LBj:,:)
255 real(r8), intent(out) :: tl_Hvom(LBi:,LBj:,:)
256
257# else
258
259# ifdef MASKING
260 real(r8), intent(in) :: umask(LBi:UBi,LBj:UBj)
261 real(r8), intent(in) :: vmask(LBi:UBi,LBj:UBj)
262# endif
263# ifdef WET_DRY_NOT_YET
264 real(r8), intent(in) :: umask_wet(LBi:UBi,LBj:UBj)
265 real(r8), intent(in) :: vmask_wet(LBi:UBi,LBj:UBj)
266# endif
267 real(r8), intent(in) :: om_v(LBi:UBi,LBj:UBj)
268 real(r8), intent(in) :: on_u(LBi:UBi,LBj:UBj)
269 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
270 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
271 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
272 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
273 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
274 real(r8), intent(in) :: Akv(LBi:UBi,LBj:UBj,0:N(ng))
275 real(r8), intent(in) :: DU_avg1(LBi:UBi,LBj:UBj)
276 real(r8), intent(in) :: DV_avg1(LBi:UBi,LBj:UBj)
277 real(r8), intent(in) :: DU_avg2(LBi:UBi,LBj:UBj)
278 real(r8), intent(in) :: DV_avg2(LBi:UBi,LBj:UBj)
279 real(r8), intent(in) :: u(LBi:UBi,LBj:UBj,N(ng),2)
280 real(r8), intent(in) :: v(LBi:UBi,LBj:UBj,N(ng),2)
281# ifdef NEARSHORE_MELLOR
282 real(r8), intent(in) :: ubar_stokes(LBi:UBi,LBj:UBj)
283 real(r8), intent(in) :: vbar_stokes(LBi:UBi,LBj:UBj)
284 real(r8), intent(in) :: tl_ubar_stokes(LBi:UBi,LBj:UBj)
285 real(r8), intent(in) :: tl_vbar_stokes(LBi:UBi,LBj:UBj)
286# endif
287
288 real(r8), intent(in) :: tl_Hz(LBi:UBi,LBj:UBj,N(ng))
289 real(r8), intent(in) :: tl_z_r(LBi:UBi,LBj:UBj,N(ng))
290 real(r8), intent(in) :: tl_z_w(LBi:UBi,LBj:UBj,0:N(ng))
291 real(r8), intent(in) :: tl_Akv(LBi:UBi,LBj:UBj,0:N(ng))
292 real(r8), intent(in) :: tl_DU_avg1(LBi:UBi,LBj:UBj)
293 real(r8), intent(in) :: tl_DV_avg1(LBi:UBi,LBj:UBj)
294 real(r8), intent(in) :: tl_DU_avg2(LBi:UBi,LBj:UBj)
295 real(r8), intent(in) :: tl_DV_avg2(LBi:UBi,LBj:UBj)
296 real(r8), intent(in) :: tl_ru(LBi:UBi,LBj:UBj,0:N(ng),2)
297 real(r8), intent(in) :: tl_rv(LBi:UBi,LBj:UBj,0:N(ng),2)
298# ifdef DIAGNOSTICS_UV
299!! real(r8), intent(inout) :: DiaU2wrk(LBi:UBi,LBj:UBj,NDM2d)
300!! real(r8), intent(inout) :: DiaV2wrk(LBi:UBi,LBj:UBj,NDM2d)
301!! real(r8), intent(inout) :: DiaU2int(LBi:UBi,LBj:UBj,NDM2d)
302!! real(r8), intent(inout) :: DiaV2int(LBi:UBi,LBj:UBj,NDM2d)
303!! real(r8), intent(inout) :: DiaU3wrk(LBi:UBi,LBj:UBj,N(ng),NDM3d)
304!! real(r8), intent(inout) :: DiaV3wrk(LBi:UBi,LBj:UBj,N(ng),NDM3d)
305!! real(r8), intent(inout) :: DiaRU(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
306!! real(r8), intent(inout) :: DiaRV(LBi:UBi,LBj:UBj,N(ng),2,NDrhs)
307# endif
308 real(r8), intent(inout) :: Huon(LBi:UBi,LBj:UBj,N(ng))
309 real(r8), intent(inout) :: Hvom(LBi:UBi,LBj:UBj,N(ng))
310 real(r8), intent(inout) :: tl_u(LBi:UBi,LBj:UBj,N(ng),2)
311 real(r8), intent(inout) :: tl_v(LBi:UBi,LBj:UBj,N(ng),2)
312# ifdef NEARSHORE_MELLOR
313 real(r8), intent(inout) :: u_stokes(LBi:UBi,LBj:UBj,N(ng))
314 real(r8), intent(inout) :: v_stokes(LBi:UBi,LBj:UBj,N(ng))
315 real(r8), intent(inout) :: tl_u_stokes(LBi:UBi,LBj:UBj,N(ng))
316 real(r8), intent(inout) :: tl_v_stokes(LBi:UBi,LBj:UBj,N(ng))
317# endif
318 real(r8), intent(out) :: tl_ubar(LBi:UBi,LBj:UBj,3)
319 real(r8), intent(out) :: tl_vbar(LBi:UBi,LBj:UBj,3)
320 real(r8), intent(out) :: tl_Huon(LBi:UBi,LBj:UBj,N(ng))
321 real(r8), intent(out) :: tl_Hvom(LBi:UBi,LBj:UBj,N(ng))
322# endif
323!
324! Local variable declarations.
325!
326 integer :: i, idiag, is, j, k
327!
328 real(r8) :: cff, cff1, cff2
329 real(r8) :: tl_cff, tl_cff1, tl_cff2
330!
331 real(r8), dimension(IminS:ImaxS) :: CF1
332 real(r8), dimension(IminS:ImaxS) :: FC1
333# ifdef NEARSHORE_MELLOR
334 real(r8), dimension(IminS:ImaxS) :: CFs1
335 real(r8), dimension(IminS:ImaxS) :: DCs1
336# endif
337 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: AK
338 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: BC
339 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: CF
340 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: DC
341 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: DC1
342 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: FC
343# ifdef NEARSHORE_MELLOR
344 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: CFs
345 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: DCs
346# endif
347 real(r8), dimension(IminS:ImaxS,N(ng)) :: Hzk
348 real(r8), dimension(IminS:ImaxS,N(ng)) :: oHz
349# ifdef DIAGNOSTICS_UV
350!! real(r8), dimension(IminS:ImaxS,0:N(ng)) :: wrk
351!! real(r8), dimension(IminS:ImaxS,1:NDM2d-1) :: Dwrk
352# endif
353 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: tl_AK
354 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: tl_BC
355 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: tl_CF
356 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: tl_DC
357 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: tl_FC
358# ifdef NEARSHORE_MELLOR
359 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: tl_CFs
360 real(r8), dimension(IminS:ImaxS,0:N(ng)) :: tl_DCs
361# endif
362 real(r8), dimension(IminS:ImaxS,N(ng)) :: tl_Hzk
363 real(r8), dimension(IminS:ImaxS,N(ng)) :: tl_oHz
364!
365# include "set_bounds.h"
366!
367!-----------------------------------------------------------------------
368! Time step momentum equation in the XI-direction.
369!-----------------------------------------------------------------------
370!
371 DO j=jstr,jend
372 DO i=istru,iend
373 ak(i,0)=0.5_r8*(akv(i-1,j,0)+ &
374 & akv(i ,j,0))
375 tl_ak(i,0)=0.5_r8*(tl_akv(i-1,j,0)+ &
376 & tl_akv(i ,j,0))
377 DO k=1,n(ng)
378 ak(i,k)=0.5_r8*(akv(i-1,j,k)+ &
379 & akv(i ,j,k))
380 tl_ak(i,k)=0.5_r8*(tl_akv(i-1,j,k)+ &
381 & tl_akv(i ,j,k))
382 hzk(i,k)=0.5_r8*(hz(i-1,j,k)+ &
383 & hz(i ,j,k))
384 tl_hzk(i,k)=0.5_r8*(tl_hz(i-1,j,k)+ &
385 & tl_hz(i ,j,k))
386# if defined SPLINES_VVISC || defined DIAGNOSTICS_UV
387 ohz(i,k)=1.0_r8/hzk(i,k)
388 tl_ohz(i,k)=-ohz(i,k)*ohz(i,k)*tl_hzk(i,k)+ &
389# ifdef TL_IOMS
390 & 2.0_r8*ohz(i,k)
391# endif
392# endif
393 END DO
394 END DO
395!
396! Time step right-hand-side terms.
397!
398 IF (iic(ng).eq.ntfirst(ng)) THEN
399 cff=0.25_r8*dt(ng)
400 ELSE IF (iic(ng).eq.(ntfirst(ng)+1)) THEN
401 cff=0.25_r8*dt(ng)*3.0_r8/2.0_r8
402 ELSE
403 cff=0.25_r8*dt(ng)*23.0_r8/12.0_r8
404 END IF
405 DO i=istru,iend
406 dc(i,0)=cff*(pm(i,j)+pm(i-1,j))*(pn(i,j)+pn(i-1,j))
407 END DO
408!
409! The BASIC STATE "u" used below must be in transport units, but "u"
410! retrived is in m/s so we multiply by "Hzk".
411!
412 DO k=1,n(ng)
413 DO i=istru,iend
414!^ u(i,j,k,nnew)=u(i,j,k,nnew)+ &
415!^ & DC(i,0)*ru(i,j,k,nrhs)
416!^
417 tl_u(i,j,k,nnew)=tl_u(i,j,k,nnew)+ &
418 & dc(i,0)*tl_ru(i,j,k,nrhs)
419# ifdef SPLINES_VVISC
420!^ u(i,j,k,nnew)=u(i,j,k,nnew)*oHz(i,k)
421!^
422 tl_u(i,j,k,nnew)=tl_u(i,j,k,nnew)*ohz(i,k)+ &
423 & (u(i,j,k,nnew)*hzk(i,k))*tl_ohz(i,k)- &
424# ifdef TL_IOMS
425 & u(i,j,k,nnew)*hzk(i,k)*ohz(i,k)
426# endif
427# endif
428# ifdef DIAGNOSTICS_UV
429!! DO idiag=1,M3pgrd
430!! DiaU3wrk(i,j,k,idiag)=(DiaU3wrk(i,j,k,idiag)+ &
431!! & DC(i,0)*DiaRU(i,j,k,nrhs,idiag))* &
432!! & oHz(i,k)
433!! END DO
434# if defined UV_VIS2 || defined UV_VIS4
435!! DiaU3wrk(i,j,k,M3xvis)=DiaU3wrk(i,j,k,M3xvis)*oHz(i,k)
436!! DiaU3wrk(i,j,k,M3yvis)=DiaU3wrk(i,j,k,M3yvis)*oHz(i,k)
437!! DiaU3wrk(i,j,k,M3hvis)=DiaU3wrk(i,j,k,M3hvis)*oHz(i,k)
438# endif
439!! DiaU3wrk(i,j,k,M3vvis)=DiaU3wrk(i,j,k,M3vvis)*oHz(i,k)
440# ifdef BODYFORCE
441!! DiaU3wrk(i,j,k,M3vvis)=DiaU3wrk(i,j,k,M3vvis)+ &
442!! & DC(i,0)*DiaRU(i,j,k,nrhs,M3vvis)* &
443!! & oHz(i,k)
444# endif
445!! DiaU3wrk(i,j,k,M3rate)=DiaU3wrk(i,j,k,M3rate)*oHz(i,k)
446# endif
447 END DO
448 END DO
449
450# ifdef SPLINES_VVISC
451!
452! Apply spline algorithim to BASIC STATE "u" which should be
453! in units of velocity (m/s). DC will be used in the tangent
454! linear spline algorithm below. Solve BASIC STATE tridiagonal
455! system.
456!
457 cff1=1.0_r8/6.0_r8
458 DO k=1,n(ng)-1
459 DO i=istru,iend
460 fc(i,k)=cff1*hzk(i,k )-dt(ng)*ak(i,k-1)*ohz(i,k )
461 cf(i,k)=cff1*hzk(i,k+1)-dt(ng)*ak(i,k+1)*ohz(i,k+1)
462 END DO
463 END DO
464 DO i=istru,iend
465 cf(i,0)=0.0_r8
466 dc(i,0)=0.0_r8
467 END DO
468!
469! LU decomposition and forward substitution.
470!
471 cff1=1.0_r8/3.0_r8
472 DO k=1,n(ng)-1
473 DO i=istru,iend
474 bc(i,k)=cff1*(hzk(i,k)+hzk(i,k+1))+ &
475 & dt(ng)*ak(i,k)*(ohz(i,k)+ohz(i,k+1))
476 cff=1.0_r8/(bc(i,k)-fc(i,k)*cf(i,k-1))
477 cf(i,k)=cff*cf(i,k)
478 dc(i,k)=cff*(u(i,j,k+1,nnew)-u(i,j,k,nnew)- &
479 & fc(i,k)*dc(i,k-1))
480 END DO
481 END DO
482!
483! Backward substitution. Save DC for the tangent linearization.
484! DC is scaled later by AK.
485!
486 DO i=istru,iend
487 dc(i,n(ng))=0.0_r8
488 END DO
489 DO k=n(ng)-1,1,-1
490 DO i=istru,iend
491 dc(i,k)=dc(i,k)-cf(i,k)*dc(i,k+1)
492 END DO
493 END DO
494!
495! Use conservative, parabolic spline reconstruction of tangent
496! linear vertical viscosity derivatives. Then, time step vertical
497! viscosity term implicitly by solving a tridiagonal system.
498# ifdef TL_IOMS
499! If tl_AK is computed, then tl_FC=tl_FC+dt*AK*oHz and
500! tl_CF=tl_CF+dt*AK*oHz below.
501# endif
502!
503 cff1=1.0_r8/6.0_r8
504 DO k=1,n(ng)-1
505 DO i=istru,iend
506 fc(i,k)=cff1*hzk(i,k )- &
507 & dt(ng)*ak(i,k-1)*ohz(i,k )
508 tl_fc(i,k)=cff1*tl_hzk(i,k )- &
509 & dt(ng)*(tl_ak(i,k-1)*ohz(i,k )+ &
510 & ak(i,k-1)*tl_ohz(i,k ))
511 cf(i,k)=cff1*hzk(i,k+1)- &
512 & dt(ng)*ak(i,k+1)*ohz(i,k+1)
513 tl_cf(i,k)=cff1*tl_hzk(i,k+1)- &
514 & dt(ng)*(tl_ak(i,k+1)*ohz(i,k+1)+ &
515 & ak(i,k+1)*tl_ohz(i,k+1))
516 END DO
517 END DO
518 DO i=istru,iend
519 cf(i,0)=0.0_r8
520 tl_cf(i,0)=0.0_r8
521 tl_dc(i,0)=0.0_r8
522 END DO
523!
524! Tangent linear LU decomposition and forward substitution.
525# ifdef TL_IOMS
526! If tl_AK is computed, then tl_BC=tl_FC-dt*AK*oHz below.
527# endif
528!
529 cff1=1.0_r8/3.0_r8
530 DO k=1,n(ng)-1
531 DO i=istru,iend
532 bc(i,k)=cff1*(hzk(i,k)+hzk(i,k+1))+ &
533 & dt(ng)*ak(i,k)*(ohz(i,k)+ohz(i,k+1))
534 tl_bc(i,k)=cff1*(tl_hzk(i,k)+tl_hzk(i,k+1))+ &
535 & dt(ng)*(tl_ak(i,k)*(ohz(i,k)+ohz(i,k+1))+ &
536 & ak(i,k)*(tl_ohz(i,k)+tl_ohz(i,k+1)))
537 cff=1.0_r8/(bc(i,k)-fc(i,k)*cf(i,k-1))
538# ifdef TL_IOMS
539 tl_dc(i,k)=cff*(tl_u(i,j,k+1,nnew)- &
540 & tl_u(i,j,k ,nnew)- &
541 & ((tl_fc(i,k)-fc(i,k))*dc(i,k-1)+ &
542 & (tl_bc(i,k)-bc(i,k))*dc(i,k )+ &
543 & (tl_cf(i,k)-cf(i,k))*dc(i,k+1))- &
544 & fc(i,k)*tl_dc(i,k-1))
545 cf(i,k)=cff*cf(i,k)
546# else
547 cf(i,k)=cff*cf(i,k)
548 tl_dc(i,k)=cff*(tl_u(i,j,k+1,nnew)- &
549 & tl_u(i,j,k ,nnew)- &
550 & (tl_fc(i,k)*dc(i,k-1)+ &
551 & tl_bc(i,k)*dc(i,k )+ &
552 & tl_cf(i,k)*dc(i,k+1))- &
553 & fc(i,k)*tl_dc(i,k-1))
554# endif
555 END DO
556 END DO
557!
558! Tangent linear backward substitution.
559!
560 DO i=istru,iend
561 tl_dc(i,n(ng))=0.0_r8
562 END DO
563 DO k=n(ng)-1,1,-1
564 DO i=istru,iend
565 tl_dc(i,k)=tl_dc(i,k)-cf(i,k)*tl_dc(i,k+1)
566 END DO
567 END DO
568!
569! Compute tl_DC before multiplying BASIC STATE spline gradients
570! DC by AK.
571!
572 DO k=1,n(ng)
573 DO i=istru,iend
574 tl_dc(i,k)=tl_dc(i,k)*ak(i,k)+dc(i,k)*tl_ak(i,k)
575# ifdef TL_IOMS
576! tl_DC(i,k)=tl_DC(i,k)- &
577! & DC(i,k)*AK(i,k) ! needed if tl_AK is computed
578# endif
579 dc(i,k)=dc(i,k)*ak(i,k)
580!^ cff=dt(ng)*oHz(i,k)*(DC(i,k)-DC(i,k-1))
581!^
582 tl_cff=dt(ng)*(tl_ohz(i,k)*(dc(i,k)-dc(i,k-1))+ &
583 & ohz(i,k)*(tl_dc(i,k)-tl_dc(i,k-1)))- &
584# ifdef TL_IOMS
585 & dt(ng)*ohz(i,k)*(dc(i,k)-dc(i,k-1))
586# endif
587!^ u(i,j,k,nnew)=u(i,j,k,nnew)+cff
588!^
589 tl_u(i,j,k,nnew)=tl_u(i,j,k,nnew)+tl_cff
590# ifdef DIAGNOSTICS_UV
591!! DiaU3wrk(i,j,k,M3vvis)=DiaU3wrk(i,j,k,M3vvis)+cff
592# endif
593 END DO
594 END DO
595# else
596!
597! Compute off-diagonal coefficients [lambda*dt*Akv/Hz] for the
598! implicit vertical viscosity term at horizontal U-points and
599! vertical W-points.
600!
601! NOTE: The original code solves the tridiagonal system A*u=r where
602! A is a matrix and u and r are vectors. We need to solve the
603! tangent linear of this system which is A*tl_u+tl_A*u=tl_r.
604! Here, tl_A*u and tl_r are known, so we must solve for tl_u
605! by inverting A*tl_u=tl_r-tl_A*u.
606!
607 cff=-lambda*dt(ng)/0.5_r8
608 DO k=1,n(ng)-1
609 DO i=istru,iend
610 cff1=1.0_r8/(z_r(i,j,k+1)+z_r(i-1,j,k+1)- &
611 & z_r(i,j,k )-z_r(i-1,j,k ))
612 tl_cff1=-cff1*cff1*(tl_z_r(i,j,k+1)+tl_z_r(i-1,j,k+1)- &
613 & tl_z_r(i,j,k )-tl_z_r(i-1,j,k ))+ &
614# ifdef TL_IOMS
615 & 2.0_r8*cff1
616# endif
617 fc(i,k)=cff*cff1*ak(i,k)
618 tl_fc(i,k)=cff*(tl_cff1*ak(i,k)+cff1*tl_ak(i,k))
619# ifdef TL_IOMS
620!
621! Uncomment if tl_Akv is computed.
622!
623! tl_FC(i,k)=tl_FC(i,k)-FC(i,k)
624# endif
625 END DO
626 END DO
627 DO i=istru,iend
628 fc(i,0)=0.0_r8
629 tl_fc(i,0)=0.0_r8
630 fc(i,n(ng))=0.0_r8
631 tl_fc(i,n(ng))=0.0_r8
632 END DO
633!
634! Solve the tangent linear tridiagonal system.
635! (DC is a tangent linear variable here).
636!
637 DO k=1,n(ng)
638 DO i=istru,iend
639 bc(i,k)=hzk(i,k)-fc(i,k)-fc(i,k-1)
640 tl_bc(i,k)=tl_hzk(i,k)-tl_fc(i,k)-tl_fc(i,k-1)
641 END DO
642 END DO
643 DO k=2,n(ng)-1
644 DO i=istru,iend
645# ifdef TL_IOMS
646 dc(i,k)=tl_u(i,j,k,nnew)- &
647 & ((tl_fc(i,k-1)-fc(i,k-1))*u(i,j,k-1,nnew)+ &
648 & (tl_bc(i,k )-bc(i,k ))*u(i,j,k ,nnew)+ &
649 & (tl_fc(i,k )-fc(i,k ))*u(i,j,k+1,nnew))
650# else
651 dc(i,k)=tl_u(i,j,k,nnew)- &
652 & (tl_fc(i,k-1)*u(i,j,k-1,nnew)+ &
653 & tl_bc(i,k )*u(i,j,k ,nnew)+ &
654 & tl_fc(i,k )*u(i,j,k+1,nnew))
655# endif
656 END DO
657 END DO
658 DO i=istru,iend
659# ifdef TL_IOMS
660 dc(i,1)=tl_u(i,j,1,nnew)- &
661 & ((tl_bc(i,1)-bc(i,1))*u(i,j,1,nnew)+ &
662 & (tl_fc(i,1)-fc(i,1))*u(i,j,2,nnew))
663 dc(i,n(ng))=tl_u(i,j,n(ng),nnew)- &
664 & ((tl_fc(i,n(ng)-1)-fc(i,n(ng)-1))* &
665 & u(i,j,n(ng)-1,nnew)+ &
666 & (tl_bc(i,n(ng) )-bc(i,n(ng) ))* &
667 & u(i,j,n(ng) ,nnew))
668# else
669 dc(i,1)=tl_u(i,j,1,nnew)- &
670 & (tl_bc(i,1)*u(i,j,1,nnew)+ &
671 & tl_fc(i,1)*u(i,j,2,nnew))
672 dc(i,n(ng))=tl_u(i,j,n(ng),nnew)- &
673 & (tl_fc(i,n(ng)-1)*u(i,j,n(ng)-1,nnew)+ &
674 & tl_bc(i,n(ng) )*u(i,j,n(ng) ,nnew))
675# endif
676 END DO
677 DO i=istru,iend
678 cff=1.0_r8/bc(i,1)
679 cf(i,1)=cff*fc(i,1)
680 dc(i,1)=cff*dc(i,1)
681 END DO
682 DO k=2,n(ng)-1
683 DO i=istru,iend
684 cff=1.0_r8/(bc(i,k)-fc(i,k-1)*cf(i,k-1))
685 cf(i,k)=cff*fc(i,k)
686 dc(i,k)=cff*(dc(i,k)-fc(i,k-1)*dc(i,k-1))
687 END DO
688 END DO
689!
690! Compute new solution by back substitution.
691! (DC is a tangent linear variable here).
692!
693 DO i=istru,iend
694# ifdef DIAGNOSTICS_UV
695!! wrk(i,N(ng))=u(i,j,N(ng),nnew)*oHz(i,N(ng))
696# endif
697 dc(i,n(ng))=(dc(i,n(ng))-fc(i,n(ng)-1)*dc(i,n(ng)-1))/ &
698 & (bc(i,n(ng))-fc(i,n(ng)-1)*cf(i,n(ng)-1))
699!^ u(i,j,N(ng),nnew)=DC(i,N(ng))
700!^
701 tl_u(i,j,n(ng),nnew)=dc(i,n(ng))
702# ifdef DIAGNOSTICS_UV
703!! DiaU3wrk(i,j,N(ng),M3vvis)=DiaU3wrk(i,j,N(ng),M3vvis)+ &
704!! & u(i,j,N(ng),nnew)-wrk(i,N(ng))
705# endif
706 END DO
707 DO k=n(ng)-1,1,-1
708 DO i=istru,iend
709# ifdef DIAGNOSTICS_UV
710!! wrk(i,k)=u(i,j,k,nnew)*oHz(i,k)
711# endif
712 dc(i,k)=dc(i,k)-cf(i,k)*dc(i,k+1)
713!^ u(i,j,k,nnew)=DC(i,k)
714!^
715 tl_u(i,j,k,nnew)=dc(i,k)
716# ifdef DIAGNOSTICS_UV
717!! DiaU3wrk(i,j,k,M3vvis)=DiaU3wrk(i,j,k,M3vvis)+ &
718!! & u(i,j,k,nnew)-wrk(i,k)
719# endif
720 END DO
721 END DO
722# endif
723!
724! Replace INTERIOR POINTS incorrect vertical mean with more accurate
725! barotropic component, ubar=DU_avg1/(D*on_u). Recall that, D=CF(:,0).
726!
727 DO i=istru,iend
728 cf(i,0)=hzk(i,1)
729 tl_cf(i,0)=tl_hzk(i,1)
730 dc(i,0)=u(i,j,1,nnew)*hzk(i,1)
731 tl_dc(i,0)=tl_u(i,j,1,nnew)*hzk(i,1)+ &
732 & u(i,j,1,nnew)*tl_hzk(i,1)- &
733# ifdef TL_IOMS
734 & dc(i,0)
735# endif
736# ifdef NEARSHORE_MELLOR
737 dcs(i,0)=u_stokes(i,j,1)*hzk(i,1)
738 tl_dcs(i,0)=tl_u_stokes(i,j,1)*hzk(i,1)+ &
739 & u_stokes(i,j,1)*tl_hzk(i,1)- &
740# ifdef TL_IOMS
741 & dcs(i,0)
742# endif
743# endif
744# ifdef DIAGNOSTICS_UV
745!! Dwrk(i,M2pgrd)=DiaU3wrk(i,j,1,M3pgrd)*Hzk(i,1)
746!! Dwrk(i,M2bstr)=DiaU3wrk(i,j,1,M3vvis)*Hzk(i,1)
747# ifdef UV_COR
748!! Dwrk(i,M2fcor)=DiaU3wrk(i,j,1,M3fcor)*Hzk(i,1)
749# endif
750# if defined UV_VIS2 || defined UV_VIS4
751!! Dwrk(i,M2xvis)=DiaU3wrk(i,j,1,M3xvis)*Hzk(i,1)
752!! Dwrk(i,M2yvis)=DiaU3wrk(i,j,1,M3yvis)*Hzk(i,1)
753!! Dwrk(i,M2hvis)=DiaU3wrk(i,j,1,M3hvis)*Hzk(i,1)
754# endif
755# ifdef UV_ADV
756!! Dwrk(i,M2xadv)=DiaU3wrk(i,j,1,M3xadv)*Hzk(i,1)
757!! Dwrk(i,M2yadv)=DiaU3wrk(i,j,1,M3yadv)*Hzk(i,1)
758!! Dwrk(i,M2hadv)=DiaU3wrk(i,j,1,M3hadv)*Hzk(i,1)
759# endif
760# ifdef NEARSHORE_MELLOR
761!! Dwrk(i,M2hrad)=DiaU3wrk(i,j,1,M3hrad)*Hzk(i,1)
762# endif
763# endif
764 END DO
765 DO k=2,n(ng)
766 DO i=istru,iend
767 cf(i,0)=cf(i,0)+hzk(i,k)
768 tl_cf(i,0)=tl_cf(i,0)+tl_hzk(i,k)
769 dc(i,0)=dc(i,0)+u(i,j,k,nnew)*hzk(i,k)
770 tl_dc(i,0)=tl_dc(i,0)+ &
771 & tl_u(i,j,k,nnew)*hzk(i,k)+ &
772 & u(i,j,k,nnew)*tl_hzk(i,k)- &
773# ifdef TL_IOMS
774 & u(i,j,k,nnew)*hzk(i,k)
775# endif
776# ifdef NEARSHORE_MELLOR
777 dcs(i,0)=dcs(i,0)+u_stokes(i,j,k)*hzk(i,k)
778 tl_dcs(i,0)=tl_dcs(i,0)+ &
779 & tl_u_stokes(i,j,k)*hzk(i,k)+ &
780 & u_stokes(i,j,k)*tl_hzk(i,k)- &
781# ifdef TL_IOMS
782 & u_stokes(i,j,k)*hzk(i,k)
783# endif
784# endif
785# ifdef DIAGNOSTICS_UV
786!! Dwrk(i,M2pgrd)=Dwrk(i,M2pgrd)+ &
787!! & DiaU3wrk(i,j,k,M3pgrd)*Hzk(i,k)
788!! Dwrk(i,M2bstr)=Dwrk(i,M2bstr)+ &
789!! & DiaU3wrk(i,j,k,M3vvis)*Hzk(i,k)
790# ifdef UV_COR
791!! Dwrk(i,M2fcor)=Dwrk(i,M2fcor)+ &
792!! & DiaU3wrk(i,j,k,M3fcor)*Hzk(i,k)
793# endif
794# if defined UV_VIS2 || defined UV_VIS4
795!! Dwrk(i,M2xvis)=Dwrk(i,M2xvis)+ &
796!! & DiaU3wrk(i,j,k,M3xvis)*Hzk(i,k)
797!! Dwrk(i,M2yvis)=Dwrk(i,M2yvis)+ &
798!! & DiaU3wrk(i,j,k,M3yvis)*Hzk(i,k)
799!! Dwrk(i,M2hvis)=Dwrk(i,M2hvis)+ &
800!! & DiaU3wrk(i,j,k,M3hvis)*Hzk(i,k)
801# endif
802# ifdef UV_ADV
803!! Dwrk(i,M2xadv)=Dwrk(i,M2xadv)+ &
804!! & DiaU3wrk(i,j,k,M3xadv)*Hzk(i,k)
805!! Dwrk(i,M2yadv)=Dwrk(i,M2yadv)+ &
806!! & DiaU3wrk(i,j,k,M3yadv)*Hzk(i,k)
807!! Dwrk(i,M2hadv)=Dwrk(i,M2hadv)+ &
808!! & DiaU3wrk(i,j,k,M3hadv)*Hzk(i,k)
809# endif
810# ifdef NEARSHORE_MELLOR
811!! Dwrk(i,M2hrad)=Dwrk(i,M2hrad)+ &
812!! & DiaU3wrk(i,j,k,M3hrad)*Hzk(i,k)
813# endif
814# endif
815 END DO
816 END DO
817 DO i=istru,iend
818 dc1(i,0)=dc(i,0) ! intermediate
819 cff1=1.0_r8/(cf(i,0)*on_u(i,j))
820 tl_cff1=-cff1*cff1*tl_cf(i,0)*on_u(i,j)+ &
821# ifdef TL_IOMS
822 & 2.0_r8*cff1
823# endif
824 dc(i,0)=(dc(i,0)*on_u(i,j)-du_avg1(i,j))*cff1 ! recursive
825 tl_dc(i,0)=(tl_dc(i,0)*on_u(i,j)-tl_du_avg1(i,j))*cff1+ &
826 & (dc1(i,0)*on_u(i,j)-du_avg1(i,j))*tl_cff1- &
827# ifdef TL_IOMS
828 & dc(i,0)
829# endif
830# ifdef NEARSHORE_MELLOR
831 dcs1(i)=dcs(i,0) ! intermediate
832 cff2=1.0_r8/cf(i,0)
833 tl_cff2=-cff2*cff2*tl_cf(i,0)+ &
834# ifdef TL_IOMS
835 & 2.0_r8*cff2
836# endif
837 dcs(i,0)=dcs(i,0)*cff2-ubar_stokes(i,j) ! recursive
838 tl_dcs(i,0)=tl_dcs(i,0)*cff2+ &
839 & dcs1(i)*tl_cff2-tl_ubar_stokes(i,j)- &
840# ifdef TL_IOMS
841 & dcs(i,0)
842# endif
843# endif
844# ifdef DIAGNOSTICS_UV
845!! DO idiag=1,M2pgrd
846!! Dwrk(i,idiag)=(Dwrk(i,idiag)*on_u(i,j)- &
847!! & DiaU2wrk(i,j,idiag))*cff1
848!! END DO
849!! Dwrk(i,M2bstr)=(Dwrk(i,M2bstr)*on_u(i,j)- &
850!! & DiaU2wrk(i,j,M2bstr)-DiaU2wrk(i,j,M2sstr))* &
851!! & cff1
852# endif
853 END DO
854!
855! Couple and update new solution.
856!
857 DO k=1,n(ng)
858 DO i=istru,iend
859!^ u(i,j,k,nnew)=u(i,j,k,nnew)-DC(i,0)
860!^
861 tl_u(i,j,k,nnew)=tl_u(i,j,k,nnew)-tl_dc(i,0)
862# ifdef MASKING
863!^ u(i,j,k,nnew)=u(i,j,k,nnew)*umask(i,j)
864!^
865 tl_u(i,j,k,nnew)=tl_u(i,j,k,nnew)*umask(i,j)
866# endif
867# ifdef WET_DRY_NOT_YET
868!^ u(i,j,k,nnew)=u(i,j,k,nnew)*umask_wet(i,j)
869!^
870 tl_u(i,j,k,nnew)=tl_u(i,j,k,nnew)*umask_wet(i,j)
871# endif
872# ifdef NEARSHORE_MELLOR
873!^ u_stokes(i,j,k)=u_stokes(i,j,k)-DCs(i,0)
874!^
875 tl_u_stokes(i,j,k)=tl_u_stokes(i,j,k)-tl_dcs(i,0)
876# ifdef MASKING
877!^ u_stokes(i,j,k)=u_stokes(i,j,k)*umask(i,j)
878!^
879 tl_u_stokes(i,j,k)=tl_u_stokes(i,j,k)*umask(i,j)
880# endif
881# ifdef WET_DRY_NOT_YET
882!^ u_stokes(i,j,k)=u_stokes(i,j,k)*umask_wet(i,j)
883!^
884 tl_u_stokes(i,j,k)=tl_u_stokes(i,j,k)*umask_wet(i,j)
885# endif
886# endif
887# ifdef DIAGNOSTICS_UV
888!! DiaU3wrk(i,j,k,M3pgrd)=DiaU3wrk(i,j,k,M3pgrd)- &
889!! & Dwrk(i,M2pgrd)
890!! DiaU3wrk(i,j,k,M3vvis)=DiaU3wrk(i,j,k,M3vvis)- &
891!! & Dwrk(i,M2bstr)
892# ifdef UV_COR
893!! DiaU3wrk(i,j,k,M3fcor)=DiaU3wrk(i,j,k,M3fcor)- &
894!! & Dwrk(i,M2fcor)
895# endif
896# if defined UV_VIS2 || defined UV_VIS4
897!! DiaU3wrk(i,j,k,M3xvis)=DiaU3wrk(i,j,k,M3xvis)- &
898!! & Dwrk(i,M2xvis)
899!! DiaU3wrk(i,j,k,M3yvis)=DiaU3wrk(i,j,k,M3yvis)- &
900!! & Dwrk(i,M2yvis)
901!! DiaU3wrk(i,j,k,M3hvis)=DiaU3wrk(i,j,k,M3hvis)- &
902!! & Dwrk(i,M2hvis)
903# endif
904# ifdef UV_ADV
905!! DiaU3wrk(i,j,k,M3xadv)=DiaU3wrk(i,j,k,M3xadv)- &
906!! & Dwrk(i,M2xadv)
907!! DiaU3wrk(i,j,k,M3yadv)=DiaU3wrk(i,j,k,M3yadv)- &
908!! & Dwrk(i,M2yadv)
909!! DiaU3wrk(i,j,k,M3hadv)=DiaU3wrk(i,j,k,M3hadv)- &
910!! & Dwrk(i,M2hadv)
911# endif
912# ifdef NEARSHORE_MELLOR
913!! DiaU3wrk(i,j,k,M3hrad)=DiaU3wrk(i,j,k,M3hrad)- &
914!! & Dwrk(i,M2hrad)
915# endif
916# endif
917 END DO
918 END DO
919
920# if defined DIAGNOSTICS_UV && defined MASKING
921!! DO k=1,N(ng)
922!! DO i=IstrU,Iend
923!! DO idiag=1,NDM3d
924!! DiaU3wrk(i,j,k,idiag)=DiaU3wrk(i,j,k,idiag)*umask(i,j)
925!! END DO
926!! END DO
927!! END DO
928# endif
929!
930!-----------------------------------------------------------------------
931! Time step momentum equation in the ETA-direction.
932!-----------------------------------------------------------------------
933!
934 IF (j.ge.jstrv) THEN
935 DO i=istr,iend
936 ak(i,0)=0.5_r8*(akv(i,j-1,0)+ &
937 & akv(i,j ,0))
938 tl_ak(i,0)=0.5_r8*(tl_akv(i,j-1,0)+ &
939 & tl_akv(i,j ,0))
940 DO k=1,n(ng)
941 ak(i,k)=0.5_r8*(akv(i,j-1,k)+ &
942 & akv(i,j ,k))
943 tl_ak(i,k)=0.5_r8*(tl_akv(i,j-1,k)+ &
944 & tl_akv(i,j ,k))
945 hzk(i,k)=0.5_r8*(hz(i,j-1,k)+ &
946 & hz(i,j ,k))
947 tl_hzk(i,k)=0.5_r8*(tl_hz(i,j-1,k)+ &
948 & tl_hz(i,j ,k))
949# if defined SPLINES_VVISC || defined DIAGNOSTICS_UV
950 ohz(i,k)=1.0_r8/hzk(i,k)
951 tl_ohz(i,k)=-ohz(i,k)*ohz(i,k)*tl_hzk(i,k)+ &
952# ifdef TL_IOMS
953 & 2.0_r8*ohz(i,k)
954# endif
955# endif
956 END DO
957 END DO
958!
959! Time step right-hand-side terms.
960!
961 IF (iic(ng).eq.ntfirst(ng)) THEN
962 cff=0.25_r8*dt(ng)
963 ELSE IF (iic(ng).eq.(ntfirst(ng)+1)) THEN
964 cff=0.25_r8*dt(ng)*3.0_r8/2.0_r8
965 ELSE
966 cff=0.25_r8*dt(ng)*23.0_r8/12.0_r8
967 END IF
968 DO i=istr,iend
969 dc(i,0)=cff*(pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))
970 END DO
971!
972! The BASIC STATE "v" used below must be in transport units, but "v"
973! retrived is in m/s so we multiply by "Hzk".
974!
975 DO k=1,n(ng)
976 DO i=istr,iend
977!^ v(i,j,k,nnew)=v(i,j,k,nnew)+DC(i,0)*rv(i,j,k,nrhs)
978!^
979 tl_v(i,j,k,nnew)=tl_v(i,j,k,nnew)+ &
980 & dc(i,0)*tl_rv(i,j,k,nrhs)
981# ifdef SPLINES_VVISC
982!^ v(i,j,k,nnew)=v(i,j,k,nnew)*oHz(i,k)
983!^
984 tl_v(i,j,k,nnew)=tl_v(i,j,k,nnew)*ohz(i,k)+ &
985 & (v(i,j,k,nnew)*hzk(i,k))*tl_ohz(i,k)- &
986# ifdef TL_IOMS
987 & v(i,j,k,nnew)*hzk(i,k)*ohz(i,k)
988# endif
989# endif
990# ifdef DIAGNOSTICS_UV
991!! DO idiag=1,M3pgrd
992!! DiaV3wrk(i,j,k,idiag)=(DiaV3wrk(i,j,k,idiag)+ &
993!! & DC(i,0)* &
994!! & DiaRV(i,j,k,nrhs,idiag))* &
995!! & oHz(i,k)
996!! END DO
997# if defined UV_VIS2 || defined UV_VIS4
998!! DiaV3wrk(i,j,k,M3xvis)=DiaV3wrk(i,j,k,M3xvis)*oHz(i,k)
999!! DiaV3wrk(i,j,k,M3yvis)=DiaV3wrk(i,j,k,M3yvis)*oHz(i,k)
1000!! DiaV3wrk(i,j,k,M3hvis)=DiaV3wrk(i,j,k,M3hvis)*oHz(i,k)
1001# endif
1002!! DiaV3wrk(i,j,k,M3vvis)=DiaV3wrk(i,j,k,M3vvis)*oHz(i,k)
1003# ifdef BODYFORCE
1004!! DiaV3wrk(i,j,k,M3vvis)=DiaV3wrk(i,j,k,M3vvis)+ &
1005!! & DC(i,0)*DiaRV(i,j,k,nrhs,M3vvis)* &
1006!! & oHz(i,k)
1007# endif
1008!! DiaV3wrk(i,j,k,M3rate)=DiaV3wrk(i,j,k,M3rate)*oHz(i,k)
1009# endif
1010 END DO
1011 END DO
1012
1013# ifdef SPLINES_VVISC
1014!
1015! Apply spline algorithim to BASIC STATE "v" which should be
1016! in units of velocity (m/s). DC will be used in the tangent
1017! linear spline algorithm below. Solve BASIC STATE tridiagonal
1018! system.
1019!
1020 cff1=1.0_r8/6.0_r8
1021 DO k=1,n(ng)-1
1022 DO i=istr,iend
1023 fc(i,k)=cff1*hzk(i,k )-dt(ng)*ak(i,k-1)*ohz(i,k )
1024 cf(i,k)=cff1*hzk(i,k+1)-dt(ng)*ak(i,k+1)*ohz(i,k+1)
1025 END DO
1026 END DO
1027 DO i=istr,iend
1028 cf(i,0)=0.0_r8
1029 dc(i,0)=0.0_r8
1030 END DO
1031!
1032! LU decomposition and forward substitution.
1033!
1034 cff1=1.0_r8/3.0_r8
1035 DO k=1,n(ng)-1
1036 DO i=istr,iend
1037 bc(i,k)=cff1*(hzk(i,k)+hzk(i,k+1))+ &
1038 & dt(ng)*ak(i,k)*(ohz(i,k)+ohz(i,k+1))
1039 cff=1.0_r8/(bc(i,k)-fc(i,k)*cf(i,k-1))
1040 cf(i,k)=cff*cf(i,k)
1041 dc(i,k)=cff*(v(i,j,k+1,nnew)-v(i,j,k,nnew)- &
1042 & fc(i,k)*dc(i,k-1))
1043 END DO
1044 END DO
1045!
1046! Backward substitution. Save DC for the tangent linearization.
1047! DC is scaled later by AK.
1048!
1049 DO i=istr,iend
1050 dc(i,n(ng))=0.0_r8
1051 END DO
1052 DO k=n(ng)-1,1,-1
1053 DO i=istr,iend
1054 dc(i,k)=dc(i,k)-cf(i,k)*dc(i,k+1)
1055 END DO
1056 END DO
1057!
1058! Use conservative, parabolic spline reconstruction of tangent
1059! linear vertical viscosity derivatives. Then, time step vertical
1060! viscosity term implicitly by solving a tridiagonal system.
1061# ifdef TL_IOMS
1062! If tl_AK is computed, then tl_FC=tl_FC+dt*AK*oHz and
1063! tl_CF=tl_CF+dt*AK*oHz below.
1064# endif
1065!
1066 cff1=1.0_r8/6.0_r8
1067 DO k=1,n(ng)-1
1068 DO i=istr,iend
1069 fc(i,k)=cff1*hzk(i,k )- &
1070 & dt(ng)*ak(i,k-1)*ohz(i,k )
1071 tl_fc(i,k)=cff1*tl_hzk(i,k )- &
1072 & dt(ng)*(tl_ak(i,k-1)*ohz(i,k )+ &
1073 & ak(i,k-1)*tl_ohz(i,k ))
1074 cf(i,k)=cff1*hzk(i,k+1)- &
1075 & dt(ng)*ak(i,k+1)*ohz(i,k+1)
1076 tl_cf(i,k)=cff1*tl_hzk(i,k+1)- &
1077 & dt(ng)*(tl_ak(i,k+1)*ohz(i,k+1)+ &
1078 & ak(i,k+1)*tl_ohz(i,k+1))
1079 END DO
1080 END DO
1081 DO i=istr,iend
1082 cf(i,0)=0.0_r8
1083 tl_cf(i,0)=0.0_r8
1084 tl_dc(i,0)=0.0_r8
1085 END DO
1086!
1087! Tangent linear LU decomposition and forward substitution.
1088# ifdef TL_IOMS
1089! If tl_AK is computed, then tl_BC=tl_BC-dt*AK*oHz below.
1090# endif
1091!
1092 cff1=1.0_r8/3.0_r8
1093 DO k=1,n(ng)-1
1094 DO i=istr,iend
1095 bc(i,k)=cff1*(hzk(i,k)+hzk(i,k+1))+ &
1096 & dt(ng)*ak(i,k)*(ohz(i,k)+ohz(i,k+1))
1097 tl_bc(i,k)=cff1*(tl_hzk(i,k)+tl_hzk(i,k+1))+ &
1098 & dt(ng)*(tl_ak(i,k)*(ohz(i,k)+ohz(i,k+1))+ &
1099 & ak(i,k)*(tl_ohz(i,k)+tl_ohz(i,k+1)))
1100 cff=1.0_r8/(bc(i,k)-fc(i,k)*cf(i,k-1))
1101# ifdef TL_IOMS
1102 tl_dc(i,k)=cff*(tl_v(i,j,k+1,nnew)- &
1103 & tl_v(i,j,k ,nnew)- &
1104 & ((tl_fc(i,k)-fc(i,k))*dc(i,k-1)+ &
1105 & (tl_bc(i,k)-bc(i,k))*dc(i,k )+ &
1106 & (tl_cf(i,k)-cf(i,k))*dc(i,k+1))- &
1107 & fc(i,k)*tl_dc(i,k-1))
1108 cf(i,k)=cff*cf(i,k)
1109# else
1110 cf(i,k)=cff*cf(i,k)
1111 tl_dc(i,k)=cff*(tl_v(i,j,k+1,nnew)- &
1112 & tl_v(i,j,k ,nnew)- &
1113 & (tl_fc(i,k)*dc(i,k-1)+ &
1114 & tl_bc(i,k)*dc(i,k )+ &
1115 & tl_cf(i,k)*dc(i,k+1))- &
1116 & fc(i,k)*tl_dc(i,k-1))
1117# endif
1118 END DO
1119 END DO
1120!
1121! Tangent linear backward substitution.
1122!
1123 DO i=istr,iend
1124 tl_dc(i,n(ng))=0.0_r8
1125 END DO
1126 DO k=n(ng)-1,1,-1
1127 DO i=istr,iend
1128 tl_dc(i,k)=tl_dc(i,k)-cf(i,k)*tl_dc(i,k+1)
1129 END DO
1130 END DO
1131!
1132! Compute tl_DC before multiplying BASIC STATE spline gradients
1133! DC by AK.
1134!
1135 DO k=1,n(ng)
1136 DO i=istr,iend
1137 tl_dc(i,k)=tl_dc(i,k)*ak(i,k)+dc(i,k)*tl_ak(i,k)
1138# ifdef TL_IOMS
1139! tl_DC(i,k)=tl_DC(i,k)- &
1140! & DC(i,k)*AK(i,k) ! needed if tl_AK is computed
1141# endif
1142 dc(i,k)=dc(i,k)*ak(i,k)
1143!^ cff=dt(ng)*oHz(i,k)*(DC(i,k)-DC(i,k-1))
1144!^
1145 tl_cff=dt(ng)*(tl_ohz(i,k)*(dc(i,k)-dc(i,k-1))+ &
1146 & ohz(i,k)*(tl_dc(i,k)-tl_dc(i,k-1)))- &
1147# ifdef TL_IOMS
1148 & dt(ng)*ohz(i,k)*(dc(i,k)-dc(i,k-1))
1149# endif
1150!^ v(i,j,k,nnew)=v(i,j,k,nnew)+cff
1151!^
1152 tl_v(i,j,k,nnew)=tl_v(i,j,k,nnew)+tl_cff
1153# ifdef DIAGNOSTICS_UV
1154!! DiaV3wrk(i,j,k,M3vvis)=DiaV3wrk(i,j,k,M3vvis)+cff
1155# endif
1156 END DO
1157 END DO
1158# else
1159!
1160! Compute off-diagonal coefficients [lambda*dt*Akv/Hz] for the
1161! implicit vertical viscosity term at horizontal V-points and
1162! vertical W-points.
1163!
1164! NOTE: The original code solves the tridiagonal system A*v=r where
1165! A is a matrix and v and r are vectors. We need to solve the
1166! tangent linear of this system which is A*tl_v+tl_A*v=tl_r.
1167! Here, tl_A*v and tl_r are known, so we must solve for tl_v
1168! by inverting A*tl_v=tl_r-tl_A*v.
1169!
1170 cff=-lambda*dt(ng)/0.5_r8
1171 DO k=1,n(ng)-1
1172 DO i=istr,iend
1173 cff1=1.0_r8/(z_r(i,j,k+1)+z_r(i,j-1,k+1)- &
1174 & z_r(i,j,k )-z_r(i,j-1,k ))
1175 tl_cff1=-cff1*cff1*(tl_z_r(i,j,k+1)+tl_z_r(i,j-1,k+1)- &
1176 & tl_z_r(i,j,k )-tl_z_r(i,j-1,k ))+ &
1177# ifdef TL_IOMS
1178 & 2.0_r8*cff1
1179# endif
1180 fc(i,k)=cff*cff1*ak(i,k)
1181 tl_fc(i,k)=cff*(tl_cff1*ak(i,k)+cff1*tl_ak(i,k))
1182# ifdef TL_IOMS
1183!
1184! Uncomment if tl_Akv is computed.
1185!
1186! tl_FC(i,k)=tl_FC(i,k)-FC(i,k)
1187# endif
1188 END DO
1189 END DO
1190 DO i=istr,iend
1191 fc(i,0)=0.0_r8
1192 tl_fc(i,0)=0.0_r8
1193 fc(i,n(ng))=0.0_r8
1194 tl_fc(i,n(ng))=0.0_r8
1195 END DO
1196!
1197! Solve the tangent linear tridiagonal system.
1198!
1199 DO k=1,n(ng)
1200 DO i=istr,iend
1201 bc(i,k)=hzk(i,k)-fc(i,k)-fc(i,k-1)
1202 tl_bc(i,k)=tl_hzk(i,k)-tl_fc(i,k)-tl_fc(i,k-1)
1203 END DO
1204 END DO
1205 DO k=2,n(ng)-1
1206 DO i=istr,iend
1207# ifdef TL_IOMS
1208 dc(i,k)=tl_v(i,j,k,nnew)- &
1209 & ((tl_fc(i,k-1)-fc(i,k-1))*v(i,j,k-1,nnew)+ &
1210 & (tl_bc(i,k )-bc(i,k ))*v(i,j,k ,nnew)+ &
1211 & (tl_fc(i,k )-fc(i,k ))*v(i,j,k+1,nnew))
1212# else
1213 dc(i,k)=tl_v(i,j,k,nnew)- &
1214 & (tl_fc(i,k-1)*v(i,j,k-1,nnew)+ &
1215 & tl_bc(i,k )*v(i,j,k ,nnew)+ &
1216 & tl_fc(i,k )*v(i,j,k+1,nnew))
1217# endif
1218 END DO
1219 END DO
1220 DO i=istr,iend
1221# ifdef TL_IOMS
1222 dc(i,1)=tl_v(i,j,1,nnew)- &
1223 & ((tl_bc(i,1)-bc(i,1))*v(i,j,1,nnew)+ &
1224 & (tl_fc(i,1)-fc(i,1))*v(i,j,2,nnew))
1225 dc(i,n(ng))=tl_v(i,j,n(ng),nnew)- &
1226 & ((tl_fc(i,n(ng)-1)-fc(i,n(ng)-1))* &
1227 & v(i,j,n(ng)-1,nnew)+ &
1228 & (tl_bc(i,n(ng) )-bc(i,n(ng) ))* &
1229 & v(i,j,n(ng) ,nnew))
1230# else
1231 dc(i,1)=tl_v(i,j,1,nnew)- &
1232 & (tl_bc(i,1)*v(i,j,1,nnew)+ &
1233 & tl_fc(i,1)*v(i,j,2,nnew))
1234 dc(i,n(ng))=tl_v(i,j,n(ng),nnew)- &
1235 & (tl_fc(i,n(ng)-1)*v(i,j,n(ng)-1,nnew)+ &
1236 & tl_bc(i,n(ng) )*v(i,j,n(ng) ,nnew))
1237# endif
1238 END DO
1239 DO i=istr,iend
1240 cff=1.0_r8/bc(i,1)
1241 cf(i,1)=cff*fc(i,1)
1242 dc(i,1)=cff*dc(i,1)
1243 END DO
1244 DO k=2,n(ng)-1
1245 DO i=istr,iend
1246 cff=1.0_r8/(bc(i,k)-fc(i,k-1)*cf(i,k-1))
1247 cf(i,k)=cff*fc(i,k)
1248 dc(i,k)=cff*(dc(i,k)-fc(i,k-1)*dc(i,k-1))
1249 END DO
1250 END DO
1251!
1252! Compute new solution by back substitution.
1253! (DC is a tangent linear variable here).
1254!
1255 DO i=istr,iend
1256# ifdef DIAGNOSTICS_UV
1257!! wrk(i,N(ng))=v(i,j,N(ng),nnew)*oHz(i,N(ng))
1258# endif
1259 dc(i,n(ng))=(dc(i,n(ng))-fc(i,n(ng)-1)*dc(i,n(ng)-1))/ &
1260 & (bc(i,n(ng))-fc(i,n(ng)-1)*cf(i,n(ng)-1))
1261!^ v(i,j,N(ng),nnew)=DC(i,N(ng))
1262!^
1263 tl_v(i,j,n(ng),nnew)=dc(i,n(ng))
1264# ifdef DIAGNOSTICS_UV
1265!! DiaV3wrk(i,j,N(ng),M3vvis)=DiaV3wrk(i,j,N(ng),M3vvis)+ &
1266!! & v(i,j,N(ng),nnew)-wrk(i,N(ng))
1267# endif
1268 END DO
1269 DO k=n(ng)-1,1,-1
1270 DO i=istr,iend
1271# ifdef DIAGNOSTICS_UV
1272!! wrk(i,k)=v(i,j,k,nnew)*oHz(i,k)
1273# endif
1274 dc(i,k)=dc(i,k)-cf(i,k)*dc(i,k+1)
1275!^ v(i,j,k,nnew)=DC(i,k)
1276!^
1277 tl_v(i,j,k,nnew)=dc(i,k)
1278# ifdef DIAGNOSTICS_UV
1279!! DiaV3wrk(i,j,k,M3vvis)=DiaV3wrk(i,j,k,M3vvis)+ &
1280!! & v(i,j,k,nnew)-wrk(i,k)
1281# endif
1282 END DO
1283 END DO
1284# endif
1285!
1286! Replace INTERIOR POINTS incorrect vertical mean with more accurate
1287! barotropic component, vbar=DV_avg1/(D*om_v). Recall that, D=CF(:,0).
1288!
1289 DO i=istr,iend
1290 cf(i,0)=hzk(i,1)
1291 tl_cf(i,0)=tl_hzk(i,1)
1292 dc(i,0)=v(i,j,1,nnew)*hzk(i,1)
1293 tl_dc(i,0)=tl_v(i,j,1,nnew)*hzk(i,1)+ &
1294 & v(i,j,1,nnew)*tl_hzk(i,1)- &
1295# ifdef TL_IOMS
1296 & dc(i,0)
1297# endif
1298# ifdef NEARSHORE_MELLOR
1299 dcs(i,0)=v_stokes(i,j,1)*hzk(i,1)
1300 tl_dcs(i,0)=tl_v_stokes(i,j,1)*hzk(i,1)+ &
1301 & v_stokes(i,j,1)*tl_hzk(i,1)- &
1302# ifdef TL_IOMS
1303 & dcs(i,0)
1304# endif
1305# endif
1306# ifdef DIAGNOSTICS_UV
1307!! Dwrk(i,M2pgrd)=DiaV3wrk(i,j,1,M3pgrd)*Hzk(i,1)
1308!! Dwrk(i,M2bstr)=DiaV3wrk(i,j,1,M3vvis)*Hzk(i,1)
1309# ifdef UV_COR
1310!! Dwrk(i,M2fcor)=DiaV3wrk(i,j,1,M3fcor)*Hzk(i,1)
1311# endif
1312# if defined UV_VIS2 || defined UV_VIS4
1313!! Dwrk(i,M2xvis)=DiaV3wrk(i,j,1,M3xvis)*Hzk(i,1)
1314!! Dwrk(i,M2yvis)=DiaV3wrk(i,j,1,M3yvis)*Hzk(i,1)
1315!! Dwrk(i,M2hvis)=DiaV3wrk(i,j,1,M3hvis)*Hzk(i,1)
1316# endif
1317# ifdef UV_ADV
1318!! Dwrk(i,M2xadv)=DiaV3wrk(i,j,1,M3xadv)*Hzk(i,1)
1319!! Dwrk(i,M2yadv)=DiaV3wrk(i,j,1,M3yadv)*Hzk(i,1)
1320!! Dwrk(i,M2hadv)=DiaV3wrk(i,j,1,M3hadv)*Hzk(i,1)
1321# endif
1322# ifdef NEARSHORE_MELLOR
1323!! Dwrk(i,M2hrad)=DiaV3wrk(i,j,1,M3hrad)*Hzk(i,1)
1324# endif
1325# endif
1326 END DO
1327 DO k=2,n(ng)
1328 DO i=istr,iend
1329 cf(i,0)=cf(i,0)+hzk(i,k)
1330 tl_cf(i,0)=tl_cf(i,0)+tl_hzk(i,k)
1331 dc(i,0)=dc(i,0)+v(i,j,k,nnew)*hzk(i,k)
1332 tl_dc(i,0)=tl_dc(i,0)+ &
1333 & tl_v(i,j,k,nnew)*hzk(i,k)+ &
1334 & v(i,j,k,nnew)*tl_hzk(i,k)- &
1335# ifdef TL_IOMS
1336 & v(i,j,k,nnew)*hzk(i,k)
1337# endif
1338# ifdef NEARSHORE_MELLOR
1339 dcs(i,0)=dcs(i,0)+v_stokes(i,j,k)*hzk(i,k)
1340 tl_dcs(i,0)=tl_dcs(i,0)+ &
1341 & tl_v_stokes(i,j,k)*hzk(i,k)+ &
1342 & v_stokes(i,j,k)*tl_hzk(i,k)- &
1343# ifdef TL_IOMS
1344 & v_stokes(i,j,k)*hzk(i,k)
1345# endif
1346# endif
1347# ifdef DIAGNOSTICS_UV
1348!! Dwrk(i,M2pgrd)=Dwrk(i,M2pgrd)+ &
1349!! & DiaV3wrk(i,j,k,M3pgrd)*Hzk(i,k)
1350!! Dwrk(i,M2bstr)=Dwrk(i,M2bstr)+ &
1351!! & DiaV3wrk(i,j,k,M3vvis)*Hzk(i,k)
1352# ifdef UV_COR
1353!! Dwrk(i,M2fcor)=Dwrk(i,M2fcor)+ &
1354!! & DiaV3wrk(i,j,k,M3fcor)*Hzk(i,k)
1355# endif
1356# if defined UV_VIS2 || defined UV_VIS4
1357!! Dwrk(i,M2xvis)=Dwrk(i,M2xvis)+ &
1358!! & DiaV3wrk(i,j,k,M3xvis)*Hzk(i,k)
1359!! Dwrk(i,M2yvis)=Dwrk(i,M2yvis)+ &
1360!! & DiaV3wrk(i,j,k,M3yvis)*Hzk(i,k)
1361!! Dwrk(i,M2hvis)=Dwrk(i,M2hvis)+ &
1362!! & DiaV3wrk(i,j,k,M3hvis)*Hzk(i,k)
1363# endif
1364# ifdef UV_ADV
1365!! Dwrk(i,M2xadv)=Dwrk(i,M2xadv)+ &
1366!! & DiaV3wrk(i,j,k,M3xadv)*Hzk(i,k)
1367!! Dwrk(i,M2yadv)=Dwrk(i,M2yadv)+ &
1368!! & DiaV3wrk(i,j,k,M3yadv)*Hzk(i,k)
1369!! Dwrk(i,M2hadv)=Dwrk(i,M2hadv)+ &
1370!! & DiaV3wrk(i,j,k,M3hadv)*Hzk(i,k)
1371# endif
1372# ifdef NEARSHORE_MELLOR
1373!! Dwrk(i,M2hrad)=Dwrk(i,M2hrad)+ &
1374!! & DiaV3wrk(i,j,k,M3hrad)*Hzk(i,k)
1375# endif
1376# endif
1377 END DO
1378 END DO
1379 DO i=istr,iend
1380 dc1(i,0)=dc(i,0) ! intermediate
1381 cff1=1.0_r8/(cf(i,0)*om_v(i,j))
1382 tl_cff1=-cff1*cff1*tl_cf(i,0)*om_v(i,j)+ &
1383# ifdef TL_IOMS
1384 & 2.0_r8*cff1
1385# endif
1386 dc(i,0)=(dc(i,0)*om_v(i,j)-dv_avg1(i,j))*cff1 ! recursive
1387 tl_dc(i,0)=(tl_dc(i,0)*om_v(i,j)-tl_dv_avg1(i,j))*cff1+ &
1388 & (dc1(i,0)*om_v(i,j)-dv_avg1(i,j))*tl_cff1- &
1389# ifdef TL_IOMS
1390 & dc(i,0)
1391# endif
1392# ifdef NEARSHORE_MELLOR
1393 dcs1(i)=dcs(i,0) ! intermediate
1394 cff2=1.0_r8/cf(i,0)
1395 tl_cff2=-cff2*cff2*tl_cf(i,0)+ &
1396# ifdef TL_IOMS
1397 & 2.0_r8*cff2
1398# endif
1399 dcs(i,0)=dcs(i,0)*cff2-vbar_stokes(i,j) ! recursive
1400 tl_dcs(i,0)=tl_dcs(i,0)*cff2+ &
1401 & dcs1(i,0)*tl_cff2-tl_vbar_stokes(i,j)- &
1402# ifdef TL_IOMS
1403 & dcs(i,0)
1404# endif
1405# endif
1406# ifdef DIAGNOSTICS_UV
1407!! DO idiag=1,M2pgrd
1408!! Dwrk(i,idiag)=(Dwrk(i,idiag)*om_v(i,j)- &
1409!! & DiaV2wrk(i,j,idiag))*cff1
1410!! END DO
1411!! Dwrk(i,M2bstr)=(Dwrk(i,M2bstr)*om_v(i,j)- &
1412!! & DiaV2wrk(i,j,M2bstr)-DiaV2wrk(i,j,M2sstr))* &
1413!! & cff1
1414# endif
1415 END DO
1416!
1417! Couple and update new solution.
1418!
1419 DO k=1,n(ng)
1420 DO i=istr,iend
1421!^ v(i,j,k,nnew)=v(i,j,k,nnew)-DC(i,0)
1422!^
1423 tl_v(i,j,k,nnew)=tl_v(i,j,k,nnew)-tl_dc(i,0)
1424# ifdef MASKING
1425!^ v(i,j,k,nnew)=v(i,j,k,nnew)*vmask(i,j)
1426!^
1427 tl_v(i,j,k,nnew)=tl_v(i,j,k,nnew)*vmask(i,j)
1428# endif
1429# ifdef WET_DRY_NOT_YET
1430!^ v(i,j,k,nnew)=v(i,j,k,nnew)*vmask_wet(i,j)
1431!^
1432 tl_v(i,j,k,nnew)=tl_v(i,j,k,nnew)*vmask_wet(i,j)
1433# endif
1434# ifdef NEARSHORE_MELLOR
1435!^ v_stokes(i,j,k)=v_stokes(i,j,k)-DCs(i,0)
1436!^
1437 tl_v_stokes(i,j,k)=tl_v_stokes(i,j,k)-tl_dcs(i,0)
1438# ifdef MASKING
1439!^ v_stokes(i,j,k)=v_stokes(i,j,k)*vmask(i,j)
1440!^
1441 tl_v_stokes(i,j,k)=tl_v_stokes(i,j,k)*vmask(i,j)
1442# endif
1443# ifdef WET_DRY_NOT_YET
1444!^ v_stokes(i,j,k)=v_stokes(i,j,k)*vmask_wet(i,j)
1445!^
1446 tl_v_stokes(i,j,k)=tl_v_stokes(i,j,k)*vmask_wet(i,j)
1447# endif
1448# endif
1449# ifdef DIAGNOSTICS_UV
1450!! DiaV3wrk(i,j,k,M3pgrd)=DiaV3wrk(i,j,k,M3pgrd)- &
1451!! & Dwrk(i,M2pgrd)
1452!! DiaV3wrk(i,j,k,M3vvis)=DiaV3wrk(i,j,k,M3vvis)- &
1453!! & Dwrk(i,M2bstr)
1454# ifdef UV_COR
1455!! DiaV3wrk(i,j,k,M3fcor)=DiaV3wrk(i,j,k,M3fcor)- &
1456!! & Dwrk(i,M2fcor)
1457# endif
1458# if defined UV_VIS2 || defined UV_VIS4
1459!! DiaV3wrk(i,j,k,M3xvis)=DiaV3wrk(i,j,k,M3xvis)- &
1460!! & Dwrk(i,M2xvis)
1461!! DiaV3wrk(i,j,k,M3yvis)=DiaV3wrk(i,j,k,M3yvis)- &
1462!! & Dwrk(i,M2yvis)
1463!! DiaV3wrk(i,j,k,M3hvis)=DiaV3wrk(i,j,k,M3hvis)- &
1464!! & Dwrk(i,M2hvis)
1465# endif
1466# ifdef UV_ADV
1467!! DiaV3wrk(i,j,k,M3xadv)=DiaV3wrk(i,j,k,M3xadv)- &
1468!! & Dwrk(i,M2xadv)
1469!! DiaV3wrk(i,j,k,M3yadv)=DiaV3wrk(i,j,k,M3yadv)- &
1470!! & Dwrk(i,M2yadv)
1471!! DiaV3wrk(i,j,k,M3hadv)=DiaV3wrk(i,j,k,M3hadv)- &
1472!! & Dwrk(i,M2hadv)
1473# endif
1474# ifdef NEARSHORE_MELLOR
1475!! DiaV3wrk(i,j,k,M3hrad)=DiaV3wrk(i,j,k,M3hrad)- &
1476!! & Dwrk(i,M2hrad)
1477# endif
1478# endif
1479 END DO
1480 END DO
1481
1482# if defined DIAGNOSTICS_UV && defined MASKING
1483!! DO k=1,N(ng)
1484!! DO i=Istr,Iend
1485!! DO idiag=1,NDM3d
1486!! DiaV3wrk(i,j,k,idiag)=DiaV3wrk(i,j,k,idiag)*vmask(i,j)
1487!! END DO
1488!! END DO
1489!! END DO
1490# endif
1491 END IF
1492 END DO
1493!
1494!-----------------------------------------------------------------------
1495! Set lateral boundary conditions.
1496!-----------------------------------------------------------------------
1497!
1498!^ CALL u3dbc_tile (ng, tile, &
1499!^ & LBi, UBi, LBj, UBj, N(ng), &
1500!^ & IminS, ImaxS, JminS, JmaxS, &
1501!^ & nstp, nnew, &
1502!^ & u)
1503!^
1504 CALL rp_u3dbc_tile (ng, tile, &
1505 & lbi, ubi, lbj, ubj, n(ng), &
1506 & imins, imaxs, jmins, jmaxs, &
1507 & nstp, nnew, &
1508 & tl_u)
1509!^ CALL v3dbc_tile (ng, tile, &
1510!^ & LBi, UBi, LBj, UBj, N(ng), &
1511!^ & IminS, ImaxS, JminS, JmaxS, &
1512!^ & nstp, nnew, &
1513!^ & v)
1514!^
1515 CALL rp_v3dbc_tile (ng, tile, &
1516 & lbi, ubi, lbj, ubj, n(ng), &
1517 & imins, imaxs, jmins, jmaxs, &
1518 & nstp, nnew, &
1519 & tl_v)
1520!
1521!-----------------------------------------------------------------------
1522! Apply momentum transport point sources (like river runoff), if any.
1523!-----------------------------------------------------------------------
1524!
1525 IF (luvsrc(ng)) THEN
1526 DO is=1,nsrc(ng)
1527 i=sources(ng)%Isrc(is)
1528 j=sources(ng)%Jsrc(is)
1529 IF (((istrr.le.i).and.(i.le.iendr)).and. &
1530 & ((jstrr.le.j).and.(j.le.jendr))) THEN
1531 IF (int(sources(ng)%Dsrc(is)).eq.0) THEN
1532 DO k=1,n(ng)
1533 cff1=1.0_r8/(on_u(i,j)* &
1534 & 0.5_r8*(z_w(i-1,j,k)-z_w(i-1,j,k-1)+ &
1535 & z_w(i ,j,k)-z_w(i ,j,k-1)))
1536 tl_cff1=-cff1*cff1*on_u(i,j)* &
1537 & 0.5_r8*(tl_z_w(i-1,j,k)-tl_z_w(i-1,j,k-1)+ &
1538 & tl_z_w(i ,j,k)-tl_z_w(i ,j,k-1))+ &
1539# ifdef TL_IOMS
1540 & 2.0_r8*cff1
1541# endif
1542!^ u(i,j,k,nnew)=SOURCES(ng)%Qsrc(is,k)*cff1
1543!^
1544 tl_u(i,j,k,nnew)=sources(ng)%tl_Qsrc(is,k)*cff1+ &
1545 & sources(ng)%Qsrc(is,k)*tl_cff1- &
1546# ifdef TL_IOMS
1547 & sources(ng)%Qsrc(is,k)*cff1
1548# endif
1549 END DO
1550 ELSE IF (int(sources(ng)%Dsrc(is)).eq.1) THEN
1551 DO k=1,n(ng)
1552 cff1=1.0_r8/(om_v(i,j)* &
1553 & 0.5_r8*(z_w(i,j-1,k)-z_w(i,j-1,k-1)+ &
1554 & z_w(i,j ,k)-z_w(i,j ,k-1)))
1555 tl_cff1=-cff1*cff1*om_v(i,j)* &
1556 & 0.5_r8*(tl_z_w(i,j-1,k)-tl_z_w(i,j-1,k-1)+ &
1557 & tl_z_w(i,j ,k)-tl_z_w(i,j ,k-1))+ &
1558# ifdef TL_IOMS
1559 & 2.0_r8*cff1
1560# endif
1561!^ v(i,j,k,nnew)=SOURCES(ng)%Qsrc(is,k)*cff1
1562!^
1563 tl_v(i,j,k,nnew)=sources(ng)%tl_Qsrc(is,k)*cff1+ &
1564 & sources(ng)%Qsrc(is,k)*tl_cff1- &
1565# ifdef TL_IOMS
1566 & sources(ng)%Qsrc(is,k)*cff1
1567# endif
1568 END DO
1569 END IF
1570 END IF
1571 END DO
1572 END IF
1573!
1574!-----------------------------------------------------------------------
1575! Couple 2D and 3D momentum equations.
1576!-----------------------------------------------------------------------
1577!
1578! Couple velocity component in the XI-direction.
1579!
1580 DO j=jstrt,jendt
1581 DO i=istrp,iendt
1582 dc(i,0)=0.0_r8
1583 tl_dc(i,0)=0.0_r8
1584 cf(i,0)=0.0_r8
1585 tl_cf(i,0)=0.0_r8
1586# ifdef NEARSHORE_MELLOR
1587 cfs(i,0)=0.0_r8
1588 tl_cfs(i,0)=0.0_r8
1589# endif
1590 fc(i,0)=0.0_r8
1591 tl_fc(i,0)=0.0_r8
1592 END DO
1593!
1594! Compute thicknesses of U-boxes DC(i,1:N), total depth of the water
1595! column DC(i,0), and incorrect vertical mean CF(i,0). Notice that
1596! barotropic component is replaced with its fast-time averaged
1597! values.
1598!
1599 DO k=1,n(ng)
1600 DO i=istrp,iendt
1601 cff=0.5_r8*on_u(i,j)
1602 dc(i,k)=cff*(hz(i,j,k)+hz(i-1,j,k))
1603 tl_dc(i,k)=cff*(tl_hz(i,j,k)+tl_hz(i-1,j,k))
1604 dc(i,0)=dc(i,0)+dc(i,k)
1605 tl_dc(i,0)=tl_dc(i,0)+tl_dc(i,k)
1606 cf(i,0)=cf(i,0)+ &
1607 & dc(i,k)*u(i,j,k,nnew)
1608 tl_cf(i,0)=tl_cf(i,0)+ &
1609 & tl_dc(i,k)*u(i,j,k,nnew)+ &
1610 & dc(i,k)*tl_u(i,j,k,nnew)- &
1611# ifdef TL_IOMS
1612 & dc(i,k)*u(i,j,k,nnew)
1613# endif
1614# ifdef NEARSHORE_MELLOR
1615 cfs(i,0)=cfs(i,0)+ &
1616 & dc(i,k)*u_stokes(i,j,k)
1617 tl_cfs(i,0)=tl_cfs(i,0)+ &
1618 & tl_dc(i,k)*u_stokes(i,j,k)+ &
1619 & dc(i,k)*tl_u_stokes(i,j,k)- &
1620# ifdef TL_IOMS
1621 & dc(i,k)*u_stokes(i,j,k)
1622# endif
1623# endif
1624 END DO
1625 END DO
1626 DO i=istrp,iendt
1627 dc1(i,0)=dc(i,0) ! intermediate
1628# ifdef NEARSHORE_MELLOR
1629 cff2=dc(i,0)*ubar_stokes(i,j)
1630 tl_cff2=tl_dc(i,0)*ubar_stokes(i,j)+ &
1631 & dc(i,0)*tl_ubar_stokes(i,j)- &
1632# ifdef TL_IOMS
1633 & cff2
1634# endif
1635# endif
1636 dc(i,0)=1.0_r8/dc(i,0) ! recursive
1637 tl_dc(i,0)=-tl_dc(i,0)/(dc1(i,0)*dc1(i,0))+ &
1638# ifdef TL_IOMS
1639 & 2.0_r8*dc(i,0)
1640# endif
1641 cf1(i)=cf(i,0) ! intermediate
1642 cf(i,0)=dc(i,0)*(cf(i,0)-du_avg1(i,j)) ! recursive
1643 tl_cf(i,0)=tl_dc(i,0)*(cf1(i)-du_avg1(i,j))+ &
1644 & dc(i,0)*(tl_cf(i,0)-tl_du_avg1(i,j))- &
1645# ifdef TL_IOMS
1646 & cf(i,0)
1647# endif
1648# ifdef NEARSHORE_MELLOR
1649 cfs1(i)=cfs(i,0) ! intermediate
1650 cfs(i,0)=dc(i,0)*(cfs(i,0)-cff2) ! recursive
1651 tl_cfs(i,0)=tl_dc(i,0)*(cfs1(i)-cff2)+ &
1652 & dc(i,0)*(tl_cfs(i,0)-tl_cff2)- &
1653# ifdef TL_IOMS
1654 & cfs(i,0)
1655# endif
1656# endif
1657!^ ubar(i,j,1)=DC(i,0)*DU_avg1(i,j)
1658!^
1659 tl_ubar(i,j,1)=tl_dc(i,0)*du_avg1(i,j)+ &
1660 & dc(i,0)*tl_du_avg1(i,j)- &
1661# ifdef TL_IOMS
1662 & dc(i,0)*du_avg1(i,j)
1663# endif
1664!^ ubar(i,j,2)=ubar(i,j,1)
1665!^
1666 tl_ubar(i,j,2)=tl_ubar(i,j,1)
1667# ifdef DIAGNOSTICS_UV
1668!! DiaU2wrk(i,j,M2rate)=ubar(i,j,1)-DiaU2int(i,j,M2rate)*DC(i,0)
1669!! DiaU2int(i,j,M2rate)=ubar(i,j,1)*DC1(i,0)
1670# endif
1671 END DO
1672# ifdef DIAGNOSTICS_UV
1673!!
1674!! Convert the units of the fast-time integrated change in mass flux
1675!! diagnostic values to change in velocity (m s-1).
1676!!
1677!! DO idiag=1,NDM2d-1
1678!! DO i=IstrP,IendT
1679!! DiaU2wrk(i,j,idiag)=DC(i,0)*DiaU2wrk(i,j,idiag)
1680# ifdef MASKING
1681!! DiaU2wrk(i,j,idiag)=DiaU2wrk(i,j,idiag)*umask(i,j)
1682# endif
1683!! END DO
1684!! END DO
1685# endif
1686!
1687! Replace only BOUNDARY POINTS incorrect vertical mean with more
1688! accurate barotropic component, ubar=DU_avg1/(D*on_u). Recall that,
1689! D=CF(:,0).
1690!
1691! NOTE: Only the BOUNDARY POINTS need to be replaced. Avoid redundant
1692! update in the interior again for computational purposes which
1693! will not affect the nonlinear code. However, the adjoint
1694! code is wrong because the interior solution is corrected
1695! twice. The replacement is avoided if the boundary edge is
1696! periodic. The J-loop is pipelined, so we need to do a special
1697! test for the southern and northern domain edges.
1698!
1699 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
1700 IF (domain(ng)%Western_Edge(tile)) THEN
1701 DO k=1,n(ng)
1702!^ u(Istr,j,k,nnew)=u(Istr,j,k,nnew)-CF(Istr,0)
1703!^
1704 tl_u(istr,j,k,nnew)=tl_u(istr,j,k,nnew)- &
1705 & tl_cf(istr,0)
1706# ifdef MASKING
1707!^ u(Istr,j,k,nnew)=u(Istr,j,k,nnew)* &
1708!^ & umask(Istr,j)
1709!^
1710 tl_u(istr,j,k,nnew)=tl_u(istr,j,k,nnew)* &
1711 & umask(istr,j)
1712# endif
1713# ifdef WET_DRY_NOT_YET
1714!^ u(Istr,j,k,nnew)=u(Istr,j,k,nnew)* &
1715!^ & umask_wet(Istr,j)
1716!^
1717 tl_u(istr,j,k,nnew)=tl_u(istr,j,k,nnew)* &
1718 & umask_wet(istr,j)
1719# endif
1720# ifdef NEARSHORE_MELLOR
1721!^ u_stokes(Istr,j,k)=u_stokes(Istr,j,k)-CFs(Istr,0)
1722!^
1723 tl_u_stokes(istr,j,k)=tl_u_stokes(istr,j,k)- &
1724 & tl_cfs(istr,0)
1725# ifdef MASKING
1726!^ u_stokes(Istr,j,k)=u_stokes(Istr,j,k)* &
1727!^ & umask(Istr,j)
1728!^
1729 tl_u_stokes(istr,j,k)=tl_u_stokes(istr,j,k)* &
1730 & umask(istr,j)
1731# endif
1732# ifdef WET_DRY_NOT_YET
1733!^ u_stokes(Istr,j,k)=u_stokes(Istr,j,k)* &
1734!^ & umask_wet(Istr,j)
1735!^
1736 tl_u_stokes(istr,j,k)=tl_u_stokes(istr,j,k)* &
1737 & umask_wet(istr,j)
1738# endif
1739# endif
1740 END DO
1741 END IF
1742 END IF
1743!
1744 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
1745 IF (domain(ng)%Eastern_Edge(tile)) THEN
1746 DO k=1,n(ng)
1747!^ u(Iend+1,j,k,nnew)=u(Iend+1,j,k,nnew)-CF(Iend+1,0)
1748!^
1749 tl_u(iend+1,j,k,nnew)=tl_u(iend+1,j,k,nnew)- &
1750 & tl_cf(iend+1,0)
1751# ifdef MASKING
1752!^ u(Iend+1,j,k,nnew)=u(Iend+1,j,k,nnew)* &
1753!^ & umask(Iend+1,j)
1754!^
1755 tl_u(iend+1,j,k,nnew)=tl_u(iend+1,j,k,nnew)* &
1756 & umask(iend+1,j)
1757# endif
1758# ifdef WET_DRY_NOT_YET
1759!^ u(Iend+1,j,k,nnew)=u(Iend+1,j,k,nnew)* &
1760!^ & umask_wet(Iend+1,j)
1761!^
1762 tl_u(iend+1,j,k,nnew)=tl_u(iend+1,j,k,nnew)* &
1763 & umask_wet(iend+1,j)
1764# endif
1765# ifdef NEARSHORE_MELLOR
1766!^ u_stokes(Iend+1,j,k)=u_stokes(Iend+1,j,k)-CFs(Iend+1,0)
1767!^
1768 tl_u_stokes(iend+1,j,k)=tl_u_stokes(iend+1,j,k)- &
1769 & tl_cfs(iend+1,0)
1770# ifdef MASKING
1771!^ u_stokes(Iend+1,j,k)=u_stokes(Iend+1,j,k)* &
1772!^ & umask(Iend+1,j)
1773!^
1774 tl_u_stokes(iend+1,j,k)=tl_u_stokes(iend+1,j,k)* &
1775 & umask(iend+1,j)
1776# endif
1777# ifdef WET_DRY_NOT_YET
1778!^ u_stokes(Iend+1,j,k)=u_stokes(Iend+1,j,k)* &
1779!^ & umask_wet(Iend+1,j)
1780!^
1781 tl_u_stokes(iend+1,j,k)=tl_u_stokes(iend+1,j,k)* &
1782 & umask_wet(iend+1,j)
1783# endif
1784# endif
1785 END DO
1786 END IF
1787 END IF
1788!
1789 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
1790 IF (j.eq.0) THEN ! southern boundary
1791 DO k=1,n(ng) ! J-loop pipelined
1792 DO i=istru,iend
1793!^ u(i,j,k,nnew)=u(i,j,k,nnew)-CF(i,0)
1794!^
1795 tl_u(i,j,k,nnew)=tl_u(i,j,k,nnew)- &
1796 & tl_cf(i,0)
1797# ifdef MASKING
1798!^ u(i,j,k,nnew)=u(i,j,k,nnew)*umask(i,j)
1799!^
1800 tl_u(i,j,k,nnew)=tl_u(i,j,k,nnew)* &
1801 & umask(i,j)
1802# endif
1803# ifdef WET_DRY_NOT_YET
1804!^ u(i,j,k,nnew)=u(i,j,k,nnew)* &
1805!^ & umask_wet(i,j)
1806!^
1807 tl_u(i,j,k,nnew)=tl_u(i,j,k,nnew)* &
1808 & umask_wet(i,j)
1809# endif
1810# ifdef NEARSHORE_MELLOR
1811!^ u_stokes(i,j,k)=u_stokes(i,j,k)-CFs(i,0)
1812!^
1813 tl_u_stokes(i,j,k)=tl_u_stokes(i,j,k)- &
1814 & tl_cfs(i,0)
1815# ifdef MASKING
1816!^ u_stokes(i,j,k)=u_stokes(i,j,k)* &
1817!^ & umask(i,j)
1818!^
1819 tl_u_stokes(i,j,k)=tl_u_stokes(i,j,k)* &
1820 & umask(i,j)
1821# endif
1822# ifdef WET_DRY_NOT_YET
1823!^ u_stokes(i,j,k)=u_stokes(i,j,k)* &
1824!^ & umask_wet(i,j)
1825!^
1826 tl_u_stokes(i,j,k)=tl_u_stokes(i,j,k)* &
1827 & umask_wet(i,j)
1828# endif
1829# endif
1830 END DO
1831 END DO
1832 END IF
1833 END IF
1834!
1835 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
1836 IF (j.eq.mm(ng)+1) THEN ! northern boundary
1837 DO k=1,n(ng) ! J-loop pipelined
1838 DO i=istru,iend
1839!^ u(i,j,k,nnew)=u(i,j,k,nnew)-CF(i,0)
1840!^
1841 tl_u(i,j,k,nnew)=tl_u(i,j,k,nnew)- &
1842 & tl_cf(i,0)
1843# ifdef MASKING
1844!^ u(i,j,k,nnew)=u(i,j,k,nnew)*umask(i,j)
1845!^
1846 tl_u(i,j,k,nnew)=tl_u(i,j,k,nnew)* &
1847 & umask(i,j)
1848# endif
1849# ifdef WET_DRY_NOT_YET
1850!^ u(i,j,k,nnew)=u(i,j,k,nnew)*umask_wet(i,j)
1851!^
1852 tl_u(i,j,k,nnew)=tl_u(i,j,k,nnew)* &
1853 & umask_wet(i,j)
1854# endif
1855# ifdef NEARSHORE_MELLOR
1856!^ u_stokes(i,j,k)=u_stokes(i,j,k)-CFs(i,0)
1857!^
1858 tl_u_stokes(i,j,k)=tl_u_stokes(i,j,k)- &
1859 & tl_cfs(i,0)
1860# ifdef MASKING
1861!^ u_stokes(i,j,k)=u_stokes(i,j,k)* &
1862!^ & umask(i,j)
1863!^
1864 tl_u_stokes(i,j,k)=tl_u_stokes(i,j,k)* &
1865 & umask(i,j)
1866# endif
1867# ifdef WET_DRY_NOT_YET
1868!^ u_stokes(i,j,k)=u_stokes(i,j,k)* &
1869!^ & umask_wet(i,j)
1870!^
1871 tl_u_stokes(i,j,k)=tl_u_stokes(i,j,k)* &
1872 & umask_wet(i,j)
1873# endif
1874# endif
1875 END DO
1876 END DO
1877 END IF
1878 END IF
1879!
1880! Compute correct mass flux, Hz*u/n.
1881!
1882 DO k=n(ng),1,-1
1883 DO i=istrp,iendt
1884 huon(i,j,k)=0.5_r8*(huon(i,j,k)+u(i,j,k,nnew)*dc(i,k))
1885 tl_huon(i,j,k)=0.5_r8*(tl_huon(i,j,k)+ &
1886 & tl_u(i,j,k,nnew)*dc(i,k)+ &
1887 & u(i,j,k,nnew)*tl_dc(i,k))- &
1888# ifdef TL_IOMS
1889 & 0.5_r8*u(i,j,k,nnew)*dc(i,k)
1890# endif
1891# ifdef NEARSHORE_MELLOR
1892 huon(i,j,k)=huon(i,j,k)+0.5_r8*u_stokes(i,j,k)*dc(i,k)
1893 tl_huon(i,j,k)=tl_huon(i,j,k)+ &
1894 & 0.5_r8*(tl_u_stokes(i,j,k)*dc(i,k)+ &
1895 & u_stokes(i,j,k)*tl_dc(i,k))- &
1896# ifdef TL_IOMS
1897 & 0.5_r8*u_stokes(i,j,k)*dc(i,k)
1898# endif
1899# endif
1900 fc(i,0)=fc(i,0)+huon(i,j,k)
1901 tl_fc(i,0)=tl_fc(i,0)+tl_huon(i,j,k)
1902# ifdef DIAGNOSTICS_UV
1903!! DiaU3wrk(i,j,k,M3rate)=u(i,j,k,nnew)-DiaU3wrk(i,j,k,M3rate)
1904# endif
1905 END DO
1906 END DO
1907 DO i=istrp,iendt
1908 fc1(i)=fc(i,0) ! intermediate
1909 fc(i,0)=dc(i,0)*(fc(i,0)-du_avg2(i,j)) ! recursive
1910 tl_fc(i,0)=tl_dc(i,0)*(fc1(i)-du_avg2(i,j))+ &
1911 & dc(i,0)*(tl_fc(i,0)-tl_du_avg2(i,j))- &
1912# ifdef TL_IOMS
1913 & fc(i,0)
1914# endif
1915 END DO
1916 DO k=1,n(ng)
1917 DO i=istrp,iendt
1918 huon(i,j,k)=huon(i,j,k)-dc(i,k)*fc(i,0)
1919 tl_huon(i,j,k)=tl_huon(i,j,k)- &
1920 & tl_dc(i,k)*fc(i,0)- &
1921 & dc(i,k)*tl_fc(i,0)+ &
1922# ifdef TL_IOMS
1923 & dc(i,k)*fc(i,0)
1924# endif
1925 END DO
1926 END DO
1927!
1928! Couple velocity component in the ETA-direction.
1929!
1930 IF (j.ge.jstr) THEN
1931 DO i=istrt,iendt
1932 dc(i,0)=0.0_r8
1933 tl_dc(i,0)=0.0_r8
1934 cf(i,0)=0.0_r8
1935 tl_cf(i,0)=0.0_r8
1936# ifdef NEARSHORE_MELLOR
1937 cfs(i,0)=0.0_r8
1938 tl_cfs(i,0)=0.0_r8
1939# endif
1940 fc(i,0)=0.0_r8
1941 tl_fc(i,0)=0.0_r8
1942 END DO
1943!
1944! Compute thicknesses of V-boxes DC(i,1:N), total depth of the water
1945! column DC(i,0), and incorrect vertical mean CF(i,0). Notice that
1946! barotropic component is replaced with its fast-time averaged
1947! values.
1948!
1949 DO k=1,n(ng)
1950 DO i=istrt,iendt
1951 cff=0.5_r8*om_v(i,j)
1952 dc(i,k)=cff*(hz(i,j,k)+hz(i,j-1,k))
1953 tl_dc(i,k)=cff*(tl_hz(i,j,k)+tl_hz(i,j-1,k))
1954 dc(i,0)=dc(i,0)+dc(i,k)
1955 tl_dc(i,0)=tl_dc(i,0)+tl_dc(i,k)
1956 cf(i,0)=cf(i,0)+ &
1957 & dc(i,k)*v(i,j,k,nnew)
1958 tl_cf(i,0)=tl_cf(i,0)+ &
1959 & tl_dc(i,k)*v(i,j,k,nnew)+ &
1960 & dc(i,k)*tl_v(i,j,k,nnew)- &
1961# ifdef TL_IOMS
1962 & dc(i,k)*v(i,j,k,nnew)
1963# endif
1964# ifdef NEARSHORE_MELLOR
1965 cfs(i,0)=cfs(i,0)+ &
1966 & dc(i,k)*v_stokes(i,j,k)
1967 tl_cfs(i,0)=tl_cfs(i,0)+ &
1968 & tl_dc(i,k)*v_stokes(i,j,k)+ &
1969 & dc(i,k)*tl_v_stokes(i,j,k)- &
1970# ifdef TL_IOMS
1971 & dc(i,k)*v_stokes(i,j,k)
1972# endif
1973# endif
1974 END DO
1975 END DO
1976 DO i=istrt,iendt
1977 dc1(i,0)=dc(i,0) ! intermediate
1978# ifdef NEARSHORE_MELLOR
1979 cff2=dc(i,0)*vbar_stokes(i,j)
1980 tl_cff2=tl_dc(i,0)*vbar_stokes(i,j)+ &
1981 & dc(i,0)*tl_vbar_stokes(i,j)- &
1982# ifdef TL_IOMS
1983 & cff2
1984# endif
1985# endif
1986 dc(i,0)=1.0_r8/dc(i,0) ! recursive
1987 tl_dc(i,0)=-tl_dc(i,0)/(dc1(i,0)*dc1(i,0))+ &
1988# ifdef TL_IOMS
1989 & 2.0_r8*dc(i,0)
1990# endif
1991 cf1(i)=cf(i,0) ! intermediate
1992 cf(i,0)=dc(i,0)*(cf(i,0)-dv_avg1(i,j)) ! recursive
1993 tl_cf(i,0)=tl_dc(i,0)*(cf1(i)-dv_avg1(i,j))+ &
1994 & dc(i,0)*(tl_cf(i,0)-tl_dv_avg1(i,j))- &
1995# ifdef TL_IOMS
1996 & cf(i,0)
1997# endif
1998# ifdef NEARSHORE_MELLOR
1999 cfs1(i)=cfs(i,0)
2000 cfs(i,0)=dc(i,0)*(cfs(i,0)-cff2) ! recursive
2001 tl_cfs(i,0)=tl_dc(i,0)*(cfs1(i)-cff2)+ &
2002 & dc(i,0)*(tl_cfs(i,0)-tl_cff2)- &
2003# ifdef TL_IOMS
2004 & cf(i,0)
2005# endif
2006# endif
2007!^ vbar(i,j,1)=DC(i,0)*DV_avg1(i,j)
2008!^
2009 tl_vbar(i,j,1)=tl_dc(i,0)*dv_avg1(i,j)+ &
2010 & dc(i,0)*tl_dv_avg1(i,j)- &
2011# ifdef TL_IOMS
2012 & dc(i,0)*dv_avg1(i,j)
2013# endif
2014!^ vbar(i,j,2)=vbar(i,j,1)
2015!^
2016 tl_vbar(i,j,2)=tl_vbar(i,j,1)
2017# ifdef DIAGNOSTICS_UV
2018!! DiaV2wrk(i,j,M2rate)=vbar(i,j,1)- &
2019!! & DiaV2int(i,j,M2rate)*DC(i,0)
2020!! DiaV2int(i,j,M2rate)=vbar(i,j,1)*DC1(i,0)
2021!! DiaV2wrk(i,j,M2rate)=vbar(i,j,1)-DiaV2int(i,j,M2rate)
2022!! DiaV2int(i,j,M2rate)=vbar(i,j,1)
2023# endif
2024 END DO
2025# ifdef DIAGNOSTICS_UV
2026!!
2027!! Convert the units of the fast-time integrated change in mass flux
2028!! diagnostic values to change in velocity (m s-1).
2029!!
2030!! DO idiag=1,NDM2d-1
2031!! DO i=IstrT,IendT
2032!! DiaV2wrk(i,j,idiag)=DC(i,0)*DiaV2wrk(i,j,idiag)
2033# ifdef MASKING
2034!! DiaV2wrk(i,j,idiag)=DiaV2wrk(i,j,idiag)*vmask(i,j)
2035# endif
2036!! END DO
2037!! END DO
2038# endif
2039!
2040! Replace only BOUNDARY POINTS incorrect vertical mean with more
2041! accurate barotropic component, vbar=DV_avg1/(D*om_v). Recall that,
2042! D=CF(:,0).
2043!
2044! NOTE: Only the BOUNDARY POINTS need to be replaced. Avoid redundant
2045! update in the interior again for computational purposes which
2046! will not affect the nonlinear code. However, the adjoint
2047! code is wrong because the interior solution is corrected
2048! twice. The replacement is avoided if the boundary edge is
2049! periodic. The J-loop is pipelined, so we need to do a special
2050! test for the southern and northern domain edges.
2051!
2052 IF (.not.(compositegrid(iwest,ng).or.ewperiodic(ng))) THEN
2053 IF (domain(ng)%Western_Edge(tile)) THEN
2054 DO k=1,n(ng)
2055!^ v(Istr-1,j,k,nnew)=v(Istr-1,j,k,nnew)-CF(Istr-1,0)
2056!^
2057 tl_v(istr-1,j,k,nnew)=tl_v(istr-1,j,k,nnew)- &
2058 & tl_cf(istr-1,0)
2059# ifdef MASKING
2060!^ v(Istr-1,j,k,nnew)=v(Istr-1,j,k,nnew)* &
2061!^ & vmask(Istr-1,j)
2062!^
2063 tl_v(istr-1,j,k,nnew)=tl_v(istr-1,j,k,nnew)* &
2064 & vmask(istr-1,j)
2065# endif
2066# ifdef WET_DRY_NOT_YET
2067!^ v(Istr-1,j,k,nnew)=v(Istr-1,j,k,nnew)* &
2068!^ & vmask_wet(Istr-1,j)
2069!^
2070 tl_v(istr-1,j,k,nnew)=tl_v(istr-1,j,k,nnew)* &
2071 & vmask_wet(istr-1,j)
2072# endif
2073# ifdef NEARSHORE_MELLOR
2074!^ v_stokes(Istr-1,j,k)=v_stokes(Istr-1,j,k)- &
2075!^ & CFs(Istr-1,0)
2076!^
2077 tl_v_stokes(istr-1,j,k)=tl_v_stokes(istr-1,j,k)- &
2078 & tl_cfs(istr-1,0)
2079# ifdef MASKING
2080!^ v_stokes(Istr-1,j,k)=v_stokes(Istr-1,j,k)* &
2081!^ & vmask(Istr-1,j)
2082!^
2083 tl_v_stokes(istr-1,j,k)=tl_v_stokes(istr-1,j,k)* &
2084 & vmask(istr-1,j)
2085# endif
2086# ifdef WET_DRY_NOT_YET
2087!^ v_stokes(Istr-1,j,k)=v_stokes(Istr-1,j,k)* &
2088!^ & vmask_wet(Istr-1,j)
2089!^
2090 tl_v_stokes(istr-1,j,k)=tl_v_stokes(istr-1,j,k)* &
2091 & vmask_wet(istr-1,j)
2092# endif
2093# endif
2094 END DO
2095 END IF
2096 END IF
2097!
2098 IF (.not.(compositegrid(ieast,ng).or.ewperiodic(ng))) THEN
2099 IF (domain(ng)%Eastern_Edge(tile)) THEN
2100 DO k=1,n(ng)
2101!^ v(Iend+1,j,k,nnew)=v(Iend+1,j,k,nnew)-CF(Iend+1,0)
2102!^
2103 tl_v(iend+1,j,k,nnew)=tl_v(iend+1,j,k,nnew)- &
2104 & tl_cf(iend+1,0)
2105# ifdef MASKING
2106!^ v(Iend+1,j,k,nnew)=v(Iend+1,j,k,nnew)* &
2107!^ & vmask(Iend+1,j)
2108!^
2109 tl_v(iend+1,j,k,nnew)=tl_v(iend+1,j,k,nnew)* &
2110 & vmask(iend+1,j)
2111# endif
2112# ifdef WET_DRY_NOT_YET
2113!^ v(Iend+1,j,k,nnew)=v(Iend+1,j,k,nnew)* &
2114!^ & vmask_wet(Iend+1,j)
2115!^
2116 tl_v(iend+1,j,k,nnew)=tl_v(iend+1,j,k,nnew)* &
2117 & vmask_wet(iend+1,j)
2118# endif
2119# ifdef NEARSHORE_MELLOR
2120!^ v_stokes(Iend+1,j,k)=v_stokes(Iend+1,j,k)- &
2121!^ & CFs(Iend+1,0)
2122!^
2123 tl_v_stokes(iend+1,j,k)=tl_v_stokes(iend+1,j,k)- &
2124 & tl_cfs(iend+1,0)
2125# ifdef MASKING
2126!^ v_stokes(Iend+1,j,k)=v_stokes(Iend+1,j,k)* &
2127!^ & vmask(Iend+1,j)
2128!^
2129 tl_v_stokes(iend+1,j,k)=tl_v_stokes(iend+1,j,k)* &
2130 & vmask(iend+1,j)
2131# endif
2132# ifdef WET_DRY_NOT_YET
2133!^ v_stokes(Iend+1,j,k)=v_stokes(Iend+1,j,k)* &
2134!^ & vmask_wet(Iend+1,j)
2135!^
2136 tl_v_stokes(iend+1,j,k)=tl_v_stokes(iend+1,j,k)* &
2137 & vmask_wet(iend+1,j)
2138# endif
2139# endif
2140 END DO
2141 END IF
2142 END IF
2143!
2144 IF (.not.(compositegrid(isouth,ng).or.nsperiodic(ng))) THEN
2145 IF (j.eq.1) THEN ! southern boundary
2146 DO k=1,n(ng) ! J-loop pipelined
2147 DO i=istr,iend
2148!^ v(i,j,k,nnew)=v(i,j,k,nnew)-CF(i,0)
2149!^
2150 tl_v(i,j,k,nnew)=tl_v(i,j,k,nnew)- &
2151 & tl_cf(i,0)
2152# ifdef MASKING
2153!^ v(i,j,k,nnew)=v(i,j,k,nnew)*
2154!^ & vmask(i,j)
2155!^
2156 tl_v(i,j,k,nnew)=tl_v(i,j,k,nnew)* &
2157 & vmask(i,j)
2158# endif
2159# ifdef WET_DRY_NOT_YET
2160!^ v(i,j,k,nnew)=v(i,j,k,nnew)* &
2161!^ & vmask_wet(i,j)
2162!^
2163 tl_v(i,j,k,nnew)=tl_v(i,j,k,nnew)*
2164 & vmask_wet(i,j)
2165# endif
2166# ifdef NEARSHORE_MELLOR
2167!^ v_stokes(i,j,k)=v_stokes(i,j,k)-CFs(i,0)
2168!^
2169 tl_v_stokes(i,j,k)=tl_v_stokes(i,j,k)- &
2170 & tl_cfs(i,0)
2171# ifdef MASKING
2172!^ v_stokes(i,j,k)=v_stokes(i,j,k)*vmask(i,j)
2173!^
2174 tl_v_stokes(i,j,k)=tl_v_stokes(i,j,k)* &
2175 & vmask(i,j)
2176# endif
2177# ifdef WET_DRY_NOT_YET
2178!^ v_stokes(i,j,k)=v_stokes(i,j,k)* &
2179!^ & vmask_wet(i,j)
2180!^
2181 tl_v_stokes(i,j,k)=tl_v_stokes(i,j,k)* &
2182 & vmask_wet(i,j)
2183# endif
2184# endif
2185 END DO
2186 END DO
2187 END IF
2188 END IF
2189!
2190 IF (.not.(compositegrid(inorth,ng).or.nsperiodic(ng))) THEN
2191 IF (j.eq.mm(ng)+1) THEN ! northern boundary
2192 DO k=1,n(ng) ! J-loop pipelined
2193 DO i=istr,iend
2194!^ v(i,j,k,nnew)=v(i,j,k,nnew)-CF(i,0)
2195!^
2196 tl_v(i,j,k,nnew)=tl_v(i,j,k,nnew)- &
2197 & tl_cf(i,0)
2198# ifdef MASKING
2199!^ v(i,j,k,nnew)=v(i,j,k,nnew)* &
2200!^ & vmask(i,j)
2201!^
2202 tl_v(i,j,k,nnew)=tl_v(i,j,k,nnew)* &
2203 & vmask(i,j)
2204# endif
2205# ifdef WET_DRY_NOT_YET
2206!^ v(i,j,k,nnew)=v(i,j,k,nnew)* &
2207!^ & vmask_wet(i,j)
2208!^
2209 tl_v(i,j,k,nnew)=tl_v(i,j,k,nnew)* &
2210 & vmask_wet(i,j)
2211# endif
2212# ifdef NEARSHORE_MELLOR
2213!^ v_stokes(i,j,k)=v_stokes(i,j,k)-CFs(i,0)
2214!^
2215 tl_v_stokes(i,j,k)=tl_v_stokes(i,j,k)- &
2216 & tl_cfs(i,0)
2217# ifdef MASKING
2218!^ v_stokes(i,j,k)=v_stokes(i,j,k)* &
2219!^ & vmask(i,j)
2220!^
2221 tl_v_stokes(i,j,k)=tl_v_stokes(i,j,k)* &
2222 & vmask(i,j)
2223# endif
2224# ifdef WET_DRY_NOT_YET
2225!^ v_stokes(i,j,k)=v_stokes(i,j,k)* &
2226!^ & vmask_wet(i,j)
2227!^
2228 tl_v_stokes(i,j,k)=tl_v_stokes(i,j,k)* &
2229 & vmask_wet(i,j)
2230# endif
2231# endif
2232 END DO
2233 END DO
2234 END IF
2235 END IF
2236!
2237! Compute correct mass flux, Hz*v/m.
2238!
2239 DO k=n(ng),1,-1
2240 DO i=istrt,iendt
2241 hvom(i,j,k)=0.5_r8*(hvom(i,j,k)+v(i,j,k,nnew)*dc(i,k))
2242 tl_hvom(i,j,k)=0.5_r8*(tl_hvom(i,j,k)+ &
2243 & tl_v(i,j,k,nnew)*dc(i,k)+ &
2244 & v(i,j,k,nnew)*tl_dc(i,k))- &
2245# ifdef TL_IOMS
2246 & 0.5_r8*v(i,j,k,nnew)*dc(i,k)
2247# endif
2248# ifdef NEARSHORE_MELLOR
2249 hvom(i,j,k)=hvom(i,j,k)+0.5_r8*v_stokes(i,j,k)*dc(i,k)
2250 tl_hvom(i,j,k)=tl_hvom(i,j,k)+ &
2251 & 0.5_r8*(tl_v_stokes(i,j,k)*dc(i,k)+ &
2252 & v_stokes(i,j,k)*tl_dc(i,k))- &
2253# ifdef TL_IOMS
2254 & 0.5_r8*v_stokes(i,j,k)*dc(i,k)
2255# endif
2256# endif
2257 fc(i,0)=fc(i,0)+hvom(i,j,k)
2258 tl_fc(i,0)=tl_fc(i,0)+tl_hvom(i,j,k)
2259# ifdef DIAGNOSTICS_UV
2260!! DiaV3wrk(i,j,k,M3rate)=v(i,j,k,nnew)- &
2261!! & DiaV3wrk(i,j,k,M3rate)
2262# endif
2263 END DO
2264 END DO
2265 DO i=istrt,iendt
2266 fc1(i)=fc(i,0) ! intermediate
2267 fc(i,0)=dc(i,0)*(fc(i,0)-dv_avg2(i,j)) ! recursive
2268 tl_fc(i,0)=tl_dc(i,0)*(fc1(i)-dv_avg2(i,j))+ &
2269 & dc(i,0)*(tl_fc(i,0)-tl_dv_avg2(i,j))- &
2270# ifdef TL_IOMS
2271 & fc(i,0)
2272# endif
2273 END DO
2274 DO k=1,n(ng)
2275 DO i=istrt,iendt
2276 hvom(i,j,k)=hvom(i,j,k)-dc(i,k)*fc(i,0)
2277 tl_hvom(i,j,k)=tl_hvom(i,j,k)- &
2278 & tl_dc(i,k)*fc(i,0)- &
2279 & dc(i,k)*tl_fc(i,0)+ &
2280# ifdef TL_IOMS
2281 & dc(i,k)*fc(i,0)
2282# endif
2283 END DO
2284 END DO
2285 END IF
2286 END DO
2287!
2288!-----------------------------------------------------------------------
2289! Exchange boundary data.
2290!-----------------------------------------------------------------------
2291!
2292 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
2293!^ CALL exchange_u3d_tile (ng, tile, &
2294!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
2295!^ & u(:,:,:,nnew))
2296!^
2297 CALL exchange_u3d_tile (ng, tile, &
2298 & lbi, ubi, lbj, ubj, 1, n(ng), &
2299 & tl_u(:,:,:,nnew))
2300!^ CALL exchange_v3d_tile (ng, tile, &
2301!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
2302!^ & v(:,:,:,nnew))
2303!^
2304 CALL exchange_v3d_tile (ng, tile, &
2305 & lbi, ubi, lbj, ubj, 1, n(ng), &
2306 & tl_v(:,:,:,nnew))
2307
2308 CALL exchange_u3d_tile (ng, tile, &
2309 & lbi, ubi, lbj, ubj, 1, n(ng), &
2310 & huon)
2311 CALL exchange_u3d_tile (ng, tile, &
2312 & lbi, ubi, lbj, ubj, 1, n(ng), &
2313 & tl_huon)
2314
2315 CALL exchange_v3d_tile (ng, tile, &
2316 & lbi, ubi, lbj, ubj, 1, n(ng), &
2317 & hvom)
2318 CALL exchange_v3d_tile (ng, tile, &
2319 & lbi, ubi, lbj, ubj, 1, n(ng), &
2320 & tl_hvom)
2321!
2322 DO k=1,2
2323!^ CALL exchange_u2d_tile (ng, tile, &
2324!^ & LBi, UBi, LBj, UBj, &
2325!^ & ubar(:,:,k))
2326!^
2327 CALL exchange_u2d_tile (ng, tile, &
2328 & lbi, ubi, lbj, ubj, &
2329 & tl_ubar(:,:,k))
2330!^ CALL exchange_v2d_tile (ng, tile, &
2331!^ & LBi, UBi, LBj, UBj, &
2332!^ & vbar(:,:,k))
2333!^
2334 CALL exchange_v2d_tile (ng, tile, &
2335 & lbi, ubi, lbj, ubj, &
2336 & tl_vbar(:,:,k))
2337 END DO
2338 END IF
2339
2340# ifdef DISTRIBUTE
2341!^ CALL mp_exchange3d (ng, tile, iNLM, 4, &
2342!^ & LBi, UBi, LBj, UBj, 1, N(ng), &
2343!^ & NghostPoints, &
2344!^ & EWperiodic(ng), NSperiodic(ng), &
2345!^ & u(:,:,:,nnew), v(:,:,:,nnew), &
2346!^ & Huon, Hvom)
2347!^
2348 CALL mp_exchange3d (ng, tile, irpm, 4, &
2349 & lbi, ubi, lbj, ubj, 1, n(ng), &
2350 & nghostpoints, &
2351 & ewperiodic(ng), nsperiodic(ng), &
2352 & tl_u(:,:,:,nnew), tl_v(:,:,:,nnew), &
2353 & tl_huon, tl_hvom)
2354 CALL mp_exchange3d (ng, tile, irpm, 2, &
2355 & lbi, ubi, lbj, ubj, 1, n(ng), &
2356 & nghostpoints, &
2357 & ewperiodic(ng), nsperiodic(ng), &
2358 & huon, hvom)
2359!
2360!^ CALL mp_exchange2d (ng, tile, iNLM, 4, &
2361!^ & LBi, UBi, LBj, UBj, &
2362!^ & NghostPoints, &
2363!^ & EWperiodic(ng), NSperiodic(ng), &
2364!^ & ubar(:,:,1), vbar(:,:,1), &
2365!^ & ubar(:,:,2), vbar(:,:,2))
2366!^
2367 CALL mp_exchange2d (ng, tile, irpm, 4, &
2368 & lbi, ubi, lbj, ubj, &
2369 & nghostpoints, &
2370 & ewperiodic(ng), nsperiodic(ng), &
2371 & tl_ubar(:,:,1), tl_vbar(:,:,1), &
2372 & tl_ubar(:,:,2), tl_vbar(:,:,2))
2373# endif
2374# ifdef UV_DESTAGGERED
2375!
2376!-----------------------------------------------------------------------
2377! Compute representer 3D momentum (tl_ua, tl_va) at RHO-points (A-Grid)
2378! for output purposes and data assimilation where the observations and
2379! state vector is located at the cell-center.
2380!-----------------------------------------------------------------------
2381!
2382!> CALL uv_C2A_grid (ng, tile, iNLM, nnew)
2383!>
2384 CALL tl_uv_c2a_grid (ng, tile, irpm, nnew)
2385# endif
2386!
2387 RETURN
2388 END SUBROUTINE rp_step3d_uv_tile
2389#endif
2390 END MODULE rp_step3d_uv_mod
subroutine exchange_u2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_v2d_tile(ng, tile, lbi, ubi, lbj, ubj, a)
subroutine exchange_u3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
subroutine exchange_v3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
type(t_coupling), dimension(:), allocatable coupling
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
integer nghostpoints
Definition mod_param.F:710
type(t_domain), dimension(:), allocatable domain
Definition mod_param.F:329
integer, dimension(:), allocatable mm
Definition mod_param.F:456
logical, dimension(:), allocatable luvsrc
integer, dimension(:), allocatable iic
real(dp), dimension(:), allocatable dt
real(r8) lambda
logical, dimension(:), allocatable ewperiodic
integer, parameter iwest
logical, dimension(:), allocatable nsperiodic
logical, dimension(:,:), allocatable compositegrid
integer, parameter isouth
integer, dimension(:), allocatable ntfirst
integer, parameter ieast
integer, parameter inorth
type(t_sources), dimension(:), allocatable sources
Definition mod_sources.F:90
integer, dimension(:), allocatable nsrc
Definition mod_sources.F:97
integer, dimension(:), allocatable nrhs
integer, dimension(:), allocatable nnew
integer, dimension(:), allocatable nstp
subroutine mp_exchange2d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine mp_exchange3d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, nghost, ew_periodic, ns_periodic, a, b, c, d)
subroutine, public rp_step3d_uv(ng, tile)
subroutine rp_step3d_uv_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, nstp, nnew, umask, vmask, umask_wet, vmask_wet, om_v, on_u, pm, pn, hz, tl_hz, z_r, tl_z_r, z_w, tl_z_w, akv, tl_akv, du_avg1, tl_du_avg1, dv_avg1, tl_dv_avg1, du_avg2, tl_du_avg2, dv_avg2, tl_dv_avg2, tl_ru, tl_rv, u, tl_u, v, tl_v, tl_ubar, tl_vbar, ubar_stokes, tl_ubar_stokes, vbar_stokes, tl_vbar_stokes, u_stokes, tl_u_stokes, v_stokes, tl_v_stokes, huon, tl_huon, hvom, tl_hvom)
subroutine, public rp_u3dbc_tile(ng, tile, lbi, ubi, lbj, ubj, ubk, imins, imaxs, jmins, jmaxs, nstp, nout, tl_u)
Definition rp_u3dbc_im.F:58
subroutine, public rp_v3dbc_tile(ng, tile, lbi, ubi, lbj, ubj, ubk, imins, imaxs, jmins, jmaxs, nstp, nout, tl_v)
Definition rp_v3dbc_im.F:59
subroutine, public tl_uv_c2a_grid(ng, tile, model, ninp)
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