Home

Context Navigation

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

Last change on this file since 1104 was 1066, checked in by hoffmann, 12 years ago
last commit documented
Property svn:keywords set to `Id`
File size: 21.5 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |