Home

Context Navigation

source: palm/trunk/SOURCE/particle_boundary_conds.f90 @ 59

Last change on this file since 59 was 58, checked in by raasch, 18 years ago
preliminary update of further changes, not running
Property svn:keywords set to `Id`
File size: 21.0 KB

Rev	Line
[58]	1	SUBROUTINE particle_boundary_conds
	2
	3	!------------------------------------------------------------------------------!
	4	! Actual revisions:
	5	! -----------------
	6	!
	7	!
	8	! Former revisions:
	9	! -----------------
	10	! $Id: particle_boundary_conds.f90 58 2007-03-09 14:27:38Z raasch $
	11	! Initial version (2007/03/09)
	12	!
	13	! Description:
	14	! ------------
	15	! Calculates the reflection of particles from vertical walls.
	16	! Routine developed by Jin Zhang (2006-2007)
	17	!------------------------------------------------------------------------------!
	18	#if defined( __particles )
	19
	20	USE control_parameters
	21	USE cpulog
	22	USE grid_variables
	23	USE indices
	24	USE interfaces
	25	USE particle_attributes
	26	USE pegrid
	27
	28	IMPLICIT NONE
	29
	30	INTEGER :: n
	31
	32
	33	CALL cpu_log( log_point_s(48), 'advec_part_refle', 'start' )
	34
	35
	36
	37	DO n = 1, number_of_particles
	38
	39	i2 = ( particles(n)%x + 0.5 * dx ) * ddx
	40	j2 = ( particles(n)%y + 0.5 * dy ) * ddy
	41	k2 = particles(n)%z / dz + 1
	42
	43	prt_x = particles(n)%x
	44	prt_y = particles(n)%y
	45	prt_z = particles(n)%z
	46
	47
	48	IF ( k2 <= nzb_s_inner(j2,i2) .AND. nzb_s_inner(j2,i2) /=0 ) THEN
	49
	50	pos_x_old = particles(n)%x - particles(n)%speed_x * dt_particle
	51	pos_y_old = particles(n)%y - particles(n)%speed_y * dt_particle
	52	pos_z_old = particles(n)%z - particles(n)%speed_z * dt_particle
	53	i1 = ( pos_x_old + 0.5 * dx ) * ddx
	54	j1 = ( pos_y_old + 0.5 * dy ) * ddy
	55	k1 = pos_z_old / dz
	56
	57	!
	58	!-- Case 1
	59	IF ( particles(n)%x > pos_x_old .AND. particles(n)%y > pos_y_old ) &
	60	THEN
	61	t_index = 1
	62
	63	DO i = i1, i2
	64	xline = i * dx + 0.5 * dx
	65	t(t_index) = ( xline - pos_x_old ) / &
	66	( particles(n)%x - pos_x_old )
	67	t_index = t_index + 1
	68	ENDDO
	69
	70	DO j = j1, j2
	71	yline = j * dy + 0.5 * dy
	72	t(t_index) = ( yline - pos_y_old ) / &
	73	( particles(n)%y - pos_y_old )
	74	t_index = t_index + 1
	75	ENDDO
	76
	77	IF ( particles(n)%z < pos_z_old ) THEN
	78	DO k = k1, k2 , -1
	79	zline = k * dz
	80	t(t_index) = ( pos_z_old - zline ) / &
	81	( pos_z_old - particles(n)%z )
	82	t_index = t_index + 1
	83	ENDDO
	84	ENDIF
	85
	86	t_index_number = t_index - 1
	87
	88	!
	89	!-- Sorting t
	90	inc = 1
	91	jr = 1
	92	DO WHILE ( inc <= t_index_number )
	93	inc = 3 * inc + 1
	94	ENDDO
	95
	96	DO WHILE ( inc > 1 )
	97	inc = inc / 3
	98	DO ir = inc+1, t_index_number
	99	tmp_t = t(ir)
	100	jr = ir
	101	DO WHILE ( t(jr-inc) > tmp_t )
	102	t(jr) = t(jr-inc)
	103	jr = jr - inc
	104	IF ( jr <= inc ) EXIT
	105	ENDDO
	106	t(jr) = tmp_t
	107	ENDDO
	108	ENDDO
	109
	110	DO t_index = 1, t_index_number
	111
	112	pos_x = pos_x_old + t(t_index) * ( prt_x - pos_x_old )
	113	pos_y = pos_y_old + t(t_index) * ( prt_y - pos_y_old )
	114	pos_z = pos_z_old + t(t_index) * ( prt_z - pos_z_old )
	115
	116	i3 = ( pos_x + 0.5 * dx ) * ddx
	117	j3 = ( pos_y + 0.5 * dy ) * ddy
	118	k3 = pos_z / dz
	119
	120	i5 = pos_x * ddx
	121	j5 = pos_y * ddy
	122	k5 = pos_z / dz
	123
	124	IF ( k5 <= nzb_s_inner(j5,i3) .AND. nzb_s_inner(j5,i3) /= 0 )&
	125	THEN
	126	IF ( pos_z == nzb_s_inner(j5,i3) * dz .AND. &
	127	k3 == nzb_s_inner(j5,i3) ) THEN
	128	particles(n)%z = 2.0 * pos_z - prt_z
	129	particles(n)%speed_z = - particles(n)%speed_z
	130
	131	IF ( use_sgs_for_particles ) THEN
	132	particles(n)%speed_z_sgs = - particles(n)%speed_z_sgs
	133	ENDIF
	134	GOTO 999
	135	ENDIF
	136	ENDIF
	137
	138	IF ( k5 <= nzb_s_inner(j3,i5) .AND. nzb_s_inner(j3,i5) /= 0 )&
	139	THEN
	140	IF ( pos_z == nzb_s_inner(j3,i5) * dz .AND. &
	141	k3 == nzb_s_inner(j3,i5) ) THEN
	142	particles(n)%z = 2.0 * pos_z - prt_z
	143	particles(n)%speed_z = - particles(n)%speed_z
	144
	145	IF ( use_sgs_for_particles ) THEN
	146	particles(n)%speed_z_sgs = - particles(n)%speed_z_sgs
	147	ENDIF
	148	GOTO 999
	149	ENDIF
	150	ENDIF
	151
	152	IF ( k3 <= nzb_s_inner(j3,i3) .AND. nzb_s_inner(j3,i3) /= 0 )&
	153	THEN
	154	IF ( pos_z == nzb_s_inner(j3,i3) * dz .AND. &
	155	k3 == nzb_s_inner(j3,i3) ) THEN
	156	particles(n)%z = 2.0 * pos_z - prt_z
	157	particles(n)%speed_z = - particles(n)%speed_z
	158
	159	IF ( use_sgs_for_particles ) THEN
	160	particles(n)%speed_z_sgs = - particles(n)%speed_z_sgs
	161	ENDIF
	162	GOTO 999
	163	ENDIF
	164
	165	IF ( pos_y == ( j3 * dy - 0.5 * dy ) .AND. &
	166	pos_z < nzb_s_inner(j3,i3) * dz ) THEN
	167	particles(n)%y = 2.0 * pos_y - prt_y
	168	particles(n)%speed_y = - particles(n)%speed_y
	169
	170	IF ( use_sgs_for_particles ) THEN
	171	particles(n)%speed_y_sgs = - particles(n)%speed_y_sgs
	172	ENDIF
	173	GOTO 999
	174	ENDIF
	175
	176	IF ( pos_x == ( i3 * dx - 0.5 * dx ) .AND. &
	177	pos_z < nzb_s_inner(j3,i3) * dz ) THEN
	178	particles(n)%x = 2.0 * pos_x - prt_x
	179	particles(n)%speed_x = - particles(n)%speed_x
	180
	181	IF ( use_sgs_for_particles ) THEN
	182	particles(n)%speed_x_sgs = - particles(n)%speed_x_sgs
	183	ENDIF
	184	GOTO 999
	185	ENDIF
	186	ENDIF
	187
	188	ENDDO
	189
	190	ENDIF
	191
	192	!
	193	!-- Case 2
	194	IF ( particles(n)%x > pos_x_old .AND. particles(n)%y < pos_y_old ) &
	195	THEN
	196	t_index = 1
	197
	198	DO i = i1, i2
	199	xline = i * dx + 0.5 * dx
	200	t(t_index) = ( xline - pos_x_old ) / &
	201	( particles(n)%x - pos_x_old )
	202	t_index = t_index + 1
	203	ENDDO
	204
	205	DO j = j1, j2, -1
	206	yline = j * dy - 0.5 * dy
	207	t(t_index) = ( pos_y_old - yline ) / &
	208	( pos_y_old - particles(n)%y )
	209	t_index = t_index + 1
	210	ENDDO
	211
	212	IF ( particles(n)%z < pos_z_old ) THEN
	213	DO k = k1, k2 , -1
	214	zline = k * dz
	215	t(t_index) = ( pos_z_old - zline ) / &
	216	( pos_z_old - particles(n)%z )
	217	t_index = t_index + 1
	218	ENDDO
	219	ENDIF
	220	t_index_number = t_index-1
	221
	222	!
	223	!-- Sorting t
	224	inc = 1
	225	jr = 1
	226	DO WHILE ( inc <= t_index_number )
	227	inc = 3 * inc + 1
	228	ENDDO
	229
	230	DO WHILE ( inc > 1 )
	231	inc = inc / 3
	232	DO ir = inc+1, t_index_number
	233	tmp_t = t(ir)
	234	jr = ir
	235	DO WHILE ( t(jr-inc) > tmp_t )
	236	t(jr) = t(jr-inc)
	237	jr = jr - inc
	238	IF ( jr <= inc ) EXIT
	239	ENDDO
	240	t(jr) = tmp_t
	241	ENDDO
	242	ENDDO
	243
	244	DO t_index = 1, t_index_number
	245
	246	pos_x = pos_x_old + t(t_index) * ( prt_x - pos_x_old )
	247	pos_y = pos_y_old + t(t_index) * ( prt_y - pos_y_old )
	248	pos_z = pos_z_old + t(t_index) * ( prt_z - pos_z_old )
	249
	250	i3 = ( pos_x + 0.5 * dx ) * ddx
	251	j3 = ( pos_y + 0.5 * dy ) * ddy
	252	k3 = pos_z / dz
	253
	254	i5 = pos_x * ddx
	255	j5 = pos_y * ddy
	256	k5 = pos_z / dz
	257
	258	IF ( k5 <= nzb_s_inner(j3,i5) .AND. nzb_s_inner(j3,i5) /= 0 )&
	259	THEN
	260	IF ( pos_z == nzb_s_inner(j3,i5) * dz .AND. &
	261	k3 == nzb_s_inner(j3,i5) ) THEN
	262	particles(n)%z = 2.0 * pos_z - prt_z
	263	particles(n)%speed_z = - particles(n)%speed_z
	264
	265	IF ( use_sgs_for_particles ) THEN
	266	particles(n)%speed_z_sgs = - particles(n)%speed_z_sgs
	267	ENDIF
	268	GOTO 999
	269	ENDIF
	270	ENDIF
	271
	272	IF ( k3 <= nzb_s_inner(j3,i3) .AND. nzb_s_inner(j3,i3) /= 0 )&
	273	THEN
	274	IF ( pos_z == nzb_s_inner(j3,i3) * dz .AND. &
	275	k3 == nzb_s_inner(j3,i3) ) THEN
	276	particles(n)%z = 2.0 * pos_z - prt_z
	277	particles(n)%speed_z = - particles(n)%speed_z
	278
	279	IF ( use_sgs_for_particles ) THEN
	280	particles(n)%speed_z_sgs = - particles(n)%speed_z_sgs
	281	ENDIF
	282	GOTO 999
	283	ENDIF
	284
	285	IF ( pos_x == ( i3 * dx - 0.5 * dx ) .AND. &
	286	pos_z < nzb_s_inner(j3,i3) * dz ) THEN
	287	particles(n)%x = 2.0 * pos_x - prt_x
	288	particles(n)%speed_x = - particles(n)%speed_x
	289
	290	IF ( use_sgs_for_particles ) THEN
	291	particles(n)%speed_x_sgs = - particles(n)%speed_x_sgs
	292	ENDIF
	293	GOTO 999
	294	ENDIF
	295	ENDIF
	296
	297	IF ( k5 <= nzb_s_inner(j5,i3) .AND. nzb_s_inner(j5,i3) /= 0 )&
	298	THEN
	299	IF ( pos_z == nzb_s_inner(j5,i3) * dz .AND. &
	300	k3 == nzb_s_inner(j5,i3) ) THEN
	301	particles(n)%z = 2.0 * pos_z - prt_z
	302	particles(n)%speed_z = - particles(n)%speed_z
	303	IF ( use_sgs_for_particles ) THEN
	304	particles(n)%speed_z_sgs = - particles(n)%speed_z_sgs
	305	ENDIF
	306	GOTO 999
	307	ENDIF
	308
	309	IF ( pos_y == ( j5 * dy + 0.5 * dy ) .AND. &
	310	pos_z < nzb_s_inner(j5,i3) * dz ) THEN
	311	particles(n)%y = 2.0 * pos_y - prt_y
	312	particles(n)%speed_y = - particles(n)%speed_y
	313
	314	IF ( use_sgs_for_particles ) THEN
	315	particles(n)%speed_y_sgs = - particles(n)%speed_y_sgs
	316	ENDIF
	317	GOTO 999
	318	ENDIF
	319	ENDIF
	320	ENDDO
	321	ENDIF
	322
	323	!
	324	!-- Case 3
	325	IF ( particles(n)%x < pos_x_old .AND. particles(n)%y > pos_y_old ) &
	326	THEN
	327	t_index = 1
	328
	329	DO i = i1, i2, -1
	330	xline = i * dx - 0.5 * dx
	331	t(t_index) = ( pos_x_old - xline ) / &
	332	( pos_x_old - particles(n)%x )
	333	t_index = t_index + 1
	334	ENDDO
	335
	336	DO j = j1, j2
	337	yline = j * dy + 0.5 * dy
	338	t(t_index) = ( yline - pos_y_old ) / &
	339	( particles(n)%y - pos_y_old )
	340	t_index = t_index + 1
	341	ENDDO
	342
	343	IF ( particles(n)%z < pos_z_old ) THEN
	344	DO k = k1, k2 , -1
	345	zline = k * dz
	346	t(t_index) = ( pos_z_old - zline ) / &
	347	( pos_z_old - particles(n)%z )
	348	t_index = t_index + 1
	349	ENDDO
	350	ENDIF
	351	t_index_number = t_index - 1
	352
	353	!
	354	!-- Sorting t
	355	inc = 1
	356	jr = 1
	357
	358	DO WHILE ( inc <= t_index_number )
	359	inc = 3 * inc + 1
	360	ENDDO
	361
	362	DO WHILE ( inc > 1 )
	363	inc = inc / 3
	364	DO ir = inc+1, t_index_number
	365	tmp_t = t(ir)
	366	jr = ir
	367	DO WHILE ( t(jr-inc) > tmp_t )
	368	t(jr) = t(jr-inc)
	369	jr = jr - inc
	370	IF ( jr <= inc ) EXIT
	371	ENDDO
	372	t(jr) = tmp_t
	373	ENDDO
	374	ENDDO
	375
	376	DO t_index = 1, t_index_number
	377
	378	pos_x = pos_x_old + t(t_index) * ( prt_x - pos_x_old )
	379	pos_y = pos_y_old + t(t_index) * ( prt_y - pos_y_old )
	380	pos_z = pos_z_old + t(t_index) * ( prt_z - pos_z_old )
	381
	382	i3 = ( pos_x + 0.5 * dx ) * ddx
	383	j3 = ( pos_y + 0.5 * dy ) * ddy
	384	k3 = pos_z / dz
	385
	386	i5 = pos_x * ddx
	387	j5 = pos_y * ddy
	388	k5 = pos_z / dz
	389
	390	IF ( k5 <= nzb_s_inner(j5,i3) .AND. nzb_s_inner(j5,i3) /= 0 )&
	391	THEN
	392	IF ( pos_z == nzb_s_inner(j5,i3) * dz .AND. &
	393	k3 == nzb_s_inner(j5,i3) ) THEN
	394	particles(n)%z = 2.0 * pos_z - prt_z
	395	particles(n)%speed_z = - particles(n)%speed_z
	396
	397	IF ( use_sgs_for_particles ) THEN
	398	particles(n)%speed_z_sgs = - particles(n)%speed_z_sgs
	399	ENDIF
	400	GOTO 999
	401	ENDIF
	402	ENDIF
	403
	404	IF ( k3 <= nzb_s_inner(j3,i3) .AND. nzb_s_inner(j3,i3) /= 0 )&
	405	THEN
	406	IF ( pos_z == nzb_s_inner(j3,i3) * dz .AND. &
	407	k3 == nzb_s_inner(j3,i3) ) THEN
	408	particles(n)%z = 2.0 * pos_z - prt_z
	409	particles(n)%speed_z = - particles(n)%speed_z
	410
	411	IF ( use_sgs_for_particles ) THEN
	412	particles(n)%speed_z_sgs = - particles(n)%speed_z_sgs
	413	ENDIF
	414	GOTO 999
	415	ENDIF
	416
	417	IF ( pos_y == ( j3 * dy - 0.5 * dy ) .AND. &
	418	pos_z < nzb_s_inner(j3,i3) * dz ) THEN
	419	particles(n)%y = 2.0 * pos_y - prt_y
	420	particles(n)%speed_y = - particles(n)%speed_y
	421
	422	IF ( use_sgs_for_particles ) THEN
	423	particles(n)%speed_y_sgs = - particles(n)%speed_y_sgs
	424	ENDIF
	425	GOTO 999
	426	ENDIF
	427	ENDIF
	428
	429	IF ( k5 <= nzb_s_inner(j3,i5) .AND. nzb_s_inner(j3,i5) /= 0 )&
	430	THEN
	431	IF ( pos_z == nzb_s_inner(j3,i5) * dz .AND. &
	432	k3 == nzb_s_inner(j3,i5) ) THEN
	433	particles(n)%z = 2.0 * pos_z - prt_z
	434	particles(n)%speed_z = - particles(n)%speed_z
	435
	436	IF ( use_sgs_for_particles ) THEN
	437	particles(n)%speed_z_sgs = - particles(n)%speed_z_sgs
	438	ENDIF
	439	GOTO 999
	440	ENDIF
	441
	442	IF ( pos_x == ( i5 * dx + 0.5 * dx ) .AND. &
	443	pos_z < nzb_s_inner(j3,i5) * dz ) THEN
	444	particles(n)%x = 2.0 * pos_x - prt_x
	445	particles(n)%speed_x = - particles(n)%speed_x
	446
	447	IF ( use_sgs_for_particles ) THEN
	448	particles(n)%speed_x_sgs = - particles(n)%speed_x_sgs
	449	ENDIF
	450	GOTO 999
	451	ENDIF
	452	ENDIF
	453	ENDDO
	454	ENDIF
	455
	456	!
	457	!-- Case 4
	458	IF ( particles(n)%x < pos_x_old .AND. particles(n)%y < pos_y_old )&
	459	THEN
	460	t_index = 1
	461
	462	DO i = i1, i2, -1
	463	xline = i * dx - 0.5 * dx
	464	t(t_index) = ( pos_x_old - xline ) / &
	465	( pos_x_old - particles(n)%x )
	466	t_index = t_index + 1
	467	ENDDO
	468
	469	DO j = j1, j2, -1
	470	yline = j * dy - 0.5 * dy
	471	t(t_index) = ( pos_y_old - yline ) / &
	472	( pos_y_old - particles(n)%y )
	473	t_index = t_index + 1
	474	ENDDO
	475
	476	IF ( particles(n)%z < pos_z_old ) THEN
	477	DO k = k1, k2 , -1
	478	zline = k * dz
	479	t(t_index) = ( pos_z_old - zline ) / &
	480	( pos_z_old-particles(n)%z )
	481	t_index = t_index + 1
	482	ENDDO
	483	ENDIF
	484	t_index_number = t_index-1
	485
	486	!
	487	!-- Sorting t
	488	inc = 1
	489	jr = 1
	490
	491	DO WHILE ( inc <= t_index_number )
	492	inc = 3 * inc + 1
	493	ENDDO
	494
	495	DO WHILE ( inc > 1 )
	496	inc = inc / 3
	497	DO ir = inc+1, t_index_number
	498	tmp_t = t(ir)
	499	jr = ir
	500	DO WHILE ( t(jr-inc) > tmp_t )
	501	t(jr) = t(jr-inc)
	502	jr = jr - inc
	503	IF ( jr <= inc ) EXIT
	504	ENDDO
	505	t(jr) = tmp_t
	506	ENDDO
	507	ENDDO
	508
	509	DO t_index = 1, t_index_number
	510
	511	pos_x = pos_x_old + t(t_index) * ( prt_x - pos_x_old )
	512	pos_y = pos_y_old + t(t_index) * ( prt_y - pos_y_old )
	513	pos_z = pos_z_old + t(t_index) * ( prt_z - pos_z_old )
	514
	515	i3 = ( pos_x + 0.5 * dx ) * ddx
	516	j3 = ( pos_y + 0.5 * dy ) * ddy
	517	k3 = pos_z / dz
	518
	519	i5 = pos_x * ddx
	520	j5 = pos_y * ddy
	521	k5 = pos_z / dz
	522
	523	IF ( k3 <= nzb_s_inner(j3,i3) .AND. nzb_s_inner(j3,i3) /= 0 )&
	524	THEN
	525	IF ( pos_z == nzb_s_inner(j3,i3) * dz .AND. &
	526	k3 == nzb_s_inner(j3,i3) ) THEN
	527	particles(n)%z = 2.0 * pos_z - prt_z
	528	particles(n)%speed_z = - particles(n)%speed_z
	529
	530	IF ( use_sgs_for_particles ) THEN
	531	particles(n)%speed_z_sgs = - particles(n)%speed_z_sgs
	532	ENDIF
	533	GOTO 999
	534	ENDIF
	535	ENDIF
	536
	537	IF ( k5 <= nzb_s_inner(j3,i5) .AND. nzb_s_inner(j3,i5) /= 0 )&
	538	THEN
	539	IF ( pos_z == nzb_s_inner(j3,i5) * dz .AND. &
	540	k3 == nzb_s_inner(j3,i5) ) THEN
	541	particles(n)%z = 2.0 * pos_z - prt_z
	542	particles(n)%speed_z = - particles(n)%speed_z
	543
	544	IF ( use_sgs_for_particles ) THEN
	545	particles(n)%speed_z_sgs = - particles(n)%speed_z_sgs
	546	ENDIF
	547	GOTO 999
	548	ENDIF
	549
	550	IF ( pos_x == ( i5 * dx + 0.5 * dx ) .AND. &
	551	nzb_s_inner(j3,i5) /=0 .AND. &
	552	pos_z < nzb_s_inner(j3,i5) * dz ) THEN
	553	particles(n)%x = 2.0 * pos_x - prt_x
	554	particles(n)%speed_x = - particles(n)%speed_x
	555
	556	IF ( use_sgs_for_particles ) THEN
	557	particles(n)%speed_x_sgs = - particles(n)%speed_x_sgs
	558	ENDIF
	559	GOTO 999
	560	ENDIF
	561	ENDIF
	562
	563	IF ( k5 <= nzb_s_inner(j5,i3) .AND. nzb_s_inner(j5,i3) /= 0 )&
	564	THEN
	565	IF ( pos_z == nzb_s_inner(j5,i3) * dz .AND. &
	566	k5 == nzb_s_inner(j5,i3) ) THEN
	567	particles(n)%z = 2.0 * pos_z - prt_z
	568	particles(n)%speed_z = - particles(n)%speed_z
	569
	570	IF ( use_sgs_for_particles ) THEN
	571	particles(n)%speed_z_sgs = - particles(n)%speed_z_sgs
	572	ENDIF
	573	GOTO 999
	574	ENDIF
	575
	576	IF ( pos_y == ( j5 * dy + 0.5 * dy ) .AND. &
	577	nzb_s_inner(j5,i3) /= 0 .AND. &
	578	pos_z < nzb_s_inner(j5,i3) * dz ) THEN
	579	particles(n)%y = 2.0 * pos_y - prt_y
	580	particles(n)%speed_y = - particles(n)%speed_y
	581
	582	IF ( use_sgs_for_particles ) THEN
	583	particles(n)%speed_y_sgs = - particles(n)%speed_y_sgs
	584	ENDIF
	585	GOTO 999
	586	ENDIF
	587	ENDIF
	588	ENDDO
	589	ENDIF
	590
	591	999 CONTINUE
	592
	593	ENDIF
	594
	595	ENDDO
	596
	597	CALL cpu_log( log_point_s(48), 'advec_part_refle', 'stop' )
	598
	599	END SUBROUTINE particle_boundary_conds

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

Download in other formats:

| Impressum | ©Leibniz Universität Hannover |