Reference documentation for deal.II version Git 1dc1051882 2021-04-22 23:57:03 +0200
\(\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\}}\)
p4est_wrappers.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2016 - 2020 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 #ifndef dealii_p4est_wrappers_h
17 #define dealii_p4est_wrappers_h
18 
19 #include <deal.II/base/config.h>
20 
22 
23 #ifdef DEAL_II_WITH_P4EST
24 # include <p4est_bits.h>
25 # include <p4est_communication.h>
26 # include <p4est_extended.h>
27 # include <p4est_ghost.h>
28 # include <p4est_iterate.h>
29 # include <p4est_vtk.h>
30 # include <p8est_bits.h>
31 # include <p8est_communication.h>
32 # include <p8est_extended.h>
33 # include <p8est_ghost.h>
34 # include <p8est_iterate.h>
35 # include <p8est_vtk.h>
36 
37 # include <map>
38 # include <set>
39 
40 
42 
43 // Forward declaration
44 # ifndef DOXYGEN
45 namespace parallel
46 {
47  namespace distributed
48  {
49  template <int dim, int spacedim>
50  class Triangulation;
51  }
52 } // namespace parallel
53 # endif
54 
55 namespace internal
56 {
57  namespace p4est
58  {
66  template <int>
67  struct types;
68 
69  // these struct mimics p4est for 1D
70  template <>
71  struct types<1>
72  {
73  // id of a quadrant is an integeger
74  using quadrant = int;
75 
76  // maximum number of children
77  static const int max_n_child_indices_bits = 27;
78 
79  // number of bits the data type of id has
80  static const int n_bits = std::numeric_limits<quadrant>::digits;
81  };
82 
83  template <>
84  struct types<2>
85  {
86  using connectivity = p4est_connectivity_t;
87  using forest = p4est_t;
88  using tree = p4est_tree_t;
89  using quadrant = p4est_quadrant_t;
90  using topidx = p4est_topidx_t;
91  using locidx = p4est_locidx_t;
92  using gloidx = p4est_gloidx_t;
93  using balance_type = p4est_connect_type_t;
94  using ghost = p4est_ghost_t;
95  using transfer_context = p4est_transfer_context_t;
96  };
97 
98  template <>
99  struct types<3>
100  {
101  using connectivity = p8est_connectivity_t;
102  using forest = p8est_t;
103  using tree = p8est_tree_t;
104  using quadrant = p8est_quadrant_t;
105  using topidx = p4est_topidx_t;
106  using locidx = p4est_locidx_t;
107  using gloidx = p4est_gloidx_t;
108  using balance_type = p8est_connect_type_t;
109  using ghost = p8est_ghost_t;
110  using transfer_context = p8est_transfer_context_t;
111  };
112 
113 
114 
122  template <int dim>
123  struct functions;
124 
125  template <>
126  struct functions<2>
127  {
128  static int (&quadrant_compare)(const void *v1, const void *v2);
129 
130  static void (&quadrant_childrenv)(const types<2>::quadrant *q,
131  types<2>::quadrant c[]);
132 
133  static int (&quadrant_overlaps_tree)(types<2>::tree * tree,
134  const types<2>::quadrant *q);
135 
136  static void (&quadrant_set_morton)(types<2>::quadrant *quadrant,
137  int level,
138  uint64_t id);
139 
141  const types<2>::quadrant *q2);
142 
143  static int (&quadrant_is_sibling)(const types<2>::quadrant *q1,
144  const types<2>::quadrant *q2);
145 
147  const types<2>::quadrant *q2);
148 
149  static int (&quadrant_ancestor_id)(const types<2>::quadrant *q,
150  int level);
151 
152  static int (&comm_find_owner)(types<2>::forest * p4est,
153  const types<2>::locidx which_tree,
154  const types<2>::quadrant *q,
155  const int guess);
156 
158  types<2>::topidx num_vertices,
159  types<2>::topidx num_trees,
160  types<2>::topidx num_corners,
161  types<2>::topidx num_vtt);
162 
164  types<2>::topidx num_vertices,
165  types<2>::topidx num_trees,
166  types<2>::topidx num_corners,
167  const double * vertices,
168  const types<2>::topidx *ttv,
169  const types<2>::topidx *ttt,
170  const int8_t * ttf,
171  const types<2>::topidx *ttc,
172  const types<2>::topidx *coff,
173  const types<2>::topidx *ctt,
174  const int8_t * ctc);
175 
176  static void (&connectivity_join_faces)(types<2>::connectivity *conn,
177  types<2>::topidx tree_left,
178  types<2>::topidx tree_right,
179  int face_left,
180  int face_right,
181  int orientation);
182 
183 
184 
185  static void (&connectivity_destroy)(p4est_connectivity_t *connectivity);
186 
188  MPI_Comm mpicomm,
189  types<2>::connectivity *connectivity,
190  types<2>::locidx min_quadrants,
191  int min_level,
192  int fill_uniform,
193  std::size_t data_size,
194  p4est_init_t init_fn,
195  void * user_pointer);
196 
198  int copy_data);
199 
200  static void (&destroy)(types<2>::forest *p4est);
201 
202  static void (&refine)(types<2>::forest *p4est,
203  int refine_recursive,
204  p4est_refine_t refine_fn,
205  p4est_init_t init_fn);
206 
207  static void (&coarsen)(types<2>::forest *p4est,
208  int coarsen_recursive,
209  p4est_coarsen_t coarsen_fn,
210  p4est_init_t init_fn);
211 
212  static void (&balance)(types<2>::forest * p4est,
214  p4est_init_t init_fn);
215 
217  int partition_for_coarsening,
218  p4est_weight_t weight_fn);
219 
220  static void (&save)(const char * filename,
221  types<2>::forest *p4est,
222  int save_data);
223 
224  static types<2>::forest *(&load_ext)(const char *filename,
225  MPI_Comm mpicomm,
226  std::size_t data_size,
227  int load_data,
228  int autopartition,
229  int broadcasthead,
230  void * user_pointer,
231  types<2>::connectivity **p4est);
232 
233  static int (&connectivity_save)(const char * filename,
234  types<2>::connectivity *connectivity);
235 
236  static int (&connectivity_is_valid)(types<2>::connectivity *connectivity);
237 
238  static types<2>::connectivity *(&connectivity_load)(const char * filename,
239  std::size_t *length);
240 
241  static unsigned int (&checksum)(types<2>::forest *p4est);
242 
243  static void (&vtk_write_file)(types<2>::forest *p4est,
244  p4est_geometry_t *,
245  const char *baseName);
246 
248  types<2>::balance_type btype);
249 
250  static void (&ghost_destroy)(types<2>::ghost *ghost);
251 
252  static void (&reset_data)(types<2>::forest *p4est,
253  std::size_t data_size,
254  p4est_init_t init_fn,
255  void * user_pointer);
256 
257  static std::size_t (&forest_memory_used)(types<2>::forest *p4est);
258 
259  static std::size_t (&connectivity_memory_used)(
260  types<2>::connectivity *p4est);
261 
262  template <int spacedim>
263  static void
264  iterate(::internal::p4est::types<2>::forest *parallel_forest,
265  ::internal::p4est::types<2>::ghost * parallel_ghost,
266  void * user_data);
267 
268  static constexpr unsigned int max_level = P4EST_MAXLEVEL;
269 
270  static void (&transfer_fixed)(const types<2>::gloidx *dest_gfq,
271  const types<2>::gloidx *src_gfq,
272  MPI_Comm mpicomm,
273  int tag,
274  void * dest_data,
275  const void * src_data,
276  std::size_t data_size);
277 
279  const types<2>::gloidx *dest_gfq,
280  const types<2>::gloidx *src_gfq,
281  MPI_Comm mpicomm,
282  int tag,
283  void * dest_data,
284  const void * src_data,
285  std::size_t data_size);
286 
287  static void (&transfer_fixed_end)(types<2>::transfer_context *tc);
288 
289  static void (&transfer_custom)(const types<2>::gloidx *dest_gfq,
290  const types<2>::gloidx *src_gfq,
291  MPI_Comm mpicomm,
292  int tag,
293  void * dest_data,
294  const int * dest_sizes,
295  const void * src_data,
296  const int * src_sizes);
297 
299  const types<2>::gloidx *dest_gfq,
300  const types<2>::gloidx *src_gfq,
301  MPI_Comm mpicomm,
302  int tag,
303  void * dest_data,
304  const int * dest_sizes,
305  const void * src_data,
306  const int * src_sizes);
307 
308  static void (&transfer_custom_end)(types<2>::transfer_context *tc);
309  };
310 
311 
312  template <>
313  struct functions<3>
314  {
315  static int (&quadrant_compare)(const void *v1, const void *v2);
316 
317  static void (&quadrant_childrenv)(const types<3>::quadrant *q,
318  types<3>::quadrant c[]);
319 
320  static int (&quadrant_overlaps_tree)(types<3>::tree * tree,
321  const types<3>::quadrant *q);
322 
323  static void (&quadrant_set_morton)(types<3>::quadrant *quadrant,
324  int level,
325  uint64_t id);
326 
328  const types<3>::quadrant *q2);
329 
330  static int (&quadrant_is_sibling)(const types<3>::quadrant *q1,
331  const types<3>::quadrant *q2);
332 
334  const types<3>::quadrant *q2);
335  static int (&quadrant_ancestor_id)(const types<3>::quadrant *q,
336  int level);
337 
338  static int (&comm_find_owner)(types<3>::forest * p4est,
339  const types<3>::locidx which_tree,
340  const types<3>::quadrant *q,
341  const int guess);
342 
344  types<3>::topidx num_vertices,
345  types<3>::topidx num_trees,
346  types<3>::topidx num_edges,
347  types<3>::topidx num_ett,
348  types<3>::topidx num_corners,
349  types<3>::topidx num_ctt);
350 
352  types<3>::topidx num_vertices,
353  types<3>::topidx num_trees,
354  types<3>::topidx num_edges,
355  types<3>::topidx num_corners,
356  const double * vertices,
357  const types<3>::topidx *ttv,
358  const types<3>::topidx *ttt,
359  const int8_t * ttf,
360  const types<3>::topidx *tte,
361  const types<3>::topidx *eoff,
362  const types<3>::topidx *ett,
363  const int8_t * ete,
364  const types<3>::topidx *ttc,
365  const types<3>::topidx *coff,
366  const types<3>::topidx *ctt,
367  const int8_t * ctc);
368 
369  static void (&connectivity_join_faces)(types<3>::connectivity *conn,
370  types<3>::topidx tree_left,
371  types<3>::topidx tree_right,
372  int face_left,
373  int face_right,
374  int orientation);
375 
376  static void (&connectivity_destroy)(p8est_connectivity_t *connectivity);
377 
379  MPI_Comm mpicomm,
380  types<3>::connectivity *connectivity,
381  types<3>::locidx min_quadrants,
382  int min_level,
383  int fill_uniform,
384  std::size_t data_size,
385  p8est_init_t init_fn,
386  void * user_pointer);
387 
389  int copy_data);
390 
391  static void (&destroy)(types<3>::forest *p8est);
392 
393  static void (&refine)(types<3>::forest *p8est,
394  int refine_recursive,
395  p8est_refine_t refine_fn,
396  p8est_init_t init_fn);
397 
398  static void (&coarsen)(types<3>::forest *p8est,
399  int coarsen_recursive,
400  p8est_coarsen_t coarsen_fn,
401  p8est_init_t init_fn);
402 
403  static void (&balance)(types<3>::forest * p8est,
405  p8est_init_t init_fn);
406 
408  int partition_for_coarsening,
409  p8est_weight_t weight_fn);
410 
411  static void (&save)(const char * filename,
412  types<3>::forest *p4est,
413  int save_data);
414 
415  static types<3>::forest *(&load_ext)(const char *filename,
416  MPI_Comm mpicomm,
417  std::size_t data_size,
418  int load_data,
419  int autopartition,
420  int broadcasthead,
421  void * user_pointer,
422  types<3>::connectivity **p4est);
423 
424  static int (&connectivity_save)(const char * filename,
425  types<3>::connectivity *connectivity);
426 
427  static int (&connectivity_is_valid)(types<3>::connectivity *connectivity);
428 
429  static types<3>::connectivity *(&connectivity_load)(const char * filename,
430  std::size_t *length);
431 
432  static unsigned int (&checksum)(types<3>::forest *p8est);
433 
434  static void (&vtk_write_file)(types<3>::forest *p8est,
435  p8est_geometry_t *,
436  const char *baseName);
438  types<3>::balance_type btype);
439 
440  static void (&ghost_destroy)(types<3>::ghost *ghost);
441 
442  static void (&reset_data)(types<3>::forest *p4est,
443  std::size_t data_size,
444  p8est_init_t init_fn,
445  void * user_pointer);
446 
447  static std::size_t (&forest_memory_used)(types<3>::forest *p4est);
448 
449  static std::size_t (&connectivity_memory_used)(
450  types<3>::connectivity *p4est);
451 
452  static constexpr unsigned int max_level = P8EST_MAXLEVEL;
453 
454  static void (&transfer_fixed)(const types<3>::gloidx *dest_gfq,
455  const types<3>::gloidx *src_gfq,
456  MPI_Comm mpicomm,
457  int tag,
458  void * dest_data,
459  const void * src_data,
460  std::size_t data_size);
461 
463  const types<3>::gloidx *dest_gfq,
464  const types<3>::gloidx *src_gfq,
465  MPI_Comm mpicomm,
466  int tag,
467  void * dest_data,
468  const void * src_data,
469  std::size_t data_size);
470 
471  static void (&transfer_fixed_end)(types<3>::transfer_context *tc);
472 
473  static void (&transfer_custom)(const types<3>::gloidx *dest_gfq,
474  const types<3>::gloidx *src_gfq,
475  MPI_Comm mpicomm,
476  int tag,
477  void * dest_data,
478  const int * dest_sizes,
479  const void * src_data,
480  const int * src_sizes);
481 
483  const types<3>::gloidx *dest_gfq,
484  const types<3>::gloidx *src_gfq,
485  MPI_Comm mpicomm,
486  int tag,
487  void * dest_data,
488  const int * dest_sizes,
489  const void * src_data,
490  const int * src_sizes);
491 
492  static void (&transfer_custom_end)(types<3>::transfer_context *tc);
493  };
494 
495 
496 
503  template <int dim>
504  struct iter;
505 
506  template <>
507  struct iter<2>
508  {
509  using corner_info = p4est_iter_corner_info_t;
510  using corner_side = p4est_iter_corner_side_t;
511  using corner_iter = p4est_iter_corner_t;
512  using face_info = p4est_iter_face_info_t;
513  using face_side = p4est_iter_face_side_t;
514  using face_iter = p4est_iter_face_t;
515  };
516 
517  template <>
518  struct iter<3>
519  {
520  using corner_info = p8est_iter_corner_info_t;
521  using corner_side = p8est_iter_corner_side_t;
522  using corner_iter = p8est_iter_corner_t;
523  using edge_info = p8est_iter_edge_info_t;
524  using edge_side = p8est_iter_edge_side_t;
525  using edge_iter = p8est_iter_edge_t;
526  using face_info = p8est_iter_face_info_t;
527  using face_side = p8est_iter_face_side_t;
528  using face_iter = p8est_iter_face_t;
529  };
530 
531 
532 
537  template <int dim>
538  void
540  const typename types<dim>::quadrant &p4est_cell,
541  typename types<dim>::quadrant (
542  &p4est_children)[::GeometryInfo<dim>::max_children_per_cell]);
543 
544 
545 
549  template <int dim>
550  void
552 
553 
554 
558  template <int dim>
559  bool
560  quadrant_is_equal(const typename types<dim>::quadrant &q1,
561  const typename types<dim>::quadrant &q2);
562 
563 
564 
568  template <int dim>
569  bool
570  quadrant_is_ancestor(const typename types<dim>::quadrant &q1,
571  const typename types<dim>::quadrant &q2);
572 
573 
574 
578  template <int dim>
579  bool
580  tree_exists_locally(const typename types<dim>::forest *parallel_forest,
581  const typename types<dim>::topidx coarse_grid_cell);
582 
583 
587  template <int dim>
588  typename types<dim>::connectivity *
589  copy_connectivity(const typename types<dim>::connectivity *connectivity);
590 
591 # ifndef DOXYGEN
592  template <>
593  typename types<2>::connectivity *
594  copy_connectivity<2>(const typename types<2>::connectivity *connectivity);
595 
596  template <>
597  typename types<3>::connectivity *
598  copy_connectivity<3>(const typename types<3>::connectivity *connectivity);
599 # endif
600  } // namespace p4est
601 } // namespace internal
602 
604 
605 #endif // DEAL_II_WITH_P4EST
606 
607 #endif // dealii_p4est_wrappers_h
static types< 2 >::forest *(&) new_forest(MPI_Comm mpicomm, types< 2 >::connectivity *connectivity, types< 2 >::locidx min_quadrants, int min_level, int fill_uniform, std::size_t data_size, p4est_init_t init_fn, void *user_pointer)
p8est_transfer_context_t transfer_context
types< 2 >::connectivity * copy_connectivity< 2 >(const typename types< 2 >::connectivity *connectivity)
p8est_iter_corner_side_t corner_side
bool tree_exists_locally(const typename types< dim >::forest *parallel_forest, const typename types< dim >::topidx coarse_grid_cell)
p8est_iter_face_t face_iter
p4est_iter_face_info_t face_info
p4est_iter_corner_info_t corner_info
static types< 2 >::connectivity *(&) connectivity_new_copy(types< 2 >::topidx num_vertices, types< 2 >::topidx num_trees, types< 2 >::topidx num_corners, const double *vertices, const types< 2 >::topidx *ttv, const types< 2 >::topidx *ttt, const int8_t *ttf, const types< 2 >::topidx *ttc, const types< 2 >::topidx *coff, const types< 2 >::topidx *ctt, const int8_t *ctc)
static types< 3 >::connectivity *(&) connectivity_new_copy(types< 3 >::topidx num_vertices, types< 3 >::topidx num_trees, types< 3 >::topidx num_edges, types< 3 >::topidx num_corners, const double *vertices, const types< 3 >::topidx *ttv, const types< 3 >::topidx *ttt, const int8_t *ttf, const types< 3 >::topidx *tte, const types< 3 >::topidx *eoff, const types< 3 >::topidx *ett, const int8_t *ete, const types< 3 >::topidx *ttc, const types< 3 >::topidx *coff, const types< 3 >::topidx *ctt, const int8_t *ctc)
p4est_connect_type_t balance_type
static types< 2 >::forest *(&) copy_forest(types< 2 >::forest *input, int copy_data)
static types< 3 >::forest *(&) new_forest(MPI_Comm mpicomm, types< 3 >::connectivity *connectivity, types< 3 >::locidx min_quadrants, int min_level, int fill_uniform, std::size_t data_size, p8est_init_t init_fn, void *user_pointer)
p4est_iter_face_t face_iter
static types< 3 >::forest *(&) load_ext(const char *filename, MPI_Comm mpicomm, std::size_t data_size, int load_data, int autopartition, int broadcasthead, void *user_pointer, types< 3 >::connectivity **p4est)
p8est_iter_face_side_t face_side
static types< 2 >::transfer_context *(&) transfer_custom_begin(const types< 2 >::gloidx *dest_gfq, const types< 2 >::gloidx *src_gfq, MPI_Comm mpicomm, int tag, void *dest_data, const int *dest_sizes, const void *src_data, const int *src_sizes)
p8est_connect_type_t balance_type
p4est_iter_face_side_t face_side
bool quadrant_is_equal(const typename types< dim >::quadrant &q1, const typename types< dim >::quadrant &q2)
p4est_connectivity_t connectivity
void init_quadrant_children(const typename types< dim >::quadrant &p4est_cell, typename types< dim >::quadrant(&p4est_children)[::GeometryInfo< dim >::max_children_per_cell])
void init_coarse_quadrant(typename types< dim >::quadrant &quad)
static types< 2 >::transfer_context *(&) transfer_fixed_begin(const types< 2 >::gloidx *dest_gfq, const types< 2 >::gloidx *src_gfq, MPI_Comm mpicomm, int tag, void *dest_data, const void *src_data, std::size_t data_size)
p8est_connectivity_t connectivity
Definition: types.h:31
void partition(const SparsityPattern &sparsity_pattern, const unsigned int n_partitions, std::vector< unsigned int > &partition_indices, const Partitioner partitioner=Partitioner::metis)
p4est_iter_corner_side_t corner_side
static types< 3 >::connectivity *(&) connectivity_new(types< 3 >::topidx num_vertices, types< 3 >::topidx num_trees, types< 3 >::topidx num_edges, types< 3 >::topidx num_ett, types< 3 >::topidx num_corners, types< 3 >::topidx num_ctt)
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:396
unsigned int level
Definition: grid_out.cc:4590
bool quadrant_is_ancestor(const typename types< dim >::quadrant &q1, const typename types< dim >::quadrant &q2)
p8est_iter_corner_info_t corner_info
Point< 3 > vertices[4]
static types< 3 >::transfer_context *(&) transfer_fixed_begin(const types< 3 >::gloidx *dest_gfq, const types< 3 >::gloidx *src_gfq, MPI_Comm mpicomm, int tag, void *dest_data, const void *src_data, std::size_t data_size)
static types< 3 >::connectivity *(&) connectivity_load(const char *filename, std::size_t *length)
static types< 2 >::forest *(&) load_ext(const char *filename, MPI_Comm mpicomm, std::size_t data_size, int load_data, int autopartition, int broadcasthead, void *user_pointer, types< 2 >::connectivity **p4est)
void coarsen(Triangulation< dim, spacedim > &tria, const Vector< Number > &criteria, const double threshold)
p8est_iter_edge_info_t edge_info
static types< 3 >::forest *(&) copy_forest(types< 3 >::forest *input, int copy_data)
static types< 3 >::ghost *(&) ghost_new(types< 3 >::forest *p4est, types< 3 >::balance_type btype)
types< 3 >::connectivity * copy_connectivity< 3 >(const typename types< 3 >::connectivity *connectivity)
p8est_iter_corner_t corner_iter
static types< 2 >::ghost *(&) ghost_new(types< 2 >::forest *p4est, types< 2 >::balance_type btype)
types< dim >::connectivity * copy_connectivity(const typename types< dim >::connectivity *connectivity)
p8est_iter_face_info_t face_info
static types< 2 >::connectivity *(&) connectivity_load(const char *filename, std::size_t *length)
p4est_transfer_context_t transfer_context
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:395
static types< 3 >::transfer_context *(&) transfer_custom_begin(const types< 3 >::gloidx *dest_gfq, const types< 3 >::gloidx *src_gfq, MPI_Comm mpicomm, int tag, void *dest_data, const int *dest_sizes, const void *src_data, const int *src_sizes)
static types< 2 >::connectivity *(&) connectivity_new(types< 2 >::topidx num_vertices, types< 2 >::topidx num_trees, types< 2 >::topidx num_corners, types< 2 >::topidx num_vtt)
void refine(Triangulation< dim, spacedim > &tria, const Vector< Number > &criteria, const double threshold, const unsigned int max_to_mark=numbers::invalid_unsigned_int)
p8est_iter_edge_side_t edge_side
p8est_iter_edge_t edge_iter
p4est_iter_corner_t corner_iter