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