Home

Context Navigation

source: palm/trunk/SOURCE/diffusion_v.f90 @ 1295

Last change on this file since 1295 was 1258, checked in by raasch, 11 years ago
last commit documented
Property svn:keywords set to `Id`
File size: 24.3 KB

Line
1	MODULE diffusion_v_mod
2
3	!--------------------------------------------------------------------------------!
4	! This file is part of PALM.
5	!
6	! PALM is free software: you can redistribute it and/or modify it under the terms
7	! of the GNU General Public License as published by the Free Software Foundation,
8	! either version 3 of the License, or (at your option) any later version.
9	!
10	! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11	! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12	! A PARTICULAR PURPOSE. See the GNU General Public License for more details.
13	!
14	! You should have received a copy of the GNU General Public License along with
15	! PALM. If not, see <http://www.gnu.org/licenses/>.
16	!
17	! Copyright 1997-2012 Leibniz University Hannover
18	!--------------------------------------------------------------------------------!
19	!
20	! Current revisions:
21	! -----------------
22	!
23	!
24	! Former revisions:
25	! -----------------
26	! $Id: diffusion_v.f90 1258 2013-11-08 16:09:09Z raasch $
27	!
28	! 1257 2013-11-08 15:18:40Z raasch
29	! openacc loop and loop vector clauses removed, declare create moved after
30	! the FORTRAN declaration statement
31	!
32	! 1128 2013-04-12 06:19:32Z raasch
33	! loop index bounds in accelerator version replaced by i_left, i_right, j_south,
34	! j_north
35	!
36	! 1036 2012-10-22 13:43:42Z raasch
37	! code put under GPL (PALM 3.9)
38	!
39	! 1015 2012-09-27 09:23:24Z raasch
40	! accelerator version (*_acc) added
41	!
42	! 1001 2012-09-13 14:08:46Z raasch
43	! arrays comunicated by module instead of parameter list
44	!
45	! 978 2012-08-09 08:28:32Z fricke
46	! outflow damping layer removed
47	! kmxm_x/_y and kmxp_x/_y change to kmxm and kmxp
48	!
49	! 667 2010-12-23 12:06:00Z suehring/gryschka
50	! nxl-1, nxr+1, nys-1, nyn+1 replaced by nxlg, nxrg, nysg, nyng
51	!
52	! 366 2009-08-25 08:06:27Z raasch
53	! bc_lr replaced by bc_lr_cyc
54	!
55	! 106 2007-08-16 14:30:26Z raasch
56	! Momentumflux at top (vswst) included as boundary condition,
57	! j loop is starting from nysv (needed for non-cyclic boundary conditions)
58	!
59	! 75 2007-03-22 09:54:05Z raasch
60	! Wall functions now include diabatic conditions, call of routine wall_fluxes,
61	! z0 removed from argument list, vynp eliminated
62	!
63	! 20 2007-02-26 00:12:32Z raasch
64	! Bugfix: ddzw dimensioned 1:nzt"+1"
65	!
66	! RCS Log replace by Id keyword, revision history cleaned up
67	!
68	! Revision 1.15 2006/02/23 10:36:00 raasch
69	! nzb_2d replaced by nzb_v_outer in horizontal diffusion and by nzb_v_inner
70	! or nzb_diff_v, respectively, in vertical diffusion,
71	! wall functions added for north and south walls, +z0 in argument list,
72	! terms containing w(k-1,..) are removed from the Prandtl-layer equation
73	! because they cause errors at the edges of topography
74	! WARNING: loops containing the MAX function are still not properly vectorized!
75	!
76	! Revision 1.1 1997/09/12 06:24:01 raasch
77	! Initial revision
78	!
79	!
80	! Description:
81	! ------------
82	! Diffusion term of the v-component
83	!------------------------------------------------------------------------------!
84
85	USE wall_fluxes_mod
86
87	PRIVATE
88	PUBLIC diffusion_v, diffusion_v_acc
89
90	INTERFACE diffusion_v
91	MODULE PROCEDURE diffusion_v
92	MODULE PROCEDURE diffusion_v_ij
93	END INTERFACE diffusion_v
94
95	INTERFACE diffusion_v_acc
96	MODULE PROCEDURE diffusion_v_acc
97	END INTERFACE diffusion_v_acc
98
99	CONTAINS
100
101
102	!------------------------------------------------------------------------------!
103	! Call for all grid points
104	!------------------------------------------------------------------------------!
105	SUBROUTINE diffusion_v
106
107	USE arrays_3d
108	USE control_parameters
109	USE grid_variables
110	USE indices
111
112	IMPLICIT NONE
113
114	INTEGER :: i, j, k
115	REAL :: kmxm, kmxp, kmzm, kmzp
116
117	REAL, DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) :: vsus
118
119	!
120	!-- First calculate horizontal momentum flux v'u' at vertical walls,
121	!-- if neccessary
122	IF ( topography /= 'flat' ) THEN
123	CALL wall_fluxes( vsus, 0.0, 1.0, 0.0, 0.0, nzb_v_inner, &
124	nzb_v_outer, wall_v )
125	ENDIF
126
127	DO i = nxl, nxr
128	DO j = nysv, nyn
129	!
130	!-- Compute horizontal diffusion
131	DO k = nzb_v_outer(j,i)+1, nzt
132	!
133	!-- Interpolate eddy diffusivities on staggered gridpoints
134	kmxp = 0.25 * &
135	( km(k,j,i)+km(k,j,i+1)+km(k,j-1,i)+km(k,j-1,i+1) )
136	kmxm = 0.25 * &
137	( km(k,j,i)+km(k,j,i-1)+km(k,j-1,i)+km(k,j-1,i-1) )
138
139	tend(k,j,i) = tend(k,j,i) &
140	& + ( kmxp * ( v(k,j,i+1) - v(k,j,i) ) * ddx &
141	& + kmxp * ( u(k,j,i+1) - u(k,j-1,i+1) ) * ddy &
142	& - kmxm * ( v(k,j,i) - v(k,j,i-1) ) * ddx &
143	& - kmxm * ( u(k,j,i) - u(k,j-1,i) ) * ddy &
144	& ) * ddx &
145	& + 2.0 * ( &
146	& km(k,j,i) * ( v(k,j+1,i) - v(k,j,i) ) &
147	& - km(k,j-1,i) * ( v(k,j,i) - v(k,j-1,i) ) &
148	& ) * ddy2
149	ENDDO
150
151	!
152	!-- Wall functions at the left and right walls, respectively
153	IF ( wall_v(j,i) /= 0.0 ) THEN
154
155	DO k = nzb_v_inner(j,i)+1, nzb_v_outer(j,i)
156	kmxp = 0.25 * &
157	( km(k,j,i)+km(k,j,i+1)+km(k,j-1,i)+km(k,j-1,i+1) )
158	kmxm = 0.25 * &
159	( km(k,j,i)+km(k,j,i-1)+km(k,j-1,i)+km(k,j-1,i-1) )
160
161	tend(k,j,i) = tend(k,j,i) &
162	+ 2.0 * ( &
163	km(k,j,i) * ( v(k,j+1,i) - v(k,j,i) ) &
164	- km(k,j-1,i) * ( v(k,j,i) - v(k,j-1,i) ) &
165	) * ddy2 &
166	+ ( fxp(j,i) * ( &
167	kmxp * ( v(k,j,i+1) - v(k,j,i) ) * ddx &
168	+ kmxp * ( u(k,j,i+1) - u(k,j-1,i+1) ) * ddy &
169	) &
170	- fxm(j,i) * ( &
171	kmxm * ( v(k,j,i) - v(k,j,i-1) ) * ddx &
172	+ kmxm * ( u(k,j,i) - u(k,j-1,i) ) * ddy &
173	) &
174	+ wall_v(j,i) * vsus(k,j,i) &
175	) * ddx
176	ENDDO
177	ENDIF
178
179	!
180	!-- Compute vertical diffusion. In case of simulating a Prandtl
181	!-- layer, index k starts at nzb_v_inner+2.
182	DO k = nzb_diff_v(j,i), nzt_diff
183	!
184	!-- Interpolate eddy diffusivities on staggered gridpoints
185	kmzp = 0.25 * &
186	( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
187	kmzm = 0.25 * &
188	( km(k,j,i)+km(k-1,j,i)+km(k,j-1,i)+km(k-1,j-1,i) )
189
190	tend(k,j,i) = tend(k,j,i) &
191	& + ( kmzp * ( ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1) &
192	& + ( w(k,j,i) - w(k,j-1,i) ) * ddy &
193	& ) &
194	& - kmzm * ( ( v(k,j,i) - v(k-1,j,i) ) * ddzu(k) &
195	& + ( w(k-1,j,i) - w(k-1,j-1,i) ) * ddy &
196	& ) &
197	& ) * ddzw(k)
198	ENDDO
199
200	!
201	!-- Vertical diffusion at the first grid point above the surface,
202	!-- if the momentum flux at the bottom is given by the Prandtl law
203	!-- or if it is prescribed by the user.
204	!-- Difference quotient of the momentum flux is not formed over
205	!-- half of the grid spacing (2.0*ddzw(k)) any more, since the
206	!-- comparison with other (LES) modell showed that the values of
207	!-- the momentum flux becomes too large in this case.
208	!-- The term containing w(k-1,..) (see above equation) is removed here
209	!-- because the vertical velocity is assumed to be zero at the surface.
210	IF ( use_surface_fluxes ) THEN
211	k = nzb_v_inner(j,i)+1
212	!
213	!-- Interpolate eddy diffusivities on staggered gridpoints
214	kmzp = 0.25 * &
215	( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
216	kmzm = 0.25 * &
217	( km(k,j,i)+km(k-1,j,i)+km(k,j-1,i)+km(k-1,j-1,i) )
218
219	tend(k,j,i) = tend(k,j,i) &
220	& + ( kmzp * ( w(k,j,i) - w(k,j-1,i) ) * ddy &
221	& ) * ddzw(k) &
222	& + ( kmzp * ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1) &
223	& + vsws(j,i) &
224	& ) * ddzw(k)
225	ENDIF
226
227	!
228	!-- Vertical diffusion at the first gridpoint below the top boundary,
229	!-- if the momentum flux at the top is prescribed by the user
230	IF ( use_top_fluxes .AND. constant_top_momentumflux ) THEN
231	k = nzt
232	!
233	!-- Interpolate eddy diffusivities on staggered gridpoints
234	kmzp = 0.25 * &
235	( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
236	kmzm = 0.25 * &
237	( km(k,j,i)+km(k-1,j,i)+km(k,j-1,i)+km(k-1,j-1,i) )
238
239	tend(k,j,i) = tend(k,j,i) &
240	& - ( kmzm * ( w(k-1,j,i) - w(k-1,j-1,i) ) * ddy &
241	& ) * ddzw(k) &
242	& + ( -vswst(j,i) &
243	& - kmzm * ( v(k,j,i) - v(k-1,j,i) ) * ddzu(k) &
244	& ) * ddzw(k)
245	ENDIF
246
247	ENDDO
248	ENDDO
249
250	END SUBROUTINE diffusion_v
251
252
253	!------------------------------------------------------------------------------!
254	! Call for all grid points - accelerator version
255	!------------------------------------------------------------------------------!
256	SUBROUTINE diffusion_v_acc
257
258	USE arrays_3d
259	USE control_parameters
260	USE grid_variables
261	USE indices
262
263	IMPLICIT NONE
264
265	INTEGER :: i, j, k
266	REAL :: kmxm, kmxp, kmzm, kmzp
267
268	REAL, DIMENSION(nzb:nzt+1,nys:nyn,nxl:nxr) :: vsus
269	!$acc declare create ( vsus )
270
271	!
272	!-- First calculate horizontal momentum flux v'u' at vertical walls,
273	!-- if neccessary
274	IF ( topography /= 'flat' ) THEN
275	CALL wall_fluxes_acc( vsus, 0.0, 1.0, 0.0, 0.0, nzb_v_inner, &
276	nzb_v_outer, wall_v )
277	ENDIF
278
279	!$acc kernels present ( u, v, w, km, tend, vsws, vswst ) &
280	!$acc present ( ddzu, ddzw, fxm, fxp, wall_v ) &
281	!$acc present ( nzb_v_inner, nzb_v_outer, nzb_diff_v )
282	DO i = i_left, i_right
283	DO j = j_south, j_north
284	!
285	!-- Compute horizontal diffusion
286	DO k = 1, nzt
287	IF ( k > nzb_v_outer(j,i) ) THEN
288	!
289	!-- Interpolate eddy diffusivities on staggered gridpoints
290	kmxp = 0.25 * &
291	( km(k,j,i)+km(k,j,i+1)+km(k,j-1,i)+km(k,j-1,i+1) )
292	kmxm = 0.25 * &
293	( km(k,j,i)+km(k,j,i-1)+km(k,j-1,i)+km(k,j-1,i-1) )
294
295	tend(k,j,i) = tend(k,j,i) &
296	& + ( kmxp * ( v(k,j,i+1) - v(k,j,i) ) * ddx &
297	& + kmxp * ( u(k,j,i+1) - u(k,j-1,i+1) ) * ddy &
298	& - kmxm * ( v(k,j,i) - v(k,j,i-1) ) * ddx &
299	& - kmxm * ( u(k,j,i) - u(k,j-1,i) ) * ddy &
300	& ) * ddx &
301	& + 2.0 * ( &
302	& km(k,j,i) * ( v(k,j+1,i) - v(k,j,i) ) &
303	& - km(k,j-1,i) * ( v(k,j,i) - v(k,j-1,i) ) &
304	& ) * ddy2
305	ENDIF
306	ENDDO
307
308	!
309	!-- Wall functions at the left and right walls, respectively
310	DO k = 1, nzt
311	IF( k > nzb_v_inner(j,i) .AND. k <= nzb_v_outer(j,i) .AND. &
312	wall_v(j,i) /= 0.0 ) THEN
313
314	kmxp = 0.25 * &
315	( km(k,j,i)+km(k,j,i+1)+km(k,j-1,i)+km(k,j-1,i+1) )
316	kmxm = 0.25 * &
317	( km(k,j,i)+km(k,j,i-1)+km(k,j-1,i)+km(k,j-1,i-1) )
318
319	tend(k,j,i) = tend(k,j,i) &
320	+ 2.0 * ( &
321	km(k,j,i) * ( v(k,j+1,i) - v(k,j,i) ) &
322	- km(k,j-1,i) * ( v(k,j,i) - v(k,j-1,i) ) &
323	) * ddy2 &
324	+ ( fxp(j,i) * ( &
325	kmxp * ( v(k,j,i+1) - v(k,j,i) ) * ddx &
326	+ kmxp * ( u(k,j,i+1) - u(k,j-1,i+1) ) * ddy &
327	) &
328	- fxm(j,i) * ( &
329	kmxm * ( v(k,j,i) - v(k,j,i-1) ) * ddx &
330	+ kmxm * ( u(k,j,i) - u(k,j-1,i) ) * ddy &
331	) &
332	+ wall_v(j,i) * vsus(k,j,i) &
333	) * ddx
334	ENDIF
335	ENDDO
336
337	!
338	!-- Compute vertical diffusion. In case of simulating a Prandtl
339	!-- layer, index k starts at nzb_v_inner+2.
340	DO k = 1, nzt_diff
341	IF ( k >= nzb_diff_v(j,i) ) THEN
342	!
343	!-- Interpolate eddy diffusivities on staggered gridpoints
344	kmzp = 0.25 * &
345	( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
346	kmzm = 0.25 * &
347	( km(k,j,i)+km(k-1,j,i)+km(k,j-1,i)+km(k-1,j-1,i) )
348
349	tend(k,j,i) = tend(k,j,i) &
350	& + ( kmzp * ( ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1)&
351	& + ( w(k,j,i) - w(k,j-1,i) ) * ddy &
352	& ) &
353	& - kmzm * ( ( v(k,j,i) - v(k-1,j,i) ) * ddzu(k)&
354	& + ( w(k-1,j,i) - w(k-1,j-1,i) ) * ddy &
355	& ) &
356	& ) * ddzw(k)
357	ENDIF
358	ENDDO
359
360	ENDDO
361	ENDDO
362
363	!
364	!-- Vertical diffusion at the first grid point above the surface,
365	!-- if the momentum flux at the bottom is given by the Prandtl law
366	!-- or if it is prescribed by the user.
367	!-- Difference quotient of the momentum flux is not formed over
368	!-- half of the grid spacing (2.0*ddzw(k)) any more, since the
369	!-- comparison with other (LES) modell showed that the values of
370	!-- the momentum flux becomes too large in this case.
371	!-- The term containing w(k-1,..) (see above equation) is removed here
372	!-- because the vertical velocity is assumed to be zero at the surface.
373	IF ( use_surface_fluxes ) THEN
374
375	DO i = i_left, i_right
376	DO j = j_south, j_north
377
378	k = nzb_v_inner(j,i)+1
379	!
380	!-- Interpolate eddy diffusivities on staggered gridpoints
381	kmzp = 0.25 * &
382	( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
383	kmzm = 0.25 * &
384	( km(k,j,i)+km(k-1,j,i)+km(k,j-1,i)+km(k-1,j-1,i) )
385
386	tend(k,j,i) = tend(k,j,i) &
387	& + ( kmzp * ( w(k,j,i) - w(k,j-1,i) ) * ddy &
388	& ) * ddzw(k) &
389	& + ( kmzp * ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1) &
390	& + vsws(j,i) &
391	& ) * ddzw(k)
392	ENDDO
393	ENDDO
394
395	ENDIF
396
397	!
398	!-- Vertical diffusion at the first gridpoint below the top boundary,
399	!-- if the momentum flux at the top is prescribed by the user
400	IF ( use_top_fluxes .AND. constant_top_momentumflux ) THEN
401
402	k = nzt
403
404	DO i = i_left, i_right
405	DO j = j_south, j_north
406
407	!
408	!-- Interpolate eddy diffusivities on staggered gridpoints
409	kmzp = 0.25 * &
410	( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
411	kmzm = 0.25 * &
412	( km(k,j,i)+km(k-1,j,i)+km(k,j-1,i)+km(k-1,j-1,i) )
413
414	tend(k,j,i) = tend(k,j,i) &
415	& - ( kmzm * ( w(k-1,j,i) - w(k-1,j-1,i) ) * ddy &
416	& ) * ddzw(k) &
417	& + ( -vswst(j,i) &
418	& - kmzm * ( v(k,j,i) - v(k-1,j,i) ) * ddzu(k) &
419	& ) * ddzw(k)
420	ENDDO
421	ENDDO
422
423	ENDIF
424	!$acc end kernels
425
426	END SUBROUTINE diffusion_v_acc
427
428
429	!------------------------------------------------------------------------------!
430	! Call for grid point i,j
431	!------------------------------------------------------------------------------!
432	SUBROUTINE diffusion_v_ij( i, j )
433
434	USE arrays_3d
435	USE control_parameters
436	USE grid_variables
437	USE indices
438
439	IMPLICIT NONE
440
441	INTEGER :: i, j, k
442	REAL :: kmxm, kmxp, kmzm, kmzp
443
444	REAL, DIMENSION(nzb:nzt+1) :: vsus
445
446	!
447	!-- Compute horizontal diffusion
448	DO k = nzb_v_outer(j,i)+1, nzt
449	!
450	!-- Interpolate eddy diffusivities on staggered gridpoints
451	kmxp = 0.25 * ( km(k,j,i)+km(k,j,i+1)+km(k,j-1,i)+km(k,j-1,i+1) )
452	kmxm = 0.25 * ( km(k,j,i)+km(k,j,i-1)+km(k,j-1,i)+km(k,j-1,i-1) )
453
454	tend(k,j,i) = tend(k,j,i) &
455	& + ( kmxp * ( v(k,j,i+1) - v(k,j,i) ) * ddx &
456	& + kmxp * ( u(k,j,i+1) - u(k,j-1,i+1) ) * ddy &
457	& - kmxm * ( v(k,j,i) - v(k,j,i-1) ) * ddx &
458	& - kmxm * ( u(k,j,i) - u(k,j-1,i) ) * ddy &
459	& ) * ddx &
460	& + 2.0 * ( &
461	& km(k,j,i) * ( v(k,j+1,i) - v(k,j,i) ) &
462	& - km(k,j-1,i) * ( v(k,j,i) - v(k,j-1,i) ) &
463	& ) * ddy2
464	ENDDO
465
466	!
467	!-- Wall functions at the left and right walls, respectively
468	IF ( wall_v(j,i) /= 0.0 ) THEN
469
470	!
471	!-- Calculate the horizontal momentum flux v'u'
472	CALL wall_fluxes( i, j, nzb_v_inner(j,i)+1, nzb_v_outer(j,i), &
473	vsus, 0.0, 1.0, 0.0, 0.0 )
474
475	DO k = nzb_v_inner(j,i)+1, nzb_v_outer(j,i)
476	kmxp = 0.25 * &
477	( km(k,j,i)+km(k,j,i+1)+km(k,j-1,i)+km(k,j-1,i+1) )
478	kmxm = 0.25 * &
479	( km(k,j,i)+km(k,j,i-1)+km(k,j-1,i)+km(k,j-1,i-1) )
480
481	tend(k,j,i) = tend(k,j,i) &
482	+ 2.0 * ( &
483	km(k,j,i) * ( v(k,j+1,i) - v(k,j,i) ) &
484	- km(k,j-1,i) * ( v(k,j,i) - v(k,j-1,i) ) &
485	) * ddy2 &
486	+ ( fxp(j,i) * ( &
487	kmxp * ( v(k,j,i+1) - v(k,j,i) ) * ddx &
488	+ kmxp * ( u(k,j,i+1) - u(k,j-1,i+1) ) * ddy &
489	) &
490	- fxm(j,i) * ( &
491	kmxm * ( v(k,j,i) - v(k,j,i-1) ) * ddx &
492	+ kmxm * ( u(k,j,i) - u(k,j-1,i) ) * ddy &
493	) &
494	+ wall_v(j,i) * vsus(k) &
495	) * ddx
496	ENDDO
497	ENDIF
498
499	!
500	!-- Compute vertical diffusion. In case of simulating a Prandtl layer,
501	!-- index k starts at nzb_v_inner+2.
502	DO k = nzb_diff_v(j,i), nzt_diff
503	!
504	!-- Interpolate eddy diffusivities on staggered gridpoints
505	kmzp = 0.25 * ( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
506	kmzm = 0.25 * ( km(k,j,i)+km(k-1,j,i)+km(k,j-1,i)+km(k-1,j-1,i) )
507
508	tend(k,j,i) = tend(k,j,i) &
509	& + ( kmzp * ( ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1) &
510	& + ( w(k,j,i) - w(k,j-1,i) ) * ddy &
511	& ) &
512	& - kmzm * ( ( v(k,j,i) - v(k-1,j,i) ) * ddzu(k) &
513	& + ( w(k-1,j,i) - w(k-1,j-1,i) ) * ddy &
514	& ) &
515	& ) * ddzw(k)
516	ENDDO
517
518	!
519	!-- Vertical diffusion at the first grid point above the surface, if the
520	!-- momentum flux at the bottom is given by the Prandtl law or if it is
521	!-- prescribed by the user.
522	!-- Difference quotient of the momentum flux is not formed over half of
523	!-- the grid spacing (2.0*ddzw(k)) any more, since the comparison with
524	!-- other (LES) modell showed that the values of the momentum flux becomes
525	!-- too large in this case.
526	!-- The term containing w(k-1,..) (see above equation) is removed here
527	!-- because the vertical velocity is assumed to be zero at the surface.
528	IF ( use_surface_fluxes ) THEN
529	k = nzb_v_inner(j,i)+1
530	!
531	!-- Interpolate eddy diffusivities on staggered gridpoints
532	kmzp = 0.25 * ( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
533	kmzm = 0.25 * ( km(k,j,i)+km(k-1,j,i)+km(k,j-1,i)+km(k-1,j-1,i) )
534
535	tend(k,j,i) = tend(k,j,i) &
536	& + ( kmzp * ( w(k,j,i) - w(k,j-1,i) ) * ddy &
537	& ) * ddzw(k) &
538	& + ( kmzp * ( v(k+1,j,i) - v(k,j,i) ) * ddzu(k+1) &
539	& + vsws(j,i) &
540	& ) * ddzw(k)
541	ENDIF
542
543	!
544	!-- Vertical diffusion at the first gridpoint below the top boundary,
545	!-- if the momentum flux at the top is prescribed by the user
546	IF ( use_top_fluxes .AND. constant_top_momentumflux ) THEN
547	k = nzt
548	!
549	!-- Interpolate eddy diffusivities on staggered gridpoints
550	kmzp = 0.25 * &
551	( km(k,j,i)+km(k+1,j,i)+km(k,j-1,i)+km(k+1,j-1,i) )
552	kmzm = 0.25 * &
553	( km(k,j,i)+km(k-1,j,i)+km(k,j-1,i)+km(k-1,j-1,i) )
554
555	tend(k,j,i) = tend(k,j,i) &
556	& - ( kmzm * ( w(k-1,j,i) - w(k-1,j-1,i) ) * ddy &
557	& ) * ddzw(k) &
558	& + ( -vswst(j,i) &
559	& - kmzm * ( v(k,j,i) - v(k-1,j,i) ) * ddzu(k) &
560	& ) * ddzw(k)
561	ENDIF
562
563	END SUBROUTINE diffusion_v_ij
564
565	END MODULE diffusion_v_mod

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |