Home

Context Navigation

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

Last change on this file since 2000 was 2000, checked in by knoop, 8 years ago
Forced header and separation lines into 80 columns
Property svn:keywords set to `Id`
File size: 25.5 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |