From: Bruno Date: Wed, 4 Nov 2020 17:21:19 +0000 (-0500) Subject: Introduction of ghost particles cache and new version to send X-Git-Tag: v9.3.0-rc1~922^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=1579a0f17e0d6a7b15e3e36090ddd2bc450976f2;p=dealii.git Introduction of ghost particles cache and new version to send information on ghost partilces without rebuilding them Added a test with asymmetric ghosts because there was an error in the previous implementation and this test enabled me to find it Added changelog --- diff --git a/doc/news/changes/major/20201110Blais b/doc/news/changes/major/20201110Blais new file mode 100644 index 0000000000..58004d83b8 --- /dev/null +++ b/doc/news/changes/major/20201110Blais @@ -0,0 +1,3 @@ +New: Added capacity to update ghost particles without rebuilding them from scratch in the particle_handler +
+(Bruno Blais, Peter Munch, 2020/11/10) diff --git a/include/deal.II/base/mpi_tags.h b/include/deal.II/base/mpi_tags.h index 37b0ff45f0..6a8e19d201 100644 --- a/include/deal.II/base/mpi_tags.h +++ b/include/deal.II/base/mpi_tags.h @@ -109,6 +109,8 @@ namespace Utilities /// ParticleHandler::send_recv_particles particle_handler_send_recv_particles_setup, /// ParticleHandler::send_recv_particles + particle_handler_send_recv_particles_cache_setup, + /// ParticleHandler::send_recv_particles particle_handler_send_recv_particles_send, /// ScaLAPACKMatrix::copy_to diff --git a/include/deal.II/particles/particle.h b/include/deal.II/particles/particle.h index 4466e79b45..a79617d223 100644 --- a/include/deal.II/particles/particle.h +++ b/include/deal.II/particles/particle.h @@ -416,6 +416,12 @@ namespace Particles void load(Archive &ar, const unsigned int version); + /** + * Free the memory of the property pool + */ + void + free_properties(); + #ifdef DOXYGEN /** * Write and read the data of this object from a stream for the purpose diff --git a/include/deal.II/particles/particle_handler.h b/include/deal.II/particles/particle_handler.h index bce703bb59..da755af1f1 100644 --- a/include/deal.II/particles/particle_handler.h +++ b/include/deal.II/particles/particle_handler.h @@ -705,7 +705,7 @@ namespace Particles * member variable. */ void - exchange_ghost_particles(); + exchange_ghost_particles(const bool enable_ghost_cache = false); /** * Update all particles that live in cells that are ghost cells to @@ -819,14 +819,6 @@ namespace Particles */ std::multimap> ghost_particles; - /** - * Set of particles that currently live in the ghost cells of the local - * domain, organized by the subdomain_id. These - * particles are equivalent to the ghost entries in distributed vectors. - */ - std::map> - ghost_particles_by_domain; - /** * This variable stores how many particles are stored globally. It is * calculated by update_cached_numbers(). @@ -927,6 +919,13 @@ namespace Particles * particle to be send in which the particle belongs. This parameter * is necessary if the cell information of the particle iterator is * outdated (e.g. after particle movement). + * + * @param [in] enable_cache Optional bool that enables updating + * the ghost particles without rebuilding them from scratch by + * building a cache of type GhostParticlePartitioner which + * stores the necessary information to update the ghost particles. + * Once this cache is built, the ghost particles can be updated + * by a call to send_recv_particles_properties_and_location(). */ void send_recv_particles( @@ -941,9 +940,111 @@ namespace Particles &new_cells_for_particles = std::map< types::subdomain_id, std::vector< - typename Triangulation::active_cell_iterator>>()); + typename Triangulation::active_cell_iterator>>(), + const bool enable_cache = false); + + /** + * Transfer particles position and properties assuming that + * the particles have not changed cells. This routine uses the + * GhostParticlePartitioner as a caching structure to update the particles. + * It inherently assumes that particles cannot have changed cell. + * All updated particles will be appended to the + * @p received_particles vector. + * + * @param [in] particles_to_send All particles for which information + * should be sent and their new subdomain_ids are in this map. + * + * @param [in,out] received_particles A map with all received + * particles. Note that it is not required nor checked that the list + * is empty, received particles are simply attached to the end of + * the vector. + * + */ + void + send_recv_particles_properties_and_location( + const std::map> + &particles_to_send, + std::multimap> + &received_particles); + + #endif + + /** + * Cache structure used to store the elements which are required to + * exchange the particle information (location and properties) accross + * processors in order to update the ghost particles. + * + * This structure should only be used when one wishes to carry out work + * using the particles without calling + * sort_particles_into_subdomain_and_cells at every iteration. This is + * useful when particle-particle interaction occur at a different time + * scale than particle-FEM interaction. + * + * This structure is similar to Utilities::MPI::Partitioner::import_targets + * when combined with neighbors. + */ + struct GhostParticlePartitioner + { + /** + * Indicates if the cache has been built to prevent updating particles + * with an invalid cache. + */ + bool valid = false; + + /** + * Vector of the subdomain id of all possible neighbors of the current + * subdomain. + */ + std::vector neighbors; + + /** + * Vector of size (neighbors.size()+1) used to store the start and the + * end point of the data that must go from the current subdomain to the + * neighbors. For neighbor i, send_pointers[i] indicates the beginning + * and send_pointers[i+1] indicates the end of the data that must be + * sent. + */ + std::vector send_pointers; + + /** + * Set of particles that currently live in the ghost cells of the local + * domain, organized by the subdomain_id. These + * particles are equivalent to the ghost entries in distributed vectors. + */ + std::map> + ghost_particles_by_domain; + + /** + * Vector of size (neighbors.size()+1) used to store the start and the + * end point of the data that must be received from neighbor[i] on + * the current subdomain. For neighbor i, recv_pointers[i] indicate the + * beggining and reicv_pointers[i+1] indicates the end of the data that + * must be received. + */ + std::vector recv_pointers; + + /** + * Vector of ghost particles in the order in which they are inserted + * in the multimap used to store particles on the triangulation. This + * information is used to update the ghost particle information + * without clearing the multimap of ghost particles, thus greatly + * reducing the cost of exchanging the ghost particles information. + */ + std::vector>::iterator> + ghost_particles_iterators; + }; + + /** + * Cache structure used to store the elements which are required to + * exchange the particle information (location and properties) accross + * processors in order to update the ghost particles. This structure + * is only used to update the ghost particles. + */ + GhostParticlePartitioner ghost_particles_cache; + /** * Called by listener functions from Triangulation for every cell * before a refinement step. All particles have to be attached to their diff --git a/source/particles/particle.cc b/source/particles/particle.cc index 69f16624cd..2befeec29c 100644 --- a/source/particles/particle.cc +++ b/source/particles/particle.cc @@ -180,6 +180,14 @@ namespace Particles property_pool->deallocate_properties_array(properties); } + template + void + Particle::free_properties() + { + if (property_pool != nullptr && properties != PropertyPool::invalid_handle) + property_pool->deallocate_properties_array(properties); + } + template diff --git a/source/particles/particle_handler.cc b/source/particles/particle_handler.cc index 56094fb7c4..783bd95091 100644 --- a/source/particles/particle_handler.cc +++ b/source/particles/particle_handler.cc @@ -90,7 +90,6 @@ namespace Particles , property_pool(std::make_unique(0)) , particles() , ghost_particles() - , ghost_particles_by_domain() , global_number_of_particles(0) , global_max_particles_per_cell(0) , next_free_particle_index(0) @@ -112,7 +111,6 @@ namespace Particles , property_pool(std::make_unique(n_properties)) , particles() , ghost_particles() - , ghost_particles_by_domain() , global_number_of_particles(0) , global_max_particles_per_cell(0) , next_free_particle_index(0) @@ -166,11 +164,13 @@ namespace Particles global_number_of_particles = particle_handler.global_number_of_particles; global_max_particles_per_cell = particle_handler.global_max_particles_per_cell; - next_free_particle_index = particle_handler.next_free_particle_index; - particles = particle_handler.particles; - ghost_particles = particle_handler.ghost_particles; - ghost_particles_by_domain = particle_handler.ghost_particles_by_domain; - handle = particle_handler.handle; + next_free_particle_index = particle_handler.next_free_particle_index; + particles = particle_handler.particles; + ghost_particles = particle_handler.ghost_particles; + + ghost_particles_cache.ghost_particles_by_domain = + particle_handler.ghost_particles_cache.ghost_particles_by_domain; + handle = particle_handler.handle; // copy dynamic properties auto from_particle = particle_handler.begin(); @@ -1160,7 +1160,8 @@ namespace Particles template void - ParticleHandler::exchange_ghost_particles() + ParticleHandler::exchange_ghost_particles( + const bool enable_cache) { // Nothing to do in serial computations const auto parallel_triangulation = @@ -1175,17 +1176,21 @@ namespace Particles else return; -#ifdef DEAL_II_WITH_MPI +#ifndef DEAL_II_WITH_MPI + (void)enable_cache; +#else // First clear the current ghost_particle information ghost_particles.clear(); - ghost_particles_by_domain.clear(); + // Clear ghost particles data structures and invalidate cache + ghost_particles_cache.ghost_particles_by_domain.clear(); + ghost_particles_cache.valid = false; const std::set ghost_owners = parallel_triangulation->ghost_owners(); for (const auto ghost_owner : ghost_owners) - ghost_particles_by_domain[ghost_owner].reserve( + ghost_particles_cache.ghost_particles_by_domain[ghost_owner].reserve( static_cast::size_type>( particles.size() * 0.25)); @@ -1215,13 +1220,21 @@ namespace Particles particle_range.begin(); particle != particle_range.end(); ++particle) - ghost_particles_by_domain[domain].push_back(particle); + ghost_particles_cache.ghost_particles_by_domain[domain] + .push_back(particle); } } } } - send_recv_particles(ghost_particles_by_domain, ghost_particles); + send_recv_particles( + ghost_particles_cache.ghost_particles_by_domain, + ghost_particles, + std::map< + types::subdomain_id, + std::vector< + typename Triangulation::active_cell_iterator>>(), + enable_cache); #endif } @@ -1240,11 +1253,18 @@ namespace Particles return; } + #ifdef DEAL_II_WITH_MPI // First clear the current ghost_particle information - ghost_particles.clear(); + // ghost_particles.clear(); + Assert( + ghost_particles_cache.valid, + ExcMessage( + "Ghost particles cannot be updated if they first have not been exchanged at least once with the cache enabled")); - send_recv_particles(ghost_particles_by_domain, ghost_particles); + + send_recv_particles_properties_and_location( + ghost_particles_cache.ghost_particles_by_domain, ghost_particles); #endif } @@ -1261,8 +1281,11 @@ namespace Particles const std::map< types::subdomain_id, std::vector::active_cell_iterator>> - &send_cells) + & send_cells, + const bool build_cache) { + ghost_particles_cache.valid = build_cache; + const auto parallel_triangulation = dynamic_cast *>( &*triangulation); @@ -1295,7 +1318,7 @@ namespace Particles Assert(ghost_owners.find(send_particles->first) != ghost_owners.end(), ExcNotImplemented()); - unsigned int n_send_particles = 0; + std::size_t n_send_particles = 0; for (auto send_particles = particles_to_send.begin(); send_particles != particles_to_send.end(); ++send_particles) @@ -1309,16 +1332,25 @@ namespace Particles std::vector send_offsets(n_neighbors, 0); std::vector send_data; + const unsigned int individual_particle_data_size = + Utilities::MPI::max(n_send_particles > 0 ? + ((begin()->serialized_size_in_bytes() + + (size_callback ? size_callback() : 0))) : + 0, + parallel_triangulation->get_communicator()); + + const unsigned int individual_total_particle_data_size = + individual_particle_data_size + cellid_size; + // Only serialize things if there are particles to be send. // We can not return early even if no particles // are send, because we might receive particles from other processes if (n_send_particles > 0) { // Allocate space for sending particle data - const unsigned int particle_size = - begin()->serialized_size_in_bytes() + cellid_size + - (size_callback ? size_callback() : 0); - send_data.resize(n_send_particles * particle_size); + send_data.resize(n_send_particles * + individual_total_particle_data_size); + void *data = static_cast(&send_data.front()); // Serialize the data sorted by receiving process @@ -1327,9 +1359,17 @@ namespace Particles send_offsets[i] = reinterpret_cast(data) - reinterpret_cast(&send_data.front()); - for (unsigned int j = 0; - j < particles_to_send.at(neighbors[i]).size(); - ++j) + const unsigned int n_particles_to_send = + particles_to_send.at(neighbors[i]).size(); + + Assert(static_cast(n_particles_to_send) * + individual_total_particle_data_size == + static_cast( + n_particles_to_send * + individual_total_particle_data_size), + ExcMessage("Overflow when trying to send particle data")); + + for (unsigned int j = 0; j < n_particles_to_send; ++j) { // If no target cells are given, use the iterator information typename Triangulation::active_cell_iterator @@ -1351,9 +1391,7 @@ namespace Particles data = store_callback(particles_to_send.at(neighbors[i])[j], data); } - n_send_data[i] = reinterpret_cast(data) - - send_offsets[i] - - reinterpret_cast(&send_data.front()); + n_send_data[i] = n_particles_to_send; } } @@ -1361,7 +1399,6 @@ namespace Particles std::vector n_recv_data(n_neighbors); std::vector recv_offsets(n_neighbors); - // Notify other processors how many particles we will send { const int mpi_tag = Utilities::MPI::internal::Tags:: particle_handler_send_recv_particles_setup; @@ -1399,7 +1436,8 @@ namespace Particles for (unsigned int neighbor_id = 0; neighbor_id < n_neighbors; ++neighbor_id) { recv_offsets[neighbor_id] = total_recv_data; - total_recv_data += n_recv_data[neighbor_id]; + total_recv_data += + n_recv_data[neighbor_id] * individual_total_particle_data_size; } // Set up the space for the received particle data @@ -1419,7 +1457,7 @@ namespace Particles { const int ierr = MPI_Irecv(&(recv_data[recv_offsets[i]]), - n_recv_data[i], + n_recv_data[i] * individual_total_particle_data_size, MPI_CHAR, neighbors[i], mpi_tag, @@ -1434,7 +1472,7 @@ namespace Particles { const int ierr = MPI_Isend(&(send_data[send_offsets[i]]), - n_send_data[i], + n_send_data[i] * individual_total_particle_data_size, MPI_CHAR, neighbors[i], mpi_tag, @@ -1452,6 +1490,33 @@ namespace Particles // triangulation const void *recv_data_it = static_cast(recv_data.data()); + // Store the particle iterators in the cache + auto &ghost_particles_iterators = + ghost_particles_cache.ghost_particles_iterators; + + if (build_cache) + { + ghost_particles_iterators.clear(); + + auto &send_pointers_particles = ghost_particles_cache.send_pointers; + send_pointers_particles.assign(n_neighbors + 1, 0); + + for (unsigned int i = 0; i < n_neighbors; ++i) + send_pointers_particles[i + 1] = + send_pointers_particles[i] + + n_send_data[i] * individual_particle_data_size; + + auto &recv_pointers_particles = ghost_particles_cache.recv_pointers; + recv_pointers_particles.assign(n_neighbors + 1, 0); + + for (unsigned int i = 0; i < n_neighbors; ++i) + recv_pointers_particles[i + 1] = + recv_pointers_particles[i] + + n_recv_data[i] * individual_particle_data_size; + + ghost_particles_cache.neighbors = neighbors; + } + while (reinterpret_cast(recv_data_it) - reinterpret_cast(recv_data.data()) < total_recv_data) @@ -1474,6 +1539,9 @@ namespace Particles recv_data_it = load_callback(particle_iterator(received_particles, recv_particle), recv_data_it); + + if (build_cache) // TODO: is this safe? + ghost_particles_iterators.push_back(recv_particle); } AssertThrow(recv_data_it == recv_data.data() + recv_data.size(), @@ -1485,6 +1553,119 @@ namespace Particles +#ifdef DEAL_II_WITH_MPI + template + void + ParticleHandler::send_recv_particles_properties_and_location( + const std::map> + &particles_to_send, + std::multimap> + &updated_particles) + { + const auto &neighbors = ghost_particles_cache.neighbors; + const auto &send_pointers = ghost_particles_cache.send_pointers; + const auto &recv_pointers = ghost_particles_cache.recv_pointers; + + const auto parallel_triangulation = + dynamic_cast *>( + &*triangulation); + Assert( + parallel_triangulation, + ExcMessage( + "This function is only implemented for parallel::TriangulationBase objects.")); + + std::vector send_data; + + // Fill data to send + if (send_pointers.back() > 0) + { + // Allocate space for sending particle data + send_data.resize(send_pointers.back()); + void *data = static_cast(&send_data.front()); + + // Serialize the data sorted by receiving process + for (const auto i : neighbors) + for (const auto &p : particles_to_send.at(i)) + { + p->write_data(data); + if (store_callback) + data = store_callback(p, data); + } + } + + // Set up the space for the received particle data + std::vector recv_data(recv_pointers.back()); + + // Exchange the particle data between domains + { + std::vector requests(2 * neighbors.size()); + unsigned int send_ops = 0; + unsigned int recv_ops = 0; + + const int mpi_tag = Utilities::MPI::internal::Tags:: + particle_handler_send_recv_particles_send; + + for (unsigned int i = 0; i < neighbors.size(); ++i) + if ((recv_pointers[i + 1] - recv_pointers[i]) > 0) + { + const int ierr = + MPI_Irecv(recv_data.data() + recv_pointers[i], + recv_pointers[i + 1] - recv_pointers[i], + MPI_CHAR, + neighbors[i], + mpi_tag, + parallel_triangulation->get_communicator(), + &(requests[send_ops])); + AssertThrowMPI(ierr); + send_ops++; + } + + for (unsigned int i = 0; i < neighbors.size(); ++i) + if ((send_pointers[i + 1] - send_pointers[i]) > 0) + { + const int ierr = + MPI_Isend(send_data.data() + send_pointers[i], + send_pointers[i + 1] - send_pointers[i], + MPI_CHAR, + neighbors[i], + mpi_tag, + parallel_triangulation->get_communicator(), + &(requests[send_ops + recv_ops])); + AssertThrowMPI(ierr); + recv_ops++; + } + const int ierr = + MPI_Waitall(send_ops + recv_ops, requests.data(), MPI_STATUSES_IGNORE); + AssertThrowMPI(ierr); + } + + // Put the received particles into the domain if they are in the + // triangulation + const void *recv_data_it = static_cast(recv_data.data()); + + // Gather ghost particle iterators from the cache + auto &ghost_particles_iterators = + ghost_particles_cache.ghost_particles_iterators; + + for (auto &recv_particle : ghost_particles_iterators) + { + recv_particle->second.free_properties(); + recv_particle->second = + Particle(recv_data_it, property_pool.get()); + + if (load_callback) + recv_data_it = + load_callback(particle_iterator(updated_particles, recv_particle), + recv_data_it); + } + + AssertThrow(recv_data_it == recv_data.data() + recv_data.size(), + ExcMessage( + "The amount of data that was read into new particles " + "does not match the amount of data sent around.")); + } +#endif + template void ParticleHandler::register_additional_store_load_functions( diff --git a/tests/particles/particle_handler_19.cc b/tests/particles/particle_handler_19.cc index a0933e5c7e..e619384bf5 100644 --- a/tests/particles/particle_handler_19.cc +++ b/tests/particles/particle_handler_19.cc @@ -43,7 +43,7 @@ test() tr.refine_global(2); MappingQ mapping(1); - Particles::ParticleHandler particle_handler(tr, mapping, 1); + Particles::ParticleHandler particle_handler(tr, mapping, 2); Point position; Point reference_position; @@ -72,18 +72,20 @@ test() { particle->get_properties()[0] = 10 + Utilities::MPI::this_mpi_process(tr.get_communicator()); + particle->get_properties()[1] = + 100 + Utilities::MPI::this_mpi_process(tr.get_communicator()); } - particle_handler.exchange_ghost_particles(); + particle_handler.exchange_ghost_particles(true); for (auto particle = particle_handler.begin(); particle != particle_handler.end(); ++particle) deallog << "Particle id : " << particle->get_id() << " location : " << particle->get_location() - << " property : " << particle->get_properties()[0] - << " is local on process : " + << " property : " << particle->get_properties()[0] << " and " + << particle->get_properties()[1] << " is local on process : " << Utilities::MPI::this_mpi_process(tr.get_communicator()) << std::endl; @@ -92,8 +94,8 @@ test() ++particle) deallog << "Particle id : " << particle->get_id() << " location : " << particle->get_location() - << " property : " << particle->get_properties()[0] - << " is ghost on process : " + << " property : " << particle->get_properties()[0] << " and " + << particle->get_properties()[1] << " is ghost on process : " << Utilities::MPI::this_mpi_process(tr.get_communicator()) << std::endl; @@ -108,17 +110,21 @@ test() auto location = particle->get_location(); location[0] += 0.1; particle->get_properties()[0] += 10; + particle->get_properties()[1] += 100; particle->set_location(location); } + + // Update the ghost particles particle_handler.update_ghost_particles(); + for (auto particle = particle_handler.begin(); particle != particle_handler.end(); ++particle) deallog << "Particle id : " << particle->get_id() << " location : " << particle->get_location() - << " property : " << particle->get_properties()[0] - << " is local on process : " + << " property : " << particle->get_properties()[0] << " and " + << particle->get_properties()[1] << " is local on process : " << Utilities::MPI::this_mpi_process(tr.get_communicator()) << std::endl; @@ -127,8 +133,8 @@ test() ++particle) deallog << "Particle id : " << particle->get_id() << " location : " << particle->get_location() - << " property : " << particle->get_properties()[0] - << " is ghost on process : " + << " property : " << particle->get_properties()[0] << " and " + << particle->get_properties()[1] << " is ghost on process : " << Utilities::MPI::this_mpi_process(tr.get_communicator()) << std::endl; } @@ -137,7 +143,6 @@ test() } - int main(int argc, char *argv[]) { diff --git a/tests/particles/particle_handler_19.with_p4est=true.mpirun=2.output b/tests/particles/particle_handler_19.with_p4est=true.mpirun=2.output index 91fc304f0f..62c15f8023 100644 --- a/tests/particles/particle_handler_19.with_p4est=true.mpirun=2.output +++ b/tests/particles/particle_handler_19.with_p4est=true.mpirun=2.output @@ -1,39 +1,39 @@ -DEAL:0:2d/2d::Particle id : 0 location : 0.475000 0.475000 property : 10.0000 is local on process : 0 -DEAL:0:2d/2d::Particle id : 1 location : 0.525000 0.525000 property : 11.0000 is ghost on process : 0 +DEAL:0:2d/2d::Particle id : 0 location : 0.475000 0.475000 property : 10.0000 and 100.000 is local on process : 0 +DEAL:0:2d/2d::Particle id : 1 location : 0.525000 0.525000 property : 11.0000 and 101.000 is ghost on process : 0 DEAL:0:2d/2d::Modifying particles positions and properties -DEAL:0:2d/2d::Particle id : 0 location : 0.575000 0.475000 property : 20.0000 is local on process : 0 -DEAL:0:2d/2d::Particle id : 1 location : 0.625000 0.525000 property : 21.0000 is ghost on process : 0 +DEAL:0:2d/2d::Particle id : 0 location : 0.575000 0.475000 property : 20.0000 and 200.000 is local on process : 0 +DEAL:0:2d/2d::Particle id : 1 location : 0.625000 0.525000 property : 21.0000 and 201.000 is ghost on process : 0 DEAL:0:2d/2d::OK -DEAL:0:2d/3d::Particle id : 0 location : 0.475000 0.475000 0.00000 property : 10.0000 is local on process : 0 -DEAL:0:2d/3d::Particle id : 1 location : 0.525000 0.525000 0.00000 property : 11.0000 is ghost on process : 0 +DEAL:0:2d/3d::Particle id : 0 location : 0.475000 0.475000 0.00000 property : 10.0000 and 100.000 is local on process : 0 +DEAL:0:2d/3d::Particle id : 1 location : 0.525000 0.525000 0.00000 property : 11.0000 and 101.000 is ghost on process : 0 DEAL:0:2d/3d::Modifying particles positions and properties -DEAL:0:2d/3d::Particle id : 0 location : 0.575000 0.475000 0.00000 property : 20.0000 is local on process : 0 -DEAL:0:2d/3d::Particle id : 1 location : 0.625000 0.525000 0.00000 property : 21.0000 is ghost on process : 0 +DEAL:0:2d/3d::Particle id : 0 location : 0.575000 0.475000 0.00000 property : 20.0000 and 200.000 is local on process : 0 +DEAL:0:2d/3d::Particle id : 1 location : 0.625000 0.525000 0.00000 property : 21.0000 and 201.000 is ghost on process : 0 DEAL:0:2d/3d::OK -DEAL:0:3d/3d::Particle id : 0 location : 0.475000 0.475000 0.475000 property : 10.0000 is local on process : 0 -DEAL:0:3d/3d::Particle id : 1 location : 0.525000 0.525000 0.525000 property : 11.0000 is ghost on process : 0 +DEAL:0:3d/3d::Particle id : 0 location : 0.475000 0.475000 0.475000 property : 10.0000 and 100.000 is local on process : 0 +DEAL:0:3d/3d::Particle id : 1 location : 0.525000 0.525000 0.525000 property : 11.0000 and 101.000 is ghost on process : 0 DEAL:0:3d/3d::Modifying particles positions and properties -DEAL:0:3d/3d::Particle id : 0 location : 0.575000 0.475000 0.475000 property : 20.0000 is local on process : 0 -DEAL:0:3d/3d::Particle id : 1 location : 0.625000 0.525000 0.525000 property : 21.0000 is ghost on process : 0 +DEAL:0:3d/3d::Particle id : 0 location : 0.575000 0.475000 0.475000 property : 20.0000 and 200.000 is local on process : 0 +DEAL:0:3d/3d::Particle id : 1 location : 0.625000 0.525000 0.525000 property : 21.0000 and 201.000 is ghost on process : 0 DEAL:0:3d/3d::OK -DEAL:1:2d/2d::Particle id : 1 location : 0.525000 0.525000 property : 11.0000 is local on process : 1 -DEAL:1:2d/2d::Particle id : 0 location : 0.475000 0.475000 property : 10.0000 is ghost on process : 1 +DEAL:1:2d/2d::Particle id : 1 location : 0.525000 0.525000 property : 11.0000 and 101.000 is local on process : 1 +DEAL:1:2d/2d::Particle id : 0 location : 0.475000 0.475000 property : 10.0000 and 100.000 is ghost on process : 1 DEAL:1:2d/2d::Modifying particles positions and properties -DEAL:1:2d/2d::Particle id : 1 location : 0.625000 0.525000 property : 21.0000 is local on process : 1 -DEAL:1:2d/2d::Particle id : 0 location : 0.575000 0.475000 property : 20.0000 is ghost on process : 1 +DEAL:1:2d/2d::Particle id : 1 location : 0.625000 0.525000 property : 21.0000 and 201.000 is local on process : 1 +DEAL:1:2d/2d::Particle id : 0 location : 0.575000 0.475000 property : 20.0000 and 200.000 is ghost on process : 1 DEAL:1:2d/2d::OK -DEAL:1:2d/3d::Particle id : 1 location : 0.525000 0.525000 0.00000 property : 11.0000 is local on process : 1 -DEAL:1:2d/3d::Particle id : 0 location : 0.475000 0.475000 0.00000 property : 10.0000 is ghost on process : 1 +DEAL:1:2d/3d::Particle id : 1 location : 0.525000 0.525000 0.00000 property : 11.0000 and 101.000 is local on process : 1 +DEAL:1:2d/3d::Particle id : 0 location : 0.475000 0.475000 0.00000 property : 10.0000 and 100.000 is ghost on process : 1 DEAL:1:2d/3d::Modifying particles positions and properties -DEAL:1:2d/3d::Particle id : 1 location : 0.625000 0.525000 0.00000 property : 21.0000 is local on process : 1 -DEAL:1:2d/3d::Particle id : 0 location : 0.575000 0.475000 0.00000 property : 20.0000 is ghost on process : 1 +DEAL:1:2d/3d::Particle id : 1 location : 0.625000 0.525000 0.00000 property : 21.0000 and 201.000 is local on process : 1 +DEAL:1:2d/3d::Particle id : 0 location : 0.575000 0.475000 0.00000 property : 20.0000 and 200.000 is ghost on process : 1 DEAL:1:2d/3d::OK -DEAL:1:3d/3d::Particle id : 1 location : 0.525000 0.525000 0.525000 property : 11.0000 is local on process : 1 -DEAL:1:3d/3d::Particle id : 0 location : 0.475000 0.475000 0.475000 property : 10.0000 is ghost on process : 1 +DEAL:1:3d/3d::Particle id : 1 location : 0.525000 0.525000 0.525000 property : 11.0000 and 101.000 is local on process : 1 +DEAL:1:3d/3d::Particle id : 0 location : 0.475000 0.475000 0.475000 property : 10.0000 and 100.000 is ghost on process : 1 DEAL:1:3d/3d::Modifying particles positions and properties -DEAL:1:3d/3d::Particle id : 1 location : 0.625000 0.525000 0.525000 property : 21.0000 is local on process : 1 -DEAL:1:3d/3d::Particle id : 0 location : 0.575000 0.475000 0.475000 property : 20.0000 is ghost on process : 1 +DEAL:1:3d/3d::Particle id : 1 location : 0.625000 0.525000 0.525000 property : 21.0000 and 201.000 is local on process : 1 +DEAL:1:3d/3d::Particle id : 0 location : 0.575000 0.475000 0.475000 property : 20.0000 and 200.000 is ghost on process : 1 DEAL:1:3d/3d::OK diff --git a/tests/particles/particle_handler_20.cc b/tests/particles/particle_handler_20.cc new file mode 100644 index 0000000000..17f97e3db6 --- /dev/null +++ b/tests/particles/particle_handler_20.cc @@ -0,0 +1,172 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2017 - 2018 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + + +// like particle_handler_06, but after exchanging the ghost particles +// modifies the location of the particles and updates them through +// the update_ghosts mechanism. This test also adds a single property +// to the particles. This property is also modified and then updated +// through the update_ghosts mechanism. + +#include + +#include + +#include +#include + +#include + +#include "../tests.h" + +template +void +test() +{ + { + parallel::distributed::Triangulation tr(MPI_COMM_WORLD); + + GridGenerator::hyper_cube(tr); + tr.refine_global(2); + MappingQ mapping(1); + + Particles::ParticleHandler particle_handler(tr, mapping, 2); + + unsigned int n_particles = 3; + Point position; + Point reference_position; + + for (unsigned int p = 0; p < n_particles; ++p) + { + if (Utilities::MPI::this_mpi_process(tr.get_communicator()) == 0) + { + for (unsigned int i = 0; i < dim; ++i) + position(i) = 0.410 + 0.01 * p; + + Particles::Particle particle( + position, + reference_position, + Utilities::MPI::this_mpi_process(tr.get_communicator()) * + n_particles + + p); + typename Triangulation::active_cell_iterator cell = + tr.begin_active(); + particle_handler.insert_particle(particle, cell); + } + } + + particle_handler.sort_particles_into_subdomains_and_cells(); + + + unsigned int counter = 0; + // Set the properties of the particle to be a unique number + for (auto particle = particle_handler.begin(); + particle != particle_handler.end(); + ++particle) + { + particle->get_properties()[0] = + 1000 + 100 * Utilities::MPI::this_mpi_process(tr.get_communicator()) + + 10 * particle->get_id(); + particle->get_properties()[1] = + 2000 + 100 * Utilities::MPI::this_mpi_process(tr.get_communicator()) + + 10 * particle->get_id(); + counter++; + } + + + particle_handler.exchange_ghost_particles(true); + + for (auto particle = particle_handler.begin(); + particle != particle_handler.end(); + ++particle) + deallog << "Particle id : " << particle->get_id() + << " location : " << particle->get_location() + << " property : " << particle->get_properties()[0] << " and " + << particle->get_properties()[1] << " is local on process : " + << Utilities::MPI::this_mpi_process(tr.get_communicator()) + << std::endl; + + for (auto particle = particle_handler.begin_ghost(); + particle != particle_handler.end_ghost(); + ++particle) + deallog << "Particle id : " << particle->get_id() + << " location : " << particle->get_location() + << " property : " << particle->get_properties()[0] << " and " + << particle->get_properties()[1] << " is ghost on process : " + << Utilities::MPI::this_mpi_process(tr.get_communicator()) + << std::endl; + + deallog << "Modifying particles positions and properties" << std::endl; + + // Modify the location of a single particle on processor 0 and update the + // ghosts + for (auto particle = particle_handler.begin(); + particle != particle_handler.end(); + ++particle) + { + auto location = particle->get_location(); + location[0] += 0.1; + particle->get_properties()[0] += 10000; + particle->get_properties()[1] += 10000; + particle->set_location(location); + } + + // Update the ghost particles + particle_handler.update_ghost_particles(); + + + for (auto particle = particle_handler.begin(); + particle != particle_handler.end(); + ++particle) + deallog << "Particle id : " << particle->get_id() + << " location : " << particle->get_location() + << " property : " << particle->get_properties()[0] << " and " + << particle->get_properties()[1] << " is local on process : " + << Utilities::MPI::this_mpi_process(tr.get_communicator()) + << std::endl; + + for (auto particle = particle_handler.begin_ghost(); + particle != particle_handler.end_ghost(); + ++particle) + deallog << "Particle id : " << particle->get_id() + << " location : " << particle->get_location() + << " property : " << particle->get_properties()[0] << " and " + << particle->get_properties()[1] << " is ghost on process : " + << Utilities::MPI::this_mpi_process(tr.get_communicator()) + << std::endl; + } + + deallog << "OK" << std::endl; +} + + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + + MPILogInitAll all; + + deallog.push("2d/2d"); + test<2, 2>(); + deallog.pop(); + deallog.push("2d/3d"); + test<2, 3>(); + deallog.pop(); + deallog.push("3d/3d"); + test<3, 3>(); + deallog.pop(); +} diff --git a/tests/particles/particle_handler_20.with_p4est=true.mpirun=2.output b/tests/particles/particle_handler_20.with_p4est=true.mpirun=2.output new file mode 100644 index 0000000000..e54f44aa47 --- /dev/null +++ b/tests/particles/particle_handler_20.with_p4est=true.mpirun=2.output @@ -0,0 +1,51 @@ + +DEAL:0:2d/2d::Particle id : 0 location : 0.410000 0.410000 property : 1000.00 and 2000.00 is local on process : 0 +DEAL:0:2d/2d::Particle id : 1 location : 0.420000 0.420000 property : 1010.00 and 2010.00 is local on process : 0 +DEAL:0:2d/2d::Particle id : 2 location : 0.430000 0.430000 property : 1020.00 and 2020.00 is local on process : 0 +DEAL:0:2d/2d::Modifying particles positions and properties +DEAL:0:2d/2d::Particle id : 0 location : 0.510000 0.410000 property : 11000.0 and 12000.0 is local on process : 0 +DEAL:0:2d/2d::Particle id : 1 location : 0.520000 0.420000 property : 11010.0 and 12010.0 is local on process : 0 +DEAL:0:2d/2d::Particle id : 2 location : 0.530000 0.430000 property : 11020.0 and 12020.0 is local on process : 0 +DEAL:0:2d/2d::OK +DEAL:0:2d/3d::Particle id : 0 location : 0.410000 0.410000 0.00000 property : 1000.00 and 2000.00 is local on process : 0 +DEAL:0:2d/3d::Particle id : 1 location : 0.420000 0.420000 0.00000 property : 1010.00 and 2010.00 is local on process : 0 +DEAL:0:2d/3d::Particle id : 2 location : 0.430000 0.430000 0.00000 property : 1020.00 and 2020.00 is local on process : 0 +DEAL:0:2d/3d::Modifying particles positions and properties +DEAL:0:2d/3d::Particle id : 0 location : 0.510000 0.410000 0.00000 property : 11000.0 and 12000.0 is local on process : 0 +DEAL:0:2d/3d::Particle id : 1 location : 0.520000 0.420000 0.00000 property : 11010.0 and 12010.0 is local on process : 0 +DEAL:0:2d/3d::Particle id : 2 location : 0.530000 0.430000 0.00000 property : 11020.0 and 12020.0 is local on process : 0 +DEAL:0:2d/3d::OK +DEAL:0:3d/3d::Particle id : 0 location : 0.410000 0.410000 0.410000 property : 1000.00 and 2000.00 is local on process : 0 +DEAL:0:3d/3d::Particle id : 1 location : 0.420000 0.420000 0.420000 property : 1010.00 and 2010.00 is local on process : 0 +DEAL:0:3d/3d::Particle id : 2 location : 0.430000 0.430000 0.430000 property : 1020.00 and 2020.00 is local on process : 0 +DEAL:0:3d/3d::Modifying particles positions and properties +DEAL:0:3d/3d::Particle id : 0 location : 0.510000 0.410000 0.410000 property : 11000.0 and 12000.0 is local on process : 0 +DEAL:0:3d/3d::Particle id : 1 location : 0.520000 0.420000 0.420000 property : 11010.0 and 12010.0 is local on process : 0 +DEAL:0:3d/3d::Particle id : 2 location : 0.530000 0.430000 0.430000 property : 11020.0 and 12020.0 is local on process : 0 +DEAL:0:3d/3d::OK + +DEAL:1:2d/2d::Particle id : 0 location : 0.410000 0.410000 property : 1000.00 and 2000.00 is ghost on process : 1 +DEAL:1:2d/2d::Particle id : 1 location : 0.420000 0.420000 property : 1010.00 and 2010.00 is ghost on process : 1 +DEAL:1:2d/2d::Particle id : 2 location : 0.430000 0.430000 property : 1020.00 and 2020.00 is ghost on process : 1 +DEAL:1:2d/2d::Modifying particles positions and properties +DEAL:1:2d/2d::Particle id : 0 location : 0.510000 0.410000 property : 11000.0 and 12000.0 is ghost on process : 1 +DEAL:1:2d/2d::Particle id : 1 location : 0.520000 0.420000 property : 11010.0 and 12010.0 is ghost on process : 1 +DEAL:1:2d/2d::Particle id : 2 location : 0.530000 0.430000 property : 11020.0 and 12020.0 is ghost on process : 1 +DEAL:1:2d/2d::OK +DEAL:1:2d/3d::Particle id : 0 location : 0.410000 0.410000 0.00000 property : 1000.00 and 2000.00 is ghost on process : 1 +DEAL:1:2d/3d::Particle id : 1 location : 0.420000 0.420000 0.00000 property : 1010.00 and 2010.00 is ghost on process : 1 +DEAL:1:2d/3d::Particle id : 2 location : 0.430000 0.430000 0.00000 property : 1020.00 and 2020.00 is ghost on process : 1 +DEAL:1:2d/3d::Modifying particles positions and properties +DEAL:1:2d/3d::Particle id : 0 location : 0.510000 0.410000 0.00000 property : 11000.0 and 12000.0 is ghost on process : 1 +DEAL:1:2d/3d::Particle id : 1 location : 0.520000 0.420000 0.00000 property : 11010.0 and 12010.0 is ghost on process : 1 +DEAL:1:2d/3d::Particle id : 2 location : 0.530000 0.430000 0.00000 property : 11020.0 and 12020.0 is ghost on process : 1 +DEAL:1:2d/3d::OK +DEAL:1:3d/3d::Particle id : 0 location : 0.410000 0.410000 0.410000 property : 1000.00 and 2000.00 is ghost on process : 1 +DEAL:1:3d/3d::Particle id : 1 location : 0.420000 0.420000 0.420000 property : 1010.00 and 2010.00 is ghost on process : 1 +DEAL:1:3d/3d::Particle id : 2 location : 0.430000 0.430000 0.430000 property : 1020.00 and 2020.00 is ghost on process : 1 +DEAL:1:3d/3d::Modifying particles positions and properties +DEAL:1:3d/3d::Particle id : 0 location : 0.510000 0.410000 0.410000 property : 11000.0 and 12000.0 is ghost on process : 1 +DEAL:1:3d/3d::Particle id : 1 location : 0.520000 0.420000 0.420000 property : 11010.0 and 12010.0 is ghost on process : 1 +DEAL:1:3d/3d::Particle id : 2 location : 0.530000 0.430000 0.430000 property : 11020.0 and 12020.0 is ghost on process : 1 +DEAL:1:3d/3d::OK +