Reference documentation for deal.II version 8.4.2
dof_handler_policy.cc
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2016 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 at
12 // the top level of the deal.II distribution.
13 //
14 // ---------------------------------------------------------------------
15 
16 
17 //TODO [TH]: renumber DoFs for multigrid is not done yet
18 
19 #include <deal.II/base/geometry_info.h>
20 #include <deal.II/base/utilities.h>
21 #include <deal.II/base/memory_consumption.h>
22 #include <deal.II/grid/tria.h>
23 #include <deal.II/grid/tria_iterator.h>
24 #include <deal.II/dofs/dof_handler.h>
25 #include <deal.II/dofs/dof_accessor.h>
26 #include <deal.II/dofs/dof_handler_policy.h>
27 #include <deal.II/fe/fe.h>
28 #include <deal.II/distributed/shared_tria.h>
29 #include <deal.II/distributed/tria.h>
30 
31 #include <set>
32 #include <algorithm>
33 #include <numeric>
34 
35 DEAL_II_NAMESPACE_OPEN
36 
37 
38 namespace internal
39 {
40  namespace DoFHandler
41  {
42  namespace Policy
43  {
44  // use class ::DoFHandler instead
45  // of namespace internal::DoFHandler in
46  // the following
47  using ::DoFHandler;
48 
49  struct Implementation
50  {
51 
52  /* -------------- distribute_dofs functionality ------------- */
53 
65  template <int spacedim>
66  static
67  types::global_dof_index
68  distribute_dofs_on_cell (const DoFHandler<1,spacedim> &dof_handler,
69  const typename DoFHandler<1,spacedim>::active_cell_iterator &cell,
70  types::global_dof_index next_free_dof)
71  {
72 
73  // distribute dofs of vertices
74  if (dof_handler.get_fe().dofs_per_vertex > 0)
75  for (unsigned int v=0; v<GeometryInfo<1>::vertices_per_cell; ++v)
76  {
77  if (cell->vertex_dof_index (v,0) ==
78  DoFHandler<1,spacedim>::invalid_dof_index)
79  for (unsigned int d=0;
80  d<dof_handler.get_fe().dofs_per_vertex; ++d)
81  {
82  Assert ((cell->vertex_dof_index (v,d) ==
83  DoFHandler<1,spacedim>::invalid_dof_index),
84  ExcInternalError());
85  cell->set_vertex_dof_index (v, d, next_free_dof++);
86  }
87  else
88  for (unsigned int d=0;
89  d<dof_handler.get_fe().dofs_per_vertex; ++d)
90  Assert ((cell->vertex_dof_index (v,d) !=
91  DoFHandler<1,spacedim>::invalid_dof_index),
92  ExcInternalError());
93  }
94 
95  // dofs of line
96  for (unsigned int d=0;
97  d<dof_handler.get_fe().dofs_per_line; ++d)
98  cell->set_dof_index (d, next_free_dof++);
99 
100  // note that this cell has been
101  // processed
102  cell->set_user_flag ();
103 
104  return next_free_dof;
105  }
106 
107 
108 
109  template <int spacedim>
110  static
111  types::global_dof_index
112  distribute_dofs_on_cell (const DoFHandler<2,spacedim> &dof_handler,
113  const typename DoFHandler<2,spacedim>::active_cell_iterator &cell,
114  types::global_dof_index next_free_dof)
115  {
116  if (dof_handler.get_fe().dofs_per_vertex > 0)
117  // number dofs on vertices
118  for (unsigned int vertex=0; vertex<GeometryInfo<2>::vertices_per_cell; ++vertex)
119  // check whether dofs for this
120  // vertex have been distributed
121  // (only check the first dof)
122  if (cell->vertex_dof_index(vertex, 0) == DoFHandler<2,spacedim>::invalid_dof_index)
123  for (unsigned int d=0; d<dof_handler.get_fe().dofs_per_vertex; ++d)
124  cell->set_vertex_dof_index (vertex, d, next_free_dof++);
125 
126  // for the four sides
127  if (dof_handler.get_fe().dofs_per_line > 0)
128  for (unsigned int side=0; side<GeometryInfo<2>::faces_per_cell; ++side)
129  {
130  const typename DoFHandler<2,spacedim>::line_iterator
131  line = cell->line(side);
132 
133  // distribute dofs if necessary:
134  // check whether line dof is already
135  // numbered (check only first dof)
136  if (line->dof_index(0) == DoFHandler<2,spacedim>::invalid_dof_index)
137  // if not: distribute dofs
138  for (unsigned int d=0; d<dof_handler.get_fe().dofs_per_line; ++d)
139  line->set_dof_index (d, next_free_dof++);
140  }
141 
142 
143  // dofs of quad
144  if (dof_handler.get_fe().dofs_per_quad > 0)
145  for (unsigned int d=0; d<dof_handler.get_fe().dofs_per_quad; ++d)
146  cell->set_dof_index (d, next_free_dof++);
147 
148 
149  // note that this cell has been processed
150  cell->set_user_flag ();
151 
152  return next_free_dof;
153  }
154 
155 
156  template <int spacedim>
157  static
158  types::global_dof_index
159  distribute_dofs_on_cell (const DoFHandler<3,spacedim> &dof_handler,
160  const typename DoFHandler<3,spacedim>::active_cell_iterator &cell,
161  types::global_dof_index next_free_dof)
162  {
163  if (dof_handler.get_fe().dofs_per_vertex > 0)
164  // number dofs on vertices
165  for (unsigned int vertex=0; vertex<GeometryInfo<3>::vertices_per_cell; ++vertex)
166  // check whether dofs for this
167  // vertex have been distributed
168  // (only check the first dof)
169  if (cell->vertex_dof_index(vertex, 0) == DoFHandler<3,spacedim>::invalid_dof_index)
170  for (unsigned int d=0; d<dof_handler.get_fe().dofs_per_vertex; ++d)
171  cell->set_vertex_dof_index (vertex, d, next_free_dof++);
172 
173  // for the lines
174  if (dof_handler.get_fe().dofs_per_line > 0)
175  for (unsigned int l=0; l<GeometryInfo<3>::lines_per_cell; ++l)
176  {
177  const typename DoFHandler<3,spacedim>::line_iterator
178  line = cell->line(l);
179 
180  // distribute dofs if necessary:
181  // check whether line dof is already
182  // numbered (check only first dof)
183  if (line->dof_index(0) == DoFHandler<3,spacedim>::invalid_dof_index)
184  // if not: distribute dofs
185  for (unsigned int d=0; d<dof_handler.get_fe().dofs_per_line; ++d)
186  line->set_dof_index (d, next_free_dof++);
187  }
188 
189  // for the quads
190  if (dof_handler.get_fe().dofs_per_quad > 0)
191  for (unsigned int q=0; q<GeometryInfo<3>::quads_per_cell; ++q)
192  {
193  const typename DoFHandler<3,spacedim>::quad_iterator
194  quad = cell->quad(q);
195 
196  // distribute dofs if necessary:
197  // check whether quad dof is already
198  // numbered (check only first dof)
199  if (quad->dof_index(0) == DoFHandler<3,spacedim>::invalid_dof_index)
200  // if not: distribute dofs
201  for (unsigned int d=0; d<dof_handler.get_fe().dofs_per_quad; ++d)
202  quad->set_dof_index (d, next_free_dof++);
203  }
204 
205 
206  // dofs of hex
207  if (dof_handler.get_fe().dofs_per_hex > 0)
208  for (unsigned int d=0; d<dof_handler.get_fe().dofs_per_hex; ++d)
209  cell->set_dof_index (d, next_free_dof++);
210 
211 
212  // note that this cell has been
213  // processed
214  cell->set_user_flag ();
215 
216  return next_free_dof;
217  }
218 
219 
226  template <int dim, int spacedim>
227  static
228  types::global_dof_index
229  distribute_dofs (const types::global_dof_index offset,
230  const types::subdomain_id subdomain_id,
231  DoFHandler<dim,spacedim> &dof_handler)
232  {
233  const ::Triangulation<dim,spacedim> &tria
234  = dof_handler.get_triangulation();
235  Assert (tria.n_levels() > 0, ExcMessage("Empty triangulation"));
236 
237  // Clear user flags because we will need them. But first we save
238  // them and make sure that we restore them later such that at the
239  // end of this function the Triangulation will be in the same state
240  // as it was at the beginning of this function.
241  std::vector<bool> user_flags;
242  tria.save_user_flags(user_flags);
243  const_cast<::Triangulation<dim,spacedim> &>(tria).clear_user_flags ();
244 
245  types::global_dof_index next_free_dof = offset;
246  typename DoFHandler<dim,spacedim>::active_cell_iterator
247  cell = dof_handler.begin_active(),
248  endc = dof_handler.end();
249 
250  for (; cell != endc; ++cell)
251  if ((subdomain_id == numbers::invalid_subdomain_id)
252  ||
253  (cell->subdomain_id() == subdomain_id))
254  next_free_dof
255  = Implementation::distribute_dofs_on_cell (dof_handler,
256  cell,
257  next_free_dof);
258 
259  // update the cache used for cell dof indices
260  for (cell = dof_handler.begin_active(); cell != endc; ++cell)
261  if (!cell->is_artificial())
262  cell->update_cell_dof_indices_cache ();
263 
264  // finally restore the user flags
265  const_cast<::Triangulation<dim,spacedim> &>(tria).load_user_flags(user_flags);
266 
267  return next_free_dof;
268  }
269 
270 
289  // These three function
290  // have an unused
291  // DoFHandler object as
292  // their first
293  // argument. Without it,
294  // the file was not
295  // compileable under gcc
296  // 4.4.5 (Debian).
297  template <int spacedim>
298  static
299  unsigned int
300  distribute_mg_dofs_on_cell (const DoFHandler<1,spacedim> &,
301  typename DoFHandler<1,spacedim>::level_cell_iterator &cell,
302  unsigned int next_free_dof)
303  {
304  const unsigned int dim = 1;
305 
306  // distribute dofs of vertices
307  if (cell->get_fe().dofs_per_vertex > 0)
308  for (unsigned int v=0; v<GeometryInfo<1>::vertices_per_cell; ++v)
309  {
310  typename DoFHandler<dim,spacedim>::level_cell_iterator neighbor = cell->neighbor(v);
311 
312  if (neighbor.state() == IteratorState::valid)
313  {
314  // has neighbor already been processed?
315  if (neighbor->user_flag_set() &&
316  (neighbor->level() == cell->level()))
317  // copy dofs if the neighbor is on
318  // the same level (only then are
319  // mg dofs the same)
320  {
321  if (v==0)
322  for (unsigned int d=0; d<cell->get_fe().dofs_per_vertex; ++d)
323  cell->set_mg_vertex_dof_index (cell->level(), 0, d,
324  neighbor->mg_vertex_dof_index (cell->level(), 1, d));
325  else
326  for (unsigned int d=0; d<cell->get_fe().dofs_per_vertex; ++d)
327  cell->set_mg_vertex_dof_index (cell->level(), 1, d,
328  neighbor->mg_vertex_dof_index (cell->level(), 0, d));
329 
330  // next neighbor
331  continue;
332  };
333  };
334 
335  // otherwise: create dofs newly
336  for (unsigned int d=0; d<cell->get_fe().dofs_per_vertex; ++d)
337  cell->set_mg_vertex_dof_index (cell->level(), v, d, next_free_dof++);
338  };
339 
340  // dofs of line
341  if (cell->get_fe().dofs_per_line > 0)
342  for (unsigned int d=0; d<cell->get_fe().dofs_per_line; ++d)
343  cell->set_mg_dof_index (cell->level(), d, next_free_dof++);
344 
345  // note that this cell has been processed
346  cell->set_user_flag ();
347 
348  return next_free_dof;
349  }
350 
351 
352  template <int spacedim>
353  static
354  unsigned int
355  distribute_mg_dofs_on_cell (const DoFHandler<2,spacedim> &,
356  typename DoFHandler<2,spacedim>::level_cell_iterator &cell,
357  unsigned int next_free_dof)
358  {
359  const unsigned int dim = 2;
360  if (cell->get_fe().dofs_per_vertex > 0)
361  // number dofs on vertices
362  for (unsigned int vertex=0; vertex<GeometryInfo<2>::vertices_per_cell; ++vertex)
363  // check whether dofs for this
364  // vertex have been distributed
365  // (only check the first dof)
366  if (cell->mg_vertex_dof_index(cell->level(), vertex, 0) == DoFHandler<2>::invalid_dof_index)
367  for (unsigned int d=0; d<cell->get_fe().dofs_per_vertex; ++d)
368  cell->set_mg_vertex_dof_index (cell->level(), vertex, d, next_free_dof++);
369 
370  // for the four sides
371  if (cell->get_fe().dofs_per_line > 0)
372  for (unsigned int side=0; side<GeometryInfo<2>::faces_per_cell; ++side)
373  {
374  typename DoFHandler<dim,spacedim>::line_iterator line = cell->line(side);
375 
376  // distribute dofs if necessary:
377  // check whether line dof is already
378  // numbered (check only first dof)
379  if (line->mg_dof_index(cell->level(), 0) == DoFHandler<2>::invalid_dof_index)
380  // if not: distribute dofs
381  for (unsigned int d=0; d<cell->get_fe().dofs_per_line; ++d)
382  line->set_mg_dof_index (cell->level(), d, next_free_dof++);
383  };
384 
385 
386  // dofs of quad
387  if (cell->get_fe().dofs_per_quad > 0)
388  for (unsigned int d=0; d<cell->get_fe().dofs_per_quad; ++d)
389  cell->set_mg_dof_index (cell->level(), d, next_free_dof++);
390 
391 
392  // note that this cell has been processed
393  cell->set_user_flag ();
394 
395  return next_free_dof;
396  }
397 
398 
399  template <int spacedim>
400  static
401  unsigned int
402  distribute_mg_dofs_on_cell (const DoFHandler<3,spacedim> &,
403  typename DoFHandler<3,spacedim>::level_cell_iterator &cell,
404  unsigned int next_free_dof)
405  {
406  const unsigned int dim = 3;
407  if (cell->get_fe().dofs_per_vertex > 0)
408  // number dofs on vertices
409  for (unsigned int vertex=0; vertex<GeometryInfo<3>::vertices_per_cell; ++vertex)
410  // check whether dofs for this
411  // vertex have been distributed
412  // (only check the first dof)
413  if (cell->mg_vertex_dof_index(cell->level(), vertex, 0) == DoFHandler<3>::invalid_dof_index)
414  for (unsigned int d=0; d<cell->get_fe().dofs_per_vertex; ++d)
415  cell->set_mg_vertex_dof_index (cell->level(), vertex, d, next_free_dof++);
416 
417  // for the lines
418  if (cell->get_fe().dofs_per_line > 0)
419  for (unsigned int l=0; l<GeometryInfo<3>::lines_per_cell; ++l)
420  {
421  typename DoFHandler<dim,spacedim>::line_iterator line = cell->line(l);
422 
423  // distribute dofs if necessary:
424  // check whether line dof is already
425  // numbered (check only first dof)
426  if (line->mg_dof_index(cell->level(), 0) == DoFHandler<3>::invalid_dof_index)
427  // if not: distribute dofs
428  for (unsigned int d=0; d<cell->get_fe().dofs_per_line; ++d)
429  line->set_mg_dof_index (cell->level(), d, next_free_dof++);
430  };
431 
432  // for the quads
433  if (cell->get_fe().dofs_per_quad > 0)
434  for (unsigned int q=0; q<GeometryInfo<3>::quads_per_cell; ++q)
435  {
436  typename DoFHandler<dim,spacedim>::quad_iterator quad = cell->quad(q);
437 
438  // distribute dofs if necessary:
439  // check whether line dof is already
440  // numbered (check only first dof)
441  if (quad->mg_dof_index(cell->level(), 0) == DoFHandler<3>::invalid_dof_index)
442  // if not: distribute dofs
443  for (unsigned int d=0; d<cell->get_fe().dofs_per_quad; ++d)
444  quad->set_mg_dof_index (cell->level(), d, next_free_dof++);
445  };
446 
447 
448  // dofs of cell
449  if (cell->get_fe().dofs_per_hex > 0)
450  for (unsigned int d=0; d<cell->get_fe().dofs_per_hex; ++d)
451  cell->set_mg_dof_index (cell->level(), d, next_free_dof++);
452 
453 
454  // note that this cell has
455  // been processed
456  cell->set_user_flag ();
457 
458  return next_free_dof;
459  }
460 
461 
462  template <int dim, int spacedim>
463  static
464  unsigned int
465  distribute_dofs_on_level (const unsigned int offset,
466  const types::subdomain_id level_subdomain_id,
467  DoFHandler<dim,spacedim> &dof_handler,
468  const unsigned int level)
469  {
470  const ::Triangulation<dim,spacedim> &tria
471  = dof_handler.get_triangulation();
472  Assert (tria.n_levels() > 0, ExcMessage("Empty triangulation"));
473  if (level>=tria.n_levels())
474  return 0; //this is allowed for multigrid
475 
476  // Clear user flags because we will
477  // need them. But first we save
478  // them and make sure that we
479  // restore them later such that at
480  // the end of this function the
481  // Triangulation will be in the
482  // same state as it was at the
483  // beginning of this function.
484  std::vector<bool> user_flags;
485  tria.save_user_flags(user_flags);
486  const_cast<::Triangulation<dim,spacedim> &>(tria).clear_user_flags ();
487 
488  unsigned int next_free_dof = offset;
489  typename DoFHandler<dim,spacedim>::level_cell_iterator
490  cell = dof_handler.begin(level),
491  endc = dof_handler.end(level);
492 
493  for (; cell != endc; ++cell)
494  if ((level_subdomain_id == numbers::invalid_subdomain_id)
495  ||
496  (cell->level_subdomain_id() == level_subdomain_id))
497  next_free_dof
498  = Implementation::distribute_mg_dofs_on_cell (dof_handler, cell, next_free_dof);
499 
500 // // update the cache used
501 // // for cell dof indices
502 // for (typename DoFHandler<dim,spacedim>::level_cell_iterator
503 // cell = dof_handler.begin(); cell != dof_handler.end(); ++cell)
504 // if (cell->subdomain_id() != numbers::artificial_subdomain_id)
505 // cell->update_cell_dof_indices_cache ();
506 
507  // finally restore the user flags
508  const_cast<::Triangulation<dim,spacedim> &>(tria).load_user_flags(user_flags);
509 
510  return next_free_dof;
511  }
512 
513 
514  /* --------------------- renumber_dofs functionality ---------------- */
515 
516 
533  template <int spacedim>
534  static
535  void
536  renumber_dofs (const std::vector<types::global_dof_index> &new_numbers,
537  const IndexSet &,
538  DoFHandler<1,spacedim> &dof_handler,
539  const bool check_validity)
540  {
541  // note that we can not use cell
542  // iterators in this function since
543  // then we would renumber the dofs on
544  // the interface of two cells more
545  // than once. Anyway, this way it's
546  // not only more correct but also
547  // faster; note, however, that dof
548  // numbers may be invalid_dof_index,
549  // namely when the appropriate
550  // vertex/line/etc is unused
551  for (std::vector<types::global_dof_index>::iterator
552  i=dof_handler.vertex_dofs.begin();
553  i!=dof_handler.vertex_dofs.end(); ++i)
554  if (*i != DoFHandler<1,spacedim>::invalid_dof_index)
555  *i = new_numbers[*i];
556  else if (check_validity)
557  // if index is
558  // invalid_dof_index:
559  // check if this one
560  // really is unused
561  Assert (dof_handler.get_triangulation()
562  .vertex_used((i-dof_handler.vertex_dofs.begin()) /
563  dof_handler.selected_fe->dofs_per_vertex)
564  == false,
565  ExcInternalError ());
566 
567  for (unsigned int level=0; level<dof_handler.levels.size(); ++level)
568  for (std::vector<types::global_dof_index>::iterator
569  i=dof_handler.levels[level]->dof_object.dofs.begin();
570  i!=dof_handler.levels[level]->dof_object.dofs.end(); ++i)
571  if (*i != DoFHandler<1,spacedim>::invalid_dof_index)
572  *i = new_numbers[*i];
573 
574  // update the cache
575  // used for cell dof
576  // indices
577  for (typename DoFHandler<1,spacedim>::level_cell_iterator
578  cell = dof_handler.begin();
579  cell != dof_handler.end(); ++cell)
580  cell->update_cell_dof_indices_cache ();
581  }
582 
583  template <int spacedim>
584  static
585  void
586  renumber_mg_dofs (const std::vector<::types::global_dof_index> &new_numbers,
587  const IndexSet &indices,
588  DoFHandler<1,spacedim> &dof_handler,
589  const unsigned int level,
590  const bool check_validity)
591  {
592  for (typename std::vector<typename DoFHandler<1,spacedim>::MGVertexDoFs>::iterator
593  i=dof_handler.mg_vertex_dofs.begin();
594  i!=dof_handler.mg_vertex_dofs.end();
595  ++i)
596  // if the present vertex lives on
597  // the current level
598  if ((i->get_coarsest_level() <= level) &&
599  (i->get_finest_level() >= level))
600  for (unsigned int d=0; d<dof_handler.get_fe().dofs_per_vertex; ++d)
601  {
602  ::types::global_dof_index idx = i->get_index (level, d);
603  if (idx != DoFHandler<1>::invalid_dof_index)
604  i->set_index (level, d,
605  (indices.n_elements() == 0)?
606  (new_numbers[idx]) :
607  (new_numbers[indices.index_within_set(idx)]));
608 
609  if (check_validity)
610  Assert(idx != DoFHandler<1>::invalid_dof_index, ExcInternalError ());
611  }
612 
613 
614  for (std::vector<types::global_dof_index>::iterator
615  i=dof_handler.mg_levels[level]->dof_object.dofs.begin();
616  i!=dof_handler.mg_levels[level]->dof_object.dofs.end();
617  ++i)
618  {
619  if (*i != DoFHandler<1>::invalid_dof_index)
620  {
621  Assert(*i<new_numbers.size(), ExcInternalError());
622  *i = (indices.n_elements() == 0)?
623  (new_numbers[*i]) :
624  (new_numbers[indices.index_within_set(*i)]);
625  }
626  }
627  }
628 
629 
630  template <int spacedim>
631  static
632  void
633  renumber_dofs (const std::vector<types::global_dof_index> &new_numbers,
634  const IndexSet &indices,
635  DoFHandler<2,spacedim> &dof_handler,
636  const bool check_validity)
637  {
638  // note that we can not use cell
639  // iterators in this function since
640  // then we would renumber the dofs on
641  // the interface of two cells more
642  // than once. Anyway, this way it's
643  // not only more correct but also
644  // faster; note, however, that dof
645  // numbers may be invalid_dof_index,
646  // namely when the appropriate
647  // vertex/line/etc is unused
648  for (std::vector<types::global_dof_index>::iterator
649  i=dof_handler.vertex_dofs.begin();
650  i!=dof_handler.vertex_dofs.end(); ++i)
651  if (*i != DoFHandler<2,spacedim>::invalid_dof_index)
652  *i = (indices.n_elements() == 0)?
653  (new_numbers[*i]) :
654  (new_numbers[indices.index_within_set(*i)]);
655  else if (check_validity)
656  // if index is invalid_dof_index:
657  // check if this one really is
658  // unused
659  Assert (dof_handler.get_triangulation()
660  .vertex_used((i-dof_handler.vertex_dofs.begin()) /
661  dof_handler.selected_fe->dofs_per_vertex)
662  == false,
663  ExcInternalError ());
664 
665  for (std::vector<types::global_dof_index>::iterator
666  i=dof_handler.faces->lines.dofs.begin();
667  i!=dof_handler.faces->lines.dofs.end(); ++i)
668  if (*i != DoFHandler<2,spacedim>::invalid_dof_index)
669  *i = ((indices.n_elements() == 0) ?
670  new_numbers[*i] :
671  new_numbers[indices.index_within_set(*i)]);
672 
673  for (unsigned int level=0; level<dof_handler.levels.size(); ++level)
674  {
675  for (std::vector<types::global_dof_index>::iterator
676  i=dof_handler.levels[level]->dof_object.dofs.begin();
677  i!=dof_handler.levels[level]->dof_object.dofs.end(); ++i)
678  if (*i != DoFHandler<2,spacedim>::invalid_dof_index)
679  *i = ((indices.n_elements() == 0) ?
680  new_numbers[*i] :
681  new_numbers[indices.index_within_set(*i)]);
682  }
683 
684  // update the cache
685  // used for cell dof
686  // indices
687  for (typename DoFHandler<2,spacedim>::level_cell_iterator
688  cell = dof_handler.begin();
689  cell != dof_handler.end(); ++cell)
690  cell->update_cell_dof_indices_cache ();
691  }
692 
693  template <int spacedim>
694  static
695  void
696  renumber_mg_dofs (const std::vector<::types::global_dof_index> &new_numbers,
697  const IndexSet &indices,
698  DoFHandler<2,spacedim> &dof_handler,
699  const unsigned int level,
700  const bool check_validity)
701  {
702  if (level>=dof_handler.get_triangulation().n_levels())
703  return;
704  for (typename std::vector<typename DoFHandler<2,spacedim>::MGVertexDoFs>::iterator i=dof_handler.mg_vertex_dofs.begin();
705  i!=dof_handler.mg_vertex_dofs.end(); ++i)
706  // if the present vertex lives on
707  // the present level
708  if ((i->get_coarsest_level() <= level) &&
709  (i->get_finest_level() >= level))
710  for (unsigned int d=0; d<dof_handler.get_fe().dofs_per_vertex; ++d)
711  {
712  ::types::global_dof_index idx =i->get_index (level, d/*, dof_handler.get_fe().dofs_per_vertex*/);
713  if (idx != DoFHandler<1>::invalid_dof_index)
714  i->set_index (level, d/*, dof_handler.get_fe().dofs_per_vertex*/,
715  ((indices.n_elements() == 0) ?
716  new_numbers[idx] :
717  new_numbers[indices.index_within_set(idx)]));
718 
719  if (check_validity)
720  Assert(idx != DoFHandler<2>::invalid_dof_index, ExcInternalError ());
721  }
722 
723  if (dof_handler.get_fe().dofs_per_line > 0)
724  {
725  // save user flags as they will be modified
726  std::vector<bool> user_flags;
727  dof_handler.get_triangulation().save_user_flags(user_flags);
728  const_cast<::Triangulation<2,spacedim> &>(dof_handler.get_triangulation()).clear_user_flags ();
729 
730  // flag all lines adjacent to cells of the current
731  // level, as those lines logically belong to the same
732  // level as the cell, at least for for isotropic
733  // refinement
734  for (typename DoFHandler<2,spacedim>::level_cell_iterator cell = dof_handler.begin(level);
735  cell != dof_handler.end(level); ++cell)
736  for (unsigned int line=0; line < GeometryInfo<2>::faces_per_cell; ++line)
737  cell->face(line)->set_user_flag();
738 
739  for (typename DoFHandler<2,spacedim>::cell_iterator cell = dof_handler.begin();
740  cell != dof_handler.end(); ++cell)
741  for (unsigned int l=0; l<GeometryInfo<2>::lines_per_cell; ++l)
742  if (cell->line(l)->user_flag_set())
743  {
744  for (unsigned int d=0; d<dof_handler.get_fe().dofs_per_line; ++d)
745  {
746  ::types::global_dof_index idx = cell->line(l)->mg_dof_index(level, d);
747  if (idx != DoFHandler<1>::invalid_dof_index)
748  cell->line(l)->set_mg_dof_index (level, d, ((indices.n_elements() == 0) ?
749  new_numbers[idx] :
750  new_numbers[indices.index_within_set(idx)]));
751  if (check_validity)
752  Assert(idx != DoFHandler<2>::invalid_dof_index, ExcInternalError ());
753  }
754  cell->line(l)->clear_user_flag();
755  }
756  // finally, restore user flags
757  const_cast<::Triangulation<2,spacedim> &>(dof_handler.get_triangulation()).load_user_flags (user_flags);
758  }
759 
760  for (std::vector<types::global_dof_index>::iterator i=dof_handler.mg_levels[level]->dof_object.dofs.begin();
761  i!=dof_handler.mg_levels[level]->dof_object.dofs.end(); ++i)
762  {
763  if (*i != DoFHandler<2>::invalid_dof_index)
764  {
765  Assert(*i<new_numbers.size(), ExcInternalError());
766  *i = ((indices.n_elements() == 0) ?
767  new_numbers[*i] :
768  new_numbers[indices.index_within_set(*i)]);
769  }
770  }
771 
772  }
773 
774 
775  template <int spacedim>
776  static
777  void
778  renumber_dofs (const std::vector<types::global_dof_index> &new_numbers,
779  const IndexSet &indices,
780  DoFHandler<3,spacedim> &dof_handler,
781  const bool check_validity)
782  {
783  // note that we can not use cell
784  // iterators in this function since
785  // then we would renumber the dofs on
786  // the interface of two cells more
787  // than once. Anyway, this way it's
788  // not only more correct but also
789  // faster; note, however, that dof
790  // numbers may be invalid_dof_index,
791  // namely when the appropriate
792  // vertex/line/etc is unused
793  for (std::vector<types::global_dof_index>::iterator
794  i=dof_handler.vertex_dofs.begin();
795  i!=dof_handler.vertex_dofs.end(); ++i)
796  if (*i != DoFHandler<3,spacedim>::invalid_dof_index)
797  *i = ((indices.n_elements() == 0) ?
798  new_numbers[*i] :
799  new_numbers[indices.index_within_set(*i)]);
800  else if (check_validity)
801  // if index is invalid_dof_index:
802  // check if this one really is
803  // unused
804  Assert (dof_handler.get_triangulation()
805  .vertex_used((i-dof_handler.vertex_dofs.begin()) /
806  dof_handler.selected_fe->dofs_per_vertex)
807  == false,
808  ExcInternalError ());
809 
810  for (std::vector<types::global_dof_index>::iterator
811  i=dof_handler.faces->lines.dofs.begin();
812  i!=dof_handler.faces->lines.dofs.end(); ++i)
813  if (*i != DoFHandler<3,spacedim>::invalid_dof_index)
814  *i = ((indices.n_elements() == 0) ?
815  new_numbers[*i] :
816  new_numbers[indices.index_within_set(*i)]);
817  for (std::vector<types::global_dof_index>::iterator
818  i=dof_handler.faces->quads.dofs.begin();
819  i!=dof_handler.faces->quads.dofs.end(); ++i)
820  if (*i != DoFHandler<3,spacedim>::invalid_dof_index)
821  *i = ((indices.n_elements() == 0) ?
822  new_numbers[*i] :
823  new_numbers[indices.index_within_set(*i)]);
824 
825  for (unsigned int level=0; level<dof_handler.levels.size(); ++level)
826  {
827  for (std::vector<types::global_dof_index>::iterator
828  i=dof_handler.levels[level]->dof_object.dofs.begin();
829  i!=dof_handler.levels[level]->dof_object.dofs.end(); ++i)
830  if (*i != DoFHandler<3,spacedim>::invalid_dof_index)
831  *i = ((indices.n_elements() == 0) ?
832  new_numbers[*i] :
833  new_numbers[indices.index_within_set(*i)]);
834  }
835 
836  // update the cache
837  // used for cell dof
838  // indices
839  for (typename DoFHandler<3,spacedim>::level_cell_iterator
840  cell = dof_handler.begin();
841  cell != dof_handler.end(); ++cell)
842  cell->update_cell_dof_indices_cache ();
843  }
844 
845  template <int spacedim>
846  static
847  void
848  renumber_mg_dofs (const std::vector<::types::global_dof_index> &new_numbers,
849  const IndexSet &indices,
850  DoFHandler<3,spacedim> &dof_handler,
851  const unsigned int level,
852  const bool check_validity)
853  {
854  if (level>=dof_handler.get_triangulation().n_levels())
855  return;
856  for (typename std::vector<typename DoFHandler<3,spacedim>::MGVertexDoFs>::iterator i=dof_handler.mg_vertex_dofs.begin();
857  i!=dof_handler.mg_vertex_dofs.end(); ++i)
858  // if the present vertex lives on
859  // the present level
860  if ((i->get_coarsest_level() <= level) &&
861  (i->get_finest_level() >= level))
862  for (unsigned int d=0; d<dof_handler.get_fe().dofs_per_vertex; ++d)
863  {
864  ::types::global_dof_index idx =i->get_index (level, d);
865  if (idx != DoFHandler<3>::invalid_dof_index)
866  i->set_index (level, d,
867  ((indices.n_elements() == 0) ?
868  new_numbers[idx] :
869  new_numbers[indices.index_within_set(idx)]));
870 
871  if (check_validity)
872  Assert(idx != DoFHandler<3>::invalid_dof_index, ExcInternalError ());
873  }
874 
875  if (dof_handler.get_fe().dofs_per_line > 0 ||
876  dof_handler.get_fe().dofs_per_quad > 0)
877  {
878  // save user flags as they will be modified
879  std::vector<bool> user_flags;
880  dof_handler.get_triangulation().save_user_flags(user_flags);
881  const_cast<::Triangulation<3,spacedim> &>(dof_handler.get_triangulation()).clear_user_flags ();
882 
883  // flag all lines adjacent to cells of the current level, as
884  // those lines logically belong to the same level as the cell,
885  // at least for isotropic refinement
886  for (typename DoFHandler<3,spacedim>::level_cell_iterator cell = dof_handler.begin(level);
887  cell != dof_handler.end(level); ++cell)
888  for (unsigned int line=0; line < GeometryInfo<3>::lines_per_cell; ++line)
889  cell->line(line)->set_user_flag();
890 
891  for (typename DoFHandler<3,spacedim>::cell_iterator cell = dof_handler.begin();
892  cell != dof_handler.end(); ++cell)
893  for (unsigned int l=0; l<GeometryInfo<3>::lines_per_cell; ++l)
894  if (cell->line(l)->user_flag_set())
895  {
896  for (unsigned int d=0; d<dof_handler.get_fe().dofs_per_line; ++d)
897  {
898  ::types::global_dof_index idx = cell->line(l)->mg_dof_index(level, d);
899  if (idx != DoFHandler<3>::invalid_dof_index)
900  cell->line(l)->set_mg_dof_index (level, d, ((indices.n_elements() == 0) ?
901  new_numbers[idx] :
902  new_numbers[indices.index_within_set(idx)]));
903  if (check_validity)
904  Assert(idx != DoFHandler<3>::invalid_dof_index, ExcInternalError ());
905  }
906  cell->line(l)->clear_user_flag();
907  }
908 
909  // flag all quads adjacent to cells of the current level, as
910  // those quads logically belong to the same level as the cell,
911  // at least for isotropic refinement
912  for (typename DoFHandler<3,spacedim>::level_cell_iterator cell = dof_handler.begin(level);
913  cell != dof_handler.end(level); ++cell)
914  for (unsigned int quad=0; quad < GeometryInfo<3>::quads_per_cell; ++quad)
915  cell->quad(quad)->set_user_flag();
916 
917  for (typename DoFHandler<3,spacedim>::cell_iterator cell = dof_handler.begin();
918  cell != dof_handler.end(); ++cell)
919  for (unsigned int l=0; l<GeometryInfo<3>::quads_per_cell; ++l)
920  if (cell->quad(l)->user_flag_set())
921  {
922  for (unsigned int d=0; d<dof_handler.get_fe().dofs_per_quad; ++d)
923  {
924  ::types::global_dof_index idx = cell->quad(l)->mg_dof_index(level, d);
925  if (idx != DoFHandler<1>::invalid_dof_index)
926  cell->quad(l)->set_mg_dof_index (level, d, ((indices.n_elements() == 0) ?
927  new_numbers[idx] :
928  new_numbers[indices.index_within_set(idx)]));
929  if (check_validity)
930  Assert(idx != DoFHandler<3>::invalid_dof_index, ExcInternalError ());
931  }
932  cell->quad(l)->clear_user_flag();
933  }
934 
935  // finally, restore user flags
936  const_cast<::Triangulation<3,spacedim> &>(dof_handler.get_triangulation()).load_user_flags (user_flags);
937  }
938 
939  for (std::vector<types::global_dof_index>::iterator i=dof_handler.mg_levels[level]->dof_object.dofs.begin();
940  i!=dof_handler.mg_levels[level]->dof_object.dofs.end(); ++i)
941  {
942  if (*i != DoFHandler<3>::invalid_dof_index)
943  {
944  Assert(*i<new_numbers.size(), ExcInternalError());
945  *i = ((indices.n_elements() == 0) ?
946  new_numbers[*i] :
947  new_numbers[indices.index_within_set(*i)]);
948  }
949  }
950  }
951 
952 
953  };
954 
955 
956 
957  /* --------------------- class PolicyBase ---------------- */
958 
959  template <int dim, int spacedim>
960  PolicyBase<dim,spacedim>::~PolicyBase ()
961  {}
962 
963 
964  /* --------------------- class Sequential ---------------- */
965 
966 
967  template <int dim, int spacedim>
968  void
969  Sequential<dim,spacedim>::
970  distribute_dofs (DoFHandler<dim,spacedim> &dof_handler,
971  NumberCache &number_cache_current ) const
972  {
973  const types::global_dof_index n_dofs =
974  Implementation::distribute_dofs (0,
975  numbers::invalid_subdomain_id,
976  dof_handler);
977 
978  // now set the elements of the
979  // number cache appropriately
980  NumberCache number_cache;
981  number_cache.n_global_dofs = n_dofs;
982  number_cache.n_locally_owned_dofs = number_cache.n_global_dofs;
983 
984  number_cache.locally_owned_dofs
985  = IndexSet (number_cache.n_global_dofs);
986  number_cache.locally_owned_dofs.add_range (0,
987  number_cache.n_global_dofs);
988  number_cache.locally_owned_dofs.compress();
989 
990  number_cache.n_locally_owned_dofs_per_processor
991  = std::vector<types::global_dof_index> (1,
992  number_cache.n_global_dofs);
993 
994  number_cache.locally_owned_dofs_per_processor
995  = std::vector<IndexSet> (1,
996  number_cache.locally_owned_dofs);
997  number_cache_current = number_cache;
998  }
999 
1000 
1001  template <int dim, int spacedim>
1002  void
1003  Sequential<dim,spacedim>::
1004  distribute_mg_dofs (DoFHandler<dim,spacedim> &dof_handler,
1005  std::vector<NumberCache> &number_caches) const
1006  {
1007  std::vector<bool> user_flags;
1008 
1009  dof_handler.get_triangulation().save_user_flags (user_flags);
1010  const_cast<::Triangulation<dim, spacedim>&>(dof_handler.get_triangulation()).clear_user_flags ();
1011 
1012  for (unsigned int level = 0; level < dof_handler.get_triangulation().n_levels(); ++level)
1013  {
1014  types::global_dof_index next_free_dof = Implementation::distribute_dofs_on_level(0, numbers::invalid_subdomain_id, dof_handler, level);
1015 
1016  number_caches[level].n_global_dofs = next_free_dof;
1017  number_caches[level].n_locally_owned_dofs = next_free_dof;
1018  number_caches[level].locally_owned_dofs = complete_index_set(next_free_dof);
1019  number_caches[level].locally_owned_dofs_per_processor.resize(1);
1020  number_caches[level].locally_owned_dofs_per_processor[0] = complete_index_set(next_free_dof);
1021  number_caches[level].n_locally_owned_dofs_per_processor.resize(1);
1022  number_caches[level].n_locally_owned_dofs_per_processor[0] = next_free_dof;
1023  }
1024  const_cast<::Triangulation<dim, spacedim>&>(dof_handler.get_triangulation()).load_user_flags (user_flags);
1025  }
1026 
1027  template <int dim, int spacedim>
1028  void
1029  Sequential<dim,spacedim>::
1030  renumber_dofs (const std::vector<types::global_dof_index> &new_numbers,
1031  ::DoFHandler<dim,spacedim> &dof_handler,
1032  NumberCache &number_cache_current) const
1033  {
1034  Implementation::renumber_dofs (new_numbers, IndexSet(0),
1035  dof_handler, true);
1036 
1037  // in the sequential case,
1038  // the number cache should
1039  // not have changed but we
1040  // have to set the elements
1041  // of the structure
1042  // appropriately anyway
1043  NumberCache number_cache;
1044  number_cache.n_global_dofs = dof_handler.n_dofs();
1045  number_cache.n_locally_owned_dofs = number_cache.n_global_dofs;
1046 
1047  number_cache.locally_owned_dofs
1048  = IndexSet (number_cache.n_global_dofs);
1049  number_cache.locally_owned_dofs.add_range (0,
1050  number_cache.n_global_dofs);
1051  number_cache.locally_owned_dofs.compress();
1052 
1053  number_cache.n_locally_owned_dofs_per_processor
1054  = std::vector<types::global_dof_index> (1,
1055  number_cache.n_global_dofs);
1056 
1057  number_cache.locally_owned_dofs_per_processor
1058  = std::vector<IndexSet> (1,
1059  number_cache.locally_owned_dofs);
1060  number_cache_current = number_cache;
1061  }
1062 
1063  /* --------------------- class ParallelShared ---------------- */
1064 
1065  template <int dim, int spacedim>
1066  void
1067  ParallelShared<dim,spacedim>::
1068  distribute_dofs (DoFHandler<dim,spacedim> &dof_handler,
1069  NumberCache &number_cache) const
1070  {
1071  // If the underlying shared::Tria allows artifical cells, we need to do
1072  // some tricks here to make Sequential algorithms play nicely.
1073  // Namely, we first restore original partition (without artificial cells)
1074  // and then turn artificial cells on at the end of this function.
1075  const parallel::shared::Triangulation<dim, spacedim> *tr =
1076  (dynamic_cast<const parallel::shared::Triangulation<dim, spacedim>*> (&dof_handler.get_triangulation()));
1077  Assert(tr != 0, ExcInternalError());
1078  typename parallel::shared::Triangulation<dim,spacedim>::active_cell_iterator
1079  cell = dof_handler.get_triangulation().begin_active(),
1080  endc = dof_handler.get_triangulation().end();
1081  std::vector<types::subdomain_id> current_subdomain_ids(tr->n_active_cells());
1082  const std::vector<types::subdomain_id> &true_subdomain_ids = tr->get_true_subdomain_ids_of_cells();
1083  if (tr->with_artificial_cells())
1084  for (unsigned int index=0; cell != endc; cell++, index++)
1085  {
1086  current_subdomain_ids[index] = cell->subdomain_id();
1087  cell->set_subdomain_id(true_subdomain_ids[index]);
1088  }
1089 
1090  Sequential<dim,spacedim>::distribute_dofs (dof_handler,number_cache);
1091  DoFRenumbering::subdomain_wise (dof_handler);
1092  // dofrenumbering will reset subdomains, this is ugly but we need to do it again:
1093  cell = tr->begin_active();
1094  if (tr->with_artificial_cells())
1095  for (unsigned int index=0; cell != endc; cell++, index++)
1096  cell->set_subdomain_id(true_subdomain_ids[index]);
1097 
1098  number_cache.locally_owned_dofs_per_processor = DoFTools::locally_owned_dofs_per_subdomain (dof_handler);
1099  number_cache.locally_owned_dofs = number_cache.locally_owned_dofs_per_processor[dof_handler.get_triangulation().locally_owned_subdomain()];
1100  number_cache.n_locally_owned_dofs_per_processor.resize (number_cache.locally_owned_dofs_per_processor.size());
1101  for (unsigned int i = 0; i < number_cache.n_locally_owned_dofs_per_processor.size(); i++)
1102  number_cache.n_locally_owned_dofs_per_processor[i] = number_cache.locally_owned_dofs_per_processor[i].n_elements();
1103  number_cache.n_locally_owned_dofs = number_cache.n_locally_owned_dofs_per_processor[dof_handler.get_triangulation().locally_owned_subdomain()];
1104 
1105  // restore current subdomain ids
1106  cell = tr->begin_active();
1107  if (tr->with_artificial_cells())
1108  for (unsigned int index=0; cell != endc; cell++, index++)
1109  cell->set_subdomain_id(current_subdomain_ids[index]);
1110  }
1111 
1112  template <int dim, int spacedim>
1113  void
1114  ParallelShared<dim,spacedim>::
1115  distribute_mg_dofs (DoFHandler<dim,spacedim> &dof_handler,
1116  std::vector<NumberCache> &number_caches) const
1117  {
1118  // first, call the sequential function to distribute dofs
1119  Sequential<dim,spacedim>:: distribute_mg_dofs (dof_handler, number_caches);
1120  // now we need to update the number cache.
1121  // This part is not yet implemented.
1122  AssertThrow(false,ExcNotImplemented());
1123  }
1124 
1125  template <int dim, int spacedim>
1126  void
1127  ParallelShared<dim,spacedim>::
1128  renumber_dofs (const std::vector<types::global_dof_index> &new_numbers,
1129  ::DoFHandler<dim,spacedim> &dof_handler,
1130  NumberCache &number_cache) const
1131  {
1132 
1133 #ifndef DEAL_II_WITH_MPI
1134  (void)new_numbers;
1135  (void)dof_handler;
1136  (void)number_cache;
1137  Assert (false, ExcNotImplemented());
1138 #else
1139  // Similar to distribute_dofs() we need to have a special treatment in
1140  // case artificial cells are present.
1141  const parallel::shared::Triangulation<dim, spacedim> *tr =
1142  (dynamic_cast<const parallel::shared::Triangulation<dim, spacedim>*> (&dof_handler.get_triangulation()));
1143  Assert(tr != 0, ExcInternalError());
1144  typename parallel::shared::Triangulation<dim,spacedim>::active_cell_iterator
1145  cell = dof_handler.get_triangulation().begin_active(),
1146  endc = dof_handler.get_triangulation().end();
1147  std::vector<types::subdomain_id> current_subdomain_ids(tr->n_active_cells());
1148  const std::vector<types::subdomain_id> &true_subdomain_ids = tr->get_true_subdomain_ids_of_cells();
1149  if (tr->with_artificial_cells())
1150  for (unsigned int index=0; cell != endc; cell++, index++)
1151  {
1152  current_subdomain_ids[index] = cell->subdomain_id();
1153  cell->set_subdomain_id(true_subdomain_ids[index]);
1154  }
1155 
1156  std::vector<types::global_dof_index> global_gathered_numbers (dof_handler.n_dofs (), 0);
1157  // as we call DoFRenumbering::subdomain_wise (dof_handler) from distribute_dofs(),
1158  // we need to support sequential-like input.
1159  // Distributed-like input from, for example, component_wise renumbering is also supported.
1160  if (new_numbers.size () == dof_handler.n_dofs ())
1161  {
1162  global_gathered_numbers = new_numbers;
1163  }
1164  else
1165  {
1166  Assert(new_numbers.size() == dof_handler.locally_owned_dofs().n_elements(),
1167  ExcInternalError());
1168  const unsigned int n_cpu = Utilities::MPI::n_mpi_processes (tr->get_communicator ());
1169  std::vector<types::global_dof_index> gathered_new_numbers (dof_handler.n_dofs (), 0);
1170  Assert(Utilities::MPI::this_mpi_process (tr->get_communicator ()) ==
1171  dof_handler.get_triangulation().locally_owned_subdomain (),
1172  ExcInternalError())
1173 
1174  //gather new numbers among processors into one vector
1175  {
1176  std::vector<types::global_dof_index> new_numbers_copy (new_numbers);
1177  // displs:
1178  // Entry i specifies the displacement (relative to recvbuf )
1179  // at which to place the incoming data from process i
1180  // rcounts:
1181  // containing the number of elements that are to be received from each process
1182  std::vector<int> displs(n_cpu),
1183  rcounts(n_cpu);
1184  types::global_dof_index shift = 0;
1185  //set rcounts based on new_numbers:
1186  int cur_count = new_numbers_copy.size ();
1187  MPI_Allgather (&cur_count, 1, MPI_INT,
1188  &rcounts[0], 1, MPI_INT,
1189  tr->get_communicator ());
1190 
1191  for (unsigned int i = 0; i < n_cpu; i++)
1192  {
1193  displs[i] = shift;
1194  shift += rcounts[i];
1195  }
1196  Assert(((int)new_numbers_copy.size()) ==
1197  rcounts[Utilities::MPI::this_mpi_process (tr->get_communicator ())],
1198  ExcInternalError());
1199  MPI_Allgatherv (&new_numbers_copy[0], new_numbers_copy.size (),
1200  DEAL_II_DOF_INDEX_MPI_TYPE,
1201  &gathered_new_numbers[0], &rcounts[0],
1202  &displs[0],
1203  DEAL_II_DOF_INDEX_MPI_TYPE,
1204  tr->get_communicator ());
1205  }
1206 
1207  // put new numbers according to the current locally_owned_dofs_per_processor IndexSets
1208  types::global_dof_index shift = 0;
1209  // flag_1 and flag_2 are
1210  // used to control that there is a
1211  // one-to-one relation between old and new DoFs.
1212  std::vector<unsigned int> flag_1 (dof_handler.n_dofs (), 0),
1213  flag_2 (dof_handler.n_dofs (), 0);
1214  for (unsigned int i = 0; i < n_cpu; i++)
1215  {
1216  const IndexSet &iset =
1217  number_cache.locally_owned_dofs_per_processor[i];
1218  for (types::global_dof_index ind = 0;
1219  ind < iset.n_elements (); ind++)
1220  {
1221  const types::global_dof_index target = iset.nth_index_in_set (ind);
1222  const types::global_dof_index value = gathered_new_numbers[shift + ind];
1223  Assert(target < dof_handler.n_dofs(), ExcInternalError());
1224  Assert(value < dof_handler.n_dofs(), ExcInternalError());
1225  global_gathered_numbers[target] = value;
1226  flag_1[target]++;
1227  flag_2[value]++;
1228  }
1229  shift += iset.n_elements ();
1230  }
1231 
1232  Assert(*std::max_element(flag_1.begin(), flag_1.end()) == 1,
1233  ExcInternalError());
1234  Assert(*std::min_element(flag_1.begin(), flag_1.end()) == 1,
1235  ExcInternalError());
1236  Assert((*std::max_element(flag_2.begin(), flag_2.end())) == 1,
1237  ExcInternalError());
1238  Assert((*std::min_element(flag_2.begin(), flag_2.end())) == 1,
1239  ExcInternalError());
1240  }
1241  Sequential<dim, spacedim>::renumber_dofs (global_gathered_numbers, dof_handler, number_cache);
1242  // correct number_cache:
1243  number_cache.locally_owned_dofs_per_processor =
1244  DoFTools::locally_owned_dofs_per_subdomain (dof_handler);
1245  number_cache.locally_owned_dofs =
1246  number_cache.locally_owned_dofs_per_processor[dof_handler.get_triangulation().locally_owned_subdomain ()];
1247  // sequential renumbering returns a vector of size 1 here,
1248  // correct this:
1249  number_cache.n_locally_owned_dofs_per_processor.resize(number_cache.locally_owned_dofs_per_processor.size());
1250  for (unsigned int i = 0;
1251  i < number_cache.n_locally_owned_dofs_per_processor.size (); i++)
1252  number_cache.n_locally_owned_dofs_per_processor[i] = number_cache.locally_owned_dofs_per_processor[i].n_elements ();
1253 
1254  number_cache.n_locally_owned_dofs =
1255  number_cache.n_locally_owned_dofs_per_processor[dof_handler.get_triangulation().locally_owned_subdomain ()];
1256 
1257  // restore artificial cells
1258  cell = tr->begin_active();
1259  if (tr->with_artificial_cells())
1260  for (unsigned int index=0; cell != endc; cell++, index++)
1261  cell->set_subdomain_id(current_subdomain_ids[index]);
1262 #endif
1263  }
1264 
1265  /* --------------------- class ParallelDistributed ---------------- */
1266 
1267 #ifdef DEAL_II_WITH_P4EST
1268 
1269  namespace
1270  {
1271  template <int dim>
1272  struct types
1273  {
1274 
1284  struct cellinfo
1285  {
1286  std::vector<unsigned int> tree_index;
1287  std::vector<typename ::internal::p4est::types<dim>::quadrant> quadrants;
1288  std::vector<::types::global_dof_index> dofs;
1289 
1290  unsigned int bytes_for_buffer () const
1291  {
1292  return (sizeof(unsigned int) +
1293  tree_index.size() * sizeof(unsigned int) +
1294  quadrants.size() * sizeof(typename ::internal::p4est
1295  ::types<dim>::quadrant) +
1296  dofs.size() * sizeof(::types::global_dof_index));
1297  }
1298 
1299  void pack_data (std::vector<char> &buffer) const
1300  {
1301  buffer.resize(bytes_for_buffer());
1302 
1303  char *ptr = &buffer[0];
1304 
1305  const unsigned int num_cells = tree_index.size();
1306  std::memcpy(ptr, &num_cells, sizeof(unsigned int));
1307  ptr += sizeof(unsigned int);
1308 
1309  std::memcpy(ptr,
1310  &tree_index[0],
1311  num_cells*sizeof(unsigned int));
1312  ptr += num_cells*sizeof(unsigned int);
1313 
1314  std::memcpy(ptr,
1315  &quadrants[0],
1316  num_cells * sizeof(typename ::internal::p4est::
1317  types<dim>::quadrant));
1318  ptr += num_cells*sizeof(typename ::internal::p4est::types<dim>::
1319  quadrant);
1320 
1321  std::memcpy(ptr,
1322  &dofs[0],
1323  dofs.size() * sizeof(::types::global_dof_index));
1324  ptr += dofs.size() * sizeof(::types::global_dof_index);
1325 
1326  Assert (ptr == &buffer[0]+buffer.size(),
1327  ExcInternalError());
1328 
1329  }
1330  };
1331  };
1332 
1333 
1334 
1335  template <int dim, int spacedim>
1336  void
1337  fill_dofindices_recursively (const typename parallel::distributed::Triangulation<dim,spacedim> &tria,
1338  const unsigned int tree_index,
1339  const typename DoFHandler<dim,spacedim>::level_cell_iterator &dealii_cell,
1340  const typename ::internal::p4est::types<dim>::quadrant &p4est_cell,
1341  const std::map<unsigned int, std::set<::types::subdomain_id> > &vertices_with_ghost_neighbors,
1342  std::map<::types::subdomain_id, typename types<dim>::cellinfo> &needs_to_get_cell)
1343  {
1344  // see if we have to
1345  // recurse...
1346  if (dealii_cell->has_children())
1347  {
1348  typename ::internal::p4est::types<dim>::quadrant
1349  p4est_child[GeometryInfo<dim>::max_children_per_cell];
1350  internal::p4est::init_quadrant_children<dim>(p4est_cell, p4est_child);
1351 
1352 
1353  for (unsigned int c=0; c<GeometryInfo<dim>::max_children_per_cell; ++c)
1354  fill_dofindices_recursively<dim,spacedim>(tria,
1355  tree_index,
1356  dealii_cell->child(c),
1357  p4est_child[c],
1358  vertices_with_ghost_neighbors,
1359  needs_to_get_cell);
1360  return;
1361  }
1362 
1363  // we're at a leaf cell. see if
1364  // the cell is flagged as
1365  // interesting. note that we
1366  // have only flagged our own
1367  // cells before
1368  if (dealii_cell->user_flag_set() && !dealii_cell->is_ghost())
1369  {
1370  Assert (!dealii_cell->is_artificial(), ExcInternalError());
1371 
1372  // check each vertex if
1373  // it is interesting and
1374  // push dofindices if yes
1375  std::set<::types::subdomain_id> send_to;
1376  for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell; ++v)
1377  {
1378  const std::map<unsigned int, std::set<::types::subdomain_id> >::const_iterator
1379  neighbor_subdomains_of_vertex
1380  = vertices_with_ghost_neighbors.find (dealii_cell->vertex_index(v));
1381 
1382  if (neighbor_subdomains_of_vertex ==
1383  vertices_with_ghost_neighbors.end())
1384  continue;
1385 
1386  Assert(neighbor_subdomains_of_vertex->second.size()!=0,
1387  ExcInternalError());
1388 
1389  send_to.insert(neighbor_subdomains_of_vertex->second.begin(),
1390  neighbor_subdomains_of_vertex->second.end());
1391  }
1392 
1393  if (send_to.size() > 0)
1394  {
1395  // this cell's dof_indices
1396  // need to be sent to
1397  // someone
1398  std::vector<::types::global_dof_index>
1399  local_dof_indices (dealii_cell->get_fe().dofs_per_cell);
1400  dealii_cell->get_dof_indices (local_dof_indices);
1401 
1402  for (std::set<::types::subdomain_id>::iterator it=send_to.begin();
1403  it!=send_to.end(); ++it)
1404  {
1405  const ::types::subdomain_id subdomain = *it;
1406 
1407  // get an iterator
1408  // to what needs to
1409  // be sent to that
1410  // subdomain (if
1411  // already exists),
1412  // or create such
1413  // an object
1414  typename std::map<::types::subdomain_id, typename types<dim>::cellinfo>::iterator
1415  p
1416  = needs_to_get_cell.insert (std::make_pair(subdomain,
1417  typename types<dim>::cellinfo()))
1418  .first;
1419 
1420  p->second.tree_index.push_back(tree_index);
1421  p->second.quadrants.push_back(p4est_cell);
1422 
1423  p->second.dofs.push_back(dealii_cell->get_fe().dofs_per_cell);
1424  p->second.dofs.insert(p->second.dofs.end(),
1425  local_dof_indices.begin(),
1426  local_dof_indices.end());
1427 
1428  }
1429  }
1430  }
1431  }
1432 
1433  template <int dim, int spacedim>
1434  void
1435  fill_mg_dofindices_recursively (const typename parallel::distributed::Triangulation<dim,spacedim> &tria,
1436  const unsigned int tree_index,
1437  const typename DoFHandler<dim,spacedim>::level_cell_iterator &dealii_cell,
1438  const typename ::internal::p4est::types<dim>::quadrant &p4est_cell,
1439  const std::map<unsigned int, std::set<::types::subdomain_id> > &vertices_with_ghost_neighbors,
1440  std::map<::types::subdomain_id, typename types<dim>::cellinfo> &needs_to_get_cell,
1441  const unsigned int level)
1442  {
1443  if (dealii_cell->level()>(int)level)
1444  return;
1445  // see if we have to
1446  // recurse...
1447  if (dealii_cell->has_children())
1448  {
1449  typename ::internal::p4est::types<dim>::quadrant
1450  p4est_child[GeometryInfo<dim>::max_children_per_cell];
1451  internal::p4est::init_quadrant_children<dim>(p4est_cell, p4est_child);
1452 
1453 
1454  for (unsigned int c=0; c<GeometryInfo<dim>::max_children_per_cell; ++c)
1455  fill_mg_dofindices_recursively<dim,spacedim>(tria,
1456  tree_index,
1457  dealii_cell->child(c),
1458  p4est_child[c],
1459  vertices_with_ghost_neighbors,
1460  needs_to_get_cell,
1461  level);
1462  }
1463 
1464  if (dealii_cell->level()<(int)level)
1465  return;
1466 
1467  // now we are on the right level!
1468  Assert(dealii_cell->level()==(int)level, ExcInternalError());
1469 
1470  // see if
1471  // the cell is flagged as
1472  // interesting. note that we
1473  // have only flagged our own
1474  // cells before
1475  if (dealii_cell->user_flag_set() && dealii_cell->level_subdomain_id() == tria.locally_owned_subdomain())
1476  {
1477  // check each vertex if
1478  // it is interesting and
1479  // push dofindices if yes
1480  std::set<::types::subdomain_id> send_to;
1481  for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell; ++v)
1482  {
1483  const std::map<unsigned int, std::set<::types::subdomain_id> >::const_iterator
1484  neighbor_subdomains_of_vertex
1485  = vertices_with_ghost_neighbors.find (dealii_cell->vertex_index(v));
1486 
1487  if (neighbor_subdomains_of_vertex ==
1488  vertices_with_ghost_neighbors.end())
1489  continue;
1490 
1491  Assert(neighbor_subdomains_of_vertex->second.size()!=0,
1492  ExcInternalError());
1493 
1494  send_to.insert(neighbor_subdomains_of_vertex->second.begin(),
1495  neighbor_subdomains_of_vertex->second.end());
1496  }
1497 
1498  // additionally, if we need to send to all our direct children (multigrid only)
1499  if (dealii_cell->has_children())
1500  {
1501  for (unsigned int c=0; c<GeometryInfo<dim>::max_children_per_cell; ++c)
1502  {
1503  //TODO: we don't know about our children if proc 0 owns all coarse cells!
1504  ::types::subdomain_id dest = dealii_cell->child(c)->level_subdomain_id();
1505  Assert(dest!=::numbers::artificial_subdomain_id && dest!=::numbers::invalid_subdomain_id, ExcInternalError());
1506  if (dest != tria.locally_owned_subdomain())
1507  send_to.insert(dest);
1508  }
1509  }
1510 
1511  //additionally (multigrid only), we can have the case that children of our neighbor
1512  //have us as a neighbor. In this case we and the children are active.
1513  if (dealii_cell->active())
1514  {
1515  for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
1516  {
1517  if (dealii_cell->at_boundary(f))
1518  continue;
1519  typename DoFHandler<dim,spacedim>::level_cell_iterator neighbor = dealii_cell->neighbor(f);
1520  if (!neighbor->has_children())
1521  continue;
1522 
1523  for (unsigned int subface=0; subface<GeometryInfo<dim>::max_children_per_face; ++subface)
1524  {
1525  typename DoFHandler<dim,spacedim>::level_cell_iterator child = dealii_cell->neighbor_child_on_subface(f,subface);
1526  ::types::subdomain_id dest = child->subdomain_id();
1527  Assert(dest != ::numbers::artificial_subdomain_id, ExcInternalError());
1528  if (dest != tria.locally_owned_subdomain())
1529  send_to.insert(dest);
1530  }
1531 
1532  }
1533 
1534  }
1535 
1536  // Finally, if we are neighboring a coarser cell, add them to
1537  // the destination list
1538  if (dealii_cell->active())
1539  {
1540  for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
1541  {
1542  if (dealii_cell->at_boundary(f))
1543  continue;
1544  typename DoFHandler<dim,spacedim>::level_cell_iterator neighbor = dealii_cell->neighbor(f);
1545  if (neighbor->level()>=dealii_cell->level())
1546  continue;
1547 
1548  ::types::subdomain_id dest = neighbor->level_subdomain_id();
1549  Assert(dest != ::numbers::artificial_subdomain_id, ExcInternalError());
1550  if (dest != tria.locally_owned_subdomain())
1551  send_to.insert(dest);
1552 
1553  }
1554  }
1555 
1556 
1557  // send if we have something to send
1558  if (send_to.size() > 0)
1559  {
1560  // this cell's dof_indices
1561  // need to be sent to
1562  // someone
1563  std::vector<::types::global_dof_index>
1564  local_dof_indices (dealii_cell->get_fe().dofs_per_cell);
1565  dealii_cell->get_mg_dof_indices (local_dof_indices);
1566 
1567  for (std::set<::types::subdomain_id>::iterator it=send_to.begin();
1568  it!=send_to.end(); ++it)
1569  {
1570  const ::types::subdomain_id subdomain = *it;
1571 
1572  // get an iterator
1573  // to what needs to
1574  // be sent to that
1575  // subdomain (if
1576  // already exists),
1577  // or create such
1578  // an object
1579  typename std::map<::types::subdomain_id, typename types<dim>::cellinfo>::iterator
1580  p
1581  = needs_to_get_cell.insert (std::make_pair(subdomain,
1582  typename types<dim>::cellinfo()))
1583  .first;
1584 
1585  p->second.tree_index.push_back(tree_index);
1586  p->second.quadrants.push_back(p4est_cell);
1587 
1588  p->second.dofs.push_back(dealii_cell->get_fe().dofs_per_cell);
1589  p->second.dofs.insert(p->second.dofs.end(),
1590  local_dof_indices.begin(),
1591  local_dof_indices.end());
1592 
1593  }
1594  }
1595  }
1596  }
1597 
1598 
1599  template <int dim, int spacedim>
1600  void
1601  set_dofindices_recursively (
1602  const parallel::distributed::Triangulation<dim,spacedim> &tria,
1603  const typename ::internal::p4est::types<dim>::quadrant &p4est_cell,
1604  const typename DoFHandler<dim,spacedim>::level_cell_iterator &dealii_cell,
1605  const typename ::internal::p4est::types<dim>::quadrant &quadrant,
1606  ::types::global_dof_index *dofs)
1607  {
1608  if (internal::p4est::quadrant_is_equal<dim>(p4est_cell, quadrant))
1609  {
1610  Assert(!dealii_cell->has_children(), ExcInternalError());
1611  Assert(dealii_cell->is_ghost(), ExcInternalError());
1612 
1613  // update dof indices of cell
1614  std::vector<::types::global_dof_index>
1615  dof_indices (dealii_cell->get_fe().dofs_per_cell);
1616  dealii_cell->update_cell_dof_indices_cache();
1617  dealii_cell->get_dof_indices(dof_indices);
1618 
1619  bool complete = true;
1620  for (unsigned int i=0; i<dof_indices.size(); ++i)
1621  if (dofs[i] != DoFHandler<dim,spacedim>::invalid_dof_index)
1622  {
1623  Assert((dof_indices[i] ==
1624  (DoFHandler<dim,spacedim>::invalid_dof_index))
1625  ||
1626  (dof_indices[i]==dofs[i]),
1627  ExcInternalError());
1628  dof_indices[i]=dofs[i];
1629  }
1630  else
1631  complete=false;
1632 
1633  if (!complete)
1634  const_cast
1635  <typename DoFHandler<dim,spacedim>::level_cell_iterator &>
1636  (dealii_cell)->set_user_flag();
1637  else
1638  const_cast
1639  <typename DoFHandler<dim,spacedim>::level_cell_iterator &>
1640  (dealii_cell)->clear_user_flag();
1641 
1642  const_cast
1643  <typename DoFHandler<dim,spacedim>::level_cell_iterator &>
1644  (dealii_cell)->set_dof_indices(dof_indices);
1645 
1646  return;
1647  }
1648 
1649  if (! dealii_cell->has_children())
1650  return;
1651 
1652  if (! internal::p4est::quadrant_is_ancestor<dim> (p4est_cell, quadrant))
1653  return;
1654 
1655  typename ::internal::p4est::types<dim>::quadrant
1656  p4est_child[GeometryInfo<dim>::max_children_per_cell];
1657  internal::p4est::init_quadrant_children<dim>(p4est_cell, p4est_child);
1658 
1659  for (unsigned int c=0; c<GeometryInfo<dim>::max_children_per_cell; ++c)
1660  set_dofindices_recursively<dim,spacedim> (tria, p4est_child[c],
1661  dealii_cell->child(c),
1662  quadrant, dofs);
1663  }
1664 
1665 
1666  template <int dim, int spacedim>
1667  void
1668  set_mg_dofindices_recursively (
1669  const parallel::distributed::Triangulation<dim,spacedim> &tria,
1670  const typename ::internal::p4est::types<dim>::quadrant &p4est_cell,
1671  const typename DoFHandler<dim,spacedim>::level_cell_iterator &dealii_cell,
1672  const typename ::internal::p4est::types<dim>::quadrant &quadrant,
1673  ::types::global_dof_index *dofs,
1674  unsigned int level)
1675  {
1676  if (internal::p4est::quadrant_is_equal<dim>(p4est_cell, quadrant))
1677  {
1678  Assert(dealii_cell->level_subdomain_id()!=::numbers::artificial_subdomain_id, ExcInternalError());
1679  Assert(dealii_cell->level()==(int)level, ExcInternalError());
1680 
1681  // update dof indices of cell
1682  std::vector<::types::global_dof_index>
1683  dof_indices (dealii_cell->get_fe().dofs_per_cell);
1684  dealii_cell->get_mg_dof_indices(dof_indices);
1685 
1686  bool complete = true;
1687  for (unsigned int i=0; i<dof_indices.size(); ++i)
1688  if (dofs[i] != DoFHandler<dim,spacedim>::invalid_dof_index)
1689  {
1690  Assert((dof_indices[i] ==
1691  (DoFHandler<dim,spacedim>::invalid_dof_index))
1692  ||
1693  (dof_indices[i]==dofs[i]),
1694  ExcInternalError());
1695  dof_indices[i]=dofs[i];
1696  }
1697  else
1698  complete=false;
1699 
1700  if (!complete)
1701  const_cast
1702  <typename DoFHandler<dim,spacedim>::level_cell_iterator &>
1703  (dealii_cell)->set_user_flag();
1704  else
1705  const_cast
1706  <typename DoFHandler<dim,spacedim>::level_cell_iterator &>
1707  (dealii_cell)->clear_user_flag();
1708 
1709  const_cast
1710  <typename DoFHandler<dim,spacedim>::level_cell_iterator &>
1711  (dealii_cell)->set_mg_dof_indices(dof_indices);
1712  return;
1713  }
1714 
1715  if (! dealii_cell->has_children())
1716  return;
1717 
1718  if (! internal::p4est::quadrant_is_ancestor<dim> (p4est_cell, quadrant))
1719  return;
1720 
1721  typename ::internal::p4est::types<dim>::quadrant
1722  p4est_child[GeometryInfo<dim>::max_children_per_cell];
1723  internal::p4est::init_quadrant_children<dim>(p4est_cell, p4est_child);
1724 
1725  for (unsigned int c=0; c<GeometryInfo<dim>::max_children_per_cell; ++c)
1726  set_mg_dofindices_recursively<dim,spacedim> (tria, p4est_child[c],
1727  dealii_cell->child(c),
1728  quadrant, dofs, level);
1729 
1730  }
1731 
1732 
1733 
1734  template <int spacedim>
1735  void
1736  communicate_dof_indices_on_marked_cells
1737  (const DoFHandler<1,spacedim> &,
1738  const std::map<unsigned int, std::set<::types::subdomain_id> > &,
1739  const std::vector<::types::global_dof_index> &,
1740  const std::vector<::types::global_dof_index> &)
1741  {
1742  Assert (false, ExcNotImplemented());
1743  }
1744 
1745 
1746 
1747  template <int dim, int spacedim>
1748  void
1749  communicate_dof_indices_on_marked_cells
1750  (const DoFHandler<dim,spacedim> &dof_handler,
1751  const std::map<unsigned int, std::set<::types::subdomain_id> > &vertices_with_ghost_neighbors,
1752  const std::vector<::types::global_dof_index> &coarse_cell_to_p4est_tree_permutation,
1753  const std::vector<::types::global_dof_index> &p4est_tree_to_coarse_cell_permutation)
1754  {
1755 #ifndef DEAL_II_WITH_P4EST
1756  (void)vertices_with_ghost_neighbors;
1757  Assert (false, ExcNotImplemented());
1758 #else
1759 
1760  const parallel::distributed::Triangulation< dim, spacedim > *tr
1761  = (dynamic_cast<const parallel::distributed::Triangulation<dim,spacedim>*>
1762  (&dof_handler.get_triangulation()));
1763  Assert (tr != 0, ExcInternalError());
1764 
1765  // now collect cells and their
1766  // dof_indices for the
1767  // interested neighbors
1768  typedef
1769  std::map<::types::subdomain_id, typename types<dim>::cellinfo>
1770  cellmap_t;
1771  cellmap_t needs_to_get_cells;
1772 
1773  for (typename DoFHandler<dim,spacedim>::level_cell_iterator
1774  cell = dof_handler.begin(0);
1775  cell != dof_handler.end(0);
1776  ++cell)
1777  {
1778  typename ::internal::p4est::types<dim>::quadrant p4est_coarse_cell;
1779  internal::p4est::init_coarse_quadrant<dim>(p4est_coarse_cell);
1780 
1781  fill_dofindices_recursively<dim,spacedim>
1782  (*tr,
1783  coarse_cell_to_p4est_tree_permutation[cell->index()],
1784  cell,
1785  p4est_coarse_cell,
1786  vertices_with_ghost_neighbors,
1787  needs_to_get_cells);
1788  }
1789 
1790 
1791  //sending
1792  std::vector<std::vector<char> > sendbuffers (needs_to_get_cells.size());
1793  std::vector<std::vector<char> >::iterator buffer = sendbuffers.begin();
1794  std::vector<MPI_Request> requests (needs_to_get_cells.size());
1795 
1796  unsigned int idx=0;
1797 
1798  for (typename cellmap_t::iterator it=needs_to_get_cells.begin();
1799  it!=needs_to_get_cells.end();
1800  ++it, ++buffer, ++idx)
1801  {
1802  const unsigned int num_cells = it->second.tree_index.size();
1803  (void)num_cells;
1804 
1805  Assert(num_cells==it->second.quadrants.size(), ExcInternalError());
1806  Assert(num_cells>0, ExcInternalError());
1807 
1808  // pack all the data into
1809  // the buffer for this
1810  // recipient and send
1811  // it. keep data around
1812  // till we can make sure
1813  // that the packet has been
1814  // received
1815  it->second.pack_data (*buffer);
1816  MPI_Isend(&(*buffer)[0], buffer->size(),
1817  MPI_BYTE, it->first,
1818  123, tr->get_communicator(), &requests[idx]);
1819  }
1820 
1821 
1822  // mark all own cells, that miss some
1823  // dof_data and collect the neighbors
1824  // that are going to send stuff to us
1825  std::set<::types::subdomain_id> senders;
1826  {
1827  std::vector<::types::global_dof_index> local_dof_indices;
1828  typename DoFHandler<dim,spacedim>::active_cell_iterator
1829  cell, endc = dof_handler.end();
1830 
1831  for (cell = dof_handler.begin_active(); cell != endc; ++cell)
1832  if (!cell->is_artificial())
1833  {
1834  if (cell->is_ghost())
1835  {
1836  if (cell->user_flag_set())
1837  senders.insert(cell->subdomain_id());
1838  }
1839  else
1840  {
1841  local_dof_indices.resize (cell->get_fe().dofs_per_cell);
1842  cell->get_dof_indices (local_dof_indices);
1843  if (local_dof_indices.end() !=
1844  std::find (local_dof_indices.begin(),
1845  local_dof_indices.end(),
1846  DoFHandler<dim,spacedim>::invalid_dof_index))
1847  cell->set_user_flag();
1848  else
1849  cell->clear_user_flag();
1850  }
1851 
1852  }
1853  }
1854 
1855 
1856  //* 5. receive ghostcelldata
1857  std::vector<char> receive;
1858  typename types<dim>::cellinfo cellinfo;
1859  for (unsigned int i=0; i<senders.size(); ++i)
1860  {
1861  MPI_Status status;
1862  int len;
1863  MPI_Probe(MPI_ANY_SOURCE, 123, tr->get_communicator(), &status);
1864  MPI_Get_count(&status, MPI_BYTE, &len);
1865  receive.resize(len);
1866 
1867  char *ptr = &receive[0];
1868  MPI_Recv(ptr, len, MPI_BYTE, status.MPI_SOURCE, status.MPI_TAG,
1869  tr->get_communicator(), &status);
1870 
1871  unsigned int cells;
1872  memcpy(&cells, ptr, sizeof(unsigned int));
1873  ptr+=sizeof(unsigned int);
1874 
1875  //TODO: reinterpret too evil?
1876  unsigned int *treeindex=reinterpret_cast<unsigned int *>(ptr);
1877  ptr+=cells*sizeof(unsigned int);
1878  typename ::internal::p4est::types<dim>::quadrant *quadrant
1879  =reinterpret_cast<typename ::internal::p4est::types<dim>::quadrant *>(ptr);
1880  ptr+=cells*sizeof(typename ::internal::p4est::types<dim>::quadrant);
1881  ::types::global_dof_index *dofs
1882  = reinterpret_cast<::types::global_dof_index *>(ptr);
1883 
1884  // the dofs pointer contains for each cell the number of dofs
1885  // on that cell (dofs[0]) followed by the dof indices itself.
1886  for (unsigned int c=0; c<cells; ++c, dofs+=1+dofs[0])
1887  {
1888  typename DoFHandler<dim,spacedim>::level_cell_iterator
1889  cell (&dof_handler.get_triangulation(),
1890  0,
1891  p4est_tree_to_coarse_cell_permutation[treeindex[c]],
1892  &dof_handler);
1893 
1894  typename ::internal::p4est::types<dim>::quadrant p4est_coarse_cell;
1895  internal::p4est::init_coarse_quadrant<dim>(p4est_coarse_cell);
1896 
1897  Assert(cell->get_fe().dofs_per_cell==dofs[0], ExcInternalError());
1898 
1899  set_dofindices_recursively<dim,spacedim> (*tr,
1900  p4est_coarse_cell,
1901  cell,
1902  quadrant[c],
1903  (dofs+1));
1904  }
1905  }
1906 
1907  // complete all sends, so that we can
1908  // safely destroy the buffers.
1909  if (requests.size() > 0)
1910  MPI_Waitall(requests.size(), &requests[0], MPI_STATUSES_IGNORE);
1911 
1912 
1913 #ifdef DEBUG
1914  {
1915  //check all msgs got sent and received
1916  unsigned int sum_send=0;
1917  unsigned int sum_recv=0;
1918  unsigned int sent=needs_to_get_cells.size();
1919  unsigned int recv=senders.size();
1920 
1921  MPI_Allreduce(&sent, &sum_send, 1, MPI_UNSIGNED, MPI_SUM, tr->get_communicator());
1922  MPI_Allreduce(&recv, &sum_recv, 1, MPI_UNSIGNED, MPI_SUM, tr->get_communicator());
1923  Assert(sum_send==sum_recv, ExcInternalError());
1924  }
1925 #endif
1926 
1927  //update dofindices
1928  {
1929  typename DoFHandler<dim,spacedim>::active_cell_iterator
1930  cell, endc = dof_handler.end();
1931 
1932  for (cell = dof_handler.begin_active(); cell != endc; ++cell)
1933  if (!cell->is_artificial())
1934  cell->update_cell_dof_indices_cache();
1935  }
1936 
1937  // important, so that sends between two
1938  // calls to this function are not mixed
1939  // up.
1940  //
1941  // this is necessary because above we
1942  // just see if there are messages and
1943  // then receive them, without
1944  // discriminating where they come from
1945  // and whether they were sent in phase
1946  // 1 or 2. the need for a global
1947  // communication step like this barrier
1948  // could be avoided by receiving
1949  // messages specifically from those
1950  // processors from which we expect
1951  // messages, and by using different
1952  // tags for phase 1 and 2
1953  MPI_Barrier(tr->get_communicator());
1954 #endif
1955  }
1956 
1957 
1958 
1959  template <int spacedim>
1960  void
1961  communicate_mg_dof_indices_on_marked_cells
1962  (const DoFHandler<1,spacedim> &,
1963  const std::map<unsigned int, std::set<::types::subdomain_id> > &,
1964  const std::vector<::types::global_dof_index> &,
1965  const std::vector<::types::global_dof_index> &,
1966  const unsigned int)
1967  {
1968  Assert (false, ExcNotImplemented());
1969  }
1970 
1971 
1972 
1973  template <int dim, int spacedim>
1974  void
1975  communicate_mg_dof_indices_on_marked_cells
1976  (const DoFHandler<dim,spacedim> &dof_handler,
1977  const std::map<unsigned int, std::set<::types::subdomain_id> > &vertices_with_ghost_neighbors,
1978  const std::vector<::types::global_dof_index> &coarse_cell_to_p4est_tree_permutation,
1979  const std::vector<::types::global_dof_index> &p4est_tree_to_coarse_cell_permutation,
1980  const unsigned int level)
1981  {
1982 #ifndef DEAL_II_WITH_P4EST
1983  (void)dof_handler;
1984  (void)vertices_with_ghost_neighbors;
1985  (void)coarse_cell_to_p4est_tree_permutation;
1986  (void)p4est_tree_to_coarse_cell_permutation;
1987  (void)level;
1988  Assert (false, ExcNotImplemented());
1989 #else
1990 
1991  const parallel::distributed::Triangulation< dim, spacedim > *tr
1992  = (dynamic_cast<const parallel::distributed::Triangulation<dim,spacedim>*>
1993  (&dof_handler.get_triangulation()));
1994  Assert (tr != 0, ExcInternalError());
1995 
1996  // now collect cells and their
1997  // dof_indices for the
1998  // interested neighbors
1999  typedef
2000  std::map<::types::subdomain_id, typename types<dim>::cellinfo>
2001  cellmap_t;
2002  cellmap_t needs_to_get_cells;
2003 
2004  for (typename DoFHandler<dim,spacedim>::level_cell_iterator
2005  cell = dof_handler.begin(0);
2006  cell != dof_handler.end(0);
2007  ++cell)
2008  {
2009  typename ::internal::p4est::types<dim>::quadrant p4est_coarse_cell;
2010  internal::p4est::init_coarse_quadrant<dim>(p4est_coarse_cell);
2011 
2012  fill_mg_dofindices_recursively<dim,spacedim>
2013  (*tr,
2014  coarse_cell_to_p4est_tree_permutation[cell->index()],
2015  cell,
2016  p4est_coarse_cell,
2017  vertices_with_ghost_neighbors,
2018  needs_to_get_cells,
2019  level);
2020  }
2021 
2022 
2023  //sending
2024  std::vector<std::vector<char> > sendbuffers (needs_to_get_cells.size());
2025  std::vector<std::vector<char> >::iterator buffer = sendbuffers.begin();
2026  std::vector<MPI_Request> requests (needs_to_get_cells.size());
2027 
2028  unsigned int idx=0;
2029 
2030  for (typename cellmap_t::iterator it=needs_to_get_cells.begin();
2031  it!=needs_to_get_cells.end();
2032  ++it, ++buffer, ++idx)
2033  {
2034  const unsigned int num_cells = it->second.tree_index.size();
2035  (void)num_cells;
2036 
2037  Assert(num_cells==it->second.quadrants.size(), ExcInternalError());
2038  Assert(num_cells>0, ExcInternalError());
2039 
2040  // pack all the data into
2041  // the buffer for this
2042  // recipient and send
2043  // it. keep data around
2044  // till we can make sure
2045  // that the packet has been
2046  // received
2047  it->second.pack_data (*buffer);
2048  MPI_Isend(&(*buffer)[0], buffer->size(),
2049  MPI_BYTE, it->first,
2050  123, tr->get_communicator(), &requests[idx]);
2051  }
2052 
2053 
2054  // mark all own cells, that miss some dof_data and collect the
2055  // neighbors that are going to send stuff to us
2056  std::set<::types::subdomain_id> senders;
2057  {
2058  std::vector<::types::global_dof_index> local_dof_indices;
2059  typename DoFHandler<dim,spacedim>::level_cell_iterator
2060  cell, endc = dof_handler.end(level);
2061 
2062  for (cell = dof_handler.begin(level); cell != endc; ++cell)
2063  {
2064  if (cell->level_subdomain_id()==::numbers::artificial_subdomain_id)
2065  {
2066  //artificial
2067  }
2068  else if (cell->level_subdomain_id()==dof_handler.get_triangulation().locally_owned_subdomain())
2069  {
2070  //own
2071  local_dof_indices.resize (cell->get_fe().dofs_per_cell);
2072  cell->get_mg_dof_indices (local_dof_indices);
2073  if (local_dof_indices.end() !=
2074  std::find (local_dof_indices.begin(),
2075  local_dof_indices.end(),
2076  DoFHandler<dim,spacedim>::invalid_dof_index))
2077  cell->set_user_flag();
2078  else
2079  cell->clear_user_flag();
2080  }
2081  else
2082  {
2083  //ghost
2084  if (cell->user_flag_set())
2085  senders.insert(cell->level_subdomain_id());
2086  }
2087  }
2088 
2089  }
2090 
2091 
2092  //* 5. receive ghostcelldata
2093  std::vector<char> receive;
2094  typename types<dim>::cellinfo cellinfo;
2095  for (unsigned int i=0; i<senders.size(); ++i)
2096  {
2097  MPI_Status status;
2098  int len;
2099  MPI_Probe(MPI_ANY_SOURCE, 123, tr->get_communicator(), &status);
2100  MPI_Get_count(&status, MPI_BYTE, &len);
2101  receive.resize(len);
2102 
2103 #ifdef DEBUG
2104  Assert(senders.find(status.MPI_SOURCE)!=senders.end(), ExcInternalError());
2105 #endif
2106 
2107  char *ptr = &receive[0];
2108  MPI_Recv(ptr, len, MPI_BYTE, status.MPI_SOURCE, status.MPI_TAG,
2109  tr->get_communicator(), &status);
2110 
2111  unsigned int cells;
2112  memcpy(&cells, ptr, sizeof(unsigned int));
2113  ptr+=sizeof(unsigned int);
2114 
2115  //reinterpret too evil?
2116  unsigned int *treeindex=reinterpret_cast<unsigned int *>(ptr);
2117  ptr+=cells*sizeof(unsigned int);
2118  typename ::internal::p4est::types<dim>::quadrant *quadrant
2119  =reinterpret_cast<typename ::internal::p4est::types<dim>::quadrant *>(ptr);
2120  ptr+=cells*sizeof(typename ::internal::p4est::types<dim>::quadrant);
2121  ::types::global_dof_index *dofs
2122  = reinterpret_cast<::types::global_dof_index *>(ptr);
2123 
2124  // the dofs pointer contains for each cell the number of dofs
2125  // on that cell (dofs[0]) followed by the dof indices itself.
2126  for (unsigned int c=0; c<cells; ++c, dofs+=1+dofs[0])
2127  {
2128  typename DoFHandler<dim,spacedim>::level_cell_iterator
2129  cell (&dof_handler.get_triangulation(),
2130  0,
2131  p4est_tree_to_coarse_cell_permutation[treeindex[c]],
2132  &dof_handler);
2133 
2134  typename ::internal::p4est::types<dim>::quadrant p4est_coarse_cell;
2135  internal::p4est::init_coarse_quadrant<dim>(p4est_coarse_cell);
2136 
2137  Assert(cell->get_fe().dofs_per_cell==dofs[0], ExcInternalError());
2138 
2139  set_mg_dofindices_recursively<dim,spacedim> (*tr,
2140  p4est_coarse_cell,
2141  cell,
2142  quadrant[c],
2143  (dofs+1),
2144  level);
2145  }
2146  }
2147 
2148  // complete all sends, so that we can
2149  // safely destroy the buffers.
2150  if (requests.size() > 0)
2151  MPI_Waitall(requests.size(), &requests[0], MPI_STATUSES_IGNORE);
2152 
2153 
2154 #ifdef DEBUG
2155  {
2156  //check all msgs got sent and received
2157  unsigned int sum_send=0;
2158  unsigned int sum_recv=0;
2159  unsigned int sent=needs_to_get_cells.size();
2160  unsigned int recv=senders.size();
2161 
2162  MPI_Allreduce(&sent, &sum_send, 1, MPI_UNSIGNED, MPI_SUM, tr->get_communicator());
2163  MPI_Allreduce(&recv, &sum_recv, 1, MPI_UNSIGNED, MPI_SUM, tr->get_communicator());
2164  Assert(sum_send==sum_recv, ExcInternalError());
2165  }
2166 #endif
2167 
2168 
2169  // important, so that sends between two
2170  // calls to this function are not mixed
2171  // up.
2172  //
2173  // this is necessary because above we
2174  // just see if there are messages and
2175  // then receive them, without
2176  // discriminating where they come from
2177  // and whether they were sent in phase
2178  // 1 or 2. the need for a global
2179  // communication step like this barrier
2180  // could be avoided by receiving
2181  // messages specifically from those
2182  // processors from which we expect
2183  // messages, and by using different
2184  // tags for phase 1 and 2
2185  MPI_Barrier(tr->get_communicator());
2186 #endif
2187  }
2188 
2189 
2190 
2191 
2192  }
2193 
2194 #endif // DEAL_II_WITH_P4EST
2195 
2196 
2197 
2198  template <int dim, int spacedim>
2199  void
2200  ParallelDistributed<dim, spacedim>::
2201  distribute_dofs (DoFHandler<dim,spacedim> &dof_handler,
2202  NumberCache &number_cache_current) const
2203  {
2204  NumberCache number_cache;
2205 
2206 #ifndef DEAL_II_WITH_P4EST
2207  (void)dof_handler;
2208  Assert (false, ExcNotImplemented());
2209 #else
2210 
2211  parallel::distributed::Triangulation< dim, spacedim > *tr
2212  = (dynamic_cast<parallel::distributed::Triangulation<dim,spacedim>*>
2213  (const_cast<::Triangulation< dim, spacedim >*>
2214  (&dof_handler.get_triangulation())));
2215  Assert (tr != 0, ExcInternalError());
2216 
2217  const unsigned int
2218  n_cpus = Utilities::MPI::n_mpi_processes (tr->get_communicator());
2219 
2220  //* 1. distribute on own
2221  //* subdomain
2222  const ::types::global_dof_index n_initial_local_dofs =
2223  Implementation::distribute_dofs (0, tr->locally_owned_subdomain(),
2224  dof_handler);
2225 
2226  //* 2. iterate over ghostcells and
2227  //kill dofs that are not owned
2228  //by us
2229  std::vector<::types::global_dof_index> renumbering(n_initial_local_dofs);
2230  for (unsigned int i=0; i<renumbering.size(); ++i)
2231  renumbering[i] = i;
2232 
2233  {
2234  std::vector<::types::global_dof_index> local_dof_indices;
2235 
2236  typename DoFHandler<dim,spacedim>::active_cell_iterator
2237  cell = dof_handler.begin_active(),
2238  endc = dof_handler.end();
2239 
2240  for (; cell != endc; ++cell)
2241  if (cell->is_ghost() &&
2242  (cell->subdomain_id() < tr->locally_owned_subdomain()))
2243  {
2244  // we found a
2245  // neighboring ghost
2246  // cell whose subdomain
2247  // is "stronger" than
2248  // our own subdomain
2249 
2250  // delete all dofs that
2251  // live there and that
2252  // we have previously
2253  // assigned a number to
2254  // (i.e. the ones on
2255  // the interface)
2256  local_dof_indices.resize (cell->get_fe().dofs_per_cell);
2257  cell->get_dof_indices (local_dof_indices);
2258  for (unsigned int i=0; i<cell->get_fe().dofs_per_cell; ++i)
2259  if (local_dof_indices[i] != DoFHandler<dim,spacedim>::invalid_dof_index)
2260  renumbering[local_dof_indices[i]]
2261  = DoFHandler<dim,spacedim>::invalid_dof_index;
2262  }
2263  }
2264 
2265 
2266  // make indices consecutive
2267  number_cache.n_locally_owned_dofs = 0;
2268  for (std::vector<::types::global_dof_index>::iterator it=renumbering.begin();
2269  it!=renumbering.end(); ++it)
2270  if (*it != DoFHandler<dim,spacedim>::invalid_dof_index)
2271  *it = number_cache.n_locally_owned_dofs++;
2272 
2273  //* 3. communicate local dofcount and
2274  //shift ids to make them unique
2275  number_cache.n_locally_owned_dofs_per_processor.resize(n_cpus);
2276 
2277  MPI_Allgather ( &number_cache.n_locally_owned_dofs,
2278  1, DEAL_II_DOF_INDEX_MPI_TYPE,
2279  &number_cache.n_locally_owned_dofs_per_processor[0],
2280  1, DEAL_II_DOF_INDEX_MPI_TYPE,
2281  tr->get_communicator());
2282 
2283  const ::types::global_dof_index
2284  shift = std::accumulate (number_cache
2285  .n_locally_owned_dofs_per_processor.begin(),
2286  number_cache
2287  .n_locally_owned_dofs_per_processor.begin()
2288  + tr->locally_owned_subdomain(),
2289  static_cast<::types::global_dof_index>(0));
2290  for (std::vector<::types::global_dof_index>::iterator it=renumbering.begin();
2291  it!=renumbering.end(); ++it)
2292  if (*it != DoFHandler<dim,spacedim>::invalid_dof_index)
2293  (*it) += shift;
2294 
2295  // now re-enumerate all dofs to
2296  // this shifted and condensed
2297  // numbering form. we renumber
2298  // some dofs as invalid, so
2299  // choose the nocheck-version.
2300  Implementation::renumber_dofs (renumbering, IndexSet(0),
2301  dof_handler, false);
2302 
2303  // now a little bit of
2304  // housekeeping
2305  number_cache.n_global_dofs
2306  = std::accumulate (number_cache
2307  .n_locally_owned_dofs_per_processor.begin(),
2308  number_cache
2309  .n_locally_owned_dofs_per_processor.end(),
2310  static_cast<::types::global_dof_index>(0));
2311 
2312  number_cache.locally_owned_dofs = IndexSet(number_cache.n_global_dofs);
2313  number_cache.locally_owned_dofs
2314  .add_range(shift,
2315  shift+number_cache.n_locally_owned_dofs);
2316  number_cache.locally_owned_dofs.compress();
2317 
2318  // fill global_dof_indexsets
2319  number_cache.locally_owned_dofs_per_processor.resize(n_cpus);
2320  {
2321  ::types::global_dof_index lshift = 0;
2322  for (unsigned int i=0; i<n_cpus; ++i)
2323  {
2324  number_cache.locally_owned_dofs_per_processor[i]
2325  = IndexSet(number_cache.n_global_dofs);
2326  number_cache.locally_owned_dofs_per_processor[i]
2327  .add_range(lshift,
2328  lshift +
2329  number_cache.n_locally_owned_dofs_per_processor[i]);
2330  lshift += number_cache.n_locally_owned_dofs_per_processor[i];
2331  }
2332  }
2333  Assert(number_cache.locally_owned_dofs_per_processor
2334  [tr->locally_owned_subdomain()].n_elements()
2335  ==
2336  number_cache.n_locally_owned_dofs,
2337  ExcInternalError());
2338  Assert(!number_cache.locally_owned_dofs_per_processor
2339  [tr->locally_owned_subdomain()].n_elements()
2340  ||
2341  number_cache.locally_owned_dofs_per_processor
2342  [tr->locally_owned_subdomain()].nth_index_in_set(0)
2343  == shift,
2344  ExcInternalError());
2345 
2346  //* 4. send dofids of cells that are
2347  //ghostcells on other machines
2348 
2349  std::vector<bool> user_flags;
2350  tr->save_user_flags(user_flags);
2351  tr->clear_user_flags ();
2352 
2353  //mark all own cells for transfer
2354  for (typename DoFHandler<dim,spacedim>::active_cell_iterator cell = dof_handler.begin_active();
2355  cell != dof_handler.end(); ++cell)
2356  if (!cell->is_artificial())
2357  cell->set_user_flag();
2358 
2359  // add each ghostcells'
2360  // subdomain to the vertex and
2361  // keep track of interesting
2362  // neighbors
2363  std::map<unsigned int, std::set<::types::subdomain_id> >
2364  vertices_with_ghost_neighbors;
2365 
2366  tr->fill_vertices_with_ghost_neighbors (vertices_with_ghost_neighbors);
2367 
2368 
2369  /* Send and receive cells. After this,
2370  only the local cells are marked,
2371  that received new data. This has to
2372  be communicated in a second
2373  communication step. */
2374  communicate_dof_indices_on_marked_cells (dof_handler,
2375  vertices_with_ghost_neighbors,
2376  tr->coarse_cell_to_p4est_tree_permutation,
2377  tr->p4est_tree_to_coarse_cell_permutation);
2378 
2379  communicate_dof_indices_on_marked_cells (dof_handler,
2380  vertices_with_ghost_neighbors,
2381  tr->coarse_cell_to_p4est_tree_permutation,
2382  tr->p4est_tree_to_coarse_cell_permutation);
2383 
2384  tr->load_user_flags(user_flags);
2385 
2386 #ifdef DEBUG
2387  //check that we are really done
2388  {
2389  std::vector<::types::global_dof_index> local_dof_indices;
2390 
2391  for (typename DoFHandler<dim,spacedim>::active_cell_iterator cell = dof_handler.begin_active();
2392  cell != dof_handler.end(); ++cell)
2393  if (!cell->is_artificial())
2394  {
2395  local_dof_indices.resize (cell->get_fe().dofs_per_cell);
2396  cell->get_dof_indices (local_dof_indices);
2397  if (local_dof_indices.end() !=
2398  std::find (local_dof_indices.begin(),
2399  local_dof_indices.end(),
2400  DoFHandler<dim,spacedim>::invalid_dof_index))
2401  {
2402  if (cell->is_ghost())
2403  {
2404  Assert(false, ExcMessage ("Not a ghost cell"));
2405  }
2406  else
2407  {
2408  Assert(false, ExcMessage ("Not one of our own cells"));
2409  }
2410  }
2411  }
2412  }
2413 #endif // DEBUG
2414 #endif // DEAL_II_WITH_P4EST
2415 
2416  number_cache_current = number_cache;
2417  }
2418 
2419 
2420 
2421  template <int dim, int spacedim>
2422  void
2423  ParallelDistributed<dim, spacedim>::
2424  distribute_mg_dofs (DoFHandler<dim,spacedim> &dof_handler,
2425  std::vector<NumberCache> &number_caches) const
2426  {
2427 #ifndef DEAL_II_WITH_P4EST
2428  (void)dof_handler;
2429  (void)number_caches;
2430  Assert (false, ExcNotImplemented());
2431 #else
2432 
2433  parallel::distributed::Triangulation< dim, spacedim > *tr
2434  = (dynamic_cast<parallel::distributed::Triangulation<dim,spacedim>*>
2435  (const_cast<::Triangulation< dim, spacedim >*>
2436  (&dof_handler.get_triangulation())));
2437  Assert (tr != 0, ExcInternalError());
2438 
2439  AssertThrow(
2440  (tr->settings & parallel::distributed::Triangulation< dim, spacedim >::construct_multigrid_hierarchy),
2441  ExcMessage("Multigrid DoFs can only be distributed on a parallel Triangulation if the flag construct_multigrid_hierarchy is set in the constructor."));
2442 
2443 
2444  const unsigned int
2445  n_cpus = Utilities::MPI::n_mpi_processes (tr->get_communicator());
2446 
2447  unsigned int n_levels = Utilities::MPI::max(dof_handler.get_triangulation().n_levels(), tr->get_communicator());
2448 
2449  for (unsigned int level = 0; level < n_levels; ++level)
2450  {
2451  NumberCache &number_cache = number_caches[level];
2452 
2453  //* 1. distribute on own
2454  //* subdomain
2455  const unsigned int n_initial_local_dofs =
2456  Implementation::distribute_dofs_on_level(0, tr->locally_owned_subdomain(), dof_handler, level);
2457 
2458  //* 2. iterate over ghostcells and
2459  //kill dofs that are not owned
2460  //by us
2461  std::vector<::types::global_dof_index> renumbering(n_initial_local_dofs);
2462  for (::types::global_dof_index i=0; i<renumbering.size(); ++i)
2463  renumbering[i] = i;
2464 
2465  if (level<tr->n_levels())
2466  {
2467  std::vector<::types::global_dof_index> local_dof_indices;
2468 
2469  typename DoFHandler<dim,spacedim>::level_cell_iterator
2470  cell = dof_handler.begin(level),
2471  endc = dof_handler.end(level);
2472 
2473  for (; cell != endc; ++cell)
2474  if (cell->level_subdomain_id()!=numbers::artificial_subdomain_id &&
2475  (cell->level_subdomain_id() < tr->locally_owned_subdomain()))
2476  {
2477  // we found a
2478  // neighboring ghost
2479  // cell whose subdomain
2480  // is "stronger" than
2481  // our own subdomain
2482 
2483  // delete all dofs that
2484  // live there and that
2485  // we have previously
2486  // assigned a number to
2487  // (i.e. the ones on
2488  // the interface)
2489  local_dof_indices.resize (cell->get_fe().dofs_per_cell);
2490  cell->get_mg_dof_indices (local_dof_indices);
2491  for (unsigned int i=0; i<cell->get_fe().dofs_per_cell; ++i)
2492  if (local_dof_indices[i] != DoFHandler<dim,spacedim>::invalid_dof_index)
2493  renumbering[local_dof_indices[i]]
2494  = DoFHandler<dim,spacedim>::invalid_dof_index;
2495  }
2496  }
2497 
2498 
2499  // make indices consecutive
2500  number_cache.n_locally_owned_dofs = 0;
2501  for (std::vector<::types::global_dof_index>::iterator it=renumbering.begin();
2502  it!=renumbering.end(); ++it)
2503  if (*it != DoFHandler<dim,spacedim>::invalid_dof_index)
2504  *it = number_cache.n_locally_owned_dofs++;
2505 
2506  //* 3. communicate local dofcount and
2507  //shift ids to make them unique
2508  number_cache.n_locally_owned_dofs_per_processor.resize(n_cpus);
2509 
2510  MPI_Allgather ( &number_cache.n_locally_owned_dofs,
2511  1, DEAL_II_DOF_INDEX_MPI_TYPE,
2512  &number_cache.n_locally_owned_dofs_per_processor[0],
2513  1, DEAL_II_DOF_INDEX_MPI_TYPE,
2514  tr->get_communicator());
2515 
2516  const ::types::global_dof_index
2517  shift = std::accumulate (number_cache
2518  .n_locally_owned_dofs_per_processor.begin(),
2519  number_cache
2520  .n_locally_owned_dofs_per_processor.begin()
2521  + tr->locally_owned_subdomain(),
2522  static_cast<::types::global_dof_index>(0));
2523  for (std::vector<::types::global_dof_index>::iterator it=renumbering.begin();
2524  it!=renumbering.end(); ++it)
2525  if (*it != DoFHandler<dim,spacedim>::invalid_dof_index)
2526  (*it) += shift;
2527 
2528  // now re-enumerate all dofs to
2529  // this shifted and condensed
2530  // numbering form. we renumber
2531  // some dofs as invalid, so
2532  // choose the nocheck-version.
2533  Implementation::renumber_mg_dofs (renumbering, IndexSet(0),
2534  dof_handler, level, false);
2535 
2536  // now a little bit of
2537  // housekeeping
2538  number_cache.n_global_dofs
2539  = std::accumulate (number_cache
2540  .n_locally_owned_dofs_per_processor.begin(),
2541  number_cache
2542  .n_locally_owned_dofs_per_processor.end(),
2543  static_cast<::types::global_dof_index>(0));
2544 
2545  number_cache.locally_owned_dofs = IndexSet(number_cache.n_global_dofs);
2546  number_cache.locally_owned_dofs
2547  .add_range(shift,
2548  shift+number_cache.n_locally_owned_dofs);
2549  number_cache.locally_owned_dofs.compress();
2550 
2551  // fill global_dof_indexsets
2552  number_cache.locally_owned_dofs_per_processor.resize(n_cpus);
2553  {
2554  ::types::global_dof_index lshift = 0;
2555  for (unsigned int i=0; i<n_cpus; ++i)
2556  {
2557  number_cache.locally_owned_dofs_per_processor[i]
2558  = IndexSet(number_cache.n_global_dofs);
2559  number_cache.locally_owned_dofs_per_processor[i]
2560  .add_range(lshift,
2561  lshift +
2562  number_cache.n_locally_owned_dofs_per_processor[i]);
2563  lshift += number_cache.n_locally_owned_dofs_per_processor[i];
2564  }
2565  }
2566  Assert(number_cache.locally_owned_dofs_per_processor
2567  [tr->locally_owned_subdomain()].n_elements()
2568  ==
2569  number_cache.n_locally_owned_dofs,
2570  ExcInternalError());
2571  Assert(!number_cache.locally_owned_dofs_per_processor
2572  [tr->locally_owned_subdomain()].n_elements()
2573  ||
2574  number_cache.locally_owned_dofs_per_processor
2575  [tr->locally_owned_subdomain()].nth_index_in_set(0)
2576  == shift,
2577  ExcInternalError());
2578 
2579  //* 4. send dofids of cells that are
2580  //ghostcells on other machines
2581  std::vector<bool> user_flags;
2582  tr->save_user_flags(user_flags);
2583  tr->clear_user_flags ();
2584 
2585  //mark all own cells for transfer
2586  if (level < tr->n_levels())
2587  {
2588  typename DoFHandler<dim,spacedim>::level_cell_iterator
2589  cell, endc = dof_handler.end(level);
2590  for (cell = dof_handler.begin(level); cell != endc; ++cell)
2591  if (cell->level_subdomain_id() != ::numbers::artificial_subdomain_id)
2592  cell->set_user_flag();
2593  }
2594 
2595  //mark the vertices we are interested in, i.e. belonging to own
2596  //and marked cells
2597  const std::vector<bool> locally_active_vertices
2598  = tr->mark_locally_active_vertices_on_level (level);
2599 
2600  // add each ghostcells' subdomain to the vertex and keep track of
2601  // interesting neighbors
2602  std::map<unsigned int, std::set<::types::subdomain_id> >
2603  vertices_with_ghost_neighbors;
2604  tr->fill_level_vertices_with_ghost_neighbors(level,
2605  vertices_with_ghost_neighbors);
2606 
2607  // Send and receive cells. After this, only the local cells are
2608  // marked, that received new data. This has to be communicated in
2609  // a second communication step.
2610 
2611  communicate_mg_dof_indices_on_marked_cells( dof_handler,
2612  vertices_with_ghost_neighbors,
2613  tr->coarse_cell_to_p4est_tree_permutation,
2614  tr->p4est_tree_to_coarse_cell_permutation,
2615  level);
2616  communicate_mg_dof_indices_on_marked_cells( dof_handler,
2617  vertices_with_ghost_neighbors,
2618  tr->coarse_cell_to_p4est_tree_permutation,
2619  tr->p4est_tree_to_coarse_cell_permutation,
2620  level);
2621 
2622  tr->load_user_flags(user_flags);
2623 
2624 #ifdef DEBUG
2625  //check that we are really done
2626  if (level < tr->n_levels())
2627  {
2628  std::vector<::types::global_dof_index> local_dof_indices;
2629  typename DoFHandler<dim,spacedim>::level_cell_iterator
2630  cell, endc = dof_handler.end(level);
2631 
2632  for (cell = dof_handler.begin(level); cell != endc; ++cell)
2633  if (cell->level_subdomain_id() != ::numbers::artificial_subdomain_id)
2634  {
2635  local_dof_indices.resize (cell->get_fe().dofs_per_cell);
2636  cell->get_mg_dof_indices (local_dof_indices);
2637  if (local_dof_indices.end() !=
2638  std::find (local_dof_indices.begin(),
2639  local_dof_indices.end(),
2640  DoFHandler<dim,spacedim>::invalid_dof_index))
2641  {
2642  Assert(false, ExcMessage ("not all DoFs got distributed!"));
2643  }
2644  }
2645  }
2646 #endif // DEBUG
2647 
2648  }
2649 
2650 #endif // DEAL_II_WITH_P4EST
2651  }
2652 
2653 
2654  template <int dim, int spacedim>
2655  void
2656  ParallelDistributed<dim, spacedim>::
2657  renumber_dofs (const std::vector<::types::global_dof_index> &new_numbers,
2658  ::DoFHandler<dim,spacedim> &dof_handler,
2659  NumberCache &number_cache_current) const
2660  {
2661  (void)new_numbers;
2662  (void)dof_handler;
2663 
2664  Assert (new_numbers.size() == dof_handler.locally_owned_dofs().n_elements(),
2665  ExcInternalError());
2666 
2667  NumberCache number_cache;
2668 
2669 #ifndef DEAL_II_WITH_P4EST
2670  Assert (false, ExcNotImplemented());
2671 #else
2672 
2673 
2674  // calculate new IndexSet. First try to find out if the new indices
2675  // are contiguous blocks. This avoids inserting each index
2676  // individually into the IndexSet, which is slow. If we own no DoFs,
2677  // we still need to go through this function, but we can skip this
2678  // calculation.
2679 
2680  number_cache.locally_owned_dofs = IndexSet (dof_handler.n_dofs());
2681  if (dof_handler.locally_owned_dofs().n_elements()>0)
2682  {
2683  std::vector<::types::global_dof_index> new_numbers_sorted (new_numbers);
2684  std::sort(new_numbers_sorted.begin(), new_numbers_sorted.end());
2685  std::vector<::types::global_dof_index>::const_iterator it = new_numbers_sorted.begin();
2686  const unsigned int n_blocks = dof_handler.get_fe().n_blocks();
2687  std::vector<std::pair<::types::global_dof_index,unsigned int> > block_indices(n_blocks);
2688  block_indices[0].first = *it++;
2689  block_indices[0].second = 1;
2690  unsigned int current_block = 0, n_filled_blocks = 1;
2691  for ( ; it != new_numbers_sorted.end(); ++it)
2692  {
2693  bool done = false;
2694 
2695  // search from the current block onwards whether the next
2696  // index is shifted by one from the previous one.
2697  for (unsigned int i=0; i<n_filled_blocks; ++i)
2698  if (*it == block_indices[current_block].first
2699  +block_indices[current_block].second)
2700  {
2701  block_indices[current_block].second++;
2702  done = true;
2703  break;
2704  }
2705  else
2706  {
2707  if (current_block == n_filled_blocks-1)
2708  current_block = 0;
2709  else
2710  ++current_block;
2711  }
2712 
2713  // could not find any contiguous range: need to add a new
2714  // block if possible. Abort otherwise, which will add all
2715  // elements individually to the IndexSet.
2716  if (done == false)
2717  {
2718  if (n_filled_blocks < n_blocks)
2719  {
2720  block_indices[n_filled_blocks].first = *it;
2721  block_indices[n_filled_blocks].second = 1;
2722  current_block = n_filled_blocks;
2723  ++n_filled_blocks;
2724  }
2725  else
2726  break;
2727  }
2728  }
2729 
2730  // check whether all indices could be assigned to blocks. If yes,
2731  // we can add the block ranges to the IndexSet, otherwise we need
2732  // to go through the indices once again and add each element
2733  // individually
2734  unsigned int sum = 0;
2735  for (unsigned int i=0; i<n_filled_blocks; ++i)
2736  sum += block_indices[i].second;
2737  if (sum == new_numbers.size())
2738  for (unsigned int i=0; i<n_filled_blocks; ++i)
2739  number_cache.locally_owned_dofs.add_range (block_indices[i].first,
2740  block_indices[i].first+
2741  block_indices[i].second);
2742  else
2743  number_cache.locally_owned_dofs.add_indices(new_numbers_sorted.begin(),
2744  new_numbers_sorted.end());
2745  }
2746 
2747 
2748  number_cache.locally_owned_dofs.compress();
2749  Assert (number_cache.locally_owned_dofs.n_elements() == new_numbers.size(),
2750  ExcInternalError());
2751  // also check with the number of locally owned degrees of freedom that
2752  // the DoFHandler object still stores
2753  Assert (number_cache.locally_owned_dofs.n_elements() ==
2754  dof_handler.n_locally_owned_dofs(),
2755  ExcInternalError());
2756 
2757  // then also set this number in our own copy
2758  number_cache.n_locally_owned_dofs = dof_handler.n_locally_owned_dofs();
2759 
2760  // mark not locally active DoFs as invalid
2761  {
2762  std::vector<::types::global_dof_index> local_dof_indices;
2763 
2764  typename DoFHandler<dim,spacedim>::active_cell_iterator
2765  cell = dof_handler.begin_active(),
2766  endc = dof_handler.end();
2767 
2768  for (; cell != endc; ++cell)
2769  if (!cell->is_artificial())
2770  {
2771  local_dof_indices.resize (cell->get_fe().dofs_per_cell);
2772  cell->get_dof_indices (local_dof_indices);
2773  for (unsigned int i=0; i<cell->get_fe().dofs_per_cell; ++i)
2774  {
2775  if (local_dof_indices[i] == DoFHandler<dim,spacedim>::invalid_dof_index)
2776  continue;
2777 
2778  if (!dof_handler.locally_owned_dofs().is_element(local_dof_indices[i]))
2779  {
2780  //this DoF is not owned by us, so set it to invalid.
2781  local_dof_indices[i]
2782  = DoFHandler<dim,spacedim>::invalid_dof_index;
2783  }
2784  }
2785 
2786  cell->set_dof_indices (local_dof_indices);
2787  }
2788  }
2789 
2790 
2791  // renumber. Skip when there is nothing to do because we own no DoF.
2792  if (dof_handler.locally_owned_dofs().n_elements() > 0)
2793  Implementation::renumber_dofs (new_numbers,
2794  dof_handler.locally_owned_dofs(),
2795  dof_handler,
2796  false);
2797 
2798  // communication
2799  {
2800  parallel::distributed::Triangulation< dim, spacedim > *tr
2801  = (dynamic_cast<parallel::distributed::Triangulation<dim,spacedim>*>
2802  (const_cast<::Triangulation< dim, spacedim >*>
2803  (&dof_handler.get_triangulation())));
2804  Assert (tr != 0, ExcInternalError());
2805 
2806  std::vector<bool> user_flags;
2807  tr->save_user_flags(user_flags);
2808  tr->clear_user_flags ();
2809 
2810  //mark all own cells for transfer
2811  typename DoFHandler<dim,spacedim>::active_cell_iterator
2812  cell, endc = dof_handler.end();
2813  for (cell = dof_handler.begin_active(); cell != endc; ++cell)
2814  if (!cell->is_artificial())
2815  cell->set_user_flag();
2816 
2817  // add each ghostcells' subdomain to the vertex and keep track of
2818  // interesting neighbors
2819  std::map<unsigned int, std::set<::types::subdomain_id> >
2820  vertices_with_ghost_neighbors;
2821 
2822  tr->fill_vertices_with_ghost_neighbors (vertices_with_ghost_neighbors);
2823 
2824  // Send and receive cells. After this, only the local cells are
2825  // marked, that received new data. This has to be communicated in a
2826  // second communication step.
2827  communicate_dof_indices_on_marked_cells (dof_handler,
2828  vertices_with_ghost_neighbors,
2829  tr->coarse_cell_to_p4est_tree_permutation,
2830  tr->p4est_tree_to_coarse_cell_permutation);
2831 
2832  communicate_dof_indices_on_marked_cells (dof_handler,
2833  vertices_with_ghost_neighbors,
2834  tr->coarse_cell_to_p4est_tree_permutation,
2835  tr->p4est_tree_to_coarse_cell_permutation);
2836 
2837 
2838  // * Create global_dof_indexsets by transferring our own owned_dofs
2839  // to every other machine.
2840  const unsigned int n_cpus = Utilities::MPI::n_mpi_processes (tr->get_communicator());
2841 
2842  // Serialize our IndexSet and determine size.
2843  std::ostringstream oss;
2844  number_cache.locally_owned_dofs.block_write(oss);
2845  std::string oss_str=oss.str();
2846  std::vector<char> my_data(oss_str.begin(), oss_str.end());
2847  unsigned int my_size = oss_str.size();
2848 
2849  // determine maximum size of IndexSet
2850  const unsigned int max_size
2851  = Utilities::MPI::max (my_size, tr->get_communicator());
2852 
2853  // as we are reading past the end, we need to increase the size of
2854  // the local buffer. This is filled with zeros.
2855  my_data.resize(max_size);
2856 
2857  std::vector<char> buffer(max_size*n_cpus);
2858  MPI_Allgather(&my_data[0], max_size, MPI_BYTE,
2859  &buffer[0], max_size, MPI_BYTE,
2860  tr->get_communicator());
2861 
2862  number_cache.locally_owned_dofs_per_processor.resize (n_cpus);
2863  number_cache.n_locally_owned_dofs_per_processor.resize (n_cpus);
2864  for (unsigned int i=0; i<n_cpus; ++i)
2865  {
2866  std::stringstream strstr;
2867  strstr.write(&buffer[i*max_size],max_size);
2868  // This does not read the whole buffer, when the size is smaller
2869  // than max_size. Therefore we need to create a new stringstream
2870  // in each iteration (resetting would be fine too).
2871  number_cache.locally_owned_dofs_per_processor[i]
2872  .block_read(strstr);
2873  number_cache.n_locally_owned_dofs_per_processor[i]
2874  = number_cache.locally_owned_dofs_per_processor[i].n_elements();
2875  }
2876 
2877  number_cache.n_global_dofs
2878  = std::accumulate (number_cache
2879  .n_locally_owned_dofs_per_processor.begin(),
2880  number_cache
2881  .n_locally_owned_dofs_per_processor.end(),
2882  static_cast<::types::global_dof_index>(0));
2883 
2884  tr->load_user_flags(user_flags);
2885  }
2886 #endif
2887 
2888  number_cache_current = number_cache;
2889  }
2890  }
2891  }
2892 }
2893 
2894 
2895 
2896 
2897 /*-------------- Explicit Instantiations -------------------------------*/
2898 #include "dof_handler_policy.inst"
2899 
2900 
2901 DEAL_II_NAMESPACE_CLOSE
Triangulation::n_active_cells
unsigned int n_active_cells() const
Definition: tria.cc:10973
DoFHandler::mg_vertex_dofs
std::vector< MGVertexDoFs > mg_vertex_dofs
Definition: dof_handler.h:1062
parallel::distributed::Triangulation::mark_locally_active_vertices_on_level
std::vector< bool > mark_locally_active_vertices_on_level(const unsigned int level) const
Definition: tria.cc:4644
DoFTools::locally_owned_dofs_per_subdomain
std::vector< IndexSet > locally_owned_dofs_per_subdomain(const DoFHandlerType &dof_handler)
Definition: dof_tools.cc:1130
parallel::distributed::Triangulation::fill_vertices_with_ghost_neighbors
void fill_vertices_with_ghost_neighbors(std::map< unsigned int, std::set<::types::subdomain_id > > &vertices_with_ghost_neighbors)
Definition: tria.cc:4574
parallel::shared::Triangulation::get_true_subdomain_ids_of_cells
const std::vector< types::subdomain_id > & get_true_subdomain_ids_of_cells() const
Definition: shared_tria.cc:108
numbers::invalid_subdomain_id
const types::subdomain_id invalid_subdomain_id
Definition: types.h:239
internal::DoFHandler::Policy::Sequential::distribute_dofs
virtual void distribute_dofs(::DoFHandler< dim, spacedim > &dof_handler, NumberCache &number_cache) const
Definition: dof_handler_policy.cc:970
Triangulation::clear_user_flags
void clear_user_flags()
Definition: tria.cc:9618
DoFHandler::begin
cell_iterator begin(const unsigned int level=0) const
Definition: dof_handler.cc:833
IndexSet::nth_index_in_set
size_type nth_index_in_set(const unsigned int local_index) const
Definition: index_set.h:1432
internal::DoFHandler::Policy::ParallelDistributed::renumber_dofs
virtual void renumber_dofs(const std::vector< types::global_dof_index > &new_numbers, ::DoFHandler< dim, spacedim > &dof_handler, NumberCache &number_cache) const
Definition: dof_handler_policy.cc:2657
Triangulation::load_user_flags
void load_user_flags(std::istream &in)
Definition: tria.cc:9675
StandardExceptions::ExcMessage
::ExceptionBase & ExcMessage(std::string arg1)
DoFHandler::end
cell_iterator end() const
Definition: dof_handler.cc:861
internal::DoFHandler::Policy::ParallelShared::distribute_dofs
virtual void distribute_dofs(::DoFHandler< dim, spacedim > &dof_handler, NumberCache &number_cache) const
Definition: dof_handler_policy.cc:1068
DoFHandler::selected_fe
SmartPointer< const FiniteElement< dim, spacedim >, DoFHandler< dim, spacedim > > selected_fe
Definition: dof_handler.h:925
DoFHandler::levels
std::vector<::internal::DoFHandler::DoFLevel< dim > * > levels
Definition: dof_handler.h:1068
DoFHandler::active_cell_iterator
ActiveSelector::active_cell_iterator active_cell_iterator
Definition: dof_handler.h:221
Triangulation::begin_active
active_cell_iterator begin_active(const unsigned int level=0) const
Definition: tria.cc:10397
AssertThrow
#define AssertThrow(cond, exc)
Definition: exceptions.h:358
parallel::distributed::Triangulation::fill_level_vertices_with_ghost_neighbors
void fill_level_vertices_with_ghost_neighbors(const unsigned int level, std::map< unsigned int, std::set<::types::subdomain_id > > &vertices_with_ghost_neighbors)
Definition: tria.cc:4620
DoFHandler::get_fe
const FiniteElement< dim, spacedim > & get_fe() const
internal::DoFHandler::NumberCache::n_locally_owned_dofs
types::global_dof_index n_locally_owned_dofs
Definition: number_cache.h:63
DoFHandler::begin_active
active_cell_iterator begin_active(const unsigned int level=0) const
Definition: dof_handler.cc:845
DoFHandler::n_locally_owned_dofs
unsigned int n_locally_owned_dofs() const
DoFHandler::cell_iterator
ActiveSelector::cell_iterator cell_iterator
Definition: dof_handler.h:249
types
Definition: types.h:30
types::global_dof_index
unsigned int global_dof_index
Definition: types.h:88
Assert
#define Assert(cond, exc)
Definition: exceptions.h:294
internal::DoFHandler::Policy::Sequential::renumber_dofs
virtual void renumber_dofs(const std::vector< types::global_dof_index > &new_numbers, ::DoFHandler< dim, spacedim > &dof_handler, NumberCache &number_cache) const
Definition: dof_handler_policy.cc:1030
IndexSet::index_within_set
size_type index_within_set(const size_type global_index) const
Definition: index_set.h:1475
internal::DoFHandler::Policy::ParallelShared::distribute_mg_dofs
virtual void distribute_mg_dofs(::DoFHandler< dim, spacedim > &dof_handler, std::vector< NumberCache > &number_caches) const
Definition: dof_handler_policy.cc:1115
parallel::Triangulation::get_communicator
virtual MPI_Comm get_communicator() const
Definition: tria_base.cc:135
DoFHandler::n_dofs
types::global_dof_index n_dofs() const
DoFHandler::faces
::internal::DoFHandler::DoFFaces< dim > * faces
Definition: dof_handler.h:1077
types::subdomain_id
unsigned int subdomain_id
Definition: types.h:42
internal::DoFHandler::Policy::ParallelDistributed::distribute_dofs
virtual void distribute_dofs(::DoFHandler< dim, spacedim > &dof_handler, NumberCache &number_cache) const
Definition: dof_handler_policy.cc:2201
Utilities::MPI::n_mpi_processes
unsigned int n_mpi_processes(const MPI_Comm &mpi_communicator)
Definition: mpi.cc:99
IndexSet::add_range
void add_range(const size_type begin, const size_type end)
Definition: index_set.cc:86
numbers::artificial_subdomain_id
const types::subdomain_id artificial_subdomain_id
Definition: types.h:255
IndexSet::compress
void compress() const
Definition: index_set.h:1256
internal::DoFHandler::Policy::Sequential::distribute_mg_dofs
virtual void distribute_mg_dofs(::DoFHandler< dim, spacedim > &dof_handler, std::vector< NumberCache > &number_caches) const
Definition: dof_handler_policy.cc:1004
internal::DoFHandler::NumberCache::n_locally_owned_dofs_per_processor
std::vector< types::global_dof_index > n_locally_owned_dofs_per_processor
Definition: number_cache.h:77
internal::DoFHandler::Policy::ParallelDistributed::distribute_mg_dofs
virtual void distribute_mg_dofs(::DoFHandler< dim, spacedim > &dof_handler, std::vector< NumberCache > &number_caches) const
Definition: dof_handler_policy.cc:2424
Utilities::MPI::this_mpi_process
unsigned int this_mpi_process(const MPI_Comm &mpi_communicator)
Definition: mpi.cc:108
Triangulation::save_user_flags
void save_user_flags(std::ostream &out) const
Definition: tria.cc:9628
DoFHandler::get_triangulation
const Triangulation< dim, spacedim > & get_triangulation() const
IteratorState::valid
Iterator points to a valid object.
Definition: tria_iterator_base.h:39
internal::DoFHandler::NumberCache::locally_owned_dofs_per_processor
std::vector< IndexSet > locally_owned_dofs_per_processor
Definition: number_cache.h:84
IndexSet::is_element
bool is_element(const size_type index) const
Definition: index_set.h:1317
internal::DoFHandler::NumberCache::n_global_dofs
types::global_dof_index n_global_dofs
Definition: number_cache.h:57
DoFRenumbering::subdomain_wise
void subdomain_wise(DoFHandlerType &dof_handler)
Definition: dof_renumbering.cc:1887
parallel::Triangulation::locally_owned_subdomain
types::subdomain_id locally_owned_subdomain() const
Definition: tria_base.cc:218
DoFHandler::locally_owned_dofs
const IndexSet & locally_owned_dofs() const
IndexSet::n_elements
size_type n_elements() const
Definition: index_set.h:1380
Utilities::MPI::max
T max(const T &t, const MPI_Comm &mpi_communicator)
Definition: mpi.h:693
internal::DoFHandler::Policy::ParallelShared::renumber_dofs
virtual void renumber_dofs(const std::vector< types::global_dof_index > &new_numbers, ::DoFHandler< dim, spacedim > &dof_handler, NumberCache &number_cache) const
Definition: dof_handler_policy.cc:1128
DoFHandler::vertex_dofs
std::vector< types::global_dof_index > vertex_dofs
Definition: dof_handler.h:1056
parallel::distributed::Triangulation::coarse_cell_to_p4est_tree_permutation
std::vector< types::global_dof_index > coarse_cell_to_p4est_tree_permutation
Definition: tria.h:879
Generated by   doxygen 1.8.12