Reference documentation for deal.II version 9.0.0
p4est_wrappers.cc
1 #include <deal.II/distributed/p4est_wrappers.h>
2 #include <deal.II/distributed/tria.h>
3 
4 DEAL_II_NAMESPACE_OPEN
5 
6 #ifdef DEAL_II_WITH_P4EST
7 
8 namespace internal
9 {
10  namespace p4est
11  {
12  namespace
13  {
14  template <int dim, int spacedim>
15  typename ::Triangulation<dim,spacedim>::cell_iterator
16  cell_from_quad
17  (const ::parallel::distributed::Triangulation<dim,spacedim> *triangulation,
18  const typename ::internal::p4est::types<dim>::topidx treeidx,
19  const typename ::internal::p4est::types<dim>::quadrant &quad)
20  {
21  int i, l = quad.level;
22  ::types::global_dof_index dealii_index =
23  triangulation->get_p4est_tree_to_coarse_cell_permutation()[treeidx];
24 
25  for (i = 0; i < l; i++)
26  {
27  typename ::Triangulation<dim,spacedim>::cell_iterator cell (triangulation, i, dealii_index);
28  const int child_id = ::internal::p4est::functions<dim>::quadrant_ancestor_id (&quad, i + 1);
29  Assert (cell->has_children (), ExcMessage ("p4est quadrant does not correspond to a cell!"));
30  dealii_index = cell->child_index(child_id);
31  }
32 
33  typename ::Triangulation<dim,spacedim>::cell_iterator out_cell (triangulation, l, dealii_index);
34 
35  return out_cell;
36  }
37 
42  template <int dim, int spacedim>
43  struct FindGhosts
44  {
45  const typename ::parallel::distributed::Triangulation<dim,spacedim> *triangulation;
46  sc_array_t *subids;
47  std::map<unsigned int, std::set<::types::subdomain_id> > *vertices_with_ghost_neighbors;
48  };
49 
50 
56  template <int dim, int spacedim>
57  void
58  find_ghosts_corner
59  (typename ::internal::p4est::iter<dim>::corner_info *info,
60  void *user_data)
61  {
62  int i, j;
63  int nsides = info->sides.elem_count;
64  typename ::internal::p4est::iter<dim>::corner_side *sides =
65  (typename ::internal::p4est::iter<dim>::corner_side *)
66  (info->sides.array);
67  FindGhosts<dim,spacedim> *fg = static_cast<FindGhosts<dim,spacedim> *>(user_data);
68  sc_array_t *subids = fg->subids;
69  const ::parallel::distributed::Triangulation<dim,spacedim> *triangulation = fg->triangulation;
70  int nsubs;
71  ::types::subdomain_id *subdomain_ids;
72  std::map<unsigned int, std::set<::types::subdomain_id> >
73  *vertices_with_ghost_neighbors = fg->vertices_with_ghost_neighbors;
74 
75  subids->elem_count = 0;
76  for (i = 0; i < nsides; i++)
77  {
78  if (sides[i].is_ghost)
79  {
80  typename ::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].quad));
81  Assert (cell->is_ghost(), ExcMessage ("ghost quad did not find ghost cell"));
82  ::types::subdomain_id *subid =
83  static_cast<::types::subdomain_id *>(sc_array_push (subids));
84  *subid = cell->subdomain_id();
85  }
86  }
87 
88  if (!subids->elem_count)
89  {
90  return;
91  }
92 
93  nsubs = (int) subids->elem_count;
94  subdomain_ids = (::types::subdomain_id *) (subids->array);
95 
96  for (i = 0; i < nsides; i++)
97  {
98  if (!sides[i].is_ghost)
99  {
100  typename ::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].quad));
101 
102  Assert (!cell->is_ghost(), ExcMessage ("local quad found ghost cell"));
103 
104  for (j = 0; j < nsubs; j++)
105  {
106  (*vertices_with_ghost_neighbors)[cell->vertex_index(sides[i].corner)]
107  .insert (subdomain_ids[j]);
108  }
109  }
110  }
111 
112  subids->elem_count = 0;
113  }
114 
118  template <int dim, int spacedim>
119  void
120  find_ghosts_edge
121  (typename ::internal::p4est::iter<dim>::edge_info *info,
122  void *user_data)
123  {
124  int i, j, k;
125  int nsides = info->sides.elem_count;
126  typename ::internal::p4est::iter<dim>::edge_side *sides =
127  (typename ::internal::p4est::iter<dim>::edge_side *)
128  (info->sides.array);
129  FindGhosts<dim,spacedim> *fg = static_cast<FindGhosts<dim,spacedim> *>(user_data);
130  sc_array_t *subids = fg->subids;
131  const ::parallel::distributed::Triangulation<dim,spacedim> *triangulation = fg->triangulation;
132  int nsubs;
133  ::types::subdomain_id *subdomain_ids;
134  std::map<unsigned int, std::set<::types::subdomain_id> >
135  *vertices_with_ghost_neighbors = fg->vertices_with_ghost_neighbors;
136 
137  subids->elem_count = 0;
138  for (i = 0; i < nsides; i++)
139  {
140  if (sides[i].is_hanging)
141  {
142  for (j = 0; j < 2; j++)
143  {
144  if (sides[i].is.hanging.is_ghost[j])
145  {
146  typename ::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].is.hanging.quad[j]));
147  ::types::subdomain_id *subid =
148  static_cast<::types::subdomain_id *>(sc_array_push (subids));
149  *subid = cell->subdomain_id();
150  }
151  }
152  }
153  }
154 
155  if (!subids->elem_count)
156  {
157  return;
158  }
159 
160  nsubs = (int) subids->elem_count;
161  subdomain_ids = (::types::subdomain_id *) (subids->array);
162 
163  for (i = 0; i < nsides; i++)
164  {
165  if (sides[i].is_hanging)
166  {
167  for (j = 0; j < 2; j++)
168  {
169  if (!sides[i].is.hanging.is_ghost[j])
170  {
171  typename ::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].is.hanging.quad[j]));
172 
173  for (k = 0; k < nsubs; k++)
174  {
175  (*vertices_with_ghost_neighbors)[cell->vertex_index(p8est_edge_corners[sides[i].edge][1^j])]
176  .insert (subdomain_ids[k]);
177  }
178  }
179  }
180  }
181  }
182 
183  subids->elem_count = 0;
184  }
185 
189  template <int dim, int spacedim>
190  void
191  find_ghosts_face
192  (typename ::internal::p4est::iter<dim>::face_info *info,
193  void *user_data)
194  {
195  int i, j, k;
196  int nsides = info->sides.elem_count;
197  typename ::internal::p4est::iter<dim>::face_side *sides =
198  (typename ::internal::p4est::iter<dim>::face_side *)
199  (info->sides.array);
200  FindGhosts<dim,spacedim> *fg = static_cast<FindGhosts<dim,spacedim> *>(user_data);
201  sc_array_t *subids = fg->subids;
202  const ::parallel::distributed::Triangulation<dim,spacedim> *triangulation = fg->triangulation;
203  int nsubs;
204  ::types::subdomain_id *subdomain_ids;
205  std::map<unsigned int, std::set<::types::subdomain_id> >
206  *vertices_with_ghost_neighbors = fg->vertices_with_ghost_neighbors;
207  int limit = (dim == 2) ? 2 : 4;
208 
209  subids->elem_count = 0;
210  for (i = 0; i < nsides; i++)
211  {
212  if (sides[i].is_hanging)
213  {
214  for (j = 0; j < limit; j++)
215  {
216  if (sides[i].is.hanging.is_ghost[j])
217  {
218  typename ::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].is.hanging.quad[j]));
219  ::types::subdomain_id *subid =
220  static_cast<::types::subdomain_id *>(sc_array_push (subids));
221  *subid = cell->subdomain_id();
222  }
223  }
224  }
225  }
226 
227  if (!subids->elem_count)
228  {
229  return;
230  }
231 
232  nsubs = (int) subids->elem_count;
233  subdomain_ids = (::types::subdomain_id *) (subids->array);
234 
235  for (i = 0; i < nsides; i++)
236  {
237  if (sides[i].is_hanging)
238  {
239  for (j = 0; j < limit; j++)
240  {
241  if (!sides[i].is.hanging.is_ghost[j])
242  {
243  typename ::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].is.hanging.quad[j]));
244 
245  for (k = 0; k < nsubs; k++)
246  {
247  if (dim == 2)
248  {
249  (*vertices_with_ghost_neighbors)[cell->vertex_index(p4est_face_corners[sides[i].face][(limit - 1)^j])]
250  .insert (subdomain_ids[k]);
251  }
252  else
253  {
254  (*vertices_with_ghost_neighbors)[cell->vertex_index(p8est_face_corners[sides[i].face][(limit - 1)^j])]
255  .insert (subdomain_ids[k]);
256  }
257  }
258  }
259  }
260  }
261  }
262 
263  subids->elem_count = 0;
264  }
265  }
266 
267 
268  int (&functions<2>::quadrant_compare) (const void *v1, const void *v2)
269  = p4est_quadrant_compare;
270 
271  void (&functions<2>::quadrant_childrenv) (const types<2>::quadrant *q,
272  types<2>::quadrant c[])
273  = p4est_quadrant_childrenv;
274 
275  int (&functions<2>::quadrant_overlaps_tree) (types<2>::tree *tree,
276  const types<2>::quadrant *q)
277  = p4est_quadrant_overlaps_tree;
278 
279  void (&functions<2>::quadrant_set_morton) (types<2>::quadrant *quadrant,
280  int level,
281  uint64_t id)
282  = p4est_quadrant_set_morton;
283 
284  int (&functions<2>::quadrant_is_equal) (const types<2>::quadrant *q1,
285  const types<2>::quadrant *q2)
286  = p4est_quadrant_is_equal;
287 
288  int (&functions<2>::quadrant_is_sibling) (const types<2>::quadrant *q1,
289  const types<2>::quadrant *q2)
290  = p4est_quadrant_is_sibling;
291 
292  int (&functions<2>::quadrant_is_ancestor) (const types<2>::quadrant *q1,
293  const types<2>::quadrant *q2)
294  = p4est_quadrant_is_ancestor;
295 
296  int (&functions<2>::quadrant_ancestor_id) (const types<2>::quadrant *q,
297  int level)
298  = p4est_quadrant_ancestor_id;
299 
300  int (&functions<2>::comm_find_owner) (types<2>::forest *p4est,
301  const types<2>::locidx which_tree,
302  const types<2>::quadrant *q,
303  const int guess)
304  = p4est_comm_find_owner;
305 
306  types<2>::connectivity *(&functions<2>::connectivity_new) (types<2>::topidx num_vertices,
307  types<2>::topidx num_trees,
308  types<2>::topidx num_corners,
309  types<2>::topidx num_vtt)
310  = p4est_connectivity_new;
311 
312 #if DEAL_II_P4EST_VERSION_GTE(0,3,4,1)
313  void (&functions<2>::connectivity_join_faces) (types<2>::connectivity *conn,
314  types<2>::topidx tree_left,
315  types<2>::topidx tree_right,
316  int face_left,
317  int face_right,
318  int orientation)
319  = p4est_connectivity_join_faces;
320 #endif
321 
322  void (&functions<2>::connectivity_destroy) (p4est_connectivity_t *connectivity)
323  = p4est_connectivity_destroy;
324 
325  types<2>::forest *(&functions<2>::new_forest) (MPI_Comm mpicomm,
326  types<2>::connectivity *connectivity,
327  types<2>::locidx min_quadrants,
328  int min_level,
329  int fill_uniform,
330  size_t data_size,
331  p4est_init_t init_fn,
332  void *user_pointer)
333  = p4est_new_ext;
334 
335  void (&functions<2>::destroy) (types<2>::forest *p4est)
336  = p4est_destroy;
337 
338  void (&functions<2>::refine) (types<2>::forest *p4est,
339  int refine_recursive,
340  p4est_refine_t refine_fn,
341  p4est_init_t init_fn)
342  = p4est_refine;
343 
344  void (&functions<2>::coarsen) (types<2>::forest *p4est,
345  int coarsen_recursive,
346  p4est_coarsen_t coarsen_fn,
347  p4est_init_t init_fn)
348  = p4est_coarsen;
349 
350  void (&functions<2>::balance) (types<2>::forest *p4est,
351  types<2>::balance_type btype,
352  p4est_init_t init_fn)
353  = p4est_balance;
354 #if DEAL_II_P4EST_VERSION_GTE(0,3,4,3)
355  p4est_gloidx_t (&functions<2>::partition) (types<2>::forest *p4est,
356  int partition_for_coarsening,
357  p4est_weight_t weight_fn)
358  = p4est_partition_ext;
359 #else
360 
361  void (&functions<2>::partition) (types<2>::forest *p4est,
362  int partition_for_coarsening,
363  p4est_weight_t weight_fn)
364  = p4est_partition_ext;
365 #endif
366 
367  void (&functions<2>::save) (const char *filename,
368  types<2>::forest *p4est,
369  int save_data)
370  = p4est_save;
371 
372 #if DEAL_II_P4EST_VERSION_GTE(0,3,4,3)
373  types<2>::forest *
374  (&functions<2>::load_ext) (const char *filename,
375  MPI_Comm mpicomm,
376  std::size_t data_size,
377  int load_data,
378  int autopartition,
379  int broadcasthead,
380  void *user_pointer,
381  types<2>::connectivity **p4est)
382  = p4est_load_ext;
383 #else
384  types<2>::forest *
385  (&functions<2>::load) (const char *filename,
386  MPI_Comm mpicomm,
387  std::size_t data_size,
388  int load_data,
389  void *user_pointer,
390  types<2>::connectivity **p4est)
391  = p4est_load;
392 #endif
393 
394 #if DEAL_II_P4EST_VERSION_GTE(0,3,4,3)
395  int (&functions<2>::connectivity_save) (const char *filename,
396  types<2>::connectivity *connectivity)
397  = p4est_connectivity_save;
398 #else
399  void (&functions<2>::connectivity_save) (const char *filename,
400  types<2>::connectivity *connectivity)
401  = p4est_connectivity_save;
402 #endif
403 
404  int (&functions<2>::connectivity_is_valid) (types<2>::connectivity
405  *connectivity)
406  = p4est_connectivity_is_valid;
407 #if DEAL_II_P4EST_VERSION_GTE(1,0,0,0)
408  types<2>::connectivity *
409  (&functions<2>::connectivity_load) (const char *filename,
410  size_t *length)
411  = p4est_connectivity_load;
412 #elif DEAL_II_P4EST_VERSION_GTE(0,3,4,3)
413  types<2>::connectivity *
414  (&functions<2>::connectivity_load) (const char *filename,
415  long unsigned *length)
416  = p4est_connectivity_load;
417 #else
418  types<2>::connectivity *
419  (&functions<2>::connectivity_load) (const char *filename,
420  long *length)
421  = p4est_connectivity_load;
422 #endif
423 
424  unsigned int (&functions<2>::checksum) (types<2>::forest *p4est)
425  = p4est_checksum;
426 
427  void (&functions<2>::vtk_write_file) (types<2>::forest *p4est,
428  p4est_geometry_t *,
429  const char *baseName)
430  = p4est_vtk_write_file;
431 
432  types<2>::ghost *(&functions<2>::ghost_new) (types<2>::forest *p4est,
433  types<2>::balance_type btype)
434  = p4est_ghost_new;
435 
436  void (&functions<2>::ghost_destroy) (types<2>::ghost *ghost)
437  = p4est_ghost_destroy;
438 
439  void (&functions<2>::reset_data) (types<2>::forest *p4est,
440  size_t data_size,
441  p4est_init_t init_fn,
442  void *user_pointer)
443  = p4est_reset_data;
444 
445  size_t (&functions<2>::forest_memory_used) (types<2>::forest *p4est)
446  = p4est_memory_used;
447 
448  size_t (&functions<2>::connectivity_memory_used) (types<2>::connectivity *p4est)
449  = p4est_connectivity_memory_used;
450 
451  template <int dim, int spacedim>
452  std::map<unsigned int, std::set<::types::subdomain_id> >
453  compute_vertices_with_ghost_neighbors(const typename ::parallel::distributed::Triangulation<dim,spacedim> &tria,
454  typename ::internal::p4est::types<dim>::forest *parallel_forest,
455  typename ::internal::p4est::types<dim>::ghost *parallel_ghost)
456  {
457  std::map<unsigned int, std::set<::types::subdomain_id> > vertices_with_ghost_neighbors;
458 
459  ::internal::p4est::FindGhosts<dim,spacedim> fg;
460  fg.subids = sc_array_new (sizeof (::types::subdomain_id));
461  fg.triangulation = &tria;
462  fg.vertices_with_ghost_neighbors = &vertices_with_ghost_neighbors;
463 
464  switch (dim)
465  {
466  case 2:
467  p4est_iterate (reinterpret_cast<::internal::p4est::types<2>::forest *>(parallel_forest),
468  reinterpret_cast<::internal::p4est::types<2>::ghost *>(parallel_ghost),
469  static_cast<void *>(&fg),
470  nullptr, find_ghosts_face<2,spacedim>, find_ghosts_corner<2,spacedim>);
471  break;
472 
473  case 3:
474  p8est_iterate (reinterpret_cast<::internal::p4est::types<3>::forest *>(parallel_forest),
475  reinterpret_cast<::internal::p4est::types<3>::ghost *>(parallel_ghost),
476  static_cast<void *>(&fg),
477  nullptr, find_ghosts_face<3,3>, find_ghosts_edge<3,3>, find_ghosts_corner<3,3>);
478  break;
479 
480  default:
481  Assert (false, ExcNotImplemented());
482  }
483 
484  sc_array_destroy (fg.subids);
485 
486  return vertices_with_ghost_neighbors;
487  }
488 
489  const unsigned int functions<2>::max_level;
490 
491 
492 
493  int (&functions<3>::quadrant_compare) (const void *v1, const void *v2)
494  = p8est_quadrant_compare;
495 
496  void (&functions<3>::quadrant_childrenv) (const types<3>::quadrant *q,
497  types<3>::quadrant c[])
498  = p8est_quadrant_childrenv;
499 
500  int (&functions<3>::quadrant_overlaps_tree) (types<3>::tree *tree,
501  const types<3>::quadrant *q)
502  = p8est_quadrant_overlaps_tree;
503 
504  void (&functions<3>::quadrant_set_morton) (types<3>::quadrant *quadrant,
505  int level,
506  uint64_t id)
507  = p8est_quadrant_set_morton;
508 
509  int (&functions<3>::quadrant_is_equal) (const types<3>::quadrant *q1,
510  const types<3>::quadrant *q2)
511  = p8est_quadrant_is_equal;
512 
513  int (&functions<3>::quadrant_is_sibling) (const types<3>::quadrant *q1,
514  const types<3>::quadrant *q2)
515  = p8est_quadrant_is_sibling;
516 
517  int (&functions<3>::quadrant_is_ancestor) (const types<3>::quadrant *q1,
518  const types<3>::quadrant *q2)
519  = p8est_quadrant_is_ancestor;
520 
521  int (&functions<3>::quadrant_ancestor_id) (const types<3>::quadrant *q,
522  int level)
523  = p8est_quadrant_ancestor_id;
524 
525  int (&functions<3>::comm_find_owner) (types<3>::forest *p4est,
526  const types<3>::locidx which_tree,
527  const types<3>::quadrant *q,
528  const int guess)
529  = p8est_comm_find_owner;
530 
531  types<3>::connectivity *(&functions<3>::connectivity_new) (types<3>::topidx num_vertices,
532  types<3>::topidx num_trees,
533  types<3>::topidx num_edges,
534  types<3>::topidx num_ett,
535  types<3>::topidx num_corners,
536  types<3>::topidx num_ctt)
537  = p8est_connectivity_new;
538 
539  void (&functions<3>::connectivity_destroy) (p8est_connectivity_t *connectivity)
540  = p8est_connectivity_destroy;
541 
542 #if DEAL_II_P4EST_VERSION_GTE(0,3,4,1)
543  void (&functions<3>::connectivity_join_faces) (types<3>::connectivity *conn,
544  types<3>::topidx tree_left,
545  types<3>::topidx tree_right,
546  int face_left,
547  int face_right,
548  int orientation)
549  = p8est_connectivity_join_faces;
550 #endif
551 
552  types<3>::forest *(&functions<3>::new_forest) (MPI_Comm mpicomm,
553  types<3>::connectivity *connectivity,
554  types<3>::locidx min_quadrants,
555  int min_level,
556  int fill_uniform,
557  size_t data_size,
558  p8est_init_t init_fn,
559  void *user_pointer)
560  = p8est_new_ext;
561 
562  void (&functions<3>::destroy) (types<3>::forest *p8est)
563  = p8est_destroy;
564 
565  void (&functions<3>::refine) (types<3>::forest *p8est,
566  int refine_recursive,
567  p8est_refine_t refine_fn,
568  p8est_init_t init_fn)
569  = p8est_refine;
570 
571  void (&functions<3>::coarsen) (types<3>::forest *p8est,
572  int coarsen_recursive,
573  p8est_coarsen_t coarsen_fn,
574  p8est_init_t init_fn)
575  = p8est_coarsen;
576 
577  void (&functions<3>::balance) (types<3>::forest *p8est,
578  types<3>::balance_type btype,
579  p8est_init_t init_fn)
580  = p8est_balance;
581 
582 #if DEAL_II_P4EST_VERSION_GTE(0,3,4,3)
583  p4est_gloidx_t (&functions<3>::partition) (types<3>::forest *p8est,
584  int partition_for_coarsening,
585  p8est_weight_t weight_fn)
586  = p8est_partition_ext;
587 #else
588  void (&functions<3>::partition) (types<3>::forest *p8est,
589  int partition_for_coarsening,
590  p8est_weight_t weight_fn)
591  = p8est_partition_ext;
592 #endif
593 
594  void (&functions<3>::save) (const char *filename,
595  types<3>::forest *p4est,
596  int save_data)
597  = p8est_save;
598 
599 #if DEAL_II_P4EST_VERSION_GTE(0,3,4,3)
600  types<3>::forest *
601  (&functions<3>::load_ext) (const char *filename,
602  MPI_Comm mpicomm,
603  std::size_t data_size,
604  int load_data,
605  int autopartition,
606  int broadcasthead,
607  void *user_pointer,
608  types<3>::connectivity **p4est)
609  = p8est_load_ext;
610 #else
611  types<3>::forest *
612  (&functions<3>::load) (const char *filename,
613  MPI_Comm mpicomm,
614  std::size_t data_size,
615  int load_data,
616  void *user_pointer,
617  types<3>::connectivity **p4est)
618  = p8est_load;
619 #endif
620 
621 #if DEAL_II_P4EST_VERSION_GTE(0,3,4,3)
622  int (&functions<3>::connectivity_save) (const char *filename,
623  types<3>::connectivity *connectivity)
624  = p8est_connectivity_save;
625 #else
626  void (&functions<3>::connectivity_save) (const char *filename,
627  types<3>::connectivity *connectivity)
628  = p8est_connectivity_save;
629 #endif
630 
631  int (&functions<3>::connectivity_is_valid) (types<3>::connectivity
632  *connectivity)
633  = p8est_connectivity_is_valid;
634 
635 #if DEAL_II_P4EST_VERSION_GTE(1,0,0,0)
636  types<3>::connectivity *
637  (&functions<3>::connectivity_load) (const char *filename,
638  size_t *length)
639  = p8est_connectivity_load;
640 #elif DEAL_II_P4EST_VERSION_GTE(0,3,4,3)
641  types<3>::connectivity *
642  (&functions<3>::connectivity_load) (const char *filename,
643  long unsigned *length)
644  = p8est_connectivity_load;
645 #else
646  types<3>::connectivity *
647  (&functions<3>::connectivity_load) (const char *filename,
648  long *length)
649  = p8est_connectivity_load;
650 #endif
651 
652  unsigned int (&functions<3>::checksum) (types<3>::forest *p8est)
653  = p8est_checksum;
654 
655  void (&functions<3>::vtk_write_file) (types<3>::forest *p8est,
656  p8est_geometry_t *,
657  const char *baseName)
658  = p8est_vtk_write_file;
659 
660  types<3>::ghost *(&functions<3>::ghost_new) (types<3>::forest *p4est,
661  types<3>::balance_type btype)
662  = p8est_ghost_new;
663 
664  void (&functions<3>::ghost_destroy) (types<3>::ghost *ghost)
665  = p8est_ghost_destroy;
666 
667  void (&functions<3>::reset_data) (types<3>::forest *p4est,
668  size_t data_size,
669  p8est_init_t init_fn,
670  void *user_pointer)
671  = p8est_reset_data;
672 
673  size_t (&functions<3>::forest_memory_used) (types<3>::forest *p4est)
674  = p8est_memory_used;
675 
676  size_t (&functions<3>::connectivity_memory_used) (types<3>::connectivity *p4est)
677  = p8est_connectivity_memory_used;
678 
679  const unsigned int functions<3>::max_level;
680 
681 
682 
683  template <int dim>
684  void
685  init_quadrant_children
686  (const typename types<dim>::quadrant &p4est_cell,
687  typename types<dim>::quadrant (&p4est_children)[::GeometryInfo<dim>::max_children_per_cell])
688  {
689 
690  for (unsigned int c=0;
691  c<::GeometryInfo<dim>::max_children_per_cell; ++c)
692  switch (dim)
693  {
694  case 2:
695  P4EST_QUADRANT_INIT(&p4est_children[c]);
696  break;
697  case 3:
698  P8EST_QUADRANT_INIT(&p4est_children[c]);
699  break;
700  default:
701  Assert (false, ExcNotImplemented());
702  }
703 
704 
705  functions<dim>::quadrant_childrenv (&p4est_cell,
706  p4est_children);
707 
708  }
709 
710  template <int dim>
711  void
712  init_coarse_quadrant(typename types<dim>::quadrant &quad)
713  {
714  switch (dim)
715  {
716  case 2:
717  P4EST_QUADRANT_INIT(&quad);
718  break;
719  case 3:
720  P8EST_QUADRANT_INIT(&quad);
721  break;
722  default:
723  Assert (false, ExcNotImplemented());
724  }
725  functions<dim>::quadrant_set_morton (&quad,
726  /*level=*/0,
727  /*index=*/0);
728  }
729 
730  template <int dim>
731  bool
732  quadrant_is_equal (const typename types<dim>::quadrant &q1,
733  const typename types<dim>::quadrant &q2)
734  {
735  return functions<dim>::quadrant_is_equal(&q1, &q2);
736  }
737 
738 
739 
740  template <int dim>
741  bool
742  quadrant_is_ancestor (const typename types<dim>::quadrant &q1,
743  const typename types<dim>::quadrant &q2)
744  {
745  return functions<dim>::quadrant_is_ancestor(&q1, &q2);
746  }
747 
748  template <int dim>
749  bool
750  tree_exists_locally (const typename types<dim>::forest *parallel_forest,
751  const typename types<dim>::topidx coarse_grid_cell)
752  {
753  Assert (coarse_grid_cell < parallel_forest->connectivity->num_trees,
754  ExcInternalError());
755  return ((coarse_grid_cell >= parallel_forest->first_local_tree)
756  &&
757  (coarse_grid_cell <= parallel_forest->last_local_tree));
758  }
759  }
760 }
761 
762 #endif //DEAL_II_WITH_P4EST
763 
764 /*-------------- Explicit Instantiations -------------------------------*/
765 #include "p4est_wrappers.inst"
766 
767 
768 DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcMessage(std::string arg1)
unsigned int global_dof_index
Definition: types.h:88
#define Assert(cond, exc)
Definition: exceptions.h:1142
unsigned int subdomain_id
Definition: types.h:42
static ::ExceptionBase & ExcNotImplemented()
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
static ::ExceptionBase & ExcInternalError()