Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |