source:
palm/trunk/SOURCE/init_particles.f90
@
402
Last change on this file since 402 was 392, checked in by raasch, 15 years ago | |
---|---|
|
|
File size: 21.8 KB |
Rev | Line | |
---|---|---|
[1] | 1 | SUBROUTINE init_particles |
2 | ||
3 | !------------------------------------------------------------------------------! | |
[254] | 4 | ! Current revisions: |
[1] | 5 | ! ----------------- |
[392] | 6 | ! |
7 | ! | |
8 | ! Former revisions: | |
9 | ! ----------------- | |
10 | ! $Id: init_particles.f90 392 2009-09-24 10:39:14Z maronga $ | |
11 | ! | |
12 | ! 336 2009-06-10 11:19:35Z raasch | |
[276] | 13 | ! Maximum number of tails is calculated from maximum number of particles and |
14 | ! skip_particles_for_tail, | |
15 | ! output of messages replaced by message handling routine | |
[229] | 16 | ! Bugfix: arrays for tails are allocated with a minimum size of 10 tails if |
17 | ! there is no tail initially | |
[1] | 18 | ! |
[198] | 19 | ! 150 2008-02-29 08:19:58Z raasch |
20 | ! Setting offset_ocean_* needed for calculating vertical indices within ocean | |
21 | ! runs | |
22 | ! | |
[139] | 23 | ! 117 2007-10-11 03:27:59Z raasch |
24 | ! Sorting of particles only in case of cloud droplets | |
25 | ! | |
[110] | 26 | ! 106 2007-08-16 14:30:26Z raasch |
27 | ! variable iran replaced by iran_part | |
28 | ! | |
[83] | 29 | ! 82 2007-04-16 15:40:52Z raasch |
30 | ! Preprocessor directives for old systems removed | |
31 | ! | |
[77] | 32 | ! 70 2007-03-18 23:46:30Z raasch |
33 | ! displacements for mpi_particle_type changed, age_m initialized, | |
34 | ! particles-package is now part of the default code | |
35 | ! | |
[39] | 36 | ! 16 2007-02-15 13:16:47Z raasch |
37 | ! Bugfix: MPI_REAL in MPI_ALLREDUCE replaced by MPI_INTEGER | |
38 | ! | |
39 | ! r4 | raasch | 2007-02-13 12:33:16 +0100 (Tue, 13 Feb 2007) | |
[3] | 40 | ! RCS Log replace by Id keyword, revision history cleaned up |
41 | ! | |
[1] | 42 | ! Revision 1.24 2007/02/11 13:00:17 raasch |
43 | ! Bugfix: allocation of tail_mask and new_tail_id in case of restart-runs | |
44 | ! Bugfix: __ was missing in a cpp-directive | |
45 | ! | |
46 | ! Revision 1.1 1999/11/25 16:22:38 raasch | |
47 | ! Initial revision | |
48 | ! | |
49 | ! | |
50 | ! Description: | |
51 | ! ------------ | |
52 | ! This routine initializes a set of particles and their attributes (position, | |
53 | ! radius, ..). Advection of these particles is carried out by advec_particles, | |
54 | ! plotting is done in data_output_dvrp. | |
55 | !------------------------------------------------------------------------------! | |
56 | ||
57 | USE arrays_3d | |
58 | USE control_parameters | |
[336] | 59 | USE dvrp_variables |
[1] | 60 | USE grid_variables |
61 | USE indices | |
62 | USE particle_attributes | |
63 | USE pegrid | |
64 | USE random_function_mod | |
65 | ||
66 | ||
67 | IMPLICIT NONE | |
68 | ||
69 | CHARACTER (LEN=10) :: particle_binary_version, version_on_file | |
70 | ||
71 | INTEGER :: i, j, n, nn | |
72 | #if defined( __parallel ) | |
73 | INTEGER, DIMENSION(3) :: blocklengths, displacements, types | |
74 | #endif | |
75 | LOGICAL :: uniform_particles_l | |
76 | REAL :: factor, pos_x, pos_y, pos_z, value | |
77 | ||
78 | ||
79 | #if defined( __parallel ) | |
80 | ! | |
81 | !-- Define MPI derived datatype for FORTRAN datatype particle_type (see module | |
[82] | 82 | !-- particle_attributes). Integer length is 4 byte, Real is 8 byte |
83 | blocklengths(1) = 19; blocklengths(2) = 4; blocklengths(3) = 1 | |
84 | displacements(1) = 0; displacements(2) = 152; displacements(3) = 168 | |
85 | ||
[1] | 86 | types(1) = MPI_REAL |
87 | types(2) = MPI_INTEGER | |
88 | types(3) = MPI_UB | |
89 | CALL MPI_TYPE_STRUCT( 3, blocklengths, displacements, types, & | |
90 | mpi_particle_type, ierr ) | |
91 | CALL MPI_TYPE_COMMIT( mpi_particle_type, ierr ) | |
92 | #endif | |
93 | ||
94 | ! | |
[150] | 95 | !-- In case of oceans runs, the vertical index calculations need an offset, |
96 | !-- because otherwise the k indices will become negative | |
97 | IF ( ocean ) THEN | |
98 | offset_ocean_nzt = nzt | |
99 | offset_ocean_nzt_m1 = nzt - 1 | |
100 | ENDIF | |
101 | ||
102 | ||
103 | ! | |
[1] | 104 | !-- Check the number of particle groups. |
105 | IF ( number_of_particle_groups > max_number_of_particle_groups ) THEN | |
[274] | 106 | WRITE( message_string, * ) 'max_number_of_particle_groups =', & |
107 | max_number_of_particle_groups , & | |
[254] | 108 | '&number_of_particle_groups reset to ', & |
109 | max_number_of_particle_groups | |
110 | CALL message( 'init_particles', 'PA0213', 0, 1, 0, 6, 0 ) | |
[1] | 111 | number_of_particle_groups = max_number_of_particle_groups |
112 | ENDIF | |
113 | ||
114 | ! | |
115 | !-- Set default start positions, if necessary | |
116 | IF ( psl(1) == 9999999.9 ) psl(1) = -0.5 * dx | |
117 | IF ( psr(1) == 9999999.9 ) psr(1) = ( nx + 0.5 ) * dx | |
118 | IF ( pss(1) == 9999999.9 ) pss(1) = -0.5 * dy | |
119 | IF ( psn(1) == 9999999.9 ) psn(1) = ( ny + 0.5 ) * dy | |
120 | IF ( psb(1) == 9999999.9 ) psb(1) = zu(nz/2) | |
121 | IF ( pst(1) == 9999999.9 ) pst(1) = psb(1) | |
122 | ||
123 | IF ( pdx(1) == 9999999.9 .OR. pdx(1) == 0.0 ) pdx(1) = dx | |
124 | IF ( pdy(1) == 9999999.9 .OR. pdy(1) == 0.0 ) pdy(1) = dy | |
125 | IF ( pdz(1) == 9999999.9 .OR. pdz(1) == 0.0 ) pdz(1) = zu(2) - zu(1) | |
126 | ||
127 | DO j = 2, number_of_particle_groups | |
128 | IF ( psl(j) == 9999999.9 ) psl(j) = psl(j-1) | |
129 | IF ( psr(j) == 9999999.9 ) psr(j) = psr(j-1) | |
130 | IF ( pss(j) == 9999999.9 ) pss(j) = pss(j-1) | |
131 | IF ( psn(j) == 9999999.9 ) psn(j) = psn(j-1) | |
132 | IF ( psb(j) == 9999999.9 ) psb(j) = psb(j-1) | |
133 | IF ( pst(j) == 9999999.9 ) pst(j) = pst(j-1) | |
134 | IF ( pdx(j) == 9999999.9 .OR. pdx(j) == 0.0 ) pdx(j) = pdx(j-1) | |
135 | IF ( pdy(j) == 9999999.9 .OR. pdy(j) == 0.0 ) pdy(j) = pdy(j-1) | |
136 | IF ( pdz(j) == 9999999.9 .OR. pdz(j) == 0.0 ) pdz(j) = pdz(j-1) | |
137 | ENDDO | |
138 | ||
139 | ! | |
140 | !-- For the first model run of a possible job chain initialize the | |
141 | !-- particles, otherwise read the particle data from file. | |
142 | IF ( TRIM( initializing_actions ) == 'read_restart_data' & | |
143 | .AND. read_particles_from_restartfile ) THEN | |
144 | ||
145 | ! | |
146 | !-- Read particle data from previous model run. | |
147 | !-- First open the input unit. | |
148 | IF ( myid_char == '' ) THEN | |
149 | OPEN ( 90, FILE='PARTICLE_RESTART_DATA_IN'//myid_char, & | |
150 | FORM='UNFORMATTED' ) | |
151 | ELSE | |
152 | OPEN ( 90, FILE='PARTICLE_RESTART_DATA_IN/'//myid_char, & | |
153 | FORM='UNFORMATTED' ) | |
154 | ENDIF | |
155 | ||
156 | ! | |
157 | !-- First compare the version numbers | |
158 | READ ( 90 ) version_on_file | |
159 | particle_binary_version = '3.0' | |
160 | IF ( TRIM( version_on_file ) /= TRIM( particle_binary_version ) ) THEN | |
[274] | 161 | message_string = 'version mismatch concerning data from prior ' // & |
162 | 'run &version on file = "' // & | |
163 | TRIM( version_on_file ) // & | |
164 | '&version in program = "' // & | |
165 | TRIM( particle_binary_version ) // '"' | |
[254] | 166 | CALL message( 'init_particles', 'PA0214', 1, 2, 0, 6, 0 ) |
[1] | 167 | ENDIF |
168 | ||
169 | ! | |
170 | !-- Read some particle parameters and the size of the particle arrays, | |
171 | !-- allocate them and read their contents. | |
172 | READ ( 90 ) bc_par_b, bc_par_lr, bc_par_ns, bc_par_t, & | |
173 | maximum_number_of_particles, maximum_number_of_tailpoints, & | |
174 | maximum_number_of_tails, number_of_initial_particles, & | |
175 | number_of_particles, number_of_particle_groups, & | |
176 | number_of_tails, particle_groups, time_prel, & | |
177 | time_write_particle_data, uniform_particles | |
178 | ||
179 | IF ( number_of_initial_particles /= 0 ) THEN | |
180 | ALLOCATE( initial_particles(1:number_of_initial_particles) ) | |
181 | READ ( 90 ) initial_particles | |
182 | ENDIF | |
183 | ||
184 | ALLOCATE( prt_count(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1), & | |
185 | prt_start_index(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1), & | |
186 | particle_mask(maximum_number_of_particles), & | |
187 | particles(maximum_number_of_particles) ) | |
188 | ||
189 | READ ( 90 ) prt_count, prt_start_index | |
190 | READ ( 90 ) particles | |
191 | ||
192 | IF ( use_particle_tails ) THEN | |
193 | ALLOCATE( particle_tail_coordinates(maximum_number_of_tailpoints,5, & | |
194 | maximum_number_of_tails), & | |
195 | new_tail_id(maximum_number_of_tails), & | |
196 | tail_mask(maximum_number_of_tails) ) | |
197 | READ ( 90 ) particle_tail_coordinates | |
198 | ENDIF | |
199 | ||
200 | CLOSE ( 90 ) | |
201 | ||
202 | ELSE | |
203 | ||
204 | ! | |
205 | !-- Allocate particle arrays and set attributes of the initial set of | |
206 | !-- particles, which can be also periodically released at later times. | |
207 | !-- Also allocate array for particle tail coordinates, if needed. | |
208 | ALLOCATE( prt_count(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1), & | |
209 | prt_start_index(nzb:nzt+1,nys-1:nyn+1,nxl-1:nxr+1), & | |
210 | particle_mask(maximum_number_of_particles), & | |
211 | particles(maximum_number_of_particles) ) | |
212 | ||
213 | ! | |
214 | !-- Initialize all particles with dummy values (otherwise errors may | |
215 | !-- occur within restart runs). The reason for this is still not clear | |
216 | !-- and may be presumably caused by errors in the respective user-interface. | |
217 | particles = particle_type( 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, & | |
218 | 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, & | |
[57] | 219 | 0.0, 0, 0, 0, 0 ) |
[1] | 220 | particle_groups = particle_groups_type( 0.0, 0.0, 0.0, 0.0 ) |
221 | ||
222 | ! | |
223 | !-- Set the default particle size used for dvrp plots | |
224 | IF ( dvrp_psize == 9999999.9 ) dvrp_psize = 0.2 * dx | |
225 | ||
226 | ! | |
227 | !-- Set values for the density ratio and radius for all particle | |
228 | !-- groups, if necessary | |
229 | IF ( density_ratio(1) == 9999999.9 ) density_ratio(1) = 0.0 | |
230 | IF ( radius(1) == 9999999.9 ) radius(1) = 0.0 | |
231 | DO i = 2, number_of_particle_groups | |
232 | IF ( density_ratio(i) == 9999999.9 ) THEN | |
233 | density_ratio(i) = density_ratio(i-1) | |
234 | ENDIF | |
235 | IF ( radius(i) == 9999999.9 ) radius(i) = radius(i-1) | |
236 | ENDDO | |
237 | ||
238 | DO i = 1, number_of_particle_groups | |
239 | IF ( density_ratio(i) /= 0.0 .AND. radius(i) == 0 ) THEN | |
[254] | 240 | WRITE( message_string, * ) 'particle group #', i, 'has a', & |
241 | 'density ratio /= 0 but radius = 0' | |
242 | CALL message( 'init_particles', 'PA0215', 1, 2, 0, 6, 0 ) | |
[1] | 243 | ENDIF |
244 | particle_groups(i)%density_ratio = density_ratio(i) | |
245 | particle_groups(i)%radius = radius(i) | |
246 | ENDDO | |
247 | ||
248 | ! | |
249 | !-- Calculate particle positions and store particle attributes, if | |
250 | !-- particle is situated on this PE | |
251 | n = 0 | |
252 | ||
253 | DO i = 1, number_of_particle_groups | |
254 | ||
255 | pos_z = psb(i) | |
256 | ||
257 | DO WHILE ( pos_z <= pst(i) ) | |
258 | ||
259 | pos_y = pss(i) | |
260 | ||
261 | DO WHILE ( pos_y <= psn(i) ) | |
262 | ||
263 | IF ( pos_y >= ( nys - 0.5 ) * dy .AND. & | |
264 | pos_y < ( nyn + 0.5 ) * dy ) THEN | |
265 | ||
266 | pos_x = psl(i) | |
267 | ||
268 | DO WHILE ( pos_x <= psr(i) ) | |
269 | ||
270 | IF ( pos_x >= ( nxl - 0.5 ) * dx .AND. & | |
271 | pos_x < ( nxr + 0.5 ) * dx ) THEN | |
272 | ||
273 | DO j = 1, particles_per_point | |
274 | ||
275 | n = n + 1 | |
276 | IF ( n > maximum_number_of_particles ) THEN | |
[254] | 277 | WRITE( message_string, * ) 'number of initial', & |
[274] | 278 | 'particles (', n, ') exceeds', & |
279 | '&maximum_number_of_particles (', & | |
280 | maximum_number_of_particles, ') on PE ', & | |
[254] | 281 | myid |
[274] | 282 | CALL message( 'init_particles', 'PA0216', & |
[277] | 283 | 2, 2, -1, 6, 1 ) |
[1] | 284 | ENDIF |
285 | particles(n)%x = pos_x | |
286 | particles(n)%y = pos_y | |
287 | particles(n)%z = pos_z | |
288 | particles(n)%age = 0.0 | |
[57] | 289 | particles(n)%age_m = 0.0 |
[1] | 290 | particles(n)%dt_sum = 0.0 |
291 | particles(n)%dvrp_psize = dvrp_psize | |
292 | particles(n)%e_m = 0.0 | |
293 | particles(n)%speed_x = 0.0 | |
294 | particles(n)%speed_x_sgs = 0.0 | |
295 | particles(n)%speed_y = 0.0 | |
296 | particles(n)%speed_y_sgs = 0.0 | |
297 | particles(n)%speed_z = 0.0 | |
298 | particles(n)%speed_z_sgs = 0.0 | |
299 | particles(n)%origin_x = pos_x | |
300 | particles(n)%origin_y = pos_y | |
301 | particles(n)%origin_z = pos_z | |
302 | particles(n)%radius = particle_groups(i)%radius | |
303 | particles(n)%weight_factor =initial_weighting_factor | |
304 | particles(n)%color = 1 | |
305 | particles(n)%group = i | |
306 | particles(n)%tailpoints = 0 | |
307 | IF ( use_particle_tails .AND. & | |
308 | MOD( n, skip_particles_for_tail ) == 0 ) THEN | |
309 | number_of_tails = number_of_tails + 1 | |
310 | ! | |
311 | !-- This is a temporary provisional setting (see | |
312 | !-- further below!) | |
313 | particles(n)%tail_id = number_of_tails | |
314 | ELSE | |
315 | particles(n)%tail_id = 0 | |
316 | ENDIF | |
317 | ||
318 | ENDDO | |
319 | ||
320 | ENDIF | |
321 | ||
322 | pos_x = pos_x + pdx(i) | |
323 | ||
324 | ENDDO | |
325 | ||
326 | ENDIF | |
327 | ||
328 | pos_y = pos_y + pdy(i) | |
329 | ||
330 | ENDDO | |
331 | ||
332 | pos_z = pos_z + pdz(i) | |
333 | ||
334 | ENDDO | |
335 | ||
336 | ENDDO | |
337 | ||
338 | number_of_initial_particles = n | |
339 | number_of_particles = n | |
340 | ||
341 | ! | |
342 | !-- Calculate the number of particles and tails of the total domain | |
343 | #if defined( __parallel ) | |
344 | CALL MPI_ALLREDUCE( number_of_particles, total_number_of_particles, 1, & | |
[16] | 345 | MPI_INTEGER, MPI_SUM, comm2d, ierr ) |
[1] | 346 | CALL MPI_ALLREDUCE( number_of_tails, total_number_of_tails, 1, & |
[16] | 347 | MPI_INTEGER, MPI_SUM, comm2d, ierr ) |
[1] | 348 | #else |
349 | total_number_of_particles = number_of_particles | |
350 | total_number_of_tails = number_of_tails | |
351 | #endif | |
352 | ||
353 | ! | |
354 | !-- Set a seed value for the random number generator to be exclusively | |
355 | !-- used for the particle code. The generated random numbers should be | |
356 | !-- different on the different PEs. | |
357 | iran_part = iran_part + myid | |
358 | ||
359 | ! | |
360 | !-- User modification of initial particles | |
361 | CALL user_init_particles | |
362 | ||
363 | ! | |
364 | !-- Store the initial set of particles for release at later times | |
365 | IF ( number_of_initial_particles /= 0 ) THEN | |
366 | ALLOCATE( initial_particles(1:number_of_initial_particles) ) | |
367 | initial_particles(1:number_of_initial_particles) = & | |
368 | particles(1:number_of_initial_particles) | |
369 | ENDIF | |
370 | ||
371 | ! | |
372 | !-- Add random fluctuation to particle positions | |
373 | IF ( random_start_position ) THEN | |
374 | ||
375 | DO n = 1, number_of_initial_particles | |
376 | IF ( psl(particles(n)%group) /= psr(particles(n)%group) ) THEN | |
377 | particles(n)%x = particles(n)%x + & | |
[106] | 378 | ( random_function( iran_part ) - 0.5 ) * & |
[1] | 379 | pdx(particles(n)%group) |
380 | IF ( particles(n)%x <= ( nxl - 0.5 ) * dx ) THEN | |
381 | particles(n)%x = ( nxl - 0.4999999999 ) * dx | |
382 | ELSEIF ( particles(n)%x >= ( nxr + 0.5 ) * dx ) THEN | |
383 | particles(n)%x = ( nxr + 0.4999999999 ) * dx | |
384 | ENDIF | |
385 | ENDIF | |
386 | IF ( pss(particles(n)%group) /= psn(particles(n)%group) ) THEN | |
387 | particles(n)%y = particles(n)%y + & | |
[106] | 388 | ( random_function( iran_part ) - 0.5 ) * & |
[1] | 389 | pdy(particles(n)%group) |
390 | IF ( particles(n)%y <= ( nys - 0.5 ) * dy ) THEN | |
391 | particles(n)%y = ( nys - 0.4999999999 ) * dy | |
392 | ELSEIF ( particles(n)%y >= ( nyn + 0.5 ) * dy ) THEN | |
393 | particles(n)%y = ( nyn + 0.4999999999 ) * dy | |
394 | ENDIF | |
395 | ENDIF | |
396 | IF ( psb(particles(n)%group) /= pst(particles(n)%group) ) THEN | |
397 | particles(n)%z = particles(n)%z + & | |
[106] | 398 | ( random_function( iran_part ) - 0.5 ) * & |
[1] | 399 | pdz(particles(n)%group) |
400 | ENDIF | |
401 | ENDDO | |
402 | ENDIF | |
403 | ||
404 | ! | |
[117] | 405 | !-- Sort particles in the sequence the gridboxes are stored in the memory. |
406 | !-- Only required if cloud droplets are used. | |
407 | IF ( cloud_droplets ) CALL sort_particles | |
[1] | 408 | |
409 | ! | |
410 | !-- Open file for statistical informations about particle conditions | |
411 | IF ( write_particle_statistics ) THEN | |
412 | CALL check_open( 80 ) | |
413 | WRITE ( 80, 8000 ) current_timestep_number, simulated_time, & | |
414 | number_of_initial_particles, & | |
415 | maximum_number_of_particles | |
416 | CALL close_file( 80 ) | |
417 | ENDIF | |
418 | ||
419 | ! | |
420 | !-- Check if particles are really uniform in color and radius (dvrp_size) | |
421 | !-- (uniform_particles is preset TRUE) | |
422 | IF ( uniform_particles ) THEN | |
423 | IF ( number_of_initial_particles == 0 ) THEN | |
424 | uniform_particles_l = .TRUE. | |
425 | ELSE | |
426 | n = number_of_initial_particles | |
427 | IF ( MINVAL( particles(1:n)%dvrp_psize ) == & | |
428 | MAXVAL( particles(1:n)%dvrp_psize ) .AND. & | |
429 | MINVAL( particles(1:n)%color ) == & | |
430 | MAXVAL( particles(1:n)%color ) ) THEN | |
431 | uniform_particles_l = .TRUE. | |
432 | ELSE | |
433 | uniform_particles_l = .FALSE. | |
434 | ENDIF | |
435 | ENDIF | |
436 | ||
437 | #if defined( __parallel ) | |
438 | CALL MPI_ALLREDUCE( uniform_particles_l, uniform_particles, 1, & | |
439 | MPI_LOGICAL, MPI_LAND, comm2d, ierr ) | |
440 | #else | |
441 | uniform_particles = uniform_particles_l | |
442 | #endif | |
443 | ||
444 | ENDIF | |
445 | ||
446 | ! | |
[336] | 447 | !-- Particles will probably become none-uniform, if their size and color |
448 | !-- will be determined by flow variables | |
449 | IF ( particle_color /= 'none' .OR. particle_dvrpsize /= 'none' ) THEN | |
450 | uniform_particles = .FALSE. | |
451 | ENDIF | |
452 | ||
453 | ! | |
[1] | 454 | !-- Set the beginning of the particle tails and their age |
455 | IF ( use_particle_tails ) THEN | |
456 | ! | |
[276] | 457 | !-- Choose the maximum number of tails with respect to the maximum number |
458 | !-- of particles and skip_particles_for_tail | |
459 | maximum_number_of_tails = maximum_number_of_particles / & | |
460 | skip_particles_for_tail | |
461 | ||
[229] | 462 | ! |
463 | !-- Create a minimum number of tails in case that there is no tail | |
464 | !-- initially (otherwise, index errors will occur when adressing the | |
465 | !-- arrays below) | |
466 | IF ( maximum_number_of_tails == 0 ) maximum_number_of_tails = 10 | |
[1] | 467 | |
468 | ALLOCATE( particle_tail_coordinates(maximum_number_of_tailpoints,5, & | |
469 | maximum_number_of_tails), & | |
470 | new_tail_id(maximum_number_of_tails), & | |
471 | tail_mask(maximum_number_of_tails) ) | |
472 | ||
473 | particle_tail_coordinates = 0.0 | |
474 | minimum_tailpoint_distance = minimum_tailpoint_distance**2 | |
475 | number_of_initial_tails = number_of_tails | |
476 | ||
477 | nn = 0 | |
478 | DO n = 1, number_of_particles | |
479 | ! | |
480 | !-- Only for those particles marked above with a provisional tail_id | |
481 | !-- tails will be created. Particles now get their final tail_id. | |
482 | IF ( particles(n)%tail_id /= 0 ) THEN | |
483 | ||
484 | nn = nn + 1 | |
485 | particles(n)%tail_id = nn | |
486 | ||
487 | particle_tail_coordinates(1,1,nn) = particles(n)%x | |
488 | particle_tail_coordinates(1,2,nn) = particles(n)%y | |
489 | particle_tail_coordinates(1,3,nn) = particles(n)%z | |
490 | particle_tail_coordinates(1,4,nn) = particles(n)%color | |
491 | particles(n)%tailpoints = 1 | |
492 | IF ( minimum_tailpoint_distance /= 0.0 ) THEN | |
493 | particle_tail_coordinates(2,1,nn) = particles(n)%x | |
494 | particle_tail_coordinates(2,2,nn) = particles(n)%y | |
495 | particle_tail_coordinates(2,3,nn) = particles(n)%z | |
496 | particle_tail_coordinates(2,4,nn) = particles(n)%color | |
497 | particle_tail_coordinates(1:2,5,nn) = 0.0 | |
498 | particles(n)%tailpoints = 2 | |
499 | ENDIF | |
500 | ||
501 | ENDIF | |
502 | ENDDO | |
503 | ENDIF | |
504 | ||
505 | ! | |
506 | !-- Plot initial positions of particles (only if particle advection is | |
507 | !-- switched on from the beginning of the simulation (t=0)) | |
508 | IF ( particle_advection_start == 0.0 ) CALL data_output_dvrp | |
509 | ||
510 | ENDIF | |
511 | ||
512 | ! | |
513 | !-- Check boundary condition and set internal variables | |
514 | SELECT CASE ( bc_par_b ) | |
515 | ||
516 | CASE ( 'absorb' ) | |
517 | ibc_par_b = 1 | |
518 | ||
519 | CASE ( 'reflect' ) | |
520 | ibc_par_b = 2 | |
521 | ||
522 | CASE DEFAULT | |
[254] | 523 | WRITE( message_string, * ) 'unknown boundary condition ', & |
524 | 'bc_par_b = "', TRIM( bc_par_b ), '"' | |
525 | CALL message( 'init_particles', 'PA0217', 1, 2, 0, 6, 0 ) | |
[1] | 526 | |
527 | END SELECT | |
528 | SELECT CASE ( bc_par_t ) | |
529 | ||
530 | CASE ( 'absorb' ) | |
531 | ibc_par_t = 1 | |
532 | ||
533 | CASE ( 'reflect' ) | |
534 | ibc_par_t = 2 | |
535 | ||
536 | CASE DEFAULT | |
[254] | 537 | WRITE( message_string, * ) 'unknown boundary condition ', & |
538 | 'bc_par_t = "', TRIM( bc_par_t ), '"' | |
539 | CALL message( 'init_particles', 'PA0218', 1, 2, 0, 6, 0 ) | |
[1] | 540 | |
541 | END SELECT | |
542 | SELECT CASE ( bc_par_lr ) | |
543 | ||
544 | CASE ( 'cyclic' ) | |
545 | ibc_par_lr = 0 | |
546 | ||
547 | CASE ( 'absorb' ) | |
548 | ibc_par_lr = 1 | |
549 | ||
550 | CASE ( 'reflect' ) | |
551 | ibc_par_lr = 2 | |
552 | ||
553 | CASE DEFAULT | |
[254] | 554 | WRITE( message_string, * ) 'unknown boundary condition ', & |
555 | 'bc_par_lr = "', TRIM( bc_par_lr ), '"' | |
556 | CALL message( 'init_particles', 'PA0219', 1, 2, 0, 6, 0 ) | |
[1] | 557 | |
558 | END SELECT | |
559 | SELECT CASE ( bc_par_ns ) | |
560 | ||
561 | CASE ( 'cyclic' ) | |
562 | ibc_par_ns = 0 | |
563 | ||
564 | CASE ( 'absorb' ) | |
565 | ibc_par_ns = 1 | |
566 | ||
567 | CASE ( 'reflect' ) | |
568 | ibc_par_ns = 2 | |
569 | ||
570 | CASE DEFAULT | |
[254] | 571 | WRITE( message_string, * ) 'unknown boundary condition ', & |
572 | 'bc_par_ns = "', TRIM( bc_par_ns ), '"' | |
573 | CALL message( 'init_particles', 'PA0220', 1, 2, 0, 6, 0 ) | |
[1] | 574 | |
575 | END SELECT | |
576 | ! | |
577 | !-- Formats | |
578 | 8000 FORMAT (I6,1X,F7.2,4X,I6,71X,I6) | |
579 | ||
580 | END SUBROUTINE init_particles |
Note: See TracBrowser
for help on using the repository browser.