ROMS
Loading...
Searching...
No Matches
uv3dmix2_s.h
Go to the documentation of this file.
1 MODULE uv3dmix2_mod
2!
3!git $Id$
4!================================================== Hernan G. Arango ===
5! Copyright (c) 2002-2025 The ROMS Group !
6! Licensed under a MIT/X style license !
7! See License_ROMS.md !
8!=======================================================================
9! !
10! This routine computes harmonic mixing of momentum, along constant !
11! S-surfaces, from the horizontal divergence of the stress tensor. !
12! A transverse isotropy is assumed so the stress tensor is split !
13! into vertical and horizontal subtensors. !
14! !
15! Reference: !
16! !
17! Wajsowicz, R.C, 1993: A consistent formulation of the !
18! anisotropic stress tensor for use in models of the !
19! large-scale ocean circulation, JCP, 105, 333-338. !
20! !
21! Sadourny, R. and K. Maynard, 1997: Formulations of !
22! lateral diffusion in geophysical fluid dynamics !
23! models, In Numerical Methods of Atmospheric and !
24! Oceanic Modelling. Lin, Laprise, and Ritchie, !
25! Eds., NRC Research Press, 547-556. !
26! !
27! Griffies, S.M. and R.W. Hallberg, 2000: Biharmonic !
28! friction with a Smagorinsky-like viscosity for !
29! use in large-scale eddy-permitting ocean models, !
30! Monthly Weather Rev., 128, 8, 2935-2946. !
31! !
32!=======================================================================
33!
34 implicit none
35!
36 PRIVATE
37 PUBLIC uv3dmix2
38!
39 CONTAINS
40!
41!***********************************************************************
42 SUBROUTINE uv3dmix2 (ng, tile)
43!***********************************************************************
44!
45 USE mod_param
46 USE mod_coupling
47#ifdef DIAGNOSTICS_UV
48 USE mod_diags
49#endif
50 USE mod_grid
51 USE mod_mixing
52 USE mod_ocean
53 USE mod_stepping
54!
55! Imported variable declarations.
56!
57 integer, intent(in) :: ng, tile
58!
59! Local variable declarations.
60!
61 character (len=*), parameter :: MyFile = &
62 & __FILE__
63!
64#include "tile.h"
65!
66#ifdef PROFILE
67 CALL wclock_on (ng, inlm, 30, __line__, myfile)
68#endif
69 CALL uv3dmix2_s_tile (ng, tile, &
70 & lbi, ubi, lbj, ubj, &
71 & imins, imaxs, jmins, jmaxs, &
72 & nrhs(ng), nnew(ng), &
73#ifdef MASKING
74 & grid(ng) % pmask, &
75#endif
76#ifdef WET_DRY
77 & grid(ng) % pmask_wet, &
78#endif
79 & grid(ng) % Hz, &
80 & grid(ng) % om_p, &
81 & grid(ng) % om_r, &
82 & grid(ng) % on_p, &
83 & grid(ng) % on_r, &
84 & grid(ng) % pm, &
85 & grid(ng) % pmon_p, &
86 & grid(ng) % pmon_r, &
87 & grid(ng) % pn, &
88 & grid(ng) % pnom_p, &
89 & grid(ng) % pnom_r, &
90#ifdef VISC_3DCOEF
91 & mixing(ng) % visc3d_r, &
92#else
93 & mixing(ng) % visc2_p, &
94 & mixing(ng) % visc2_r, &
95#endif
96#ifdef DIAGNOSTICS_UV
97 & diags(ng) % DiaRUfrc, &
98 & diags(ng) % DiaRVfrc, &
99 & diags(ng) % DiaU3wrk, &
100 & diags(ng) % DiaV3wrk, &
101#endif
102 & coupling(ng) % rufrc, &
103 & coupling(ng) % rvfrc, &
104 & ocean(ng) % u, &
105 & ocean(ng) % v)
106#ifdef PROFILE
107 CALL wclock_off (ng, inlm, 30, __line__, myfile)
108#endif
109!
110 RETURN
111 END SUBROUTINE uv3dmix2
112!
113!***********************************************************************
114 SUBROUTINE uv3dmix2_s_tile (ng, tile, &
115 & LBi, UBi, LBj, UBj, &
116 & IminS, ImaxS, JminS, JmaxS, &
117 & nrhs, nnew, &
118#ifdef MASKING
119 & pmask, &
120#endif
121#ifdef WET_DRY
122 & pmask_wet, &
123#endif
124 & Hz, &
125 & om_p, om_r, on_p, on_r, &
126 & pm, pmon_p, pmon_r, &
127 & pn, pnom_p, pnom_r, &
128#ifdef VISC_3DCOEF
129 & visc3d_r, &
130#else
131 & visc2_p, visc2_r, &
132#endif
133#ifdef DIAGNOSTICS_UV
134 & DiaRUfrc, DiaRVfrc, &
135 & DiaU3wrk, DiaV3wrk, &
136#endif
137 & rufrc, rvfrc, u, v)
138!***********************************************************************
139!
140 USE mod_param
141 USE mod_scalars
142!
143! Imported variable declarations.
144!
145 integer, intent(in) :: ng, tile
146 integer, intent(in) :: LBi, UBi, LBj, UBj
147 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
148 integer, intent(in) :: nrhs, nnew
149
150#ifdef ASSUMED_SHAPE
151# ifdef MASKING
152 real(r8), intent(in) :: pmask(LBi:,LBj:)
153# endif
154# ifdef WET_DRY
155 real(r8), intent(in) :: pmask_wet(LBi:,LBj:)
156# endif
157 real(r8), intent(in) :: Hz(LBi:,LBj:,:)
158 real(r8), intent(in) :: om_p(LBi:,LBj:)
159 real(r8), intent(in) :: om_r(LBi:,LBj:)
160 real(r8), intent(in) :: on_p(LBi:,LBj:)
161 real(r8), intent(in) :: on_r(LBi:,LBj:)
162 real(r8), intent(in) :: pm(LBi:,LBj:)
163 real(r8), intent(in) :: pmon_p(LBi:,LBj:)
164 real(r8), intent(in) :: pmon_r(LBi:,LBj:)
165 real(r8), intent(in) :: pn(LBi:,LBj:)
166 real(r8), intent(in) :: pnom_p(LBi:,LBj:)
167 real(r8), intent(in) :: pnom_r(LBi:,LBj:)
168# ifdef VISC_3DCOEF
169 real(r8), intent(in) :: visc3d_r(LBi:,LBj:,:)
170# else
171 real(r8), intent(in) :: visc2_p(LBi:,LBj:)
172 real(r8), intent(in) :: visc2_r(LBi:,LBj:)
173# endif
174# ifdef DIAGNOSTICS_UV
175 real(r8), intent(inout) :: DiaRUfrc(LBi:,LBj:,:,:)
176 real(r8), intent(inout) :: DiaRVfrc(LBi:,LBj:,:,:)
177 real(r8), intent(inout) :: DiaU3wrk(LBi:,LBj:,:,:)
178 real(r8), intent(inout) :: DiaV3wrk(LBi:,LBj:,:,:)
179# endif
180 real(r8), intent(inout) :: rufrc(LBi:,LBj:)
181 real(r8), intent(inout) :: rvfrc(LBi:,LBj:)
182 real(r8), intent(inout) :: u(LBi:,LBj:,:,:)
183 real(r8), intent(inout) :: v(LBi:,LBj:,:,:)
184#else
185# ifdef MASKING
186 real(r8), intent(in) :: pmask(LBi:UBi,LBj:UBj)
187# endif
188# ifdef WET_DRY
189 real(r8), intent(in) :: pmask_wet(LBi:UBi,LBj:UBj)
190# endif
191 real(r8), intent(in) :: Hz(LBi:UBi,LBj:UBj,N(ng))
192 real(r8), intent(in) :: om_p(LBi:UBi,LBj:UBj)
193 real(r8), intent(in) :: om_r(LBi:UBi,LBj:UBj)
194 real(r8), intent(in) :: on_p(LBi:UBi,LBj:UBj)
195 real(r8), intent(in) :: on_r(LBi:UBi,LBj:UBj)
196 real(r8), intent(in) :: pm(LBi:UBi,LBj:UBj)
197 real(r8), intent(in) :: pmon_p(LBi:UBi,LBj:UBj)
198 real(r8), intent(in) :: pmon_r(LBi:UBi,LBj:UBj)
199 real(r8), intent(in) :: pn(LBi:UBi,LBj:UBj)
200 real(r8), intent(in) :: pnom_p(LBi:UBi,LBj:UBj)
201 real(r8), intent(in) :: pnom_r(LBi:UBi,LBj:UBj)
202# ifdef VISC_3DCOEF
203 real(r8), intent(in) :: visc3d_r(LBi:UBi,LBj:UBj,N(ng))
204# else
205 real(r8), intent(in) :: visc2_p(LBi:UBi,LBj:UBj)
206 real(r8), intent(in) :: visc2_r(LBi:UBi,LBj:UBj)
207# endif
208# ifdef DIAGNOSTICS_UV
209 real(r8), intent(inout) :: DiaRUfrc(LBi:UBi,LBj:UBj,3,NDM2d-1)
210 real(r8), intent(inout) :: DiaRVfrc(LBi:UBi,LBj:UBj,3,NDM2d-1)
211 real(r8), intent(inout) :: DiaU3wrk(LBi:UBi,LBj:UBj,N(ng),NDM3d)
212 real(r8), intent(inout) :: DiaV3wrk(LBi:UBi,LBj:UBj,N(ng),NDM3d)
213# endif
214 real(r8), intent(inout) :: rufrc(LBi:UBi,LBj:UBj)
215 real(r8), intent(inout) :: rvfrc(LBi:UBi,LBj:UBj)
216 real(r8), intent(inout) :: u(LBi:UBi,LBj:UBj,N(ng),2)
217 real(r8), intent(inout) :: v(LBi:UBi,LBj:UBj,N(ng),2)
218#endif
219!
220! Local variable declarations.
221!
222 integer :: i, j, k
223
224 real(r8) :: cff, cff1, cff2, cff3
225#ifdef VISC_3DCOEF
226 real(r8) :: visc_p
227#endif
228 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: UFe
229 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: VFe
230 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: UFx
231 real(r8), dimension(IminS:ImaxS,JminS:JmaxS) :: VFx
232
233#include "set_bounds.h"
234!
235!-----------------------------------------------------------------------
236! Compute horizontal harmonic viscosity along constant S-surfaces.
237!-----------------------------------------------------------------------
238!
239 k_loop : DO k=1,n(ng)
240!
241! Compute flux-components of the horizontal divergence of the stress
242! tensor (m5/s2) in XI- and ETA-directions.
243!
244 DO j=jstrv-1,jend
245 DO i=istru-1,iend
246 cff=hz(i,j,k)*0.5_r8* &
247 & (pmon_r(i,j)* &
248 & ((pn(i ,j)+pn(i+1,j))*u(i+1,j,k,nrhs)- &
249 & (pn(i-1,j)+pn(i ,j))*u(i ,j,k,nrhs))- &
250 & pnom_r(i,j)* &
251 & ((pm(i,j )+pm(i,j+1))*v(i,j+1,k,nrhs)- &
252 & (pm(i,j-1)+pm(i,j ))*v(i,j ,k,nrhs)))
253#ifdef VISC_3DCOEF
254 ufx(i,j)=on_r(i,j)*on_r(i,j)*visc3d_r(i,j,k)*cff
255 vfe(i,j)=om_r(i,j)*om_r(i,j)*visc3d_r(i,j,k)*cff
256#else
257 ufx(i,j)=on_r(i,j)*on_r(i,j)*visc2_r(i,j)*cff
258 vfe(i,j)=om_r(i,j)*om_r(i,j)*visc2_r(i,j)*cff
259#endif
260 END DO
261 END DO
262 DO j=jstr,jend+1
263 DO i=istr,iend+1
264 cff=0.125_r8*(hz(i-1,j ,k)+hz(i,j ,k)+ &
265 & hz(i-1,j-1,k)+hz(i,j-1,k))* &
266 & (pmon_p(i,j)* &
267 & ((pn(i ,j-1)+pn(i ,j))*v(i ,j,k,nrhs)- &
268 & (pn(i-1,j-1)+pn(i-1,j))*v(i-1,j,k,nrhs))+ &
269 & pnom_p(i,j)* &
270 & ((pm(i-1,j )+pm(i,j ))*u(i,j ,k,nrhs)- &
271 & (pm(i-1,j-1)+pm(i,j-1))*u(i,j-1,k,nrhs)))
272#ifdef MASKING
273 cff=cff*pmask(i,j)
274#endif
275#ifdef WET_DRY
276 cff=cff*pmask_wet(i,j)
277#endif
278#ifdef VISC_3DCOEF
279 visc_p=0.25_r8*(visc3d_r(i-1,j-1,k)+visc3d_r(i-1,j,k)+ &
280 & visc3d_r(i ,j-1,k)+visc3d_r(i ,j,k))
281 ufe(i,j)=om_p(i,j)*om_p(i,j)*visc_p*cff
282 vfx(i,j)=on_p(i,j)*on_p(i,j)*visc_p*cff
283#else
284 ufe(i,j)=om_p(i,j)*om_p(i,j)*visc2_p(i,j)*cff
285 vfx(i,j)=on_p(i,j)*on_p(i,j)*visc2_p(i,j)*cff
286#endif
287 END DO
288 END DO
289!
290! Time-step harmonic, S-surfaces viscosity term. Notice that momentum
291! at this stage is HzU and HzV and has m2/s units. Add contribution for
292! barotropic forcing terms.
293!
294 DO j=jstr,jend
295 DO i=istru,iend
296 cff=dt(ng)*0.25_r8*(pm(i-1,j)+pm(i,j))*(pn(i-1,j)+pn(i,j))
297 cff1=0.5_r8*(pn(i-1,j)+pn(i,j))*(ufx(i,j )-ufx(i-1,j))
298 cff2=0.5_r8*(pm(i-1,j)+pm(i,j))*(ufe(i,j+1)-ufe(i ,j))
299 cff3=cff*(cff1+cff2)
300 rufrc(i,j)=rufrc(i,j)+cff1+cff2
301 u(i,j,k,nnew)=u(i,j,k,nnew)+cff3
302#ifdef DIAGNOSTICS_UV
303 diarufrc(i,j,3,m2hvis)=diarufrc(i,j,3,m2hvis)+cff1+cff2
304 diarufrc(i,j,3,m2xvis)=diarufrc(i,j,3,m2xvis)+cff1
305 diarufrc(i,j,3,m2yvis)=diarufrc(i,j,3,m2yvis)+cff2
306 diau3wrk(i,j,k,m3hvis)=cff3
307 diau3wrk(i,j,k,m3xvis)=cff*cff1
308 diau3wrk(i,j,k,m3yvis)=cff*cff2
309#endif
310 END DO
311 END DO
312 DO j=jstrv,jend
313 DO i=istr,iend
314 cff=dt(ng)*0.25_r8*(pm(i,j)+pm(i,j-1))*(pn(i,j)+pn(i,j-1))
315 cff1=0.5_r8*(pn(i,j-1)+pn(i,j))*(vfx(i+1,j)-vfx(i,j ))
316 cff2=0.5_r8*(pm(i,j-1)+pm(i,j))*(vfe(i ,j)-vfe(i,j-1))
317 cff3=cff*(cff1-cff2)
318 rvfrc(i,j)=rvfrc(i,j)+cff1-cff2
319 v(i,j,k,nnew)=v(i,j,k,nnew)+cff3
320#ifdef DIAGNOSTICS_UV
321 diarvfrc(i,j,3,m2hvis)=diarvfrc(i,j,3,m2hvis)+cff1-cff2
322 diarvfrc(i,j,3,m2xvis)=diarvfrc(i,j,3,m2xvis)+cff1
323 diarvfrc(i,j,3,m2yvis)=diarvfrc(i,j,3,m2yvis)-cff2
324 diav3wrk(i,j,k,m3hvis)=cff3
325 diav3wrk(i,j,k,m3xvis)= cff*cff1
326 diav3wrk(i,j,k,m3yvis)=-cff*cff2
327#endif
328 END DO
329 END DO
330 END DO k_loop
331!
332 RETURN
333 END SUBROUTINE uv3dmix2_s_tile
334
335 END MODULE uv3dmix2_mod
type(t_coupling), dimension(:), allocatable coupling
type(t_diags), dimension(:), allocatable diags
Definition mod_diags.F:100
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_mixing), dimension(:), allocatable mixing
Definition mod_mixing.F:399
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
integer, parameter inlm
Definition mod_param.F:662
real(dp), dimension(:), allocatable dt
integer m3xvis
integer m2hvis
integer m3hvis
integer m3yvis
integer m2yvis
integer m2xvis
integer, dimension(:), allocatable nrhs
integer, dimension(:), allocatable nnew
subroutine uv3dmix2_s_tile(ng, tile, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, nrhs, nnew, pmask, pmask_wet, hz, om_p, om_r, on_p, on_r, pm, pmon_p, pmon_r, pn, pnom_p, pnom_r, visc3d_r, visc2_p, visc2_r, diarufrc, diarvfrc, diau3wrk, diav3wrk, rufrc, rvfrc, u, v)
Definition uv3dmix2_s.h:138
subroutine, public uv3dmix2(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