Home

Context Navigation

diffusion_e_mod.f90 @ 1867

Last change on this file since 1867 was 1851, checked in by maronga, 9 years ago
last commit documented
Property svn:keywords set to `Id`
File size: 25.4 KB

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

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

Context Navigation

source: palm/trunk/SOURCE/diffusion_e_mod.f90 @ 1867

Download in other formats: