Home

Context Navigation

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

Last change on this file since 2016 was 2001, checked in by knoop, 8 years ago
last commit documented
Property svn:keywords set to `Id`
File size: 25.5 KB

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |