Home

Context Navigation

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |