Ignore:
Timestamp:
Apr 11, 2014 5:15:14 PM (8 years ago)
Author:
hoffmann
Message:

new Lagrangian particle structure integrated

File:
1 edited

Legend:

Unmodified
Added
Removed
  • palm/trunk/SOURCE/lpm_exchange_horiz.f90

    r1329 r1359  
    1  SUBROUTINE lpm_exchange_horiz
     1 MODULE lpm_exchange_horiz_mod
    22
    33!--------------------------------------------------------------------------------!
     
    2020! Current revisions:
    2121! ------------------
    22 !
     22! New particle structure integrated.
     23! Kind definition added to all floating point numbers.
    2324!
    2425! Former revisions:
     
    5657
    5758    USE control_parameters,                                                    &
    58         ONLY:  message_string, netcdf_data_format
     59        ONLY:  dz, message_string, netcdf_data_format, simulated_time
    5960
    6061    USE cpulog,                                                                &
     
    6566
    6667    USE indices,                                                               &
    67         ONLY:  nx, nxl, nxr, ny, nyn, nys
     68        ONLY:  nx, nxl, nxr, ny, nyn, nys, nzb, nzt
    6869
    6970    USE kinds
    7071
     72    USE lpm_pack_arrays_mod,                                                   &
     73        ONLY:  lpm_pack_arrays
     74
    7175    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
     76        ONLY:  alloc_factor, deleted_particles, deleted_tails, grid_particles, &
     77               ibc_par_lr, ibc_par_ns, maximum_number_of_tails,                &
     78               maximum_number_of_tailpoints, min_nr_particle,                  &
     79               mpi_particle_type, number_of_tails, number_of_particles,        &
     80               offset_ocean_nzt, offset_ocean_nzt_m1, particles,               &
     81               particle_tail_coordinates, particle_type, prt_count,            &
     82               tail_mask, trlp_count_sum,                                      &
     83               trlp_count_recv_sum, trnp_count_sum, trnp_count_recv_sum,       &
     84               trrp_count_sum, trrp_count_recv_sum, trsp_count_sum,            &
     85               trsp_count_recv_sum, use_particle_tails, zero_particle
    8086
    8187    USE pegrid
     
    8389    IMPLICIT NONE
    8490
     91    INTEGER(iwp), PARAMETER ::  NR_2_direction_move = 10000 !:
     92    INTEGER(iwp)            ::  nr_move_north               !:
     93    INTEGER(iwp)            ::  nr_move_south               !:
     94
     95    TYPE(particle_type), DIMENSION(NR_2_direction_move) ::  move_also_north
     96    TYPE(particle_type), DIMENSION(NR_2_direction_move) ::  move_also_south
     97
     98    SAVE
     99
     100    PRIVATE
     101    PUBLIC lpm_exchange_horiz, lpm_move_particle, realloc_particles_array
     102
     103    INTERFACE lpm_exchange_horiz
     104       MODULE PROCEDURE lpm_exchange_horiz
     105    END INTERFACE lpm_exchange_horiz
     106
     107    INTERFACE lpm_move_particle
     108       MODULE PROCEDURE lpm_move_particle
     109    END INTERFACE lpm_move_particle
     110
     111    INTERFACE realloc_particles_array
     112       MODULE PROCEDURE realloc_particles_array
     113    END INTERFACE realloc_particles_array
     114
     115CONTAINS
     116
     117 SUBROUTINE lpm_exchange_horiz
     118
     119    IMPLICIT NONE
     120
    85121    INTEGER(iwp) ::  i                                       !:
     122    INTEGER(iwp) ::  ip                                      !:
    86123    INTEGER(iwp) ::  j                                       !:
     124    INTEGER(iwp) ::  jp                                      !:
     125    INTEGER(iwp) ::  k                                       !:
     126    INTEGER(iwp) ::  kp                                      !:
    87127    INTEGER(iwp) ::  n                                       !:
    88128    INTEGER(iwp) ::  nn                                      !:
     
    110150    REAL(wp), DIMENSION(:,:,:), ALLOCATABLE ::  trspt        !:
    111151
    112 
     152    TYPE(particle_type), DIMENSION(:), ALLOCATABLE ::  rvlp  !:
     153    TYPE(particle_type), DIMENSION(:), ALLOCATABLE ::  rvnp  !:
     154    TYPE(particle_type), DIMENSION(:), ALLOCATABLE ::  rvrp  !:
     155    TYPE(particle_type), DIMENSION(:), ALLOCATABLE ::  rvsp  !:
    113156    TYPE(particle_type), DIMENSION(:), ALLOCATABLE ::  trlp  !:
    114157    TYPE(particle_type), DIMENSION(:), ALLOCATABLE ::  trnp  !:
    115158    TYPE(particle_type), DIMENSION(:), ALLOCATABLE ::  trrp  !:
    116159    TYPE(particle_type), DIMENSION(:), ALLOCATABLE ::  trsp  !:
     160
     161    CALL cpu_log( log_point_s(23), 'lpm_exchange_horiz', 'start' )
    117162
    118163#if defined( __parallel )
     
    141186    IF ( pdims(1) /= 1 )  THEN
    142187!
    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
     188!--    First calculate the storage necessary for sending and receiving the data.
     189!--    Compute only first (nxl) and last (nxr) loop iterration.
     190       DO  ip = nxl, nxr, nxr - nxl
     191          DO  jp = nys, nyn
     192             DO  kp = nzb+1, nzt
     193
     194                number_of_particles = prt_count(kp,jp,ip)
     195                IF ( number_of_particles <= 0 )  CYCLE
     196                particles => grid_particles(kp,jp,ip)%particles(1:number_of_particles)
     197                DO  n = 1, number_of_particles
     198                   IF ( particles(n)%particle_mask )  THEN
     199                      i = ( particles(n)%x + 0.5_wp * dx ) * ddx
     200   !
     201   !--                Above calculation does not work for indices less than zero
     202                      IF ( particles(n)%x < -0.5_wp * dx )  i = -1
     203
     204                      IF ( i < nxl )  THEN
     205                         trlp_count = trlp_count + 1
     206                         IF ( particles(n)%tail_id /= 0 )  trlpt_count = trlpt_count + 1
     207                      ELSEIF ( i > nxr )  THEN
     208                         trrp_count = trrp_count + 1
     209                         IF ( particles(n)%tail_id /= 0 )  trrpt_count = trrpt_count + 1
     210                      ENDIF
     211                   ENDIF
     212                ENDDO
     213
     214             ENDDO
     215          ENDDO
    157216       ENDDO
    158217
     
    164223       ALLOCATE( trlp(trlp_count), trrp(trrp_count) )
    165224
    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 )
     225       trlp = zero_particle
     226       trrp = zero_particle
    172227
    173228       IF ( use_particle_tails )  THEN
     
    183238
    184239    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
     240!
     241!-- Compute only first (nxl) and last (nxr) loop iterration
     242    DO  ip = nxl, nxr,nxr-nxl
     243       DO  jp = nys, nyn
     244          DO  kp = nzb+1, nzt
     245             number_of_particles = prt_count(kp,jp,ip)
     246             IF ( number_of_particles <= 0 ) CYCLE
     247             particles => grid_particles(kp,jp,ip)%particles(1:number_of_particles)
     248             DO  n = 1, number_of_particles
     249
     250                nn = particles(n)%tail_id
     251!
     252!--             Only those particles that have not been marked as 'deleted' may
     253!--             be moved.
     254                IF ( particles(n)%particle_mask )  THEN
     255
     256                   i = ( particles(n)%x + 0.5_wp * dx ) * ddx
     257   !
     258   !--             Above calculation does not work for indices less than zero
     259                   IF ( particles(n)%x < - 0.5_wp * dx )  i = -1
     260
     261                   IF ( i <  nxl )  THEN
     262                      IF ( i < 0 )  THEN
     263   !
     264   !--                   Apply boundary condition along x
     265                         IF ( ibc_par_lr == 0 )  THEN
     266   !
     267   !--                      Cyclic condition
     268                            IF ( pdims(1) == 1 )  THEN
     269                               particles(n)%x        = ( nx + 1 ) * dx + particles(n)%x
     270                               particles(n)%origin_x = ( nx + 1 ) * dx + &
     271                               particles(n)%origin_x
     272                               IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
     273                                  i  = particles(n)%tailpoints
     274                                  particle_tail_coordinates(1:i,1,nn) = ( nx + 1 ) * dx &
     275                                  + particle_tail_coordinates(1:i,1,nn)
     276                               ENDIF
     277                            ELSE
     278                               trlp_count = trlp_count + 1
     279                               trlp(trlp_count)   = particles(n)
     280                               trlp(trlp_count)%x = ( nx + 1 ) * dx + trlp(trlp_count)%x
     281                               trlp(trlp_count)%origin_x = trlp(trlp_count)%origin_x + &
     282                               ( nx + 1 ) * dx
     283                               particles(n)%particle_mask  = .FALSE.
     284                               deleted_particles = deleted_particles + 1
     285
     286                               IF ( trlp(trlp_count)%x >= (nx + 0.5_wp)* dx - 1.0E-12_wp )  THEN
     287                                  trlp(trlp_count)%x = trlp(trlp_count)%x - 1.0E-10_wp
     288                                  !++ why is 1 subtracted in next statement???
     289                                  trlp(trlp_count)%origin_x = trlp(trlp_count)%origin_x - 1
     290                               ENDIF
     291
     292                               IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
     293                                  trlpt_count = trlpt_count + 1
     294                                  trlpt(:,:,trlpt_count) = particle_tail_coordinates(:,:,nn)
     295                                  trlpt(:,1,trlpt_count) = ( nx + 1 ) * dx + &
     296                                  trlpt(:,1,trlpt_count)
     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                            particles(n)%particle_mask = .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       = -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
     322   !--                   send to the neighbouring PE
     323                         trlp_count = trlp_count + 1
     324                         trlp(trlp_count) = particles(n)
     325                         particles(n)%particle_mask = .FALSE.
     326                         deleted_particles = deleted_particles + 1
     327
     328                         IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
     329                            trlpt_count = trlpt_count + 1
     330                            trlpt(:,:,trlpt_count) = particle_tail_coordinates(:,:,nn)
     331                            tail_mask(nn) = .FALSE.
     332                            deleted_tails = deleted_tails + 1
     333                         ENDIF
     334                      ENDIF
     335
     336                   ELSEIF ( i > nxr )  THEN
     337                      IF ( i > nx )  THEN
     338   !
     339   !--                   Apply boundary condition along x
     340                         IF ( ibc_par_lr == 0 )  THEN
     341   !
     342   !--                      Cyclic condition
     343                            IF ( pdims(1) == 1 )  THEN
     344                               particles(n)%x = particles(n)%x - ( nx + 1 ) * dx
     345                               particles(n)%origin_x = particles(n)%origin_x - &
     346                               ( nx + 1 ) * dx
     347                               IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
     348                                  i = particles(n)%tailpoints
     349                                  particle_tail_coordinates(1:i,1,nn) = - ( nx+1 ) * dx &
     350                                  + particle_tail_coordinates(1:i,1,nn)
     351                               ENDIF
     352                            ELSE
     353                               trrp_count = trrp_count + 1
     354                               trrp(trrp_count) = particles(n)
     355                               trrp(trrp_count)%x = trrp(trrp_count)%x - ( nx + 1 ) * dx
     356                               trrp(trrp_count)%origin_x = trrp(trrp_count)%origin_x - &
     357                               ( nx + 1 ) * dx
     358                               particles(n)%particle_mask = .FALSE.
     359                               deleted_particles = deleted_particles + 1
     360
     361                               IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
     362                                  trrpt_count = trrpt_count + 1
     363                                  trrpt(:,:,trrpt_count) = particle_tail_coordinates(:,:,nn)
     364                                  trrpt(:,1,trrpt_count) = trrpt(:,1,trrpt_count) - &
     365                                  ( nx + 1 ) * dx
     366                                  tail_mask(nn) = .FALSE.
     367                                  deleted_tails = deleted_tails + 1
     368                               ENDIF
     369                            ENDIF
     370
     371                         ELSEIF ( ibc_par_lr == 1 )  THEN
     372   !
     373   !--                      Particle absorption
     374                            particles(n)%particle_mask = .FALSE.
     375                            deleted_particles = deleted_particles + 1
     376                            IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
     377                               tail_mask(nn) = .FALSE.
     378                               deleted_tails = deleted_tails + 1
     379                            ENDIF
     380
     381                         ELSEIF ( ibc_par_lr == 2 )  THEN
     382   !
     383   !--                      Particle reflection
     384                            particles(n)%x       = 2 * ( nx * dx ) - particles(n)%x
     385                            particles(n)%speed_x = -particles(n)%speed_x
     386
     387                         ENDIF
     388                      ELSE
     389   !
     390   !--                   Store particle data in the transfer array, which will be send
     391   !--                   to the neighbouring PE
     392                         trrp_count = trrp_count + 1
     393                         trrp(trrp_count) = particles(n)
     394                         particles(n)%particle_mask = .FALSE.
     395                         deleted_particles = deleted_particles + 1
     396
     397                         IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
     398                            trrpt_count = trrpt_count + 1
     399                            trrpt(:,:,trrpt_count) = particle_tail_coordinates(:,:,nn)
     400                            tail_mask(nn) = .FALSE.
     401                            deleted_tails = deleted_tails + 1
     402                         ENDIF
     403                      ENDIF
     404
    233405                   ENDIF
    234406                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
     407             ENDDO
     408          ENDDO
     409       ENDDO
    340410    ENDDO
    341411
     
    345415    IF ( pdims(1) /= 1 )  THEN
    346416
    347        CALL cpu_log( log_point_s(23), 'sendrcv_particles', 'start' )
    348417       CALL MPI_SENDRECV( trlp_count,      1, MPI_INTEGER, pleft,  0, &
    349418                          trrp_count_recv, 1, MPI_INTEGER, pright, 0, &
    350419                          comm2d, status, ierr )
    351420
    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,  &
     421       ALLOCATE(rvrp(MAX(1,trrp_count_recv)))
     422
     423       CALL MPI_SENDRECV( trlp(1)%radius, max(1,trlp_count), mpi_particle_type,&
     424                          pleft, 1, rvrp(1)%radius,                            &
     425                          max(1,trrp_count_recv), mpi_particle_type, pright, 1,&
    369426                          comm2d, status, ierr )
     427
     428       IF ( trrp_count_recv > 0 )  CALL Add_particles_to_gridcell(rvrp(1:trrp_count_recv))
     429
     430       DEALLOCATE(rvrp)
    370431
    371432       IF ( use_particle_tails )  THEN
     
    405466       ENDIF
    406467
    407        number_of_particles = number_of_particles + trrp_count_recv
    408        number_of_tails     = number_of_tails     + trrpt_count_recv
    409 
    410468!
    411469!--    Send right boundary, receive left boundary
     
    414472                          comm2d, status, ierr )
    415473
    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,    &
     474       ALLOCATE(rvlp(MAX(1,trlp_count_recv)))
     475
     476       CALL MPI_SENDRECV( trrp(1)%radius, max(1,trrp_count), mpi_particle_type,&
     477                          pright, 1, rvlp(1)%radius,                           &
     478                          max(1,trlp_count_recv), mpi_particle_type, pleft, 1, &
    433479                          comm2d, status, ierr )
     480
     481       IF ( trlp_count_recv > 0 )  CALL Add_particles_to_gridcell(rvlp(1:trlp_count_recv))
     482
     483       DEALLOCATE(rvlp)
    434484
    435485       IF ( use_particle_tails )  THEN
     
    469519       ENDIF
    470520
    471        number_of_particles = number_of_particles + trlp_count_recv
    472        number_of_tails     = number_of_tails     + trlpt_count_recv
     521!       number_of_particles = number_of_particles + trlp_count_recv
     522!       number_of_tails     = number_of_tails     + trlpt_count_recv
    473523
    474524       IF ( use_particle_tails )  THEN
    475525          DEALLOCATE( trlpt, trrpt )
    476526       ENDIF
    477        DEALLOCATE( trlp, trrp )
    478 
    479        CALL cpu_log( log_point_s(23), 'sendrcv_particles', 'pause' )
     527       DEALLOCATE( trlp, trrp )
    480528
    481529    ENDIF
     
    490538!-- For a one-dimensional decomposition along x, no transfer is necessary,
    491539!-- because the particle remains on the PE.
    492     trsp_count  = 0
     540    trsp_count  = nr_move_south
    493541    trspt_count = 0
    494     trnp_count  = 0
     542    trnp_count  = nr_move_north
    495543    trnpt_count = 0
    496544
     
    504552!--    First calculate the storage necessary for sending and receiving the
    505553!--    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
     554       DO  ip = nxl, nxr
     555          DO  jp = nys, nyn, nyn-nys                                 !compute only first (nys) and last (nyn) loop iterration
     556             DO  kp = nzb+1, nzt
     557                number_of_particles = prt_count(kp,jp,ip)
     558                IF ( number_of_particles <= 0 )  CYCLE
     559                particles => grid_particles(kp,jp,ip)%particles(1:number_of_particles)
     560                DO  n = 1, number_of_particles
     561                   IF ( particles(n)%particle_mask )  THEN
     562                      j = ( particles(n)%y + 0.5_wp * dy ) * ddy
     563!
     564!--                   Above calculation does not work for indices less than zero
     565                      IF ( particles(n)%y < -0.5_wp * dy )  j = -1
     566
     567                      IF ( j < nys )  THEN
     568                         trsp_count = trsp_count + 1
     569                         IF ( particles(n)%tail_id /= 0 )  trspt_count = trspt_count + 1
     570                      ELSEIF ( j > nyn )  THEN
     571                         trnp_count = trnp_count + 1
     572                         IF ( particles(n)%tail_id /= 0 )  trnpt_count = trnpt_count + 1
     573                      ENDIF
     574                   ENDIF
     575                ENDDO
     576             ENDDO
     577          ENDDO
    521578       ENDDO
    522579
     
    528585       ALLOCATE( trsp(trsp_count), trnp(trnp_count) )
    529586
    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 )
     587       trsp = zero_particle
     588       trnp = zero_particle
    536589
    537590       IF ( use_particle_tails )  THEN
     
    541594       ENDIF
    542595
    543        trsp_count  = 0
     596       trsp_count  = nr_move_south
    544597       trspt_count = 0
    545        trnp_count  = 0
     598       trnp_count  = nr_move_north
    546599       trnpt_count = 0
    547600
     601       trsp(1:nr_move_south) = move_also_south(1:nr_move_south)
     602       trnp(1:nr_move_north) = move_also_north(1:nr_move_north)
     603
    548604    ENDIF
    549605
    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)
     606    DO  ip = nxl, nxr
     607       DO  jp = nys, nyn, nyn-nys ! compute only first (nys) and last (nyn) loop iterration
     608          DO  kp = nzb+1, nzt
     609             number_of_particles = prt_count(kp,jp,ip)
     610             IF ( number_of_particles <= 0 )  CYCLE
     611             particles => grid_particles(kp,jp,ip)%particles(1:number_of_particles)
     612             DO  n = 1, number_of_particles
     613
     614                nn = particles(n)%tail_id
     615!
     616!--             Only those particles that have not been marked as 'deleted' may
     617!--             be moved.
     618                IF ( particles(n)%particle_mask )  THEN
     619                   j = ( particles(n)%y + 0.5_wp * dy ) * ddy
     620!
     621!--                Above calculation does not work for indices less than zero
     622                   IF ( particles(n)%y < -0.5_wp * dy )  j = -1
     623
     624                   IF ( j < nys )  THEN
     625                      IF ( j < 0 )  THEN
     626!
     627!--                      Apply boundary condition along y
     628                         IF ( ibc_par_ns == 0 )  THEN
     629!
     630!--                         Cyclic condition
     631                            IF ( pdims(2) == 1 )  THEN
     632                               particles(n)%y = ( ny + 1 ) * dy + particles(n)%y
     633                               particles(n)%origin_y = ( ny + 1 ) * dy + &
     634                               particles(n)%origin_y
     635                               IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
     636                                  i = particles(n)%tailpoints
     637                                  particle_tail_coordinates(1:i,2,nn) =        &
     638                                     ( ny+1 ) * dy + particle_tail_coordinates(1:i,2,nn)
     639                               ENDIF
     640                            ELSE
     641                               trsp_count = trsp_count + 1
     642                               trsp(trsp_count) = particles(n)
     643                               trsp(trsp_count)%y = ( ny + 1 ) * dy + &
     644                               trsp(trsp_count)%y
     645                               trsp(trsp_count)%origin_y = trsp(trsp_count)%origin_y &
     646                               + ( ny + 1 ) * dy
     647                               particles(n)%particle_mask = .FALSE.
     648                               deleted_particles = deleted_particles + 1
     649
     650                               IF ( trsp(trsp_count)%y >= (ny+0.5_wp)* dy - 1.0E-12_wp )  THEN
     651                                  trsp(trsp_count)%y = trsp(trsp_count)%y - 1.0E-10_wp
     652                                  !++ why is 1 subtracted in next statement???
     653                                  trsp(trsp_count)%origin_y =                        &
     654                                  trsp(trsp_count)%origin_y - 1
     655                               ENDIF
     656
     657                               IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
     658                                  trspt_count = trspt_count + 1
     659                                  trspt(:,:,trspt_count) = &
     660                                  particle_tail_coordinates(:,:,nn)
     661                                  trspt(:,2,trspt_count) = ( ny + 1 ) * dy + &
     662                                  trspt(:,2,trspt_count)
     663                                  tail_mask(nn) = .FALSE.
     664                                  deleted_tails = deleted_tails + 1
     665                               ENDIF
     666                            ENDIF
     667
     668                         ELSEIF ( ibc_par_ns == 1 )  THEN
     669!
     670!--                         Particle absorption
     671                            particles(n)%particle_mask = .FALSE.
     672                            deleted_particles = deleted_particles + 1
     673                            IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
     674                               tail_mask(nn) = .FALSE.
     675                               deleted_tails = deleted_tails + 1
     676                            ENDIF
     677
     678                         ELSEIF ( ibc_par_ns == 2 )  THEN
     679!
     680!--                         Particle reflection
     681                            particles(n)%y       = -particles(n)%y
     682                            particles(n)%speed_y = -particles(n)%speed_y
     683
     684                         ENDIF
     685                      ELSE
     686!
     687!--                      Store particle data in the transfer array, which will
     688!--                      be send to the neighbouring PE
     689                         trsp_count = trsp_count + 1
     690                         trsp(trsp_count) = particles(n)
     691                         particles(n)%particle_mask = .FALSE.
     692                         deleted_particles = deleted_particles + 1
     693
     694                         IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
     695                            trspt_count = trspt_count + 1
     696                            trspt(:,:,trspt_count) = particle_tail_coordinates(:,:,nn)
     697                            tail_mask(nn) = .FALSE.
     698                            deleted_tails = deleted_tails + 1
     699                         ENDIF
    577700                      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
     701
     702                   ELSEIF ( j > nyn )  THEN
     703                      IF ( j > ny )  THEN
     704!
     705!--                       Apply boundary condition along x
     706                         IF ( ibc_par_ns == 0 )  THEN
     707!
     708!--                         Cyclic condition
     709                            IF ( pdims(2) == 1 )  THEN
     710                               particles(n)%y = particles(n)%y - ( ny + 1 ) * dy
     711                               particles(n)%origin_y =                         &
     712                                  particles(n)%origin_y - ( ny + 1 ) * dy
     713                               IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
     714                                  i = particles(n)%tailpoints
     715                                  particle_tail_coordinates(1:i,2,nn) =        &
     716                                     - (ny+1) * dy + particle_tail_coordinates(1:i,2,nn)
     717                               ENDIF
     718                            ELSE
     719                               trnp_count = trnp_count + 1
     720                               trnp(trnp_count) = particles(n)
     721                               trnp(trnp_count)%y =                            &
     722                                  trnp(trnp_count)%y - ( ny + 1 ) * dy
     723                               trnp(trnp_count)%origin_y =                     &
     724                                  trnp(trnp_count)%origin_y - ( ny + 1 ) * dy
     725                               particles(n)%particle_mask = .FALSE.
     726                               deleted_particles = deleted_particles + 1
     727                               IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
     728                                  trnpt_count = trnpt_count + 1
     729                                  trnpt(:,:,trnpt_count) =                     &
     730                                     particle_tail_coordinates(:,:,nn)
     731                                  trnpt(:,2,trnpt_count) =                     &
     732                                     trnpt(:,2,trnpt_count) - ( ny + 1 ) * dy
     733                                  tail_mask(nn) = .FALSE.
     734                                  deleted_tails = deleted_tails + 1
     735                               ENDIF
     736                            ENDIF
     737
     738                         ELSEIF ( ibc_par_ns == 1 )  THEN
     739!
     740!--                         Particle absorption
     741                            particles(n)%particle_mask = .FALSE.
     742                            deleted_particles = deleted_particles + 1
     743                            IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
     744                               tail_mask(nn) = .FALSE.
     745                               deleted_tails = deleted_tails + 1
     746                            ENDIF
     747
     748                         ELSEIF ( ibc_par_ns == 2 )  THEN
     749!
     750!--                         Particle reflection
     751                            particles(n)%y       = 2 * ( ny * dy ) - particles(n)%y
     752                            particles(n)%speed_y = -particles(n)%speed_y
     753
     754                         ENDIF
     755                      ELSE
     756!
     757!--                      Store particle data in the transfer array, which will
     758!--                      be send to the neighbouring PE
     759                         trnp_count = trnp_count + 1
     760                         trnp(trnp_count) = particles(n)
     761                         particles(n)%particle_mask = .FALSE.
     762                         deleted_particles = deleted_particles + 1
     763
     764                         IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
     765                            trnpt_count = trnpt_count + 1
     766                            trnpt(:,:,trnpt_count) = particle_tail_coordinates(:,:,nn)
     767                            tail_mask(nn) = .FALSE.
     768                            deleted_tails = deleted_tails + 1
     769                         ENDIF
    593770                      ENDIF
    594771
    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
    604772                   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 
    622773                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
     774             ENDDO
     775          ENDDO
     776       ENDDO
    713777    ENDDO
    714778
     
    718782    IF ( pdims(2) /= 1 )  THEN
    719783
    720        CALL cpu_log( log_point_s(23), 'sendrcv_particles', 'continue' )
    721784       CALL MPI_SENDRECV( trsp_count,      1, MPI_INTEGER, psouth, 0, &
    722785                          trnp_count_recv, 1, MPI_INTEGER, pnorth, 0, &
    723786                          comm2d, status, ierr )
    724787
    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, &
     788       ALLOCATE(rvnp(MAX(1,trnp_count_recv)))
     789
     790       CALL MPI_SENDRECV( trsp(1)%radius, trsp_count, mpi_particle_type,      &
     791                          psouth, 1, rvnp(1)%radius,                             &
    741792                          trnp_count_recv, mpi_particle_type, pnorth, 1,   &
    742793                          comm2d, status, ierr )
     794
     795       IF ( trnp_count_recv  > 0 )  CALL Add_particles_to_gridcell(rvnp(1:trnp_count_recv))
     796
     797       DEALLOCATE(rvnp)
    743798
    744799       IF ( use_particle_tails )  THEN
     
    779834       ENDIF
    780835
    781        number_of_particles = number_of_particles + trnp_count_recv
    782        number_of_tails     = number_of_tails     + trnpt_count_recv
     836!       number_of_particles = number_of_particles + trnp_count_recv
     837!       number_of_tails     = number_of_tails     + trnpt_count_recv
    783838
    784839!
     
    788843                          comm2d, status, ierr )
    789844
    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, &
     845       ALLOCATE(rvsp(MAX(1,trsp_count_recv)))
     846
     847       CALL MPI_SENDRECV( trnp(1)%radius, trnp_count, mpi_particle_type,      &
     848                          pnorth, 1, rvsp(1)%radius,                          &
    806849                          trsp_count_recv, mpi_particle_type, psouth, 1,   &
    807850                          comm2d, status, ierr )
     851
     852       IF ( trsp_count_recv > 0 )  CALL Add_particles_to_gridcell(rvsp(1:trsp_count_recv))
     853
     854       DEALLOCATE(rvsp)
    808855
    809856       IF ( use_particle_tails )  THEN
     
    851898       DEALLOCATE( trsp, trnp )
    852899
    853        CALL cpu_log( log_point_s(23), 'sendrcv_particles', 'stop' )
    854 
    855900    ENDIF
    856901
     
    863908       nn = particles(n)%tail_id
    864909
    865        IF ( particles(n)%x < -0.5 * dx )  THEN
     910       IF ( particles(n)%x < -0.5_wp * dx )  THEN
    866911
    867912          IF ( ibc_par_lr == 0 )  THEN
     
    877922!
    878923!--          Particle absorption
    879              particle_mask(n) = .FALSE.
     924             particles(n)%particle_mask = .FALSE.
    880925             deleted_particles = deleted_particles + 1
    881926             IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
     
    890935          ENDIF
    891936
    892        ELSEIF ( particles(n)%x >= ( nx + 0.5 ) * dx )  THEN
     937       ELSEIF ( particles(n)%x >= ( nx + 0.5_wp ) * dx )  THEN
    893938
    894939          IF ( ibc_par_lr == 0 )  THEN
     
    904949!
    905950!--          Particle absorption
    906              particle_mask(n) = .FALSE.
     951             particles(n)%particle_mask = .FALSE.
    907952             deleted_particles = deleted_particles + 1
    908953             IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
     
    919964       ENDIF
    920965
    921        IF ( particles(n)%y < -0.5 * dy )  THEN
     966       IF ( particles(n)%y < -0.5_wp * dy )  THEN
    922967
    923968          IF ( ibc_par_ns == 0 )  THEN
     
    933978!
    934979!--          Particle absorption
    935              particle_mask(n) = .FALSE.
     980             particles(n)%particle_mask = .FALSE.
    936981             deleted_particles = deleted_particles + 1
    937982             IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
     
    946991          ENDIF
    947992
    948        ELSEIF ( particles(n)%y >= ( ny + 0.5 ) * dy )  THEN
     993       ELSEIF ( particles(n)%y >= ( ny + 0.5_wp ) * dy )  THEN
    949994
    950995          IF ( ibc_par_ns == 0 )  THEN
     
    9601005!
    9611006!--          Particle absorption
    962              particle_mask(n) = .FALSE.
     1007             particles(n)%particle_mask = .FALSE.
    9631008             deleted_particles = deleted_particles + 1
    9641009             IF ( use_particle_tails  .AND.  nn /= 0 )  THEN
     
    9911036#endif
    9921037
     1038    CALL cpu_log( log_point_s(23), 'lpm_exchange_horiz', 'stop' )
    9931039
    9941040 END SUBROUTINE lpm_exchange_horiz
     1041
     1042 SUBROUTINE Add_particles_to_gridcell (particle_array)
     1043
     1044!
     1045!-- If a particle moves from one processor to another, this subroutine moves
     1046!-- the corresponding elements from the particle arrays of the old grid cells
     1047!-- to the particle arrays of the new grid cells.
     1048
     1049    IMPLICIT NONE
     1050
     1051    INTEGER(iwp)        ::  ip        !:
     1052    INTEGER(iwp)        ::  jp        !:
     1053    INTEGER(iwp)        ::  kp        !:
     1054    INTEGER(iwp)        ::  n         !:
     1055    INTEGER(iwp)        ::  pindex    !:
     1056
     1057    LOGICAL             ::  pack_done !:
     1058
     1059    TYPE(particle_type), DIMENSION(:), INTENT(IN)       :: particle_array
     1060
     1061    pack_done     = .FALSE.
     1062
     1063    nr_move_north = 0
     1064    nr_move_south = 0
     1065
     1066    DO n = 1, SIZE(particle_array)
     1067       ip = ( particle_array(n)%x + 0.5_wp * dx ) * ddx
     1068       jp = ( particle_array(n)%y + 0.5_wp * dy ) * ddy
     1069       kp = particle_array(n)%z / dz + 1 + offset_ocean_nzt_m1
     1070
     1071       IF ( ip >= nxl  .AND.  ip <= nxr  .AND.  jp >= nys  .AND.  jp <= nyn    &
     1072            .AND.  kp >= nzb+1  .AND.  kp <= nzt)  THEN ! particle stays on processor
     1073          number_of_particles = prt_count(kp,jp,ip)
     1074          particles => grid_particles(kp,jp,ip)%particles(1:number_of_particles)
     1075
     1076          pindex = prt_count(kp,jp,ip)+1
     1077          IF( pindex > SIZE(grid_particles(kp,jp,ip)%particles) )  THEN
     1078             IF ( pack_done )  THEN
     1079                CALL realloc_particles_array (ip,jp,kp)
     1080             ELSE
     1081                CALL lpm_pack_arrays
     1082                prt_count(kp,jp,ip) = number_of_particles
     1083                pindex = prt_count(kp,jp,ip)+1
     1084                IF ( pindex > SIZE(grid_particles(kp,jp,ip)%particles) )  THEN
     1085                   CALL realloc_particles_array (ip,jp,kp)
     1086                ENDIF
     1087                pack_done = .TRUE.
     1088             ENDIF
     1089          ENDIF
     1090          grid_particles(kp,jp,ip)%particles(pindex) = particle_array(n)
     1091          prt_count(kp,jp,ip) = pindex
     1092       ELSE
     1093          IF ( jp == nys - 1 )  THEN
     1094             nr_move_south = nr_move_south+1
     1095             move_also_south(nr_move_south) = particle_array(n)
     1096             IF ( jp == -1 )  THEN
     1097                move_also_south(nr_move_south)%y =                             &
     1098                   move_also_south(nr_move_south)%y + ( ny + 1 ) * dy
     1099                move_also_south(nr_move_south)%origin_y =                      &
     1100                   move_also_south(nr_move_south)%origin_y + ( ny + 1 ) * dy
     1101             ENDIF
     1102          ELSEIF ( jp == nyn+1 )  THEN
     1103             nr_move_north = nr_move_north+1
     1104             move_also_north(nr_move_north) = particle_array(n)
     1105             IF ( jp == ny+1 )  THEN
     1106                move_also_north(nr_move_north)%y =                             &
     1107                   move_also_north(nr_move_north)%y - ( ny + 1 ) * dy
     1108                move_also_north(nr_move_north)%origin_y =                      &
     1109                   move_also_north(nr_move_north)%origin_y - ( ny + 1 ) * dy
     1110             ENDIF
     1111          ELSE
     1112             WRITE(0,'(a,8i7)') 'particle out of range ',myid,ip,jp,kp,nxl,nxr,nys,nyn
     1113          ENDIF
     1114       ENDIF
     1115    ENDDO
     1116
     1117    RETURN
     1118
     1119 END SUBROUTINE Add_particles_to_gridcell
     1120
     1121
     1122
     1123
     1124 SUBROUTINE lpm_move_particle
     1125
     1126!
     1127!-- If a particle moves from one grid cell to another (on the current
     1128!-- processor!), this subroutine moves the corresponding element from the
     1129!-- particle array of the old grid cell to the particle array of the new grid
     1130!-- cell.
     1131
     1132    IMPLICIT NONE
     1133
     1134    INTEGER(iwp)        ::  i           !:
     1135    INTEGER(iwp)        ::  ip          !:
     1136    INTEGER(iwp)        ::  j           !:
     1137    INTEGER(iwp)        ::  jp          !:
     1138    INTEGER(iwp)        ::  k           !:
     1139    INTEGER(iwp)        ::  kp          !:
     1140    INTEGER(iwp)        ::  n           !:
     1141    INTEGER(iwp)        ::  np_old_cell !:
     1142    INTEGER(iwp)        ::  n_start     !:
     1143    INTEGER(iwp)        ::  pindex      !:
     1144
     1145    LOGICAL             ::  pack_done   !:
     1146
     1147    TYPE(particle_type), DIMENSION(:), POINTER  ::  particles_old_cell !:
     1148
     1149    CALL cpu_log( log_point_s(41), 'lpm_move_particle', 'start' )
     1150
     1151    DO  ip = nxl, nxr
     1152       DO  jp = nys, nyn
     1153          DO  kp = nzb+1, nzt
     1154
     1155             np_old_cell = prt_count(kp,jp,ip)
     1156             IF ( np_old_cell <= 0 )  CYCLE
     1157             particles_old_cell => grid_particles(kp,jp,ip)%particles(1:np_old_cell)
     1158             n_start = -1
     1159             
     1160             DO  n = 1, np_old_cell
     1161                i = ( particles_old_cell(n)%x + 0.5_wp * dx ) * ddx
     1162                j = ( particles_old_cell(n)%y + 0.5_wp * dy ) * ddy
     1163                k = particles_old_cell(n)%z / dz + 1 + offset_ocean_nzt
     1164!
     1165!--             Check, if particle has moved to another grid cell.
     1166                IF ( i /= ip  .OR.  j /= jp  .OR.  k /= kp )  THEN
     1167!
     1168!--                The particle has moved to another grid cell. Now check, if
     1169!--                particle stays on the same processor.
     1170                   IF ( i >= nxl  .AND.  i <= nxr  .AND.  j >= nys  .AND.      &
     1171                        j <= nyn  .AND.  k >= nzb+1  .AND.  k <= nzt)  THEN
     1172!
     1173!--                   If the particle stays on the same processor, the particle
     1174!--                   will be added to the particle array of the new processor.
     1175                      number_of_particles = prt_count(k,j,i)
     1176                      particles => grid_particles(k,j,i)%particles(1:number_of_particles)
     1177
     1178                      pindex = prt_count(k,j,i)+1
     1179                      IF (  pindex > SIZE(grid_particles(k,j,i)%particles)  )  &
     1180                      THEN
     1181                         n_start = n
     1182                         EXIT
     1183                      ENDIF
     1184
     1185                      grid_particles(k,j,i)%particles(pindex) = particles_old_cell(n)
     1186                      prt_count(k,j,i) = pindex
     1187
     1188                      particles_old_cell(n)%particle_mask = .FALSE.
     1189                   ENDIF
     1190                ENDIF
     1191             ENDDO
     1192
     1193             IF ( n_start .GE. 0 )  THEN
     1194                pack_done = .FALSE.
     1195                DO  n = n_start, np_old_cell
     1196                   i = ( particles_old_cell(n)%x + 0.5_wp * dx ) * ddx
     1197                   j = ( particles_old_cell(n)%y + 0.5_wp * dy ) * ddy
     1198                   k = particles_old_cell(n)%z / dz + 1 + offset_ocean_nzt
     1199                   IF ( i /= ip  .OR.  j /= jp  .OR.  k /= kp )  THEN
     1200!
     1201!--                   Particle is in different box
     1202                      IF ( i >= nxl  .AND.  i <= nxr  .AND.  j >= nys  .AND.   &
     1203                           j <= nyn  .AND.  k >= nzb+1  .AND.  k <= nzt)  THEN
     1204!
     1205!--                      Particle stays on processor
     1206                         number_of_particles = prt_count(k,j,i)
     1207                         particles => grid_particles(k,j,i)%particles(1:number_of_particles)
     1208
     1209                         pindex = prt_count(k,j,i)+1
     1210                         IF ( pindex > SIZE(grid_particles(k,j,i)%particles) ) &
     1211                         THEN
     1212                            IF ( pack_done )  THEN
     1213                               CALL realloc_particles_array(i,j,k)
     1214                               pindex = prt_count(k,j,i)+1
     1215                            ELSE
     1216                               CALL lpm_pack_arrays
     1217                               prt_count(k,j,i) = number_of_particles
     1218
     1219                               pindex = prt_count(k,j,i)+1
     1220!
     1221!--                            If number of particles in the new grid box
     1222!--                            exceeds its allocated memory, the particle array
     1223!--                            will be reallocated
     1224                               IF ( pindex > SIZE(grid_particles(k,j,i)%particles) )  THEN
     1225                                  CALL realloc_particles_array(i,j,k)
     1226                                  pindex = prt_count(k,j,i)+1
     1227                               ENDIF
     1228
     1229                               pack_done = .TRUE.
     1230                            ENDIF
     1231                         ENDIF
     1232
     1233                         grid_particles(k,j,i)%particles(pindex) = particles_old_cell(n)
     1234                         prt_count(k,j,i) = pindex
     1235
     1236                         particles_old_cell(n)%particle_mask = .FALSE.
     1237                      ENDIF
     1238                   ENDIF
     1239                ENDDO
     1240             ENDIF
     1241          ENDDO
     1242       ENDDO
     1243    ENDDO
     1244
     1245    CALL cpu_log( log_point_s(41), 'lpm_move_particle', 'stop' )
     1246
     1247    RETURN
     1248
     1249 END SUBROUTINE lpm_move_particle
     1250
     1251 SUBROUTINE realloc_particles_array (i,j,k,size_in)
     1252
     1253    IMPLICIT NONE
     1254
     1255    INTEGER(iwp), INTENT(in)                       ::  i              !:
     1256    INTEGER(iwp), INTENT(in)                       ::  j              !:
     1257    INTEGER(iwp), INTENT(in)                       ::  k              !:
     1258    INTEGER(iwp), INTENT(in), optional             ::  size_in        !:
     1259
     1260    INTEGER(iwp)                                   :: old_size        !:
     1261    INTEGER(iwp)                                   :: new_size        !:
     1262    TYPE(particle_type), DIMENSION(:), ALLOCATABLE :: tmp_particles_d !:
     1263    TYPE(particle_type), DIMENSION(500)            :: tmp_particles_s !:
     1264
     1265
     1266    old_size = SIZE(grid_particles(k,j,i)%particles)
     1267
     1268    IF ( PRESENT(size_in) )   THEN
     1269       new_size = size_in
     1270    ELSE
     1271       new_size = old_size * ( 1.0 + alloc_factor / 100.0 )
     1272    ENDIF
     1273
     1274    new_size = MAX( new_size, min_nr_particle )
     1275
     1276    IF ( old_size <= 500 )  THEN
     1277
     1278       tmp_particles_s(1:old_size) = grid_particles(k,j,i)%particles(1:old_size)
     1279
     1280       DEALLOCATE(grid_particles(k,j,i)%particles)
     1281       ALLOCATE(grid_particles(k,j,i)%particles(new_size))
     1282
     1283       grid_particles(k,j,i)%particles(1:old_size)          = tmp_particles_s(1:old_size)
     1284       grid_particles(k,j,i)%particles(old_size+1:new_size) = zero_particle
     1285
     1286    ELSE
     1287
     1288       ALLOCATE(tmp_particles_d(new_size))
     1289       tmp_particles_d(1:old_size) = grid_particles(k,j,i)%particles
     1290
     1291       DEALLOCATE(grid_particles(k,j,i)%particles)
     1292       ALLOCATE(grid_particles(k,j,i)%particles(new_size))
     1293
     1294       grid_particles(k,j,i)%particles(1:old_size)          = tmp_particles_d(1:old_size)
     1295       grid_particles(k,j,i)%particles(old_size+1:new_size) = zero_particle
     1296
     1297       DEALLOCATE(tmp_particles_d)
     1298
     1299    ENDIF
     1300    particles => grid_particles(k,j,i)%particles(1:number_of_particles)
     1301
     1302    RETURN
     1303 END SUBROUTINE realloc_particles_array
     1304
     1305END MODULE lpm_exchange_horiz_mod
Note: See TracChangeset for help on using the changeset viewer.