Home

Context Navigation

source: palm/trunk/SOURCE/diffusion_e.f90 @ 1517

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

Line
1	MODULE diffusion_e_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-2014 Leibniz Universitaet Hannover
18	!--------------------------------------------------------------------------------!
19	!
20	! Current revisions:
21	! -----------------
22	!
23	!
24	! Former revisions:
25	! -----------------
26	! $Id: diffusion_e.f90 1375 2014-04-25 13:07:08Z hoffmann $
27	!
28	! 1374 2014-04-25 12:55:07Z raasch
29	! missing variables added to ONLY list
30	! rif removed from acc-present-list
31	!
32	! 1340 2014-03-25 19:45:13Z kanani
33	! REAL constants defined as wp-kind
34	!
35	! 1320 2014-03-20 08:40:49Z raasch
36	! ONLY-attribute added to USE-statements,
37	! kind-parameters added to all INTEGER and REAL declaration statements,
38	! kinds are defined in new module kinds,
39	! revision history before 2012 removed,
40	! comment fields (!:) to be used for variable explanations added to
41	! all variable declaration statements
42	!
43	! 1257 2013-11-08 15:18:40Z raasch
44	! openacc loop and loop vector clauses removed
45	!
46	! 1179 2013-06-14 05:57:58Z raasch
47	! use_reference renamed use_single_reference_value
48	!
49	! 1171 2013-05-30 11:27:45Z raasch
50	! use_reference-case activated in accelerator version
51	!
52	! 1128 2013-04-12 06:19:32Z raasch
53	! loop index bounds in accelerator version replaced by i_left, i_right, j_south,
54	! j_north
55	!
56	! 1065 2012-11-22 17:42:36Z hoffmann
57	! Enabled the claculation of diss in case of turbulence = .TRUE. (parameterized
58	! effects of turbulence on autoconversion and accretion in two-moments cloud
59	! physics scheme).
60	!
61	! 1036 2012-10-22 13:43:42Z raasch
62	! code put under GPL (PALM 3.9)
63	!
64	! 1015 2012-09-27 09:23:24Z raasch
65	! accelerator version (*_acc) added,
66	! adjustment of mixing length to the Prandtl mixing length at first grid point
67	! above ground removed
68	!
69	! 1010 2012-09-20 07:59:54Z raasch
70	! cpp switch __nopointer added for pointer free version
71	!
72	! 1001 2012-09-13 14:08:46Z raasch
73	! most arrays comunicated by module instead of parameter list
74	!
75	! 825 2012-02-19 03:03:44Z raasch
76	! wang_collision_kernel renamed wang_kernel
77	!
78	! Revision 1.1 1997/09/19 07:40:24 raasch
79	! Initial revision
80	!
81	!
82	! Description:
83	! ------------
84	! Diffusion- and dissipation terms for the TKE
85	!------------------------------------------------------------------------------!
86
87	PRIVATE
88	PUBLIC diffusion_e, diffusion_e_acc
89
90
91	INTERFACE diffusion_e
92	MODULE PROCEDURE diffusion_e
93	MODULE PROCEDURE diffusion_e_ij
94	END INTERFACE diffusion_e
95
96	INTERFACE diffusion_e_acc
97	MODULE PROCEDURE diffusion_e_acc
98	END INTERFACE diffusion_e_acc
99
100	CONTAINS
101
102
103	!------------------------------------------------------------------------------!
104	! Call for all grid points
105	!------------------------------------------------------------------------------!
106	SUBROUTINE diffusion_e( var, var_reference )
107
108	USE arrays_3d, &
109	ONLY: dd2zu, ddzu, ddzw, diss, e, km, l_grid, tend, zu, zw
110
111	USE control_parameters, &
112	ONLY: atmos_ocean_sign, g, turbulence, use_single_reference_value, &
113	wall_adjustment, wall_adjustment_factor
114
115	USE grid_variables, &
116	ONLY: ddx2, ddy2
117
118	USE indices, &
119	ONLY: nxl, nxlg, nxr, nxrg, nyn, nyng, nys, nysg, nzb, nzb_s_inner,&
120	nzt
121
122	USE kinds
123
124	USE particle_attributes, &
125	ONLY: use_sgs_for_particles, wang_kernel
126
127	IMPLICIT NONE
128
129	INTEGER(iwp) :: i !:
130	INTEGER(iwp) :: j !:
131	INTEGER(iwp) :: k !:
132	REAL(wp) :: dvar_dz !:
133	REAL(wp) :: l_stable !:
134	REAL(wp) :: var_reference !:
135
136	#if defined( __nopointer )
137	REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) :: var !:
138	#else
139	REAL(wp), DIMENSION(:,:,:), POINTER :: var !:
140	#endif
141	REAL(wp), DIMENSION(nzb+1:nzt,nys:nyn) :: dissipation !:
142	REAL(wp), DIMENSION(nzb+1:nzt,nys:nyn) :: l !:
143	REAL(wp), DIMENSION(nzb+1:nzt,nys:nyn) :: ll !:
144
145
146	!
147	!-- This if clause must be outside the k-loop because otherwise
148	!-- runtime errors occur with -C hopt on NEC
149	IF ( use_single_reference_value ) THEN
150
151	DO i = nxl, nxr
152	DO j = nys, nyn
153	DO k = nzb_s_inner(j,i)+1, nzt
154	!
155	!-- Calculate the mixing length (for dissipation)
156	dvar_dz = atmos_ocean_sign * &
157	( var(k+1,j,i) - var(k-1,j,i) ) * dd2zu(k)
158	IF ( dvar_dz > 0.0_wp ) THEN
159	l_stable = 0.76_wp * SQRT( e(k,j,i) ) / &
160	SQRT( g / var_reference * dvar_dz ) + 1E-5_wp
161	ELSE
162	l_stable = l_grid(k)
163	ENDIF
164	!
165	!-- Adjustment of the mixing length
166	IF ( wall_adjustment ) THEN
167	l(k,j) = MIN( wall_adjustment_factor * &
168	( zu(k) - zw(nzb_s_inner(j,i)) ), &
169	l_grid(k), l_stable )
170	ll(k,j) = MIN( wall_adjustment_factor * &
171	( zu(k) - zw(nzb_s_inner(j,i)) ), &
172	l_grid(k) )
173	ELSE
174	l(k,j) = MIN( l_grid(k), l_stable )
175	ll(k,j) = l_grid(k)
176	ENDIF
177
178	ENDDO
179	ENDDO
180
181	!
182	!-- Calculate the tendency terms
183	DO j = nys, nyn
184	DO k = nzb_s_inner(j,i)+1, nzt
185
186	dissipation(k,j) = ( 0.19_wp + 0.74_wp * l(k,j) / ll(k,j) ) * &
187	e(k,j,i) * SQRT( e(k,j,i) ) / l(k,j)
188
189	tend(k,j,i) = tend(k,j,i) &
190	+ ( &
191	( km(k,j,i)+km(k,j,i+1) ) * ( e(k,j,i+1)-e(k,j,i) ) &
192	- ( km(k,j,i)+km(k,j,i-1) ) * ( e(k,j,i)-e(k,j,i-1) ) &
193	) * ddx2 &
194	+ ( &
195	( km(k,j,i)+km(k,j+1,i) ) * ( e(k,j+1,i)-e(k,j,i) ) &
196	- ( km(k,j,i)+km(k,j-1,i) ) * ( e(k,j,i)-e(k,j-1,i) ) &
197	) * ddy2 &
198	+ ( &
199	( km(k,j,i)+km(k+1,j,i) ) * ( e(k+1,j,i)-e(k,j,i) ) * ddzu(k+1) &
200	- ( km(k,j,i)+km(k-1,j,i) ) * ( e(k,j,i)-e(k-1,j,i) ) * ddzu(k) &
201	) * ddzw(k) &
202	- dissipation(k,j)
203
204	ENDDO
205	ENDDO
206
207	!
208	!-- Store dissipation if needed for calculating the sgs particle
209	!-- velocities
210	IF ( use_sgs_for_particles .OR. wang_kernel .OR. &
211	turbulence ) THEN
212	DO j = nys, nyn
213	DO k = nzb_s_inner(j,i)+1, nzt
214	diss(k,j,i) = dissipation(k,j)
215	ENDDO
216	ENDDO
217	ENDIF
218
219	ENDDO
220
221	ELSE
222
223	DO i = nxl, nxr
224	DO j = nys, nyn
225	DO k = nzb_s_inner(j,i)+1, nzt
226	!
227	!-- Calculate the mixing length (for dissipation)
228	dvar_dz = atmos_ocean_sign * &
229	( var(k+1,j,i) - var(k-1,j,i) ) * dd2zu(k)
230	IF ( dvar_dz > 0.0_wp ) THEN
231	l_stable = 0.76_wp * SQRT( e(k,j,i) ) / &
232	SQRT( g / var(k,j,i) * dvar_dz ) + 1E-5_wp
233	ELSE
234	l_stable = l_grid(k)
235	ENDIF
236	!
237	!-- Adjustment of the mixing length
238	IF ( wall_adjustment ) THEN
239	l(k,j) = MIN( wall_adjustment_factor * &
240	( zu(k) - zw(nzb_s_inner(j,i)) ), &
241	l_grid(k), l_stable )
242	ll(k,j) = MIN( wall_adjustment_factor * &
243	( zu(k) - zw(nzb_s_inner(j,i)) ), &
244	l_grid(k) )
245	ELSE
246	l(k,j) = MIN( l_grid(k), l_stable )
247	ll(k,j) = l_grid(k)
248	ENDIF
249
250	ENDDO
251	ENDDO
252
253	!
254	!-- Calculate the tendency terms
255	DO j = nys, nyn
256	DO k = nzb_s_inner(j,i)+1, nzt
257
258	dissipation(k,j) = ( 0.19_wp + 0.74_wp * l(k,j) / ll(k,j) ) * &
259	e(k,j,i) * SQRT( e(k,j,i) ) / l(k,j)
260
261	tend(k,j,i) = tend(k,j,i) &
262	+ ( &
263	( km(k,j,i)+km(k,j,i+1) ) * ( e(k,j,i+1)-e(k,j,i) ) &
264	- ( km(k,j,i)+km(k,j,i-1) ) * ( e(k,j,i)-e(k,j,i-1) ) &
265	) * ddx2 &
266	+ ( &
267	( km(k,j,i)+km(k,j+1,i) ) * ( e(k,j+1,i)-e(k,j,i) ) &
268	- ( km(k,j,i)+km(k,j-1,i) ) * ( e(k,j,i)-e(k,j-1,i) ) &
269	) * ddy2 &
270	+ ( &
271	( km(k,j,i)+km(k+1,j,i) ) * ( e(k+1,j,i)-e(k,j,i) ) * ddzu(k+1) &
272	- ( km(k,j,i)+km(k-1,j,i) ) * ( e(k,j,i)-e(k-1,j,i) ) * ddzu(k) &
273	) * ddzw(k) &
274	- dissipation(k,j)
275
276	ENDDO
277	ENDDO
278
279	!
280	!-- Store dissipation if needed for calculating the sgs particle
281	!-- velocities
282	IF ( use_sgs_for_particles .OR. wang_kernel .OR. &
283	turbulence ) THEN
284	DO j = nys, nyn
285	DO k = nzb_s_inner(j,i)+1, nzt
286	diss(k,j,i) = dissipation(k,j)
287	ENDDO
288	ENDDO
289	ENDIF
290
291	ENDDO
292
293	ENDIF
294
295	!
296	!-- Boundary condition for dissipation
297	IF ( use_sgs_for_particles .OR. wang_kernel .OR. turbulence ) THEN
298	DO i = nxl, nxr
299	DO j = nys, nyn
300	diss(nzb_s_inner(j,i),j,i) = diss(nzb_s_inner(j,i)+1,j,i)
301	ENDDO
302	ENDDO
303	ENDIF
304
305	END SUBROUTINE diffusion_e
306
307
308	!------------------------------------------------------------------------------!
309	! Call for all grid points - accelerator version
310	!------------------------------------------------------------------------------!
311	SUBROUTINE diffusion_e_acc( var, var_reference )
312
313	USE arrays_3d, &
314	ONLY: dd2zu, ddzu, ddzw, diss, e, km, l_grid, tend, zu, zw
315
316	USE control_parameters, &
317	ONLY: atmos_ocean_sign, g, turbulence, use_single_reference_value, &
318	wall_adjustment, wall_adjustment_factor
319
320	USE grid_variables, &
321	ONLY: ddx2, ddy2
322
323	USE indices, &
324	ONLY: i_left, i_right, j_north, j_south, nxlg, nxrg, nyng, nysg, &
325	nzb, nzb_s_inner, nzt
326
327	USE kinds
328
329	USE particle_attributes, &
330	ONLY: use_sgs_for_particles, wang_kernel
331
332	IMPLICIT NONE
333
334	INTEGER(iwp) :: i !:
335	INTEGER(iwp) :: j !:
336	INTEGER(iwp) :: k !:
337	REAL(wp) :: dissipation !:
338	REAL(wp) :: dvar_dz !:
339	REAL(wp) :: l !:
340	REAL(wp) :: ll !:
341	REAL(wp) :: l_stable !:
342	REAL(wp) :: var_reference !:
343
344	#if defined( __nopointer )
345	REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) :: var !:
346	#else
347	REAL(wp), DIMENSION(:,:,:), POINTER :: var !:
348	#endif
349
350
351	!
352	!-- This if clause must be outside the k-loop because otherwise
353	!-- runtime errors occur with -C hopt on NEC
354	IF ( use_single_reference_value ) THEN
355
356	!$acc kernels present( ddzu, ddzw, dd2zu, diss, e, km, l_grid ) &
357	!$acc present( nzb_s_inner, tend, var, zu, zw )
358	DO i = i_left, i_right
359	DO j = j_south, j_north
360	DO k = 1, nzt
361
362	IF ( k > nzb_s_inner(j,i) ) THEN
363	!
364	!-- Calculate the mixing length (for dissipation)
365	dvar_dz = atmos_ocean_sign * &
366	( var(k+1,j,i) - var(k-1,j,i) ) * dd2zu(k)
367	IF ( dvar_dz > 0.0_wp ) THEN
368	l_stable = 0.76_wp * SQRT( e(k,j,i) ) / &
369	SQRT( g / var_reference * dvar_dz ) + 1E-5_wp
370	ELSE
371	l_stable = l_grid(k)
372	ENDIF
373	!
374	!-- Adjustment of the mixing length
375	IF ( wall_adjustment ) THEN
376	l = MIN( wall_adjustment_factor * &
377	( zu(k) - zw(nzb_s_inner(j,i)) ), &
378	l_grid(k), l_stable )
379	ll = MIN( wall_adjustment_factor * &
380	( zu(k) - zw(nzb_s_inner(j,i)) ), &
381	l_grid(k) )
382	ELSE
383	l = MIN( l_grid(k), l_stable )
384	ll = l_grid(k)
385	ENDIF
386	!
387	!-- Calculate the tendency terms
388	dissipation = ( 0.19_wp + 0.74_wp * l / ll ) * &
389	e(k,j,i) * SQRT( e(k,j,i) ) / l
390
391	tend(k,j,i) = tend(k,j,i) &
392	+ ( &
393	( km(k,j,i)+km(k,j,i+1) ) * ( e(k,j,i+1)-e(k,j,i) ) &
394	- ( km(k,j,i)+km(k,j,i-1) ) * ( e(k,j,i)-e(k,j,i-1) ) &
395	) * ddx2 &
396	+ ( &
397	( km(k,j,i)+km(k,j+1,i) ) * ( e(k,j+1,i)-e(k,j,i) ) &
398	- ( km(k,j,i)+km(k,j-1,i) ) * ( e(k,j,i)-e(k,j-1,i) ) &
399	) * ddy2 &
400	+ ( &
401	( km(k,j,i)+km(k+1,j,i) ) * ( e(k+1,j,i)-e(k,j,i) ) * ddzu(k+1) &
402	- ( km(k,j,i)+km(k-1,j,i) ) * ( e(k,j,i)-e(k-1,j,i) ) * ddzu(k) &
403	) * ddzw(k) &
404	- dissipation
405
406	!
407	!-- Store dissipation if needed for calculating the sgs particle
408	!-- velocities
409	IF ( use_sgs_for_particles .OR. wang_kernel .OR. &
410	turbulence ) THEN
411	diss(k,j,i) = dissipation
412	ENDIF
413
414	ENDIF
415
416	ENDDO
417	ENDDO
418	ENDDO
419	!$acc end kernels
420
421	ELSE
422
423	!$acc kernels present( ddzu, ddzw, dd2zu, diss, e, km, l_grid ) &
424	!$acc present( nzb_s_inner, tend, var, zu, zw )
425	DO i = i_left, i_right
426	DO j = j_south, j_north
427	DO k = 1, nzt
428
429	IF ( k > nzb_s_inner(j,i) ) THEN
430	!
431	!-- Calculate the mixing length (for dissipation)
432	dvar_dz = atmos_ocean_sign * &
433	( var(k+1,j,i) - var(k-1,j,i) ) * dd2zu(k)
434	IF ( dvar_dz > 0.0_wp ) THEN
435	l_stable = 0.76_wp * SQRT( e(k,j,i) ) / &
436	SQRT( g / var(k,j,i) * dvar_dz ) + 1E-5_wp
437	ELSE
438	l_stable = l_grid(k)
439	ENDIF
440	!
441	!-- Adjustment of the mixing length
442	IF ( wall_adjustment ) THEN
443	l = MIN( wall_adjustment_factor * &
444	( zu(k) - zw(nzb_s_inner(j,i)) ), &
445	l_grid(k), l_stable )
446	ll = MIN( wall_adjustment_factor * &
447	( zu(k) - zw(nzb_s_inner(j,i)) ), &
448	l_grid(k) )
449	ELSE
450	l = MIN( l_grid(k), l_stable )
451	ll = l_grid(k)
452	ENDIF
453	!
454	!-- Calculate the tendency terms
455	dissipation = ( 0.19_wp + 0.74_wp * l / ll ) * &
456	e(k,j,i) * SQRT( e(k,j,i) ) / l
457
458	tend(k,j,i) = tend(k,j,i) &
459	+ ( &
460	( km(k,j,i)+km(k,j,i+1) ) * ( e(k,j,i+1)-e(k,j,i) ) &
461	- ( km(k,j,i)+km(k,j,i-1) ) * ( e(k,j,i)-e(k,j,i-1) ) &
462	) * ddx2 &
463	+ ( &
464	( km(k,j,i)+km(k,j+1,i) ) * ( e(k,j+1,i)-e(k,j,i) ) &
465	- ( km(k,j,i)+km(k,j-1,i) ) * ( e(k,j,i)-e(k,j-1,i) ) &
466	) * ddy2 &
467	+ ( &
468	( km(k,j,i)+km(k+1,j,i) ) * ( e(k+1,j,i)-e(k,j,i) ) * ddzu(k+1) &
469	- ( km(k,j,i)+km(k-1,j,i) ) * ( e(k,j,i)-e(k-1,j,i) ) * ddzu(k) &
470	) * ddzw(k) &
471	- dissipation
472
473	!
474	!-- Store dissipation if needed for calculating the sgs
475	!-- particle velocities
476	IF ( use_sgs_for_particles .OR. wang_kernel .OR. &
477	turbulence ) THEN
478	diss(k,j,i) = dissipation
479	ENDIF
480
481	ENDIF
482
483	ENDDO
484	ENDDO
485	ENDDO
486	!$acc end kernels
487
488	ENDIF
489
490	!
491	!-- Boundary condition for dissipation
492	IF ( use_sgs_for_particles .OR. wang_kernel .OR. turbulence ) THEN
493	!$acc kernels present( diss, nzb_s_inner )
494	DO i = i_left, i_right
495	DO j = j_south, j_north
496	diss(nzb_s_inner(j,i),j,i) = diss(nzb_s_inner(j,i)+1,j,i)
497	ENDDO
498	ENDDO
499	!$acc end kernels
500	ENDIF
501
502	END SUBROUTINE diffusion_e_acc
503
504
505	!------------------------------------------------------------------------------!
506	! Call for grid point i,j
507	!------------------------------------------------------------------------------!
508	SUBROUTINE diffusion_e_ij( i, j, var, var_reference )
509
510	USE arrays_3d, &
511	ONLY: dd2zu, ddzu, ddzw, diss, e, km, l_grid, tend, zu, zw
512
513	USE control_parameters, &
514	ONLY: atmos_ocean_sign, g, turbulence, use_single_reference_value, &
515	wall_adjustment, wall_adjustment_factor
516
517	USE grid_variables, &
518	ONLY: ddx2, ddy2
519
520	USE indices, &
521	ONLY: nxlg, nxrg, nyng, nysg, nzb, nzb_s_inner, nzt
522
523	USE kinds
524
525	USE particle_attributes, &
526	ONLY: use_sgs_for_particles, wang_kernel
527
528	IMPLICIT NONE
529
530	INTEGER(iwp) :: i !:
531	INTEGER(iwp) :: j !:
532	INTEGER(iwp) :: k !:
533	REAL(wp) :: dvar_dz !:
534	REAL(wp) :: l_stable !:
535	REAL(wp) :: var_reference !:
536
537	#if defined( __nopointer )
538	REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,nxlg:nxrg) :: var !:
539	#else
540	REAL(wp), DIMENSION(:,:,:), POINTER :: var !:
541	#endif
542	REAL(wp), DIMENSION(nzb+1:nzt) :: dissipation !:
543	REAL(wp), DIMENSION(nzb+1:nzt) :: l !:
544	REAL(wp), DIMENSION(nzb+1:nzt) :: ll !:
545
546
547	!
548	!-- Calculate the mixing length (for dissipation)
549	DO k = nzb_s_inner(j,i)+1, nzt
550	dvar_dz = atmos_ocean_sign * &
551	( var(k+1,j,i) - var(k-1,j,i) ) * dd2zu(k)
552	IF ( dvar_dz > 0.0_wp ) THEN
553	IF ( use_single_reference_value ) THEN
554	l_stable = 0.76_wp * SQRT( e(k,j,i) ) / &
555	SQRT( g / var_reference * dvar_dz ) + 1E-5_wp
556	ELSE
557	l_stable = 0.76_wp * SQRT( e(k,j,i) ) / &
558	SQRT( g / var(k,j,i) * dvar_dz ) + 1E-5_wp
559	ENDIF
560	ELSE
561	l_stable = l_grid(k)
562	ENDIF
563	!
564	!-- Adjustment of the mixing length
565	IF ( wall_adjustment ) THEN
566	l(k) = MIN( wall_adjustment_factor * &
567	( zu(k) - zw(nzb_s_inner(j,i)) ), l_grid(k), &
568	l_stable )
569	ll(k) = MIN( wall_adjustment_factor * &
570	( zu(k) - zw(nzb_s_inner(j,i)) ), l_grid(k) )
571	ELSE
572	l(k) = MIN( l_grid(k), l_stable )
573	ll(k) = l_grid(k)
574	ENDIF
575	!
576	!-- Calculate the tendency term
577	dissipation(k) = ( 0.19_wp + 0.74_wp * l(k) / ll(k) ) * e(k,j,i) * &
578	SQRT( e(k,j,i) ) / l(k)
579
580	tend(k,j,i) = tend(k,j,i) &
581	+ ( &
582	( km(k,j,i)+km(k,j,i+1) ) * ( e(k,j,i+1)-e(k,j,i) ) &
583	- ( km(k,j,i)+km(k,j,i-1) ) * ( e(k,j,i)-e(k,j,i-1) ) &
584	) * ddx2 &
585	+ ( &
586	( km(k,j,i)+km(k,j+1,i) ) * ( e(k,j+1,i)-e(k,j,i) ) &
587	- ( km(k,j,i)+km(k,j-1,i) ) * ( e(k,j,i)-e(k,j-1,i) ) &
588	) * ddy2 &
589	+ ( &
590	( km(k,j,i)+km(k+1,j,i) ) * ( e(k+1,j,i)-e(k,j,i) ) * ddzu(k+1) &
591	- ( km(k,j,i)+km(k-1,j,i) ) * ( e(k,j,i)-e(k-1,j,i) ) * ddzu(k) &
592	) * ddzw(k) &
593	- dissipation(k)
594
595	ENDDO
596
597	!
598	!-- Store dissipation if needed for calculating the sgs particle velocities
599	IF ( use_sgs_for_particles .OR. wang_kernel .OR. turbulence ) THEN
600	DO k = nzb_s_inner(j,i)+1, nzt
601	diss(k,j,i) = dissipation(k)
602	ENDDO
603	!
604	!-- Boundary condition for dissipation
605	diss(nzb_s_inner(j,i),j,i) = diss(nzb_s_inner(j,i)+1,j,i)
606	ENDIF
607
608	END SUBROUTINE diffusion_e_ij
609
610	END MODULE diffusion_e_mod

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |