source: palm/tags/release-3.4a/SOURCE/init_particles.f90 @ 710

Last change on this file since 710 was 139, checked in by raasch, 16 years ago

New:
---

Plant canopy model of Watanabe (2004,BLM 112,307-341) added.
It can be switched on by the inipar parameter plant_canopy.
The inipar parameter canopy_mode can be used to prescribe a
plant canopy type. The default case is a homogeneous plant
canopy. Heterogeneous distributions of the leaf area
density and the canopy drag coefficient can be defined in the
new routine user_init_plant_canopy (user_interface).
The inipar parameters lad_surface, lad_vertical_gradient and
lad_vertical_gradient_level can be used in order to
prescribe the vertical profile of leaf area density. The
inipar parameter drag_coefficient determines the canopy
drag coefficient.
Finally, the inipar parameter pch_index determines the
index of the upper boundary of the plant canopy.

Allow new case bc_uv_t = 'dirichlet_0' for channel flow.

For unknown variables (CASE DEFAULT) call new subroutine user_data_output_dvrp

Pressure boundary conditions for vertical walls added to the multigrid solver.
They are applied using new wall flag arrays (wall_flags_..) which are defined
for each grid level. New argument gls added to routine user_init_grid
(user_interface).

Frequence of sorting particles can be controlled with new particles_par
parameter dt_sort_particles. Sorting is moved from the SGS timestep loop in
advec_particles after the end of this loop.

advec_particles, check_parameters, data_output_dvrp, header, init_3d_model, init_grid, init_particles, init_pegrid, modules, package_parin, parin, plant_canopy_model, read_var_list, read_3d_binary, user_interface, write_var_list, write_3d_binary

Changed:


Redefine initial nzb_local as the actual total size of topography (later the
extent of topography in nzb_local is reduced by 1dx at the E topography walls
and by 1dy at the N topography walls to form the basis for nzb_s_inner);
for consistency redefine 'single_building' case.

Vertical profiles now based on nzb_s_inner; they are divided by
ngp_2dh_s_inner (scalars, procucts of scalars) and ngp_2dh (staggered velocity
components and their products, procucts of scalars and velocity components),
respectively.

Allow two instead of one digit to specify isosurface and slicer variables.

Status of 3D-volume NetCDF data file only depends on switch netcdf_64bit_3d (check_open)

prognostic_equations include the respective wall_*flux in the parameter list of
calls of diffusion_s. Same as before, only the values of wall_heatflux(0:4)
can be assigned. At present, wall_humidityflux, wall_qflux, wall_salinityflux,
and wall_scalarflux are kept zero. diffusion_s uses the respective wall_*flux
instead of wall_heatflux. This update serves two purposes:

  • it avoids errors in calculations with humidity/scalar/salinity and prescribed

non-zero wall_heatflux,

  • it prepares PALM for a possible assignment of wall fluxes of

humidity/scalar/salinity in a future release.

buoyancy, check_open, data_output_dvrp, diffusion_s, diffusivities, flow_statistics, header, init_3d_model, init_dvrp, init_grid, modules, prognostic_equations

Errors:


Bugfix: summation of sums_l_l in diffusivities.

Several bugfixes in the ocean part: Initial density rho is calculated
(init_ocean). Error in initializing u_init and v_init removed
(check_parameters). Calculation of density flux now starts from
nzb+1 (production_e).

Bugfix: pleft/pright changed to pnorth/psouth in sendrecv of particle tail
numbers along y, small bugfixes in the SGS part (advec_particles)

Bugfix: model_string needed a default value (combine_plot_fields)

Bugfix: wavenumber calculation for even nx in routines maketri (poisfft)

Bugfix: assignment of fluxes at walls

Bugfix: absolute value of f must be used when calculating the Blackadar mixing length (init_1d_model)

advec_particles, check_parameters, combine_plot_fields, diffusion_s, diffusivities, init_ocean, init_1d_model, poisfft, production_e

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