Home

Context Navigation

source: palm/trunk/SOURCE/lpm_init_sgs_tke.f90 @ 2270

Last change on this file since 2270 was 2233, checked in by suehring, 7 years ago
last commit documented
Property svn:keywords set to `Id`
File size: 11.5 KB

Rev	Line
[1682]	1	!> @file lpm_init_sgs_tke.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	!
[849]	20	! Current revisions:
	21	! ------------------
[1360]	22	!
[2233]	23	!
[1321]	24	! Former revisions:
	25	! -----------------
	26	! $Id: lpm_init_sgs_tke.f90 2233 2017-05-30 18:08:54Z maronga $
	27	!
[2233]	28	! 2232 2017-05-30 17:47:52Z suehring
	29	! Adjustments according to new topography realization
	30	!
[2001]	31	! 2000 2016-08-20 18:09:15Z knoop
	32	! Forced header and separation lines into 80 columns
	33	!
[1930]	34	! 1929 2016-06-09 16:25:25Z suehring
	35	! sgs_wfu_par, sgs_wfv_par and sgs_wfw_par are replaced by sgs_wf_par
	36	!
[1823]	37	! 1822 2016-04-07 07:49:42Z hoffmann
	38	! Unused variables removed.
	39	!
[1683]	40	! 1682 2015-10-07 23:56:08Z knoop
	41	! Code annotations made doxygen readable
	42	!
[1360]	43	! 1359 2014-04-11 17:15:14Z hoffmann
	44	! New particle structure integrated.
	45	! Kind definition added to all floating point numbers.
	46	!
[1321]	47	! 1320 2014-03-20 08:40:49Z raasch
[1320]	48	! ONLY-attribute added to USE-statements,
	49	! kind-parameters added to all INTEGER and REAL declaration statements,
	50	! kinds are defined in new module kinds,
	51	! comment fields (!:) to be used for variable explanations added to
	52	! all variable declaration statements
[849]	53	!
[1037]	54	! 1036 2012-10-22 13:43:42Z raasch
	55	! code put under GPL (PALM 3.9)
	56	!
[850]	57	! 849 2012-03-15 10:35:09Z raasch
	58	! initial revision (former part of advec_particles)
[849]	59	!
[850]	60	!
[849]	61	! Description:
	62	! ------------
[1682]	63	!> Calculates quantities required for considering the SGS velocity fluctuations
	64	!> in the particle transport by a stochastic approach. The respective
	65	!> quantities are: SGS-TKE gradients and horizontally averaged profiles of the
	66	!> SGS TKE and the resolved-scale velocity variances.
[849]	67	!------------------------------------------------------------------------------!
[1682]	68	SUBROUTINE lpm_init_sgs_tke
	69
[849]	70
[1320]	71	USE arrays_3d, &
	72	ONLY: de_dx, de_dy, de_dz, diss, e, u, v, w, zu
	73
	74	USE grid_variables, &
	75	ONLY: ddx, ddy
	76
	77	USE indices, &
[2232]	78	ONLY: nbgp, ngp_2dh_outer, nxl, nxr, nyn, nys, nzb, &
	79	nzb_s_outer, nzt, wall_flags_0
[1320]	80
	81	USE kinds
	82
	83	USE particle_attributes, &
[1929]	84	ONLY: sgs_wf_part
[1320]	85
[849]	86	USE pegrid
	87
[1320]	88	USE statistics, &
	89	ONLY: flow_statistics_called, hom, sums, sums_l
	90
[2232]	91	USE surface_mod, &
	92	ONLY: bc_h
	93
[849]	94	IMPLICIT NONE
	95
[2232]	96	INTEGER(iwp) :: i !< index variable along x
	97	INTEGER(iwp) :: j !< index variable along y
	98	INTEGER(iwp) :: k !< index variable along z
	99	INTEGER(iwp) :: m !< running index for the surface elements
[849]	100
[2232]	101	REAL(wp) :: flag1 !< flag to mask topography
[849]	102
	103	!
	104	!-- TKE gradient along x and y
	105	DO i = nxl, nxr
	106	DO j = nys, nyn
	107	DO k = nzb, nzt+1
	108
[2232]	109	IF ( .NOT. BTEST( wall_flags_0(k,j,i-1), 0 ) .AND. &
	110	BTEST( wall_flags_0(k,j,i), 0 ) .AND. &
	111	BTEST( wall_flags_0(k,j,i+1), 0 ) ) &
[849]	112	THEN
[1929]	113	de_dx(k,j,i) = 2.0_wp * sgs_wf_part * &
[1359]	114	( e(k,j,i+1) - e(k,j,i) ) * ddx
[2232]	115	ELSEIF ( BTEST( wall_flags_0(k,j,i-1), 0 ) .AND. &
	116	BTEST( wall_flags_0(k,j,i), 0 ) .AND. &
	117	.NOT. BTEST( wall_flags_0(k,j,i+1), 0 ) ) &
[849]	118	THEN
[1929]	119	de_dx(k,j,i) = 2.0_wp * sgs_wf_part * &
[1359]	120	( e(k,j,i) - e(k,j,i-1) ) * ddx
[2232]	121	ELSEIF ( .NOT. BTEST( wall_flags_0(k,j,i), 22 ) .AND. &
	122	.NOT. BTEST( wall_flags_0(k,j,i+1), 22 ) ) &
[849]	123	THEN
[1359]	124	de_dx(k,j,i) = 0.0_wp
[2232]	125	ELSEIF ( .NOT. BTEST( wall_flags_0(k,j,i-1), 22 ) .AND. &
	126	.NOT. BTEST( wall_flags_0(k,j,i), 22 ) ) &
[849]	127	THEN
[1359]	128	de_dx(k,j,i) = 0.0_wp
[849]	129	ELSE
[1929]	130	de_dx(k,j,i) = sgs_wf_part * ( e(k,j,i+1) - e(k,j,i-1) ) * ddx
[849]	131	ENDIF
	132
[2232]	133	IF ( .NOT. BTEST( wall_flags_0(k,j-1,i), 0 ) .AND. &
	134	BTEST( wall_flags_0(k,j,i), 0 ) .AND. &
	135	BTEST( wall_flags_0(k,j+1,i), 0 ) ) &
[849]	136	THEN
[1929]	137	de_dy(k,j,i) = 2.0_wp * sgs_wf_part * &
[1359]	138	( e(k,j+1,i) - e(k,j,i) ) * ddy
[2232]	139	ELSEIF ( BTEST( wall_flags_0(k,j-1,i), 0 ) .AND. &
	140	BTEST( wall_flags_0(k,j,i), 0 ) .AND. &
	141	.NOT. BTEST( wall_flags_0(k,j+1,i), 0 ) ) &
[849]	142	THEN
[1929]	143	de_dy(k,j,i) = 2.0_wp * sgs_wf_part * &
[1359]	144	( e(k,j,i) - e(k,j-1,i) ) * ddy
[2232]	145	ELSEIF ( .NOT. BTEST( wall_flags_0(k,j,i), 22 ) .AND. &
	146	.NOT. BTEST( wall_flags_0(k,j+1,i), 22 ) ) &
[849]	147	THEN
[1359]	148	de_dy(k,j,i) = 0.0_wp
[2232]	149	ELSEIF ( .NOT. BTEST( wall_flags_0(k,j-1,i), 22 ) .AND. &
	150	.NOT. BTEST( wall_flags_0(k,j,i), 22 ) ) &
[849]	151	THEN
[1359]	152	de_dy(k,j,i) = 0.0_wp
[849]	153	ELSE
[1929]	154	de_dy(k,j,i) = sgs_wf_part * ( e(k,j+1,i) - e(k,j-1,i) ) * ddy
[849]	155	ENDIF
	156
	157	ENDDO
	158	ENDDO
	159	ENDDO
	160
	161	!
[2232]	162	!-- TKE gradient along z at topograhy and including bottom and top boundary conditions
[849]	163	DO i = nxl, nxr
	164	DO j = nys, nyn
[2232]	165	DO k = nzb+1, nzt-1
	166	!
	167	!-- Flag to mask topography
	168	flag1 = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_0(k,j,i), 0 ) )
[849]	169
[2232]	170	de_dz(k,j,i) = 2.0_wp * sgs_wf_part * &
	171	( e(k+1,j,i) - e(k-1,j,i) ) / ( zu(k+1) - zu(k-1) ) &
	172	* flag1
[849]	173	ENDDO
[2232]	174	!
	175	!-- upward-facing surfaces
	176	DO m = bc_h(0)%start_index(j,i), bc_h(0)%end_index(j,i)
	177	k = bc_h(0)%k(m)
	178	de_dz(k,j,i) = 2.0_wp * sgs_wf_part * &
	179	( e(k+1,j,i) - e(k,j,i) ) / ( zu(k+1) - zu(k) )
	180	ENDDO
	181	!
	182	!-- downward-facing surfaces
	183	DO m = bc_h(1)%start_index(j,i), bc_h(1)%end_index(j,i)
	184	k = bc_h(1)%k(m)
	185	de_dz(k,j,i) = 2.0_wp * sgs_wf_part * &
	186	( e(k,j,i) - e(k-1,j,i) ) / ( zu(k) - zu(k-1) )
	187	ENDDO
[849]	188
[2232]	189	de_dz(nzb,j,i) = 0.0_wp
[1359]	190	de_dz(nzt,j,i) = 0.0_wp
	191	de_dz(nzt+1,j,i) = 0.0_wp
[849]	192	ENDDO
	193	ENDDO
	194	!
[2232]	195	!-- Ghost point exchange
[849]	196	CALL exchange_horiz( de_dx, nbgp )
	197	CALL exchange_horiz( de_dy, nbgp )
	198	CALL exchange_horiz( de_dz, nbgp )
	199	CALL exchange_horiz( diss, nbgp )
	200
	201
	202	!
	203	!-- Calculate the horizontally averaged profiles of SGS TKE and resolved
	204	!-- velocity variances (they may have been already calculated in routine
	205	!-- flow_statistics).
	206	IF ( .NOT. flow_statistics_called ) THEN
	207
	208	!
	209	!-- First calculate horizontally averaged profiles of the horizontal
	210	!-- velocities.
[1359]	211	sums_l(:,1,0) = 0.0_wp
	212	sums_l(:,2,0) = 0.0_wp
[849]	213
	214	DO i = nxl, nxr
	215	DO j = nys, nyn
[2232]	216	DO k = nzb, nzt+1
	217	!
	218	!-- Flag indicate nzb_s_outer
	219	flag1 = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_0(k,j,i), 24 ) )
	220
	221	sums_l(k,1,0) = sums_l(k,1,0) + u(k,j,i) * flag1
	222	sums_l(k,2,0) = sums_l(k,2,0) + v(k,j,i) * flag1
[849]	223	ENDDO
	224	ENDDO
	225	ENDDO
	226
	227	#if defined( __parallel )
	228	!
	229	!-- Compute total sum from local sums
	230	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
	231	CALL MPI_ALLREDUCE( sums_l(nzb,1,0), sums(nzb,1), nzt+2-nzb, &
	232	MPI_REAL, MPI_SUM, comm2d, ierr )
	233	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
	234	CALL MPI_ALLREDUCE( sums_l(nzb,2,0), sums(nzb,2), nzt+2-nzb, &
	235	MPI_REAL, MPI_SUM, comm2d, ierr )
	236	#else
	237	sums(:,1) = sums_l(:,1,0)
	238	sums(:,2) = sums_l(:,2,0)
	239	#endif
	240
	241	!
	242	!-- Final values are obtained by division by the total number of grid
	243	!-- points used for the summation.
	244	hom(:,1,1,0) = sums(:,1) / ngp_2dh_outer(:,0) ! u
	245	hom(:,1,2,0) = sums(:,2) / ngp_2dh_outer(:,0) ! v
	246
	247	!
	248	!-- Now calculate the profiles of SGS TKE and the resolved-scale
	249	!-- velocity variances
[1359]	250	sums_l(:,8,0) = 0.0_wp
	251	sums_l(:,30,0) = 0.0_wp
	252	sums_l(:,31,0) = 0.0_wp
	253	sums_l(:,32,0) = 0.0_wp
[849]	254	DO i = nxl, nxr
	255	DO j = nys, nyn
[2232]	256	DO k = nzb, nzt+1
	257	!
	258	!-- Flag indicate nzb_s_outer
	259	flag1 = MERGE( 1.0_wp, 0.0_wp, BTEST( wall_flags_0(k,j,i), 24 ) )
	260
	261	sums_l(k,8,0) = sums_l(k,8,0) + e(k,j,i) * flag1
	262	sums_l(k,30,0) = sums_l(k,30,0) + ( u(k,j,i) - hom(k,1,1,0) )*2 flag1
	263	sums_l(k,31,0) = sums_l(k,31,0) + ( v(k,j,i) - hom(k,1,2,0) )*2 flag1
	264	sums_l(k,32,0) = sums_l(k,32,0) + w(k,j,i)*2 flag1
[849]	265	ENDDO
	266	ENDDO
	267	ENDDO
	268
	269	#if defined( __parallel )
	270	!
	271	!-- Compute total sum from local sums
	272	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
	273	CALL MPI_ALLREDUCE( sums_l(nzb,8,0), sums(nzb,8), nzt+2-nzb, &
	274	MPI_REAL, MPI_SUM, comm2d, ierr )
	275	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
	276	CALL MPI_ALLREDUCE( sums_l(nzb,30,0), sums(nzb,30), nzt+2-nzb, &
	277	MPI_REAL, MPI_SUM, comm2d, ierr )
	278	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
	279	CALL MPI_ALLREDUCE( sums_l(nzb,31,0), sums(nzb,31), nzt+2-nzb, &
	280	MPI_REAL, MPI_SUM, comm2d, ierr )
	281	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
	282	CALL MPI_ALLREDUCE( sums_l(nzb,32,0), sums(nzb,32), nzt+2-nzb, &
	283	MPI_REAL, MPI_SUM, comm2d, ierr )
	284
	285	#else
	286	sums(:,8) = sums_l(:,8,0)
	287	sums(:,30) = sums_l(:,30,0)
	288	sums(:,31) = sums_l(:,31,0)
	289	sums(:,32) = sums_l(:,32,0)
	290	#endif
	291
	292	!
	293	!-- Final values are obtained by division by the total number of grid
	294	!-- points used for the summation.
	295	hom(:,1,8,0) = sums(:,8) / ngp_2dh_outer(:,0) ! e
	296	hom(:,1,30,0) = sums(:,30) / ngp_2dh_outer(:,0) ! u*2
	297	hom(:,1,31,0) = sums(:,31) / ngp_2dh_outer(:,0) ! v*2
	298	hom(:,1,32,0) = sums(:,32) / ngp_2dh_outer(:,0) ! w*2
	299
	300	ENDIF
	301
	302	END SUBROUTINE lpm_init_sgs_tke

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |