17 #include <deal.II/distributed/p4est_wrappers.h> 18 #include <deal.II/distributed/tria.h> 20 DEAL_II_NAMESPACE_OPEN
22 #ifdef DEAL_II_WITH_P4EST 30 template <
int dim,
int spacedim>
31 typename ::Triangulation<dim, spacedim>::cell_iterator
33 const ::parallel::distributed::Triangulation<dim, spacedim>
35 const typename ::internal::p4est::types<dim>::topidx treeidx,
36 const typename ::internal::p4est::types<dim>::quadrant &quad)
38 int i,
l = quad.level;
40 triangulation->get_p4est_tree_to_coarse_cell_permutation()[treeidx];
42 for (i = 0; i <
l; i++)
44 typename ::Triangulation<dim, spacedim>::cell_iterator cell(
45 triangulation, i, dealii_index);
49 Assert(cell->has_children(),
50 ExcMessage(
"p4est quadrant does not correspond to a cell!"));
51 dealii_index = cell->child_index(child_id);
54 typename ::Triangulation<dim, spacedim>::cell_iterator out_cell(
55 triangulation, l, dealii_index);
64 template <
int dim,
int spacedim>
67 const typename ::parallel::distributed::Triangulation<dim,
71 std::map<unsigned int, std::set<::types::subdomain_id>>
72 *vertices_with_ghost_neighbors;
81 template <
int dim,
int spacedim>
84 typename ::internal::p4est::iter<dim>::corner_info *info,
88 int nsides = info->sides.elem_count;
89 auto *sides =
reinterpret_cast< 90 typename ::internal::p4est::iter<dim>::corner_side *
>(
92 FindGhosts<dim, spacedim> *fg =
93 static_cast<FindGhosts<dim, spacedim> *
>(user_data);
94 sc_array_t *subids = fg->subids;
95 const ::parallel::distributed::Triangulation<dim, spacedim>
96 * triangulation = fg->triangulation;
99 std::map<unsigned int, std::set<::types::subdomain_id>>
100 *vertices_with_ghost_neighbors = fg->vertices_with_ghost_neighbors;
102 subids->elem_count = 0;
103 for (i = 0; i < nsides; i++)
105 if (sides[i].is_ghost)
107 typename ::parallel::distributed::
108 Triangulation<dim, spacedim>::cell_iterator cell =
109 cell_from_quad(triangulation,
113 ExcMessage(
"ghost quad did not find ghost cell"));
116 sc_array_push(subids));
117 *subid = cell->subdomain_id();
121 if (!subids->elem_count)
126 nsubs =
static_cast<int>(subids->elem_count);
130 for (i = 0; i < nsides; i++)
132 if (!sides[i].is_ghost)
134 typename ::parallel::distributed::
135 Triangulation<dim, spacedim>::cell_iterator cell =
136 cell_from_quad(triangulation,
143 for (j = 0; j < nsubs; j++)
145 (*vertices_with_ghost_neighbors)[cell->vertex_index(
147 .insert(subdomain_ids[j]);
152 subids->elem_count = 0;
158 template <
int dim,
int spacedim>
161 typename ::internal::p4est::iter<dim>::edge_info *info,
165 int nsides = info->sides.elem_count;
166 auto *sides =
reinterpret_cast< 167 typename ::internal::p4est::iter<dim>::edge_side *
>(
169 auto * fg =
static_cast<FindGhosts<dim, spacedim> *
>(user_data);
170 sc_array_t *subids = fg->subids;
171 const ::parallel::distributed::Triangulation<dim, spacedim>
172 * triangulation = fg->triangulation;
175 std::map<unsigned int, std::set<::types::subdomain_id>>
176 *vertices_with_ghost_neighbors = fg->vertices_with_ghost_neighbors;
178 subids->elem_count = 0;
179 for (i = 0; i < nsides; i++)
181 if (sides[i].is_hanging)
183 for (j = 0; j < 2; j++)
185 if (sides[i].is.hanging.is_ghost[j])
187 typename ::parallel::distributed::
188 Triangulation<dim, spacedim>::cell_iterator cell =
189 cell_from_quad(triangulation,
191 *(sides[i].is.hanging.quad[j]));
194 sc_array_push(subids));
195 *subid = cell->subdomain_id();
201 if (!subids->elem_count)
206 nsubs =
static_cast<int>(subids->elem_count);
210 for (i = 0; i < nsides; i++)
212 if (sides[i].is_hanging)
214 for (j = 0; j < 2; j++)
216 if (!sides[i].is.hanging.is_ghost[j])
218 typename ::parallel::distributed::
219 Triangulation<dim, spacedim>::cell_iterator cell =
220 cell_from_quad(triangulation,
222 *(sides[i].is.hanging.quad[j]));
224 for (k = 0; k < nsubs; k++)
226 (*vertices_with_ghost_neighbors)
228 p8est_edge_corners[sides[i].edge][1 ^ j])]
229 .insert(subdomain_ids[k]);
236 subids->elem_count = 0;
242 template <
int dim,
int spacedim>
245 typename ::internal::p4est::iter<dim>::face_info *info,
249 int nsides = info->sides.elem_count;
250 auto *sides =
reinterpret_cast< 251 typename ::internal::p4est::iter<dim>::face_side *
>(
253 FindGhosts<dim, spacedim> *fg =
254 static_cast<FindGhosts<dim, spacedim> *
>(user_data);
255 sc_array_t *subids = fg->subids;
256 const ::parallel::distributed::Triangulation<dim, spacedim>
257 * triangulation = fg->triangulation;
260 std::map<unsigned int, std::set<::types::subdomain_id>>
261 * vertices_with_ghost_neighbors = fg->vertices_with_ghost_neighbors;
262 int limit = (dim == 2) ? 2 : 4;
264 subids->elem_count = 0;
265 for (i = 0; i < nsides; i++)
267 if (sides[i].is_hanging)
269 for (j = 0; j < limit; j++)
271 if (sides[i].is.hanging.is_ghost[j])
273 typename ::parallel::distributed::
274 Triangulation<dim, spacedim>::cell_iterator cell =
275 cell_from_quad(triangulation,
277 *(sides[i].is.hanging.quad[j]));
280 sc_array_push(subids));
281 *subid = cell->subdomain_id();
287 if (!subids->elem_count)
292 nsubs =
static_cast<int>(subids->elem_count);
296 for (i = 0; i < nsides; i++)
298 if (sides[i].is_hanging)
300 for (j = 0; j < limit; j++)
302 if (!sides[i].is.hanging.is_ghost[j])
304 typename ::parallel::distributed::
305 Triangulation<dim, spacedim>::cell_iterator cell =
306 cell_from_quad(triangulation,
308 *(sides[i].is.hanging.quad[j]));
310 for (k = 0; k < nsubs; k++)
314 (*vertices_with_ghost_neighbors)
316 p4est_face_corners[sides[i].face]
318 .insert(subdomain_ids[k]);
322 (*vertices_with_ghost_neighbors)
324 p8est_face_corners[sides[i].face]
326 .insert(subdomain_ids[k]);
334 subids->elem_count = 0;
339 int (&functions<2>::quadrant_compare)(
const void *v1,
const void *v2) =
340 p4est_quadrant_compare;
342 void (&functions<2>::quadrant_childrenv)(
const types<2>::quadrant *q,
343 types<2>::quadrant c[]) =
344 p4est_quadrant_childrenv;
346 int (&functions<2>::quadrant_overlaps_tree)(types<2>::tree * tree,
347 const types<2>::quadrant *q) =
348 p4est_quadrant_overlaps_tree;
350 void (&functions<2>::quadrant_set_morton)(types<2>::quadrant *quadrant,
353 p4est_quadrant_set_morton;
355 int (&functions<2>::quadrant_is_equal)(
const types<2>::quadrant *q1,
356 const types<2>::quadrant *q2) =
357 p4est_quadrant_is_equal;
359 int (&functions<2>::quadrant_is_sibling)(
const types<2>::quadrant *q1,
360 const types<2>::quadrant *q2) =
361 p4est_quadrant_is_sibling;
363 int (&functions<2>::quadrant_is_ancestor)(
const types<2>::quadrant *q1,
364 const types<2>::quadrant *q2) =
365 p4est_quadrant_is_ancestor;
367 int (&functions<2>::quadrant_ancestor_id)(
const types<2>::quadrant *q,
369 p4est_quadrant_ancestor_id;
371 int (&functions<2>::comm_find_owner)(types<2>::forest * p4est,
372 const types<2>::locidx which_tree,
373 const types<2>::quadrant *q,
375 p4est_comm_find_owner;
377 types<2>::connectivity *(&functions<2>::connectivity_new)(
378 types<2>::topidx num_vertices,
379 types<2>::topidx num_trees,
380 types<2>::topidx num_corners,
381 types<2>::topidx num_vtt) = p4est_connectivity_new;
383 void (&functions<2>::connectivity_join_faces)(types<2>::connectivity *conn,
384 types<2>::topidx tree_left,
385 types<2>::topidx tree_right,
389 p4est_connectivity_join_faces;
391 void (&functions<2>::connectivity_destroy)(
392 p4est_connectivity_t *connectivity) = p4est_connectivity_destroy;
394 types<2>::forest *(&functions<2>::new_forest)(
396 types<2>::connectivity *connectivity,
397 types<2>::locidx min_quadrants,
400 std::size_t data_size,
401 p4est_init_t init_fn,
402 void * user_pointer) = p4est_new_ext;
404 void (&functions<2>::destroy)(types<2>::forest *p4est) = p4est_destroy;
406 void (&functions<2>::refine)(types<2>::forest *p4est,
407 int refine_recursive,
408 p4est_refine_t refine_fn,
409 p4est_init_t init_fn) = p4est_refine;
411 void (&functions<2>::coarsen)(types<2>::forest *p4est,
412 int coarsen_recursive,
413 p4est_coarsen_t coarsen_fn,
414 p4est_init_t init_fn) = p4est_coarsen;
416 void (&functions<2>::balance)(types<2>::forest * p4est,
417 types<2>::balance_type btype,
418 p4est_init_t init_fn) = p4est_balance;
420 types<2>::gloidx (&functions<2>::partition)(types<2>::forest *p4est,
421 int partition_for_coarsening,
422 p4est_weight_t weight_fn) =
425 void (&functions<2>::save)(
const char * filename,
426 types<2>::forest *p4est,
427 int save_data) = p4est_save;
429 types<2>::forest *(&functions<2>::load_ext)(
430 const char * filename,
432 std::size_t data_size,
437 types<2>::connectivity **p4est) = p4est_load_ext;
439 int (&functions<2>::connectivity_save)(
440 const char * filename,
441 types<2>::connectivity *connectivity) = p4est_connectivity_save;
443 int (&functions<2>::connectivity_is_valid)(
444 types<2>::connectivity *connectivity) = p4est_connectivity_is_valid;
446 types<2>::connectivity *(&functions<2>::connectivity_load)(
447 const char * filename,
448 std::size_t *length) = p4est_connectivity_load;
450 unsigned int (&functions<2>::checksum)(types<2>::forest *p4est) =
453 void (&functions<2>::vtk_write_file)(types<2>::forest *p4est,
455 const char *baseName) =
456 p4est_vtk_write_file;
458 types<2>::ghost *(&functions<2>::ghost_new)(types<2>::forest * p4est,
459 types<2>::balance_type btype) =
462 void (&functions<2>::ghost_destroy)(types<2>::ghost *ghost) =
465 void (&functions<2>::reset_data)(types<2>::forest *p4est,
466 std::size_t data_size,
467 p4est_init_t init_fn,
468 void *user_pointer) = p4est_reset_data;
470 std::size_t (&functions<2>::forest_memory_used)(types<2>::forest *p4est) =
473 std::size_t (&functions<2>::connectivity_memory_used)(
474 types<2>::connectivity *p4est) = p4est_connectivity_memory_used;
476 template <
int dim,
int spacedim>
477 std::map<unsigned int, std::set<::types::subdomain_id>>
478 compute_vertices_with_ghost_neighbors(
479 const typename ::parallel::distributed::Triangulation<dim, spacedim>
481 typename ::internal::p4est::types<dim>::forest *parallel_forest,
482 typename ::internal::p4est::types<dim>::ghost * parallel_ghost)
484 std::map<unsigned int, std::set<::types::subdomain_id>>
485 vertices_with_ghost_neighbors;
487 ::internal::p4est::FindGhosts<dim, spacedim> fg;
489 fg.triangulation = &tria;
490 fg.vertices_with_ghost_neighbors = &vertices_with_ghost_neighbors;
500 static_cast<void *>(&fg),
502 find_ghosts_face<2, spacedim>,
503 find_ghosts_corner<2, spacedim>);
512 static_cast<void *>(&fg),
514 find_ghosts_face<3, 3>,
515 find_ghosts_edge<3, 3>,
516 find_ghosts_corner<3, 3>);
523 sc_array_destroy(fg.subids);
525 return vertices_with_ghost_neighbors;
528 constexpr
unsigned int functions<2>::max_level;
530 void (&functions<2>::transfer_fixed)(
const types<2>::gloidx *dest_gfq,
531 const types<2>::gloidx *src_gfq,
535 const void * src_data,
536 std::size_t data_size) =
537 p4est_transfer_fixed;
539 types<2>::transfer_context *(&functions<2>::transfer_fixed_begin)(
540 const types<2>::gloidx *dest_gfq,
541 const types<2>::gloidx *src_gfq,
545 const void * src_data,
546 std::size_t data_size) = p4est_transfer_fixed_begin;
548 void (&functions<2>::transfer_fixed_end)(types<2>::transfer_context *tc) =
549 p4est_transfer_fixed_end;
551 void (&functions<2>::transfer_custom)(
const types<2>::gloidx *dest_gfq,
552 const types<2>::gloidx *src_gfq,
556 const int * dest_sizes,
557 const void * src_data,
558 const int * src_sizes) =
559 p4est_transfer_custom;
561 types<2>::transfer_context *(&functions<2>::transfer_custom_begin)(
562 const types<2>::gloidx *dest_gfq,
563 const types<2>::gloidx *src_gfq,
567 const int * dest_sizes,
568 const void * src_data,
569 const int * src_sizes) = p4est_transfer_custom_begin;
571 void (&functions<2>::transfer_custom_end)(types<2>::transfer_context *tc) =
572 p4est_transfer_custom_end;
576 int (&functions<3>::quadrant_compare)(
const void *v1,
const void *v2) =
577 p8est_quadrant_compare;
579 void (&functions<3>::quadrant_childrenv)(
const types<3>::quadrant *q,
580 types<3>::quadrant c[]) =
581 p8est_quadrant_childrenv;
583 int (&functions<3>::quadrant_overlaps_tree)(types<3>::tree * tree,
584 const types<3>::quadrant *q) =
585 p8est_quadrant_overlaps_tree;
587 void (&functions<3>::quadrant_set_morton)(types<3>::quadrant *quadrant,
590 p8est_quadrant_set_morton;
592 int (&functions<3>::quadrant_is_equal)(
const types<3>::quadrant *q1,
593 const types<3>::quadrant *q2) =
594 p8est_quadrant_is_equal;
596 int (&functions<3>::quadrant_is_sibling)(
const types<3>::quadrant *q1,
597 const types<3>::quadrant *q2) =
598 p8est_quadrant_is_sibling;
600 int (&functions<3>::quadrant_is_ancestor)(
const types<3>::quadrant *q1,
601 const types<3>::quadrant *q2) =
602 p8est_quadrant_is_ancestor;
604 int (&functions<3>::quadrant_ancestor_id)(
const types<3>::quadrant *q,
606 p8est_quadrant_ancestor_id;
608 int (&functions<3>::comm_find_owner)(types<3>::forest * p4est,
609 const types<3>::locidx which_tree,
610 const types<3>::quadrant *q,
612 p8est_comm_find_owner;
614 types<3>::connectivity *(&functions<3>::connectivity_new)(
615 types<3>::topidx num_vertices,
616 types<3>::topidx num_trees,
617 types<3>::topidx num_edges,
618 types<3>::topidx num_ett,
619 types<3>::topidx num_corners,
620 types<3>::topidx num_ctt) = p8est_connectivity_new;
622 void (&functions<3>::connectivity_destroy)(
623 p8est_connectivity_t *connectivity) = p8est_connectivity_destroy;
625 void (&functions<3>::connectivity_join_faces)(types<3>::connectivity *conn,
626 types<3>::topidx tree_left,
627 types<3>::topidx tree_right,
631 p8est_connectivity_join_faces;
633 types<3>::forest *(&functions<3>::new_forest)(
635 types<3>::connectivity *connectivity,
636 types<3>::locidx min_quadrants,
639 std::size_t data_size,
640 p8est_init_t init_fn,
641 void * user_pointer) = p8est_new_ext;
643 void (&functions<3>::destroy)(types<3>::forest *p8est) = p8est_destroy;
645 void (&functions<3>::refine)(types<3>::forest *p8est,
646 int refine_recursive,
647 p8est_refine_t refine_fn,
648 p8est_init_t init_fn) = p8est_refine;
650 void (&functions<3>::coarsen)(types<3>::forest *p8est,
651 int coarsen_recursive,
652 p8est_coarsen_t coarsen_fn,
653 p8est_init_t init_fn) = p8est_coarsen;
655 void (&functions<3>::balance)(types<3>::forest * p8est,
656 types<3>::balance_type btype,
657 p8est_init_t init_fn) = p8est_balance;
659 types<3>::gloidx (&functions<3>::partition)(types<3>::forest *p8est,
660 int partition_for_coarsening,
661 p8est_weight_t weight_fn) =
664 void (&functions<3>::save)(
const char * filename,
665 types<3>::forest *p4est,
666 int save_data) = p8est_save;
668 types<3>::forest *(&functions<3>::load_ext)(
669 const char * filename,
671 std::size_t data_size,
676 types<3>::connectivity **p4est) = p8est_load_ext;
678 int (&functions<3>::connectivity_save)(
679 const char * filename,
680 types<3>::connectivity *connectivity) = p8est_connectivity_save;
682 int (&functions<3>::connectivity_is_valid)(
683 types<3>::connectivity *connectivity) = p8est_connectivity_is_valid;
685 types<3>::connectivity *(&functions<3>::connectivity_load)(
686 const char * filename,
687 std::size_t *length) = p8est_connectivity_load;
689 unsigned int (&functions<3>::checksum)(types<3>::forest *p8est) =
692 void (&functions<3>::vtk_write_file)(types<3>::forest *p8est,
694 const char *baseName) =
695 p8est_vtk_write_file;
697 types<3>::ghost *(&functions<3>::ghost_new)(types<3>::forest * p4est,
698 types<3>::balance_type btype) =
701 void (&functions<3>::ghost_destroy)(types<3>::ghost *ghost) =
704 void (&functions<3>::reset_data)(types<3>::forest *p4est,
705 std::size_t data_size,
706 p8est_init_t init_fn,
707 void *user_pointer) = p8est_reset_data;
709 std::size_t (&functions<3>::forest_memory_used)(types<3>::forest *p4est) =
712 std::size_t (&functions<3>::connectivity_memory_used)(
713 types<3>::connectivity *p4est) = p8est_connectivity_memory_used;
715 constexpr
unsigned int functions<3>::max_level;
717 void (&functions<3>::transfer_fixed)(
const types<3>::gloidx *dest_gfq,
718 const types<3>::gloidx *src_gfq,
722 const void * src_data,
723 std::size_t data_size) =
724 p8est_transfer_fixed;
726 types<3>::transfer_context *(&functions<3>::transfer_fixed_begin)(
727 const types<3>::gloidx *dest_gfq,
728 const types<3>::gloidx *src_gfq,
732 const void * src_data,
733 std::size_t data_size) = p8est_transfer_fixed_begin;
735 void (&functions<3>::transfer_fixed_end)(types<3>::transfer_context *tc) =
736 p8est_transfer_fixed_end;
738 void (&functions<3>::transfer_custom)(
const types<3>::gloidx *dest_gfq,
739 const types<3>::gloidx *src_gfq,
743 const int * dest_sizes,
744 const void * src_data,
745 const int * src_sizes) =
746 p8est_transfer_custom;
748 types<3>::transfer_context *(&functions<3>::transfer_custom_begin)(
749 const types<3>::gloidx *dest_gfq,
750 const types<3>::gloidx *src_gfq,
754 const int * dest_sizes,
755 const void * src_data,
756 const int * src_sizes) = p8est_transfer_custom_begin;
758 void (&functions<3>::transfer_custom_end)(types<3>::transfer_context *tc) =
759 p8est_transfer_custom_end;
765 init_quadrant_children(
766 const typename types<dim>::quadrant &p4est_cell,
767 typename types<dim>::quadrant (
770 for (
unsigned int c = 0;
771 c < ::GeometryInfo<dim>::max_children_per_cell;
776 P4EST_QUADRANT_INIT(&p4est_children[c]);
779 P8EST_QUADRANT_INIT(&p4est_children[c]);
786 functions<dim>::quadrant_childrenv(&p4est_cell, p4est_children);
791 init_coarse_quadrant(
typename types<dim>::quadrant &quad)
796 P4EST_QUADRANT_INIT(&quad);
799 P8EST_QUADRANT_INIT(&quad);
804 functions<dim>::quadrant_set_morton(&quad,
811 quadrant_is_equal(
const typename types<dim>::quadrant &q1,
812 const typename types<dim>::quadrant &q2)
814 return functions<dim>::quadrant_is_equal(&q1, &q2);
821 quadrant_is_ancestor(
const typename types<dim>::quadrant &q1,
822 const typename types<dim>::quadrant &q2)
824 return functions<dim>::quadrant_is_ancestor(&q1, &q2);
829 tree_exists_locally(
const typename types<dim>::forest *parallel_forest,
830 const typename types<dim>::topidx coarse_grid_cell)
832 Assert(coarse_grid_cell < parallel_forest->connectivity->num_trees,
834 return ((coarse_grid_cell >= parallel_forest->first_local_tree) &&
835 (coarse_grid_cell <= parallel_forest->last_local_tree));
840 #endif // DEAL_II_WITH_P4EST 843 #include "p4est_wrappers.inst" 846 DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcMessage(std::string arg1)
unsigned int subdomain_id
#define Assert(cond, exc)
unsigned int global_dof_index
static ::ExceptionBase & ExcNotImplemented()
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
static ::ExceptionBase & ExcInternalError()