Home

Context Navigation

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

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

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |