source: palm/trunk/SOURCE/lpm_exchange_horiz.f90 @ 1322

Last change on this file since 1322 was 1321, checked in by raasch, 10 years ago

last commit documented

  • Property svn:keywords set to Id
File size: 39.5 KB
Line 
1 SUBROUTINE lpm_exchange_horiz
2
3!--------------------------------------------------------------------------------!
4! This file is part of PALM.
5!
6! PALM is free software: you can redistribute it and/or modify it under the terms
7! of the GNU General Public License as published by the Free Software Foundation,
8! either version 3 of the License, or (at your option) any later version.
9!
10! PALM is distributed in the hope that it will be useful, but WITHOUT ANY
11! WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
12! A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
13!
14! You should have received a copy of the GNU General Public License along with
15! PALM. If not, see <http://www.gnu.org/licenses/>.
16!
17! Copyright 1997-2014 Leibniz Universitaet Hannover
18!--------------------------------------------------------------------------------!
19!
20! Current revisions:
21! ------------------
22!
23!
24! Former revisions:
25! -----------------
26! $Id: lpm_exchange_horiz.f90 1321 2014-03-20 09:40:40Z raasch $
27!
28! 1320 2014-03-20 08:40:49Z raasch
29! ONLY-attribute added to USE-statements,
30! kind-parameters added to all INTEGER and REAL declaration statements,
31! kinds are defined in new module kinds,
32! comment fields (!:) to be used for variable explanations added to
33! all variable declaration statements
34!
35! 1318 2014-03-17 13:35:16Z raasch
36! module interfaces removed
37!
38! 1036 2012-10-22 13:43:42Z raasch
39! code put under GPL (PALM 3.9)
40!
41! 851 2012-03-15 14:32:58Z raasch
42! Bugfix: resetting of particle_mask and tail mask moved from end of this
43! routine to lpm
44!
45! 849 2012-03-15 10:35:09Z raasch
46! initial revision (former part of advec_particles)
47!
48!
49! Description:
50! ------------
51! Exchange of particles (and tails) between the subdomains.
52!------------------------------------------------------------------------------!
53
54    USE control_parameters,                                                    &
55        ONLY:  message_string, netcdf_output, netcdf_data_format
56
57    USE cpulog,                                                                &
58        ONLY:  cpu_log, log_point_s
59
60    USE grid_variables,                                                        &
61        ONLY:  ddx, ddy, dx, dy
62
63    USE indices,                                                               &
64        ONLY:  nx, nxl, nxr, ny, nyn, nys
65
66    USE kinds
67
68    USE particle_attributes,                                                   &
69        ONLY:  deleted_particles, deleted_tails, ibc_par_lr, ibc_par_ns,       &
70               maximum_number_of_particles, maximum_number_of_tails,           &
71               maximum_number_of_tailpoints, mpi_particle_type,                &
72               number_of_tails, number_of_particles, particles, particle_mask, &
73               particle_tail_coordinates, particle_type, tail_mask,            &
74               trlp_count_sum, trlp_count_recv_sum, trnp_count_sum,            &
75               trnp_count_recv_sum, trrp_count_sum, trrp_count_recv_sum,       &
76               trsp_count_sum, trsp_count_recv_sum, use_particle_tails
77
78    USE pegrid
79
80    IMPLICIT NONE
81
82    INTEGER(iwp) ::  i                                       !:
83    INTEGER(iwp) ::  j                                       !:
84    INTEGER(iwp) ::  n                                       !:
85    INTEGER(iwp) ::  nn                                      !:
86    INTEGER(iwp) ::  tlength                                 !:
87    INTEGER(iwp) ::  trlp_count                              !:
88    INTEGER(iwp) ::  trlp_count_recv                         !:
89    INTEGER(iwp) ::  trlpt_count                             !:
90    INTEGER(iwp) ::  trlpt_count_recv                        !:
91    INTEGER(iwp) ::  trnp_count                              !:
92    INTEGER(iwp) ::  trnp_count_recv                         !:
93    INTEGER(iwp) ::  trnpt_count                             !:
94    INTEGER(iwp) ::  trnpt_count_recv                        !:
95    INTEGER(iwp) ::  trrp_count                              !:
96    INTEGER(iwp) ::  trrp_count_recv                         !:
97    INTEGER(iwp) ::  trrpt_count                             !:
98    INTEGER(iwp) ::  trrpt_count_recv                        !:
99    INTEGER(iwp) ::  trsp_count                              !:
100    INTEGER(iwp) ::  trsp_count_recv                         !:
101    INTEGER(iwp) ::  trspt_count                             !:
102    INTEGER(iwp) ::  trspt_count_recv                        !:
103
104    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  trlpt        !:
105    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  trnpt        !:
106    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  trrpt        !:
107    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  trspt        !:
108
109
110    TYPE(particle_type), DIMENSION(:), ALLOCATABLE ::  trlp  !:
111    TYPE(particle_type), DIMENSION(:), ALLOCATABLE ::  trnp  !:
112    TYPE(particle_type), DIMENSION(:), ALLOCATABLE ::  trrp  !:
113    TYPE(particle_type), DIMENSION(:), ALLOCATABLE ::  trsp  !:
114
115#if defined( __parallel )
116
117!
118!-- Exchange between subdomains.
119!-- As soon as one particle has moved beyond the boundary of the domain, it
120!-- is included in the relevant transfer arrays and marked for subsequent
121!-- deletion on this PE.
122!-- First sweep for crossings in x direction. Find out first the number of
123!-- particles to be transferred and allocate temporary arrays needed to store
124!-- them.
125!-- For a one-dimensional decomposition along y, no transfer is necessary,
126!-- because the particle remains on the PE, but the particle coordinate has to
127!-- be adjusted.
128    trlp_count  = 0
129    trlpt_count = 0
130    trrp_count  = 0
131    trrpt_count = 0
132
133    trlp_count_recv   = 0
134    trlpt_count_recv  = 0
135    trrp_count_recv   = 0
136    trrpt_count_recv  = 0
137
138    IF ( pdims(1) /= 1 )  THEN
139!
140!--    First calculate the storage necessary for sending and receiving the data
141       DO  n = 1, number_of_particles
142          i = ( particles(n)%x + 0.5 * dx ) * ddx
143!
144!--       Above calculation does not work for indices less than zero
145          IF ( particles(n)%x < -0.5 * dx )  i = -1
146
147          IF ( i < nxl )  THEN
148             trlp_count = trlp_count + 1
149             IF ( particles(n)%tail_id /= 0 )  trlpt_count = trlpt_count + 1
150          ELSEIF ( i > nxr )  THEN
151             trrp_count = trrp_count + 1
152             IF ( particles(n)%tail_id /= 0 )  trrpt_count = trrpt_count + 1
153          ENDIF
154       ENDDO
155
156       IF ( trlp_count  == 0 )  trlp_count  = 1
157       IF ( trlpt_count == 0 )  trlpt_count = 1
158       IF ( trrp_count  == 0 )  trrp_count  = 1
159       IF ( trrpt_count == 0 )  trrpt_count = 1
160
161       ALLOCATE( trlp(trlp_count), trrp(trrp_count) )
162
163       trlp = particle_type( 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, &
164                             0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, &
165                             0.0, 0, 0, 0, 0 )
166       trrp = particle_type( 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, &
167                             0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, &
168                             0.0, 0, 0, 0, 0 )
169
170       IF ( use_particle_tails )  THEN
171          ALLOCATE( trlpt(maximum_number_of_tailpoints,5,trlpt_count), &
172                    trrpt(maximum_number_of_tailpoints,5,trrpt_count) )
173          tlength = maximum_number_of_tailpoints * 5
174       ENDIF
175
176       trlp_count  = 0
177       trlpt_count = 0
178       trrp_count  = 0
179       trrpt_count = 0
180
181    ENDIF
182
183    DO  n = 1, number_of_particles
184
185       nn = particles(n)%tail_id
186
187       i = ( particles(n)%x + 0.5 * dx ) * ddx
188!
189!--    Above calculation does not work for indices less than zero
190       IF ( particles(n)%x < - 0.5 * dx )  i = -1
191
192       IF ( i <  nxl )  THEN
193          IF ( i < 0 )  THEN
194!
195!--          Apply boundary condition along x
196             IF ( ibc_par_lr == 0 )  THEN
197!
198!--             Cyclic condition
199                IF ( pdims(1) == 1 )  THEN
200                   particles(n)%x        = ( nx + 1 ) * dx + particles(n)%x
201                   particles(n)%origin_x = ( nx + 1 ) * dx + &
202                                           particles(n)%origin_x
203                   IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
204                      i  = particles(n)%tailpoints
205                      particle_tail_coordinates(1:i,1,nn) = ( nx + 1 ) * dx &
206                                        + particle_tail_coordinates(1:i,1,nn)
207                   ENDIF
208                ELSE
209                   trlp_count = trlp_count + 1
210                   trlp(trlp_count)   = particles(n)
211                   trlp(trlp_count)%x = ( nx + 1 ) * dx + trlp(trlp_count)%x
212                   trlp(trlp_count)%origin_x = trlp(trlp_count)%origin_x + &
213                                               ( nx + 1 ) * dx
214                   particle_mask(n)  = .FALSE.
215                   deleted_particles = deleted_particles + 1
216
217                   IF ( trlp(trlp_count)%x >= (nx + 0.5)* dx - 1.0E-12 ) THEN
218                      trlp(trlp_count)%x = trlp(trlp_count)%x - 1.0E-10
219!++ why is 1 subtracted in next statement???
220                      trlp(trlp_count)%origin_x = trlp(trlp_count)%origin_x - 1
221                   ENDIF
222
223                   IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
224                      trlpt_count = trlpt_count + 1
225                      trlpt(:,:,trlpt_count) = particle_tail_coordinates(:,:,nn)
226                      trlpt(:,1,trlpt_count) = ( nx + 1 ) * dx + &
227                                               trlpt(:,1,trlpt_count)
228                      tail_mask(nn) = .FALSE.
229                      deleted_tails = deleted_tails + 1
230                   ENDIF
231                ENDIF
232
233             ELSEIF ( ibc_par_lr == 1 )  THEN
234!
235!--             Particle absorption
236                particle_mask(n) = .FALSE.
237                deleted_particles = deleted_particles + 1
238                IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
239                   tail_mask(nn) = .FALSE.
240                   deleted_tails = deleted_tails + 1
241                ENDIF
242
243             ELSEIF ( ibc_par_lr == 2 )  THEN
244!
245!--             Particle reflection
246                particles(n)%x       = -particles(n)%x
247                particles(n)%speed_x = -particles(n)%speed_x
248
249             ENDIF
250          ELSE
251!
252!--          Store particle data in the transfer array, which will be send
253!--          to the neighbouring PE
254             trlp_count = trlp_count + 1
255             trlp(trlp_count) = particles(n)
256             particle_mask(n) = .FALSE.
257             deleted_particles = deleted_particles + 1
258
259             IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
260                trlpt_count = trlpt_count + 1
261                trlpt(:,:,trlpt_count) = particle_tail_coordinates(:,:,nn)
262                tail_mask(nn) = .FALSE.
263                deleted_tails = deleted_tails + 1
264             ENDIF
265         ENDIF
266
267       ELSEIF ( i > nxr )  THEN
268          IF ( i > nx )  THEN
269!
270!--          Apply boundary condition along x
271             IF ( ibc_par_lr == 0 )  THEN
272!
273!--             Cyclic condition
274                IF ( pdims(1) == 1 )  THEN
275                   particles(n)%x = particles(n)%x - ( nx + 1 ) * dx
276                   particles(n)%origin_x = particles(n)%origin_x - &
277                                           ( nx + 1 ) * dx
278                   IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
279                      i = particles(n)%tailpoints
280                      particle_tail_coordinates(1:i,1,nn) = - ( nx+1 ) * dx &
281                                           + particle_tail_coordinates(1:i,1,nn)
282                   ENDIF
283                ELSE
284                   trrp_count = trrp_count + 1
285                   trrp(trrp_count) = particles(n)
286                   trrp(trrp_count)%x = trrp(trrp_count)%x - ( nx + 1 ) * dx
287                   trrp(trrp_count)%origin_x = trrp(trrp_count)%origin_x - &
288                                               ( nx + 1 ) * dx
289                   particle_mask(n) = .FALSE.
290                   deleted_particles = deleted_particles + 1
291
292                   IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
293                      trrpt_count = trrpt_count + 1
294                      trrpt(:,:,trrpt_count) = particle_tail_coordinates(:,:,nn)
295                      trrpt(:,1,trrpt_count) = trrpt(:,1,trrpt_count) - &
296                                               ( nx + 1 ) * dx
297                      tail_mask(nn) = .FALSE.
298                      deleted_tails = deleted_tails + 1
299                   ENDIF
300                ENDIF
301
302             ELSEIF ( ibc_par_lr == 1 )  THEN
303!
304!--             Particle absorption
305                particle_mask(n) = .FALSE.
306                deleted_particles = deleted_particles + 1
307                IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
308                   tail_mask(nn) = .FALSE.
309                   deleted_tails = deleted_tails + 1
310                ENDIF
311
312             ELSEIF ( ibc_par_lr == 2 )  THEN
313!
314!--             Particle reflection
315                particles(n)%x       = 2 * ( nx * dx ) - particles(n)%x
316                particles(n)%speed_x = -particles(n)%speed_x
317
318             ENDIF
319          ELSE
320!
321!--          Store particle data in the transfer array, which will be send
322!--          to the neighbouring PE
323             trrp_count = trrp_count + 1
324             trrp(trrp_count) = particles(n)
325             particle_mask(n) = .FALSE.
326             deleted_particles = deleted_particles + 1
327
328             IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
329                trrpt_count = trrpt_count + 1
330                trrpt(:,:,trrpt_count) = particle_tail_coordinates(:,:,nn)
331                tail_mask(nn) = .FALSE.
332                deleted_tails = deleted_tails + 1
333             ENDIF
334          ENDIF
335
336       ENDIF
337    ENDDO
338
339!
340!-- Send left boundary, receive right boundary (but first exchange how many
341!-- and check, if particle storage must be extended)
342    IF ( pdims(1) /= 1 )  THEN
343
344       CALL cpu_log( log_point_s(23), 'sendrcv_particles', 'start' )
345       CALL MPI_SENDRECV( trlp_count,      1, MPI_INTEGER, pleft,  0, &
346                          trrp_count_recv, 1, MPI_INTEGER, pright, 0, &
347                          comm2d, status, ierr )
348
349       IF ( number_of_particles + trrp_count_recv > &
350            maximum_number_of_particles )           &
351       THEN
352          IF ( netcdf_output  .AND.  netcdf_data_format < 3 )  THEN
353              message_string = 'maximum_number_of_particles ' //    &
354                               'needs to be increased ' //          &
355                               '&but this is not allowed with ' //  &
356                               'netcdf-data_format < 3'
357             CALL message( 'lpm_exch_horiz', 'PA0146', 2, 2, -1, 6, 1 )
358          ELSE
359             CALL lpm_extend_particle_array( trrp_count_recv )
360          ENDIF
361       ENDIF
362
363       CALL MPI_SENDRECV( trlp(1)%age, trlp_count, mpi_particle_type,     &
364                          pleft, 1, particles(number_of_particles+1)%age, &
365                          trrp_count_recv, mpi_particle_type, pright, 1,  &
366                          comm2d, status, ierr )
367
368       IF ( use_particle_tails )  THEN
369
370          CALL MPI_SENDRECV( trlpt_count,      1, MPI_INTEGER, pleft,  0, &
371                             trrpt_count_recv, 1, MPI_INTEGER, pright, 0, &
372                             comm2d, status, ierr )
373
374          IF ( number_of_tails+trrpt_count_recv > maximum_number_of_tails ) &
375          THEN
376             IF ( netcdf_output  .AND.  netcdf_data_format < 3 )  THEN
377                message_string = 'maximum_number_of_tails ' //   &
378                                 'needs to be increased ' //     &
379                                 '&but this is not allowed wi'// &
380                                 'th netcdf_data_format < 3'
381                CALL message( 'lpm_exch_horiz', 'PA0147', 2, 2, -1, 6, 1 )
382             ELSE
383                CALL lpm_extend_tail_array( trrpt_count_recv )
384             ENDIF
385          ENDIF
386
387          CALL MPI_SENDRECV( trlpt(1,1,1), trlpt_count*tlength, MPI_REAL,      &
388                             pleft, 1,                                         &
389                             particle_tail_coordinates(1,1,number_of_tails+1), &
390                             trrpt_count_recv*tlength, MPI_REAL, pright, 1,    &
391                             comm2d, status, ierr )
392!
393!--       Update the tail ids for the transferred particles
394          nn = number_of_tails
395          DO  n = number_of_particles+1, number_of_particles+trrp_count_recv
396             IF ( particles(n)%tail_id /= 0 )  THEN
397                nn = nn + 1
398                particles(n)%tail_id = nn
399             ENDIF
400          ENDDO
401
402       ENDIF
403
404       number_of_particles = number_of_particles + trrp_count_recv
405       number_of_tails     = number_of_tails     + trrpt_count_recv
406
407!
408!--    Send right boundary, receive left boundary
409       CALL MPI_SENDRECV( trrp_count,      1, MPI_INTEGER, pright, 0, &
410                          trlp_count_recv, 1, MPI_INTEGER, pleft,  0, &
411                          comm2d, status, ierr )
412
413       IF ( number_of_particles + trlp_count_recv > &
414            maximum_number_of_particles )           &
415       THEN
416          IF ( netcdf_output  .AND.  netcdf_data_format < 3 )  THEN
417             message_string = 'maximum_number_of_particles ' //  &
418                              'needs to be increased ' //        &
419                              '&but this is not allowed with '// &
420                              'netcdf_data_format < 3'
421             CALL message( 'lpm_exch_horiz', 'PA0146', 2, 2, -1, 6, 1 )
422          ELSE
423             CALL lpm_extend_particle_array( trlp_count_recv )
424          ENDIF
425       ENDIF
426
427       CALL MPI_SENDRECV( trrp(1)%age, trrp_count, mpi_particle_type,      &
428                          pright, 1, particles(number_of_particles+1)%age, &
429                          trlp_count_recv, mpi_particle_type, pleft, 1,    &
430                          comm2d, status, ierr )
431
432       IF ( use_particle_tails )  THEN
433
434          CALL MPI_SENDRECV( trrpt_count,      1, MPI_INTEGER, pright, 0, &
435                             trlpt_count_recv, 1, MPI_INTEGER, pleft,  0, &
436                             comm2d, status, ierr )
437
438          IF ( number_of_tails+trlpt_count_recv > maximum_number_of_tails ) &
439          THEN
440             IF ( netcdf_output  .AND.  netcdf_data_format < 3 )  THEN
441                message_string = 'maximum_number_of_tails ' //   &
442                                 'needs to be increased ' //     &
443                                 '&but this is not allowed wi'// &
444                                 'th netcdf_data_format < 3'
445                CALL message( 'lpm_exch_horiz', 'PA0147', 2, 2, -1, 6, 1 ) 
446             ELSE
447                CALL lpm_extend_tail_array( trlpt_count_recv )
448             ENDIF
449          ENDIF
450
451          CALL MPI_SENDRECV( trrpt(1,1,1), trrpt_count*tlength, MPI_REAL,      &
452                             pright, 1,                                        &
453                             particle_tail_coordinates(1,1,number_of_tails+1), &
454                             trlpt_count_recv*tlength, MPI_REAL, pleft, 1,     &
455                             comm2d, status, ierr )
456!
457!--       Update the tail ids for the transferred particles
458          nn = number_of_tails
459          DO  n = number_of_particles+1, number_of_particles+trlp_count_recv
460             IF ( particles(n)%tail_id /= 0 )  THEN
461                nn = nn + 1
462                particles(n)%tail_id = nn
463             ENDIF
464          ENDDO
465
466       ENDIF
467
468       number_of_particles = number_of_particles + trlp_count_recv
469       number_of_tails     = number_of_tails     + trlpt_count_recv
470
471       IF ( use_particle_tails )  THEN
472          DEALLOCATE( trlpt, trrpt )
473       ENDIF
474       DEALLOCATE( trlp, trrp ) 
475
476       CALL cpu_log( log_point_s(23), 'sendrcv_particles', 'pause' )
477
478    ENDIF
479
480
481!
482!-- Check whether particles have crossed the boundaries in y direction. Note
483!-- that this case can also apply to particles that have just been received
484!-- from the adjacent right or left PE.
485!-- Find out first the number of particles to be transferred and allocate
486!-- temporary arrays needed to store them.
487!-- For a one-dimensional decomposition along x, no transfer is necessary,
488!-- because the particle remains on the PE.
489    trsp_count  = 0
490    trspt_count = 0
491    trnp_count  = 0
492    trnpt_count = 0
493
494    trsp_count_recv   = 0
495    trspt_count_recv  = 0
496    trnp_count_recv   = 0
497    trnpt_count_recv  = 0
498
499    IF ( pdims(2) /= 1 )  THEN
500!
501!--    First calculate the storage necessary for sending and receiving the
502!--    data
503       DO  n = 1, number_of_particles
504          IF ( particle_mask(n) )  THEN
505             j = ( particles(n)%y + 0.5 * dy ) * ddy
506!
507!--          Above calculation does not work for indices less than zero
508             IF ( particles(n)%y < -0.5 * dy )  j = -1
509
510             IF ( j < nys )  THEN
511                trsp_count = trsp_count + 1
512                IF ( particles(n)%tail_id /= 0 )  trspt_count = trspt_count+1
513             ELSEIF ( j > nyn )  THEN
514                trnp_count = trnp_count + 1
515                IF ( particles(n)%tail_id /= 0 )  trnpt_count = trnpt_count+1
516             ENDIF
517          ENDIF
518       ENDDO
519
520       IF ( trsp_count  == 0 )  trsp_count  = 1
521       IF ( trspt_count == 0 )  trspt_count = 1
522       IF ( trnp_count  == 0 )  trnp_count  = 1
523       IF ( trnpt_count == 0 )  trnpt_count = 1
524
525       ALLOCATE( trsp(trsp_count), trnp(trnp_count) )
526
527       trsp = particle_type( 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, &
528                             0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, &
529                             0.0, 0, 0, 0, 0 )
530       trnp = particle_type( 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, &
531                             0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, &
532                             0.0, 0, 0, 0, 0 )
533
534       IF ( use_particle_tails )  THEN
535          ALLOCATE( trspt(maximum_number_of_tailpoints,5,trspt_count), &
536                    trnpt(maximum_number_of_tailpoints,5,trnpt_count) )
537          tlength = maximum_number_of_tailpoints * 5
538       ENDIF
539
540       trsp_count  = 0
541       trspt_count = 0
542       trnp_count  = 0
543       trnpt_count = 0
544
545    ENDIF
546
547    DO  n = 1, number_of_particles
548
549       nn = particles(n)%tail_id
550!
551!--    Only those particles that have not been marked as 'deleted' may be
552!--    moved.
553       IF ( particle_mask(n) )  THEN
554          j = ( particles(n)%y + 0.5 * dy ) * ddy
555!
556!--       Above calculation does not work for indices less than zero
557          IF ( particles(n)%y < -0.5 * dy )  j = -1
558
559          IF ( j < nys )  THEN
560             IF ( j < 0 )  THEN
561!
562!--             Apply boundary condition along y
563                IF ( ibc_par_ns == 0 )  THEN
564!
565!--                Cyclic condition
566                   IF ( pdims(2) == 1 )  THEN
567                      particles(n)%y = ( ny + 1 ) * dy + particles(n)%y
568                      particles(n)%origin_y = ( ny + 1 ) * dy + &
569                                              particles(n)%origin_y
570                      IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
571                         i = particles(n)%tailpoints
572                         particle_tail_coordinates(1:i,2,nn) = ( ny+1 ) * dy&
573                                        + particle_tail_coordinates(1:i,2,nn)
574                      ENDIF
575                   ELSE
576                      trsp_count = trsp_count + 1
577                      trsp(trsp_count) = particles(n)
578                      trsp(trsp_count)%y = ( ny + 1 ) * dy + &
579                                           trsp(trsp_count)%y
580                      trsp(trsp_count)%origin_y = trsp(trsp_count)%origin_y &
581                                                  + ( ny + 1 ) * dy
582                      particle_mask(n) = .FALSE.
583                      deleted_particles = deleted_particles + 1
584
585                      IF ( trsp(trsp_count)%y >= (ny+0.5)* dy - 1.0E-12 ) THEN
586                         trsp(trsp_count)%y = trsp(trsp_count)%y - 1.0E-10
587!++ why is 1 subtracted in next statement???
588                         trsp(trsp_count)%origin_y =                        &
589                                                   trsp(trsp_count)%origin_y - 1
590                      ENDIF
591
592                      IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
593                         trspt_count = trspt_count + 1
594                         trspt(:,:,trspt_count) = &
595                                               particle_tail_coordinates(:,:,nn)
596                         trspt(:,2,trspt_count) = ( ny + 1 ) * dy + &
597                                                  trspt(:,2,trspt_count)
598                         tail_mask(nn) = .FALSE.
599                         deleted_tails = deleted_tails + 1
600                      ENDIF
601                   ENDIF
602
603                ELSEIF ( ibc_par_ns == 1 )  THEN
604!
605!--                Particle absorption
606                   particle_mask(n) = .FALSE.
607                   deleted_particles = deleted_particles + 1
608                   IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
609                      tail_mask(nn) = .FALSE.
610                      deleted_tails = deleted_tails + 1
611                   ENDIF
612
613                ELSEIF ( ibc_par_ns == 2 )  THEN
614!
615!--                Particle reflection
616                   particles(n)%y       = -particles(n)%y
617                   particles(n)%speed_y = -particles(n)%speed_y
618
619                ENDIF
620             ELSE
621!
622!--             Store particle data in the transfer array, which will be send
623!--             to the neighbouring PE
624                trsp_count = trsp_count + 1
625                trsp(trsp_count) = particles(n)
626                particle_mask(n) = .FALSE.
627                deleted_particles = deleted_particles + 1
628
629                IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
630                   trspt_count = trspt_count + 1
631                   trspt(:,:,trspt_count) = particle_tail_coordinates(:,:,nn)
632                   tail_mask(nn) = .FALSE.
633                   deleted_tails = deleted_tails + 1
634                ENDIF
635             ENDIF
636
637          ELSEIF ( j > nyn )  THEN
638             IF ( j > ny )  THEN
639!
640!--             Apply boundary condition along x
641                IF ( ibc_par_ns == 0 )  THEN
642!
643!--                Cyclic condition
644                   IF ( pdims(2) == 1 )  THEN
645                      particles(n)%y = particles(n)%y - ( ny + 1 ) * dy
646                      particles(n)%origin_y = particles(n)%origin_y - &
647                                              ( ny + 1 ) * dy
648                      IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
649                         i = particles(n)%tailpoints
650                         particle_tail_coordinates(1:i,2,nn) = - (ny+1) * dy &
651                                           + particle_tail_coordinates(1:i,2,nn)
652                      ENDIF
653                   ELSE
654                      trnp_count = trnp_count + 1
655                      trnp(trnp_count) = particles(n)
656                      trnp(trnp_count)%y = trnp(trnp_count)%y - &
657                                           ( ny + 1 ) * dy
658                      trnp(trnp_count)%origin_y = trnp(trnp_count)%origin_y &
659                                                  - ( ny + 1 ) * dy
660                      particle_mask(n) = .FALSE.
661                      deleted_particles = deleted_particles + 1
662
663                      IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
664                         trnpt_count = trnpt_count + 1
665                         trnpt(:,:,trnpt_count) = &
666                                               particle_tail_coordinates(:,:,nn)
667                         trnpt(:,2,trnpt_count) = trnpt(:,2,trnpt_count) - &
668                                                  ( ny + 1 ) * dy
669                         tail_mask(nn) = .FALSE.
670                         deleted_tails = deleted_tails + 1
671                      ENDIF
672                   ENDIF
673
674                ELSEIF ( ibc_par_ns == 1 )  THEN
675!
676!--                Particle absorption
677                   particle_mask(n) = .FALSE.
678                   deleted_particles = deleted_particles + 1
679                   IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
680                      tail_mask(nn) = .FALSE.
681                      deleted_tails = deleted_tails + 1
682                   ENDIF
683
684                ELSEIF ( ibc_par_ns == 2 )  THEN
685!
686!--                Particle reflection
687                   particles(n)%y       = 2 * ( ny * dy ) - particles(n)%y
688                   particles(n)%speed_y = -particles(n)%speed_y
689
690                ENDIF
691             ELSE
692!
693!--             Store particle data in the transfer array, which will be send
694!--             to the neighbouring PE
695                trnp_count = trnp_count + 1
696                trnp(trnp_count) = particles(n)
697                particle_mask(n) = .FALSE.
698                deleted_particles = deleted_particles + 1
699
700                IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
701                   trnpt_count = trnpt_count + 1
702                   trnpt(:,:,trnpt_count) = particle_tail_coordinates(:,:,nn)
703                   tail_mask(nn) = .FALSE.
704                   deleted_tails = deleted_tails + 1
705                ENDIF
706             ENDIF
707
708          ENDIF
709       ENDIF
710    ENDDO
711
712!
713!-- Send front boundary, receive back boundary (but first exchange how many
714!-- and check, if particle storage must be extended)
715    IF ( pdims(2) /= 1 )  THEN
716
717       CALL cpu_log( log_point_s(23), 'sendrcv_particles', 'continue' )
718       CALL MPI_SENDRECV( trsp_count,      1, MPI_INTEGER, psouth, 0, &
719                          trnp_count_recv, 1, MPI_INTEGER, pnorth, 0, &
720                          comm2d, status, ierr )
721
722       IF ( number_of_particles + trnp_count_recv > &
723            maximum_number_of_particles )           &
724       THEN
725          IF ( netcdf_output  .AND.  netcdf_data_format < 3 )  THEN
726             message_string = 'maximum_number_of_particles ' //  &
727                              'needs to be increased ' //        &
728                              '&but this is not allowed with '// &
729                              'netcdf_data_format < 3'
730             CALL message( 'lpm_exch_horiz', 'PA0146', 2, 2, -1, 6, 1 ) 
731          ELSE
732             CALL lpm_extend_particle_array( trnp_count_recv )
733          ENDIF
734       ENDIF
735
736       CALL MPI_SENDRECV( trsp(1)%age, trsp_count, mpi_particle_type,      &
737                          psouth, 1, particles(number_of_particles+1)%age, &
738                          trnp_count_recv, mpi_particle_type, pnorth, 1,   &
739                          comm2d, status, ierr )
740
741       IF ( use_particle_tails )  THEN
742
743          CALL MPI_SENDRECV( trspt_count,      1, MPI_INTEGER, psouth, 0, &
744                             trnpt_count_recv, 1, MPI_INTEGER, pnorth, 0, &
745                             comm2d, status, ierr )
746
747          IF ( number_of_tails+trnpt_count_recv > maximum_number_of_tails ) &
748          THEN
749             IF ( netcdf_output  .AND.  netcdf_data_format < 3 )  THEN
750                message_string = 'maximum_number_of_tails ' //    &
751                                 'needs to be increased ' //      &
752                                 '&but this is not allowed wi' // &
753                                 'th netcdf_data_format < 3'
754                CALL message( 'lpm_exch_horiz', 'PA0147', 2, 2, -1, 6, 1 ) 
755             ELSE
756                CALL lpm_extend_tail_array( trnpt_count_recv )
757             ENDIF
758          ENDIF
759
760          CALL MPI_SENDRECV( trspt(1,1,1), trspt_count*tlength, MPI_REAL,      &
761                             psouth, 1,                                        &
762                             particle_tail_coordinates(1,1,number_of_tails+1), &
763                             trnpt_count_recv*tlength, MPI_REAL, pnorth, 1,    &
764                             comm2d, status, ierr )
765
766!
767!--       Update the tail ids for the transferred particles
768          nn = number_of_tails
769          DO  n = number_of_particles+1, number_of_particles+trnp_count_recv
770             IF ( particles(n)%tail_id /= 0 )  THEN
771                nn = nn + 1
772                particles(n)%tail_id = nn
773             ENDIF
774          ENDDO
775
776       ENDIF
777
778       number_of_particles = number_of_particles + trnp_count_recv
779       number_of_tails     = number_of_tails     + trnpt_count_recv
780
781!
782!--    Send back boundary, receive front boundary
783       CALL MPI_SENDRECV( trnp_count,      1, MPI_INTEGER, pnorth, 0, &
784                          trsp_count_recv, 1, MPI_INTEGER, psouth, 0, &
785                          comm2d, status, ierr )
786
787       IF ( number_of_particles + trsp_count_recv > &
788            maximum_number_of_particles )           &
789       THEN
790          IF ( netcdf_output  .AND.  netcdf_data_format < 3 )  THEN
791             message_string = 'maximum_number_of_particles ' //   &
792                              'needs to be increased ' //         &
793                              '&but this is not allowed with ' // &
794                              'netcdf_data_format < 3'
795            CALL message( 'lpm_exch_horiz', 'PA0146', 2, 2, -1, 6, 1 ) 
796          ELSE
797             CALL lpm_extend_particle_array( trsp_count_recv )
798          ENDIF
799       ENDIF
800
801       CALL MPI_SENDRECV( trnp(1)%age, trnp_count, mpi_particle_type,      &
802                          pnorth, 1, particles(number_of_particles+1)%age, &
803                          trsp_count_recv, mpi_particle_type, psouth, 1,   &
804                          comm2d, status, ierr )
805
806       IF ( use_particle_tails )  THEN
807
808          CALL MPI_SENDRECV( trnpt_count,      1, MPI_INTEGER, pnorth, 0, &
809                             trspt_count_recv, 1, MPI_INTEGER, psouth, 0, &
810                             comm2d, status, ierr )
811
812          IF ( number_of_tails+trspt_count_recv > maximum_number_of_tails ) &
813          THEN
814             IF ( netcdf_output  .AND.  netcdf_data_format < 3 )  THEN
815                message_string = 'maximum_number_of_tails ' //   &
816                                 'needs to be increased ' //     &
817                                 '&but this is not allowed wi'// &
818                                 'th NetCDF output switched on'
819                CALL message( 'lpm_exch_horiz', 'PA0147', 2, 2, -1, 6, 1 )
820             ELSE
821                CALL lpm_extend_tail_array( trspt_count_recv )
822             ENDIF
823          ENDIF
824
825          CALL MPI_SENDRECV( trnpt(1,1,1), trnpt_count*tlength, MPI_REAL,      &
826                             pnorth, 1,                                        &
827                             particle_tail_coordinates(1,1,number_of_tails+1), &
828                             trspt_count_recv*tlength, MPI_REAL, psouth, 1,    &
829                             comm2d, status, ierr )
830!
831!--       Update the tail ids for the transferred particles
832          nn = number_of_tails
833          DO  n = number_of_particles+1, number_of_particles+trsp_count_recv
834             IF ( particles(n)%tail_id /= 0 )  THEN
835                nn = nn + 1
836                particles(n)%tail_id = nn
837             ENDIF
838          ENDDO
839
840       ENDIF
841
842       number_of_particles = number_of_particles + trsp_count_recv
843       number_of_tails     = number_of_tails     + trspt_count_recv
844
845       IF ( use_particle_tails )  THEN
846          DEALLOCATE( trspt, trnpt )
847       ENDIF
848       DEALLOCATE( trsp, trnp )
849
850       CALL cpu_log( log_point_s(23), 'sendrcv_particles', 'stop' )
851
852    ENDIF
853
854#else
855
856!
857!-- Apply boundary conditions
858    DO  n = 1, number_of_particles
859
860       nn = particles(n)%tail_id
861
862       IF ( particles(n)%x < -0.5 * dx )  THEN
863
864          IF ( ibc_par_lr == 0 )  THEN
865!
866!--          Cyclic boundary. Relevant coordinate has to be changed.
867             particles(n)%x = ( nx + 1 ) * dx + particles(n)%x
868             IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
869                i = particles(n)%tailpoints
870                particle_tail_coordinates(1:i,1,nn) = ( nx + 1 ) * dx + &
871                                             particle_tail_coordinates(1:i,1,nn)
872             ENDIF
873          ELSEIF ( ibc_par_lr == 1 )  THEN
874!
875!--          Particle absorption
876             particle_mask(n) = .FALSE.
877             deleted_particles = deleted_particles + 1
878             IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
879                tail_mask(nn) = .FALSE.
880                deleted_tails = deleted_tails + 1
881             ENDIF
882          ELSEIF ( ibc_par_lr == 2 )  THEN
883!
884!--          Particle reflection
885             particles(n)%x       = -dx - particles(n)%x
886             particles(n)%speed_x = -particles(n)%speed_x
887          ENDIF
888
889       ELSEIF ( particles(n)%x >= ( nx + 0.5 ) * dx )  THEN
890
891          IF ( ibc_par_lr == 0 )  THEN
892!
893!--          Cyclic boundary. Relevant coordinate has to be changed.
894             particles(n)%x = particles(n)%x - ( nx + 1 ) * dx
895             IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
896                i = particles(n)%tailpoints
897                particle_tail_coordinates(1:i,1,nn) = - ( nx + 1 ) * dx + &
898                                             particle_tail_coordinates(1:i,1,nn)
899             ENDIF
900          ELSEIF ( ibc_par_lr == 1 )  THEN
901!
902!--          Particle absorption
903             particle_mask(n) = .FALSE.
904             deleted_particles = deleted_particles + 1
905             IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
906                tail_mask(nn) = .FALSE.
907                deleted_tails = deleted_tails + 1
908             ENDIF
909          ELSEIF ( ibc_par_lr == 2 )  THEN
910!
911!--          Particle reflection
912             particles(n)%x       = ( nx + 1 ) * dx - particles(n)%x
913             particles(n)%speed_x = -particles(n)%speed_x
914          ENDIF
915
916       ENDIF
917
918       IF ( particles(n)%y < -0.5 * dy )  THEN
919
920          IF ( ibc_par_ns == 0 )  THEN
921!
922!--          Cyclic boundary. Relevant coordinate has to be changed.
923             particles(n)%y = ( ny + 1 ) * dy + particles(n)%y
924             IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
925                i = particles(n)%tailpoints
926                particle_tail_coordinates(1:i,2,nn) = ( ny + 1 ) * dy + &
927                                             particle_tail_coordinates(1:i,2,nn)
928             ENDIF
929          ELSEIF ( ibc_par_ns == 1 )  THEN
930!
931!--          Particle absorption
932             particle_mask(n) = .FALSE.
933             deleted_particles = deleted_particles + 1
934             IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
935                tail_mask(nn) = .FALSE.
936                deleted_tails = deleted_tails + 1
937             ENDIF
938          ELSEIF ( ibc_par_ns == 2 )  THEN
939!
940!--          Particle reflection
941             particles(n)%y       = -dy - particles(n)%y
942             particles(n)%speed_y = -particles(n)%speed_y
943          ENDIF
944
945       ELSEIF ( particles(n)%y >= ( ny + 0.5 ) * dy )  THEN
946
947          IF ( ibc_par_ns == 0 )  THEN
948!
949!--          Cyclic boundary. Relevant coordinate has to be changed.
950             particles(n)%y = particles(n)%y - ( ny + 1 ) * dy
951             IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
952                i = particles(n)%tailpoints
953                particle_tail_coordinates(1:i,2,nn) = - ( ny + 1 ) * dy + &
954                                             particle_tail_coordinates(1:i,2,nn)
955             ENDIF
956          ELSEIF ( ibc_par_ns == 1 )  THEN
957!
958!--          Particle absorption
959             particle_mask(n) = .FALSE.
960             deleted_particles = deleted_particles + 1
961             IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
962                tail_mask(nn) = .FALSE.
963                deleted_tails = deleted_tails + 1
964             ENDIF
965          ELSEIF ( ibc_par_ns == 2 )  THEN
966!
967!--          Particle reflection
968             particles(n)%y       = ( ny + 1 ) * dy - particles(n)%y
969             particles(n)%speed_y = -particles(n)%speed_y
970          ENDIF
971
972       ENDIF
973    ENDDO
974
975#endif
976
977!
978!-- Accumulate the number of particles transferred between the subdomains
979#if defined( __parallel )
980    trlp_count_sum      = trlp_count_sum      + trlp_count
981    trlp_count_recv_sum = trlp_count_recv_sum + trlp_count_recv
982    trrp_count_sum      = trrp_count_sum      + trrp_count
983    trrp_count_recv_sum = trrp_count_recv_sum + trrp_count_recv
984    trsp_count_sum      = trsp_count_sum      + trsp_count
985    trsp_count_recv_sum = trsp_count_recv_sum + trsp_count_recv
986    trnp_count_sum      = trnp_count_sum      + trnp_count
987    trnp_count_recv_sum = trnp_count_recv_sum + trnp_count_recv
988#endif
989
990
991 END SUBROUTINE lpm_exchange_horiz
Note: See TracBrowser for help on using the repository browser.