DOLFIN-X
DOLFIN-X C++ interface
Partitioning.h
1 // Copyright (C) 2020 Garth N. Wells
2 //
3 // This file is part of DOLFINX (https://www.fenicsproject.org)
4 //
5 // SPDX-License-Identifier: LGPL-3.0-or-later
6 
7 #pragma once
8 
9 #include <cstdint>
10 #include <dolfinx/common/IndexMap.h>
11 #include <dolfinx/common/Timer.h>
12 #include <dolfinx/graph/AdjacencyList.h>
13 #include <mpi.h>
14 #include <set>
15 #include <utility>
16 #include <vector>
17 
18 namespace dolfinx::graph
19 {
20 
24 
26 {
27 public:
49  static std::tuple<std::vector<std::int32_t>, std::vector<std::int64_t>,
50  std::vector<int>>
51  reorder_global_indices(MPI_Comm comm,
52  const std::vector<std::int64_t>& global_indices,
53  const std::vector<bool>& shared_indices);
54 
63  static std::pair<graph::AdjacencyList<std::int32_t>,
64  std::vector<std::int64_t>>
66 
77  static std::tuple<graph::AdjacencyList<std::int32_t>, common::IndexMap>
79  MPI_Comm comm, const graph::AdjacencyList<std::int32_t>& list_local,
80  const std::vector<std::int64_t>& local_to_global_links,
81  const std::vector<bool>& shared_links);
82 
94  static std::tuple<graph::AdjacencyList<std::int64_t>, std::vector<int>,
95  std::vector<std::int64_t>, std::vector<int>>
96  distribute(MPI_Comm comm, const graph::AdjacencyList<std::int64_t>& list,
97  const graph::AdjacencyList<std::int32_t>& destinations);
98 
106  static std::vector<std::int64_t>
107  compute_ghost_indices(MPI_Comm comm,
108  const std::vector<std::int64_t>& global_indices,
109  const std::vector<int>& ghost_owners);
110 
120  template <typename T>
121  static Eigen::Array<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>
122  distribute_data(MPI_Comm comm, const std::vector<std::int64_t>& indices,
123  const Eigen::Ref<const Eigen::Array<
124  T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>>& x);
125 
137  static std::vector<std::int64_t> compute_local_to_global_links(
140 
149  static std::vector<std::int32_t>
150  compute_local_to_local(const std::vector<std::int64_t>& local0_to_global,
151  const std::vector<std::int64_t>& local1_to_global);
152 };
153 
154 //---------------------------------------------------------------------------
155 // Implementation
156 //---------------------------------------------------------------------------
157 template <typename T>
158 Eigen::Array<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>
160  MPI_Comm comm, const std::vector<std::int64_t>& indices,
161  const Eigen::Ref<const Eigen::Array<T, Eigen::Dynamic, Eigen::Dynamic,
162  Eigen::RowMajor>>& x)
163 {
164  common::Timer timer("Fetch float data from remote processes");
165 
166  const std::int64_t num_points_local = x.rows();
167  const int size = dolfinx::MPI::size(comm);
168  const int rank = dolfinx::MPI::rank(comm);
169  std::vector<std::int64_t> global_sizes(size);
170  MPI_Allgather(&num_points_local, 1, MPI_INT64_T, global_sizes.data(), 1,
171  MPI_INT64_T, comm);
172  std::vector<std::int64_t> global_offsets(size + 1, 0);
173  std::partial_sum(global_sizes.begin(), global_sizes.end(),
174  global_offsets.begin() + 1);
175 
176  // Build index data requests
177  std::vector<int> number_index_send(size, 0);
178  std::vector<int> index_owner(indices.size());
179  std::vector<int> index_order(indices.size());
180  std::iota(index_order.begin(), index_order.end(), 0);
181  std::sort(index_order.begin(), index_order.end(),
182  [&indices](int a, int b) { return (indices[a] < indices[b]); });
183 
184  int p = 0;
185  for (std::size_t i = 0; i < index_order.size(); ++i)
186  {
187  int j = index_order[i];
188  while (indices[j] >= global_offsets[p + 1])
189  ++p;
190  index_owner[j] = p;
191  number_index_send[p]++;
192  }
193 
194  // Compute send displacements
195  std::vector<int> disp_index_send(size + 1, 0);
196  std::partial_sum(number_index_send.begin(), number_index_send.end(),
197  disp_index_send.begin() + 1);
198 
199  // Pack global index send data
200  std::vector<std::int64_t> indices_send(disp_index_send.back());
201  std::vector<int> disp_tmp = disp_index_send;
202  for (std::size_t i = 0; i < indices.size(); ++i)
203  {
204  const int owner = index_owner[i];
205  indices_send[disp_tmp[owner]++] = indices[i];
206  }
207 
208  // Send/receive number of indices to communicate to each process
209  std::vector<int> number_index_recv(size);
210  MPI_Alltoall(number_index_send.data(), 1, MPI_INT, number_index_recv.data(),
211  1, MPI_INT, comm);
212 
213  // Compute receive displacements
214  std::vector<int> disp_index_recv(size + 1, 0);
215  std::partial_sum(number_index_recv.begin(), number_index_recv.end(),
216  disp_index_recv.begin() + 1);
217 
218  // Send/receive global indices
219  std::vector<std::int64_t> indices_recv(disp_index_recv.back());
220  MPI_Alltoallv(indices_send.data(), number_index_send.data(),
221  disp_index_send.data(), MPI_INT64_T, indices_recv.data(),
222  number_index_recv.data(), disp_index_recv.data(), MPI_INT64_T,
223  comm);
224 
225  const int item_size = x.cols();
226  assert(item_size != 0);
227  // Pack point data to send back (transpose)
228  Eigen::Array<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> x_return(
229  indices_recv.size(), item_size);
230  for (int p = 0; p < size; ++p)
231  {
232  for (int i = disp_index_recv[p]; i < disp_index_recv[p + 1]; ++i)
233  {
234  const std::int32_t index_local = indices_recv[i] - global_offsets[rank];
235  assert(index_local >= 0);
236  x_return.row(i) = x.row(index_local);
237  }
238  }
239 
240  MPI_Datatype compound_type;
241  MPI_Type_contiguous(item_size, dolfinx::MPI::mpi_type<T>(), &compound_type);
242  MPI_Type_commit(&compound_type);
243 
244  // Send back point data
245  Eigen::Array<T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor> my_x(
246  disp_index_send.back(), item_size);
247  MPI_Alltoallv(x_return.data(), number_index_recv.data(),
248  disp_index_recv.data(), compound_type, my_x.data(),
249  number_index_send.data(), disp_index_send.data(), compound_type,
250  comm);
251 
252  return my_x;
253 }
254 
255 } // namespace dolfinx::graph
dolfinx::common::IndexMap
This class represents the distribution index arrays across processes. An index array is a contiguous ...
Definition: IndexMap.h:44
dolfinx::graph::Partitioning::reorder_global_indices
static std::tuple< std::vector< std::int32_t >, std::vector< std::int64_t >, std::vector< int > > reorder_global_indices(MPI_Comm comm, const std::vector< std::int64_t > &global_indices, const std::vector< bool > &shared_indices)
Definition: Partitioning.cpp:22
dolfinx::graph
Graph data structures and algorithms.
Definition: AdjacencyList.h:17
dolfinx::graph::Partitioning::distribute
static std::tuple< graph::AdjacencyList< std::int64_t >, std::vector< int >, std::vector< std::int64_t >, std::vector< int > > distribute(MPI_Comm comm, const graph::AdjacencyList< std::int64_t > &list, const graph::AdjacencyList< std::int32_t > &destinations)
Distribute adjacency list nodes to destination ranks. The global index of each node is assumed to be ...
Definition: Partitioning.cpp:420
dolfinx::graph::AdjacencyList
This class provides a static adjacency list data structure. It is commonly used to store directed gra...
Definition: AdjacencyList.h:27
dolfinx::MPI::rank
static int rank(MPI_Comm comm)
Return process rank for the communicator.
Definition: MPI.cpp:79
dolfinx::MPI::size
static int size(MPI_Comm comm)
Return size of the group (number of processes) associated with the communicator.
Definition: MPI.cpp:87
dolfinx::graph::Partitioning::create_local_adjacency_list
static std::pair< graph::AdjacencyList< std::int32_t >, std::vector< std::int64_t > > create_local_adjacency_list(const graph::AdjacencyList< std::int64_t > &list)
Compute a local AdjacencyList list with contiguous indices from an AdjacencyList that may have non-co...
Definition: Partitioning.cpp:359
dolfinx::graph::Partitioning::create_distributed_adjacency_list
static std::tuple< graph::AdjacencyList< std::int32_t >, common::IndexMap > create_distributed_adjacency_list(MPI_Comm comm, const graph::AdjacencyList< std::int32_t > &list_local, const std::vector< std::int64_t > &local_to_global_links, const std::vector< bool > &shared_links)
Build a distributed AdjacencyList list with re-numbered links from an AdjacencyList that may have non...
Definition: Partitioning.cpp:391
dolfinx::graph::Partitioning::compute_ghost_indices
static std::vector< std::int64_t > compute_ghost_indices(MPI_Comm comm, const std::vector< std::int64_t > &global_indices, const std::vector< int > &ghost_owners)
Compute ghost indices in a global IndexMap space, from a list of arbitrary global indices,...
Definition: Partitioning.cpp:535
dolfinx::graph::Partitioning::compute_local_to_global_links
static std::vector< std::int64_t > compute_local_to_global_links(const graph::AdjacencyList< std::int64_t > &global, const graph::AdjacencyList< std::int32_t > &local)
Given an adjacency list with global, possibly non-contiguous, link indices and a local adjacency list...
Definition: Partitioning.cpp:652
dolfinx::graph::Partitioning
Tools for distributed graphs.
Definition: Partitioning.h:25
dolfinx::graph::Partitioning::distribute_data
static Eigen::Array< T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > distribute_data(MPI_Comm comm, const std::vector< std::int64_t > &indices, const Eigen::Ref< const Eigen::Array< T, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor >> &x)
Distribute data to process ranks where it it required.
Definition: Partitioning.h:159
dolfinx::common::Timer
A timer can be used for timing tasks. The basic usage is.
Definition: Timer.h:30
dolfinx::graph::Partitioning::compute_local_to_local
static std::vector< std::int32_t > compute_local_to_local(const std::vector< std::int64_t > &local0_to_global, const std::vector< std::int64_t > &local1_to_global)
Compute a local0-to-local1 map from two local-to-global maps with common global indices.
Definition: Partitioning.cpp:686