Loading [MathJax]/extensions/TeX/newcommand.js
 Reference documentation for deal.II version 9.1.1
\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
p4est_wrappers.cc
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2008 - 2019 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 
17 #include <deal.II/distributed/p4est_wrappers.h>
18 #include <deal.II/distributed/tria.h>
19 
20 DEAL_II_NAMESPACE_OPEN
21 
22 #ifdef DEAL_II_WITH_P4EST
23 
24 namespace internal
25 {
26  namespace p4est
27  {
28  namespace
29  {
30  template <int dim, int spacedim>
31  typename ::Triangulation<dim, spacedim>::cell_iterator
32  cell_from_quad(
33  const ::parallel::distributed::Triangulation<dim, spacedim>
34  *triangulation,
35  const typename ::internal::p4est::types<dim>::topidx treeidx,
36  const typename ::internal::p4est::types<dim>::quadrant &quad)
37  {
38  int i, l = quad.level;
39  ::types::global_dof_index dealii_index =
40  triangulation->get_p4est_tree_to_coarse_cell_permutation()[treeidx];
41 
42  for (i = 0; i < l; i++)
43  {
44  typename ::Triangulation<dim, spacedim>::cell_iterator cell(
45  triangulation, i, dealii_index);
46  const int child_id =
48  &quad, i + 1);
49  Assert(cell->has_children(),
50  ExcMessage("p4est quadrant does not correspond to a cell!"));
51  dealii_index = cell->child_index(child_id);
52  }
53 
54  typename ::Triangulation<dim, spacedim>::cell_iterator out_cell(
55  triangulation, l, dealii_index);
56 
57  return out_cell;
58  }
59 
64  template <int dim, int spacedim>
65  struct FindGhosts
66  {
67  const typename ::parallel::distributed::Triangulation<dim,
68  spacedim>
69  * triangulation;
70  sc_array_t *subids;
71  std::map<unsigned int, std::set<::types::subdomain_id>>
72  *vertices_with_ghost_neighbors;
73  };
74 
75 
81  template <int dim, int spacedim>
82  void
83  find_ghosts_corner(
84  typename ::internal::p4est::iter<dim>::corner_info *info,
85  void * user_data)
86  {
87  int i, j;
88  int nsides = info->sides.elem_count;
89  auto *sides = reinterpret_cast<
90  typename ::internal::p4est::iter<dim>::corner_side *>(
91  info->sides.array);
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;
97  int nsubs;
98  ::types::subdomain_id *subdomain_ids;
99  std::map<unsigned int, std::set<::types::subdomain_id>>
100  *vertices_with_ghost_neighbors = fg->vertices_with_ghost_neighbors;
101 
102  subids->elem_count = 0;
103  for (i = 0; i < nsides; i++)
104  {
105  if (sides[i].is_ghost)
106  {
107  typename ::parallel::distributed::
108  Triangulation<dim, spacedim>::cell_iterator cell =
109  cell_from_quad(triangulation,
110  sides[i].treeid,
111  *(sides[i].quad));
112  Assert(cell->is_ghost(),
113  ExcMessage("ghost quad did not find ghost cell"));
114  ::types::subdomain_id *subid =
115  static_cast<::types::subdomain_id *>(
116  sc_array_push(subids));
117  *subid = cell->subdomain_id();
118  }
119  }
120 
121  if (!subids->elem_count)
122  {
123  return;
124  }
125 
126  nsubs = static_cast<int>(subids->elem_count);
127  subdomain_ids =
128  reinterpret_cast<::types::subdomain_id *>(subids->array);
129 
130  for (i = 0; i < nsides; i++)
131  {
132  if (!sides[i].is_ghost)
133  {
134  typename ::parallel::distributed::
135  Triangulation<dim, spacedim>::cell_iterator cell =
136  cell_from_quad(triangulation,
137  sides[i].treeid,
138  *(sides[i].quad));
139 
140  Assert(!cell->is_ghost(),
141  ExcMessage("local quad found ghost cell"));
142 
143  for (j = 0; j < nsubs; j++)
144  {
145  (*vertices_with_ghost_neighbors)[cell->vertex_index(
146  sides[i].corner)]
147  .insert(subdomain_ids[j]);
148  }
149  }
150  }
151 
152  subids->elem_count = 0;
153  }
154 
158  template <int dim, int spacedim>
159  void
160  find_ghosts_edge(
161  typename ::internal::p4est::iter<dim>::edge_info *info,
162  void * user_data)
163  {
164  int i, j, k;
165  int nsides = info->sides.elem_count;
166  auto *sides = reinterpret_cast<
167  typename ::internal::p4est::iter<dim>::edge_side *>(
168  info->sides.array);
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;
173  int nsubs;
174  ::types::subdomain_id *subdomain_ids;
175  std::map<unsigned int, std::set<::types::subdomain_id>>
176  *vertices_with_ghost_neighbors = fg->vertices_with_ghost_neighbors;
177 
178  subids->elem_count = 0;
179  for (i = 0; i < nsides; i++)
180  {
181  if (sides[i].is_hanging)
182  {
183  for (j = 0; j < 2; j++)
184  {
185  if (sides[i].is.hanging.is_ghost[j])
186  {
187  typename ::parallel::distributed::
188  Triangulation<dim, spacedim>::cell_iterator cell =
189  cell_from_quad(triangulation,
190  sides[i].treeid,
191  *(sides[i].is.hanging.quad[j]));
192  ::types::subdomain_id *subid =
193  static_cast<::types::subdomain_id *>(
194  sc_array_push(subids));
195  *subid = cell->subdomain_id();
196  }
197  }
198  }
199  }
200 
201  if (!subids->elem_count)
202  {
203  return;
204  }
205 
206  nsubs = static_cast<int>(subids->elem_count);
207  subdomain_ids =
208  reinterpret_cast<::types::subdomain_id *>(subids->array);
209 
210  for (i = 0; i < nsides; i++)
211  {
212  if (sides[i].is_hanging)
213  {
214  for (j = 0; j < 2; j++)
215  {
216  if (!sides[i].is.hanging.is_ghost[j])
217  {
218  typename ::parallel::distributed::
219  Triangulation<dim, spacedim>::cell_iterator cell =
220  cell_from_quad(triangulation,
221  sides[i].treeid,
222  *(sides[i].is.hanging.quad[j]));
223 
224  for (k = 0; k < nsubs; k++)
225  {
226  (*vertices_with_ghost_neighbors)
227  [cell->vertex_index(
228  p8est_edge_corners[sides[i].edge][1 ^ j])]
229  .insert(subdomain_ids[k]);
230  }
231  }
232  }
233  }
234  }
235 
236  subids->elem_count = 0;
237  }
238 
242  template <int dim, int spacedim>
243  void
244  find_ghosts_face(
245  typename ::internal::p4est::iter<dim>::face_info *info,
246  void * user_data)
247  {
248  int i, j, k;
249  int nsides = info->sides.elem_count;
250  auto *sides = reinterpret_cast<
251  typename ::internal::p4est::iter<dim>::face_side *>(
252  info->sides.array);
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;
258  int nsubs;
259  ::types::subdomain_id *subdomain_ids;
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;
263 
264  subids->elem_count = 0;
265  for (i = 0; i < nsides; i++)
266  {
267  if (sides[i].is_hanging)
268  {
269  for (j = 0; j < limit; j++)
270  {
271  if (sides[i].is.hanging.is_ghost[j])
272  {
273  typename ::parallel::distributed::
274  Triangulation<dim, spacedim>::cell_iterator cell =
275  cell_from_quad(triangulation,
276  sides[i].treeid,
277  *(sides[i].is.hanging.quad[j]));
278  ::types::subdomain_id *subid =
279  static_cast<::types::subdomain_id *>(
280  sc_array_push(subids));
281  *subid = cell->subdomain_id();
282  }
283  }
284  }
285  }
286 
287  if (!subids->elem_count)
288  {
289  return;
290  }
291 
292  nsubs = static_cast<int>(subids->elem_count);
293  subdomain_ids =
294  reinterpret_cast<::types::subdomain_id *>(subids->array);
295 
296  for (i = 0; i < nsides; i++)
297  {
298  if (sides[i].is_hanging)
299  {
300  for (j = 0; j < limit; j++)
301  {
302  if (!sides[i].is.hanging.is_ghost[j])
303  {
304  typename ::parallel::distributed::
305  Triangulation<dim, spacedim>::cell_iterator cell =
306  cell_from_quad(triangulation,
307  sides[i].treeid,
308  *(sides[i].is.hanging.quad[j]));
309 
310  for (k = 0; k < nsubs; k++)
311  {
312  if (dim == 2)
313  {
314  (*vertices_with_ghost_neighbors)
315  [cell->vertex_index(
316  p4est_face_corners[sides[i].face]
317  [(limit - 1) ^ j])]
318  .insert(subdomain_ids[k]);
319  }
320  else
321  {
322  (*vertices_with_ghost_neighbors)
323  [cell->vertex_index(
324  p8est_face_corners[sides[i].face]
325  [(limit - 1) ^ j])]
326  .insert(subdomain_ids[k]);
327  }
328  }
329  }
330  }
331  }
332  }
333 
334  subids->elem_count = 0;
335  }
336  } // namespace
337 
338 
339  int (&functions<2>::quadrant_compare)(const void *v1, const void *v2) =
340  p4est_quadrant_compare;
341 
342  void (&functions<2>::quadrant_childrenv)(const types<2>::quadrant *q,
343  types<2>::quadrant c[]) =
344  p4est_quadrant_childrenv;
345 
346  int (&functions<2>::quadrant_overlaps_tree)(types<2>::tree * tree,
347  const types<2>::quadrant *q) =
348  p4est_quadrant_overlaps_tree;
349 
350  void (&functions<2>::quadrant_set_morton)(types<2>::quadrant *quadrant,
351  int level,
352  uint64_t id) =
353  p4est_quadrant_set_morton;
354 
355  int (&functions<2>::quadrant_is_equal)(const types<2>::quadrant *q1,
356  const types<2>::quadrant *q2) =
357  p4est_quadrant_is_equal;
358 
359  int (&functions<2>::quadrant_is_sibling)(const types<2>::quadrant *q1,
360  const types<2>::quadrant *q2) =
361  p4est_quadrant_is_sibling;
362 
363  int (&functions<2>::quadrant_is_ancestor)(const types<2>::quadrant *q1,
364  const types<2>::quadrant *q2) =
365  p4est_quadrant_is_ancestor;
366 
367  int (&functions<2>::quadrant_ancestor_id)(const types<2>::quadrant *q,
368  int level) =
369  p4est_quadrant_ancestor_id;
370 
371  int (&functions<2>::comm_find_owner)(types<2>::forest * p4est,
372  const types<2>::locidx which_tree,
373  const types<2>::quadrant *q,
374  const int guess) =
375  p4est_comm_find_owner;
376 
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;
382 
383  void (&functions<2>::connectivity_join_faces)(types<2>::connectivity *conn,
384  types<2>::topidx tree_left,
385  types<2>::topidx tree_right,
386  int face_left,
387  int face_right,
388  int orientation) =
389  p4est_connectivity_join_faces;
390 
391  void (&functions<2>::connectivity_destroy)(
392  p4est_connectivity_t *connectivity) = p4est_connectivity_destroy;
393 
394  types<2>::forest *(&functions<2>::new_forest)(
395  MPI_Comm mpicomm,
396  types<2>::connectivity *connectivity,
397  types<2>::locidx min_quadrants,
398  int min_level,
399  int fill_uniform,
400  std::size_t data_size,
401  p4est_init_t init_fn,
402  void * user_pointer) = p4est_new_ext;
403 
404  void (&functions<2>::destroy)(types<2>::forest *p4est) = p4est_destroy;
405 
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;
410 
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;
415 
416  void (&functions<2>::balance)(types<2>::forest * p4est,
417  types<2>::balance_type btype,
418  p4est_init_t init_fn) = p4est_balance;
419 
420  types<2>::gloidx (&functions<2>::partition)(types<2>::forest *p4est,
421  int partition_for_coarsening,
422  p4est_weight_t weight_fn) =
423  p4est_partition_ext;
424 
425  void (&functions<2>::save)(const char * filename,
426  types<2>::forest *p4est,
427  int save_data) = p4est_save;
428 
429  types<2>::forest *(&functions<2>::load_ext)(
430  const char * filename,
431  MPI_Comm mpicomm,
432  std::size_t data_size,
433  int load_data,
434  int autopartition,
435  int broadcasthead,
436  void * user_pointer,
437  types<2>::connectivity **p4est) = p4est_load_ext;
438 
439  int (&functions<2>::connectivity_save)(
440  const char * filename,
441  types<2>::connectivity *connectivity) = p4est_connectivity_save;
442 
443  int (&functions<2>::connectivity_is_valid)(
444  types<2>::connectivity *connectivity) = p4est_connectivity_is_valid;
445 
446  types<2>::connectivity *(&functions<2>::connectivity_load)(
447  const char * filename,
448  std::size_t *length) = p4est_connectivity_load;
449 
450  unsigned int (&functions<2>::checksum)(types<2>::forest *p4est) =
451  p4est_checksum;
452 
453  void (&functions<2>::vtk_write_file)(types<2>::forest *p4est,
454  p4est_geometry_t *,
455  const char *baseName) =
456  p4est_vtk_write_file;
457 
458  types<2>::ghost *(&functions<2>::ghost_new)(types<2>::forest * p4est,
459  types<2>::balance_type btype) =
460  p4est_ghost_new;
461 
462  void (&functions<2>::ghost_destroy)(types<2>::ghost *ghost) =
463  p4est_ghost_destroy;
464 
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;
469 
470  std::size_t (&functions<2>::forest_memory_used)(types<2>::forest *p4est) =
471  p4est_memory_used;
472 
473  std::size_t (&functions<2>::connectivity_memory_used)(
474  types<2>::connectivity *p4est) = p4est_connectivity_memory_used;
475 
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>
480  & tria,
481  typename ::internal::p4est::types<dim>::forest *parallel_forest,
482  typename ::internal::p4est::types<dim>::ghost * parallel_ghost)
483  {
484  std::map<unsigned int, std::set<::types::subdomain_id>>
485  vertices_with_ghost_neighbors;
486 
487  ::internal::p4est::FindGhosts<dim, spacedim> fg;
488  fg.subids = sc_array_new(sizeof(::types::subdomain_id));
489  fg.triangulation = &tria;
490  fg.vertices_with_ghost_neighbors = &vertices_with_ghost_neighbors;
491 
492  switch (dim)
493  {
494  case 2:
495  p4est_iterate(
496  reinterpret_cast<::internal::p4est::types<2>::forest *>(
497  parallel_forest),
498  reinterpret_cast<::internal::p4est::types<2>::ghost *>(
499  parallel_ghost),
500  static_cast<void *>(&fg),
501  nullptr,
502  find_ghosts_face<2, spacedim>,
503  find_ghosts_corner<2, spacedim>);
504  break;
505 
506  case 3:
507  p8est_iterate(
508  reinterpret_cast<::internal::p4est::types<3>::forest *>(
509  parallel_forest),
510  reinterpret_cast<::internal::p4est::types<3>::ghost *>(
511  parallel_ghost),
512  static_cast<void *>(&fg),
513  nullptr,
514  find_ghosts_face<3, 3>,
515  find_ghosts_edge<3, 3>,
516  find_ghosts_corner<3, 3>);
517  break;
518 
519  default:
520  Assert(false, ExcNotImplemented());
521  }
522 
523  sc_array_destroy(fg.subids);
524 
525  return vertices_with_ghost_neighbors;
526  }
527 
528  constexpr unsigned int functions<2>::max_level;
529 
530  void (&functions<2>::transfer_fixed)(const types<2>::gloidx *dest_gfq,
531  const types<2>::gloidx *src_gfq,
532  MPI_Comm mpicomm,
533  int tag,
534  void * dest_data,
535  const void * src_data,
536  std::size_t data_size) =
537  p4est_transfer_fixed;
538 
539  types<2>::transfer_context *(&functions<2>::transfer_fixed_begin)(
540  const types<2>::gloidx *dest_gfq,
541  const types<2>::gloidx *src_gfq,
542  MPI_Comm mpicomm,
543  int tag,
544  void * dest_data,
545  const void * src_data,
546  std::size_t data_size) = p4est_transfer_fixed_begin;
547 
548  void (&functions<2>::transfer_fixed_end)(types<2>::transfer_context *tc) =
549  p4est_transfer_fixed_end;
550 
551  void (&functions<2>::transfer_custom)(const types<2>::gloidx *dest_gfq,
552  const types<2>::gloidx *src_gfq,
553  MPI_Comm mpicomm,
554  int tag,
555  void * dest_data,
556  const int * dest_sizes,
557  const void * src_data,
558  const int * src_sizes) =
559  p4est_transfer_custom;
560 
561  types<2>::transfer_context *(&functions<2>::transfer_custom_begin)(
562  const types<2>::gloidx *dest_gfq,
563  const types<2>::gloidx *src_gfq,
564  MPI_Comm mpicomm,
565  int tag,
566  void * dest_data,
567  const int * dest_sizes,
568  const void * src_data,
569  const int * src_sizes) = p4est_transfer_custom_begin;
570 
571  void (&functions<2>::transfer_custom_end)(types<2>::transfer_context *tc) =
572  p4est_transfer_custom_end;
573 
574 
575 
576  int (&functions<3>::quadrant_compare)(const void *v1, const void *v2) =
577  p8est_quadrant_compare;
578 
579  void (&functions<3>::quadrant_childrenv)(const types<3>::quadrant *q,
580  types<3>::quadrant c[]) =
581  p8est_quadrant_childrenv;
582 
583  int (&functions<3>::quadrant_overlaps_tree)(types<3>::tree * tree,
584  const types<3>::quadrant *q) =
585  p8est_quadrant_overlaps_tree;
586 
587  void (&functions<3>::quadrant_set_morton)(types<3>::quadrant *quadrant,
588  int level,
589  uint64_t id) =
590  p8est_quadrant_set_morton;
591 
592  int (&functions<3>::quadrant_is_equal)(const types<3>::quadrant *q1,
593  const types<3>::quadrant *q2) =
594  p8est_quadrant_is_equal;
595 
596  int (&functions<3>::quadrant_is_sibling)(const types<3>::quadrant *q1,
597  const types<3>::quadrant *q2) =
598  p8est_quadrant_is_sibling;
599 
600  int (&functions<3>::quadrant_is_ancestor)(const types<3>::quadrant *q1,
601  const types<3>::quadrant *q2) =
602  p8est_quadrant_is_ancestor;
603 
604  int (&functions<3>::quadrant_ancestor_id)(const types<3>::quadrant *q,
605  int level) =
606  p8est_quadrant_ancestor_id;
607 
608  int (&functions<3>::comm_find_owner)(types<3>::forest * p4est,
609  const types<3>::locidx which_tree,
610  const types<3>::quadrant *q,
611  const int guess) =
612  p8est_comm_find_owner;
613 
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;
621 
622  void (&functions<3>::connectivity_destroy)(
623  p8est_connectivity_t *connectivity) = p8est_connectivity_destroy;
624 
625  void (&functions<3>::connectivity_join_faces)(types<3>::connectivity *conn,
626  types<3>::topidx tree_left,
627  types<3>::topidx tree_right,
628  int face_left,
629  int face_right,
630  int orientation) =
631  p8est_connectivity_join_faces;
632 
633  types<3>::forest *(&functions<3>::new_forest)(
634  MPI_Comm mpicomm,
635  types<3>::connectivity *connectivity,
636  types<3>::locidx min_quadrants,
637  int min_level,
638  int fill_uniform,
639  std::size_t data_size,
640  p8est_init_t init_fn,
641  void * user_pointer) = p8est_new_ext;
642 
643  void (&functions<3>::destroy)(types<3>::forest *p8est) = p8est_destroy;
644 
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;
649 
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;
654 
655  void (&functions<3>::balance)(types<3>::forest * p8est,
656  types<3>::balance_type btype,
657  p8est_init_t init_fn) = p8est_balance;
658 
659  types<3>::gloidx (&functions<3>::partition)(types<3>::forest *p8est,
660  int partition_for_coarsening,
661  p8est_weight_t weight_fn) =
662  p8est_partition_ext;
663 
664  void (&functions<3>::save)(const char * filename,
665  types<3>::forest *p4est,
666  int save_data) = p8est_save;
667 
668  types<3>::forest *(&functions<3>::load_ext)(
669  const char * filename,
670  MPI_Comm mpicomm,
671  std::size_t data_size,
672  int load_data,
673  int autopartition,
674  int broadcasthead,
675  void * user_pointer,
676  types<3>::connectivity **p4est) = p8est_load_ext;
677 
678  int (&functions<3>::connectivity_save)(
679  const char * filename,
680  types<3>::connectivity *connectivity) = p8est_connectivity_save;
681 
682  int (&functions<3>::connectivity_is_valid)(
683  types<3>::connectivity *connectivity) = p8est_connectivity_is_valid;
684 
685  types<3>::connectivity *(&functions<3>::connectivity_load)(
686  const char * filename,
687  std::size_t *length) = p8est_connectivity_load;
688 
689  unsigned int (&functions<3>::checksum)(types<3>::forest *p8est) =
690  p8est_checksum;
691 
692  void (&functions<3>::vtk_write_file)(types<3>::forest *p8est,
693  p8est_geometry_t *,
694  const char *baseName) =
695  p8est_vtk_write_file;
696 
697  types<3>::ghost *(&functions<3>::ghost_new)(types<3>::forest * p4est,
698  types<3>::balance_type btype) =
699  p8est_ghost_new;
700 
701  void (&functions<3>::ghost_destroy)(types<3>::ghost *ghost) =
702  p8est_ghost_destroy;
703 
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;
708 
709  std::size_t (&functions<3>::forest_memory_used)(types<3>::forest *p4est) =
710  p8est_memory_used;
711 
712  std::size_t (&functions<3>::connectivity_memory_used)(
713  types<3>::connectivity *p4est) = p8est_connectivity_memory_used;
714 
715  constexpr unsigned int functions<3>::max_level;
716 
717  void (&functions<3>::transfer_fixed)(const types<3>::gloidx *dest_gfq,
718  const types<3>::gloidx *src_gfq,
719  MPI_Comm mpicomm,
720  int tag,
721  void * dest_data,
722  const void * src_data,
723  std::size_t data_size) =
724  p8est_transfer_fixed;
725 
726  types<3>::transfer_context *(&functions<3>::transfer_fixed_begin)(
727  const types<3>::gloidx *dest_gfq,
728  const types<3>::gloidx *src_gfq,
729  MPI_Comm mpicomm,
730  int tag,
731  void * dest_data,
732  const void * src_data,
733  std::size_t data_size) = p8est_transfer_fixed_begin;
734 
735  void (&functions<3>::transfer_fixed_end)(types<3>::transfer_context *tc) =
736  p8est_transfer_fixed_end;
737 
738  void (&functions<3>::transfer_custom)(const types<3>::gloidx *dest_gfq,
739  const types<3>::gloidx *src_gfq,
740  MPI_Comm mpicomm,
741  int tag,
742  void * dest_data,
743  const int * dest_sizes,
744  const void * src_data,
745  const int * src_sizes) =
746  p8est_transfer_custom;
747 
748  types<3>::transfer_context *(&functions<3>::transfer_custom_begin)(
749  const types<3>::gloidx *dest_gfq,
750  const types<3>::gloidx *src_gfq,
751  MPI_Comm mpicomm,
752  int tag,
753  void * dest_data,
754  const int * dest_sizes,
755  const void * src_data,
756  const int * src_sizes) = p8est_transfer_custom_begin;
757 
758  void (&functions<3>::transfer_custom_end)(types<3>::transfer_context *tc) =
759  p8est_transfer_custom_end;
760 
761 
762 
763  template <int dim>
764  void
765  init_quadrant_children(
766  const typename types<dim>::quadrant &p4est_cell,
767  typename types<dim>::quadrant (
768  &p4est_children)[::GeometryInfo<dim>::max_children_per_cell])
769  {
770  for (unsigned int c = 0;
771  c < ::GeometryInfo<dim>::max_children_per_cell;
772  ++c)
773  switch (dim)
774  {
775  case 2:
776  P4EST_QUADRANT_INIT(&p4est_children[c]);
777  break;
778  case 3:
779  P8EST_QUADRANT_INIT(&p4est_children[c]);
780  break;
781  default:
782  Assert(false, ExcNotImplemented());
783  }
784 
785 
786  functions<dim>::quadrant_childrenv(&p4est_cell, p4est_children);
787  }
788 
789  template <int dim>
790  void
791  init_coarse_quadrant(typename types<dim>::quadrant &quad)
792  {
793  switch (dim)
794  {
795  case 2:
796  P4EST_QUADRANT_INIT(&quad);
797  break;
798  case 3:
799  P8EST_QUADRANT_INIT(&quad);
800  break;
801  default:
802  Assert(false, ExcNotImplemented());
803  }
804  functions<dim>::quadrant_set_morton(&quad,
805  /*level=*/0,
806  /*index=*/0);
807  }
808 
809  template <int dim>
810  bool
811  quadrant_is_equal(const typename types<dim>::quadrant &q1,
812  const typename types<dim>::quadrant &q2)
813  {
814  return functions<dim>::quadrant_is_equal(&q1, &q2);
815  }
816 
817 
818 
819  template <int dim>
820  bool
821  quadrant_is_ancestor(const typename types<dim>::quadrant &q1,
822  const typename types<dim>::quadrant &q2)
823  {
824  return functions<dim>::quadrant_is_ancestor(&q1, &q2);
825  }
826 
827  template <int dim>
828  bool
829  tree_exists_locally(const typename types<dim>::forest *parallel_forest,
830  const typename types<dim>::topidx coarse_grid_cell)
831  {
832  Assert(coarse_grid_cell < parallel_forest->connectivity->num_trees,
833  ExcInternalError());
834  return ((coarse_grid_cell >= parallel_forest->first_local_tree) &&
835  (coarse_grid_cell <= parallel_forest->last_local_tree));
836  }
837  } // namespace p4est
838 } // namespace internal
839 
840 #endif // DEAL_II_WITH_P4EST
841 
842 /*-------------- Explicit Instantiations -------------------------------*/
843 #include "p4est_wrappers.inst"
844 
845 
846 DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcMessage(std::string arg1)
unsigned int subdomain_id
Definition: types.h:43
#define Assert(cond, exc)
Definition: exceptions.h:1407
unsigned int global_dof_index
Definition: types.h:89
static ::ExceptionBase & ExcNotImplemented()
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
static ::ExceptionBase & ExcInternalError()