Reference documentation for deal.II version 9.3.3
\(\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 - 2021 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
45namespace parallel
46{
47 namespace distributed
48 {
49 template <int dim, int spacedim>
50 class Triangulation;
51 }
52} // namespace parallel
53# endif
54
55namespace 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,
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
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)(
261
262 template <int spacedim>
263 static void
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,
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);
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)(
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 (
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
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
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:402
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:403
Point< 3 > vertices[4]
unsigned int level
Definition: grid_out.cc:4590
const unsigned int v1
Definition: grid_tools.cc:963
void refine(Triangulation< dim, spacedim > &tria, const Vector< Number > &criteria, const double threshold, const unsigned int max_to_mark=numbers::invalid_unsigned_int)
void coarsen(Triangulation< dim, spacedim > &tria, const Vector< Number > &criteria, const double threshold)
void partition(const SparsityPattern &sparsity_pattern, const unsigned int n_partitions, std::vector< unsigned int > &partition_indices, const Partitioner partitioner=Partitioner::metis)
void init_coarse_quadrant(typename types< dim >::quadrant &quad)
types< 2 >::connectivity * copy_connectivity< 2 >(const typename types< 2 >::connectivity *connectivity)
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)
bool tree_exists_locally(const typename types< dim >::forest *parallel_forest, const typename types< dim >::topidx coarse_grid_cell)
types< dim >::connectivity * copy_connectivity(const typename types< dim >::connectivity *connectivity)
Definition: types.h:32
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)
static types< 2 >::forest *(&) copy_forest(types< 2 >::forest *input, int copy_data)
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)
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)
static types< 2 >::connectivity *(&) connectivity_load(const char *filename, std::size_t *length)
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 void iterate(::internal::p4est::types< 2 >::forest *parallel_forest, ::internal::p4est::types< 2 >::ghost *parallel_ghost, void *user_data)
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)
static types< 2 >::ghost *(&) ghost_new(types< 2 >::forest *p4est, types< 2 >::balance_type btype)
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)
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)
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)
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< 3 >::ghost *(&) ghost_new(types< 3 >::forest *p4est, types< 3 >::balance_type btype)
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)
static types< 3 >::connectivity *(&) connectivity_load(const char *filename, std::size_t *length)
static types< 3 >::forest *(&) copy_forest(types< 3 >::forest *input, int copy_data)
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 >::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_corner_t corner_iter
p4est_iter_corner_info_t corner_info
p4est_iter_face_info_t face_info
p4est_iter_face_side_t face_side
p4est_iter_corner_side_t corner_side
p4est_iter_face_t face_iter
p8est_iter_corner_t corner_iter
p8est_iter_corner_info_t corner_info
p8est_iter_face_info_t face_info
p8est_iter_edge_side_t edge_side
p8est_iter_edge_t edge_iter
p8est_iter_edge_info_t edge_info
p8est_iter_face_side_t face_side
p8est_iter_face_t face_iter
p8est_iter_corner_side_t corner_side
p4est_connectivity_t connectivity
p4est_connect_type_t balance_type
p4est_transfer_context_t transfer_context
p8est_connectivity_t connectivity
p8est_transfer_context_t transfer_context
p8est_connect_type_t balance_type