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

Last change on this file since 1332 was 1329, checked in by raasch, 11 years ago

last commit documented

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