source: palm/trunk/SOURCE/pmc_particle_interface.f90 @ 4828

Last change on this file since 4828 was 4828, checked in by Giersch, 3 years ago

Copyright updated to year 2021, interface pmc_sort removed to accelarate the nesting code

  • Property svn:keywords set to Id
File size: 51.6 KB
RevLine 
[2801]1MODULE pmc_particle_interface
[4649]2!--------------------------------------------------------------------------------------------------!
3! This file is part of the PALM model system.
[2801]4!
[4649]5! PALM is free software: you can redistribute it and/or modify it under the terms of the GNU General
6! Public License as published by the Free Software Foundation, either version 3 of the License, or
7! (at your option) any later version.
[2801]8!
[4649]9! PALM is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the
10! implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General
11! Public License for more details.
[2801]12!
[4649]13! You should have received a copy of the GNU General Public License along with PALM. If not, see
14! <http://www.gnu.org/licenses/>.
[2801]15!
[4828]16! Copyright 1997-2021 Leibniz Universitaet Hannover
[4649]17!--------------------------------------------------------------------------------------------------!
[2801]18!
[4649]19!
[2801]20! Current revisions:
[4649]21! -----------------
[4629]22!
23!
[2801]24! Former revisions:
[4629]25! -----------------
[2809]26! $Id: pmc_particle_interface.f90 4828 2021-01-05 11:21:41Z Giersch $
[4649]27! File re-formatted to follow the PALM coding standard
28!
29!
30! 4629 2020-07-29 09:37:56Z raasch
31! Support for MPI Fortran77 interface (mpif.h) removed
32!
[4629]33! 4444 2020-03-05 15:59:50Z raasch
[4649]34! Bugfix: preprocessor directives for serial mode added
35!
[4444]36! 4360 2020-01-07 11:25:50Z suehring
[4182]37! Corrected "Former revisions" section
[4649]38!
[4182]39! 4043 2019-06-18 16:59:00Z schwenkel
[4043]40! Remove min_nr_particle
[4649]41!
[4043]42! 4017 2019-06-06 12:16:46Z schwenkel
[4649]43! Coarse bound renamed as parent_bound and icl, icr, jcs, jcn as ipl, ipr, jps, jpn.
44!
[3948]45! 3883 2019-04-10 12:51:50Z hellstea
[4649]46! Function get_number_of_childs renamed to get_number_of_children and cg renamed to pg according to
47! their definitions in pmc_interface_mod
48!
[3883]49! 3655 2019-01-07 16:51:22Z knoop
[4649]50! Unused variables removed
[4182]51!
52! Initial Version (by Klaus Ketelsen)
53!
[4649]54!
55!--------------------------------------------------------------------------------------------------!
[2801]56! Description:
57! ------------
[4649]58! Introduce particle transfer in nested models.
59! Interface to palm lpm model to handel particle transfer between parent and child model.
60!--------------------------------------------------------------------------------------------------!
[4444]61#if defined( __parallel )
[2801]62
63   USE, INTRINSIC ::  ISO_C_BINDING
[2841]64
[2801]65   USE MPI
66
67   USE kinds
68
[4649]69   USE pegrid,                                                                                     &
70       ONLY: myidx,                                                                                &
71             myidy
[2801]72
[4649]73   USE indices,                                                                                    &
74       ONLY: nbgp,                                                                                 &
75             nx,                                                                                   &
76             nxl,                                                                                  &
77             nxr,                                                                                  &
78             nxlg,                                                                                 &
79             nxrg,                                                                                 &
80             ny,                                                                                   &
81             nys,                                                                                  &
82             nyn,                                                                                  &
83             nysg,                                                                                 &
84             nyng,                                                                                 &
85             nzb,                                                                                  &
86             nzt
[2801]87
88
[4649]89   USE grid_variables,                                                                             &
90       ONLY:  dx,                                                                                  &
91              dy
92
93   USE arrays_3d,                                                                                  &
[2801]94        ONLY: zw
95
[4649]96   USE control_parameters,                                                                         &
[2801]97       ONLY:  message_string
98
[4649]99   USE particle_attributes,                                                                        &
100       ONLY:  alloc_factor,                                                                        &
101              grid_particles,                                                                      &
102              ibc_par_lr,                                                                          &
103              ibc_par_ns,                                                                          &
104              ibc_par_t,                                                                           &
105              particles,                                                                           &
106              particle_type,                                                                       &
107              prt_count,                                                                           &
108              number_of_particles,                                                                 &
109              zero_particle
[2801]110
[4649]111
112
113
114
[4017]115!    USE lpm_pack_and_sort_mod
[2801]116
[2967]117#if defined( __parallel )
[4649]118   USE pmc_general,                                                                                &
[2801]119       ONLY: pedef
120
[4649]121   USE pmc_parent,                                                                                 &
122       ONLY: children,                                                                             &
123             pmc_s_fillbuffer,                                                                     &
124             pmc_s_getdata_from_buffer,                                                            &
[2801]125             pmc_s_get_child_npes
126
[4649]127   USE pmc_child,                                                                                  &
128       ONLY:  me,                                                                                  &
129              pmc_c_getbuffer,                                                                     &
130              pmc_c_putbuffer
[2801]131
[4649]132   USE pmc_interface,                                                                              &
133       ONLY:  coord_x,                                                                             &
134              coord_y,                                                                             &
135              cpl_id,                                                                              &
136              get_childid,                                                                         &
137              get_child_edges,                                                                     &
138              get_child_gridspacing,                                                               &
139              get_number_of_children,                                                              &
140              lower_left_coord_x,                                                                  &
141              lower_left_coord_y,                                                                  &
142              nr_part,                                                                             &
143              nr_partc,                                                                            &
144              parent_bound,                                                                        &
145              part_adr,                                                                            &
146              part_adrc,                                                                           &
147              pg,                                                                                  &
148              nested_run
[2801]149
[4649]150
151    USE pmc_handle_communicator,                                                                   &
[2801]152        ONLY:  pmc_parent_for_child
153
[4649]154   USE pmc_mpi_wrapper,                                                                            &
155       ONLY:   pmc_recv_from_child,                                                                &
156               pmc_send_to_parent
[2801]157
[2967]158#endif
[2801]159
160   IMPLICIT NONE
[2841]161
162
[2801]163   PRIVATE
164   SAVE
165
166   TYPE  coarse_particle_def
[4649]167      INTEGER(iwp) ::  nr_particle  !<
168
169      TYPE(particle_type),ALLOCATABLE,DIMENSION(:) ::  parent_particles  !<
[2801]170   END TYPE  coarse_particle_def
171
[4649]172   INTEGER(iwp),PARAMETER ::  max_nr_particle_in_rma_win = 100000  !<
173   INTEGER(iwp),PARAMETER ::  min_particles_per_column   = 100     !<
[2801]174
175
[4649]176   INTEGER(iwp) ::  nr_fine_in_coarse   !< Number of fine grid cells in coarse grid (one direction)
177   INTEGER(iwp) ::  particle_win_child  !<
[2801]178
[4649]179   INTEGER(iwp),ALLOCATABLE,DIMENSION(:) ::  particle_win_parent  !<
[2801]180
[4649]181   TYPE(C_PTR), ALLOCATABLE,DIMENSION(:) ::  buf_ptr  !<
[2801]182
[4649]183   TYPE(particle_type), DIMENSION(:),POINTER ::  particle_in_win  !<
184
185   TYPE(coarse_particle_def),ALLOCATABLE,DIMENSION(:,:) ::  coarse_particles  !<
186
187
188
[2801]189   INTERFACE pmcp_g_init
190      MODULE PROCEDURE pmcp_g_init
191   END  INTERFACE pmcp_g_init
192
193   INTERFACE pmcp_g_alloc_win
194      MODULE PROCEDURE pmcp_g_alloc_win
195   END  INTERFACE pmcp_g_alloc_win
196
197   INTERFACE pmcp_c_get_particle_from_parent
198      MODULE PROCEDURE pmcp_c_get_particle_from_parent
199   END  INTERFACE pmcp_c_get_particle_from_parent
200
201   INTERFACE pmcp_c_send_particle_to_parent
202      MODULE PROCEDURE pmcp_c_send_particle_to_parent
203   END  INTERFACE pmcp_c_send_particle_to_parent
204
205   INTERFACE pmcp_p_fill_particle_win
206      MODULE PROCEDURE pmcp_p_fill_particle_win
207   END  INTERFACE pmcp_p_fill_particle_win
208
209   INTERFACE pmcp_p_empty_particle_win
210      MODULE PROCEDURE pmcp_p_empty_particle_win
211   END  INTERFACE pmcp_p_empty_particle_win
212
213   INTERFACE pmcp_g_print_number_of_particles
214      MODULE PROCEDURE pmcp_g_print_number_of_particles
215   END  INTERFACE pmcp_g_print_number_of_particles
216
217   INTERFACE pmcp_p_delete_particles_in_fine_grid_area
218      MODULE PROCEDURE pmcp_p_delete_particles_in_fine_grid_area
219   END  INTERFACE pmcp_p_delete_particles_in_fine_grid_area
220
221   PUBLIC pmcp_g_init, pmcp_g_alloc_win, pmcp_c_get_particle_from_parent
222   PUBLIC pmcp_c_send_particle_to_parent, pmcp_p_fill_particle_win, pmcp_g_print_number_of_particles
223   PUBLIC pmcp_p_empty_particle_win, pmcp_p_delete_particles_in_fine_grid_area
224
225 CONTAINS
226
[4649]227!--------------------------------------------------------------------------------------------------!
[2801]228! Description:
229! ------------
[4649]230!> General routine:
231!> Initializing actions of the particle interface check particle boundary conditions for the child
232!> models
233!--------------------------------------------------------------------------------------------------!
[2801]234 SUBROUTINE pmcp_g_init
[4649]235
[2801]236    IMPLICIT NONE
237
[4649]238    INTEGER(iwp) ::  nr_childs  !< Number of child models of the current model
239
[2967]240#if defined( __parallel )
241
[3883]242    nr_childs = get_number_of_children()
[2801]243!
244!-- Check if the current model has child models
[4649]245    IF ( nr_childs > 0 )  THEN
[2801]246       ALLOCATE( nr_part(nysg:nyng, nxlg:nxrg) )
247       ALLOCATE( part_adr(nysg:nyng, nxlg:nxrg) )
[4649]248       nr_part  = 0
249       part_adr = 0
[2801]250    ENDIF
251
252!
253!-- Set the boundary conditions to nested for all non root (i.e child) models
[4649]254    IF ( cpl_id > 1 )  THEN
255
256       IF ( ibc_par_t /= 3 )  THEN
[2801]257          ibc_par_t  = 3
[3046]258          message_string = 'In Child model:  ibc_par_t is automatically set to nested '
259          CALL message( 'pmcp_g_init ', 'PA0477', 0, 1, 0, 6, 0 )
[2801]260       ENDIF
[4649]261
262       IF ( ibc_par_lr /= 3 )  THEN
[2801]263          ibc_par_lr = 3
[3046]264          message_string = 'In Child model:  ibc_par_lr is automatically set to nested '
265          CALL message( 'pmcp_g_init ', 'PA0478', 0, 1, 0, 6, 0 )
[2801]266       ENDIF
[4649]267
268       IF ( ibc_par_ns /= 3 )  THEN
[2801]269          ibc_par_ns = 3
[3046]270          message_string = 'In Child model:  ibc_par_ns is automatically set to nested '
271          CALL message( 'pmcp_g_init ', 'PA0479', 0, 1, 0, 6, 0 )
[2801]272       ENDIF
[4649]273
[2801]274    ENDIF
275
[2967]276#endif
[2801]277 END SUBROUTINE pmcp_g_init
[4649]278!--------------------------------------------------------------------------------------------------!
[2801]279! Description:
280! ------------
[4649]281!> General routine:
282!> Allocate the MPI windows
283!--------------------------------------------------------------------------------------------------!
[2801]284 SUBROUTINE pmcp_g_alloc_win
[4649]285
[2801]286    IMPLICIT NONE
287
[4649]288    INTEGER(iwp) ::  child_id   !< Id of a child model
289    INTEGER(iwp) ::  ierr       !< error code
290    INTEGER(iwp) ::  ipl        !< left boundary in coarse(parent) index space
291    INTEGER(iwp) ::  ipr        !< right boundary in coarse(parent) index space
292    INTEGER(iwp) ::  jps        !< south boundary in coarse(parent) index space
293    INTEGER(iwp) ::  jpn        !< north boundary in coarse(parent) index space
294    INTEGER(iwp) ::  m          !< loop index
295    INTEGER(iwp) ::  nr_childs  !< Number of child models of the current model
296
297    INTEGER ::  parsize  !<
298
299    INTEGER(iwp),DIMENSION(1) ::  buf_shape  !<
300
301    TYPE(C_PTR), SAVE ::  ptr  !<
302
303    TYPE(particle_type),DIMENSION(:),POINTER ::  win_buffer  !<
304
305
[2967]306#if defined( __parallel )
[4649]307    INTEGER(KIND=MPI_ADDRESS_KIND) ::  parsize_mpi_address_kind  !<
308    INTEGER(KIND=MPI_ADDRESS_KIND) ::  winsize                   !<
[2967]309
[2801]310!
311!-- If the model has a parent model prepare the structures for transfer
[4649]312    IF ( cpl_id > 1 )  THEN
[2809]313
[4649]314       parsize_mpi_address_kind = STORAGE_SIZE(zero_particle) / 8
[2809]315
[2884]316       CALL MPI_ALLOC_MEM( parsize_mpi_address_kind , MPI_INFO_NULL, ptr, ierr )
317       parsize = parsize_mpi_address_kind
[2801]318       buf_shape(1) = 1
319       CALL C_F_POINTER( ptr, win_buffer, buf_shape )
[4649]320       CALL MPI_WIN_CREATE( win_buffer, parsize_mpi_address_kind, parsize, MPI_INFO_NULL,          &
321                            me%intra_comm, particle_win_child, ierr )
[2801]322
323!
324!--    Child domain boundaries in the parent index space
[3948]325       ipl = parent_bound(1)
326       ipr = parent_bound(2)
327       jps = parent_bound(3)
328       jpn = parent_bound(4)
[2801]329
[3948]330       ALLOCATE( coarse_particles(jps:jpn,ipl:ipr) )
[2801]331
332       coarse_particles(:,:)%nr_particle = 0
333    ENDIF
334
335!
336!-- If the model has child models prepare the structures for transfer
[3883]337    nr_childs = get_number_of_children()
[2801]338    IF ( nr_childs > 0 )   THEN
339       ALLOCATE( particle_win_parent(nr_childs) )
340       ALLOCATE( buf_ptr(nr_childs) )
341       DO  m = 1, nr_childs
342          child_id = get_childid(m)
[4649]343          parsize_mpi_address_kind = STORAGE_SIZE(zero_particle) / 8
[2884]344          parsize = parsize_mpi_address_kind
[2809]345
[2884]346          winsize = max_nr_particle_in_rma_win * parsize_mpi_address_kind
[2801]347          CALL MPI_ALLOC_MEM( winsize , MPI_INFO_NULL, buf_ptr(m), ierr )
348          buf_shape(1) = max_nr_particle_in_rma_win
349          CALL C_F_POINTER( buf_ptr(m), win_buffer, buf_shape )
[4649]350          CALL MPI_WIN_CREATE( win_buffer, winsize, parsize, MPI_INFO_NULL,                        &
351                               children(child_id)%intra_comm, particle_win_parent(m), ierr )
[2801]352          ENDDO
353    ENDIF
354
[2967]355#endif
[2801]356 END SUBROUTINE pmcp_g_alloc_win
[2967]357
358
[4649]359!--------------------------------------------------------------------------------------------------!
[2801]360! Description:
361! ------------
[4649]362!> Child routine:
363!> Read/get particles out of the parent MPI window
364!--------------------------------------------------------------------------------------------------!
[2801]365 SUBROUTINE pmcp_c_get_particle_from_parent
[4649]366
[2801]367    IMPLICIT NONE
368
[4649]369    INTEGER(iwp) ::  i     !< x grid index
370    INTEGER(iwp) ::  ierr  !< error code
371    INTEGER(iwp) ::  ij    !< combined xy index for the buffer array
372    INTEGER(iwp) ::  ip    !< loop index (child PEs)
373    INTEGER(iwp) ::  ipl   !< left boundary in coarse(parent) index space
374    INTEGER(iwp) ::  j     !< y grid index
375    INTEGER(iwp) ::  jps   !< south boundary in coarse(parent) index space
376    INTEGER(iwp) ::  nr    !< number of particles to receive from a parent box
377
[2884]378    INTEGER ::  parsize !<
[2967]379
380#if defined( __parallel )
[4649]381    TYPE(pedef), POINTER ::  ape  !< TO_DO Klaus: give a description and better name of the variable
[2967]382
[4649]383    INTEGER(KIND=MPI_ADDRESS_KIND) ::  parsize_mpi_address_kind  !<
384    INTEGER(KIND=MPI_ADDRESS_KIND) ::  target_disp               !<
[2801]385
386    IF ( cpl_id > 1 )  THEN
[4649]387
[2801]388       CALL pmc_c_getbuffer( particle_transfer = .TRUE. ) !Get number of particle/column and offset in RMA window xx
389
390!
391!--    Wait for buffer to fill.
392!
[4649]393!--    The parent side (in pmc_s_fillbuffer) is filling the buffer in the MPI RMA window. When the
394!--    filling is complete, a MPI_BARRIER is called. The child is not allowd to access the
395!--    parent-buffer before it is completely filled. Synchronization is done implicitely in
396!--    pmc_c_getbuffer and pmc_s_fillbuffer on the parent side.
[2801]397
[3948]398       ipl = parent_bound(1)
399       jps = parent_bound(3)
[2801]400
401       DO  ip = 1, me%inter_npes
402
403          ape => me%pes(ip)
404
405          DO  ij = 1, ape%nrele
[4649]406              j  = ape%locind(ij)%j + jps - 1
407              i  = ape%locind(ij)%i + ipl - 1
[2801]408              nr = nr_partc(j,i)
409              IF ( nr > 0 )  THEN
410
411                 CALL check_and_alloc_coarse_particle (i, j, nr)
[4649]412                 parsize_mpi_address_kind = STORAGE_SIZE(zero_particle) / 8
[2884]413                 parsize = parsize_mpi_address_kind
[2801]414                 target_disp = part_adrc(j,i) - 1
[4649]415                 CALL MPI_WIN_LOCK( MPI_LOCK_SHARED , ip - 1, 0, particle_win_child, ierr )
416                 CALL MPI_GET( coarse_particles(j,i)%parent_particles, nr * parsize, MPI_BYTE,     &
417                               ip - 1, target_disp, nr * parsize, MPI_BYTE, particle_win_child,    &
418                               ierr )
[2801]419                 CALL MPI_WIN_UNLOCK( ip - 1, particle_win_child, ierr )
420             ENDIF
421             coarse_particles(j,i)%nr_particle = nr
422          ENDDO
423       ENDDO
424
425       CALL c_copy_particle_to_child_grid
426    ENDIF
427
[2967]428#endif
[2801]429 END SUBROUTINE pmcp_c_get_particle_from_parent
[2967]430
431
[4649]432!--------------------------------------------------------------------------------------------------!
[2801]433! Description:
434! ------------
[4649]435!> Child routine:
[2801]436!> Write/put particles into the parent MPI window
[4649]437!--------------------------------------------------------------------------------------------------!
438 SUBROUTINE pmcp_c_send_particle_to_parent
439
[2801]440    IMPLICIT NONE
441
[4649]442    INTEGER(iwp) ::  disp_offset             !<
443    INTEGER(iwp) ::  i                       !< x loop index
444    INTEGER(iwp) ::  ierr                    !< error code
445    INTEGER(iwp) ::  ij                      !< combined xy index for the buffer array
446    INTEGER(iwp) ::  ip                      !< loop index (child PEs)
447    INTEGER(iwp) ::  ipl                     !< left boundary in coarse(parent) index space
448    INTEGER(iwp) ::  ipr                     !< right boundary in coarse(parent) index space
449    INTEGER(iwp) ::  j                       !< y loop index
450    INTEGER(iwp) ::  jpn                     !< north boundary in coarse(parent) index space
451    INTEGER(iwp) ::  jps                     !< south boundary in coarse(parent) index space
452    INTEGER(iwp) ::  max_nr_particle_per_pe  !< maximum number of particles per PE (depending on grid apect ratio)
453    INTEGER(iwp) ::  n                       !< shorter variable name for nr_fine_in_coarse
454    INTEGER(iwp) ::  nr                      !< shorter variabel name for nr_partc
455    INTEGER(iwp) ::  pe_offset               !< offset index of the current PE
456
457    INTEGER ::  parsize  !<
458
459    REAL(wp) ::  eps=0.00001  !< used in calculations to avoid rounding errors
460    REAL(wp) ::  xx           !< number of fine grid cells inside a coarse grid cell in x-direction
461    REAL(wp) ::  yy           !< number of fine grid cells inside a coarse grid cell in y-direction
462
[2967]463 !   TYPE(particle_type) ::  dummy_part !< dummy particle (needed for size calculations)
464
465#if defined( __parallel )
[4649]466    TYPE(pedef), POINTER ::  ape  !< TO_DO Klaus: give a description and better name of the variable
[2967]467
[4649]468    INTEGER(KIND=MPI_ADDRESS_KIND) ::  parsize_mpi_address_kind  !<
469    INTEGER(KIND=MPI_ADDRESS_KIND) ::  target_disp               !<
470
471
[2801]472    IF ( cpl_id > 1 )  THEN
473       CALL c_copy_particle_to_coarse_grid
474
475!
476!--    Child domain boundaries in the parent index space
477
[3948]478       ipl = parent_bound(1)
479       ipr = parent_bound(2)
480       jps = parent_bound(3)
481       jpn = parent_bound(4)
[2801]482
483       nr_partc = 0
[4649]484
[3948]485       DO i = ipl, ipr
486          DO j = jps, jpn
[2801]487             nr_partc(j,i) = coarse_particles(j,i)%nr_particle
488          ENDDO
489       ENDDO
490       part_adrc = 0
491
492!
[4649]493!--    Compute number of fine grid cells in coarse grid (one direction)
[3883]494       xx = ( pg%dx + eps ) / dx ! +eps to avoid rounding error
495       yy = ( pg%dy + eps ) / dy
[4649]496       nr_fine_in_coarse = MAX( INT( xx ), INT( yy ) )
[2801]497
[4649]498       IF ( MOD( coord_x(0), pg%dx ) /= 0.0 .OR. MOD( coord_y(0), pg%dy ) /= 0.0 )  THEN
[2801]499          nr_fine_in_coarse = nr_fine_in_coarse + 1
500       ENDIF
501
502!
503!--    Assign a number to my child PE to select different areas in the RMA window on server side
504!--    With this number a square of child PEs is defined which share the same coarse grid cells
505
[4649]506       n           = nr_fine_in_coarse ! Local variable n to make folloing statements shorter
[2801]507       pe_offset   = MOD( myidx, n ) * n + MOD( myidy, n)
508       max_nr_particle_per_pe = max_nr_particle_in_rma_win / ( n * n )
509       disp_offset            = pe_offset * max_nr_particle_per_pe
[4649]510       parsize_mpi_address_kind = STORAGE_SIZE(zero_particle) /8
[2884]511       parsize = parsize_mpi_address_kind
[2801]512       DO  ip = 1, me%inter_npes
513
514          ape => me%pes(ip)
515
516          target_disp = disp_offset
[4649]517          DO  ij = 1, ape%nrele
[3948]518             j  = ape%locind(ij)%j + jps - 1
519             i  = ape%locind(ij)%i + ipl - 1
[2801]520             nr = nr_partc(j,i)
[4649]521             IF( nr > 0 )  THEN
522                IF ( target_disp + nr - disp_offset >= max_nr_particle_per_pe )  THEN
523                   WRITE( 9, * ) 'RMA window too small on child ',                                 &
524                                 target_disp + nr - disp_offset, max_nr_particle_per_pe,           &
525                                 max_nr_particle_in_rma_win
[2801]526                   message_string = 'RMA window too small on child'
[3046]527                   CALL message( 'pmci_create_child_arrays', 'PA0480', 3, 2, 0, 6, 0 )
[2801]528                ENDIF
[4649]529                CALL MPI_WIN_LOCK( MPI_LOCK_SHARED , ip - 1, 0, particle_win_child, ierr )
530                CALL MPI_PUT( coarse_particles(j,i)%parent_particles, nr * parsize, MPI_BYTE,      &
531                              ip - 1, target_disp, nr * parsize, MPI_BYTE,  particle_win_child,    &
532                              ierr )
[2801]533                CALL MPI_WIN_UNLOCK( ip - 1, particle_win_child, ierr )
534                part_adrc(j,i) = target_disp + 1
535                target_disp    = target_disp + nr
536             ENDIF
537          ENDDO
538       ENDDO
539
540       CALL pmc_c_putbuffer ( particle_transfer = .TRUE. )   !Send new number of particle/column and offset to parent
541
542    ENDIF
543
[2967]544#endif
[2801]545 END SUBROUTINE pmcp_c_send_particle_to_parent
[2967]546
547
[4649]548!--------------------------------------------------------------------------------------------------!
[2801]549! Description:
550! ------------
[4649]551!> Parent routine:
552!> Write particles into the MPI window
553!--------------------------------------------------------------------------------------------------!
[2801]554 SUBROUTINE pmcp_p_fill_particle_win
[4649]555
[2801]556    IMPLICIT NONE
557
[4649]558    INTEGER(iwp) ::  child_id             !< id of the child model
559    INTEGER(iwp) ::  i                    !< x grid box index
560    INTEGER(iwp) ::  ij                   !< combined xy index for the buffer array
561    INTEGER(iwp) ::  ip                   !< loop index (child PEs)
562    INTEGER(iwp) ::  j                    !< y grid box index
563    INTEGER(iwp) ::  k                    !< z grid box index
564    INTEGER(iwp) ::  m                    !< loop index (number of childs)
565    INTEGER(iwp) ::  n                    !< loop index (number of particles)
566    INTEGER(iwp) ::  nr_part_col          !< Number of particles to transfer per column
567    INTEGER(iwp) ::  number_of_particles  !<
568    INTEGER(iwp) ::  pindex               !<
569    INTEGER(iwp) ::  tot_particle_count   !< Total number of particles per child
570
571    INTEGER(iwp),DIMENSION(1) ::  buf_shape  !<
572
573    LOGICAL      ::  active_particle  !< Particles located in the fine/child grid area are marked as active (to be transferred)
574    LOGICAL,SAVE ::  lfirst = .TRUE.  !<
575
576    REAL(wp) ::  dx_child    !< child grid spacing
577    REAL(wp) ::  dy_child    !< child grid spacing
578    REAL(wp) ::  dz_child    !< child grid spacing
579    REAL(wp) ::  ny_coord    !< north coordinate of child grid
580    REAL(wp) ::  ny_coord_b  !< north coordinate of child grid boundary
581    REAL(wp) ::  lx_coord    !< left coordinate of child grid
582    REAL(wp) ::  lx_coord_b  !< left coordinate of child grid boundary
583    REAL(wp) ::  rx_coord    !< right coordinate of child grid
584    REAL(wp) ::  rx_coord_b  !< right coordinate of child grid boundary
585    REAL(wp) ::  sy_coord    !< south coordinate of child grid
586    REAL(wp) ::  sy_coord_b  !< south coordinate of child grid boundary
587    REAL(wp) ::  uz_coord    !< upper coordinate of child grid
588    REAL(wp) ::  uz_coord_b  !< upper coordinate of child grid boundary
589    REAL(wp) ::  x           !< particle position
590    REAL(wp) ::  y           !< particle position
591    REAL(wp) ::  z           !< particle position
592
593
[2967]594#if defined( __parallel )
[4649]595    TYPE(pedef), POINTER ::  ape  !< TO_DO Klaus: give a description and better name of the variable
[2801]596
[3883]597    DO  m = 1, get_number_of_children()
[2801]598
599       child_id = get_childid(m)
600
[4649]601       CALL get_child_edges( m, lx_coord, lx_coord_b, rx_coord, rx_coord_b, sy_coord, sy_coord_b,  &
602                             ny_coord, ny_coord_b, uz_coord, uz_coord_b)
[2801]603
604       CALL get_child_gridspacing( m, dx_child, dy_child, dz_child )
605
606       IF ( lfirst )   THEN
[4649]607          WRITE( 9, '(a,5f10.2)' ) 'edges          ', lx_coord,rx_coord, sy_coord, ny_coord,       &
608                                                      uz_coord
609          WRITE( 9, '(a,5f10.2)' ) 'edges boundary ', lx_coord_b, rx_coord_b, sy_coord_b,          &
610                                                      ny_coord_b, uz_coord_b
611          WRITE( 9, '(a,5f10.2)' ) 'child spacing  ', dx_child, dy_child, dz_child,                &
612                                                      lower_left_coord_x,lower_left_coord_y
613     ENDIF
[2801]614!
[4649]615!--    Reset values for every child
616       tot_particle_count = 0
617       nr_part  = 0
[2801]618       part_adr = 0
619       pindex   = 1
620
621       buf_shape(1) = max_nr_particle_in_rma_win
622       CALL C_F_POINTER( buf_ptr(m), particle_in_win , buf_shape )
623
624       DO  ip = 1, children(child_id)%inter_npes
625
626          ape => children(child_id)%pes(ip)
627
[4649]628          nr_part_col   = 0
629
[2801]630          DO  ij = 1, ape%nrele
[4649]631
[2801]632!
633!--          Inside the PMC adressing of 3d arrays starts with 1
634             i = ape%locind(ij)%i + nxl - nbgp - 1
635             j = ape%locind(ij)%j + nys - nbgp - 1
[3241]636             nr_part_col = 0   ! Number of particles to transfer per column
[2801]637             part_adr(j,i) = pindex
[4649]638
[2801]639             DO  k = nzb + 1, nzt
640                number_of_particles = prt_count(k,j,i)
[4649]641
[2801]642                IF ( number_of_particles <= 0 )  CYCLE
[4649]643
[2801]644                particles => grid_particles(k,j,i)%particles(1:number_of_particles)
645
646                ! Select particles within boundary area
647
[4649]648                DO  n = 1, number_of_particles
[2801]649                   x = particles(n)%x
650                   y = particles(n)%y
651                   z = particles(n)%z
652!
[4649]653!--                Check if the particle is located in the fine grid area
654                   active_particle = ( (x > lx_coord .AND. x < rx_coord) .AND.                     &
655                                       (y > sy_coord .AND. y < ny_coord) .AND.                     &
656                                       (z > 0.000000001 .AND. z < uz_coord) )
[2801]657                   IF ( active_particle .AND. particles(n)%particle_mask )  THEN
[4649]658
[2801]659                      particle_in_win(pindex) = particles(n)
660!
661!--                   Change particle positions and origin relative to global origin
662                      particle_in_win(pindex)%x = particle_in_win(pindex)%x + lower_left_coord_x
663                      particle_in_win(pindex)%y = particle_in_win(pindex)%y + lower_left_coord_y
[4649]664                      particle_in_win(pindex)%origin_x = particle_in_win(pindex)%origin_x          &
665                                                         + lower_left_coord_x
666                      particle_in_win(pindex)%origin_y = particle_in_win(pindex)%origin_y          &
667                                                         + lower_left_coord_y
[2801]668
669                      tot_particle_count = tot_particle_count + 1
670                      nr_part_col        = nr_part_col + 1
671                      pindex             = pindex + 1
[4649]672                      IF ( pindex > max_nr_particle_in_rma_win )  THEN
673                         WRITE( 9, * ) 'RMA window too small on parent ', pindex,                  &
674                                                                          max_nr_particle_in_rma_win
[2801]675                         message_string = 'RMA window too small on parent'
676                         CALL message( 'pmci_create_child_arrays', 'PA0481', 3, 2, 0, 6, 0 )   ! PA number has to be adjusted
677                     ENDIF
678                   END IF
679                ENDDO
680             ENDDO
681             nr_part(j,i) = nr_part_col
682          ENDDO
683       ENDDO
684
685       CALL pmc_s_fillbuffer( child_id, particle_transfer = .TRUE. )
686    ENDDO
687
688    lfirst = .FALSE.
689
[2967]690#endif
[2801]691 END SUBROUTINE pmcp_p_fill_particle_win
[2967]692
[4649]693
694!--------------------------------------------------------------------------------------------------!
[2801]695! Description:
696! ------------
[4649]697!> Parent routine:
698!> Delete particles from the MPI window
699!--------------------------------------------------------------------------------------------------!
[2801]700 SUBROUTINE pmcp_p_empty_particle_win
[2967]701
[2801]702    IMPLICIT NONE
703
[4649]704    INTEGER(iwp) ::  child_id  !< model id of the child
705    INTEGER(iwp) ::  ip        !< loop index (child PEs)
706    INTEGER(iwp) ::  m         !< loop index (number of childs)
[2801]707
[4649]708    INTEGER(iwp),DIMENSION(1) ::  buf_shape  !<
[2801]709
[2967]710#if defined( __parallel )
[3883]711    DO  m = 1, get_number_of_children()
[2801]712
713       child_id = get_childid(m)
714
715       buf_shape(1) = max_nr_particle_in_rma_win
716       CALL C_F_POINTER( buf_ptr(m), particle_in_win , buf_shape )
717
718!
[4649]719!--    In some cells of the coarse grid, there are contributions from more than one child process.
720!--    Therefore p_copy_particle_to_org_grid is done for one child process per call
721       DO  ip = 1, pmc_s_get_child_npes( child_id )
[2801]722
723          nr_part  = 0
724          part_adr = 0
725
[4649]726          CALL pmc_s_getdata_from_buffer( child_id, particle_transfer = .TRUE.,                    &
[3241]727                                          child_process_nr = ip )
728          CALL p_copy_particle_to_org_grid( m )
[2801]729       ENDDO
730
731    ENDDO
732
[2967]733#endif
[2801]734 END SUBROUTINE pmcp_p_empty_particle_win
[2967]735
736
[4649]737!--------------------------------------------------------------------------------------------------!
[2801]738! Description:
739! ------------
[4649]740!> Parent routine:
741!> After the transfer mark all parent particles that are still inside on of the child areas for
742!> deletion.
743!--------------------------------------------------------------------------------------------------!
[2801]744 SUBROUTINE pmcp_p_delete_particles_in_fine_grid_area
745
746    IMPLICIT NONE
747
[4649]748    INTEGER(iwp) ::  i  !< loop index (x grid)
749    INTEGER(iwp) ::  j  !< loop index (y grid)
750    INTEGER(iwp) ::  k  !< loop index (z grid)
751    INTEGER(iwp) ::  m  !< loop index (number of particles)
752    INTEGER(iwp) ::  n  !< loop index (number of childs)
753
754    LOGICAL ::  to_delete  !< particles outside of model domain are marked as to_delete
755
756    REAL(wp) ::  dx_child    !< child grid spacing
757    REAL(wp) ::  dy_child    !< child grid spacing
758    REAL(wp) ::  dz_child    !< child grid spacing
759    REAL(wp) ::  ny_coord    !< north coordinate of child grid
760    REAL(wp) ::  ny_coord_b  !< north coordinate of child grid boundary
761    REAL(wp) ::  lx_coord    !< left coordinate of child grid
762    REAL(wp) ::  lx_coord_b  !< left coordinate of child grid boundary
763    REAL(wp) ::  rx_coord    !< right coordinate of child grid
764    REAL(wp) ::  rx_coord_b  !< right coordinate of child grid boundary
765    REAL(wp) ::  sy_coord    !< south coordinate of child grid
766    REAL(wp) ::  sy_coord_b  !< south coordinate of child grid boundary
767    REAL(wp) ::  uz_coord    !< upper coordinate of child grid
768    REAL(wp) ::  uz_coord_b  !< upper coordinate of child grid boundary
769    REAL(wp) ::  x           !< particle position
770    REAL(wp) ::  y           !< particle position
771    REAL(wp) ::  z           !< particle position
772
[2967]773#if defined( __parallel )
[3883]774    DO  m = 1, get_number_of_children()
[4649]775       CALL get_child_edges( m, lx_coord, lx_coord_b, rx_coord, rx_coord_b, sy_coord, sy_coord_b,  &
776                             ny_coord, ny_coord_b, uz_coord, uz_coord_b )
[2801]777
778
[4649]779       CALL get_child_gridspacing( m, dx_child, dy_child, dz_child )
[2801]780
[4649]781       DO  i = nxl, nxr
782          DO  j = nys, nyn
783             DO  k = nzb, nzt
[2801]784                number_of_particles = prt_count(k,j,i)
785
[4649]786                IF ( number_of_particles == 0 )  CYCLE
[2801]787
788                particles => grid_particles(k,j,i)%particles(1:number_of_particles)
789
[4649]790                DO  n = 1, number_of_particles
[2801]791                   x = particles(n)%x
792                   y = particles(n)%y
793                   z = particles(n)%z
794
[4649]795                   to_delete = ( (x > lx_coord .AND. x < rx_coord) .AND.                           &
796                                 (y > sy_coord .AND. y < ny_coord) .AND.                           &
797                                 (z > 0.000000001 .AND. z < uz_coord) )
[2801]798
[4649]799                   IF ( to_delete )  THEN
[2801]800                      particles(n)%particle_mask = .FALSE.
801                   ENDIF
802                ENDDO
803             ENDDO
804          ENDDO
805       ENDDO
806
807    ENDDO
808
[2967]809#endif
[2801]810 END SUBROUTINE pmcp_p_delete_particles_in_fine_grid_area
[2967]811
812
[4649]813!--------------------------------------------------------------------------------------------------!
[2801]814! Description:
815! ------------
[4649]816!> General routine:
817!> Print the total number of of the current model and its child models into a file
818!--------------------------------------------------------------------------------------------------!
819 SUBROUTINE pmcp_g_print_number_of_particles( my_time, local_nr_particles )
820
821    USE pegrid,                                                                                    &
[2801]822        ONLY: myid
[4649]823
[2801]824    IMPLICIT NONE
825
[4649]826    CHARACTER(LEN=16) ::  fname  !< filename
827
828    INTEGER(iwp) ::  child_id            !<
829    INTEGER(iwp) ::  child_nr_particles  !< total number of particles in all child models
830    INTEGER(iwp) ::  ierr                !< error code
831    INTEGER(iwp) ::  m                   !< loop index (number of childs
832
833    INTEGER(iwp),INTENT(IN) ::  local_nr_particles  !<
834
835    INTEGER(iwp),DIMENSION(2) ::  ivalr  !< integer value to be received
836    INTEGER(iwp),DIMENSION(2) ::  ivals  !< integer value to be send
837
838    LOGICAL, SAVE ::  is_file_open=.FALSE.  !<
839
840    REAL(wp),INTENT(IN) ::  my_time  !<
841
842
[2967]843#if defined( __parallel )
[2801]844    child_nr_particles = 0
845    IF ( myid == 0 )  THEN
[4649]846       IF ( cpl_id > 1 )  THEN
[2801]847          ivals(1) = local_nr_particles
848          CALL pmc_send_to_parent( ivals, 1, 0, 400, ierr )
849       ENDIF
850       DO  m = 1, SIZE( pmc_parent_for_child ) - 1
851
852          child_id = pmc_parent_for_child(m)
853          CALL pmc_recv_from_child( child_id, ivalr, 1, 0, 400, ierr )
854          child_nr_particles = child_nr_particles + ivalr(1)
855       ENDDO
856
[4649]857       IF ( SIZE( pmc_parent_for_child ) > 1 )  THEN
858          IF ( .NOT. is_file_open )  THEN  !kk muss noch auf file_open umgestellt werden
859             WRITE( fname, '(a,i2.2)' ) 'nr_particles_', cpl_id
860             OPEN ( 333, FILE = fname )
[2801]861             is_file_open = .true.
862          ENDIF
[4649]863          WRITE( 333, '(f12.2,3i10)' ) my_time, local_nr_particles + child_nr_particles,           &
864                                       local_nr_particles, child_nr_particles
[2801]865       ENDIF
866    ENDIF
867
[2967]868#endif
[2801]869 END SUBROUTINE pmcp_g_print_number_of_particles
870
[2967]871
[4649]872!--------------------------------------------------------------------------------------------------!
873!--------------------------------------------------------------------------------------------------!
[2801]874! Private subroutines
[4649]875!--------------------------------------------------------------------------------------------------!
876!--------------------------------------------------------------------------------------------------!
[2801]877
[4649]878!--------------------------------------------------------------------------------------------------!
[2801]879! Description:
880! ------------
[4649]881!> Child routine
882!> Update the size of the local buffer (coarse_particles)
883!--------------------------------------------------------------------------------------------------!
884 SUBROUTINE check_and_alloc_coarse_particle( ic, jc, nr, with_copy )
885
[2801]886    IMPLICIT NONE
887
[4649]888    INTEGER(iwp),INTENT(IN) ::  ic  !< coarse x grid index
889    INTEGER(iwp),INTENT(IN) ::  jc  !< coarse y grid index
890    INTEGER(iwp),INTENT(IN) ::  nr  !< number of particles to be transferred/stored in local buffer
[2801]891
[4649]892    INTEGER(iwp) ::  new_size  !< new size of the local buffer
893    INTEGER(iwp) ::  old_size  !< old size of the local buffer
894
895    LOGICAL :: with_copy_lo  !< local variable of with copy
896
897    LOGICAL,INTENT(IN),OPTIONAL ::  with_copy  !< copy particles in buffer? or reallocate empty buffer
898
899    TYPE(particle_type), DIMENSION(:), ALLOCATABLE ::  tmp_particles_d  !<
900
[2967]901#if defined( __parallel )
[2801]902    with_copy_lo = .FALSE.
[4649]903    IF ( PRESENT( with_copy ) )  with_copy_lo = with_copy
[2801]904
[4649]905    IF ( .NOT. ALLOCATED( coarse_particles(jc,ic)%parent_particles ) )  THEN
906       new_size = MAX( nr, min_particles_per_column )
[2801]907       ALLOCATE( coarse_particles(jc,ic)%parent_particles(new_size) )
[4649]908    ELSEIF ( nr > SIZE( coarse_particles(jc,ic)%parent_particles ) )  THEN
[2801]909
910       old_size = SIZE( coarse_particles(jc,ic)%parent_particles )
911       new_size = old_size * ( 1.0_wp + alloc_factor / 100.0_wp )
[4649]912       new_size = MAX( nr, new_size, old_size + min_particles_per_column )
[2801]913
914!
915!--    Copy existing particles to new particle buffer
[4649]916       IF ( with_copy_lo )  THEN
[2801]917          ALLOCATE( tmp_particles_d(old_size) )
918          tmp_particles_d(1:old_size) = coarse_particles(jc,ic)%parent_particles
919
920          DEALLOCATE( coarse_particles(jc,ic)%parent_particles )
921          ALLOCATE( coarse_particles(jc,ic)%parent_particles(new_size) )
922
923          coarse_particles(jc,ic)%parent_particles(1:old_size)          = tmp_particles_d(1:old_size)
924          coarse_particles(jc,ic)%parent_particles(old_size+1:new_size) = zero_particle
925
926          DEALLOCATE( tmp_particles_d )
927!
928!--    allocate or reallocate an empty buffer
[4649]929       ELSE
[2801]930          DEALLOCATE( coarse_particles(jc,ic)%parent_particles )
931          ALLOCATE( coarse_particles(jc,ic)%parent_particles(new_size) )
932       ENDIF
933    ENDIF
934
[2967]935#endif
[2801]936 END SUBROUTINE check_and_alloc_coarse_particle
937
[2967]938
[4649]939!--------------------------------------------------------------------------------------------------!
[2801]940! Description:
941! ------------
[4649]942!> Child routine:
943!> Copy/sort particles out of the local buffer into the respective grid boxes
944!--------------------------------------------------------------------------------------------------!
[2801]945 SUBROUTINE c_copy_particle_to_child_grid
[4649]946
[2801]947    IMPLICIT NONE
948
[4649]949    INTEGER(iwp) ::  ic   !< coarse x grid index
950    INTEGER(iwp) ::  ipl  !< left boundary in coarse(parent) index space
951    INTEGER(iwp) ::  ipr  !< right boundary in coarse(parent) index space
952    INTEGER(iwp) ::  ip   !< x grid index
953    INTEGER(iwp) ::  jc   !< coarse y grid index
954    INTEGER(iwp) ::  jpn  !< north boundary in coarse(parent) index space
955    INTEGER(iwp) ::  jps  !< south boundary in coarse(parent) index space
956    INTEGER(iwp) ::  jp   !< y grid index
957    INTEGER(iwp) ::  kp   !< z grid index
958    INTEGER(iwp) ::  n    !< loop index (number of particles)
959    INTEGER(iwp) ::  nr   !< short variable for name number or particles
960
961    REAL(wp) ::  xc   !< child x coordinate
962    REAL(wp) ::  xoc  !< child x origin
963    REAL(wp) ::  yc   !< child y coordinate
964    REAL(wp) ::  yoc  !< child y origin
965    REAL(wp) ::  zc   !< child z coordinate
966
[2967]967#if defined( __parallel )
[2801]968!
969!-- Child domain boundaries in the parent index space
[3948]970    ipl = parent_bound(1)
971    ipr = parent_bound(2)
972    jps = parent_bound(3)
973    jpn = parent_bound(4)
[2801]974
[4649]975    DO  ic = ipl, ipr
976       DO  jc = jps, jpn
[2801]977          nr = coarse_particles(jc,ic)%nr_particle
978
[4649]979          IF ( nr > 0 )  THEN
980             DO  n = 1, nr
[2801]981                xc =  coarse_particles(jc,ic)%parent_particles(n)%x-lower_left_coord_x
982                yc =  coarse_particles(jc,ic)%parent_particles(n)%y-lower_left_coord_y
983                zc =  coarse_particles(jc,ic)%parent_particles(n)%z
984                xoc = coarse_particles(jc,ic)%parent_particles(n)%origin_x-lower_left_coord_x
985                yoc = coarse_particles(jc,ic)%parent_particles(n)%origin_y-lower_left_coord_y
986                ip = xc / dx
987                jp = yc / dy
988                kp = nzt
[4649]989                DO WHILE ( zw(kp-1) > zc .AND. kp > nzb + 1 )  ! kk search loop has to be optimzed !!!
[2801]990                   kp = kp - 1
[4649]991                ENDDO
[2801]992
[4649]993                IF ( ip >= nxl  .AND.  ip <= nxr  .AND.  jp >= nys  .AND.  jp <= nyn )  THEN
[2801]994                   prt_count(kp,jp,ip) = prt_count(kp,jp,ip) + 1
[4649]995                   IF ( prt_count(kp,jp,ip) > SIZE( grid_particles(kp,jp,ip)%particles ) )  THEN
[4017]996                      CALL pmc_realloc_particles_array( ip, jp, kp )
[2801]997                   ENDIF
[4649]998                   coarse_particles(jc,ic)%parent_particles(n)%x = xc         ! Adjust coordinates to child grid
[3241]999                   coarse_particles(jc,ic)%parent_particles(n)%y = yc
[4649]1000                   coarse_particles(jc,ic)%parent_particles(n)%origin_x = xoc ! Adjust origins to child grid
[2801]1001                   coarse_particles(jc,ic)%parent_particles(n)%origin_y = yoc
[4649]1002                   grid_particles(kp,jp,ip)%particles(prt_count(kp,jp,ip))                         &
1003                                                       = coarse_particles(jc,ic)%parent_particles(n)
[2801]1004                ENDIF
1005             ENDDO
1006          ENDIF
1007       ENDDO
1008    ENDDO
1009
[2967]1010#endif
[2801]1011 END SUBROUTINE c_copy_particle_to_child_grid
[2967]1012
1013
[4649]1014!--------------------------------------------------------------------------------------------------!
[2801]1015! Description:
1016! ------------
[4649]1017!> Child routine:
1018!> Copy particles which left the model area into the local buffer
1019!--------------------------------------------------------------------------------------------------!
[2801]1020 SUBROUTINE c_copy_particle_to_coarse_grid
[4649]1021
[2801]1022    IMPLICIT NONE
1023
[4649]1024    INTEGER(iwp) ::  i     !< loop index (x grid)
1025    INTEGER(iwp) ::  ic    !< loop index (coarse x grid)
1026    INTEGER(iwp) ::  ipl   !< left boundary in coarse(parent) index space
1027    INTEGER(iwp) ::  ipr   !< left boundary in coarse(parent) index space
1028    INTEGER(iwp) ::  ierr  !< error code
1029    INTEGER(iwp) ::  j     !< loop index (y grid)
1030    INTEGER(iwp) ::  jc    !< loop index (coarse y grid)
1031    INTEGER(iwp) ::  jps   !< south boundary in coarse(parent) index space
1032    INTEGER(iwp) ::  jpn   !< north boundary in coarse(parent) index space
1033    INTEGER(iwp) ::  k     !< loop index (z grid)
1034    INTEGER(iwp) ::  n     !< loop index (number of particles)
1035
1036    LOGICAL ::  boundary_particle  !<
1037
1038    REAL(wp) ::  x        !< x coordinate
1039    REAL(wp) ::  xo       !< x origin
1040    REAL(wp) ::  x_left   !< absolute left boundary
1041    REAL(wp) ::  x_right  !< absolute right boundary
1042    REAL(wp) ::  y        !< left boundary in coarse(parent) index space
1043    REAL(wp) ::  yo       !< y origin
1044    REAL(wp) ::  y_north  !< absolute north boundary
1045    REAL(wp) ::  y_south  !< absoulte south boundary
1046    REAL(wp) ::  z        !< z coordinate
1047
[2967]1048#if defined( __parallel )
[2801]1049!
1050!-- Child domain boundaries in the parent index space
1051
[3948]1052    ipl = parent_bound(1)
1053    ipr = parent_bound(2)
1054    jps = parent_bound(3)
1055    jpn = parent_bound(4)
[2801]1056
1057!
[4649]1058!-- Absolute coordinates
[2801]1059    x_left  = coord_x(0)
1060    x_right = coord_x(nx) + dx
1061    y_south = coord_y(0)
1062    y_north = coord_y(ny) + dy
[4649]1063!
1064!-- Clear Particle Buffer
1065    DO  ic = ipl, ipr
1066       DO  jc = jps, jpn
[2801]1067          coarse_particles(jc,ic)%nr_particle = 0
1068       ENDDO
1069    ENDDO
1070
1071!
[4649]1072!-- Determine particles which leave the inner area in east or west dirextion.
[2801]1073!-- Compute only first (nxl) and last (nxr) loop iterration.
1074    DO  i = nxl, nxr
1075       DO  j = nys, nyn
1076          DO  k = nzb + 1, nzt
1077             number_of_particles = prt_count(k,j,i)
1078             IF ( number_of_particles <= 0 )  CYCLE
1079             particles => grid_particles(k,j,i)%particles(1:number_of_particles)
1080             DO  n = 1, number_of_particles
1081                IF ( particles(n)%particle_mask )  THEN
[4649]1082                   x  = particles(n)%x+ lower_left_coord_x
1083                   y  = particles(n)%y+ lower_left_coord_y
[2801]1084                   xo = particles(n)%origin_x + lower_left_coord_x
1085                   yo = particles(n)%origin_y + lower_left_coord_y
[4649]1086                   z  = particles(n)%z
1087
[2801]1088                   boundary_particle = .FALSE.
[4649]1089                   boundary_particle = boundary_particle .OR. ( x < x_left  .OR. x > x_right )
1090                   boundary_particle = boundary_particle .OR. ( y < y_south .OR. y > y_north )
1091                   boundary_particle = boundary_particle .OR. ( z > zw(nzt) )
[2801]1092
[4649]1093                   IF ( boundary_particle )  THEN
[3883]1094                      ic = x / pg%dx                     !TODO anpassen auf Mehrfachnesting
1095                      jc = y / pg%dy
[2801]1096
[4649]1097                      IF ( ic >= ipl  .AND.  ic <= ipr  .AND.  jc >= jps  .AND.  jc <= jpn )  THEN
1098                         coarse_particles(jc,ic)%nr_particle = coarse_particles(jc,ic)%nr_particle &
1099                                                               + 1
1100                         CALL check_and_alloc_coarse_particle( ic, jc,                             &
1101                                                               coarse_particles(jc,ic)%nr_particle,&
1102                                                               with_copy = .TRUE. )
[2801]1103
[4649]1104                         coarse_particles(jc,ic)%parent_particles(                                 &
1105                         coarse_particles(jc,ic)%nr_particle) = particles(n)
1106                         coarse_particles(jc,ic)%parent_particles(                                 &
1107                         coarse_particles(jc,ic)%nr_particle)%x = x   !Adapt to global coordinates
1108                         coarse_particles(jc,ic)%parent_particles(                                 &
1109                         coarse_particles(jc,ic)%nr_particle)%y = y
1110                         coarse_particles(jc,ic)%parent_particles(                                 &
1111                         coarse_particles(jc,ic)%nr_particle)%origin_x = xo
1112                         coarse_particles(jc,ic)%parent_particles(                                 &
1113                         coarse_particles(jc,ic)%nr_particle)%origin_y = yo
1114!
[2801]1115!--                      Mark particle as deleted after copying it to the transfer buffer
1116                         grid_particles(k,j,i)%particles(n)%particle_mask = .FALSE.
1117                      ELSE
[4649]1118                         WRITE( 9, '(a,10i6)' ) 'This should not happen ', i, j, k, ic, jc, ipl,   &
1119                                                                           ipr, jps, jpn
[2801]1120                         CALL MPI_Abort( MPI_COMM_WORLD, 9999, ierr )
1121                      ENDIF
1122                   ENDIF
1123                ENDIF
1124             ENDDO
1125          ENDDO
1126       ENDDO
1127    ENDDO
1128
1129!
[4649]1130!-- Pack particles (eliminate those marked for deletion), determine new number of particles
[4017]1131!    CALL lpm_sort_in_subboxes
[2801]1132
[2967]1133#endif
[2801]1134 END SUBROUTINE c_copy_particle_to_coarse_grid
[2967]1135
1136
[4649]1137!--------------------------------------------------------------------------------------------------!
[2801]1138! Description:
1139! ------------
[4649]1140!> Parent routine:
1141!> Copy/sort particles from the MPI window into the respective grid boxes
1142!--------------------------------------------------------------------------------------------------!
[3241]1143 SUBROUTINE p_copy_particle_to_org_grid( m )
[2801]1144
1145    IMPLICIT NONE
1146
[4649]1147    INTEGER(iwp) ::  i       !< loop index (x grid)
1148    INTEGER(iwp) ::  j       !< loop index (y grid)
1149    INTEGER(iwp) ::  k       !< loop index (z grid)
1150    INTEGER(iwp) ::  n       !< loop index (nr part)
1151    INTEGER(iwp) ::  nr      !< short variable name for nr_part
1152    INTEGER(iwp) ::  pindex  !< short variable name part_adr
[2801]1153
[4649]1154    INTEGER(iwp),INTENT(IN) ::  m  !<
[2801]1155
[4649]1156    INTEGER(iwp),DIMENSION(1) ::  buf_shape  !<
1157
1158    REAL(wp) ::  x   !< x coordinate
1159    REAL(wp) ::  xo  !< x origin
1160    REAL(wp) ::  y   !< y coordinate
1161    REAL(wp) ::  yo  !< y origin
1162    REAL(wp) ::  z   !< z coordinate
1163
1164
[2967]1165#if defined( __parallel )
[2801]1166    buf_shape(1) = max_nr_particle_in_rma_win
1167    CALL C_F_POINTER( buf_ptr(m), particle_in_win , buf_shape )
1168
[4649]1169    DO  i = nxl, nxr
1170       DO  j = nys, nyn
[2801]1171          nr = nr_part(j,i)
[4649]1172          IF ( nr > 0 )  THEN
[2801]1173             pindex = part_adr(j,i)
[4649]1174             DO  n = 1, nr
[2801]1175                x = particle_in_win(pindex)%x-lower_left_coord_x
1176                y = particle_in_win(pindex)%y-lower_left_coord_y
1177                z = particle_in_win(pindex)%z
1178                xo = particle_in_win(pindex)%origin_x-lower_left_coord_x
1179                yo = particle_in_win(pindex)%origin_y-lower_left_coord_y
1180                k = nzt + 1
[4649]1181                DO WHILE ( zw(k-1) > z .AND. k > nzb + 1 )  ! kk search loop has to be optimzed !!!
[2801]1182                   k = k - 1
1183                END DO
1184
1185                prt_count(k,j,i) = prt_count(k,j,i) + 1
[4649]1186                IF ( prt_count(k,j,i) > SIZE( grid_particles(k,j,i)%particles ) )  THEN
[4017]1187                   CALL pmc_realloc_particles_array( i, j, k )
[2801]1188                ENDIF
1189                grid_particles(k,j,i)%particles(prt_count(k,j,i)) = particle_in_win(pindex)
[4649]1190
[2801]1191!
1192!--             Update particle positions and origins relative to parent domain
1193                grid_particles(k,j,i)%particles(prt_count(k,j,i))%x = x
1194                grid_particles(k,j,i)%particles(prt_count(k,j,i))%y = y
1195                grid_particles(k,j,i)%particles(prt_count(k,j,i))%origin_x = xo
1196                grid_particles(k,j,i)%particles(prt_count(k,j,i))%origin_y = yo
1197                pindex = pindex + 1
1198             ENDDO
1199          ENDIF
1200       ENDDO
1201    ENDDO
[2967]1202
1203#endif
[2801]1204 END SUBROUTINE p_copy_particle_to_org_grid
[4649]1205
1206!--------------------------------------------------------------------------------------------------!
[4017]1207! Description:
1208! ------------
[4649]1209!> If the allocated memory for the particle array do not suffice to add arriving particles from
1210!> neighbour grid cells, this subrouting reallocates the particle array to assure enough memory is
1211!> available.
1212!--------------------------------------------------------------------------------------------------!
1213 SUBROUTINE pmc_realloc_particles_array( i, j, k, size_in )
[2801]1214
1215
[4649]1216    INTEGER(iwp) ::  new_size  !<
1217    INTEGER(iwp) ::  old_size  !<
[4017]1218
[4649]1219    INTEGER(iwp), INTENT(IN) ::  i  !<
1220    INTEGER(iwp), INTENT(IN) ::  j  !<
1221    INTEGER(iwp), INTENT(IN) ::  k  !<
[4017]1222
[4649]1223    INTEGER(iwp), INTENT(IN), OPTIONAL ::  size_in  !<
1224
1225
1226    TYPE(particle_type), DIMENSION(:), ALLOCATABLE ::  tmp_particles_d  !<
1227
1228    TYPE(particle_type), DIMENSION(500)            ::  tmp_particles_s  !<
1229
1230    old_size = SIZE( grid_particles(k,j,i)%particles )
1231
1232    IF ( PRESENT( size_in ) )  THEN
[4017]1233       new_size = size_in
1234    ELSE
1235       new_size = old_size * ( 1.0_wp + alloc_factor / 100.0_wp )
1236    ENDIF
1237
[4043]1238    new_size = MAX( new_size, 1, old_size + 1 )
[4017]1239
1240    IF ( old_size <= 500 )  THEN
1241
1242       tmp_particles_s(1:old_size) = grid_particles(k,j,i)%particles(1:old_size)
1243
[4649]1244       DEALLOCATE( grid_particles(k,j,i)%particles )
1245       ALLOCATE( grid_particles(k,j,i)%particles(new_size) )
[4017]1246
1247       grid_particles(k,j,i)%particles(1:old_size)          = tmp_particles_s(1:old_size)
1248       grid_particles(k,j,i)%particles(old_size+1:new_size) = zero_particle
1249
1250    ELSE
1251
[4649]1252       ALLOCATE( tmp_particles_d(new_size) )
[4017]1253       tmp_particles_d(1:old_size) = grid_particles(k,j,i)%particles
1254
[4649]1255       DEALLOCATE( grid_particles(k,j,i)%particles )
1256       ALLOCATE( grid_particles(k,j,i)%particles(new_size) )
[4017]1257
1258       grid_particles(k,j,i)%particles(1:old_size)          = tmp_particles_d(1:old_size)
1259       grid_particles(k,j,i)%particles(old_size+1:new_size) = zero_particle
1260
[4649]1261       DEALLOCATE( tmp_particles_d )
[4017]1262
1263    ENDIF
1264    particles => grid_particles(k,j,i)%particles(1:new_size)
1265
1266    RETURN
[4649]1267
[4017]1268 END SUBROUTINE pmc_realloc_particles_array
1269
[4444]1270#endif
[2801]1271END MODULE pmc_particle_interface
Note: See TracBrowser for help on using the repository browser.