Home

Context Navigation

source: palm/trunk/SOURCE/lpm_boundary_conds.f90 @ 1822

Last change on this file since 1822 was 1822, checked in by hoffmann, 8 years ago
changes in LPM and bulk cloud microphysics
Property svn:keywords set to `Id`
File size: 25.9 KB

Line
1	!> @file lpm_boundary_conds.f90
2	!--------------------------------------------------------------------------------!
3	! This file is part of PALM.
4	!
5	! PALM is free software: you can redistribute it and/or modify it under the terms
6	! of the GNU General Public License as published by the Free Software Foundation,
7	! either version 3 of the License, or (at your option) any later version.
8	!
9	! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
10	! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
11	! A PARTICULAR PURPOSE. See the GNU General Public License for more details.
12	!
13	! You should have received a copy of the GNU General Public License along with
14	! PALM. If not, see <http://www.gnu.org/licenses/>.
15	!
16	! Copyright 1997-2016 Leibniz Universitaet Hannover
17	!--------------------------------------------------------------------------------!
18	!
19	! Current revisions:
20	! -----------------
21	! Tails removed. Unused variables removed.
22	!
23	! Former revisions:
24	! -----------------
25	! $Id: lpm_boundary_conds.f90 1822 2016-04-07 07:49:42Z hoffmann $
26	!
27	! 1682 2015-10-07 23:56:08Z knoop
28	! Code annotations made doxygen readable
29	!
30	! 1359 2014-04-11 17:15:14Z hoffmann
31	! New particle structure integrated.
32	! Kind definition added to all floating point numbers.
33	!
34	! 1320 2014-03-20 08:40:49Z raasch
35	! ONLY-attribute added to USE-statements,
36	! kind-parameters added to all INTEGER and REAL declaration statements,
37	! kinds are defined in new module kinds,
38	! revision history before 2012 removed,
39	! comment fields (!:) to be used for variable explanations added to
40	! all variable declaration statements
41	!
42	! 1036 2012-10-22 13:43:42Z raasch
43	! code put under GPL (PALM 3.9)
44	!
45	! 849 2012-03-15 10:35:09Z raasch
46	! routine renamed lpm_boundary_conds, bottom and top boundary conditions
47	! included (former part of advec_particles)
48	!
49	! 824 2012-02-17 09:09:57Z raasch
50	! particle attributes speed_x\|y\|z_sgs renamed rvar1\|2\|3
51	!
52	! Initial version (2007/03/09)
53	!
54	! Description:
55	! ------------
56	!> Boundary conditions for the Lagrangian particles.
57	!> The routine consists of two different parts. One handles the bottom (flat)
58	!> and top boundary. In this part, also particles which exceeded their lifetime
59	!> are deleted.
60	!> The other part handles the reflection of particles from vertical walls.
61	!> This part was developed by Jin Zhang during 2006-2007.
62	!>
63	!> To do: Code structure for finding the t_index values and for checking the
64	!> ----- reflection conditions is basically the same for all four cases, so it
65	!> should be possible to further simplify/shorten it.
66	!>
67	!> THE WALLS PART OF THIS ROUTINE HAS NOT BEEN TESTED FOR OCEAN RUNS SO FAR!!!!
68	!> (see offset_ocean_*)
69	!------------------------------------------------------------------------------!
70	SUBROUTINE lpm_boundary_conds( range )
71
72
73	USE arrays_3d, &
74	ONLY: zu, zw
75
76	USE control_parameters, &
77	ONLY: dz, message_string, particle_maximum_age
78
79	USE cpulog, &
80	ONLY: cpu_log, log_point_s
81
82	USE grid_variables, &
83	ONLY: ddx, dx, ddy, dy
84
85	USE indices, &
86	ONLY: nxl, nxr, nyn, nys, nz, nzb_s_inner
87
88	USE kinds
89
90	USE particle_attributes, &
91	ONLY: deleted_particles, ibc_par_b, ibc_par_t, number_of_particles, &
92	particles, particle_type, offset_ocean_nzt_m1, &
93	use_sgs_for_particles
94
95	USE pegrid
96
97	IMPLICIT NONE
98
99	CHARACTER (LEN=*) :: range !<
100
101	INTEGER(iwp) :: i !<
102	INTEGER(iwp) :: inc !<
103	INTEGER(iwp) :: ir !<
104	INTEGER(iwp) :: i1 !<
105	INTEGER(iwp) :: i2 !<
106	INTEGER(iwp) :: i3 !<
107	INTEGER(iwp) :: i5 !<
108	INTEGER(iwp) :: j !<
109	INTEGER(iwp) :: jr !<
110	INTEGER(iwp) :: j1 !<
111	INTEGER(iwp) :: j2 !<
112	INTEGER(iwp) :: j3 !<
113	INTEGER(iwp) :: j5 !<
114	INTEGER(iwp) :: k !<
115	INTEGER(iwp) :: k1 !<
116	INTEGER(iwp) :: k2 !<
117	INTEGER(iwp) :: k3 !<
118	INTEGER(iwp) :: k5 !<
119	INTEGER(iwp) :: n !<
120	INTEGER(iwp) :: t_index !<
121	INTEGER(iwp) :: t_index_number !<
122
123	LOGICAL :: reflect_x !<
124	LOGICAL :: reflect_y !<
125	LOGICAL :: reflect_z !<
126
127	REAL(wp) :: dt_particle !<
128	REAL(wp) :: pos_x !<
129	REAL(wp) :: pos_x_old !<
130	REAL(wp) :: pos_y !<
131	REAL(wp) :: pos_y_old !<
132	REAL(wp) :: pos_z !<
133	REAL(wp) :: pos_z_old !<
134	REAL(wp) :: prt_x !<
135	REAL(wp) :: prt_y !<
136	REAL(wp) :: prt_z !<
137	REAL(wp) :: t(1:200) !<
138	REAL(wp) :: tmp_t !<
139	REAL(wp) :: xline !<
140	REAL(wp) :: yline !<
141	REAL(wp) :: zline !<
142
143	IF ( range == 'bottom/top' ) THEN
144
145	!
146	!-- Apply boundary conditions to those particles that have crossed the top or
147	!-- bottom boundary and delete those particles, which are older than allowed
148	DO n = 1, number_of_particles
149
150	!
151	!-- Stop if particles have moved further than the length of one
152	!-- PE subdomain (newly released particles have age = age_m!)
153	IF ( particles(n)%age /= particles(n)%age_m ) THEN
154	IF ( ABS(particles(n)%speed_x) > &
155	((nxr-nxl+2)*dx)/(particles(n)%age-particles(n)%age_m) .OR. &
156	ABS(particles(n)%speed_y) > &
157	((nyn-nys+2)*dy)/(particles(n)%age-particles(n)%age_m) ) THEN
158
159	WRITE( message_string, * ) 'particle too fast. n = ', n
160	CALL message( 'lpm_boundary_conds', 'PA0148', 2, 2, -1, 6, 1 )
161	ENDIF
162	ENDIF
163
164	IF ( particles(n)%age > particle_maximum_age .AND. &
165	particles(n)%particle_mask ) &
166	THEN
167	particles(n)%particle_mask = .FALSE.
168	deleted_particles = deleted_particles + 1
169	ENDIF
170
171	IF ( particles(n)%z >= zu(nz) .AND. particles(n)%particle_mask ) THEN
172	IF ( ibc_par_t == 1 ) THEN
173	!
174	!-- Particle absorption
175	particles(n)%particle_mask = .FALSE.
176	deleted_particles = deleted_particles + 1
177	ELSEIF ( ibc_par_t == 2 ) THEN
178	!
179	!-- Particle reflection
180	particles(n)%z = 2.0_wp * zu(nz) - particles(n)%z
181	particles(n)%speed_z = -particles(n)%speed_z
182	IF ( use_sgs_for_particles .AND. &
183	particles(n)%rvar3 > 0.0_wp ) THEN
184	particles(n)%rvar3 = -particles(n)%rvar3
185	ENDIF
186	ENDIF
187	ENDIF
188
189	IF ( particles(n)%z < zw(0) .AND. particles(n)%particle_mask ) THEN
190	IF ( ibc_par_b == 1 ) THEN
191	!
192	!-- Particle absorption
193	particles(n)%particle_mask = .FALSE.
194	deleted_particles = deleted_particles + 1
195	ELSEIF ( ibc_par_b == 2 ) THEN
196	!
197	!-- Particle reflection
198	particles(n)%z = 2.0_wp * zw(0) - particles(n)%z
199	particles(n)%speed_z = -particles(n)%speed_z
200	IF ( use_sgs_for_particles .AND. &
201	particles(n)%rvar3 < 0.0_wp ) THEN
202	particles(n)%rvar3 = -particles(n)%rvar3
203	ENDIF
204	ENDIF
205	ENDIF
206	ENDDO
207
208	ELSEIF ( range == 'walls' ) THEN
209
210	CALL cpu_log( log_point_s(48), 'lpm_wall_reflect', 'start' )
211
212	reflect_x = .FALSE.
213	reflect_y = .FALSE.
214	reflect_z = .FALSE.
215
216	DO n = 1, number_of_particles
217
218	dt_particle = particles(n)%age - particles(n)%age_m
219
220	i2 = ( particles(n)%x + 0.5_wp * dx ) * ddx
221	j2 = ( particles(n)%y + 0.5_wp * dy ) * ddy
222	k2 = particles(n)%z / dz + 1 + offset_ocean_nzt_m1
223
224	prt_x = particles(n)%x
225	prt_y = particles(n)%y
226	prt_z = particles(n)%z
227
228	!
229	!-- If the particle position is below the surface, it has to be reflected
230	IF ( k2 <= nzb_s_inner(j2,i2) .AND. nzb_s_inner(j2,i2) /=0 ) THEN
231
232	pos_x_old = particles(n)%x - particles(n)%speed_x * dt_particle
233	pos_y_old = particles(n)%y - particles(n)%speed_y * dt_particle
234	pos_z_old = particles(n)%z - particles(n)%speed_z * dt_particle
235	i1 = ( pos_x_old + 0.5_wp * dx ) * ddx
236	j1 = ( pos_y_old + 0.5_wp * dy ) * ddy
237	k1 = pos_z_old / dz + offset_ocean_nzt_m1
238
239	!
240	!-- Case 1
241	IF ( particles(n)%x > pos_x_old .AND. particles(n)%y > pos_y_old )&
242	THEN
243	t_index = 1
244
245	DO i = i1, i2
246	xline = i * dx + 0.5_wp * dx
247	t(t_index) = ( xline - pos_x_old ) / &
248	( particles(n)%x - pos_x_old )
249	t_index = t_index + 1
250	ENDDO
251
252	DO j = j1, j2
253	yline = j * dy + 0.5_wp * dy
254	t(t_index) = ( yline - pos_y_old ) / &
255	( particles(n)%y - pos_y_old )
256	t_index = t_index + 1
257	ENDDO
258
259	IF ( particles(n)%z < pos_z_old ) THEN
260	DO k = k1, k2 , -1
261	zline = k * dz
262	t(t_index) = ( pos_z_old - zline ) / &
263	( pos_z_old - particles(n)%z )
264	t_index = t_index + 1
265	ENDDO
266	ENDIF
267
268	t_index_number = t_index - 1
269
270	!
271	!-- Sorting t
272	inc = 1
273	jr = 1
274	DO WHILE ( inc <= t_index_number )
275	inc = 3 * inc + 1
276	ENDDO
277
278	DO WHILE ( inc > 1 )
279	inc = inc / 3
280	DO ir = inc+1, t_index_number
281	tmp_t = t(ir)
282	jr = ir
283	DO WHILE ( t(jr-inc) > tmp_t )
284	t(jr) = t(jr-inc)
285	jr = jr - inc
286	IF ( jr <= inc ) EXIT
287	ENDDO
288	t(jr) = tmp_t
289	ENDDO
290	ENDDO
291
292	case1: DO t_index = 1, t_index_number
293
294	pos_x = pos_x_old + t(t_index) * ( prt_x - pos_x_old )
295	pos_y = pos_y_old + t(t_index) * ( prt_y - pos_y_old )
296	pos_z = pos_z_old + t(t_index) * ( prt_z - pos_z_old )
297
298	i3 = ( pos_x + 0.5_wp * dx ) * ddx
299	j3 = ( pos_y + 0.5_wp * dy ) * ddy
300	k3 = pos_z / dz + offset_ocean_nzt_m1
301
302	i5 = pos_x * ddx
303	j5 = pos_y * ddy
304	k5 = pos_z / dz + offset_ocean_nzt_m1
305
306	IF ( k5 <= nzb_s_inner(j5,i3) .AND. &
307	nzb_s_inner(j5,i3) /= 0 ) THEN
308
309	IF ( pos_z == nzb_s_inner(j5,i3) * dz .AND. &
310	k3 == nzb_s_inner(j5,i3) ) THEN
311	reflect_z = .TRUE.
312	EXIT case1
313	ENDIF
314
315	ENDIF
316
317	IF ( k5 <= nzb_s_inner(j3,i5) .AND. &
318	nzb_s_inner(j3,i5) /= 0 ) THEN
319
320	IF ( pos_z == nzb_s_inner(j3,i5) * dz .AND. &
321	k3 == nzb_s_inner(j3,i5) ) THEN
322	reflect_z = .TRUE.
323	EXIT case1
324	ENDIF
325
326	ENDIF
327
328	IF ( k3 <= nzb_s_inner(j3,i3) .AND. &
329	nzb_s_inner(j3,i3) /= 0 ) THEN
330
331	IF ( pos_z == nzb_s_inner(j3,i3) * dz .AND. &
332	k3 == nzb_s_inner(j3,i3) ) THEN
333	reflect_z = .TRUE.
334	EXIT case1
335	ENDIF
336
337	IF ( pos_y == ( j3 * dy - 0.5_wp * dy ) .AND. &
338	pos_z < nzb_s_inner(j3,i3) * dz ) THEN
339	reflect_y = .TRUE.
340	EXIT case1
341	ENDIF
342
343	IF ( pos_x == ( i3 * dx - 0.5_wp * dx ) .AND. &
344	pos_z < nzb_s_inner(j3,i3) * dz ) THEN
345	reflect_x = .TRUE.
346	EXIT case1
347	ENDIF
348
349	ENDIF
350
351	ENDDO case1
352
353	!
354	!-- Case 2
355	ELSEIF ( particles(n)%x > pos_x_old .AND. &
356	particles(n)%y < pos_y_old ) THEN
357
358	t_index = 1
359
360	DO i = i1, i2
361	xline = i * dx + 0.5_wp * dx
362	t(t_index) = ( xline - pos_x_old ) / &
363	( particles(n)%x - pos_x_old )
364	t_index = t_index + 1
365	ENDDO
366
367	DO j = j1, j2, -1
368	yline = j * dy - 0.5_wp * dy
369	t(t_index) = ( pos_y_old - yline ) / &
370	( pos_y_old - particles(n)%y )
371	t_index = t_index + 1
372	ENDDO
373
374	IF ( particles(n)%z < pos_z_old ) THEN
375	DO k = k1, k2 , -1
376	zline = k * dz
377	t(t_index) = ( pos_z_old - zline ) / &
378	( pos_z_old - particles(n)%z )
379	t_index = t_index + 1
380	ENDDO
381	ENDIF
382	t_index_number = t_index-1
383
384	!
385	!-- Sorting t
386	inc = 1
387	jr = 1
388	DO WHILE ( inc <= t_index_number )
389	inc = 3 * inc + 1
390	ENDDO
391
392	DO WHILE ( inc > 1 )
393	inc = inc / 3
394	DO ir = inc+1, t_index_number
395	tmp_t = t(ir)
396	jr = ir
397	DO WHILE ( t(jr-inc) > tmp_t )
398	t(jr) = t(jr-inc)
399	jr = jr - inc
400	IF ( jr <= inc ) EXIT
401	ENDDO
402	t(jr) = tmp_t
403	ENDDO
404	ENDDO
405
406	case2: DO t_index = 1, t_index_number
407
408	pos_x = pos_x_old + t(t_index) * ( prt_x - pos_x_old )
409	pos_y = pos_y_old + t(t_index) * ( prt_y - pos_y_old )
410	pos_z = pos_z_old + t(t_index) * ( prt_z - pos_z_old )
411
412	i3 = ( pos_x + 0.5_wp * dx ) * ddx
413	j3 = ( pos_y + 0.5_wp * dy ) * ddy
414	k3 = pos_z / dz + offset_ocean_nzt_m1
415
416	i5 = pos_x * ddx
417	j5 = pos_y * ddy
418	k5 = pos_z / dz + offset_ocean_nzt_m1
419
420	IF ( k5 <= nzb_s_inner(j3,i5) .AND. &
421	nzb_s_inner(j3,i5) /= 0 ) THEN
422
423	IF ( pos_z == nzb_s_inner(j3,i5) * dz .AND. &
424	k3 == nzb_s_inner(j3,i5) ) THEN
425	reflect_z = .TRUE.
426	EXIT case2
427	ENDIF
428
429	ENDIF
430
431	IF ( k3 <= nzb_s_inner(j3,i3) .AND. &
432	nzb_s_inner(j3,i3) /= 0 ) THEN
433
434	IF ( pos_z == nzb_s_inner(j3,i3) * dz .AND. &
435	k3 == nzb_s_inner(j3,i3) ) THEN
436	reflect_z = .TRUE.
437	EXIT case2
438	ENDIF
439
440	IF ( pos_x == ( i3 * dx - 0.5_wp * dx ) .AND. &
441	pos_z < nzb_s_inner(j3,i3) * dz ) THEN
442	reflect_x = .TRUE.
443	EXIT case2
444	ENDIF
445
446	ENDIF
447
448	IF ( k5 <= nzb_s_inner(j5,i3) .AND. &
449	nzb_s_inner(j5,i3) /= 0 ) THEN
450
451	IF ( pos_z == nzb_s_inner(j5,i3) * dz .AND. &
452	k3 == nzb_s_inner(j5,i3) ) THEN
453	reflect_z = .TRUE.
454	EXIT case2
455	ENDIF
456
457	IF ( pos_y == ( j5 * dy + 0.5_wp * dy ) .AND. &
458	pos_z < nzb_s_inner(j5,i3) * dz ) THEN
459	reflect_y = .TRUE.
460	EXIT case2
461	ENDIF
462
463	ENDIF
464
465	ENDDO case2
466
467	!
468	!-- Case 3
469	ELSEIF ( particles(n)%x < pos_x_old .AND. &
470	particles(n)%y > pos_y_old ) THEN
471
472	t_index = 1
473
474	DO i = i1, i2, -1
475	xline = i * dx - 0.5_wp * dx
476	t(t_index) = ( pos_x_old - xline ) / &
477	( pos_x_old - particles(n)%x )
478	t_index = t_index + 1
479	ENDDO
480
481	DO j = j1, j2
482	yline = j * dy + 0.5_wp * dy
483	t(t_index) = ( yline - pos_y_old ) / &
484	( particles(n)%y - pos_y_old )
485	t_index = t_index + 1
486	ENDDO
487
488	IF ( particles(n)%z < pos_z_old ) THEN
489	DO k = k1, k2 , -1
490	zline = k * dz
491	t(t_index) = ( pos_z_old - zline ) / &
492	( pos_z_old - particles(n)%z )
493	t_index = t_index + 1
494	ENDDO
495	ENDIF
496	t_index_number = t_index - 1
497
498	!
499	!-- Sorting t
500	inc = 1
501	jr = 1
502
503	DO WHILE ( inc <= t_index_number )
504	inc = 3 * inc + 1
505	ENDDO
506
507	DO WHILE ( inc > 1 )
508	inc = inc / 3
509	DO ir = inc+1, t_index_number
510	tmp_t = t(ir)
511	jr = ir
512	DO WHILE ( t(jr-inc) > tmp_t )
513	t(jr) = t(jr-inc)
514	jr = jr - inc
515	IF ( jr <= inc ) EXIT
516	ENDDO
517	t(jr) = tmp_t
518	ENDDO
519	ENDDO
520
521	case3: DO t_index = 1, t_index_number
522
523	pos_x = pos_x_old + t(t_index) * ( prt_x - pos_x_old )
524	pos_y = pos_y_old + t(t_index) * ( prt_y - pos_y_old )
525	pos_z = pos_z_old + t(t_index) * ( prt_z - pos_z_old )
526
527	i3 = ( pos_x + 0.5_wp * dx ) * ddx
528	j3 = ( pos_y + 0.5_wp * dy ) * ddy
529	k3 = pos_z / dz + offset_ocean_nzt_m1
530
531	i5 = pos_x * ddx
532	j5 = pos_y * ddy
533	k5 = pos_z / dz + offset_ocean_nzt_m1
534
535	IF ( k5 <= nzb_s_inner(j5,i3) .AND. &
536	nzb_s_inner(j5,i3) /= 0 ) THEN
537
538	IF ( pos_z == nzb_s_inner(j5,i3) * dz .AND. &
539	k3 == nzb_s_inner(j5,i3) ) THEN
540	reflect_z = .TRUE.
541	EXIT case3
542	ENDIF
543
544	ENDIF
545
546	IF ( k3 <= nzb_s_inner(j3,i3) .AND. &
547	nzb_s_inner(j3,i3) /= 0 ) THEN
548
549	IF ( pos_z == nzb_s_inner(j3,i3) * dz .AND. &
550	k3 == nzb_s_inner(j3,i3) ) THEN
551	reflect_z = .TRUE.
552	EXIT case3
553	ENDIF
554
555	IF ( pos_y == ( j3 * dy - 0.5_wp * dy ) .AND. &
556	pos_z < nzb_s_inner(j3,i3) * dz ) THEN
557	reflect_y = .TRUE.
558	EXIT case3
559	ENDIF
560
561	ENDIF
562
563	IF ( k5 <= nzb_s_inner(j3,i5) .AND. &
564	nzb_s_inner(j3,i5) /= 0 ) THEN
565
566	IF ( pos_z == nzb_s_inner(j3,i5) * dz .AND. &
567	k3 == nzb_s_inner(j3,i5) ) THEN
568	reflect_z = .TRUE.
569	EXIT case3
570	ENDIF
571
572	IF ( pos_x == ( i5 * dx + 0.5_wp * dx ) .AND. &
573	pos_z < nzb_s_inner(j3,i5) * dz ) THEN
574	reflect_x = .TRUE.
575	EXIT case3
576	ENDIF
577
578	ENDIF
579
580	ENDDO case3
581
582	!
583	!-- Case 4
584	ELSEIF ( particles(n)%x < pos_x_old .AND. &
585	particles(n)%y < pos_y_old ) THEN
586
587	t_index = 1
588
589	DO i = i1, i2, -1
590	xline = i * dx - 0.5_wp * dx
591	t(t_index) = ( pos_x_old - xline ) / &
592	( pos_x_old - particles(n)%x )
593	t_index = t_index + 1
594	ENDDO
595
596	DO j = j1, j2, -1
597	yline = j * dy - 0.5_wp * dy
598	t(t_index) = ( pos_y_old - yline ) / &
599	( pos_y_old - particles(n)%y )
600	t_index = t_index + 1
601	ENDDO
602
603	IF ( particles(n)%z < pos_z_old ) THEN
604	DO k = k1, k2 , -1
605	zline = k * dz
606	t(t_index) = ( pos_z_old - zline ) / &
607	( pos_z_old-particles(n)%z )
608	t_index = t_index + 1
609	ENDDO
610	ENDIF
611	t_index_number = t_index-1
612
613	!
614	!-- Sorting t
615	inc = 1
616	jr = 1
617
618	DO WHILE ( inc <= t_index_number )
619	inc = 3 * inc + 1
620	ENDDO
621
622	DO WHILE ( inc > 1 )
623	inc = inc / 3
624	DO ir = inc+1, t_index_number
625	tmp_t = t(ir)
626	jr = ir
627	DO WHILE ( t(jr-inc) > tmp_t )
628	t(jr) = t(jr-inc)
629	jr = jr - inc
630	IF ( jr <= inc ) EXIT
631	ENDDO
632	t(jr) = tmp_t
633	ENDDO
634	ENDDO
635
636	case4: DO t_index = 1, t_index_number
637
638	pos_x = pos_x_old + t(t_index) * ( prt_x - pos_x_old )
639	pos_y = pos_y_old + t(t_index) * ( prt_y - pos_y_old )
640	pos_z = pos_z_old + t(t_index) * ( prt_z - pos_z_old )
641
642	i3 = ( pos_x + 0.5_wp * dx ) * ddx
643	j3 = ( pos_y + 0.5_wp * dy ) * ddy
644	k3 = pos_z / dz + offset_ocean_nzt_m1
645
646	i5 = pos_x * ddx
647	j5 = pos_y * ddy
648	k5 = pos_z / dz + offset_ocean_nzt_m1
649
650	IF ( k3 <= nzb_s_inner(j3,i3) .AND. &
651	nzb_s_inner(j3,i3) /= 0 ) THEN
652
653	IF ( pos_z == nzb_s_inner(j3,i3) * dz .AND. &
654	k3 == nzb_s_inner(j3,i3) ) THEN
655	reflect_z = .TRUE.
656	EXIT case4
657	ENDIF
658
659	ENDIF
660
661	IF ( k5 <= nzb_s_inner(j3,i5) .AND. &
662	nzb_s_inner(j3,i5) /= 0 ) THEN
663
664	IF ( pos_z == nzb_s_inner(j3,i5) * dz .AND. &
665	k3 == nzb_s_inner(j3,i5) ) THEN
666	reflect_z = .TRUE.
667	EXIT case4
668	ENDIF
669
670	IF ( pos_x == ( i5 * dx + 0.5_wp * dx ) .AND. &
671	nzb_s_inner(j3,i5) /=0 .AND. &
672	pos_z < nzb_s_inner(j3,i5) * dz ) THEN
673	reflect_x = .TRUE.
674	EXIT case4
675	ENDIF
676
677	ENDIF
678
679	IF ( k5 <= nzb_s_inner(j5,i3) .AND. &
680	nzb_s_inner(j5,i3) /= 0 ) THEN
681
682	IF ( pos_z == nzb_s_inner(j5,i3) * dz .AND. &
683	k5 == nzb_s_inner(j5,i3) ) THEN
684	reflect_z = .TRUE.
685	EXIT case4
686	ENDIF
687
688	IF ( pos_y == ( j5 * dy + 0.5_wp * dy ) .AND. &
689	nzb_s_inner(j5,i3) /= 0 .AND. &
690	pos_z < nzb_s_inner(j5,i3) * dz ) THEN
691	reflect_y = .TRUE.
692	EXIT case4
693	ENDIF
694
695	ENDIF
696
697	ENDDO case4
698
699	ENDIF
700
701	ENDIF ! Check, if particle position is below the surface
702
703	!
704	!-- Do the respective reflection, in case that one of the above conditions
705	!-- is found to be true
706	IF ( reflect_z ) THEN
707
708	particles(n)%z = 2.0_wp * pos_z - prt_z
709	particles(n)%speed_z = - particles(n)%speed_z
710
711	IF ( use_sgs_for_particles ) THEN
712	particles(n)%rvar3 = - particles(n)%rvar3
713	ENDIF
714	reflect_z = .FALSE.
715
716	ELSEIF ( reflect_y ) THEN
717
718	particles(n)%y = 2.0_wp * pos_y - prt_y
719	particles(n)%speed_y = - particles(n)%speed_y
720
721	IF ( use_sgs_for_particles ) THEN
722	particles(n)%rvar2 = - particles(n)%rvar2
723	ENDIF
724	reflect_y = .FALSE.
725
726	ELSEIF ( reflect_x ) THEN
727
728	particles(n)%x = 2.0_wp * pos_x - prt_x
729	particles(n)%speed_x = - particles(n)%speed_x
730
731	IF ( use_sgs_for_particles ) THEN
732	particles(n)%rvar1 = - particles(n)%rvar1
733	ENDIF
734	reflect_x = .FALSE.
735
736	ENDIF
737
738	ENDDO
739
740	CALL cpu_log( log_point_s(48), 'lpm_wall_reflect', 'stop' )
741
742	ENDIF
743
744	END SUBROUTINE lpm_boundary_conds

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |