11 #include <dolfinx/common/MPI.h>
12 #include <dolfinx/common/Timer.h>
13 #include <dolfinx/graph/AdjacencyList.h>
33 using partition_fn = std::function<graph::AdjacencyList<std::int32_t>(
34 MPI_Comm comm,
int nparts,
const AdjacencyList<std::int64_t>& local_graph,
35 std::int32_t num_ghost_nodes,
bool ghosting)>;
51 std::int32_t num_ghost_nodes,
bool ghosting);
69 std::tuple<graph::AdjacencyList<std::int64_t>, std::vector<int>,
70 std::vector<std::int64_t>, std::vector<int>>
81 std::vector<std::int64_t>
83 const std::vector<std::int64_t>& global_indices,
84 const std::vector<int>& ghost_owners);
97 const std::vector<std::int64_t>& indices,
98 const xt::xtensor<T, 2>& x);
111 std::vector<std::int64_t>
123 std::vector<std::int32_t>
125 const std::vector<std::int64_t>& local1_to_global);
131 template <
typename T>
134 const xt::xtensor<T, 2>& x)
136 common::Timer timer(
"Fetch float data from remote processes");
138 const std::int64_t num_points_local = x.shape(0);
141 std::vector<std::int64_t> global_sizes(size);
142 MPI_Allgather(&num_points_local, 1, MPI_INT64_T, global_sizes.data(), 1,
144 std::vector<std::int64_t> global_offsets(size + 1, 0);
145 std::partial_sum(global_sizes.begin(), global_sizes.end(),
146 global_offsets.begin() + 1);
149 std::vector<int> number_index_send(size, 0);
150 std::vector<int> index_owner(indices.size());
151 std::vector<int> index_order(indices.size());
152 std::iota(index_order.begin(), index_order.end(), 0);
153 std::sort(index_order.begin(), index_order.end(),
154 [&indices](
int a,
int b) { return (indices[a] < indices[b]); });
157 for (std::size_t i = 0; i < index_order.size(); ++i)
159 int j = index_order[i];
160 while (indices[j] >= global_offsets[p + 1])
163 number_index_send[p]++;
167 std::vector<int> disp_index_send(size + 1, 0);
168 std::partial_sum(number_index_send.begin(), number_index_send.end(),
169 disp_index_send.begin() + 1);
172 std::vector<std::int64_t> indices_send(disp_index_send.back());
173 std::vector<int> disp_tmp = disp_index_send;
174 for (std::size_t i = 0; i < indices.size(); ++i)
176 const int owner = index_owner[i];
177 indices_send[disp_tmp[owner]++] = indices[i];
181 std::vector<int> number_index_recv(size);
182 MPI_Alltoall(number_index_send.data(), 1, MPI_INT, number_index_recv.data(),
186 std::vector<int> disp_index_recv(size + 1, 0);
187 std::partial_sum(number_index_recv.begin(), number_index_recv.end(),
188 disp_index_recv.begin() + 1);
191 std::vector<std::int64_t> indices_recv(disp_index_recv.back());
192 MPI_Alltoallv(indices_send.data(), number_index_send.data(),
193 disp_index_send.data(), MPI_INT64_T, indices_recv.data(),
194 number_index_recv.data(), disp_index_recv.data(), MPI_INT64_T,
197 assert(x.shape(1) != 0);
199 xt::xtensor<T, 2> x_return({indices_recv.size(), x.shape(1)});
200 for (
int p = 0; p < size; ++p)
202 for (
int i = disp_index_recv[p]; i < disp_index_recv[p + 1]; ++i)
204 const std::int32_t index_local = indices_recv[i] - global_offsets[rank];
205 assert(index_local >= 0);
206 for (std::size_t j = 0; j < x.shape(1); ++j)
207 x_return(i, j) = x(index_local, j);
211 MPI_Datatype compound_type;
212 MPI_Type_contiguous(x.shape(1), dolfinx::MPI::mpi_type<T>(), &compound_type);
213 MPI_Type_commit(&compound_type);
216 xt::xtensor<T, 2> my_x(
217 {
static_cast<std::size_t
>(disp_index_send.back()), x.shape(1)});
218 MPI_Alltoallv(x_return.data(), number_index_recv.data(),
219 disp_index_recv.data(), compound_type, my_x.data(),
220 number_index_send.data(), disp_index_send.data(), compound_type,