Reference documentation for deal.II version 9.5.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
vector_data_exchange.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2020 - 2023 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
17#include <deal.II/base/mpi.h>
18#include <deal.II/base/mpi.templates.h>
22#include <deal.II/base/timer.h>
23
25
26#include <boost/serialization/utility.hpp>
27
28#include <map>
29#include <vector>
30
31
33
34namespace internal
35{
36 namespace MatrixFreeFunctions
37 {
38 namespace VectorDataExchange
39 {
41 const std::shared_ptr<const Utilities::MPI::Partitioner> &partitioner)
42 : partitioner(partitioner)
43 {}
44
45
46
47 unsigned int
49 {
50 return partitioner->locally_owned_size();
51 }
52
53
54
55 unsigned int
57 {
58 return partitioner->n_ghost_indices();
59 }
60
61
62
63 unsigned int
65 {
66 return partitioner->n_import_indices();
67 }
68
69
70
71 unsigned int
73 {
74 return 0;
75 }
76
77
78
81 {
82 return partitioner->size();
83 }
84
85
86
87 void
89 const unsigned int communication_channel,
90 const ArrayView<const double> & locally_owned_array,
91 const std::vector<ArrayView<const double>> &shared_arrays,
92 const ArrayView<double> & ghost_array,
93 const ArrayView<double> & temporary_storage,
94 std::vector<MPI_Request> & requests) const
95 {
96 (void)shared_arrays;
97#ifndef DEAL_II_WITH_MPI
98 (void)communication_channel;
99 (void)locally_owned_array;
100 (void)ghost_array;
101 (void)temporary_storage;
102 (void)requests;
103#else
104 partitioner->export_to_ghosted_array_start(communication_channel,
105 locally_owned_array,
106 temporary_storage,
107 ghost_array,
108 requests);
109#endif
110 }
111
112
113
114 void
116 const ArrayView<const double> & locally_owned_array,
117 const std::vector<ArrayView<const double>> &shared_arrays,
118 const ArrayView<double> & ghost_array,
119 std::vector<MPI_Request> & requests) const
120 {
121 (void)locally_owned_array;
122 (void)shared_arrays;
123#ifndef DEAL_II_WITH_MPI
124 (void)ghost_array;
125 (void)requests;
126#else
127 partitioner->export_to_ghosted_array_finish(ghost_array, requests);
128#endif
129 }
130
131
132
133 void
135 const VectorOperation::values vector_operation,
136 const unsigned int communication_channel,
137 const ArrayView<const double> & locally_owned_array,
138 const std::vector<ArrayView<const double>> &shared_arrays,
139 const ArrayView<double> & ghost_array,
140 const ArrayView<double> & temporary_storage,
141 std::vector<MPI_Request> & requests) const
142 {
143 (void)locally_owned_array;
144 (void)shared_arrays;
145#ifndef DEAL_II_WITH_MPI
146 (void)vector_operation;
147 (void)communication_channel;
148 (void)ghost_array;
149 (void)temporary_storage;
150 (void)requests;
151#else
152 partitioner->import_from_ghosted_array_start(vector_operation,
153 communication_channel,
154 ghost_array,
155 temporary_storage,
156 requests);
157#endif
158 }
159
160
161
162 void
164 const VectorOperation::values vector_operation,
165 const ArrayView<double> & locally_owned_storage,
166 const std::vector<ArrayView<const double>> &shared_arrays,
167 const ArrayView<double> & ghost_array,
168 const ArrayView<const double> & temporary_storage,
169 std::vector<MPI_Request> & requests) const
170 {
171 (void)shared_arrays;
172#ifndef DEAL_II_WITH_MPI
173 (void)vector_operation;
174 (void)locally_owned_storage;
175 (void)ghost_array;
176 (void)temporary_storage;
177 (void)requests;
178#else
179 partitioner->import_from_ghosted_array_finish(vector_operation,
180 temporary_storage,
181 locally_owned_storage,
182 ghost_array,
183 requests);
184#endif
185 }
186
187
188
189 void
191 const ArrayView<double> &ghost_array) const
192 {
193 reset_ghost_values_impl(ghost_array);
194 }
195
196
197
198 void
200 const unsigned int communication_channel,
201 const ArrayView<const float> & locally_owned_array,
202 const std::vector<ArrayView<const float>> &shared_arrays,
203 const ArrayView<float> & ghost_array,
204 const ArrayView<float> & temporary_storage,
205 std::vector<MPI_Request> & requests) const
206 {
207 (void)shared_arrays;
208#ifndef DEAL_II_WITH_MPI
209 (void)communication_channel;
210 (void)locally_owned_array;
211 (void)ghost_array;
212 (void)temporary_storage;
213 (void)requests;
214#else
215 partitioner->export_to_ghosted_array_start(communication_channel,
216 locally_owned_array,
217 temporary_storage,
218 ghost_array,
219 requests);
220#endif
221 }
222
223
224
225 void
227 const ArrayView<const float> & locally_owned_array,
228 const std::vector<ArrayView<const float>> &shared_arrays,
229 const ArrayView<float> & ghost_array,
230 std::vector<MPI_Request> & requests) const
231 {
232 (void)locally_owned_array;
233 (void)shared_arrays;
234#ifndef DEAL_II_WITH_MPI
235 (void)ghost_array;
236 (void)requests;
237#else
238 partitioner->export_to_ghosted_array_finish(ghost_array, requests);
239#endif
240 }
241
242
243
244 void
246 const VectorOperation::values vector_operation,
247 const unsigned int communication_channel,
248 const ArrayView<const float> & locally_owned_array,
249 const std::vector<ArrayView<const float>> &shared_arrays,
250 const ArrayView<float> & ghost_array,
251 const ArrayView<float> & temporary_storage,
252 std::vector<MPI_Request> & requests) const
253 {
254 (void)locally_owned_array;
255 (void)shared_arrays;
256#ifndef DEAL_II_WITH_MPI
257 (void)vector_operation;
258 (void)communication_channel;
259 (void)ghost_array;
260 (void)temporary_storage;
261 (void)requests;
262#else
263 partitioner->import_from_ghosted_array_start(vector_operation,
264 communication_channel,
265 ghost_array,
266 temporary_storage,
267 requests);
268#endif
269 }
270
271
272
273 void
275 const VectorOperation::values vector_operation,
276 const ArrayView<float> & locally_owned_storage,
277 const std::vector<ArrayView<const float>> &shared_arrays,
278 const ArrayView<float> & ghost_array,
279 const ArrayView<const float> & temporary_storage,
280 std::vector<MPI_Request> & requests) const
281 {
282 (void)shared_arrays;
283#ifndef DEAL_II_WITH_MPI
284 (void)vector_operation;
285 (void)locally_owned_storage;
286 (void)ghost_array;
287 (void)temporary_storage;
288 (void)requests;
289#else
290 partitioner->import_from_ghosted_array_finish(vector_operation,
291 temporary_storage,
292 locally_owned_storage,
293 ghost_array,
294 requests);
295#endif
296 }
297
298
299
300 void
302 const ArrayView<float> &ghost_array) const
303 {
304 reset_ghost_values_impl(ghost_array);
305 }
306
307
308
309 template <typename Number>
310 void
312 const ArrayView<Number> &ghost_array) const
313 {
314 for (const auto &my_ghosts :
315 partitioner->ghost_indices_within_larger_ghost_set())
316 for (unsigned int j = my_ghosts.first; j < my_ghosts.second; ++j)
317 ghost_array[j] = 0.;
318 }
319
320
321
322 namespace internal
323 {
324 std::pair<std::vector<unsigned int>,
325 std::vector<std::pair<unsigned int, unsigned int>>>
327 const std::vector<unsigned int> &sm_export_ptr,
328 const std::vector<unsigned int> &sm_export_indices)
329 {
330 std::vector<unsigned int> recv_ptr = {0};
331 std::vector<unsigned int> recv_indices;
332 std::vector<unsigned int> recv_len;
333
334 for (unsigned int i = 0; i + 1 < sm_export_ptr.size(); ++i)
335 {
336 if (sm_export_ptr[i] != sm_export_ptr[i + 1])
337 {
338 recv_indices.push_back(sm_export_indices[sm_export_ptr[i]]);
339 recv_len.push_back(1);
340
341 for (unsigned int j = sm_export_ptr[i] + 1;
342 j < sm_export_ptr[i + 1];
343 j++)
344 if (recv_indices.back() + recv_len.back() !=
345 sm_export_indices[j])
346 {
347 recv_indices.push_back(sm_export_indices[j]);
348 recv_len.push_back(1);
349 }
350 else
351 recv_len.back()++;
352 }
353 recv_ptr.push_back(recv_indices.size());
354 }
355
356 std::pair<std::vector<unsigned int>,
357 std::vector<std::pair<unsigned int, unsigned int>>>
358 result;
359
360 result.first = recv_ptr;
361
362 for (unsigned int i = 0; i < recv_indices.size(); ++i)
363 result.second.emplace_back(recv_indices[i], recv_len[i]);
364
365 return result;
366 }
367
368 } // namespace internal
369
370
371
373 const std::shared_ptr<const Utilities::MPI::Partitioner> &partitioner,
374 const MPI_Comm communicator_sm)
375 : comm(partitioner->get_mpi_communicator())
376 , comm_sm(communicator_sm)
377 , n_local_elements(partitioner->locally_owned_range().n_elements())
378 , n_ghost_elements(partitioner->ghost_indices().n_elements())
379 , n_global_elements(partitioner->locally_owned_range().size())
380 {
381#ifndef DEAL_II_WITH_MPI
382 Assert(false, ExcNeedsMPI());
383#else
385 return; // nothing to do in serial case
386
387 const auto &is_locally_owned = partitioner->locally_owned_range();
388 const auto &is_locally_ghost = partitioner->ghost_indices();
389 const auto &ghost_indices_within_larger_ghost_set =
390 partitioner->ghost_indices_within_larger_ghost_set();
391
392 // temporal data strucutures
393 std::vector<unsigned int> n_ghost_indices_in_larger_set_by_remote_rank;
394
395 std::vector<std::array<unsigned int, 3>> ghost_targets_data;
396
397 std::vector<std::array<unsigned int, 3>> import_targets_data;
398
399 std::vector<unsigned int> sm_ghost_ranks;
400
401 std::vector<unsigned int> sm_import_ranks;
402
403 // temporary uncompressed data structures for ghost_indices_subset_data
404 std::vector<unsigned int> ghost_indices_subset_data_ptr = {0};
405 std::vector<unsigned int> ghost_indices_subset_data_indices;
406
407 // ... for import_indices_data
408 std::vector<unsigned int> import_indices_data_ptr = {0};
409 std::vector<unsigned int> import_indices_data_indices;
410
411 // ... for sm_export_data
412 std::vector<unsigned int> sm_export_data_ptr = {0};
413 std::vector<unsigned int> sm_export_data_indices;
414
415 // ... for sm_export_data_this
416 std::vector<unsigned int> sm_export_data_this_ptr = {0};
417 std::vector<unsigned int> sm_export_data_this_indices;
418
419 // ... for sm_import_data
420 std::vector<unsigned int> sm_import_data_ptr = {};
421 std::vector<unsigned int> sm_import_data_indices;
422
423 // ... for sm_import_data_this
424 std::vector<unsigned int> sm_import_data_this_ptr = {0};
425 std::vector<unsigned int> sm_import_data_this_indices;
426
427 // collect ranks of processes of shared-memory domain
428 const auto sm_ranks =
430
431 // determine owners of ghost indices and determine requesters
432 std::vector<unsigned int> owning_ranks_of_ghosts(
433 is_locally_ghost.n_elements());
434
436 process(is_locally_owned,
437 is_locally_ghost,
438 comm,
439 owning_ranks_of_ghosts,
440 /*track_index_requests = */ true);
441
443 std::vector<
444 std::pair<types::global_dof_index, types::global_dof_index>>,
445 std::vector<unsigned int>>
446 consensus_algorithm;
447 consensus_algorithm.run(process, comm);
448
449 // decompress ghost_indices_within_larger_ghost_set for simpler
450 // data access during setup
451 std::vector<unsigned int> shifts_indices;
452 for (const auto &pair : ghost_indices_within_larger_ghost_set)
453 for (unsigned int k = pair.first; k < pair.second; ++k)
454 shifts_indices.push_back(k);
455
456 // process ghost indices
457 {
458 // collect ghost indices according to owning rank
459 std::map<unsigned int, std::vector<types::global_dof_index>>
460 rank_to_local_indices;
461
462 for (unsigned int i = 0; i < owning_ranks_of_ghosts.size(); ++i)
463 rank_to_local_indices[owning_ranks_of_ghosts[i]].push_back(i);
464
465 unsigned int compressed_offset = 0;
466
467 for (const auto &rank_and_local_indices : rank_to_local_indices)
468 {
469 const auto sm_ranks_ptr = std::find(sm_ranks.begin(),
470 sm_ranks.end(),
471 rank_and_local_indices.first);
472
473 if (sm_ranks_ptr == sm_ranks.end()) // remote process
474 {
475 ghost_targets_data.emplace_back(std::array<unsigned int, 3>{{
476 rank_and_local_indices.first, // rank
477 shifts_indices[compressed_offset], // offset
478 static_cast<unsigned int>(
479 rank_and_local_indices.second.size()) // length
480 }});
481
482 for (unsigned int i = 0;
483 i < rank_and_local_indices.second.size();
484 ++i)
485 ghost_indices_subset_data_indices.push_back(
486 shifts_indices[i + compressed_offset]);
487
488 ghost_indices_subset_data_ptr.push_back(
489 ghost_indices_subset_data_indices.size());
490
491 ghost_indices_subset_data.first.push_back(compressed_offset);
492
493 unsigned int i =
495
497 (shifts_indices[ghost_indices_subset_data.first[i] +
498 (ghost_targets_data[i][2] - 1)] -
499 shifts_indices[ghost_indices_subset_data.first[i]]) +
500 1);
501 }
502 else // shared process
503 {
504 sm_ghost_ranks.push_back(
505 std::distance(sm_ranks.begin(), sm_ranks_ptr));
506
507 sm_export_data_ptr.push_back(
508 sm_export_data_ptr.back() +
509 rank_and_local_indices.second.size());
510
511 for (unsigned int i = compressed_offset;
512 i <
513 rank_and_local_indices.second.size() + compressed_offset;
514 ++i)
515 sm_export_data_this_indices.push_back(
516 shifts_indices[i] + is_locally_owned.n_elements());
517
518 sm_export_data_this_ptr.push_back(
519 sm_export_data_this_indices.size());
520 }
521 compressed_offset += rank_and_local_indices.second.size();
522 }
523
524 sm_export_data_indices.resize(sm_export_data_ptr.back());
525 }
526
527 // process requesters
528 {
529 const auto rank_to_global_indices = process.get_requesters();
530
531 for (const auto &rank_and_global_indices : rank_to_global_indices)
532 {
533 const auto sm_ranks_ptr =
534 std::find(sm_ranks.begin(),
535 sm_ranks.end(),
536 rank_and_global_indices.first);
537
538 if (sm_ranks_ptr == sm_ranks.end()) // remote process
539 {
540 import_targets_data.emplace_back(std::array<unsigned int, 3>{{
541 rank_and_global_indices.first, // rank
542 static_cast<unsigned int>(
543 import_indices_data_indices.size()), // offset
544 static_cast<unsigned int>(
545 rank_and_global_indices.second.n_elements()) // length
546 }});
547
548 for (const auto i : rank_and_global_indices.second)
549 import_indices_data_indices.push_back(
550 is_locally_owned.index_within_set(i));
551
552 import_indices_data_ptr.push_back(
553 import_indices_data_indices.size());
554 }
555 else // shared process
556 {
557 sm_import_ranks.push_back(
558 std::distance(sm_ranks.begin(), sm_ranks_ptr));
559
560 for (const auto i : rank_and_global_indices.second)
561 sm_import_data_this_indices.push_back(
562 is_locally_owned.index_within_set(i));
563
564 sm_import_data_this_ptr.push_back(
565 sm_import_data_this_indices.size());
566 }
567 }
568
569 sm_import_data_ptr = sm_import_data_this_ptr;
570 sm_import_data_indices.resize(sm_import_data_this_ptr.back());
571 }
572
573 // send sm_export_data_this to sm-neighbor -> sm_import_data
574 {
575 std::vector<MPI_Request> requests(sm_ghost_ranks.size() +
576 sm_import_ranks.size());
577
578 for (unsigned int i = 0; i < sm_ghost_ranks.size(); ++i)
579 {
580 const int ierr = MPI_Isend(sm_export_data_this_indices.data() +
581 sm_export_data_this_ptr[i],
582 sm_export_data_this_ptr[i + 1] -
583 sm_export_data_this_ptr[i],
584 MPI_UNSIGNED,
586 4,
587 comm_sm,
588 requests.data() + i);
589 AssertThrowMPI(ierr);
590 }
591
592 for (unsigned int i = 0; i < sm_import_ranks.size(); ++i)
593 {
594 const int ierr =
595 MPI_Irecv(sm_import_data_indices.data() + sm_import_data_ptr[i],
596 sm_import_data_ptr[i + 1] - sm_import_data_ptr[i],
597 MPI_UNSIGNED,
599 4,
600 comm_sm,
601 requests.data() + sm_ghost_ranks.size() + i);
602 AssertThrowMPI(ierr);
603 }
604
605 const int ierr =
606 MPI_Waitall(requests.size(), requests.data(), MPI_STATUSES_IGNORE);
607 AssertThrowMPI(ierr);
608 }
609
610 // send sm_import_data_this to sm-neighbor -> sm_export_data_indices
611 {
612 std::vector<MPI_Request> requests(sm_import_ranks.size() +
613 sm_ghost_ranks.size());
614
615 for (unsigned int i = 0; i < sm_import_ranks.size(); ++i)
616 {
617 const int ierr = MPI_Isend(sm_import_data_this_indices.data() +
618 sm_import_data_this_ptr[i],
619 sm_import_data_this_ptr[i + 1] -
620 sm_import_data_this_ptr[i],
621 MPI_UNSIGNED,
623 2,
624 comm_sm,
625 requests.data() + i);
626 AssertThrowMPI(ierr);
627 }
628
629 for (unsigned int i = 0; i < sm_ghost_ranks.size(); ++i)
630 {
631 const int ierr =
632 MPI_Irecv(sm_export_data_indices.data() + sm_export_data_ptr[i],
633 sm_export_data_ptr[i + 1] - sm_export_data_ptr[i],
634 MPI_UNSIGNED,
636 2,
637 comm_sm,
638 requests.data() + sm_import_ranks.size() + i);
639 AssertThrowMPI(ierr);
640 }
641
642 const int ierr =
643 MPI_Waitall(requests.size(), requests.data(), MPI_STATUSES_IGNORE);
644 AssertThrowMPI(ierr);
645 }
646
647 // store data structures and, if needed, compress them
648 this->n_ghost_indices_in_larger_set_by_remote_rank =
650
653 ghost_indices_subset_data_ptr, ghost_indices_subset_data_indices);
654
655 this->ghost_targets_data = ghost_targets_data;
656
657 this->import_targets_data = import_targets_data;
658
659 this->import_indices_data =
660 internal::compress_to_contiguous_ranges(import_indices_data_ptr,
661 import_indices_data_indices);
662
663 this->sm_ghost_ranks = sm_ghost_ranks;
664
665 this->sm_export_data =
667 sm_export_data_indices);
668
669 this->sm_export_data_this =
670 internal::compress_to_contiguous_ranges(sm_export_data_this_ptr,
671 sm_export_data_this_indices);
672
673 this->sm_import_ranks = sm_import_ranks;
674
675 this->sm_import_data =
677 sm_import_data_indices);
678
679 this->sm_import_data_this =
680 internal::compress_to_contiguous_ranges(sm_import_data_this_ptr,
681 sm_import_data_this_indices);
682
683#endif
684 }
685
686
687
688 void
690 const unsigned int communication_channel,
691 const ArrayView<const double> & locally_owned_array,
692 const std::vector<ArrayView<const double>> &shared_arrays,
693 const ArrayView<double> & ghost_array,
694 const ArrayView<double> & temporary_storage,
695 std::vector<MPI_Request> & requests) const
696 {
697 export_to_ghosted_array_start_impl(communication_channel,
698 locally_owned_array,
699 shared_arrays,
700 ghost_array,
701 temporary_storage,
702 requests);
703 }
704
705
706
707 void
709 const ArrayView<const double> & locally_owned_array,
710 const std::vector<ArrayView<const double>> &shared_arrays,
711 const ArrayView<double> & ghost_array,
712 std::vector<MPI_Request> & requests) const
713 {
714 export_to_ghosted_array_finish_impl(locally_owned_array,
715 shared_arrays,
716 ghost_array,
717 requests);
718 }
719
720
721
722 void
724 const VectorOperation::values vector_operation,
725 const unsigned int communication_channel,
726 const ArrayView<const double> & locally_owned_array,
727 const std::vector<ArrayView<const double>> &shared_arrays,
728 const ArrayView<double> & ghost_array,
729 const ArrayView<double> & temporary_storage,
730 std::vector<MPI_Request> & requests) const
731 {
733 communication_channel,
734 locally_owned_array,
735 shared_arrays,
736 ghost_array,
737 temporary_storage,
738 requests);
739 }
740
741
742
743 void
745 const VectorOperation::values vector_operation,
746 const ArrayView<double> & locally_owned_storage,
747 const std::vector<ArrayView<const double>> &shared_arrays,
748 const ArrayView<double> & ghost_array,
749 const ArrayView<const double> & temporary_storage,
750 std::vector<MPI_Request> & requests) const
751 {
753 locally_owned_storage,
754 shared_arrays,
755 ghost_array,
756 temporary_storage,
757 requests);
758 }
759
760
761
762 void
764 const unsigned int communication_channel,
765 const ArrayView<const float> & locally_owned_array,
766 const std::vector<ArrayView<const float>> &shared_arrays,
767 const ArrayView<float> & ghost_array,
768 const ArrayView<float> & temporary_storage,
769 std::vector<MPI_Request> & requests) const
770 {
771 export_to_ghosted_array_start_impl(communication_channel,
772 locally_owned_array,
773 shared_arrays,
774 ghost_array,
775 temporary_storage,
776 requests);
777 }
778
779
780
781 void
783 const ArrayView<const float> & locally_owned_array,
784 const std::vector<ArrayView<const float>> &shared_arrays,
785 const ArrayView<float> & ghost_array,
786 std::vector<MPI_Request> & requests) const
787 {
788 export_to_ghosted_array_finish_impl(locally_owned_array,
789 shared_arrays,
790 ghost_array,
791 requests);
792 }
793
794
795
796 void
798 const VectorOperation::values vector_operation,
799 const unsigned int communication_channel,
800 const ArrayView<const float> & locally_owned_array,
801 const std::vector<ArrayView<const float>> &shared_arrays,
802 const ArrayView<float> & ghost_array,
803 const ArrayView<float> & temporary_storage,
804 std::vector<MPI_Request> & requests) const
805 {
807 communication_channel,
808 locally_owned_array,
809 shared_arrays,
810 ghost_array,
811 temporary_storage,
812 requests);
813 }
814
815
816
817 void
819 const VectorOperation::values vector_operation,
820 const ArrayView<float> & locally_owned_storage,
821 const std::vector<ArrayView<const float>> &shared_arrays,
822 const ArrayView<float> & ghost_array,
823 const ArrayView<const float> & temporary_storage,
824 std::vector<MPI_Request> & requests) const
825 {
827 locally_owned_storage,
828 shared_arrays,
829 ghost_array,
830 temporary_storage,
831 requests);
832 }
833
834
835
836 template <typename Number>
837 void
839 const unsigned int communication_channel,
840 const ArrayView<const Number> & data_this,
841 const std::vector<ArrayView<const Number>> &data_others,
842 const ArrayView<Number> & buffer,
843 const ArrayView<Number> & temporary_storage,
844 std::vector<MPI_Request> & requests) const
845 {
846#ifndef DEAL_II_WITH_MPI
847 Assert(false, ExcNeedsMPI());
848
849 (void)communication_channel;
850 (void)data_this;
851 (void)data_others;
852 (void)buffer;
853 (void)temporary_storage;
854 (void)requests;
855#else
856 (void)data_others;
857
858 requests.resize(sm_import_ranks.size() + sm_ghost_ranks.size() +
860
861 int dummy;
862 // receive a signal that relevant sm neighbors are ready
863 for (unsigned int i = 0; i < sm_ghost_ranks.size(); ++i)
864 {
865 const int ierr =
866 MPI_Irecv(&dummy,
867 0,
868 MPI_INT,
870 communication_channel + 0,
871 comm_sm,
872 requests.data() + sm_import_ranks.size() + i);
873 AssertThrowMPI(ierr);
874 }
875
876 // signal to all relevant sm neighbors that this process is ready
877 for (unsigned int i = 0; i < sm_import_ranks.size(); ++i)
878 {
879 const int ierr = MPI_Isend(&dummy,
880 0,
881 MPI_INT,
883 communication_channel + 0,
884 comm_sm,
885 requests.data() + i);
886 AssertThrowMPI(ierr);
887 }
888
889 // receive data from remote processes
890 for (unsigned int i = 0; i < ghost_targets_data.size(); ++i)
891 {
892 const unsigned int offset =
894 ghost_targets_data[i][2];
895
896 const int ierr = MPI_Irecv(
897 buffer.data() + ghost_targets_data[i][1] + offset,
898 ghost_targets_data[i][2],
899 Utilities::MPI::mpi_type_id_for_type<decltype(*buffer.data())>,
900 ghost_targets_data[i][0],
901 communication_channel + 1,
902 comm,
903 requests.data() + sm_import_ranks.size() + sm_ghost_ranks.size() +
904 i);
905 AssertThrowMPI(ierr);
906 }
907
908 // send data to remote processes
909 for (unsigned int i = 0, k = 0; i < import_targets_data.size(); ++i)
910 {
911 for (unsigned int j = import_indices_data.first[i];
912 j < import_indices_data.first[i + 1];
913 j++)
914 for (unsigned int l = 0; l < import_indices_data.second[j].second;
915 l++, k++)
916 temporary_storage[k] =
917 data_this[import_indices_data.second[j].first + l];
918
919 // send data away
920 const int ierr = MPI_Isend(
921 temporary_storage.data() + import_targets_data[i][1],
923 Utilities::MPI::mpi_type_id_for_type<decltype(*data_this.data())>,
925 communication_channel + 1,
926 comm,
927 requests.data() + sm_import_ranks.size() + sm_ghost_ranks.size() +
928 ghost_targets_data.size() + i);
929 AssertThrowMPI(ierr);
930 }
931#endif
932 }
933
934
935
936 template <typename Number>
937 void
939 const ArrayView<const Number> & data_this,
940 const std::vector<ArrayView<const Number>> &data_others,
941 const ArrayView<Number> & ghost_array,
942 std::vector<MPI_Request> & requests) const
943 {
944 (void)data_this;
945
946#ifndef DEAL_II_WITH_MPI
947 Assert(false, ExcNeedsMPI());
948
949 (void)data_others;
950 (void)ghost_array;
951 (void)requests;
952#else
953
954 AssertDimension(requests.size(),
955 sm_import_ranks.size() + sm_ghost_ranks.size() +
956 ghost_targets_data.size() +
957 import_targets_data.size());
958
959 const auto split =
960 [&](const unsigned int i) -> std::pair<unsigned int, unsigned int> {
962 (sm_ghost_ranks.size() + ghost_targets_data.size()));
963
964 if (i < sm_ghost_ranks.size())
965 return {0, i};
966 else
967 return {1, i - sm_ghost_ranks.size()};
968 };
969
970 for (unsigned int c = 0;
971 c < sm_ghost_ranks.size() + ghost_targets_data.size();
972 c++)
973 {
974 int i;
975 const int ierr =
976 MPI_Waitany(sm_ghost_ranks.size() + ghost_targets_data.size(),
977 requests.data() + sm_import_ranks.size(),
978 &i,
979 MPI_STATUS_IGNORE);
980 AssertThrowMPI(ierr);
981
982 const auto s = split(i);
983 i = s.second;
984
985 if (s.first == 0)
986 {
987 const Number *DEAL_II_RESTRICT data_others_ptr =
988 data_others[sm_ghost_ranks[i]].data();
989 Number *DEAL_II_RESTRICT data_this_ptr = ghost_array.data();
990
991 for (unsigned int lo = sm_export_data.first[i],
992 ko = sm_export_data_this.first[i],
993 li = 0,
994 ki = 0;
995 (lo < sm_export_data.first[i + 1]) &&
996 (ko < sm_export_data_this.first[i + 1]);)
997 {
998 for (; (li < sm_export_data.second[lo].second) &&
999 (ki < sm_export_data_this.second[ko].second);
1000 ++li, ++ki)
1001 data_this_ptr[sm_export_data_this.second[ko].first + ki -
1003 data_others_ptr[sm_export_data.second[lo].first + li];
1004
1005 if (li == sm_export_data.second[lo].second)
1006 {
1007 lo++; // increment outer counter
1008 li = 0; // reset inner counter
1009 }
1010
1011 if (ki == sm_export_data_this.second[ko].second)
1012 {
1013 ko++; // increment outer counter
1014 ki = 0; // reset inner counter
1015 }
1016 }
1017 }
1018 else /*if(s.second == 1)*/
1019 {
1020 const unsigned int offset =
1022 ghost_targets_data[i][2];
1023
1024 for (unsigned int c = 0,
1025 ko = ghost_indices_subset_data.first[i],
1026 ki = 0;
1027 c < ghost_targets_data[i][2];
1028 ++c)
1029 {
1031 ghost_indices_subset_data.second.size());
1032
1033 const unsigned int idx_1 =
1034 ghost_indices_subset_data.second[ko].first + ki;
1035 const unsigned int idx_2 =
1036 ghost_targets_data[i][1] + c + offset;
1037
1038 AssertIndexRange(idx_1, ghost_array.size());
1039 AssertIndexRange(idx_2, ghost_array.size());
1040
1041 if (idx_1 == idx_2)
1042 {
1043 // noting to do
1044 }
1045 else if (idx_1 < idx_2)
1046 {
1047 ghost_array[idx_1] = ghost_array[idx_2];
1048 ghost_array[idx_2] = 0.0;
1049 }
1050 else
1051 {
1052 Assert(false, ExcNotImplemented());
1053 }
1054
1055 ++ki;
1056
1057 if (ki == ghost_indices_subset_data.second[ko].second)
1058 {
1059 ko++; // increment outer counter
1060 ki = 0; // reset inner counter
1061 }
1062 }
1063 }
1064 }
1065
1066 const int ierr =
1067 MPI_Waitall(requests.size(), requests.data(), MPI_STATUSES_IGNORE);
1068 AssertThrowMPI(ierr);
1069
1070#endif
1071 }
1072
1073
1074
1075 template <typename Number>
1076 void
1078 const VectorOperation::values operation,
1079 const unsigned int communication_channel,
1080 const ArrayView<const Number> & data_this,
1081 const std::vector<ArrayView<const Number>> &data_others,
1082 const ArrayView<Number> & buffer,
1083 const ArrayView<Number> & temporary_storage,
1084 std::vector<MPI_Request> & requests) const
1085 {
1086 (void)data_this;
1087
1088#ifndef DEAL_II_WITH_MPI
1089 Assert(false, ExcNeedsMPI());
1090
1091 (void)operation;
1092 (void)communication_channel;
1093 (void)data_others;
1094 (void)buffer;
1095 (void)temporary_storage;
1096 (void)requests;
1097#else
1098 // return;
1099
1100 (void)data_others;
1101 (void)operation;
1102
1104
1105 requests.resize(sm_ghost_ranks.size() + sm_import_ranks.size() +
1106 ghost_targets_data.size() + import_targets_data.size());
1107
1108 int dummy;
1109 for (unsigned int i = 0; i < sm_ghost_ranks.size(); ++i)
1110 {
1111 const int ierr = MPI_Isend(&dummy,
1112 0,
1113 MPI_INT,
1114 sm_ghost_ranks[i],
1115 communication_channel + 1,
1116 comm_sm,
1117 requests.data() + i);
1118 AssertThrowMPI(ierr);
1119 }
1120
1121 for (unsigned int i = 0; i < sm_import_ranks.size(); ++i)
1122 {
1123 const int ierr =
1124 MPI_Irecv(&dummy,
1125 0,
1126 MPI_INT,
1127 sm_import_ranks[i],
1128 communication_channel + 1,
1129 comm_sm,
1130 requests.data() + sm_ghost_ranks.size() + i);
1131 AssertThrowMPI(ierr);
1132 }
1133
1134 for (unsigned int i = 0; i < ghost_targets_data.size(); ++i)
1135 {
1136 for (unsigned int c = 0,
1137 ko = ghost_indices_subset_data.first[i],
1138 ki = 0;
1139 c < ghost_targets_data[i][2];
1140 ++c)
1141 {
1143
1144 const unsigned int idx_1 =
1145 ghost_indices_subset_data.second[ko].first + ki;
1146 const unsigned int idx_2 = ghost_targets_data[i][1] + c;
1147
1148 AssertIndexRange(idx_1, buffer.size());
1149 AssertIndexRange(idx_2, buffer.size());
1150
1151 if (idx_1 == idx_2)
1152 {
1153 // nothing to do
1154 }
1155 else if (idx_2 < idx_1)
1156 {
1157 buffer[idx_2] = buffer[idx_1];
1158 buffer[idx_1] = 0.0;
1159 }
1160 else
1161 {
1162 Assert(false, ExcNotImplemented());
1163 }
1164
1165 if (++ki == ghost_indices_subset_data.second[ko].second)
1166 {
1167 ko++; // increment outer counter
1168 ki = 0; // reset inner counter
1169 }
1170 }
1171
1172 const int ierr = MPI_Isend(
1173 buffer.data() + ghost_targets_data[i][1],
1174 ghost_targets_data[i][2],
1175 Utilities::MPI::mpi_type_id_for_type<decltype(*buffer.data())>,
1176 ghost_targets_data[i][0],
1177 communication_channel + 0,
1178 comm,
1179 requests.data() + sm_ghost_ranks.size() + sm_import_ranks.size() +
1180 i);
1181 AssertThrowMPI(ierr);
1182 }
1183
1184 for (unsigned int i = 0; i < import_targets_data.size(); ++i)
1185 {
1186 const int ierr =
1187 MPI_Irecv(temporary_storage.data() + import_targets_data[i][1],
1188 import_targets_data[i][2],
1190 *temporary_storage.data())>,
1191 import_targets_data[i][0],
1192 communication_channel + 0,
1193 comm,
1194 requests.data() + sm_ghost_ranks.size() +
1195 sm_import_ranks.size() + ghost_targets_data.size() +
1196 i);
1197 AssertThrowMPI(ierr);
1198 }
1199#endif
1200 }
1201
1202
1203
1204 template <typename Number>
1205 void
1207 const VectorOperation::values operation,
1208 const ArrayView<Number> & data_this,
1209 const std::vector<ArrayView<const Number>> &data_others,
1210 const ArrayView<Number> & buffer,
1211 const ArrayView<const Number> & temporary_storage,
1212 std::vector<MPI_Request> & requests) const
1213 {
1214#ifndef DEAL_II_WITH_MPI
1215 Assert(false, ExcNeedsMPI());
1216
1217 (void)operation;
1218 (void)data_this;
1219 (void)data_others;
1220 (void)buffer;
1221 (void)temporary_storage;
1222 (void)requests;
1223#else
1224
1225 (void)operation;
1226
1228
1229 AssertDimension(requests.size(),
1230 sm_ghost_ranks.size() + sm_import_ranks.size() +
1231 ghost_targets_data.size() +
1232 import_targets_data.size());
1233
1234 const auto split =
1235 [&](const unsigned int i) -> std::pair<unsigned int, unsigned int> {
1237 (sm_import_ranks.size() + ghost_targets_data.size() +
1238 import_targets_data.size()));
1239
1240 if (i < sm_import_ranks.size())
1241 return {0, i};
1242 else if (i < (sm_import_ranks.size() + ghost_targets_data.size()))
1243 return {2, i - sm_import_ranks.size()};
1244 else
1245 return {1, i - sm_import_ranks.size() - ghost_targets_data.size()};
1246 };
1247
1248 for (unsigned int c = 0;
1249 c < sm_import_ranks.size() + import_targets_data.size() +
1250 ghost_targets_data.size();
1251 c++)
1252 {
1253 int i;
1254 const int ierr =
1255 MPI_Waitany(sm_import_ranks.size() + import_targets_data.size() +
1256 ghost_targets_data.size(),
1257 requests.data() + sm_ghost_ranks.size(),
1258 &i,
1259 MPI_STATUS_IGNORE);
1260 AssertThrowMPI(ierr);
1261
1262 const auto &s = split(i);
1263 i = s.second;
1264
1265 if (s.first == 0)
1266 {
1267 Number *DEAL_II_RESTRICT data_others_ptr =
1268 const_cast<Number *>(data_others[sm_import_ranks[i]].data());
1269 Number *DEAL_II_RESTRICT data_this_ptr = data_this.data();
1270
1271 for (unsigned int lo = sm_import_data_this.first[i],
1272 ko = sm_import_data.first[i],
1273 li = 0,
1274 ki = 0;
1275 (lo < sm_import_data_this.first[i + 1]) &&
1276 (ko < sm_import_data.first[i + 1]);)
1277 {
1278 for (; (li < sm_import_data_this.second[lo].second) &&
1279 (ki < sm_import_data.second[ko].second);
1280 ++li, ++ki)
1281 {
1282 data_this_ptr[sm_import_data_this.second[lo].first +
1283 li] +=
1284 data_others_ptr[sm_import_data.second[ko].first + ki];
1285 data_others_ptr[sm_import_data.second[ko].first + ki] =
1286 0.0;
1287 }
1288
1289 if (li == sm_import_data_this.second[lo].second)
1290 {
1291 lo++; // increment outer counter
1292 li = 0; // reset inner counter
1293 }
1294 if (ki == sm_import_data.second[ko].second)
1295 {
1296 ko++; // increment outer counter
1297 ki = 0; // reset inner counter
1298 }
1299 }
1300 }
1301 else if (s.first == 1)
1302 {
1303 for (unsigned int j = import_indices_data.first[i],
1304 k = import_targets_data[i][1];
1305 j < import_indices_data.first[i + 1];
1306 j++)
1307 for (unsigned int l = 0;
1308 l < import_indices_data.second[j].second;
1309 l++)
1310 data_this[import_indices_data.second[j].first + l] +=
1311 temporary_storage[k++];
1312 }
1313 else /*if (s.first == 2)*/
1314 {
1315 std::memset(buffer.data() + ghost_targets_data[i][1],
1316 0.0,
1317 (ghost_targets_data[i][2]) * sizeof(Number));
1318 }
1319 }
1320
1321 const int ierr =
1322 MPI_Waitall(requests.size(), requests.data(), MPI_STATUSES_IGNORE);
1323 AssertThrowMPI(ierr);
1324#endif
1325 }
1326
1327
1328
1329 unsigned int
1331 {
1332 return n_local_elements;
1333 }
1334
1335
1336
1337 unsigned int
1339 {
1340 return n_ghost_elements;
1341 }
1342
1343
1344
1345 unsigned int
1347 {
1348 if (import_targets_data.size() == 0)
1349 return 0;
1350 return import_targets_data.back()[1] + import_targets_data.back()[2];
1351 }
1352
1353
1354
1355 unsigned int
1357 {
1358 return sm_import_ranks.size() + sm_ghost_ranks.size(); // TODO
1359 }
1360
1361
1362
1365 {
1366 return n_global_elements;
1367 }
1368
1369
1370
1371 MPI_Comm
1373 {
1374 return this->comm_sm;
1375 }
1376
1377
1378
1379 void
1381 {
1382 reset_ghost_values_impl(ghost_array);
1383 }
1384
1385
1386
1387 void
1389 {
1390 reset_ghost_values_impl(ghost_array);
1391 }
1392
1393
1394
1395 template <typename Number>
1396 void
1398 {
1399 // reset ghost values coming from shared-memory neighbors
1400 // TODO: only needed if values are buffered
1401 for (const auto &i : sm_export_data_this.second)
1402 std::memset(ghost_array.data() + (i.first - n_local_elements),
1403 0,
1404 sizeof(Number) * i.second);
1405
1406 // reset ghost values coming from remote neighbors
1407 for (const auto &i : ghost_indices_subset_data.second)
1408 std::memset(ghost_array.data() + i.first,
1409 0,
1410 sizeof(Number) * i.second);
1411 }
1412
1413
1414
1415 } // namespace VectorDataExchange
1416 } // namespace MatrixFreeFunctions
1417} // namespace internal
1418
1419
value_type * data() const noexcept
Definition array_view.h:553
std::size_t size() const
Definition array_view.h:576
virtual std::vector< unsigned int > run(const std::vector< unsigned int > &targets, const std::function< RequestType(const unsigned int)> &create_request, const std::function< AnswerType(const unsigned int, const RequestType &)> &answer_request, const std::function< void(const unsigned int, const AnswerType &)> &process_answer, const MPI_Comm comm) override
void export_to_ghosted_array_start(const unsigned int communication_channel, const ArrayView< const double > &locally_owned_array, const std::vector< ArrayView< const double > > &shared_arrays, const ArrayView< double > &ghost_array, const ArrayView< double > &temporary_storage, std::vector< MPI_Request > &requests) const override
void reset_ghost_values(const ArrayView< double > &ghost_array) const override
std::vector< std::array< unsigned int, 3 > > ghost_targets_data
std::pair< std::vector< unsigned int >, std::vector< std::pair< unsigned int, unsigned int > > > sm_export_data_this
void export_to_ghosted_array_finish(const ArrayView< const double > &locally_owned_array, const std::vector< ArrayView< const double > > &shared_arrays, const ArrayView< double > &ghost_array, std::vector< MPI_Request > &requests) const override
void import_from_ghosted_array_start(const VectorOperation::values vector_operation, const unsigned int communication_channel, const ArrayView< const double > &locally_owned_array, const std::vector< ArrayView< const double > > &shared_arrays, const ArrayView< double > &ghost_array, const ArrayView< double > &temporary_storage, std::vector< MPI_Request > &requests) const override
void export_to_ghosted_array_start_impl(const unsigned int communication_channel, const ArrayView< const Number > &locally_owned_array, const std::vector< ArrayView< const Number > > &shared_arrays, const ArrayView< Number > &ghost_array, const ArrayView< Number > &temporary_storage, std::vector< MPI_Request > &requests) const
void import_from_ghosted_array_finish_impl(const VectorOperation::values vector_operation, const ArrayView< Number > &locally_owned_storage, const std::vector< ArrayView< const Number > > &shared_arrays, const ArrayView< Number > &ghost_array, const ArrayView< const Number > &temporary_storage, std::vector< MPI_Request > &requests) const
std::pair< std::vector< unsigned int >, std::vector< std::pair< unsigned int, unsigned int > > > sm_export_data
std::vector< std::array< unsigned int, 3 > > import_targets_data
std::pair< std::vector< unsigned int >, std::vector< std::pair< unsigned int, unsigned int > > > ghost_indices_subset_data
void import_from_ghosted_array_finish(const VectorOperation::values vector_operation, const ArrayView< double > &locally_owned_storage, const std::vector< ArrayView< const double > > &shared_arrays, const ArrayView< double > &ghost_array, const ArrayView< const double > &temporary_storage, std::vector< MPI_Request > &requests) const override
void export_to_ghosted_array_finish_impl(const ArrayView< const Number > &locally_owned_array, const std::vector< ArrayView< const Number > > &shared_arrays, const ArrayView< Number > &ghost_array, std::vector< MPI_Request > &requests) const
std::pair< std::vector< unsigned int >, std::vector< std::pair< unsigned int, unsigned int > > > sm_import_data
void import_from_ghosted_array_start_impl(const VectorOperation::values vector_operation, const unsigned int communication_channel, const ArrayView< const Number > &locally_owned_array, const std::vector< ArrayView< const Number > > &shared_arrays, const ArrayView< Number > &ghost_array, const ArrayView< Number > &temporary_storage, std::vector< MPI_Request > &requests) const
Full(const std::shared_ptr< const Utilities::MPI::Partitioner > &partitioner, const MPI_Comm communicator_sm)
void reset_ghost_values_impl(const ArrayView< Number > &ghost_array) const
std::pair< std::vector< unsigned int >, std::vector< std::pair< unsigned int, unsigned int > > > sm_import_data_this
virtual types::global_dof_index size() const override
std::pair< std::vector< unsigned int >, std::vector< std::pair< unsigned int, unsigned int > > > import_indices_data
const std::shared_ptr< const Utilities::MPI::Partitioner > partitioner
void export_to_ghosted_array_finish(const ArrayView< const double > &locally_owned_array, const std::vector< ArrayView< const double > > &shared_arrays, const ArrayView< double > &ghost_array, std::vector< MPI_Request > &requests) const override
PartitionerWrapper(const std::shared_ptr< const Utilities::MPI::Partitioner > &partitioner)
void export_to_ghosted_array_start(const unsigned int communication_channel, const ArrayView< const double > &locally_owned_array, const std::vector< ArrayView< const double > > &shared_arrays, const ArrayView< double > &ghost_array, const ArrayView< double > &temporary_storage, std::vector< MPI_Request > &requests) const override
void reset_ghost_values_impl(const ArrayView< Number > &ghost_array) const
void import_from_ghosted_array_finish(const VectorOperation::values vector_operation, const ArrayView< double > &locally_owned_storage, const std::vector< ArrayView< const double > > &shared_arrays, const ArrayView< double > &ghost_array, const ArrayView< const double > &temporary_storage, std::vector< MPI_Request > &requests) const override
void reset_ghost_values(const ArrayView< double > &ghost_array) const override
void import_from_ghosted_array_start(const VectorOperation::values vector_operation, const unsigned int communication_channel, const ArrayView< const double > &locally_owned_array, const std::vector< ArrayView< const double > > &shared_arrays, const ArrayView< double > &ghost_array, const ArrayView< double > &temporary_storage, std::vector< MPI_Request > &requests) const override
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_RESTRICT
Definition config.h:107
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
static ::ExceptionBase & ExcNotImplemented()
static ::ExceptionBase & ExcNeedsMPI()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
#define AssertThrowMPI(error_code)
#define AssertIndexRange(index, range)
const std::vector< unsigned int > mpi_processes_within_communicator(const MPI_Comm comm_large, const MPI_Comm comm_small)
Definition mpi.cc:173
bool job_supports_mpi()
Definition mpi.cc:1057
const MPI_Datatype mpi_type_id_for_type
Definition mpi.h:1732
std::pair< std::vector< unsigned int >, std::vector< std::pair< unsigned int, unsigned int > > > compress_to_contiguous_ranges(const std::vector< unsigned int > &sm_export_ptr, const std::vector< unsigned int > &sm_export_indices)
const MPI_Comm comm