comm_particles Subroutine

public subroutine comm_particles()

Uses

  • proc~~comm_particles~~UsesGraph proc~comm_particles comm_particles module~mod_global_parameters mod_global_parameters proc~comm_particles->module~mod_global_parameters module~mod_basic_types mod_basic_types module~mod_global_parameters->module~mod_basic_types module~mod_connectivity mod_connectivity module~mod_global_parameters->module~mod_connectivity module~mod_constants mod_constants module~mod_global_parameters->module~mod_constants module~mod_physicaldata mod_physicaldata module~mod_global_parameters->module~mod_physicaldata module~mod_variables mod_variables module~mod_global_parameters->module~mod_variables mpi mpi module~mod_global_parameters->mpi module~mod_variables->module~mod_basic_types

Arguments

None

Calls

proc~~comm_particles~~CallsGraph proc~comm_particles comm_particles mpi_waitall mpi_waitall proc~comm_particles->mpi_waitall proc~apply_periodb apply_periodB proc~comm_particles->proc~apply_periodb proc~destroy_particles destroy_particles proc~comm_particles->proc~destroy_particles proc~find_particle_ipe find_particle_ipe proc~comm_particles->proc~find_particle_ipe proc~mpi_irecv_wrapper mpi_irecv_wrapper proc~comm_particles->proc~mpi_irecv_wrapper proc~mpi_isend_wrapper mpi_isend_wrapper proc~comm_particles->proc~mpi_isend_wrapper proc~mpistop mpistop proc~comm_particles->proc~mpistop proc~particle_in_igrid particle_in_igrid proc~comm_particles->proc~particle_in_igrid proc~pull_particle_from_particles_on_mype pull_particle_from_particles_on_mype proc~comm_particles->proc~pull_particle_from_particles_on_mype proc~push_particle_into_particles_on_mype push_particle_into_particles_on_mype proc~comm_particles->proc~push_particle_into_particles_on_mype proc~apply_periodb->proc~find_particle_ipe proc~destroy_particles->proc~pull_particle_from_particles_on_mype mpi_allreduce mpi_allreduce proc~destroy_particles->mpi_allreduce proc~output_ensemble output_ensemble proc~destroy_particles->proc~output_ensemble proc~get_igslice get_igslice proc~find_particle_ipe->proc~get_igslice proc~particle_in_domain particle_in_domain proc~find_particle_ipe->proc~particle_in_domain mpi_irecv mpi_irecv proc~mpi_irecv_wrapper->mpi_irecv mpi_isend mpi_isend proc~mpi_isend_wrapper->mpi_isend mpi_abort mpi_abort proc~mpistop->mpi_abort proc~get_igslice->proc~mpistop mpi_allgather mpi_allgather proc~output_ensemble->mpi_allgather proc~mpi_recv_wrapper mpi_recv_wrapper proc~output_ensemble->proc~mpi_recv_wrapper proc~mpi_send_wrapper mpi_send_wrapper proc~output_ensemble->proc~mpi_send_wrapper proc~output_particle output_particle proc~output_ensemble->proc~output_particle mpi_recv mpi_recv proc~mpi_recv_wrapper->mpi_recv mpi_send mpi_send proc~mpi_send_wrapper->mpi_send proc~partcoord_to_cartesian partcoord_to_cartesian proc~output_particle->proc~partcoord_to_cartesian

Called by

proc~~comm_particles~~CalledByGraph proc~comm_particles comm_particles proc~advance_particles advance_particles proc~advance_particles->proc~comm_particles proc~handle_particles handle_particles proc~handle_particles->proc~advance_particles proc~advance advance proc~advance->proc~handle_particles proc~main main proc~main->proc~handle_particles proc~timeintegration timeintegration proc~main->proc~timeintegration proc~timeintegration->proc~advance program~agile agile program~agile->proc~main

Source Code

  subroutine comm_particles()

    use mod_global_parameters

    integer                         :: ipart, iipart, igrid_particle,&
        ipe_particle, ipe, iipe
    integer                         :: index
    integer                         :: tag_send, tag_receive, send_buff,&
        rcv_buff
    integer                         :: status(MPI_STATUS_SIZE)
    integer, dimension(0:npe-1)     :: send_n_particles_to_ipe
    integer, dimension(0:npe-1)    :: receive_n_particles_from_ipe
    type(particle_t), allocatable, dimension(:,:)  :: send_particles
    type(particle_t), allocatable, dimension(:,:)  :: receive_particles
    double precision, allocatable, dimension(:,:,:)  :: send_payload
    double precision, allocatable, dimension(:,:,:)  :: receive_payload
    integer, allocatable, dimension(:,:)      :: &
       particle_index_to_be_sent_to_ipe
    integer, dimension(nparticles_per_cpu_hi) :: &
       particle_index_to_be_destroyed
    integer                                   :: destroy_n_particles_mype
    logical                                   :: BC_applied
    integer, allocatable, dimension(:)        :: sndrqst, rcvrqst
    integer, allocatable, dimension(:)        :: sndrqst_payload,&
        rcvrqst_payload
    integer                                   :: isnd, ircv
    integer, parameter                        :: maxneighbors=56 !maximum case: coarse-fine in 3D
    !-----------------------------------------------------------------------------

    send_n_particles_to_ipe(:)      = 0
    receive_n_particles_from_ipe(:) = 0
    destroy_n_particles_mype        = 0

    allocate( particle_index_to_be_sent_to_ipe(nparticles_per_cpu_hi,&
       0:npe-1) )
    allocate( sndrqst(1:npe-1), rcvrqst(1:npe-1) )
    allocate( sndrqst_payload(1:npe-1), rcvrqst_payload(1:npe-1) )
    sndrqst = MPI_REQUEST_NULL; rcvrqst = MPI_REQUEST_NULL;
    sndrqst_payload = MPI_REQUEST_NULL; rcvrqst_payload = MPI_REQUEST_NULL;

    ! check if and where to send each particle, destroy if necessary
    do iipart=1,nparticles_on_mype;ipart=particles_on_mype(iipart);

       ! first check if the particle should be destroyed (user defined criterion)
       if (associated(usr_destroy_particle)) then
         if (usr_destroy_particle(particle(ipart))) then
           destroy_n_particles_mype  = destroy_n_particles_mype + 1
           particle_index_to_be_destroyed(destroy_n_particles_mype) = ipart
           cycle
         end if
       end if

       ! is my particle still in the same igrid?
       if (.not.particle_in_igrid(ipart,particle(ipart)%igrid)) then
          call find_particle_ipe(particle(ipart)%self%x,igrid_particle,&
             ipe_particle)

          ! destroy particle if out of domain (signalled by return value of -1)
          if (igrid_particle == -1 )then
             call apply_periodB(particle(ipart)%self,igrid_particle,&
                ipe_particle,BC_applied)
             if (.not. BC_applied .or. igrid_particle == -1) then
                destroy_n_particles_mype  = destroy_n_particles_mype + 1
                particle_index_to_be_destroyed(destroy_n_particles_mype) = &
                   ipart
                cycle
             end if
          end if

          ! particle still present
          particle(ipart)%igrid = igrid_particle
          particle(ipart)%ipe = ipe_particle

          ! if we have more than one core, is it on another cpu?
          if (npe .gt. 1 .and. particle(ipart)%ipe .ne. mype) then
             send_n_particles_to_ipe(ipe_particle) = &
                send_n_particles_to_ipe(ipe_particle) + 1
             if (send_n_particles_to_ipe(ipe_particle) .gt. &
                nparticles_per_cpu_hi) call &
                mpistop('too many particles, increase nparticles_per_cpu_hi')
             particle_index_to_be_sent_to_ipe(send_n_particles_to_ipe(&
                ipe_particle),ipe_particle) = ipart
          end if ! ipe_particle

       end if ! particle_in_grid

    end do ! ipart

    call destroy_particles(destroy_n_particles_mype,&
       particle_index_to_be_destroyed(1:destroy_n_particles_mype))

    ! get out when only one core:
    if (npe == 1) return

    ! communicate amount of particles to be sent/received
    isnd = 0; ircv = 0;
    do iipe=1,npe_neighbors;ipe=ipe_neighbor(iipe);
       tag_send    = mype * npe + ipe
       tag_receive = ipe * npe + mype
       isnd = isnd + 1; ircv = ircv + 1
       call mpi_isend_wrapper(send_n_particles_to_ipe(ipe),1,MPI_INTEGER,ipe,tag_send,&
          icomm,sndrqst(isnd),ierrmpi)
       call mpi_irecv_wrapper(receive_n_particles_from_ipe(ipe),1,MPI_INTEGER,ipe,&
          tag_receive,icomm,rcvrqst(ircv),ierrmpi)
    end do

    call MPI_WAITALL(isnd,sndrqst,MPI_STATUSES_IGNORE,ierrmpi)
    call MPI_WAITALL(ircv,rcvrqst,MPI_STATUSES_IGNORE,ierrmpi)
    sndrqst = MPI_REQUEST_NULL; rcvrqst = MPI_REQUEST_NULL;

    ! allocate the send and receive buffers,
    ! If this is still not good enough, consider optimizing using lists of lists:
    allocate(send_particles(maxval(send_n_particles_to_ipe),1:min(maxneighbors,&
       npe-1)))
    allocate(receive_particles(maxval(receive_n_particles_from_ipe),&
       1:min(maxneighbors,npe-1)))
    allocate(send_payload(npayload,maxval(send_n_particles_to_ipe),&
       1:min(maxneighbors,npe-1)))
    allocate(receive_payload(npayload,maxval(receive_n_particles_from_ipe),&
       1:min(maxneighbors,npe-1)))

    isnd = 0; ircv = 0;
    ! send and receive the data of the particles
    do iipe=1,npe_neighbors;ipe=ipe_neighbor(iipe);
       tag_send    = mype * npe + ipe
       tag_receive = ipe * npe + mype

       ! should i send some particles to ipe?
       if (send_n_particles_to_ipe(ipe) .gt. 0) then

          ! create the send buffer
          do ipart = 1, send_n_particles_to_ipe(ipe)
             send_particles(ipart,iipe) = &
                particle(particle_index_to_be_sent_to_ipe(ipart,ipe))%self
             send_payload(1:npayload,ipart,&
                iipe) = particle(particle_index_to_be_sent_to_ipe(ipart,&
                ipe))%payload(1:npayload)
          end do ! ipart
          isnd = isnd + 1
          call mpi_isend_wrapper(send_particles(:,iipe),send_n_particles_to_ipe(ipe),&
             type_particle,ipe,tag_send,icomm,sndrqst(isnd),ierrmpi)
          call mpi_isend_wrapper(send_payload(:,:,iipe),&
             npayload*send_n_particles_to_ipe(ipe),MPI_DOUBLE_PRECISION,ipe,&
             tag_send,icomm,sndrqst_payload(isnd),ierrmpi)

       end if ! send .gt. 0

       ! should i receive some particles from ipe?
       if (receive_n_particles_from_ipe(ipe) .gt. 0) then

          ircv = ircv + 1
          call mpi_irecv_wrapper(receive_particles(:,iipe),&
             receive_n_particles_from_ipe(ipe),type_particle,ipe,tag_receive,&
             icomm,rcvrqst(ircv),ierrmpi)
          call mpi_irecv_wrapper(receive_payload(:,:,iipe),&
             npayload*receive_n_particles_from_ipe(ipe),MPI_DOUBLE_PRECISION,&
             ipe,tag_receive,icomm,rcvrqst_payload(ircv),ierrmpi)

       end if ! receive .gt. 0
    end do ! ipe

    call MPI_WAITALL(isnd,sndrqst,MPI_STATUSES_IGNORE,ierrmpi)
    call MPI_WAITALL(ircv,rcvrqst,MPI_STATUSES_IGNORE,ierrmpi)
    call MPI_WAITALL(isnd,sndrqst_payload,MPI_STATUSES_IGNORE,ierrmpi)
    call MPI_WAITALL(ircv,rcvrqst_payload,MPI_STATUSES_IGNORE,ierrmpi)

    do iipe=1,npe_neighbors;ipe=ipe_neighbor(iipe);
       do ipart = 1, send_n_particles_to_ipe(ipe)
          deallocate(particle(particle_index_to_be_sent_to_ipe(ipart,&
             ipe))%self)
          deallocate(particle(particle_index_to_be_sent_to_ipe(ipart,&
             ipe))%payload)
          call pull_particle_from_particles_on_mype(&
             particle_index_to_be_sent_to_ipe(ipart,ipe))
       end do ! ipart
    end do

    do iipe=1,npe_neighbors;ipe=ipe_neighbor(iipe);
       if (receive_n_particles_from_ipe(ipe) .gt. 0) then
          do ipart = 1, receive_n_particles_from_ipe(ipe)

             index = receive_particles(ipart,iipe)%index
             call push_particle_into_particles_on_mype(index)
             !if (.not. allocated(particle(index)%self)) &
               allocate(particle(index)%self)
             particle(index)%self = receive_particles(ipart,iipe)
             !if (.not. allocated(particle(index)%payload)) &
               allocate(particle(index)%payload(npayload))
             particle(index)%payload(1:npayload) = receive_payload(1:npayload,&
                ipart,iipe)

             ! since we don't send the igrid, need to re-locate it
             call find_particle_ipe(particle(index)%self%x,igrid_particle,&
                ipe_particle)
             particle(index)%igrid = igrid_particle
             particle(index)%ipe = ipe_particle

          end do ! ipart

       end if ! receive .gt. 0
    end do ! ipe

  end subroutine comm_particles