Reference documentation for deal.II version GIT relicensing-216-gcda843ca70 2024-03-28 12:20:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
p4est_wrappers.cc
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2016 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15
18
20
21#ifdef DEAL_II_WITH_P4EST
22
23namespace internal
24{
25 namespace p4est
26 {
27 namespace
28 {
29 template <int dim, int spacedim>
30 typename ::Triangulation<dim, spacedim>::cell_iterator
31 cell_from_quad(
32 const ::parallel::distributed::Triangulation<dim, spacedim>
34 const typename ::internal::p4est::types<dim>::topidx treeidx,
35 const typename ::internal::p4est::types<dim>::quadrant &quad)
36 {
37 int i, l = quad.level;
38 ::types::global_dof_index dealii_index =
39 triangulation->get_p4est_tree_to_coarse_cell_permutation()[treeidx];
40
41 for (i = 0; i < l; ++i)
42 {
43 typename ::Triangulation<dim, spacedim>::cell_iterator cell(
44 triangulation, i, dealii_index);
45 const int child_id =
47 &quad, i + 1);
48 Assert(cell->has_children(),
49 ExcMessage("p4est quadrant does not correspond to a cell!"));
50 dealii_index = cell->child_index(child_id);
51 }
52
53 typename ::Triangulation<dim, spacedim>::cell_iterator out_cell(
54 triangulation, l, dealii_index);
55
56 return out_cell;
57 }
58
63 template <int dim, int spacedim>
64 struct FindGhosts
65 {
66 const typename ::parallel::distributed::Triangulation<dim,
67 spacedim>
68 *triangulation;
69 sc_array_t *subids;
70 std::map<unsigned int, std::set<::types::subdomain_id>>
72 };
73
74
80 template <int dim, int spacedim>
81 void
82 find_ghosts_corner(
83 typename ::internal::p4est::iter<dim>::corner_info *info,
84 void *user_data)
85 {
86 int i, j;
87 int nsides = info->sides.elem_count;
88 auto *sides = reinterpret_cast<
89 typename ::internal::p4est::iter<dim>::corner_side *>(
90 info->sides.array);
91 FindGhosts<dim, spacedim> *fg =
92 static_cast<FindGhosts<dim, spacedim> *>(user_data);
93 sc_array_t *subids = fg->subids;
94 const ::parallel::distributed::Triangulation<dim, spacedim>
95 *triangulation = fg->triangulation;
96 int nsubs;
97 ::types::subdomain_id *subdomain_ids;
98 std::map<unsigned int, std::set<::types::subdomain_id>>
99 *vertices_with_ghost_neighbors = fg->vertices_with_ghost_neighbors;
100
101 subids->elem_count = 0;
102 for (i = 0; i < nsides; ++i)
103 {
104 if (sides[i].is_ghost)
105 {
106 typename ::parallel::distributed::
107 Triangulation<dim, spacedim>::cell_iterator cell =
108 cell_from_quad(triangulation,
109 sides[i].treeid,
110 *(sides[i].quad));
111 Assert(cell->is_ghost(),
112 ExcMessage("ghost quad did not find ghost cell"));
113 ::types::subdomain_id *subid =
114 static_cast<::types::subdomain_id *>(
115 sc_array_push(subids));
116 *subid = cell->subdomain_id();
117 }
118 }
119
120 if (!subids->elem_count)
121 {
122 return;
123 }
124
125 nsubs = static_cast<int>(subids->elem_count);
126 subdomain_ids =
127 reinterpret_cast<::types::subdomain_id *>(subids->array);
128
129 for (i = 0; i < nsides; ++i)
130 {
131 if (!sides[i].is_ghost)
132 {
133 typename ::parallel::distributed::
134 Triangulation<dim, spacedim>::cell_iterator cell =
135 cell_from_quad(triangulation,
136 sides[i].treeid,
137 *(sides[i].quad));
138
139 Assert(!cell->is_ghost(),
140 ExcMessage("local quad found ghost cell"));
141
142 for (j = 0; j < nsubs; ++j)
143 {
144 (*vertices_with_ghost_neighbors)[cell->vertex_index(
145 sides[i].corner)]
146 .insert(subdomain_ids[j]);
147 }
148 }
149 }
150
151 subids->elem_count = 0;
152 }
153
157 template <int dim, int spacedim>
158 void
159 find_ghosts_edge(
160 typename ::internal::p4est::iter<dim>::edge_info *info,
161 void *user_data)
162 {
163 int i, j, k;
164 int nsides = info->sides.elem_count;
165 auto *sides = reinterpret_cast<
166 typename ::internal::p4est::iter<dim>::edge_side *>(
167 info->sides.array);
168 auto *fg = static_cast<FindGhosts<dim, spacedim> *>(user_data);
169 sc_array_t *subids = fg->subids;
170 const ::parallel::distributed::Triangulation<dim, spacedim>
171 *triangulation = fg->triangulation;
172 int nsubs;
173 ::types::subdomain_id *subdomain_ids;
174 std::map<unsigned int, std::set<::types::subdomain_id>>
175 *vertices_with_ghost_neighbors = fg->vertices_with_ghost_neighbors;
176
177 subids->elem_count = 0;
178 for (i = 0; i < nsides; ++i)
179 {
180 if (sides[i].is_hanging)
181 {
182 for (j = 0; j < 2; ++j)
183 {
184 if (sides[i].is.hanging.is_ghost[j])
185 {
186 typename ::parallel::distributed::
187 Triangulation<dim, spacedim>::cell_iterator cell =
188 cell_from_quad(triangulation,
189 sides[i].treeid,
190 *(sides[i].is.hanging.quad[j]));
191 ::types::subdomain_id *subid =
192 static_cast<::types::subdomain_id *>(
193 sc_array_push(subids));
194 *subid = cell->subdomain_id();
195 }
196 }
197 }
198 }
199
200 if (!subids->elem_count)
201 {
202 return;
203 }
204
205 nsubs = static_cast<int>(subids->elem_count);
206 subdomain_ids =
207 reinterpret_cast<::types::subdomain_id *>(subids->array);
208
209 for (i = 0; i < nsides; ++i)
210 {
211 if (sides[i].is_hanging)
212 {
213 for (j = 0; j < 2; ++j)
214 {
215 if (!sides[i].is.hanging.is_ghost[j])
216 {
217 typename ::parallel::distributed::
218 Triangulation<dim, spacedim>::cell_iterator cell =
219 cell_from_quad(triangulation,
220 sides[i].treeid,
221 *(sides[i].is.hanging.quad[j]));
222
223 for (k = 0; k < nsubs; ++k)
224 {
225 (*vertices_with_ghost_neighbors)
226 [cell->vertex_index(
227 p8est_edge_corners[sides[i].edge][1 ^ j])]
228 .insert(subdomain_ids[k]);
229 }
230 }
231 }
232 }
233 }
234
235 subids->elem_count = 0;
236 }
237
241 template <int dim, int spacedim>
242 void
243 find_ghosts_face(
244 typename ::internal::p4est::iter<dim>::face_info *info,
245 void *user_data)
246 {
247 int i, j, k;
248 int nsides = info->sides.elem_count;
249 auto *sides = reinterpret_cast<
250 typename ::internal::p4est::iter<dim>::face_side *>(
251 info->sides.array);
252 FindGhosts<dim, spacedim> *fg =
253 static_cast<FindGhosts<dim, spacedim> *>(user_data);
254 sc_array_t *subids = fg->subids;
255 const ::parallel::distributed::Triangulation<dim, spacedim>
256 *triangulation = fg->triangulation;
257 int nsubs;
258 ::types::subdomain_id *subdomain_ids;
259 std::map<unsigned int, std::set<::types::subdomain_id>>
260 *vertices_with_ghost_neighbors = fg->vertices_with_ghost_neighbors;
261 int limit = (dim == 2) ? 2 : 4;
262
263 subids->elem_count = 0;
264 for (i = 0; i < nsides; ++i)
265 {
266 if (sides[i].is_hanging)
267 {
268 for (j = 0; j < limit; ++j)
269 {
270 if (sides[i].is.hanging.is_ghost[j])
271 {
272 typename ::parallel::distributed::
273 Triangulation<dim, spacedim>::cell_iterator cell =
274 cell_from_quad(triangulation,
275 sides[i].treeid,
276 *(sides[i].is.hanging.quad[j]));
277 ::types::subdomain_id *subid =
278 static_cast<::types::subdomain_id *>(
279 sc_array_push(subids));
280 *subid = cell->subdomain_id();
281 }
282 }
283 }
284 }
285
286 if (!subids->elem_count)
287 {
288 return;
289 }
290
291 nsubs = static_cast<int>(subids->elem_count);
292 subdomain_ids =
293 reinterpret_cast<::types::subdomain_id *>(subids->array);
294
295 for (i = 0; i < nsides; ++i)
296 {
297 if (sides[i].is_hanging)
298 {
299 for (j = 0; j < limit; ++j)
300 {
301 if (!sides[i].is.hanging.is_ghost[j])
302 {
303 typename ::parallel::distributed::
304 Triangulation<dim, spacedim>::cell_iterator cell =
305 cell_from_quad(triangulation,
306 sides[i].treeid,
307 *(sides[i].is.hanging.quad[j]));
308
309 for (k = 0; k < nsubs; ++k)
310 {
311 if (dim == 2)
312 {
313 (*vertices_with_ghost_neighbors)
314 [cell->vertex_index(
315 p4est_face_corners[sides[i].face]
316 [(limit - 1) ^ j])]
317 .insert(subdomain_ids[k]);
318 }
319 else
320 {
321 (*vertices_with_ghost_neighbors)
322 [cell->vertex_index(
323 p8est_face_corners[sides[i].face]
324 [(limit - 1) ^ j])]
325 .insert(subdomain_ids[k]);
326 }
327 }
328 }
329 }
330 }
331 }
332
333 subids->elem_count = 0;
334 }
335 } // namespace
336
337
338 int (&functions<2>::quadrant_compare)(const void *v1, const void *v2) =
339 p4est_quadrant_compare;
340
342 types<2>::quadrant c[]) =
343 p4est_quadrant_childrenv;
344
346 const types<2>::quadrant *q) =
347 p4est_quadrant_overlaps_tree;
348
350 int level,
351 std::uint64_t id) =
352 p4est_quadrant_set_morton;
353
355 const types<2>::quadrant *q2) =
356 p4est_quadrant_is_equal;
357
359 const types<2>::quadrant *q2) =
360 p4est_quadrant_is_sibling;
361
363 const types<2>::quadrant *q2) =
364 p4est_quadrant_is_ancestor;
365
367 int level) =
368 p4est_quadrant_ancestor_id;
369
371 const types<2>::locidx which_tree,
372 const types<2>::quadrant *q,
373 const int guess) =
374 p4est_comm_find_owner;
375
377 types<2>::topidx num_vertices,
378 types<2>::topidx num_trees,
379 types<2>::topidx num_corners,
380 types<2>::topidx num_vtt) = p4est_connectivity_new;
381
383 types<2>::topidx num_vertices,
384 types<2>::topidx num_trees,
385 types<2>::topidx num_corners,
386 const double *vertices,
387 const types<2>::topidx *ttv,
388 const types<2>::topidx *ttt,
389 const int8_t *ttf,
390 const types<2>::topidx *ttc,
391 const types<2>::topidx *coff,
392 const types<2>::topidx *ctt,
393 const int8_t *ctc) = p4est_connectivity_new_copy;
394
396 types<2>::topidx tree_left,
397 types<2>::topidx tree_right,
398 int face_left,
399 int face_right,
400 int orientation) =
401 p4est_connectivity_join_faces;
402
404 p4est_connectivity_t *connectivity) = p4est_connectivity_destroy;
405
407 MPI_Comm mpicomm,
408 types<2>::connectivity *connectivity,
409 types<2>::locidx min_quadrants,
410 int min_level,
411 int fill_uniform,
412 std::size_t data_size,
413 p4est_init_t init_fn,
414 void *user_pointer) = p4est_new_ext;
415
417 int copy_data) = p4est_copy;
418
419 void (&functions<2>::destroy)(types<2>::forest *p4est) = p4est_destroy;
420
422 int refine_recursive,
423 p4est_refine_t refine_fn,
424 p4est_init_t init_fn) = p4est_refine;
425
427 int coarsen_recursive,
428 p4est_coarsen_t coarsen_fn,
429 p4est_init_t init_fn) = p4est_coarsen;
430
433 p4est_init_t init_fn) = p4est_balance;
434
436 int partition_for_coarsening,
437 p4est_weight_t weight_fn) =
438 p4est_partition_ext;
439
440 void (&functions<2>::save)(const char *filename,
441 types<2>::forest *p4est,
442 int save_data) = p4est_save;
443
445 const char *filename,
446 MPI_Comm mpicomm,
447 std::size_t data_size,
448 int load_data,
449 int autopartition,
450 int broadcasthead,
451 void *user_pointer,
452 types<2>::connectivity **p4est) = p4est_load_ext;
453
455 const char *filename,
456 types<2>::connectivity *connectivity) = p4est_connectivity_save;
457
459 types<2>::connectivity *connectivity) = p4est_connectivity_is_valid;
460
462 const char *filename,
463 std::size_t *length) = p4est_connectivity_load;
464
465 unsigned int (&functions<2>::checksum)(types<2>::forest *p4est) =
466 p4est_checksum;
467
469 p4est_geometry_t *,
470 const char *baseName) =
471 p4est_vtk_write_file;
472
475 p4est_ghost_new;
476
478 p4est_ghost_destroy;
479
481 std::size_t data_size,
482 p4est_init_t init_fn,
483 void *user_pointer) = p4est_reset_data;
484
486 p4est_memory_used;
487
489 types<2>::connectivity *p4est) = p4est_connectivity_memory_used;
490
491 constexpr unsigned int functions<2>::max_level;
492
493 void (&functions<2>::transfer_fixed)(const types<2>::gloidx *dest_gfq,
494 const types<2>::gloidx *src_gfq,
495 MPI_Comm mpicomm,
496 int tag,
497 void *dest_data,
498 const void *src_data,
499 std::size_t data_size) =
500 p4est_transfer_fixed;
501
503 const types<2>::gloidx *dest_gfq,
504 const types<2>::gloidx *src_gfq,
505 MPI_Comm mpicomm,
506 int tag,
507 void *dest_data,
508 const void *src_data,
509 std::size_t data_size) = p4est_transfer_fixed_begin;
510
512 p4est_transfer_fixed_end;
513
514 void (&functions<2>::transfer_custom)(const types<2>::gloidx *dest_gfq,
515 const types<2>::gloidx *src_gfq,
516 MPI_Comm mpicomm,
517 int tag,
518 void *dest_data,
519 const int *dest_sizes,
520 const void *src_data,
521 const int *src_sizes) =
522 p4est_transfer_custom;
523
525 const types<2>::gloidx *dest_gfq,
526 const types<2>::gloidx *src_gfq,
527 MPI_Comm mpicomm,
528 int tag,
529 void *dest_data,
530 const int *dest_sizes,
531 const void *src_data,
532 const int *src_sizes) = p4est_transfer_custom_begin;
533
535 p4est_transfer_custom_end;
536
537# ifdef P4EST_SEARCH_LOCAL
539 types<2>::forest *p4est,
540 int call_post,
543 sc_array_t *points) = p4est_search_partition;
544# endif
545
547 types<2>::connectivity *connectivity,
548 types<2>::topidx treeid,
551 double vxyz[3]) = p4est_qcoord_to_vertex;
552
553 int (&functions<3>::quadrant_compare)(const void *v1, const void *v2) =
554 p8est_quadrant_compare;
555
557 types<3>::quadrant c[]) =
558 p8est_quadrant_childrenv;
559
561 const types<3>::quadrant *q) =
562 p8est_quadrant_overlaps_tree;
563
565 int level,
566 std::uint64_t id) =
567 p8est_quadrant_set_morton;
568
570 const types<3>::quadrant *q2) =
571 p8est_quadrant_is_equal;
572
574 const types<3>::quadrant *q2) =
575 p8est_quadrant_is_sibling;
576
578 const types<3>::quadrant *q2) =
579 p8est_quadrant_is_ancestor;
580
582 int level) =
583 p8est_quadrant_ancestor_id;
584
586 const types<3>::locidx which_tree,
587 const types<3>::quadrant *q,
588 const int guess) =
589 p8est_comm_find_owner;
590
592 types<3>::topidx num_vertices,
593 types<3>::topidx num_trees,
594 types<3>::topidx num_edges,
595 types<3>::topidx num_ett,
596 types<3>::topidx num_corners,
597 types<3>::topidx num_ctt) = p8est_connectivity_new;
598
600 types<3>::topidx num_vertices,
601 types<3>::topidx num_trees,
602 types<3>::topidx num_edges,
603 types<3>::topidx num_corners,
604 const double *vertices,
605 const types<3>::topidx *ttv,
606 const types<3>::topidx *ttt,
607 const int8_t *ttf,
608 const types<3>::topidx *tte,
609 const types<3>::topidx *eoff,
610 const types<3>::topidx *ett,
611 const int8_t *ete,
612 const types<3>::topidx *ttc,
613 const types<3>::topidx *coff,
614 const types<3>::topidx *ctt,
615 const int8_t *ctc) = p8est_connectivity_new_copy;
616
618 p8est_connectivity_t *connectivity) = p8est_connectivity_destroy;
619
621 types<3>::topidx tree_left,
622 types<3>::topidx tree_right,
623 int face_left,
624 int face_right,
625 int orientation) =
626 p8est_connectivity_join_faces;
627
629 MPI_Comm mpicomm,
630 types<3>::connectivity *connectivity,
631 types<3>::locidx min_quadrants,
632 int min_level,
633 int fill_uniform,
634 std::size_t data_size,
635 p8est_init_t init_fn,
636 void *user_pointer) = p8est_new_ext;
637
639 int copy_data) = p8est_copy;
640
641 void (&functions<3>::destroy)(types<3>::forest *p8est) = p8est_destroy;
642
644 int refine_recursive,
645 p8est_refine_t refine_fn,
646 p8est_init_t init_fn) = p8est_refine;
647
649 int coarsen_recursive,
650 p8est_coarsen_t coarsen_fn,
651 p8est_init_t init_fn) = p8est_coarsen;
652
655 p8est_init_t init_fn) = p8est_balance;
656
658 int partition_for_coarsening,
659 p8est_weight_t weight_fn) =
660 p8est_partition_ext;
661
662 void (&functions<3>::save)(const char *filename,
663 types<3>::forest *p4est,
664 int save_data) = p8est_save;
665
667 const char *filename,
668 MPI_Comm mpicomm,
669 std::size_t data_size,
670 int load_data,
671 int autopartition,
672 int broadcasthead,
673 void *user_pointer,
674 types<3>::connectivity **p4est) = p8est_load_ext;
675
677 const char *filename,
678 types<3>::connectivity *connectivity) = p8est_connectivity_save;
679
681 types<3>::connectivity *connectivity) = p8est_connectivity_is_valid;
682
684 const char *filename,
685 std::size_t *length) = p8est_connectivity_load;
686
687 unsigned int (&functions<3>::checksum)(types<3>::forest *p8est) =
688 p8est_checksum;
689
691 p8est_geometry_t *,
692 const char *baseName) =
693 p8est_vtk_write_file;
694
697 p8est_ghost_new;
698
700 p8est_ghost_destroy;
701
703 std::size_t data_size,
704 p8est_init_t init_fn,
705 void *user_pointer) = p8est_reset_data;
706
708 p8est_memory_used;
709
711 types<3>::connectivity *p4est) = p8est_connectivity_memory_used;
712
713 constexpr unsigned int functions<3>::max_level;
714
715 void (&functions<3>::transfer_fixed)(const types<3>::gloidx *dest_gfq,
716 const types<3>::gloidx *src_gfq,
717 MPI_Comm mpicomm,
718 int tag,
719 void *dest_data,
720 const void *src_data,
721 std::size_t data_size) =
722 p8est_transfer_fixed;
723
725 const types<3>::gloidx *dest_gfq,
726 const types<3>::gloidx *src_gfq,
727 MPI_Comm mpicomm,
728 int tag,
729 void *dest_data,
730 const void *src_data,
731 std::size_t data_size) = p8est_transfer_fixed_begin;
732
734 p8est_transfer_fixed_end;
735
736 void (&functions<3>::transfer_custom)(const types<3>::gloidx *dest_gfq,
737 const types<3>::gloidx *src_gfq,
738 MPI_Comm mpicomm,
739 int tag,
740 void *dest_data,
741 const int *dest_sizes,
742 const void *src_data,
743 const int *src_sizes) =
744 p8est_transfer_custom;
745
747 const types<3>::gloidx *dest_gfq,
748 const types<3>::gloidx *src_gfq,
749 MPI_Comm mpicomm,
750 int tag,
751 void *dest_data,
752 const int *dest_sizes,
753 const void *src_data,
754 const int *src_sizes) = p8est_transfer_custom_begin;
755
757 p8est_transfer_custom_end;
758
759# ifdef P4EST_SEARCH_LOCAL
761 types<3>::forest *p4est,
762 int call_post,
765 sc_array_t *points) = p8est_search_partition;
766# endif
767
769 types<3>::connectivity *connectivity,
770 types<3>::topidx treeid,
774 double vxyz[3]) = p8est_qcoord_to_vertex;
775
776 template <int dim>
777 void
779 const typename types<dim>::quadrant &p4est_cell,
780 typename types<dim>::quadrant (
782 {
783 for (unsigned int c = 0;
784 c < ::GeometryInfo<dim>::max_children_per_cell;
785 ++c)
786 switch (dim)
787 {
788 case 2:
789 P4EST_QUADRANT_INIT(&p4est_children[c]);
790 break;
791 case 3:
792 P8EST_QUADRANT_INIT(&p4est_children[c]);
793 break;
794 default:
796 }
797
798
799 functions<dim>::quadrant_childrenv(&p4est_cell, p4est_children);
800 }
801
802 template <int dim>
803 void
805 {
806 switch (dim)
807 {
808 case 2:
809 P4EST_QUADRANT_INIT(&quad);
810 break;
811 case 3:
812 P8EST_QUADRANT_INIT(&quad);
813 break;
814 default:
816 }
818 /*level=*/0,
819 /*index=*/0);
820 }
821
822 template <int dim>
823 bool
825 const typename types<dim>::quadrant &q2)
826 {
827 return functions<dim>::quadrant_is_equal(&q1, &q2);
828 }
829
830
831
832 template <int dim>
833 bool
835 const typename types<dim>::quadrant &q2)
836 {
838 }
839
840 template <int dim>
841 bool
842 tree_exists_locally(const typename types<dim>::forest *parallel_forest,
843 const typename types<dim>::topidx coarse_grid_cell)
844 {
845 Assert(coarse_grid_cell < parallel_forest->connectivity->num_trees,
847 return ((coarse_grid_cell >= parallel_forest->first_local_tree) &&
848 (coarse_grid_cell <= parallel_forest->last_local_tree));
849 }
850
851
852
853 // template specializations
854
855 template <>
857 copy_connectivity<2>(const typename types<2>::connectivity *connectivity)
858 {
860 connectivity->num_vertices,
861 connectivity->num_trees,
862 connectivity->num_corners,
863 connectivity->vertices,
864 connectivity->tree_to_vertex,
865 connectivity->tree_to_tree,
866 connectivity->tree_to_face,
867 connectivity->tree_to_corner,
868 connectivity->ctt_offset,
869 connectivity->corner_to_tree,
870 connectivity->corner_to_corner);
871 }
872
873 template <>
875 copy_connectivity<3>(const typename types<3>::connectivity *connectivity)
876 {
878 connectivity->num_vertices,
879 connectivity->num_trees,
880 connectivity->num_edges,
881 connectivity->num_corners,
882 connectivity->vertices,
883 connectivity->tree_to_vertex,
884 connectivity->tree_to_tree,
885 connectivity->tree_to_face,
886 connectivity->tree_to_edge,
887 connectivity->ett_offset,
888 connectivity->edge_to_tree,
889 connectivity->edge_to_edge,
890 connectivity->tree_to_corner,
891 connectivity->ctt_offset,
892 connectivity->corner_to_tree,
893 connectivity->corner_to_corner);
894 }
895
896
897
898 template <>
899 bool
900 quadrant_is_equal<1>(const typename types<1>::quadrant &q1,
901 const typename types<1>::quadrant &q2)
902 {
903 return q1 == q2;
904 }
905
906
907
908 template <>
909 bool
911 const types<1>::quadrant &q2)
912 {
913 // determine level of quadrants
914 const int level_1 = (q1 << types<1>::max_n_child_indices_bits) >>
916 const int level_2 = (q2 << types<1>::max_n_child_indices_bits) >>
918
919 // q1 can be an ancestor of q2 if q1's level is smaller
920 if (level_1 >= level_2)
921 return false;
922
923 // extract path of quadrants up to level of possible ancestor q1
924 const int truncated_id_1 = (q1 >> (types<1>::n_bits - 1 - level_1))
925 << (types<1>::n_bits - 1 - level_1);
926 const int truncated_id_2 = (q2 >> (types<1>::n_bits - 1 - level_1))
927 << (types<1>::n_bits - 1 - level_1);
928
929 // compare paths
930 return truncated_id_1 == truncated_id_2;
931 }
932
933
934
935 template <>
936 void
938 const typename types<1>::quadrant &q,
939 typename types<1>::quadrant (
941 {
942 // determine the current level of quadrant
943 const int level_parent = (q << types<1>::max_n_child_indices_bits) >>
945 const int level_child = level_parent + 1;
946
947 // left child: only n_child_indices has to be incremented
948 p4est_children[0] = (q + 1);
949
950 // right child: increment and set a bit to 1 indicating that it is a right
951 // child
952 p4est_children[1] = (q + 1) | (1 << (types<1>::n_bits - 1 - level_child));
953 }
954
955
956
957 template <>
958 void
960 {
961 quad = 0;
962 }
963
964 } // namespace p4est
965} // namespace internal
966
967#endif // DEAL_II_WITH_P4EST
968
969/*-------------- Explicit Instantiations -------------------------------*/
970#include "p4est_wrappers.inst"
971
972
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
Point< 3 > vertices[4]
unsigned int level
Definition grid_out.cc:4616
const unsigned int v1
#define Assert(cond, exc)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define DEAL_II_NOT_IMPLEMENTED()
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
void init_quadrant_children< 1 >(const typename types< 1 >::quadrant &q, typename types< 1 >::quadrant(&p4est_children)[::GeometryInfo< 1 >::max_children_per_cell])
void init_coarse_quadrant(typename types< dim >::quadrant &quad)
bool quadrant_is_equal< 1 >(const typename types< 1 >::quadrant &q1, const typename types< 1 >::quadrant &q2)
types< 2 >::connectivity * copy_connectivity< 2 >(const typename types< 2 >::connectivity *connectivity)
bool quadrant_is_ancestor< 1 >(const types< 1 >::quadrant &q1, const types< 1 >::quadrant &q2)
types< 3 >::connectivity * copy_connectivity< 3 >(const typename types< 3 >::connectivity *connectivity)
bool quadrant_is_equal(const typename types< dim >::quadrant &q1, const typename types< dim >::quadrant &q2)
void init_quadrant_children(const typename types< dim >::quadrant &p4est_cell, typename types< dim >::quadrant(&p4est_children)[::GeometryInfo< dim >::max_children_per_cell])
bool quadrant_is_ancestor(const typename types< dim >::quadrant &q1, const typename types< dim >::quadrant &q2)
void init_coarse_quadrant< 1 >(typename types< 1 >::quadrant &quad)
bool tree_exists_locally(const typename types< dim >::forest *parallel_forest, const typename types< dim >::topidx coarse_grid_cell)
sc_array_t * subids
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation
std::map< unsigned int, std::set<::types::subdomain_id > > * vertices_with_ghost_neighbors