Ticket #33: ana_vmix.h

File ana_vmix.h, 10.8 KB (added by m.hadfield, 17 years ago)
Line 
1 SUBROUTINE ana_vmix (ng, tile, model)
2!
3!! svn $Id: ana_vmix.h 34 2007-04-27 04:40:21Z arango $
4!!======================================================================
5!! Copyright (c) 2002-2007 The ROMS/TOMS Group !
6!! Licensed under a MIT/X style license !
7!! See License_ROMS.txt !
8!! !
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#include "tile.h"
28!
29 CALL ana_vmix_tile (ng, model, Istr, Iend, Jstr, Jend, &
30 & LBi, UBi, LBj, UBj, &
31 & knew(ng), &
32 & GRID(ng) % h, &
33 & GRID(ng) % z_r, &
34 & GRID(ng) % z_w, &
35 & OCEAN(ng) % zeta, &
36 & MIXING(ng) % Akv, &
37 & MIXING(ng) % Akt)
38!
39! Set analytical header file name used.
40!
41 IF (Lanafile) THEN
42 ANANAME(35)='ROMS/Functionals/ana_vmix.h'
43 END IF
44
45 RETURN
46 END SUBROUTINE ana_vmix
47!
48!***********************************************************************
49 SUBROUTINE ana_vmix_tile (ng, model, Istr, Iend, Jstr, Jend, &
50 & LBi, UBi, LBj, UBj, &
51 & knew, &
52 & h, z_r, z_w, zeta, Akv, Akt)
53!***********************************************************************
54!
55 USE mod_param
56 USE mod_scalars
57!
58#if defined EW_PERIODIC || defined NS_PERIODIC
59 USE exchange_3d_mod, ONLY : exchange_w3d_tile
60#endif
61#ifdef DISTRIBUTE
62 USE mp_exchange_mod, ONLY : mp_exchange3d, mp_exchange4d
63#endif
64!
65! Imported variable declarations.
66!
67 integer, intent(in) :: ng, model, Iend, Istr, Jend, Jstr
68 integer, intent(in) :: LBi, UBi, LBj, UBj
69 integer, intent(in) :: knew
70!
71#ifdef ASSUMED_SHAPE
72 real(r8), intent(in) :: h(LBi:,LBj:)
73 real(r8), intent(in) :: z_r(LBi:,LBj:,:)
74 real(r8), intent(in) :: z_w(LBi:,LBj:,0:)
75 real(r8), intent(in) :: zeta(LBi:,LBj:,:)
76 real(r8), intent(out) :: Akv(LBi:,LBj:,0:)
77 real(r8), intent(out) :: Akt(LBi:,LBj:,0:,:)
78#else
79 real(r8), intent(in) :: h(LBi:UBi,LBj:UBj)
80 real(r8), intent(in) :: z_r(LBi:UBi,LBj:UBj,N(ng))
81 real(r8), intent(in) :: z_w(LBi:UBi,LBj:UBj,0:N(ng))
82 real(r8), intent(in) :: zeta(LBi:UBi,LBj:UBj,3)
83 real(r8), intent(out) :: Akv(LBi:UBi,LBj:UBj,0:N(ng))
84 real(r8), intent(out) :: Akt(LBi:UBi,LBj:UBj,0:N(ng),NAT)
85#endif
86!
87! Local variable declarations.
88!
89#ifdef DISTRIBUTE
90# ifdef EW_PERIODIC
91 logical :: EWperiodic=.TRUE.
92# else
93 logical :: EWperiodic=.FALSE.
94# endif
95# ifdef NS_PERIODIC
96 logical :: NSperiodic=.TRUE.
97# else
98 logical :: NSperiodic=.FALSE.
99# endif
100#endif
101 integer :: IstrR, IendR, JstrR, JendR, IstrU, JstrV
102 integer :: i, itrc, j, k
103
104#include "set_bounds.h"
105!
106!-----------------------------------------------------------------------
107! Set vertical viscosity coefficient (m2/s).
108!-----------------------------------------------------------------------
109!
110#if defined CANYON
111 DO k=1,N(ng)-1
112 DO j=JstrR,JendR
113 DO i=IstrR,IendR
114 Akv(i,j,k)=1.0E-03_r8+95.0E-04_r8*EXP(z_w(i,j,k)/50.0_r8)+ &
115 & 95.0E-04_r8*EXP(-(z_w(i,j,k)+h(i,j))/50.0_r8)
116 END DO
117 END DO
118 END DO
119#elif defined CHANNEL_NECK
120 DO k=1,N(ng)-1
121 DO j=JstrR,JendR
122 DO i=IstrR,IendR
123 Akv(i,j,k)=2.0E-04_r8+8.0E-04_r8*EXP(z_w(i,j,k)/5.0_r8)
124 END DO
125 END DO
126 END DO
127#elif defined COUPLING_TEST
128 DO k=1,N(ng)-1
129 DO j=JstrR,JendR
130 DO i=IstrR,IendR
131 Akv(i,j,k)=2.0E-03_r8+8.0E-03_r8*EXP(z_w(i,j,k)/1500.0_r8)
132 END DO
133 END DO
134 END DO
135#elif defined ESTUARY_TEST
136 DO k=1,N(ng)-1
137 DO j=JstrR,JendR
138 DO i=IstrR,IendR
139 Akv(i,j,k)=0.002_r8
140 END DO
141 END DO
142 END DO
143#elif defined LAKE_SIGNELL
144 DO k=1,N(ng)-1
145 DO j=JstrR,JendR
146 DO i=IstrR,IendR
147 Akv(i,j,k)=0.0005_r8
148 END DO
149 END DO
150 END DO
151#elif defined NJ_BIGHT
152 DO k=1,N(ng)-1
153 DO j=JstrR,JendR
154 DO i=IstrR,IendR
155 Akv(i,j,k)=1.0E-03_r8+2.0E-04_r8*EXP(z_r(i,j,k)/10.0_r8)
156 END DO
157 END DO
158 END DO
159#elif defined SED_TEST1
160 DO k=1,N(ng)-1 ! vonkar*ustar*z*(1-z/D)
161 DO j=JstrR,JendR
162 DO i=IstrR,IendR
163 Akv(i,j,k)=0.025_r8*(h(i,j)+z_w(i,j,k))* &
164 & (1.0_r8-(h(i,j)+z_w(i,j,k))/ &
165 & (h(i,j)+zeta(i,j,knew)))
166 Akt(i,j,k,itemp)=Akv(i,j,k)*0.49_r8/0.39_r8
167 Akt(i,j,k,isalt)=Akt(i,j,k,itemp)
168 END DO
169 END DO
170 END DO
171#elif defined SED_TOY
172 DO k=1,N(ng)-1 ! vonkar*ustar*z*(1-z/D)
173 DO j=JstrR,JendR
174 DO i=IstrR,IendR
175 Akv(i,j,k)=0.41_r8*0.01_r8*(h(i,j)+z_w(i,j,k))* &
176 & (1.0_r8-(h(i,j)+z_w(i,j,k))/ &
177 & (h(i,j)+zeta(i,j,knew)))
178 END DO
179 END DO
180 END DO
181#elif defined SHOREFACE
182 DO k=1,N(ng)-1
183 DO j=JstrR,JendR
184 DO i=IstrR,IendR
185 Akv(i,j,k)=0.025_r8*(h(i,j)+z_w(i,j,k))* &
186 & (1.0_r8-(h(i,j)+z_w(i,j,k))/ &
187 & (h(i,j)+zeta(i,j,knew)))
188 END DO
189 END DO
190 END DO
191#elif defined TEST_CHAN
192 DO k=1,N(ng)-1 ! vonkar*ustar*z*(1-z/D)
193 DO j=JstrR,JendR
194 DO i=IstrR,IendR
195 Akv(i,j,k)=0.41_r8*0.0625_r8*(h(i,j)+z_w(i,j,k))* &
196 & (1.0_r8-(h(i,j)+z_w(i,j,k))/ &
197 & (h(i,j)+zeta(i,j,knew)))
198 END DO
199 END DO
200 END DO
201#elif defined UPWELLING
202 DO k=1,N(ng)-1
203 DO j=JstrR,JendR
204 DO i=IstrR,IendR
205 Akv(i,j,k)=2.0E-03_r8+8.0E-03_r8*EXP(z_w(i,j,k)/150.0_r8)
206 END DO
207 END DO
208 END DO
209#else
210 ana_vmix.h: no values provided for Akv.
211#endif
212#if defined EW_PERIODIC || defined NS_PERIODIC
213 CALL exchange_w3d_tile (ng, Istr, Iend, Jstr, Jend, &
214 & LBi, UBi, LBj, UBj, 0, N(ng), &
215 & Akv)
216#endif
217#ifdef DISTRIBUTE
218 CALL mp_exchange3d (ng, model, 1, Istr, Iend, Jstr, Jend, &
219 & LBi, UBi, LBj, UBj, 0, N(ng), &
220 & NghostPoints, EWperiodic, NSperiodic, &
221 & Akv)
222#endif
223!
224!-----------------------------------------------------------------------
225! Set vertical diffusion coefficient (m2/s).
226!-----------------------------------------------------------------------
227!
228#if defined CANYON
229 DO k=1,N(ng)-1
230 DO j=JstrR,JendR
231 DO i=IstrR,IendR
232 Akt(i,j,k,itemp)=Akt_bak(itemp,ng)
233 END DO
234 END DO
235 END DO
236#elif defined CHANNEL_NECK
237 DO k=1,N(ng)-1
238 DO j=JstrR,JendR
239 DO i=IstrR,IendR
240 Akt(i,j,k,itemp)=2.0E-06_r8+ &
241 & 8.0E-06_r8*EXP(z_w(i,j,k)/5.0_r8)
242 END DO
243 END DO
244 END DO
245#elif defined COUPLING_TEST
246 DO k=1,N(ng)-1
247 DO j=JstrR,JendR
248 DO i=IstrR,IendR
249 Akt(i,j,k,itemp)=Akt_bak(itemp,ng)
250 Akt(i,j,k,isalt)=Akt_bak(isalt,ng)
251 END DO
252 END DO
253 END DO
254#elif defined ESTUARY_TEST
255 DO k=1,N(ng)-1
256 DO j=JstrR,JendR
257 DO i=IstrR,IendR
258 Akt(i,j,k,itemp)=Akv(i,j,k)
259 Akt(i,j,k,isalt)=Akv(i,j,k)
260 END DO
261 END DO
262 END DO
263#elif defined LAKE_SIGNELL
264 DO k=1,N(ng)-1
265 DO j=JstrR,JendR
266 DO i=IstrR,IendR
267 Akt(i,j,k,itemp)=Akt_bak(itemp,ng)
268 Akt(i,j,k,isalt)=Akt_bak(isalt,ng)
269 END DO
270 END DO
271 END DO
272#elif defined NJ_BIGHT
273 DO k=1,N(ng)-1
274 DO j=JstrR,JendR
275 DO i=IstrR,IendR
276 Akt(i,j,k,itemp)=1.0E-05_r8+ &
277 & 2.0E-06_r8*EXP(z_r(i,j,k)/10.0_r8)
278 Akt(i,j,k,isalt)=Akt(i,j,k,itemp)
279 END DO
280 END DO
281 END DO
282#elif defined SED_TOY
283 DO k=1,N(ng)-1 ! vonkar*ustar*z*(1-z/D)
284 DO j=JstrR,JendR
285 DO i=IstrR,IendR
286 Akt(i,j,k,itemp)=Akv(i,j,k)
287 Akt(i,j,k,isalt)=Akv(i,j,k)
288 END DO
289 END DO
290 END DO
291#elif defined SHOREFACE
292 DO k=1,N(ng)-1
293 DO j=JstrR,JendR
294 DO i=IstrR,IendR
295 Akt(i,j,k,itemp)=Akv(i,j,k)
296 Akt(i,j,k,isalt)=Akv(i,j,k)
297 END DO
298 END DO
299 END DO
300#elif defined TEST_CHAN
301 DO k=1,N(ng)-1
302 DO j=JstrR,JendR
303 DO i=IstrR,IendR
304 Akt(i,j,k,itemp)=Akv(i,j,k)*0.49_r8/0.39_r8
305 Akt(i,j,k,isalt)=Akt(i,j,k,itemp)
306 END DO
307 END DO
308 END DO
309#elif defined UPWELLING
310 DO k=1,N(ng)-1
311 DO j=JstrR,JendR
312 DO i=IstrR,IendR
313 Akt(i,j,k,itemp)=Akt_bak(itemp,ng)
314 Akt(i,j,k,isalt)=Akt_bak(isalt,ng)
315 END DO
316 END DO
317 END DO
318#else
319 ana_vmix.h: no values provided for Akt.
320#endif
321#if defined EW_PERIODIC || defined NS_PERIODIC
322 DO itrc=1,NAT
323 CALL exchange_w3d_tile (ng, Istr, Iend, Jstr, Jend, &
324 & LBi, UBi, LBj, UBj, 0, N(ng), &
325 & Akt(:,:,:,itrc))
326 END DO
327#endif
328#ifdef DISTRIBUTE
329 CALL mp_exchange4d (ng, model, 1, Istr, Iend, Jstr, Jend, &
330 & LBi, UBi, LBj, UBj, 0, N(ng), 1, NAT, &
331 & NghostPoints, EWperiodic, NSperiodic, &
332 & Akt)
333#endif
334 RETURN
335 END SUBROUTINE ana_vmix_tile