Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |