ROMS
Loading...
Searching...
No Matches
ana_vmix.h
Go to the documentation of this file.
1!!
2 SUBROUTINE ana_vmix (ng, tile, model)
3!
4!! git $Id$
5!!======================================================================
6!! Copyright (c) 2002-2025 The ROMS Group !
7!! Licensed under a MIT/X style license !
8!! See License_ROMS.md !
9!=======================================================================
10! !
11! This routine sets vertical mixing coefficients for momentum "Akv" !
12! and tracers "Akt" (m2/s) using analytical expressions. !
13! !
14!=======================================================================
15!
16 USE mod_param
17 USE mod_grid
18 USE mod_mixing
19 USE mod_ncparam
20 USE mod_ocean
21 USE mod_stepping
22!
23! Imported variable declarations.
24!
25 integer, intent(in) :: ng, tile, model
26!
27! Local variable declarations.
28!
29 character (len=*), parameter :: MyFile = &
30 & __FILE__
31!
32#include "tile.h"
33!
34 CALL ana_vmix_tile (ng, tile, model, &
35 & lbi, ubi, lbj, ubj, &
36 & imins, imaxs, jmins, jmaxs, &
37 & knew(ng), &
38 & grid(ng) % h, &
39 & grid(ng) % z_r, &
40 & grid(ng) % z_w, &
41 & ocean(ng) % zeta, &
42 & mixing(ng) % Akv, &
43 & mixing(ng) % Akt)
44!
45! Set analytical header file name used.
46!
47#ifdef DISTRIBUTE
48 IF (lanafile) THEN
49#else
50 IF (lanafile.and.(tile.eq.0)) THEN
51#endif
52 ananame(35)=myfile
53 END IF
54!
55 RETURN
56 END SUBROUTINE ana_vmix
57!
58!***********************************************************************
59 SUBROUTINE ana_vmix_tile (ng, tile, model, &
60 & LBi, UBi, LBj, UBj, &
61 & IminS, ImaxS, JminS, JmaxS, &
62 & knew, &
63 & h, z_r, z_w, zeta, Akv, Akt)
64!***********************************************************************
65!
66 USE mod_param
67 USE mod_scalars
68!
70#ifdef DISTRIBUTE
72#endif
73!
74! Imported variable declarations.
75!
76 integer, intent(in) :: ng, tile, model
77 integer, intent(in) :: LBi, UBi, LBj, UBj
78 integer, intent(in) :: IminS, ImaxS, JminS, JmaxS
79 integer, intent(in) :: knew
80!
81#ifdef ASSUMED_SHAPE
82 real(r8), intent(in) :: h(LBi:,LBj:)
83 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
84 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
85 real(r8), intent(in) :: zeta(LBi:,LBj:,:)
86 real(r8), intent(out) :: Akv(LBi:,LBj:,0:)
87 real(r8), intent(out) :: Akt(LBi:,LBj:,0:,:)
88#else
89 real(r8), intent(in) :: h(LBi:UBi,LBj:UBj)
90 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
91 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
92 real(r8), intent(in) :: zeta(LBi:UBi,LBj:UBj,:)
93 real(r8), intent(out) :: Akv(LBi:UBi,LBj:UBj,0:N(ng))
94 real(r8), intent(out) :: Akt(LBi:UBi,LBj:UBj,0:N(ng),NAT)
95#endif
96!
97! Local variable declarations.
98!
99 integer :: i, itrc, j, k
100
101#include "set_bounds.h"
102!
103!-----------------------------------------------------------------------
104! Set vertical viscosity coefficient (m2/s).
105!-----------------------------------------------------------------------
106!
107#if defined CANYON
108 DO k=1,n(ng)-1
109 DO j=jstrt,jendt
110 DO i=istrt,iendt
111 akv(i,j,k)=1.0e-03_r8+95.0e-04_r8*exp(z_w(i,j,k)/50.0_r8)+ &
112 & 95.0e-04_r8*exp(-(z_w(i,j,k)+h(i,j))/50.0_r8)
113 END DO
114 END DO
115 END DO
116#elif defined CHANNEL_NECK
117 DO k=1,n(ng)-1
118 DO j=jstrt,jendt
119 DO i=istrt,iendt
120 akv(i,j,k)=2.0e-04_r8+8.0e-04_r8*exp(z_w(i,j,k)/5.0_r8)
121 END DO
122 END DO
123 END DO
124#elif defined COUPLING_TEST
125 DO k=1,n(ng)-1
126 DO j=jstrt,jendt
127 DO i=istrt,iendt
128 akv(i,j,k)=2.0e-03_r8+8.0e-03_r8*exp(z_w(i,j,k)/1500.0_r8)
129 END DO
130 END DO
131 END DO
132#elif defined ESTUARY_TEST
133 DO k=1,n(ng)-1
134 DO j=jstrt,jendt
135 DO i=istrt,iendt
136 akv(i,j,k)=0.002_r8
137 END DO
138 END DO
139 END DO
140#elif defined LAKE_SIGNELL
141 DO k=1,n(ng)-1
142 DO j=jstrt,jendt
143 DO i=istrt,iendt
144 akv(i,j,k)=0.0005_r8
145 END DO
146 END DO
147 END DO
148#elif defined NJ_BIGHT
149 DO k=1,n(ng)-1
150 DO j=jstrt,jendt
151 DO i=istrt,iendt
152 akv(i,j,k)=1.0e-03_r8+2.0e-04_r8*exp(z_r(i,j,k)/10.0_r8)
153 END DO
154 END DO
155 END DO
156#elif defined SED_TEST1
157 DO k=1,n(ng)-1 ! vonkar*ustar*z*(1-z/D)
158 DO j=jstrt,jendt
159 DO i=istrt,iendt
160 akv(i,j,k)=0.025_r8*(h(i,j)+z_w(i,j,k))* &
161 & (1.0_r8-(h(i,j)+z_w(i,j,k))/ &
162 & (h(i,j)+zeta(i,j,knew)))
163 akt(i,j,k,itemp)=akv(i,j,k)*0.49_r8/0.39_r8
164# ifdef SALINITY
165 akt(i,j,k,isalt)=akt(i,j,k,itemp)
166# endif
167 END DO
168 END DO
169 END DO
170#elif defined SED_TOY
171 DO k=1,n(ng)-1 ! vonkar*ustar*z*(1-z/D)
172 DO j=jstrt,jendt
173 DO i=istrt,iendt
174 akv(i,j,k)=0.41_r8*0.01_r8*(h(i,j)+z_w(i,j,k))* &
175 & (1.0_r8-(h(i,j)+z_w(i,j,k))/ &
176 & (h(i,j)+zeta(i,j,knew)))
177 END DO
178 END DO
179 END DO
180#elif defined SHOREFACE
181 DO k=1,n(ng)-1
182 DO j=jstrt,jendt
183 DO i=istrt,iendt
184 akv(i,j,k)=0.025_r8*(h(i,j)+z_w(i,j,k))* &
185 & (1.0_r8-(h(i,j)+z_w(i,j,k))/ &
186 & (h(i,j)+zeta(i,j,knew)))
187 END DO
188 END DO
189 END DO
190#elif defined TEST_CHAN
191 DO k=1,n(ng)-1 ! vonkar*ustar*z*(1-z/D)
192 DO j=jstrt,jendt
193 DO i=istrt,iendt
194 akv(i,j,k)=0.41_r8*0.0625_r8*(h(i,j)+z_w(i,j,k))* &
195 & (1.0_r8-(h(i,j)+z_w(i,j,k))/ &
196 & (h(i,j)+zeta(i,j,knew)))
197 END DO
198 END DO
199 END DO
200#elif defined UPWELLING
201 DO k=1,n(ng)-1
202 DO j=jstrt,jendt
203 DO i=istrt,iendt
204 akv(i,j,k)=2.0e-03_r8+8.0e-03_r8*exp(z_w(i,j,k)/150.0_r8)
205 END DO
206 END DO
207 END DO
208#else
209 ana_vmix.h: no values provided for akv.
210#endif
211!
212! Exchange boundary data.
213!
214 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
215 CALL exchange_w3d_tile (ng, tile, &
216 & lbi, ubi, lbj, ubj, 0, n(ng), &
217 & akv)
218 END IF
219
220#ifdef DISTRIBUTE
221 CALL mp_exchange3d (ng, tile, model, 1, &
222 & lbi, ubi, lbj, ubj, 0, n(ng), &
223 & nghostpoints, &
224 & ewperiodic(ng), nsperiodic(ng), &
225 & akv)
226#endif
227!
228!-----------------------------------------------------------------------
229! Set vertical diffusion coefficient (m2/s).
230!-----------------------------------------------------------------------
231!
232#if defined CANYON
233 DO k=1,n(ng)-1
234 DO j=jstrt,jendt
235 DO i=istrt,iendt
236 akt(i,j,k,itemp)=akt_bak(itemp,ng)
237 END DO
238 END DO
239 END DO
240#elif defined CHANNEL_NECK
241 DO k=1,n(ng)-1
242 DO j=jstrt,jendt
243 DO i=istrt,iendt
244 akt(i,j,k,itemp)=2.0e-06_r8+ &
245 & 8.0e-06_r8*exp(z_w(i,j,k)/5.0_r8)
246 END DO
247 END DO
248 END DO
249#elif defined COUPLING_TEST
250 DO k=1,n(ng)-1
251 DO j=jstrt,jendt
252 DO i=istrt,iendt
253 akt(i,j,k,itemp)=akt_bak(itemp,ng)
254# ifdef SALINITY
255 akt(i,j,k,isalt)=akt_bak(isalt,ng)
256# endif
257 END DO
258 END DO
259 END DO
260#elif defined ESTUARY_TEST
261 DO k=1,n(ng)-1
262 DO j=jstrt,jendt
263 DO i=istrt,iendt
264 akt(i,j,k,itemp)=akv(i,j,k)
265# ifdef SALINITY
266 akt(i,j,k,isalt)=akv(i,j,k)
267# endif
268 END DO
269 END DO
270 END DO
271#elif defined LAKE_SIGNELL
272 DO k=1,n(ng)-1
273 DO j=jstrt,jendt
274 DO i=istrt,iendt
275 akt(i,j,k,itemp)=akt_bak(itemp,ng)
276# ifdef SALINITY
277 akt(i,j,k,isalt)=akt_bak(isalt,ng)
278# endif
279 END DO
280 END DO
281 END DO
282#elif defined NJ_BIGHT
283 DO k=1,n(ng)-1
284 DO j=jstrt,jendt
285 DO i=istrt,iendt
286 akt(i,j,k,itemp)=1.0e-05_r8+ &
287 & 2.0e-06_r8*exp(z_r(i,j,k)/10.0_r8)
288# ifdef SALINITY
289 akt(i,j,k,isalt)=akt(i,j,k,itemp)
290# endif
291 END DO
292 END DO
293 END DO
294#elif defined SED_TOY
295 DO k=1,n(ng)-1 ! vonkar*ustar*z*(1-z/D)
296 DO j=jstrt,jendt
297 DO i=istrt,iendt
298 akt(i,j,k,itemp)=akv(i,j,k)
299# ifdef SALINITY
300 akt(i,j,k,isalt)=akv(i,j,k)
301# endif
302 END DO
303 END DO
304 END DO
305#elif defined SHOREFACE
306 DO k=1,n(ng)-1
307 DO j=jstrt,jendt
308 DO i=istrt,iendt
309 akt(i,j,k,itemp)=akv(i,j,k)
310# ifdef SALINITY
311 akt(i,j,k,isalt)=akv(i,j,k)
312# endif
313 END DO
314 END DO
315 END DO
316#elif defined TEST_CHAN
317 DO k=1,n(ng)-1
318 DO j=jstrt,jendt
319 DO i=istrt,iendt
320 akt(i,j,k,itemp)=akv(i,j,k)*0.49_r8/0.39_r8
321# ifdef SALINITY
322 akt(i,j,k,isalt)=akt(i,j,k,itemp)
323# endif
324 END DO
325 END DO
326 END DO
327#elif defined UPWELLING
328 DO k=1,n(ng)-1
329 DO j=jstrt,jendt
330 DO i=istrt,iendt
331 akt(i,j,k,itemp)=akt_bak(itemp,ng)
332# ifdef SALINITY
333 akt(i,j,k,isalt)=akt_bak(isalt,ng)
334# endif
335 END DO
336 END DO
337 END DO
338#else
339 ana_vmix.h: no values provided for akt.
340#endif
341!
342! Exchange boundary data.
343!
344 IF (ewperiodic(ng).or.nsperiodic(ng)) THEN
345 DO itrc=1,nat
346 CALL exchange_w3d_tile (ng, tile, &
347 & lbi, ubi, lbj, ubj, 0, n(ng), &
348 & akt(:,:,:,itrc))
349 END DO
350 END IF
351
352#ifdef DISTRIBUTE
353 CALL mp_exchange4d (ng, tile, model, 1, &
354 & lbi, ubi, lbj, ubj, 0, n(ng), 1, nat, &
355 & nghostpoints, &
356 & ewperiodic(ng), nsperiodic(ng), &
357 & akt)
358#endif
359!
360 RETURN
361 END SUBROUTINE ana_vmix_tile
subroutine ana_vmix(ng, tile, model)
Definition ana_vmix.h:3
subroutine ana_vmix_tile(ng, tile, model, lbi, ubi, lbj, ubj, imins, imaxs, jmins, jmaxs, knew, h, z_r, z_w, zeta, akv, akt)
Definition ana_vmix.h:64
subroutine exchange_w3d_tile(ng, tile, lbi, ubi, lbj, ubj, lbk, ubk, a)
type(t_grid), dimension(:), allocatable grid
Definition mod_grid.F:365
type(t_mixing), dimension(:), allocatable mixing
Definition mod_mixing.F:399
logical lanafile
character(len=256), dimension(39) ananame
type(t_ocean), dimension(:), allocatable ocean
Definition mod_ocean.F:351
integer nghostpoints
Definition mod_param.F:710
logical, dimension(:), allocatable ewperiodic
logical, dimension(:), allocatable nsperiodic
integer isalt
real(r8), dimension(:,:), allocatable akt_bak
integer itemp
integer, dimension(:), allocatable knew
subroutine mp_exchange4d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, lbt, ubt, nghost, ew_periodic, ns_periodic, a, b, c)
subroutine mp_exchange3d(ng, tile, model, nvar, lbi, ubi, lbj, ubj, lbk, ubk, nghost, ew_periodic, ns_periodic, a, b, c, d)