Reference documentation for deal.II version 9.2.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\}}\)
task_info.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2018 - 2020 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 
16 
19 #include <deal.II/base/mpi.h>
21 #include <deal.II/base/parallel.h>
22 #include <deal.II/base/utilities.h>
23 
25 
26 
27 #ifdef DEAL_II_WITH_THREADS
28 # include <tbb/blocked_range.h>
29 # include <tbb/parallel_for.h>
30 # define TBB_SUPPRESS_DEPRECATED_MESSAGES 1
31 # include <tbb/task.h>
32 # include <tbb/task_scheduler_init.h>
33 # undef TBB_SUPPRESS_DEPRECATED_MESSAGES
34 #endif
35 
36 #include <iostream>
37 #include <set>
38 
40 
41 
42 
43 /*-------------------- Implementation of the matrix-free loop --------------*/
44 namespace internal
45 {
46  namespace MatrixFreeFunctions
47  {
48 #ifdef DEAL_II_WITH_THREADS
49 
50  // This defines the TBB data structures that are needed to schedule the
51  // partition-partition variant
52 
53  namespace partition
54  {
56  {
57  public:
59  const unsigned int partition,
60  const TaskInfo & task_info)
61  : worker(nullptr)
65  {}
66 
68  const unsigned int partition,
69  const TaskInfo & task_info)
70  : worker(&worker)
71  , worker_pointer(nullptr)
74  {}
75 
76  void
77  operator()() const
78  {
79  MFWorkerInterface *used_worker =
80  worker != nullptr ? worker : *worker_pointer;
81  Assert(used_worker != nullptr, ExcInternalError());
82  used_worker->cell(
83  std::make_pair(task_info.cell_partition_data[partition],
85 
86  if (task_info.face_partition_data.empty() == false)
87  {
88  used_worker->face(
89  std::make_pair(task_info.face_partition_data[partition],
91 
92  used_worker->boundary(std::make_pair(
95  }
96  }
97 
98  private:
101  const unsigned int partition;
103  };
104 
105  class CellWork : public tbb::task
106  {
107  public:
109  const unsigned int partition,
110  const TaskInfo & task_info,
111  const bool is_blocked)
112  : dummy(nullptr)
113  , work(worker, partition, task_info)
115  {}
116 
117  tbb::task *
118  execute() override
119  {
120  work();
121 
122  if (is_blocked == true)
123  tbb::empty_task::spawn(*dummy);
124  return nullptr;
125  }
126 
127  tbb::empty_task *dummy;
128 
129  private:
131  const bool is_blocked;
132  };
133 
134 
135 
136  class PartitionWork : public tbb::task
137  {
138  public:
140  const unsigned int partition_in,
141  const TaskInfo & task_info_in,
142  const bool is_blocked_in = false)
143  : dummy(nullptr)
144  , function(function_in)
145  , partition(partition_in)
146  , task_info(task_info_in)
147  , is_blocked(is_blocked_in)
148  {}
149 
150  tbb::task *
151  execute() override
152  {
153  tbb::empty_task *root =
154  new (tbb::task::allocate_root()) tbb::empty_task;
155  const unsigned int evens = task_info.partition_evens[partition];
156  const unsigned int odds = task_info.partition_odds[partition];
157  const unsigned int n_blocked_workers =
159  const unsigned int n_workers =
161  std::vector<CellWork *> worker(n_workers);
162  std::vector<CellWork *> blocked_worker(n_blocked_workers);
163 
164  root->set_ref_count(evens + 1);
165  for (unsigned int j = 0; j < evens; j++)
166  {
167  worker[j] = new (root->allocate_child())
168  CellWork(function,
170  task_info,
171  false);
172  if (j > 0)
173  {
174  worker[j]->set_ref_count(2);
175  blocked_worker[j - 1]->dummy =
176  new (worker[j]->allocate_child()) tbb::empty_task;
177  tbb::task::spawn(*blocked_worker[j - 1]);
178  }
179  else
180  worker[j]->set_ref_count(1);
181  if (j < evens - 1)
182  {
183  blocked_worker[j] = new (worker[j]->allocate_child())
184  CellWork(function,
186  1,
187  task_info,
188  true);
189  }
190  else
191  {
192  if (odds == evens)
193  {
194  worker[evens] = new (worker[j]->allocate_child())
195  CellWork(function,
197  2 * j + 1,
198  task_info,
199  false);
200  tbb::task::spawn(*worker[evens]);
201  }
202  else
203  {
204  tbb::empty_task *child =
205  new (worker[j]->allocate_child()) tbb::empty_task();
206  tbb::task::spawn(*child);
207  }
208  }
209  }
210 
211  root->wait_for_all();
212  root->destroy(*root);
213  if (is_blocked == true)
214  tbb::empty_task::spawn(*dummy);
215  return nullptr;
216  }
217 
218  tbb::empty_task *dummy;
219 
220  private:
221  MFWorkerInterface &function;
222  const unsigned int partition;
224  const bool is_blocked;
225  };
226 
227  } // end of namespace partition
228 
229 
230 
231  namespace color
232  {
233  class CellWork
234  {
235  public:
237  const TaskInfo & task_info_in,
238  const unsigned int partition_in)
239  : worker(worker_in)
240  , task_info(task_info_in)
241  , partition(partition_in)
242  {}
243 
244  void
245  operator()(const tbb::blocked_range<unsigned int> &r) const
246  {
247  const unsigned int start_index =
249  task_info.block_size * r.begin();
250  const unsigned int end_index =
251  std::min(start_index + task_info.block_size * (r.end() - r.begin()),
253  worker.cell(std::make_pair(start_index, end_index));
254 
255  if (task_info.face_partition_data.empty() == false)
256  {
257  AssertThrow(false, ExcNotImplemented());
258  }
259  }
260 
261  private:
264  const unsigned int partition;
265  };
266 
267 
268 
269  class PartitionWork : public tbb::task
270  {
271  public:
273  const unsigned int partition_in,
274  const TaskInfo & task_info_in,
275  const bool is_blocked_in)
276  : dummy(nullptr)
277  , worker(worker_in)
278  , partition(partition_in)
279  , task_info(task_info_in)
280  , is_blocked(is_blocked_in)
281  {}
282 
283  tbb::task *
284  execute() override
285  {
286  const unsigned int n_chunks =
289  1) /
291  parallel_for(tbb::blocked_range<unsigned int>(0, n_chunks, 1),
293  if (is_blocked == true)
294  tbb::empty_task::spawn(*dummy);
295  return nullptr;
296  }
297 
298  tbb::empty_task *dummy;
299 
300  private:
302  const unsigned int partition;
304  const bool is_blocked;
305  };
306 
307  } // end of namespace color
308 
309 
310 
311  class MPICommunication : public tbb::task
312  {
313  public:
315  : worker(worker_in)
317  {}
318 
319  tbb::task *
320  execute() override
321  {
322  if (do_compress == false)
324  else
326  return nullptr;
327  }
328 
329  private:
331  const bool do_compress;
332  };
333 
334 #endif // DEAL_II_WITH_THREADS
335 
336 
337 
338  void
340  {
341  // If we use thread parallelism, we do not currently support to schedule
342  // pieces of updates within the loop, so this index will collect all
343  // calls in that case and work like a single complete loop over all
344  // cells
345  if (scheme != none)
347  else
348  funct.cell_loop_pre_range(
350 
352 
353 #ifdef DEAL_II_WITH_THREADS
354 
355  if (scheme != none)
356  {
359  {
360  tbb::empty_task *root =
361  new (tbb::task::allocate_root()) tbb::empty_task;
362  root->set_ref_count(evens + 1);
363  std::vector<partition::PartitionWork *> worker(n_workers);
364  std::vector<partition::PartitionWork *> blocked_worker(
366  MPICommunication *worker_compr =
367  new (root->allocate_child()) MPICommunication(funct, true);
368  worker_compr->set_ref_count(1);
369  for (unsigned int j = 0; j < evens; j++)
370  {
371  if (j > 0)
372  {
373  worker[j] = new (root->allocate_child())
374  partition::PartitionWork(funct, 2 * j, *this, false);
375  worker[j]->set_ref_count(2);
376  blocked_worker[j - 1]->dummy =
377  new (worker[j]->allocate_child()) tbb::empty_task;
378  tbb::task::spawn(*blocked_worker[j - 1]);
379  }
380  else
381  {
382  worker[j] = new (worker_compr->allocate_child())
383  partition::PartitionWork(funct, 2 * j, *this, false);
384  worker[j]->set_ref_count(2);
385  MPICommunication *worker_dist =
386  new (worker[j]->allocate_child())
387  MPICommunication(funct, false);
388  tbb::task::spawn(*worker_dist);
389  }
390  if (j < evens - 1)
391  {
392  blocked_worker[j] = new (worker[j]->allocate_child())
393  partition::PartitionWork(funct, 2 * j + 1, *this, true);
394  }
395  else
396  {
397  if (odds == evens)
398  {
399  worker[evens] = new (worker[j]->allocate_child())
401  2 * j + 1,
402  *this,
403  false);
404  tbb::task::spawn(*worker[evens]);
405  }
406  else
407  {
408  tbb::empty_task *child =
409  new (worker[j]->allocate_child()) tbb::empty_task();
410  tbb::task::spawn(*child);
411  }
412  }
413  }
414 
415  root->wait_for_all();
416  root->destroy(*root);
417  }
418  else // end of partition-partition, start of partition-color
419  {
420  // check whether there is only one partition. if not, build up the
421  // tree of partitions
422  if (odds > 0)
423  {
424  tbb::empty_task *root =
425  new (tbb::task::allocate_root()) tbb::empty_task;
426  root->set_ref_count(evens + 1);
427  const unsigned int n_blocked_workers =
428  odds - (odds + evens + 1) % 2;
429  const unsigned int n_workers =
431  std::vector<color::PartitionWork *> worker(n_workers);
432  std::vector<color::PartitionWork *> blocked_worker(
434  unsigned int worker_index = 0, slice_index = 0;
435  int spawn_index_child = -2;
436  MPICommunication *worker_compr =
437  new (root->allocate_child()) MPICommunication(funct, true);
438  worker_compr->set_ref_count(1);
439  for (unsigned int part = 0;
440  part < partition_row_index.size() - 1;
441  part++)
442  {
443  if (part == 0)
444  worker[worker_index] =
445  new (worker_compr->allocate_child())
446  color::PartitionWork(funct,
447  slice_index,
448  *this,
449  false);
450  else
451  worker[worker_index] = new (root->allocate_child())
452  color::PartitionWork(funct,
453  slice_index,
454  *this,
455  false);
456  slice_index++;
457  for (; slice_index < partition_row_index[part + 1];
458  slice_index++)
459  {
460  worker[worker_index]->set_ref_count(1);
461  worker_index++;
462  worker[worker_index] =
463  new (worker[worker_index - 1]->allocate_child())
464  color::PartitionWork(funct,
465  slice_index,
466  *this,
467  false);
468  }
469  worker[worker_index]->set_ref_count(2);
470  if (part > 0)
471  {
472  blocked_worker[(part - 1) / 2]->dummy =
473  new (worker[worker_index]->allocate_child())
474  tbb::empty_task;
475  worker_index++;
476  if (spawn_index_child == -1)
477  tbb::task::spawn(*blocked_worker[(part - 1) / 2]);
478  else
479  {
480  Assert(spawn_index_child >= 0,
481  ExcInternalError());
482  tbb::task::spawn(*worker[spawn_index_child]);
483  }
484  spawn_index_child = -2;
485  }
486  else
487  {
488  MPICommunication *worker_dist =
489  new (worker[worker_index]->allocate_child())
490  MPICommunication(funct, false);
491  tbb::task::spawn(*worker_dist);
492  worker_index++;
493  }
494  part += 1;
495  if (part < partition_row_index.size() - 1)
496  {
497  if (part < partition_row_index.size() - 2)
498  {
499  blocked_worker[part / 2] =
500  new (worker[worker_index - 1]->allocate_child())
501  color::PartitionWork(funct,
502  slice_index,
503  *this,
504  true);
505  slice_index++;
506  if (slice_index < partition_row_index[part + 1])
507  {
508  blocked_worker[part / 2]->set_ref_count(1);
509  worker[worker_index] = new (
510  blocked_worker[part / 2]->allocate_child())
511  color::PartitionWork(funct,
512  slice_index,
513  *this,
514  false);
515  slice_index++;
516  }
517  else
518  {
519  spawn_index_child = -1;
520  continue;
521  }
522  }
523  for (; slice_index < partition_row_index[part + 1];
524  slice_index++)
525  {
526  if (slice_index > partition_row_index[part])
527  {
528  worker[worker_index]->set_ref_count(1);
529  worker_index++;
530  }
531  worker[worker_index] =
532  new (worker[worker_index - 1]->allocate_child())
533  color::PartitionWork(funct,
534  slice_index,
535  *this,
536  false);
537  }
538  spawn_index_child = worker_index;
539  worker_index++;
540  }
541  else
542  {
543  tbb::empty_task *final =
544  new (worker[worker_index - 1]->allocate_child())
545  tbb::empty_task;
546  tbb::task::spawn(*final);
547  spawn_index_child = worker_index - 1;
548  }
549  }
550  if (evens == odds)
551  {
552  Assert(spawn_index_child >= 0, ExcInternalError());
553  tbb::task::spawn(*worker[spawn_index_child]);
554  }
555  root->wait_for_all();
556  root->destroy(*root);
557  }
558  // case when we only have one partition: this is the usual
559  // coloring scheme, and we just schedule a parallel for loop for
560  // each color
561  else
562  {
563  Assert(evens <= 1, ExcInternalError());
565 
566  for (unsigned int color = 0; color < partition_row_index[1];
567  ++color)
568  {
569  tbb::empty_task *root =
570  new (tbb::task::allocate_root()) tbb::empty_task;
571  root->set_ref_count(2);
572  color::PartitionWork *worker =
573  new (root->allocate_child())
574  color::PartitionWork(funct, color, *this, false);
575  tbb::empty_task::spawn(*worker);
576  root->wait_for_all();
577  root->destroy(*root);
578  }
579 
580  funct.vector_compress_start();
581  }
582  }
583  }
584  else
585 #endif
586  // serial loop, go through up to three times and do the MPI transfer at
587  // the beginning/end of the second part
588  {
589  for (unsigned int part = 0; part < partition_row_index.size() - 2;
590  ++part)
591  {
592  if (part == 1)
594 
595  for (unsigned int i = partition_row_index[part];
596  i < partition_row_index[part + 1];
597  ++i)
598  {
599  funct.cell_loop_pre_range(i);
600  funct.zero_dst_vector_range(i);
601  AssertIndexRange(i + 1, cell_partition_data.size());
603  {
604  funct.cell(std::make_pair(cell_partition_data[i],
605  cell_partition_data[i + 1]));
606  }
607 
608  if (face_partition_data.empty() == false)
609  {
611  funct.face(std::make_pair(face_partition_data[i],
612  face_partition_data[i + 1]));
613  if (boundary_partition_data[i + 1] >
615  funct.boundary(
616  std::make_pair(boundary_partition_data[i],
617  boundary_partition_data[i + 1]));
618  }
619  funct.cell_loop_post_range(i);
620  }
621 
622  if (part == 1)
623  funct.vector_compress_start();
624  }
625  }
626  funct.vector_compress_finish();
627 
628  if (scheme != none)
630  else
631  funct.cell_loop_post_range(
633  }
634 
635 
636 
638  {
639  clear();
640  }
641 
642 
643 
644  void
646  {
647  n_active_cells = 0;
648  n_ghost_cells = 0;
650  block_size = 0;
651  n_blocks = 0;
652  scheme = none;
653  partition_row_index.clear();
654  partition_row_index.resize(2);
655  cell_partition_data.clear();
656  face_partition_data.clear();
657  boundary_partition_data.clear();
658  evens = 0;
659  odds = 0;
660  n_blocked_workers = 0;
661  n_workers = 0;
662  partition_evens.clear();
663  partition_odds.clear();
665  partition_n_workers.clear();
666  communicator = MPI_COMM_SELF;
667  my_pid = 0;
668  n_procs = 1;
669  }
670 
671 
672 
673  template <typename StreamType>
674  void
676  const std::size_t data_length) const
677  {
678  Utilities::MPI::MinMaxAvg memory_c =
679  Utilities::MPI::min_max_avg(1e-6 * data_length, communicator);
680  if (n_procs < 2)
681  out << memory_c.min;
682  else
683  out << memory_c.min << "/" << memory_c.avg << "/" << memory_c.max;
684  out << " MB" << std::endl;
685  }
686 
687 
688 
689  std::size_t
691  {
692  return (
693  sizeof(*this) +
702  }
703 
704 
705 
706  void
708  std::vector<unsigned int> &boundary_cells)
709  {
710  // try to make the number of boundary cells divisible by the number of
711  // vectors in vectorization
712  unsigned int fillup_needed =
713  (vectorization_length - boundary_cells.size() % vectorization_length) %
715  if (fillup_needed > 0 && boundary_cells.size() < n_active_cells)
716  {
717  // fill additional cells into the list of boundary cells to get a
718  // balanced number. Go through the indices successively until we
719  // found enough indices
720  std::vector<unsigned int> new_boundary_cells;
721  new_boundary_cells.reserve(boundary_cells.size());
722 
723  unsigned int next_free_slot = 0, bound_index = 0;
724  while (fillup_needed > 0 && bound_index < boundary_cells.size())
725  {
726  if (next_free_slot < boundary_cells[bound_index])
727  {
728  // check if there are enough cells to fill with in the
729  // current slot
730  if (next_free_slot + fillup_needed <=
731  boundary_cells[bound_index])
732  {
733  for (unsigned int j =
734  boundary_cells[bound_index] - fillup_needed;
735  j < boundary_cells[bound_index];
736  ++j)
737  new_boundary_cells.push_back(j);
738  fillup_needed = 0;
739  }
740  // ok, not enough indices, so just take them all up to the
741  // next boundary cell
742  else
743  {
744  for (unsigned int j = next_free_slot;
745  j < boundary_cells[bound_index];
746  ++j)
747  new_boundary_cells.push_back(j);
748  fillup_needed -=
749  boundary_cells[bound_index] - next_free_slot;
750  }
751  }
752  new_boundary_cells.push_back(boundary_cells[bound_index]);
753  next_free_slot = boundary_cells[bound_index] + 1;
754  ++bound_index;
755  }
756  while (fillup_needed > 0 &&
757  (new_boundary_cells.size() == 0 ||
758  new_boundary_cells.back() < n_active_cells - 1))
759  new_boundary_cells.push_back(new_boundary_cells.back() + 1);
760  while (bound_index < boundary_cells.size())
761  new_boundary_cells.push_back(boundary_cells[bound_index++]);
762 
763  boundary_cells.swap(new_boundary_cells);
764  }
765 
766  // set the number of cells
767  std::sort(boundary_cells.begin(), boundary_cells.end());
768 
769  // check that number of boundary cells is divisible by
770  // vectorization_length or that it contains all cells
771  Assert(boundary_cells.size() % vectorization_length == 0 ||
772  boundary_cells.size() == n_active_cells,
773  ExcInternalError());
774  }
775 
776 
777 
778  void
780  const std::vector<unsigned int> &cells_with_comm,
781  const unsigned int dofs_per_cell,
782  const bool categories_are_hp,
783  const std::vector<unsigned int> &cell_vectorization_categories,
784  const bool cell_vectorization_categories_strict,
785  const std::vector<unsigned int> &parent_relation,
786  std::vector<unsigned int> & renumbering,
787  std::vector<unsigned char> & incompletely_filled_vectorization)
788  {
789  // This function is decomposed into several steps to determine a good
790  // ordering that satisfies the following constraints:
791  // a. Only cells belonging to the same category (or next higher if the
792  // cell_vectorization_categories_strict is false) can be grouped into
793  // the same SIMD batch
794  // b. hp adaptive computations must form contiguous ranges for the same
795  // degree (category) in cell_partition_data
796  // c. We want to group the cells with the same parent in the same SIMD
797  // lane if possible
798  // d. The cell order should be similar to the initial one
799  // e. Form sets without MPI communication and those with to overlap
800  // communication with computation
801  //
802  // These constraints are satisfied by first grouping by the categories
803  // and, within the groups, to distinguish between cells with a parent
804  // and those without. All of this is set up with batches of cells (with
805  // padding if the size does not match). Then we define a vector of
806  // arrays where we define sorting criteria for the cell batches to
807  // satisfy the items b and d together, split by different parts to
808  // satisfy item e.
809 
810  // Give the compiler a chance to detect that vectorization_length is a
811  // power of two, which allows it to replace integer divisions by shifts
812  unsigned int vectorization_length_bits = 0;
813  unsigned int my_length = vectorization_length;
814  while (my_length >>= 1)
815  ++vectorization_length_bits;
816  const unsigned int n_lanes = 1 << vectorization_length_bits;
817 
818  // Step 1: find tight map of categories for not taking exceeding amounts
819  // of memory below. Sort the new categories by the numbers in the
820  // old one to ensure we respect the given rules
821  unsigned int n_categories = 1;
822  std::vector<unsigned int> tight_category_map(n_active_cells, 0);
823  if (cell_vectorization_categories.empty() == false)
824  {
825  AssertDimension(cell_vectorization_categories.size(),
827 
828  std::set<unsigned int> used_categories;
829  for (unsigned int i = 0; i < n_active_cells; ++i)
830  used_categories.insert(cell_vectorization_categories[i]);
831  std::vector<unsigned int> used_categories_vector(
832  used_categories.size());
833  n_categories = 0;
834  for (const auto &it : used_categories)
835  used_categories_vector[n_categories++] = it;
836  for (unsigned int i = 0; i < n_active_cells; ++i)
837  {
838  const unsigned int index =
839  std::lower_bound(used_categories_vector.begin(),
840  used_categories_vector.end(),
841  cell_vectorization_categories[i]) -
842  used_categories_vector.begin();
843  AssertIndexRange(index, used_categories_vector.size());
844  tight_category_map[i] = index;
845  }
846  }
847 
848  // Step 2: Sort the cells by the category. If we want to fill up the
849  // ranges in vectorization, promote some of the cells to a higher
850  // category
851  std::vector<std::vector<unsigned int>> renumbering_category(n_categories);
852  for (unsigned int i = 0; i < n_active_cells; ++i)
853  renumbering_category[tight_category_map[i]].push_back(i);
854 
855  if (cell_vectorization_categories_strict == false && n_categories > 1)
856  for (unsigned int j = n_categories - 1; j > 0; --j)
857  {
858  unsigned int lower_index = j - 1;
859  while (renumbering_category[j].size() % n_lanes)
860  {
861  while (renumbering_category[j].size() % n_lanes &&
862  !renumbering_category[lower_index].empty())
863  {
864  renumbering_category[j].push_back(
865  renumbering_category[lower_index].back());
866  renumbering_category[lower_index].pop_back();
867  }
868  if (lower_index == 0)
869  break;
870  else
871  --lower_index;
872  }
873  }
874 
875  // Step 3: Use the parent relation to find a good grouping of cells. To
876  // do this, we first put cells of each category defined above into two
877  // bins, those which we know can be grouped together by the given parent
878  // relation and those which cannot
879  std::vector<unsigned int> temporary_numbering;
880  temporary_numbering.reserve(n_active_cells +
881  (n_lanes - 1) * n_categories);
882  const unsigned int n_cells_per_parent =
883  std::count(parent_relation.begin(), parent_relation.end(), 0);
884  std::vector<unsigned int> category_size;
885  for (unsigned int j = 0; j < n_categories; ++j)
886  {
887  std::vector<std::pair<unsigned int, unsigned int>> grouped_cells;
888  std::vector<unsigned int> other_cells;
889  for (const unsigned int cell : renumbering_category[j])
890  if (parent_relation.empty() ||
891  parent_relation[cell] == numbers::invalid_unsigned_int)
892  other_cells.push_back(cell);
893  else
894  grouped_cells.emplace_back(parent_relation[cell], cell);
895 
896  // Compute the number of cells per group
897  std::sort(grouped_cells.begin(), grouped_cells.end());
898  std::vector<unsigned int> n_cells_per_group;
899  unsigned int length = 0;
900  for (unsigned int i = 0; i < grouped_cells.size(); ++i, ++length)
901  if (i > 0 && grouped_cells[i].first != grouped_cells[i - 1].first)
902  {
903  n_cells_per_group.push_back(length);
904  length = 0;
905  }
906  if (length > 0)
907  n_cells_per_group.push_back(length);
908 
909  // Move groups that do not have the complete size (due to
910  // categories) to the 'other_cells'. The cells with correct group
911  // size are immediately appended to the temporary cell numbering
912  auto group_it = grouped_cells.begin();
913  for (unsigned int length : n_cells_per_group)
914  if (length < n_cells_per_parent)
915  for (unsigned int j = 0; j < length; ++j)
916  other_cells.push_back((group_it++)->second);
917  else
918  {
919  // we should not have more cells in a group than in the first
920  // check we did above
921  AssertDimension(length, n_cells_per_parent);
922  for (unsigned int j = 0; j < length; ++j)
923  temporary_numbering.push_back((group_it++)->second);
924  }
925 
926  // Sort the remaining cells and append them as well
927  std::sort(other_cells.begin(), other_cells.end());
928  temporary_numbering.insert(temporary_numbering.end(),
929  other_cells.begin(),
930  other_cells.end());
931 
932  while (temporary_numbering.size() % n_lanes != 0)
933  temporary_numbering.push_back(numbers::invalid_unsigned_int);
934 
935  category_size.push_back(temporary_numbering.size());
936  }
937 
938  // Step 4: Identify the batches with cells marked as "comm"
939  std::vector<bool> batch_with_comm(temporary_numbering.size() / n_lanes,
940  false);
941  std::vector<unsigned int> temporary_numbering_inverse(n_active_cells);
942  for (unsigned int i = 0; i < temporary_numbering.size(); ++i)
943  if (temporary_numbering[i] != numbers::invalid_unsigned_int)
944  temporary_numbering_inverse[temporary_numbering[i]] = i;
945  for (const unsigned int cell : cells_with_comm)
946  batch_with_comm[temporary_numbering_inverse[cell] / n_lanes] = true;
947 
948  // Step 5: Sort the batches of cells by their last cell index to get
949  // good locality, assuming that the initial cell order is of good
950  // locality. In case we have hp calculations with categories, we need to
951  // sort also by the category.
952  std::vector<std::array<unsigned int, 3>> batch_order;
953  std::vector<std::array<unsigned int, 3>> batch_order_comm;
954  for (unsigned int i = 0; i < temporary_numbering.size(); i += n_lanes)
955  {
956  unsigned int max_index = 0;
957  for (unsigned int j = 0; j < n_lanes; ++j)
958  if (temporary_numbering[i + j] < numbers::invalid_unsigned_int)
959  max_index = std::max(temporary_numbering[i + j], max_index);
960  const unsigned int category_hp =
961  categories_are_hp ?
962  std::upper_bound(category_size.begin(), category_size.end(), i) -
963  category_size.begin() :
964  0;
965  const std::array<unsigned int, 3> next{{category_hp, max_index, i}};
966  if (batch_with_comm[i / n_lanes])
967  batch_order_comm.emplace_back(next);
968  else
969  batch_order.emplace_back(next);
970  }
971 
972  std::sort(batch_order.begin(), batch_order.end());
973  std::sort(batch_order_comm.begin(), batch_order_comm.end());
974 
975  // Step 6: Put the cells with communication in the middle of the cell
976  // range. For the MPI case, we need three groups to enable overlap for
977  // communication and computation (part before comm, part with comm, part
978  // after comm), whereas we need one for the other case. And in each
979  // case, we allow for a slot of "ghosted" cells.
980  std::vector<unsigned int> blocks;
981  if (n_procs == 1)
982  {
983  if (batch_order.empty())
984  std::swap(batch_order_comm, batch_order);
985  Assert(batch_order_comm.empty(), ExcInternalError());
986  partition_row_index.resize(3);
987  blocks = {0, static_cast<unsigned int>(batch_order.size())};
988  }
989  else
990  {
991  partition_row_index.resize(5);
992  const unsigned int comm_begin = batch_order.size() / 2;
993  batch_order.insert(batch_order.begin() + comm_begin,
994  batch_order_comm.begin(),
995  batch_order_comm.end());
996  const unsigned int comm_end = comm_begin + batch_order_comm.size();
997  const unsigned int end = batch_order.size();
998  blocks = {0, comm_begin, comm_end, end};
999  }
1000 
1001  // Step 7: Fill in the data by batches for the locally owned cells.
1002  const unsigned int n_cell_batches = batch_order.size();
1003  const unsigned int n_ghost_batches =
1004  (n_ghost_cells + n_lanes - 1) / n_lanes;
1005  incompletely_filled_vectorization.resize(n_cell_batches +
1006  n_ghost_batches);
1007 
1008  cell_partition_data.clear();
1009  cell_partition_data.resize(1, 0);
1010 
1011  renumbering.clear();
1012  renumbering.resize(n_active_cells + n_ghost_cells,
1014 
1015  unsigned int counter = 0;
1016  for (unsigned int block = 0; block < blocks.size() - 1; ++block)
1017  {
1018  const unsigned int grain_size =
1019  std::max((2048U / dofs_per_cell) / 8 * 4, 2U);
1020  for (unsigned int k = blocks[block]; k < blocks[block + 1];
1021  k += grain_size)
1022  cell_partition_data.push_back(
1023  std::min(k + grain_size, blocks[block + 1]));
1024  partition_row_index[block + 1] = cell_partition_data.size() - 1;
1025 
1026  // Set the numbering according to the reordered temporary one
1027  for (unsigned int k = blocks[block]; k < blocks[block + 1]; ++k)
1028  {
1029  const unsigned int pos = batch_order[k][2];
1030  unsigned int j = 0;
1031  for (; j < n_lanes && temporary_numbering[pos + j] !=
1033  ++j)
1034  renumbering[counter++] = temporary_numbering[pos + j];
1035  if (j < n_lanes)
1036  incompletely_filled_vectorization[k] = j;
1037  }
1038  }
1039  AssertDimension(counter, n_active_cells);
1040 
1041  // Step 8: Treat the ghost cells
1042  for (unsigned int cell = n_active_cells;
1043  cell < n_active_cells + n_ghost_cells;
1044  ++cell)
1045  {
1046  if (!cell_vectorization_categories.empty())
1047  AssertDimension(cell_vectorization_categories[cell],
1048  cell_vectorization_categories[n_active_cells]);
1049  renumbering[cell] = cell;
1050  }
1051  if (n_ghost_cells % n_lanes)
1052  incompletely_filled_vectorization.back() = n_ghost_cells % n_lanes;
1053  cell_partition_data.push_back(n_cell_batches + n_ghost_batches);
1054  partition_row_index.back() = cell_partition_data.size() - 1;
1055 
1056 #ifdef DEBUG
1057  std::vector<unsigned int> renumber_cpy(renumbering);
1058  std::sort(renumber_cpy.begin(), renumber_cpy.end());
1059  for (unsigned int i = 0; i < renumber_cpy.size(); ++i)
1060  AssertDimension(i, renumber_cpy[i]);
1061 #endif
1062  }
1063 
1064 
1065 
1066  void
1068  const std::vector<unsigned int> &boundary_cells,
1069  std::vector<unsigned int> & renumbering,
1070  std::vector<unsigned char> & incompletely_filled_vectorization)
1071  {
1072  const unsigned int n_macro_cells =
1074  const unsigned int n_ghost_slots =
1076  incompletely_filled_vectorization.resize(n_macro_cells + n_ghost_slots);
1077  if (n_macro_cells * vectorization_length > n_active_cells)
1078  incompletely_filled_vectorization[n_macro_cells - 1] =
1080  (n_macro_cells * vectorization_length - n_active_cells);
1081  if (n_ghost_slots * vectorization_length > n_ghost_cells)
1082  incompletely_filled_vectorization[n_macro_cells + n_ghost_slots - 1] =
1084  (n_ghost_slots * vectorization_length - n_ghost_cells);
1085 
1086  std::vector<unsigned int> reverse_numbering(
1088  for (unsigned int j = 0; j < boundary_cells.size(); ++j)
1089  reverse_numbering[boundary_cells[j]] = j;
1090  unsigned int counter = boundary_cells.size();
1091  for (unsigned int j = 0; j < n_active_cells; ++j)
1092  if (reverse_numbering[j] == numbers::invalid_unsigned_int)
1093  reverse_numbering[j] = counter++;
1094 
1095  AssertDimension(counter, n_active_cells);
1096  renumbering = Utilities::invert_permutation(reverse_numbering);
1097 
1098  for (unsigned int j = n_active_cells; j < n_active_cells + n_ghost_cells;
1099  ++j)
1100  renumbering.push_back(j);
1101 
1102  // TODO: might be able to simplify this code by not relying on the cell
1103  // partition data while computing the thread graph
1104  cell_partition_data.clear();
1105  cell_partition_data.push_back(0);
1106  if (n_procs > 1)
1107  {
1108  const unsigned int n_macro_boundary_cells =
1109  (boundary_cells.size() + vectorization_length - 1) /
1111  cell_partition_data.push_back(
1112  (n_macro_cells - n_macro_boundary_cells) / 2);
1114  n_macro_boundary_cells);
1115  }
1116  else
1117  AssertDimension(boundary_cells.size(), 0);
1118  cell_partition_data.push_back(n_macro_cells);
1119  cell_partition_data.push_back(cell_partition_data.back() + n_ghost_slots);
1120  partition_row_index.resize(n_procs > 1 ? 4 : 2);
1121  partition_row_index[0] = 0;
1122  partition_row_index[1] = 1;
1123  if (n_procs > 1)
1124  {
1125  partition_row_index[2] = 2;
1126  partition_row_index[3] = 3;
1127  }
1128  }
1129 
1130 
1131 
1132  void
1133  TaskInfo::guess_block_size(const unsigned int dofs_per_cell)
1134  {
1135  // user did not say a positive number, so we have to guess
1136  if (block_size == 0)
1137  {
1138  // we would like to have enough work to do, so as first guess, try
1139  // to get 16 times as many chunks as we have threads on the system.
1142 
1143  // if there are too few degrees of freedom per cell, need to
1144  // increase the block size
1145  const unsigned int minimum_parallel_grain_size = 200;
1146  if (dofs_per_cell * block_size < minimum_parallel_grain_size)
1147  block_size = (minimum_parallel_grain_size / dofs_per_cell + 1);
1148  if (dofs_per_cell * block_size > 10000)
1149  block_size /= 4;
1150 
1151  block_size =
1152  1 << static_cast<unsigned int>(std::log2(block_size + 1));
1153  }
1154  if (block_size > n_active_cells)
1156  }
1157 
1158 
1159 
1160  void
1162  DynamicSparsityPattern & connectivity_large,
1163  std::vector<unsigned int> & renumbering,
1164  std::vector<unsigned char> &irregular_cells,
1165  const bool)
1166  {
1167  const unsigned int n_macro_cells = *(cell_partition_data.end() - 2);
1168  if (n_macro_cells == 0)
1169  return;
1170 
1172 
1173  unsigned int partition = 0, counter = 0;
1174 
1175  // Create connectivity graph for blocks based on connectivity graph for
1176  // cells.
1177  DynamicSparsityPattern connectivity(n_blocks, n_blocks);
1178  make_connectivity_cells_to_blocks(irregular_cells,
1179  connectivity_large,
1180  connectivity);
1181 
1182  // Create cell-block partitioning.
1183 
1184  // For each block of cells, this variable saves to which partitions the
1185  // block belongs. Initialize all to -1 to mark them as not yet assigned
1186  // a partition.
1187  std::vector<unsigned int> cell_partition(n_blocks,
1189 
1190  // In element j of this variable, one puts the old number of the block
1191  // that should be the jth block in the new numeration.
1192  std::vector<unsigned int> partition_list(n_blocks, 0);
1193  std::vector<unsigned int> partition_color_list(n_blocks, 0);
1194 
1195  // This vector points to the start of each partition.
1196  std::vector<unsigned int> partition_size(2, 0);
1197 
1198  // blocking_connectivity = true;
1199 
1200  // The cluster_size in make_partitioning defines that the no. of cells
1201  // in each partition should be a multiple of cluster_size.
1202  unsigned int cluster_size = 1;
1203 
1204  // Make the partitioning of the first layer of the blocks of cells.
1205  make_partitioning(connectivity,
1206  cluster_size,
1207  cell_partition,
1208  partition_list,
1209  partition_size,
1210  partition);
1211 
1212  // Color the cells within each partition
1214  partition,
1215  cell_partition,
1216  partition_list,
1217  partition_size,
1218  partition_color_list);
1219 
1220  partition_list = renumbering;
1221 
1222 #ifdef DEBUG
1223  // in debug mode, check that the partition color list is one-to-one
1224  {
1225  std::vector<unsigned int> sorted_pc_list(partition_color_list);
1226  std::sort(sorted_pc_list.begin(), sorted_pc_list.end());
1227  for (unsigned int i = 0; i < sorted_pc_list.size(); ++i)
1228  Assert(sorted_pc_list[i] == i, ExcInternalError());
1229  }
1230 #endif
1231 
1232  // set the start list for each block and compute the renumbering of
1233  // cells
1234  std::vector<unsigned int> block_start(n_macro_cells + 1);
1235  std::vector<unsigned char> irregular(n_macro_cells);
1236 
1237  unsigned int mcell_start = 0;
1238  block_start[0] = 0;
1239  for (unsigned int block = 0; block < n_blocks; block++)
1240  {
1241  block_start[block + 1] = block_start[block];
1242  for (unsigned int mcell = mcell_start;
1243  mcell < std::min(mcell_start + block_size, n_macro_cells);
1244  ++mcell)
1245  {
1246  unsigned int n_comp = (irregular_cells[mcell] > 0) ?
1247  irregular_cells[mcell] :
1249  block_start[block + 1] += n_comp;
1250  ++counter;
1251  }
1252  mcell_start += block_size;
1253  }
1254  counter = 0;
1255  unsigned int counter_macro = 0;
1256  unsigned int block_size_last =
1257  n_macro_cells - block_size * (n_blocks - 1);
1258  if (block_size_last == 0)
1259  block_size_last = block_size;
1260 
1261  unsigned int tick = 0;
1262  for (unsigned int block = 0; block < n_blocks; block++)
1263  {
1264  unsigned int present_block = partition_color_list[block];
1265  for (unsigned int cell = block_start[present_block];
1266  cell < block_start[present_block + 1];
1267  ++cell)
1268  renumbering[counter++] = partition_list[cell];
1269  unsigned int this_block_size =
1270  (present_block == n_blocks - 1) ? block_size_last : block_size;
1271 
1272  // Also re-compute the content of cell_partition_data to
1273  // contain the numbers of cells, not blocks
1274  if (cell_partition_data[tick] == block)
1275  cell_partition_data[tick++] = counter_macro;
1276 
1277  for (unsigned int j = 0; j < this_block_size; j++)
1278  irregular[counter_macro++] =
1279  irregular_cells[present_block * block_size + j];
1280  }
1281  AssertDimension(tick + 1, cell_partition_data.size());
1282  cell_partition_data.back() = counter_macro;
1283 
1284  irregular_cells.swap(irregular);
1285  AssertDimension(counter, n_active_cells);
1286  AssertDimension(counter_macro, n_macro_cells);
1287 
1288  // check that the renumbering is one-to-one
1289 #ifdef DEBUG
1290  {
1291  std::vector<unsigned int> sorted_renumbering(renumbering);
1292  std::sort(sorted_renumbering.begin(), sorted_renumbering.end());
1293  for (unsigned int i = 0; i < sorted_renumbering.size(); ++i)
1294  Assert(sorted_renumbering[i] == i, ExcInternalError());
1295  }
1296 #endif
1297 
1298 
1300  partition); // Actually sets too much for partition color case
1301 
1302  AssertDimension(cell_partition_data.back(), n_macro_cells);
1303  }
1304 
1305 
1306 
1307  void
1309  const std::vector<unsigned int> &cell_active_fe_index,
1310  DynamicSparsityPattern & connectivity,
1311  std::vector<unsigned int> & renumbering,
1312  std::vector<unsigned char> & irregular_cells,
1313  const bool hp_bool)
1314  {
1315  const unsigned int n_macro_cells = *(cell_partition_data.end() - 2);
1316  if (n_macro_cells == 0)
1317  return;
1318 
1320 
1321  // if we want to block before partitioning, create connectivity graph
1322  // for blocks based on connectivity graph for cells.
1323  DynamicSparsityPattern connectivity_blocks(n_blocks, n_blocks);
1324  make_connectivity_cells_to_blocks(irregular_cells,
1325  connectivity,
1326  connectivity_blocks);
1327 
1328  unsigned int n_blocks = 0;
1329  if (scheme == partition_color ||
1330  scheme == color) // blocking_connectivity == true
1331  n_blocks = this->n_blocks;
1332  else
1333  n_blocks = n_active_cells;
1334 
1335  // For each block of cells, this variable saves to which partitions the
1336  // block belongs. Initialize all to -1 to mark them as not yet assigned
1337  // a partition.
1338  std::vector<unsigned int> cell_partition(n_blocks,
1340 
1341  // In element j of this variable, one puts the old number (but after
1342  // renumbering according to the input renumbering) of the block that
1343  // should be the jth block in the new numeration.
1344  std::vector<unsigned int> partition_list(n_blocks, 0);
1345  std::vector<unsigned int> partition_2layers_list(n_blocks, 0);
1346 
1347  // This vector points to the start of each partition.
1348  std::vector<unsigned int> partition_size(2, 0);
1349 
1350  unsigned int partition = 0;
1351 
1352  // Within the partitions we want to be able to block for the case that
1353  // we do not block already in the connectivity. The cluster_size in
1354  // make_partitioning defines that the no. of cells in each partition
1355  // should be a multiple of cluster_size.
1356  unsigned int cluster_size = 1;
1357  if (scheme == partition_partition)
1358  cluster_size = block_size * vectorization_length;
1359 
1360  // Make the partitioning of the first layer of the blocks of cells.
1361  if (scheme == partition_color || scheme == color)
1362  make_partitioning(connectivity_blocks,
1363  cluster_size,
1364  cell_partition,
1365  partition_list,
1366  partition_size,
1367  partition);
1368  else
1369  make_partitioning(connectivity,
1370  cluster_size,
1371  cell_partition,
1372  partition_list,
1373  partition_size,
1374  partition);
1375 
1376  // Partition or color second layer
1377  if (scheme == partition_partition)
1378 
1379  {
1380  // Partition within partitions.
1382  connectivity,
1383  cell_active_fe_index,
1384  partition,
1385  cluster_size,
1386  hp_bool,
1387  cell_partition,
1388  partition_list,
1389  partition_size,
1390  partition_2layers_list,
1391  irregular_cells);
1392  }
1393  else if (scheme == partition_color || scheme == color)
1394  {
1395  make_coloring_within_partitions_pre_blocked(connectivity_blocks,
1396  partition,
1397  cell_partition,
1398  partition_list,
1399  partition_size,
1400  partition_2layers_list);
1401  }
1402 
1403  // in debug mode, check that the partition_2layers_list is one-to-one
1404 #ifdef DEBUG
1405  {
1406  std::vector<unsigned int> sorted_pc_list(partition_2layers_list);
1407  std::sort(sorted_pc_list.begin(), sorted_pc_list.end());
1408  for (unsigned int i = 0; i < sorted_pc_list.size(); ++i)
1409  Assert(sorted_pc_list[i] == i, ExcInternalError());
1410  }
1411 #endif
1412 
1413  // Set the new renumbering
1414  std::vector<unsigned int> renumbering_in(n_active_cells, 0);
1415  renumbering_in.swap(renumbering);
1416  if (scheme == partition_partition) // blocking_connectivity == false
1417  {
1418  // This is the simple case. The renumbering is just a combination of
1419  // the renumbering that we were given as an input and the
1420  // renumbering of partition/coloring given in partition_2layers_list
1421  for (unsigned int j = 0; j < renumbering.size(); j++)
1422  renumbering[j] = renumbering_in[partition_2layers_list[j]];
1423  // Account for the ghost cells, finally.
1424  for (unsigned int i = 0; i < n_ghost_cells; ++i)
1425  renumbering.push_back(i + n_active_cells);
1426  }
1427  else
1428  {
1429  // set the start list for each block and compute the renumbering of
1430  // cells
1431  std::vector<unsigned int> block_start(n_macro_cells + 1);
1432  std::vector<unsigned char> irregular(n_macro_cells);
1433 
1434  unsigned int counter = 0;
1435  unsigned int mcell_start = 0;
1436  block_start[0] = 0;
1437  for (unsigned int block = 0; block < n_blocks; block++)
1438  {
1439  block_start[block + 1] = block_start[block];
1440  for (unsigned int mcell = mcell_start;
1441  mcell < std::min(mcell_start + block_size, n_macro_cells);
1442  ++mcell)
1443  {
1444  unsigned int n_comp = (irregular_cells[mcell] > 0) ?
1445  irregular_cells[mcell] :
1447  block_start[block + 1] += n_comp;
1448  ++counter;
1449  }
1450  mcell_start += block_size;
1451  }
1452  counter = 0;
1453  unsigned int counter_macro = 0;
1454  unsigned int block_size_last =
1455  n_macro_cells - block_size * (n_blocks - 1);
1456  if (block_size_last == 0)
1457  block_size_last = block_size;
1458 
1459  unsigned int tick = 0;
1460  for (unsigned int block = 0; block < n_blocks; block++)
1461  {
1462  unsigned int present_block = partition_2layers_list[block];
1463  for (unsigned int cell = block_start[present_block];
1464  cell < block_start[present_block + 1];
1465  ++cell)
1466  renumbering[counter++] = renumbering_in[cell];
1467  unsigned int this_block_size =
1468  (present_block == n_blocks - 1) ? block_size_last : block_size;
1469 
1470  // Also re-compute the content of cell_partition_data to
1471  // contain the numbers of cells, not blocks
1472  if (cell_partition_data[tick] == block)
1473  cell_partition_data[tick++] = counter_macro;
1474 
1475  for (unsigned int j = 0; j < this_block_size; j++)
1476  irregular[counter_macro++] =
1477  irregular_cells[present_block * block_size + j];
1478  }
1479  AssertDimension(tick + 1, cell_partition_data.size());
1480  cell_partition_data.back() = counter_macro;
1481 
1482  irregular_cells.swap(irregular);
1483  AssertDimension(counter, n_active_cells);
1484  AssertDimension(counter_macro, n_macro_cells);
1485  // check that the renumbering is one-to-one
1486 #ifdef DEBUG
1487  {
1488  std::vector<unsigned int> sorted_renumbering(renumbering);
1489  std::sort(sorted_renumbering.begin(), sorted_renumbering.end());
1490  for (unsigned int i = 0; i < sorted_renumbering.size(); ++i)
1491  Assert(sorted_renumbering[i] == i, ExcInternalError());
1492  }
1493 #endif
1494  }
1495 
1496  // Update the task_info with the more information for the thread graph.
1498  }
1499 
1500 
1501 
1502  void
1504  const std::vector<unsigned int> &cell_active_fe_index,
1505  DynamicSparsityPattern & connectivity,
1506  std::vector<unsigned int> & renumbering,
1507  std::vector<unsigned char> & irregular_cells,
1508  const bool hp_bool)
1509  {
1510  const unsigned int n_macro_cells = *(cell_partition_data.end() - 2);
1511  if (n_macro_cells == 0)
1512  return;
1513 
1514  const unsigned int cluster_size = block_size * vectorization_length;
1515 
1516  // Create cell-block partitioning.
1517 
1518  // For each block of cells, this variable saves to which partitions the
1519  // block belongs. Initialize all to n_macro_cells to mark them as not
1520  // yet assigned a partition.
1521  std::vector<unsigned int> cell_partition(n_active_cells,
1523 
1524 
1525  // In element j of this variable, one puts the old number of the block
1526  // that should be the jth block in the new numeration.
1527  std::vector<unsigned int> partition_list(n_active_cells, 0);
1528  std::vector<unsigned int> partition_partition_list(n_active_cells, 0);
1529 
1530  // This vector points to the start of each partition.
1531  std::vector<unsigned int> partition_size(2, 0);
1532 
1533  unsigned int partition = 0;
1534  // Here, we do not block inside the connectivity graph
1535  // blocking_connectivity = false;
1536 
1537  // Make the partitioning of the first layer of the blocks of cells.
1538  make_partitioning(connectivity,
1539  cluster_size,
1540  cell_partition,
1541  partition_list,
1542  partition_size,
1543  partition);
1544 
1545  // Partition within partitions.
1547  cell_active_fe_index,
1548  partition,
1549  cluster_size,
1550  hp_bool,
1551  cell_partition,
1552  partition_list,
1553  partition_size,
1554  partition_partition_list,
1555  irregular_cells);
1556 
1557  partition_list.swap(renumbering);
1558 
1559  for (unsigned int j = 0; j < renumbering.size(); j++)
1560  renumbering[j] = partition_list[partition_partition_list[j]];
1561 
1562  for (unsigned int i = 0; i < n_ghost_cells; ++i)
1563  renumbering.push_back(i + n_active_cells);
1564 
1566  }
1567 
1568 
1569 
1570  void
1572  const std::vector<unsigned char> &irregular_cells,
1573  const DynamicSparsityPattern & connectivity_cells,
1574  DynamicSparsityPattern & connectivity_blocks) const
1575  {
1576  std::vector<std::vector<unsigned int>> cell_blocks(n_blocks);
1577  std::vector<unsigned int> touched_cells(n_active_cells);
1578  unsigned int cell = 0;
1579  for (unsigned int i = 0, mcell = 0; i < n_blocks; ++i)
1580  {
1581  for (unsigned int c = 0;
1582  c < block_size && mcell < *(cell_partition_data.end() - 2);
1583  ++c, ++mcell)
1584  {
1585  unsigned int ncomp = (irregular_cells[mcell] > 0) ?
1586  irregular_cells[mcell] :
1588  for (unsigned int c = 0; c < ncomp; ++c, ++cell)
1589  {
1590  cell_blocks[i].push_back(cell);
1591  touched_cells[cell] = i;
1592  }
1593  }
1594  }
1596  for (unsigned int i = 0; i < cell_blocks.size(); ++i)
1597  for (unsigned int col = 0; col < cell_blocks[i].size(); ++col)
1598  {
1600  connectivity_cells.begin(cell_blocks[i][col]);
1601  it != connectivity_cells.end(cell_blocks[i][col]);
1602  ++it)
1603  {
1604  if (touched_cells[it->column()] != i)
1605  connectivity_blocks.add(i, touched_cells[it->column()]);
1606  }
1607  }
1608  }
1609 
1610 
1611 
1612  // Function to create partitioning on the second layer within each
1613  // partition. Version without preblocking.
1614  void
1616  const DynamicSparsityPattern & connectivity,
1617  const std::vector<unsigned int> &cell_active_fe_index,
1618  const unsigned int partition,
1619  const unsigned int cluster_size,
1620  const bool hp_bool,
1621  const std::vector<unsigned int> &cell_partition,
1622  const std::vector<unsigned int> &partition_list,
1623  const std::vector<unsigned int> &partition_size,
1624  std::vector<unsigned int> & partition_partition_list,
1625  std::vector<unsigned char> & irregular_cells)
1626  {
1627  const unsigned int n_macro_cells = *(cell_partition_data.end() - 2);
1628  const unsigned int n_ghost_slots =
1629  *(cell_partition_data.end() - 1) - n_macro_cells;
1630 
1631  // List of cells in previous partition
1632  std::vector<unsigned int> neighbor_list;
1633  // List of cells in current partition for use as neighbors in next
1634  // partition
1635  std::vector<unsigned int> neighbor_neighbor_list;
1636 
1637  std::vector<unsigned int> renumbering(n_active_cells);
1638 
1639  irregular_cells.back() = 0;
1640  irregular_cells.resize(n_active_cells + n_ghost_slots);
1641 
1642  unsigned int max_fe_index = 0;
1643  for (const unsigned int fe_index : cell_active_fe_index)
1644  max_fe_index = std::max(fe_index, max_fe_index);
1645 
1646  Assert(!hp_bool || cell_active_fe_index.size() == n_active_cells,
1647  ExcInternalError());
1648 
1649  {
1650  unsigned int n_macro_cells_before = 0;
1651  // Create partitioning within partitions.
1652 
1653  // For each block of cells, this variable saves to which partitions
1654  // the block belongs. Initialize all to n_macro_cells to mark them as
1655  // not yet assigned a partition.
1656  std::vector<unsigned int> cell_partition_l2(
1658  partition_row_index.clear();
1659  partition_row_index.resize(partition + 1, 0);
1660  cell_partition_data.resize(1, 0);
1661 
1662  unsigned int counter = 0;
1663  unsigned int missing_macros;
1664  for (unsigned int part = 0; part < partition; ++part)
1665  {
1666  neighbor_neighbor_list.resize(0);
1667  neighbor_list.resize(0);
1668  bool work = true;
1669  unsigned int partition_l2 = 0;
1670  unsigned int start_up = partition_size[part];
1671  unsigned int partition_counter = 0;
1672  while (work)
1673  {
1674  if (neighbor_list.size() == 0)
1675  {
1676  work = false;
1677  partition_counter = 0;
1678  for (unsigned int j = start_up;
1679  j < partition_size[part + 1];
1680  ++j)
1681  if (cell_partition[partition_list[j]] == part &&
1682  cell_partition_l2[partition_list[j]] ==
1684  {
1685  start_up = j;
1686  work = true;
1687  partition_counter = 1;
1688  // To start up, set the start_up cell to partition
1689  // and list all its neighbors.
1690  AssertIndexRange(start_up, partition_size[part + 1]);
1691  cell_partition_l2[partition_list[start_up]] =
1692  partition_l2;
1693  neighbor_neighbor_list.push_back(
1694  partition_list[start_up]);
1695  partition_partition_list[counter++] =
1696  partition_list[start_up];
1697  start_up++;
1698  break;
1699  }
1700  }
1701  else
1702  {
1703  partition_counter = 0;
1704  for (const unsigned int neighbor : neighbor_list)
1705  {
1706  Assert(cell_partition[neighbor] == part,
1707  ExcInternalError());
1708  Assert(cell_partition_l2[neighbor] == partition_l2 - 1,
1709  ExcInternalError());
1710  auto neighbor_it = connectivity.begin(neighbor);
1711  const auto end_it = connectivity.end(neighbor);
1712  for (; neighbor_it != end_it; ++neighbor_it)
1713  {
1714  if (cell_partition[neighbor_it->column()] == part &&
1715  cell_partition_l2[neighbor_it->column()] ==
1717  {
1718  cell_partition_l2[neighbor_it->column()] =
1719  partition_l2;
1720  neighbor_neighbor_list.push_back(
1721  neighbor_it->column());
1722  partition_partition_list[counter++] =
1723  neighbor_it->column();
1724  partition_counter++;
1725  }
1726  }
1727  }
1728  }
1729  if (partition_counter > 0)
1730  {
1731  int index_before = neighbor_neighbor_list.size(),
1732  index = index_before;
1733  {
1734  // put the cells into separate lists for each FE index
1735  // within one partition-partition
1736  missing_macros = 0;
1737  std::vector<unsigned int> remaining_per_macro_cell(
1738  max_fe_index + 1);
1739  std::vector<std::vector<unsigned int>>
1740  renumbering_fe_index;
1741  unsigned int cell;
1742  bool filled = true;
1743  if (hp_bool == true)
1744  {
1745  renumbering_fe_index.resize(max_fe_index + 1);
1746  for (cell = counter - partition_counter;
1747  cell < counter;
1748  ++cell)
1749  {
1750  renumbering_fe_index
1751  [cell_active_fe_index.empty() ?
1752  0 :
1753  cell_active_fe_index
1754  [partition_partition_list[cell]]]
1755  .push_back(partition_partition_list[cell]);
1756  }
1757  // check how many more cells are needed in the lists
1758  for (unsigned int j = 0; j < max_fe_index + 1; j++)
1759  {
1760  remaining_per_macro_cell[j] =
1761  renumbering_fe_index[j].size() %
1763  if (remaining_per_macro_cell[j] != 0)
1764  filled = false;
1765  missing_macros +=
1766  ((renumbering_fe_index[j].size() +
1767  vectorization_length - 1) /
1769  }
1770  }
1771  else
1772  {
1773  remaining_per_macro_cell.resize(1);
1774  remaining_per_macro_cell[0] =
1775  partition_counter % vectorization_length;
1776  missing_macros =
1777  partition_counter / vectorization_length;
1778  if (remaining_per_macro_cell[0] != 0)
1779  {
1780  filled = false;
1781  missing_macros++;
1782  }
1783  }
1784  missing_macros =
1785  cluster_size - (missing_macros % cluster_size);
1786 
1787  // now we realized that there are some cells missing.
1788  while (missing_macros > 0 || filled == false)
1789  {
1790  if (index == 0)
1791  {
1792  index = neighbor_neighbor_list.size();
1793  if (index == index_before)
1794  {
1795  if (missing_macros != 0)
1796  {
1797  neighbor_neighbor_list.resize(0);
1798  }
1799  start_up--;
1800  break; // not connected - start again
1801  }
1802  index_before = index;
1803  }
1804  index--;
1805  unsigned int additional =
1806  neighbor_neighbor_list[index];
1807 
1808  // go through the neighbors of the last cell in the
1809  // current partition and check if we find some to
1810  // fill up with.
1812  connectivity.begin(
1813  additional),
1814  end =
1815  connectivity.end(
1816  additional);
1817  for (; neighbor != end; ++neighbor)
1818  {
1819  if (cell_partition[neighbor->column()] == part &&
1820  cell_partition_l2[neighbor->column()] ==
1822  {
1823  unsigned int this_index = 0;
1824  if (hp_bool == true)
1825  this_index =
1826  cell_active_fe_index.empty() ?
1827  0 :
1828  cell_active_fe_index[neighbor
1829  ->column()];
1830 
1831  // Only add this cell if we need more macro
1832  // cells in the current block or if there is
1833  // a macro cell with the FE index that is
1834  // not yet fully populated
1835  if (missing_macros > 0 ||
1836  remaining_per_macro_cell[this_index] > 0)
1837  {
1838  cell_partition_l2[neighbor->column()] =
1839  partition_l2;
1840  neighbor_neighbor_list.push_back(
1841  neighbor->column());
1842  if (hp_bool == true)
1843  renumbering_fe_index[this_index]
1844  .push_back(neighbor->column());
1845  partition_partition_list[counter] =
1846  neighbor->column();
1847  counter++;
1848  partition_counter++;
1849  if (remaining_per_macro_cell
1850  [this_index] == 0 &&
1851  missing_macros > 0)
1852  missing_macros--;
1853  remaining_per_macro_cell[this_index]++;
1854  if (remaining_per_macro_cell
1855  [this_index] ==
1857  {
1858  remaining_per_macro_cell[this_index] =
1859  0;
1860  }
1861  if (missing_macros == 0)
1862  {
1863  filled = true;
1864  for (unsigned int fe_ind = 0;
1865  fe_ind < max_fe_index + 1;
1866  ++fe_ind)
1867  if (remaining_per_macro_cell
1868  [fe_ind] != 0)
1869  filled = false;
1870  }
1871  if (filled == true)
1872  break;
1873  }
1874  }
1875  }
1876  }
1877  if (hp_bool == true)
1878  {
1879  // set the renumbering according to their active FE
1880  // index within one partition-partition which was
1881  // implicitly assumed above
1882  cell = counter - partition_counter;
1883  for (unsigned int j = 0; j < max_fe_index + 1; j++)
1884  {
1885  for (const unsigned int jj :
1886  renumbering_fe_index[j])
1887  renumbering[cell++] = jj;
1888  if (renumbering_fe_index[j].size() %
1890  0)
1891  irregular_cells[renumbering_fe_index[j].size() /
1893  n_macro_cells_before] =
1894  renumbering_fe_index[j].size() %
1896  n_macro_cells_before +=
1897  (renumbering_fe_index[j].size() +
1898  vectorization_length - 1) /
1900  renumbering_fe_index[j].resize(0);
1901  }
1902  }
1903  else
1904  {
1905  n_macro_cells_before +=
1906  partition_counter / vectorization_length;
1907  if (partition_counter % vectorization_length != 0)
1908  {
1909  irregular_cells[n_macro_cells_before] =
1910  partition_counter % vectorization_length;
1911  n_macro_cells_before++;
1912  }
1913  }
1914  }
1915  cell_partition_data.push_back(n_macro_cells_before);
1916  partition_l2++;
1917  }
1918  neighbor_list = neighbor_neighbor_list;
1919  neighbor_neighbor_list.resize(0);
1920  }
1921  partition_row_index[part + 1] =
1922  partition_row_index[part] + partition_l2;
1923  }
1924  }
1925  if (hp_bool == true)
1926  {
1927  partition_partition_list.swap(renumbering);
1928  }
1929  }
1930 
1931 
1932 
1933  // Function to create coloring on the second layer within each partition.
1934  // Version assumes preblocking.
1935  void
1937  const DynamicSparsityPattern & connectivity,
1938  const unsigned int partition,
1939  const std::vector<unsigned int> &cell_partition,
1940  const std::vector<unsigned int> &partition_list,
1941  const std::vector<unsigned int> &partition_size,
1942  std::vector<unsigned int> & partition_color_list)
1943  {
1944  const unsigned int n_macro_cells = *(cell_partition_data.end() - 2);
1945  std::vector<unsigned int> cell_color(n_blocks, n_macro_cells);
1946  std::vector<bool> color_finder;
1947 
1948  partition_row_index.resize(partition + 1);
1949  cell_partition_data.clear();
1950  unsigned int color_counter = 0, index_counter = 0;
1951  for (unsigned int part = 0; part < partition; part++)
1952  {
1953  partition_row_index[part] = index_counter;
1954  unsigned int max_color = 0;
1955  for (unsigned int k = partition_size[part];
1956  k < partition_size[part + 1];
1957  k++)
1958  {
1959  unsigned int cell = partition_list[k];
1960  unsigned int n_neighbors = connectivity.row_length(cell);
1961 
1962  // In the worst case, each neighbor has a different color. So we
1963  // find at least one available color between 0 and n_neighbors.
1964  color_finder.resize(n_neighbors + 1);
1965  for (unsigned int j = 0; j <= n_neighbors; ++j)
1966  color_finder[j] = true;
1968  connectivity.begin(cell),
1969  end = connectivity.end(cell);
1970  for (; neighbor != end; ++neighbor)
1971  {
1972  // Mark the color that a neighbor within the partition has
1973  // as taken
1974  if (cell_partition[neighbor->column()] == part &&
1975  cell_color[neighbor->column()] <= n_neighbors)
1976  color_finder[cell_color[neighbor->column()]] = false;
1977  }
1978  // Choose the smallest color that is not taken for the block
1979  cell_color[cell] = 0;
1980  while (color_finder[cell_color[cell]] == false)
1981  cell_color[cell]++;
1982  if (cell_color[cell] > max_color)
1983  max_color = cell_color[cell];
1984  }
1985  // Reorder within partition: First, all blocks that belong the 0 and
1986  // then so on until those with color max (Note that the smaller the
1987  // number the larger the partition)
1988  for (unsigned int color = 0; color <= max_color; color++)
1989  {
1990  cell_partition_data.push_back(color_counter);
1991  index_counter++;
1992  for (unsigned int k = partition_size[part];
1993  k < partition_size[part + 1];
1994  k++)
1995  {
1996  unsigned int cell = partition_list[k];
1997  if (cell_color[cell] == color)
1998  {
1999  partition_color_list[color_counter++] = cell;
2000  }
2001  }
2002  }
2003  }
2004  cell_partition_data.push_back(n_blocks);
2005  partition_row_index[partition] = index_counter;
2006  AssertDimension(color_counter, n_blocks);
2007  }
2008 
2009 
2010  // Function to create partitioning on the first layer.
2011  void
2013  const unsigned int cluster_size,
2014  std::vector<unsigned int> & cell_partition,
2015  std::vector<unsigned int> & partition_list,
2016  std::vector<unsigned int> & partition_size,
2017  unsigned int & partition) const
2018 
2019  {
2020  // For each block of cells, this variable saves to which partitions the
2021  // block belongs. Initialize all to n_macro_cells to mark them as not
2022  // yet assigned a partition.
2023  // std::vector<unsigned int> cell_partition (n_active_cells,
2024  // numbers::invalid_unsigned_int);
2025  // List of cells in previous partition
2026  std::vector<unsigned int> neighbor_list;
2027  // List of cells in current partition for use as neighbors in next
2028  // partition
2029  std::vector<unsigned int> neighbor_neighbor_list;
2030 
2031  // In element j of this variable, one puts the old number of the block
2032  // that should be the jth block in the new numeration.
2033  // std::vector<unsigned int> partition_list(n_active_cells,0);
2034 
2035  // This vector points to the start of each partition.
2036  // std::vector<unsigned int> partition_size(2,0);
2037 
2038  partition = 0;
2039  unsigned int counter = 0;
2040  unsigned int start_nonboundary =
2041  cell_partition_data.size() == 5 ?
2044  0;
2045 
2046  const unsigned int n_macro_cells = *(cell_partition_data.end() - 2);
2047  if (n_macro_cells == 0)
2048  return;
2049  if (scheme == color)
2050  start_nonboundary = n_macro_cells;
2051  if (scheme == partition_color ||
2052  scheme == color) // blocking_connectivity == true
2053  start_nonboundary = ((start_nonboundary + block_size - 1) / block_size);
2054  unsigned int n_blocks;
2055  if (scheme == partition_color ||
2056  scheme == color) // blocking_connectivity == true
2057  n_blocks = this->n_blocks;
2058  else
2059  n_blocks = n_active_cells;
2060 
2061  if (start_nonboundary > n_blocks)
2062  start_nonboundary = n_blocks;
2063 
2064 
2065  unsigned int start_up = 0;
2066  bool work = true;
2067  unsigned int remainder = cluster_size;
2068 
2069  // this performs a classical breath-first search in the connectivity
2070  // graph of the cells under the restriction that the size of the
2071  // partitions should be a multiple of the given block size
2072  while (work)
2073  {
2074  // put the cells with neighbors on remote MPI processes up front
2075  if (start_nonboundary > 0)
2076  {
2077  for (unsigned int cell = 0; cell < start_nonboundary; ++cell)
2078  {
2079  const unsigned int cell_nn = cell;
2080  cell_partition[cell_nn] = partition;
2081  neighbor_list.push_back(cell_nn);
2082  partition_list[counter++] = cell_nn;
2083  partition_size.back()++;
2084  }
2085  start_nonboundary = 0;
2086  remainder -= (start_nonboundary % cluster_size);
2087  if (remainder == cluster_size)
2088  remainder = 0;
2089  }
2090  else
2091  {
2092  // To start up, set the start_up cell to partition and list all
2093  // its neighbors.
2094  cell_partition[start_up] = partition;
2095  neighbor_list.push_back(start_up);
2096  partition_list[counter++] = start_up;
2097  partition_size.back()++;
2098  start_up++;
2099  remainder--;
2100  if (remainder == cluster_size)
2101  remainder = 0;
2102  }
2103  int index_before = neighbor_list.size(), index = index_before,
2104  index_stop = 0;
2105  while (remainder > 0)
2106  {
2107  if (index == index_stop)
2108  {
2109  index = neighbor_list.size();
2110  if (index == index_before)
2111  {
2112  neighbor_list.resize(0);
2113  goto not_connect;
2114  }
2115  index_stop = index_before;
2116  index_before = index;
2117  }
2118  index--;
2119  unsigned int additional = neighbor_list[index];
2121  connectivity.begin(additional),
2122  end =
2123  connectivity.end(additional);
2124  for (; neighbor != end; ++neighbor)
2125  {
2126  if (cell_partition[neighbor->column()] ==
2128  {
2129  partition_size.back()++;
2130  cell_partition[neighbor->column()] = partition;
2131  neighbor_list.push_back(neighbor->column());
2132  partition_list[counter++] = neighbor->column();
2133  remainder--;
2134  if (remainder == 0)
2135  break;
2136  }
2137  }
2138  }
2139 
2140  while (neighbor_list.size() > 0)
2141  {
2142  partition++;
2143 
2144  // counter for number of cells so far in current partition
2145  unsigned int partition_counter = 0;
2146 
2147  // Mark the start of the new partition
2148  partition_size.push_back(partition_size.back());
2149 
2150  // Loop through the list of cells in previous partition and put
2151  // all their neighbors in current partition
2152  for (const unsigned int cell : neighbor_list)
2153  {
2154  Assert(cell_partition[cell] == partition - 1,
2155  ExcInternalError());
2156  auto neighbor = connectivity.begin(cell);
2157  const auto end = connectivity.end(cell);
2158  for (; neighbor != end; ++neighbor)
2159  {
2160  if (cell_partition[neighbor->column()] ==
2162  {
2163  partition_size.back()++;
2164  cell_partition[neighbor->column()] = partition;
2165 
2166  // collect the cells of the current partition for
2167  // use as neighbors in next partition
2168  neighbor_neighbor_list.push_back(neighbor->column());
2169  partition_list[counter++] = neighbor->column();
2170  partition_counter++;
2171  }
2172  }
2173  }
2174  remainder = cluster_size - (partition_counter % cluster_size);
2175  if (remainder == cluster_size)
2176  remainder = 0;
2177  int index_stop = 0;
2178  int index_before = neighbor_neighbor_list.size(),
2179  index = index_before;
2180  while (remainder > 0)
2181  {
2182  if (index == index_stop)
2183  {
2184  index = neighbor_neighbor_list.size();
2185  if (index == index_before)
2186  {
2187  neighbor_neighbor_list.resize(0);
2188  break;
2189  }
2190  index_stop = index_before;
2191  index_before = index;
2192  }
2193  index--;
2194  unsigned int additional = neighbor_neighbor_list[index];
2196  connectivity.begin(
2197  additional),
2198  end = connectivity.end(
2199  additional);
2200  for (; neighbor != end; ++neighbor)
2201  {
2202  if (cell_partition[neighbor->column()] ==
2204  {
2205  partition_size.back()++;
2206  cell_partition[neighbor->column()] = partition;
2207  neighbor_neighbor_list.push_back(neighbor->column());
2208  partition_list[counter++] = neighbor->column();
2209  remainder--;
2210  if (remainder == 0)
2211  break;
2212  }
2213  }
2214  }
2215 
2216  neighbor_list = neighbor_neighbor_list;
2217  neighbor_neighbor_list.resize(0);
2218  }
2219  not_connect:
2220  // One has to check if the graph is not connected so we have to find
2221  // another partition.
2222  work = false;
2223  for (unsigned int j = start_up; j < n_blocks; ++j)
2224  if (cell_partition[j] == numbers::invalid_unsigned_int)
2225  {
2226  start_up = j;
2227  work = true;
2228  if (remainder == 0)
2229  remainder = cluster_size;
2230  break;
2231  }
2232  }
2233  if (remainder != 0)
2234  partition++;
2235 
2236  AssertDimension(partition_size[partition], n_blocks);
2237  }
2238 
2239 
2240  void
2242  {
2243  evens = (partition + 1) / 2;
2244  odds = partition / 2;
2245  n_blocked_workers = odds - (odds + evens + 1) % 2;
2247  // From here only used for partition partition option.
2248  partition_evens.resize(partition);
2249  partition_odds.resize(partition);
2252  for (unsigned int part = 0; part < partition; part++)
2253  {
2254  partition_evens[part] =
2255  (partition_row_index[part + 1] - partition_row_index[part] + 1) / 2;
2256  partition_odds[part] =
2257  (partition_row_index[part + 1] - partition_row_index[part]) / 2;
2259  partition_odds[part] -
2260  (partition_odds[part] + partition_evens[part] + 1) % 2;
2261  partition_n_workers[part] = partition_evens[part] +
2262  partition_odds[part] -
2264  }
2265  }
2266  } // namespace MatrixFreeFunctions
2267 } // namespace internal
2268 
2269 
2270 
2271 // explicit instantiations of template functions
2272 template void
2273 internal::MatrixFreeFunctions::TaskInfo::print_memory_statistics<std::ostream>(
2274  std::ostream &,
2275  const std::size_t) const;
2276 template void
2278  ConditionalOStream>(ConditionalOStream &, const std::size_t) const;
2279 
2280 
internal::MFWorkerInterface::face
virtual void face(const std::pair< unsigned int, unsigned int > &face_range)=0
internal::MatrixFreeFunctions::TaskInfo::make_thread_graph_partition_color
void make_thread_graph_partition_color(DynamicSparsityPattern &connectivity, std::vector< unsigned int > &renumbering, std::vector< unsigned char > &irregular_cells, const bool hp_bool)
Definition: task_info.cc:1161
internal::MFWorkerInterface::vector_compress_start
virtual void vector_compress_start()=0
Starts the communication for the vector compress operation.
internal::MatrixFreeFunctions::TaskInfo::TaskInfo
TaskInfo()
Definition: task_info.cc:637
internal::MatrixFreeFunctions::TaskInfo::loop
void loop(MFWorkerInterface &worker) const
Definition: task_info.cc:339
internal::MatrixFreeFunctions::partition::PartitionWork::is_blocked
const bool is_blocked
Definition: task_info.cc:224
internal::MatrixFreeFunctions::TaskInfo::make_partitioning
void make_partitioning(const DynamicSparsityPattern &connectivity, const unsigned int cluster_size, std::vector< unsigned int > &cell_partition, std::vector< unsigned int > &partition_list, std::vector< unsigned int > &partition_size, unsigned int &partition) const
Definition: task_info.cc:2012
MultithreadInfo::n_threads
static unsigned int n_threads()
Definition: multithread_info.cc:169
internal::MatrixFreeFunctions::partition::ActualCellWork::task_info
const TaskInfo & task_info
Definition: task_info.cc:102
internal::MatrixFreeFunctions::TaskInfo::boundary_partition_data
std::vector< unsigned int > boundary_partition_data
Definition: task_info.h:491
internal::MatrixFreeFunctions::color::PartitionWork
Definition: task_info.cc:269
internal::MatrixFreeFunctions::partition::CellWork
Definition: task_info.cc:105
internal::MatrixFreeFunctions::TaskInfo::partition_n_blocked_workers
std::vector< unsigned int > partition_n_blocked_workers
Definition: task_info.h:530
internal::MatrixFreeFunctions::TaskInfo::cell_partition_data
std::vector< unsigned int > cell_partition_data
Definition: task_info.h:473
StandardExceptions::ExcNotImplemented
static ::ExceptionBase & ExcNotImplemented()
internal::MatrixFreeFunctions::TaskInfo::communicator
MPI_Comm communicator
Definition: task_info.h:570
DynamicSparsityPattern::begin
iterator begin() const
Definition: dynamic_sparsity_pattern.h:1105
memory_consumption.h
utilities.h
internal::MatrixFreeFunctions::partition::CellWork::execute
tbb::task * execute() override
Definition: task_info.cc:118
internal::MatrixFreeFunctions::TaskInfo::none
@ none
Definition: task_info.h:113
LAPACKSupport::U
static const char U
Definition: lapack_support.h:167
internal::MatrixFreeFunctions::TaskInfo::n_blocked_workers
unsigned int n_blocked_workers
Definition: task_info.h:554
AssertIndexRange
#define AssertIndexRange(index, range)
Definition: exceptions.h:1649
internal::MatrixFreeFunctions::TaskInfo::n_procs
unsigned int n_procs
Definition: task_info.h:580
internal::MatrixFreeFunctions::TaskInfo::make_partitioning_within_partitions_post_blocked
void make_partitioning_within_partitions_post_blocked(const DynamicSparsityPattern &connectivity, const std::vector< unsigned int > &cell_active_fe_index, const unsigned int partition, const unsigned int cluster_size, const bool hp_bool, const std::vector< unsigned int > &cell_partition, const std::vector< unsigned int > &partition_list, const std::vector< unsigned int > &partition_size, std::vector< unsigned int > &partition_partition_list, std::vector< unsigned char > &irregular_cells)
Definition: task_info.cc:1615
internal::MatrixFreeFunctions::color::PartitionWork::partition
const unsigned int partition
Definition: task_info.cc:302
multithread_info.h
Physics::Elasticity::Kinematics::e
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
internal::MatrixFreeFunctions::TaskInfo::my_pid
unsigned int my_pid
Definition: task_info.h:575
internal::MatrixFreeFunctions::partition::CellWork::is_blocked
const bool is_blocked
Definition: task_info.cc:131
internal::MFWorkerInterface
Definition: task_info.h:46
internal::MatrixFreeFunctions::partition::PartitionWork::function
MFWorkerInterface & function
Definition: task_info.cc:221
SparsityTools::partition
void partition(const SparsityPattern &sparsity_pattern, const unsigned int n_partitions, std::vector< unsigned int > &partition_indices, const Partitioner partitioner=Partitioner::metis)
Definition: sparsity_tools.cc:400
internal::MatrixFreeFunctions::partition::ActualCellWork::operator()
void operator()() const
Definition: task_info.cc:77
Utilities::MPI::MinMaxAvg::max
double max
Definition: mpi.h:705
internal::MatrixFreeFunctions::TaskInfo
Definition: task_info.h:105
internal::MatrixFreeFunctions::partition::PartitionWork::execute
tbb::task * execute() override
Definition: task_info.cc:151
internal::MatrixFreeFunctions::TaskInfo::partition_odds
std::vector< unsigned int > partition_odds
Definition: task_info.h:524
internal::MatrixFreeFunctions::MPICommunication
Definition: task_info.cc:311
Utilities::MPI::MinMaxAvg
Definition: mpi.h:687
DynamicSparsityPatternIterators::Accessor::column
size_type column() const
Definition: dynamic_sparsity_pattern.h:812
conditional_ostream.h
internal::MatrixFreeFunctions::TaskInfo::guess_block_size
void guess_block_size(const unsigned int dofs_per_cell)
Definition: task_info.cc:1133
internal::MatrixFreeFunctions::color::CellWork::partition
const unsigned int partition
Definition: task_info.cc:264
internal::MatrixFreeFunctions::TaskInfo::n_ghost_cells
unsigned int n_ghost_cells
Definition: task_info.h:437
DynamicSparsityPatternIterators::Iterator
Definition: dynamic_sparsity_pattern.h:191
DynamicSparsityPattern::row_length
size_type row_length(const size_type row) const
Definition: dynamic_sparsity_pattern.h:1072
internal::MatrixFreeFunctions::TaskInfo::initial_setup_blocks_tasks
void initial_setup_blocks_tasks(const std::vector< unsigned int > &boundary_cells, std::vector< unsigned int > &renumbering, std::vector< unsigned char > &incompletely_filled_vectorization)
Definition: task_info.cc:1067
internal::MatrixFreeFunctions::partition::ActualCellWork
Definition: task_info.cc:55
internal::MatrixFreeFunctions::partition::ActualCellWork::ActualCellWork
ActualCellWork(MFWorkerInterface &worker, const unsigned int partition, const TaskInfo &task_info)
Definition: task_info.cc:67
internal::MatrixFreeFunctions::TaskInfo::evens
unsigned int evens
Definition: task_info.h:542
internal::MatrixFreeFunctions::color::PartitionWork::is_blocked
const bool is_blocked
Definition: task_info.cc:304
internal::MatrixFreeFunctions::TaskInfo::partition_partition
@ partition_partition
Definition: task_info.h:114
internal::MatrixFreeFunctions::TaskInfo::update_task_info
void update_task_info(const unsigned int partition)
Definition: task_info.cc:2241
mpi.h
internal::MatrixFreeFunctions::color::PartitionWork::dummy
tbb::empty_task * dummy
Definition: task_info.cc:298
internal::MatrixFreeFunctions::color::CellWork::operator()
void operator()(const tbb::blocked_range< unsigned int > &r) const
Definition: task_info.cc:245
internal::MatrixFreeFunctions::partition::ActualCellWork::worker
MFWorkerInterface * worker
Definition: task_info.cc:99
Utilities::lower_bound
Iterator lower_bound(Iterator first, Iterator last, const T &val)
Definition: utilities.h:1102
internal::MatrixFreeFunctions::TaskInfo::make_thread_graph
void make_thread_graph(const std::vector< unsigned int > &cell_active_fe_index, DynamicSparsityPattern &connectivity, std::vector< unsigned int > &renumbering, std::vector< unsigned char > &irregular_cells, const bool hp_bool)
Definition: task_info.cc:1308
DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:358
internal::MatrixFreeFunctions::TaskInfo::n_active_cells
unsigned int n_active_cells
Definition: task_info.h:431
internal::MatrixFreeFunctions::TaskInfo::make_thread_graph_partition_partition
void make_thread_graph_partition_partition(const std::vector< unsigned int > &cell_active_fe_index, DynamicSparsityPattern &connectivity, std::vector< unsigned int > &renumbering, std::vector< unsigned char > &irregular_cells, const bool hp_bool)
Definition: task_info.cc:1503
MemoryConsumption::memory_consumption
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
Definition: memory_consumption.h:268
internal::MatrixFreeFunctions::partition::PartitionWork::task_info
const TaskInfo & task_info
Definition: task_info.cc:223
internal::MatrixFreeFunctions::color::PartitionWork::task_info
const TaskInfo & task_info
Definition: task_info.cc:303
internal::MatrixFreeFunctions::TaskInfo::memory_consumption
std::size_t memory_consumption() const
Definition: task_info.cc:690
internal::MatrixFreeFunctions::TaskInfo::print_memory_statistics
void print_memory_statistics(StreamType &out, std::size_t data_length) const
Definition: task_info.cc:675
DynamicSparsityPattern
Definition: dynamic_sparsity_pattern.h:323
internal::VectorImplementation::minimum_parallel_grain_size
unsigned int minimum_parallel_grain_size
Definition: parallel.cc:34
internal::MatrixFreeFunctions::TaskInfo::n_blocks
unsigned int n_blocks
Definition: task_info.h:452
TrilinosWrappers::internal::end
VectorType::value_type * end(VectorType &V)
Definition: trilinos_sparse_matrix.cc:65
internal::MFWorkerInterface::boundary
virtual void boundary(const std::pair< unsigned int, unsigned int > &face_range)=0
internal::MatrixFreeFunctions::MPICommunication::execute
tbb::task * execute() override
Definition: task_info.cc:320
AssertDimension
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1579
MemorySpace::swap
void swap(MemorySpaceData< Number, MemorySpace > &, MemorySpaceData< Number, MemorySpace > &)
Definition: memory_space.h:103
internal::MatrixFreeFunctions::color::PartitionWork::execute
tbb::task * execute() override
Definition: task_info.cc:284
internal::MatrixFreeFunctions::partition::CellWork::CellWork
CellWork(MFWorkerInterface &worker, const unsigned int partition, const TaskInfo &task_info, const bool is_blocked)
Definition: task_info.cc:108
internal::MatrixFreeFunctions::color::PartitionWork::PartitionWork
PartitionWork(MFWorkerInterface &worker_in, const unsigned int partition_in, const TaskInfo &task_info_in, const bool is_blocked_in)
Definition: task_info.cc:272
task_info.h
internal::MatrixFreeFunctions::partition::CellWork::dummy
tbb::empty_task * dummy
Definition: task_info.cc:127
internal::MatrixFreeFunctions::TaskInfo::n_workers
unsigned int n_workers
Definition: task_info.h:559
internal::MatrixFreeFunctions::partition::PartitionWork
Definition: task_info.cc:136
internal::MFWorkerInterface::vector_update_ghosts_finish
virtual void vector_update_ghosts_finish()=0
Finishes the communication for the update ghost values operation.
internal::MatrixFreeFunctions::color::CellWork::CellWork
CellWork(MFWorkerInterface &worker_in, const TaskInfo &task_info_in, const unsigned int partition_in)
Definition: task_info.cc:236
internal::MatrixFreeFunctions::partition::CellWork::work
ActualCellWork work
Definition: task_info.cc:130
internal::MatrixFreeFunctions::MPICommunication::worker
MFWorkerInterface & worker
Definition: task_info.cc:330
internal::MatrixFreeFunctions::color::CellWork::worker
MFWorkerInterface & worker
Definition: task_info.cc:262
internal::dummy
const types::global_dof_index * dummy()
Definition: dof_handler.cc:59
internal::MatrixFreeFunctions::partition::PartitionWork::dummy
tbb::empty_task * dummy
Definition: task_info.cc:218
StandardExceptions::ExcInternalError
static ::ExceptionBase & ExcInternalError()
internal::MFWorkerInterface::cell_loop_pre_range
virtual void cell_loop_pre_range(const unsigned int range_index)=0
internal::MatrixFreeFunctions::partition::PartitionWork::PartitionWork
PartitionWork(MFWorkerInterface &function_in, const unsigned int partition_in, const TaskInfo &task_info_in, const bool is_blocked_in=false)
Definition: task_info.cc:139
Assert
#define Assert(cond, exc)
Definition: exceptions.h:1419
internal::MatrixFreeFunctions::partition::PartitionWork::partition
const unsigned int partition
Definition: task_info.cc:222
Utilities::MPI::MinMaxAvg::avg
double avg
Definition: mpi.h:731
internal::MFWorkerInterface::zero_dst_vector_range
virtual void zero_dst_vector_range(const unsigned int range_index)=0
internal::MFWorkerInterface::cell_loop_post_range
virtual void cell_loop_post_range(const unsigned int range_index)=0
internal::MFWorkerInterface::cell
virtual void cell(const std::pair< unsigned int, unsigned int > &cell_range)=0
internal::MatrixFreeFunctions::TaskInfo::face_partition_data
std::vector< unsigned int > face_partition_data
Definition: task_info.h:482
internal::MatrixFreeFunctions::TaskInfo::odds
unsigned int odds
Definition: task_info.h:548
DynamicSparsityPattern::end
iterator end() const
Definition: dynamic_sparsity_pattern.h:1115
internal::MatrixFreeFunctions::TaskInfo::partition_row_index
std::vector< unsigned int > partition_row_index
Definition: task_info.h:465
ConditionalOStream
Definition: conditional_ostream.h:81
numbers::invalid_unsigned_int
static const unsigned int invalid_unsigned_int
Definition: types.h:191
internal::MFWorkerInterface::vector_compress_finish
virtual void vector_compress_finish()=0
Finishes the communication for the vector compress operation.
Utilities::MPI::min_max_avg
MinMaxAvg min_max_avg(const double my_value, const MPI_Comm &mpi_communicator)
Definition: mpi.cc:91
Utilities::MPI::min
T min(const T &t, const MPI_Comm &mpi_communicator)
internal::MatrixFreeFunctions::color::CellWork
Definition: task_info.cc:233
internal::MatrixFreeFunctions::TaskInfo::vectorization_length
unsigned int vectorization_length
Definition: task_info.h:442
internal::MatrixFreeFunctions::TaskInfo::make_coloring_within_partitions_pre_blocked
void make_coloring_within_partitions_pre_blocked(const DynamicSparsityPattern &connectivity, const unsigned int partition, const std::vector< unsigned int > &cell_partition, const std::vector< unsigned int > &partition_list, const std::vector< unsigned int > &partition_size, std::vector< unsigned int > &partition_color_list)
Definition: task_info.cc:1936
internal::MatrixFreeFunctions::TaskInfo::partition_color
@ partition_color
Definition: task_info.h:115
internal::MatrixFreeFunctions::MPICommunication::do_compress
const bool do_compress
Definition: task_info.cc:331
internal::MatrixFreeFunctions::color::PartitionWork::worker
MFWorkerInterface & worker
Definition: task_info.cc:301
internal::MatrixFreeFunctions::partition::ActualCellWork::partition
const unsigned int partition
Definition: task_info.cc:101
internal
Definition: aligned_vector.h:369
internal::MatrixFreeFunctions::color::CellWork::task_info
const TaskInfo & task_info
Definition: task_info.cc:263
internal::MatrixFreeFunctions::TaskInfo::color
@ color
Definition: task_info.h:116
internal::MatrixFreeFunctions::partition::ActualCellWork::ActualCellWork
ActualCellWork(MFWorkerInterface **worker_pointer, const unsigned int partition, const TaskInfo &task_info)
Definition: task_info.cc:58
parallel.h
Utilities::invert_permutation
std::vector< Integer > invert_permutation(const std::vector< Integer > &permutation)
Definition: utilities.h:1446
DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:359
Utilities::MPI::MinMaxAvg::min
double min
Definition: mpi.h:699
internal::MatrixFreeFunctions::TaskInfo::block_size
unsigned int block_size
Definition: task_info.h:447
first
Point< 2 > first
Definition: grid_out.cc:4352
DynamicSparsityPattern::add
void add(const size_type i, const size_type j)
Definition: dynamic_sparsity_pattern.h:1032
internal::MatrixFreeFunctions::TaskInfo::clear
void clear()
Definition: task_info.cc:645
internal::MatrixFreeFunctions::TaskInfo::create_blocks_serial
void create_blocks_serial(const std::vector< unsigned int > &cells_with_comm, const unsigned int dofs_per_cell, const bool categories_are_hp, const std::vector< unsigned int > &cell_vectorization_categories, const bool cell_vectorization_categories_strict, const std::vector< unsigned int > &parent_relation, std::vector< unsigned int > &renumbering, std::vector< unsigned char > &incompletely_filled_vectorization)
Definition: task_info.cc:779
internal::MatrixFreeFunctions::TaskInfo::make_boundary_cells_divisible
void make_boundary_cells_divisible(std::vector< unsigned int > &boundary_cells)
Definition: task_info.cc:707
AssertThrow
#define AssertThrow(cond, exc)
Definition: exceptions.h:1531
internal::MatrixFreeFunctions::TaskInfo::partition_evens
std::vector< unsigned int > partition_evens
Definition: task_info.h:518
internal::MatrixFreeFunctions::TaskInfo::scheme
TasksParallelScheme scheme
Definition: task_info.h:457
internal::MatrixFreeFunctions::partition::ActualCellWork::worker_pointer
MFWorkerInterface ** worker_pointer
Definition: task_info.cc:100
internal::MatrixFreeFunctions::TaskInfo::partition_n_workers
std::vector< unsigned int > partition_n_workers
Definition: task_info.h:536
internal::MatrixFreeFunctions::MPICommunication::MPICommunication
MPICommunication(MFWorkerInterface &worker_in, const bool do_compress)
Definition: task_info.cc:314
internal::MFWorkerInterface::vector_update_ghosts_start
virtual void vector_update_ghosts_start()=0
Starts the communication for the update ghost values operation.
Utilities::MPI::max
T max(const T &t, const MPI_Comm &mpi_communicator)
internal::MatrixFreeFunctions::TaskInfo::make_connectivity_cells_to_blocks
void make_connectivity_cells_to_blocks(const std::vector< unsigned char > &irregular_cells, const DynamicSparsityPattern &connectivity_cells, DynamicSparsityPattern &connectivity_blocks) const
Definition: task_info.cc:1571
parallel::internal::parallel_for
void parallel_for(Iterator x_begin, Iterator x_end, const Functor &functor, const unsigned int grainsize)
Definition: parallel.h:156