Home

Context Navigation

source: palm/trunk/SOURCE/inflow_turbulence.f90 @ 1560

Last change on this file since 1560 was 1560, checked in by keck, 9 years ago
implemented possibility of adding a y shift to the recycled inflow turbulence
Property svn:keywords set to `Id`
File size: 9.7 KB

Line
1	SUBROUTINE inflow_turbulence
2
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	!
17	! Copyright 1997-2014 Leibniz Universitaet Hannover
18	!--------------------------------------------------------------------------------!
19	!
20	! Current revisions:
21	! -----------------
22	! Option recycling_yshift added. If this option is switched on, the turbulence
23	! data, which is mapped from the recycling plane to the inflow, is shifted in
24	! y direction (by ny * dy / 2 )
25	!
26	! Former revisions:
27	! -----------------
28	! $Id: inflow_turbulence.f90 1560 2015-03-06 10:48:54Z keck $
29	!
30	! 1353 2014-04-08 15:21:23Z heinze
31	! REAL constants provided with KIND-attribute
32	!
33	! 1346 2014-03-27 13:18:20Z heinze
34	! Bugfix: REAL constants provided with KIND-attribute especially in call of
35	! intrinsic function like MAX, MIN, SIGN
36	!
37	! 1320 2014-03-20 08:40:49Z raasch
38	! ONLY-attribute added to USE-statements,
39	! kind-parameters added to all INTEGER and REAL declaration statements,
40	! kinds are defined in new module kinds,
41	! revision history before 2012 removed,
42	! comment fields (!:) to be used for variable explanations added to
43	! all variable declaration statements
44	!
45	! 1092 2013-02-02 11:24:22Z raasch
46	! unused variables removed
47	!
48	! 1036 2012-10-22 13:43:42Z raasch
49	! code put under GPL (PALM 3.9)
50	!
51	! Initial version (2008/03/07)
52	!
53	! Description:
54	! ------------
55	! Imposing turbulence at the respective inflow using the turbulence
56	! recycling method of Kataoka and Mizuno (2002).
57	!------------------------------------------------------------------------------!
58
59	USE arrays_3d, &
60	ONLY: e, inflow_damping_factor, mean_inflow_profiles, pt, u, v, w
61
62	USE control_parameters, &
63	ONLY: recycling_plane, recycling_yshift
64
65	USE cpulog, &
66	ONLY: cpu_log, log_point
67
68	USE grid_variables, &
69	ONLY:
70
71	USE indices, &
72	ONLY: nbgp, nxl, ny, nyn, nys, nyng, nysg, nzb, nzt
73
74	USE kinds
75
76	USE pegrid
77
78
79	IMPLICIT NONE
80
81	INTEGER(iwp) :: i !:
82	INTEGER(iwp) :: j !:
83	INTEGER(iwp) :: k !:
84	INTEGER(iwp) :: l !:
85	INTEGER(iwp) :: next !:
86	INTEGER(iwp) :: ngp_ifd !:
87	INTEGER(iwp) :: ngp_pr !:
88	INTEGER(iwp) :: prev !:
89
90	REAL(wp), DIMENSION(nzb:nzt+1,5,nbgp) :: &
91	avpr, avpr_l !:
92	REAL(wp), DIMENSION(nzb:nzt+1,nysg:nyng,5,nbgp) :: &
93	inflow_dist, local_inflow_dist !:
94
95	CALL cpu_log( log_point(40), 'inflow_turbulence', 'start' )
96
97	!
98	!-- Carry out spanwise averaging in the recycling plane
99	avpr_l = 0.0_wp
100	ngp_pr = ( nzt - nzb + 2 ) * 5 * nbgp
101	ngp_ifd = ngp_pr * ( nyn - nys + 1 + 2 * nbgp )
102
103	!
104	!-- First, local averaging within the recycling domain
105	i = recycling_plane
106
107	#if defined( __parallel )
108	IF ( myidx == id_recycling ) THEN
109
110	DO l = 1, nbgp
111	DO j = nys, nyn
112	DO k = nzb, nzt + 1
113
114	avpr_l(k,1,l) = avpr_l(k,1,l) + u(k,j,i)
115	avpr_l(k,2,l) = avpr_l(k,2,l) + v(k,j,i)
116	avpr_l(k,3,l) = avpr_l(k,3,l) + w(k,j,i)
117	avpr_l(k,4,l) = avpr_l(k,4,l) + pt(k,j,i)
118	avpr_l(k,5,l) = avpr_l(k,5,l) + e(k,j,i)
119
120	ENDDO
121	ENDDO
122	i = i + 1
123	ENDDO
124
125	ENDIF
126	!
127	!-- Now, averaging over all PEs
128	IF ( collective_wait ) CALL MPI_BARRIER( comm2d, ierr )
129	CALL MPI_ALLREDUCE( avpr_l(nzb,1,1), avpr(nzb,1,1), ngp_pr, MPI_REAL, &
130	MPI_SUM, comm2d, ierr )
131
132	#else
133	DO l = 1, nbgp
134	DO j = nys, nyn
135	DO k = nzb, nzt + 1
136
137	avpr_l(k,1,l) = avpr_l(k,1,l) + u(k,j,i)
138	avpr_l(k,2,l) = avpr_l(k,2,l) + v(k,j,i)
139	avpr_l(k,3,l) = avpr_l(k,3,l) + w(k,j,i)
140	avpr_l(k,4,l) = avpr_l(k,4,l) + pt(k,j,i)
141	avpr_l(k,5,l) = avpr_l(k,5,l) + e(k,j,i)
142
143	ENDDO
144	ENDDO
145	i = i + 1
146	ENDDO
147
148	avpr = avpr_l
149	#endif
150
151	avpr = avpr / ( ny + 1 )
152	!
153	!-- Calculate the disturbances at the recycling plane
154	i = recycling_plane
155
156	#if defined( __parallel )
157	IF ( myidx == id_recycling ) THEN
158	DO l = 1, nbgp
159	DO j = nysg, nyng
160	DO k = nzb, nzt + 1
161
162	inflow_dist(k,j,1,l) = u(k,j,i+1) - avpr(k,1,l)
163	inflow_dist(k,j,2,l) = v(k,j,i) - avpr(k,2,l)
164	inflow_dist(k,j,3,l) = w(k,j,i) - avpr(k,3,l)
165	inflow_dist(k,j,4,l) = pt(k,j,i) - avpr(k,4,l)
166	inflow_dist(k,j,5,l) = e(k,j,i) - avpr(k,5,l)
167
168	ENDDO
169	ENDDO
170	i = i + 1
171	ENDDO
172
173	ENDIF
174	#else
175	DO l = 1, nbgp
176	DO j = nysg, nyng
177	DO k = nzb, nzt+1
178
179	inflow_dist(k,j,1,l) = u(k,j,i+1) - avpr(k,1,l)
180	inflow_dist(k,j,2,l) = v(k,j,i) - avpr(k,2,l)
181	inflow_dist(k,j,3,l) = w(k,j,i) - avpr(k,3,l)
182	inflow_dist(k,j,4,l) = pt(k,j,i) - avpr(k,4,l)
183	inflow_dist(k,j,5,l) = e(k,j,i) - avpr(k,5,l)
184
185	ENDDO
186	ENDDO
187	i = i + 1
188	ENDDO
189	#endif
190
191	!
192	!-- For parallel runs, send the disturbances to the respective inflow PE
193	#if defined( __parallel )
194	IF ( myidx == id_recycling .AND. myidx /= id_inflow ) THEN
195
196	CALL MPI_SEND( inflow_dist(nzb,nysg,1,1), ngp_ifd, MPI_REAL, &
197	id_inflow, 1, comm1dx, ierr )
198
199	ELSEIF ( myidx /= id_recycling .AND. myidx == id_inflow ) THEN
200
201	inflow_dist = 0.0_wp
202	CALL MPI_RECV( inflow_dist(nzb,nysg,1,1), ngp_ifd, MPI_REAL, &
203	id_recycling, 1, comm1dx, status, ierr )
204
205	ENDIF
206
207
208	IF ( recycling_yshift .AND. myidx == id_inflow ) THEN
209
210	IF ( pdims(2) >= 2 ) THEN
211
212	IF ( myidy >= INT( pdims(2) / 2 ) ) THEN
213	prev = myidy - INT( pdims(2) / 2 )
214	ELSE
215	prev = pdims(2) - ( INT( pdims(2) / 2 ) - myidy )
216	ENDIF
217
218	IF ( myidy < pdims(2) - INT( pdims(2) / 2 ) ) THEN
219	next = myidy + INT( pdims(2) / 2 )
220	ELSE
221	next = INT( pdims(2) / 2 ) - ( pdims(2) - myidy )
222	ENDIF
223
224	ENDIF
225
226	local_inflow_dist = 0.0_wp
227
228	CALL MPI_SENDRECV( inflow_dist(nzb,nysg,1,1), ngp_ifd, MPI_REAL, &
229	next, 1, local_inflow_dist(nzb,nysg,1,1), ngp_ifd, &
230	MPI_REAL, prev, 1, comm1dy, status, ierr )
231
232	ENDIF
233
234	#endif
235
236	!
237	!-- Add the disturbance at the inflow
238	IF ( nxl == 0 ) THEN
239
240	IF ( recycling_yshift ) THEN
241
242	DO j = nysg, nyng
243	DO k = nzb, nzt + 1
244
245	u(k,j,-nbgp+1:0) = mean_inflow_profiles(k,1) + &
246	local_inflow_dist(k,j,1,1:nbgp) * &
247	inflow_damping_factor(k)
248	v(k,j,-nbgp:-1) = mean_inflow_profiles(k,2) + &
249	local_inflow_dist(k,j,2,1:nbgp) * &
250	inflow_damping_factor(k)
251	w(k,j,-nbgp:-1) = &
252	local_inflow_dist(k,j,3,1:nbgp) * &
253	inflow_damping_factor(k)
254	pt(k,j,-nbgp:-1) = mean_inflow_profiles(k,4) + &
255	local_inflow_dist(k,j,4,1:nbgp) * &
256	inflow_damping_factor(k)
257	e(k,j,-nbgp:-1) = mean_inflow_profiles(k,5) + &
258	local_inflow_dist(k,j,5,1:nbgp) * &
259	inflow_damping_factor(k)
260	e(k,j,-nbgp:-1) = MAX( e(k,j,-nbgp:-1), 0.0_wp )
261
262	ENDDO
263	ENDDO
264
265	ELSE
266
267	DO j = nysg, nyng
268	DO k = nzb, nzt + 1
269
270	u(k,j,-nbgp+1:0) = mean_inflow_profiles(k,1) + &
271	inflow_dist(k,j,1,1:nbgp) * inflow_damping_factor(k)
272	v(k,j,-nbgp:-1) = mean_inflow_profiles(k,2) + &
273	inflow_dist(k,j,2,1:nbgp) * inflow_damping_factor(k)
274	w(k,j,-nbgp:-1) = &
275	inflow_dist(k,j,3,1:nbgp) * inflow_damping_factor(k)
276	pt(k,j,-nbgp:-1) = mean_inflow_profiles(k,4) + &
277	inflow_dist(k,j,4,1:nbgp) * inflow_damping_factor(k)
278	e(k,j,-nbgp:-1) = mean_inflow_profiles(k,5) + &
279	inflow_dist(k,j,5,1:nbgp) * inflow_damping_factor(k)
280	e(k,j,-nbgp:-1) = MAX( e(k,j,-nbgp:-1), 0.0_wp )
281
282	ENDDO
283	ENDDO
284
285	ENDIF
286
287	ENDIF
288
289
290	CALL cpu_log( log_point(40), 'inflow_turbulence', 'stop' )
291
292
293	END SUBROUTINE inflow_turbulence

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |