1 #include <deal.II/distributed/p4est_wrappers.h> 2 #include <deal.II/distributed/tria.h> 6 #ifdef DEAL_II_WITH_P4EST 14 template <
int dim,
int spacedim>
15 typename ::Triangulation<dim,spacedim>::cell_iterator
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)
21 int i,
l = quad.level;
23 triangulation->get_p4est_tree_to_coarse_cell_permutation()[treeidx];
25 for (i = 0; i <
l; i++)
27 typename ::Triangulation<dim,spacedim>::cell_iterator cell (triangulation, i, dealii_index);
29 Assert (cell->has_children (),
ExcMessage (
"p4est quadrant does not correspond to a cell!"));
30 dealii_index = cell->child_index(child_id);
33 typename ::Triangulation<dim,spacedim>::cell_iterator out_cell (triangulation, l, dealii_index);
42 template <
int dim,
int spacedim>
45 const typename ::parallel::distributed::Triangulation<dim,spacedim> *triangulation;
47 std::map<unsigned int, std::set<::types::subdomain_id> > *vertices_with_ghost_neighbors;
56 template <
int dim,
int spacedim>
59 (typename ::internal::p4est::iter<dim>::corner_info *info,
63 int nsides = info->sides.elem_count;
64 typename ::internal::p4est::iter<dim>::corner_side *sides =
65 (typename ::internal::p4est::iter<dim>::corner_side *)
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;
72 std::map<unsigned int, std::set<::types::subdomain_id> >
73 *vertices_with_ghost_neighbors = fg->vertices_with_ghost_neighbors;
75 subids->elem_count = 0;
76 for (i = 0; i < nsides; i++)
78 if (sides[i].is_ghost)
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"));
84 *subid = cell->subdomain_id();
88 if (!subids->elem_count)
93 nsubs = (int) subids->elem_count;
96 for (i = 0; i < nsides; i++)
98 if (!sides[i].is_ghost)
100 typename ::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].quad));
104 for (j = 0; j < nsubs; j++)
106 (*vertices_with_ghost_neighbors)[cell->vertex_index(sides[i].corner)]
107 .insert (subdomain_ids[j]);
112 subids->elem_count = 0;
118 template <
int dim,
int spacedim>
121 (typename ::internal::p4est::iter<dim>::edge_info *info,
125 int nsides = info->sides.elem_count;
126 typename ::internal::p4est::iter<dim>::edge_side *sides =
127 (typename ::internal::p4est::iter<dim>::edge_side *)
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;
134 std::map<unsigned int, std::set<::types::subdomain_id> >
135 *vertices_with_ghost_neighbors = fg->vertices_with_ghost_neighbors;
137 subids->elem_count = 0;
138 for (i = 0; i < nsides; i++)
140 if (sides[i].is_hanging)
142 for (j = 0; j < 2; j++)
144 if (sides[i].is.hanging.is_ghost[j])
146 typename ::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].is.hanging.quad[j]));
149 *subid = cell->subdomain_id();
155 if (!subids->elem_count)
160 nsubs = (int) subids->elem_count;
163 for (i = 0; i < nsides; i++)
165 if (sides[i].is_hanging)
167 for (j = 0; j < 2; j++)
169 if (!sides[i].is.hanging.is_ghost[j])
171 typename ::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].is.hanging.quad[j]));
173 for (k = 0; k < nsubs; k++)
175 (*vertices_with_ghost_neighbors)[cell->vertex_index(p8est_edge_corners[sides[i].edge][1^j])]
176 .insert (subdomain_ids[k]);
183 subids->elem_count = 0;
189 template <
int dim,
int spacedim>
192 (typename ::internal::p4est::iter<dim>::face_info *info,
196 int nsides = info->sides.elem_count;
197 typename ::internal::p4est::iter<dim>::face_side *sides =
198 (typename ::internal::p4est::iter<dim>::face_side *)
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;
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;
209 subids->elem_count = 0;
210 for (i = 0; i < nsides; i++)
212 if (sides[i].is_hanging)
214 for (j = 0; j < limit; j++)
216 if (sides[i].is.hanging.is_ghost[j])
218 typename ::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].is.hanging.quad[j]));
221 *subid = cell->subdomain_id();
227 if (!subids->elem_count)
232 nsubs = (int) subids->elem_count;
235 for (i = 0; i < nsides; i++)
237 if (sides[i].is_hanging)
239 for (j = 0; j < limit; j++)
241 if (!sides[i].is.hanging.is_ghost[j])
243 typename ::parallel::distributed::Triangulation<dim,spacedim>::cell_iterator cell = cell_from_quad (triangulation, sides[i].treeid, *(sides[i].is.hanging.quad[j]));
245 for (k = 0; k < nsubs; k++)
249 (*vertices_with_ghost_neighbors)[cell->vertex_index(p4est_face_corners[sides[i].face][(limit - 1)^j])]
250 .insert (subdomain_ids[k]);
254 (*vertices_with_ghost_neighbors)[cell->vertex_index(p8est_face_corners[sides[i].face][(limit - 1)^j])]
255 .insert (subdomain_ids[k]);
263 subids->elem_count = 0;
268 int (&functions<2>::quadrant_compare) (
const void *v1,
const void *v2)
269 = p4est_quadrant_compare;
271 void (&functions<2>::quadrant_childrenv) (
const types<2>::quadrant *q,
272 types<2>::quadrant c[])
273 = p4est_quadrant_childrenv;
275 int (&functions<2>::quadrant_overlaps_tree) (types<2>::tree *tree,
276 const types<2>::quadrant *q)
277 = p4est_quadrant_overlaps_tree;
279 void (&functions<2>::quadrant_set_morton) (types<2>::quadrant *quadrant,
282 = p4est_quadrant_set_morton;
284 int (&functions<2>::quadrant_is_equal) (
const types<2>::quadrant *q1,
285 const types<2>::quadrant *q2)
286 = p4est_quadrant_is_equal;
288 int (&functions<2>::quadrant_is_sibling) (
const types<2>::quadrant *q1,
289 const types<2>::quadrant *q2)
290 = p4est_quadrant_is_sibling;
292 int (&functions<2>::quadrant_is_ancestor) (
const types<2>::quadrant *q1,
293 const types<2>::quadrant *q2)
294 = p4est_quadrant_is_ancestor;
296 int (&functions<2>::quadrant_ancestor_id) (
const types<2>::quadrant *q,
298 = p4est_quadrant_ancestor_id;
300 int (&functions<2>::comm_find_owner) (types<2>::forest *p4est,
301 const types<2>::locidx which_tree,
302 const types<2>::quadrant *q,
304 = p4est_comm_find_owner;
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;
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,
319 = p4est_connectivity_join_faces;
322 void (&functions<2>::connectivity_destroy) (p4est_connectivity_t *connectivity)
323 = p4est_connectivity_destroy;
325 types<2>::forest *(&functions<2>::new_forest) (MPI_Comm mpicomm,
326 types<2>::connectivity *connectivity,
327 types<2>::locidx min_quadrants,
331 p4est_init_t init_fn,
335 void (&functions<2>::destroy) (types<2>::forest *p4est)
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)
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)
350 void (&functions<2>::balance) (types<2>::forest *p4est,
351 types<2>::balance_type btype,
352 p4est_init_t init_fn)
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;
361 void (&functions<2>::partition) (types<2>::forest *p4est,
362 int partition_for_coarsening,
363 p4est_weight_t weight_fn)
364 = p4est_partition_ext;
367 void (&functions<2>::save) (
const char *filename,
368 types<2>::forest *p4est,
372 #if DEAL_II_P4EST_VERSION_GTE(0,3,4,3) 374 (&functions<2>::load_ext) (
const char *filename,
376 std::size_t data_size,
381 types<2>::connectivity **p4est)
385 (&functions<2>::load) (
const char *filename,
387 std::size_t data_size,
390 types<2>::connectivity **p4est)
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;
399 void (&functions<2>::connectivity_save) (
const char *filename,
400 types<2>::connectivity *connectivity)
401 = p4est_connectivity_save;
404 int (&functions<2>::connectivity_is_valid) (types<2>::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,
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;
418 types<2>::connectivity *
419 (&functions<2>::connectivity_load) (
const char *filename,
421 = p4est_connectivity_load;
424 unsigned int (&functions<2>::checksum) (types<2>::forest *p4est)
427 void (&functions<2>::vtk_write_file) (types<2>::forest *p4est,
429 const char *baseName)
430 = p4est_vtk_write_file;
432 types<2>::ghost *(&functions<2>::ghost_new) (types<2>::forest *p4est,
433 types<2>::balance_type btype)
436 void (&functions<2>::ghost_destroy) (types<2>::ghost *ghost)
437 = p4est_ghost_destroy;
439 void (&functions<2>::reset_data) (types<2>::forest *p4est,
441 p4est_init_t init_fn,
445 size_t (&functions<2>::forest_memory_used) (types<2>::forest *p4est)
448 size_t (&functions<2>::connectivity_memory_used) (types<2>::connectivity *p4est)
449 = p4est_connectivity_memory_used;
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)
457 std::map<unsigned int, std::set<::types::subdomain_id> > vertices_with_ghost_neighbors;
459 ::internal::p4est::FindGhosts<dim,spacedim> fg;
461 fg.triangulation = &tria;
462 fg.vertices_with_ghost_neighbors = &vertices_with_ghost_neighbors;
469 static_cast<void *>(&fg),
470 nullptr, find_ghosts_face<2,spacedim>, find_ghosts_corner<2,spacedim>);
476 static_cast<void *>(&fg),
477 nullptr, find_ghosts_face<3,3>, find_ghosts_edge<3,3>, find_ghosts_corner<3,3>);
484 sc_array_destroy (fg.subids);
486 return vertices_with_ghost_neighbors;
489 const unsigned int functions<2>::max_level;
493 int (&functions<3>::quadrant_compare) (
const void *v1,
const void *v2)
494 = p8est_quadrant_compare;
496 void (&functions<3>::quadrant_childrenv) (
const types<3>::quadrant *q,
497 types<3>::quadrant c[])
498 = p8est_quadrant_childrenv;
500 int (&functions<3>::quadrant_overlaps_tree) (types<3>::tree *tree,
501 const types<3>::quadrant *q)
502 = p8est_quadrant_overlaps_tree;
504 void (&functions<3>::quadrant_set_morton) (types<3>::quadrant *quadrant,
507 = p8est_quadrant_set_morton;
509 int (&functions<3>::quadrant_is_equal) (
const types<3>::quadrant *q1,
510 const types<3>::quadrant *q2)
511 = p8est_quadrant_is_equal;
513 int (&functions<3>::quadrant_is_sibling) (
const types<3>::quadrant *q1,
514 const types<3>::quadrant *q2)
515 = p8est_quadrant_is_sibling;
517 int (&functions<3>::quadrant_is_ancestor) (
const types<3>::quadrant *q1,
518 const types<3>::quadrant *q2)
519 = p8est_quadrant_is_ancestor;
521 int (&functions<3>::quadrant_ancestor_id) (
const types<3>::quadrant *q,
523 = p8est_quadrant_ancestor_id;
525 int (&functions<3>::comm_find_owner) (types<3>::forest *p4est,
526 const types<3>::locidx which_tree,
527 const types<3>::quadrant *q,
529 = p8est_comm_find_owner;
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;
539 void (&functions<3>::connectivity_destroy) (p8est_connectivity_t *connectivity)
540 = p8est_connectivity_destroy;
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,
549 = p8est_connectivity_join_faces;
552 types<3>::forest *(&functions<3>::new_forest) (MPI_Comm mpicomm,
553 types<3>::connectivity *connectivity,
554 types<3>::locidx min_quadrants,
558 p8est_init_t init_fn,
562 void (&functions<3>::destroy) (types<3>::forest *p8est)
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)
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)
577 void (&functions<3>::balance) (types<3>::forest *p8est,
578 types<3>::balance_type btype,
579 p8est_init_t init_fn)
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;
588 void (&functions<3>::partition) (types<3>::forest *p8est,
589 int partition_for_coarsening,
590 p8est_weight_t weight_fn)
591 = p8est_partition_ext;
594 void (&functions<3>::save) (
const char *filename,
595 types<3>::forest *p4est,
599 #if DEAL_II_P4EST_VERSION_GTE(0,3,4,3) 601 (&functions<3>::load_ext) (
const char *filename,
603 std::size_t data_size,
608 types<3>::connectivity **p4est)
612 (&functions<3>::load) (
const char *filename,
614 std::size_t data_size,
617 types<3>::connectivity **p4est)
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;
626 void (&functions<3>::connectivity_save) (
const char *filename,
627 types<3>::connectivity *connectivity)
628 = p8est_connectivity_save;
631 int (&functions<3>::connectivity_is_valid) (types<3>::connectivity
633 = p8est_connectivity_is_valid;
635 #if DEAL_II_P4EST_VERSION_GTE(1,0,0,0) 636 types<3>::connectivity *
637 (&functions<3>::connectivity_load) (
const char *filename,
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;
646 types<3>::connectivity *
647 (&functions<3>::connectivity_load) (
const char *filename,
649 = p8est_connectivity_load;
652 unsigned int (&functions<3>::checksum) (types<3>::forest *p8est)
655 void (&functions<3>::vtk_write_file) (types<3>::forest *p8est,
657 const char *baseName)
658 = p8est_vtk_write_file;
660 types<3>::ghost *(&functions<3>::ghost_new) (types<3>::forest *p4est,
661 types<3>::balance_type btype)
664 void (&functions<3>::ghost_destroy) (types<3>::ghost *ghost)
665 = p8est_ghost_destroy;
667 void (&functions<3>::reset_data) (types<3>::forest *p4est,
669 p8est_init_t init_fn,
673 size_t (&functions<3>::forest_memory_used) (types<3>::forest *p4est)
676 size_t (&functions<3>::connectivity_memory_used) (types<3>::connectivity *p4est)
677 = p8est_connectivity_memory_used;
679 const unsigned int functions<3>::max_level;
685 init_quadrant_children
686 (
const typename types<dim>::quadrant &p4est_cell,
690 for (
unsigned int c=0;
691 c<::GeometryInfo<dim>::max_children_per_cell; ++c)
695 P4EST_QUADRANT_INIT(&p4est_children[c]);
698 P8EST_QUADRANT_INIT(&p4est_children[c]);
705 functions<dim>::quadrant_childrenv (&p4est_cell,
712 init_coarse_quadrant(
typename types<dim>::quadrant &quad)
717 P4EST_QUADRANT_INIT(&quad);
720 P8EST_QUADRANT_INIT(&quad);
725 functions<dim>::quadrant_set_morton (&quad,
732 quadrant_is_equal (
const typename types<dim>::quadrant &q1,
733 const typename types<dim>::quadrant &q2)
735 return functions<dim>::quadrant_is_equal(&q1, &q2);
742 quadrant_is_ancestor (
const typename types<dim>::quadrant &q1,
743 const typename types<dim>::quadrant &q2)
745 return functions<dim>::quadrant_is_ancestor(&q1, &q2);
750 tree_exists_locally (
const typename types<dim>::forest *parallel_forest,
751 const typename types<dim>::topidx coarse_grid_cell)
753 Assert (coarse_grid_cell < parallel_forest->connectivity->num_trees,
755 return ((coarse_grid_cell >= parallel_forest->first_local_tree)
757 (coarse_grid_cell <= parallel_forest->last_local_tree));
762 #endif //DEAL_II_WITH_P4EST 765 #include "p4est_wrappers.inst" 768 DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcMessage(std::string arg1)
unsigned int global_dof_index
#define Assert(cond, exc)
unsigned int subdomain_id
static ::ExceptionBase & ExcNotImplemented()
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
static ::ExceptionBase & ExcInternalError()