p4est 2.8.6
p4est is a software library for parallel adaptive mesh refinement.
Data Structures | Macros | Typedefs | Enumerations | Functions | Variables
p8est_connectivity.h File Reference

The connectivity defines the coarse topology of the forest. More...

#include <sc_io.h>
#include <p4est_base.h>
Include dependency graph for p8est_connectivity.h:
This graph shows which files directly or indirectly include this file:

Go to the source code of this file.

Data Structures

struct  p8est_connectivity
 This structure holds the 3D inter-tree connectivity information. More...
 
struct  p8est_edge_transform_t
 Generic interface for transformations between a tree and any of its edge. More...
 
struct  p8est_edge_info_t
 Information about the neighbors of an edge. More...
 
struct  p8est_corner_transform_t
 Generic interface for transformations between a tree and any of its corner. More...
 
struct  p8est_corner_info_t
 Information about the neighbors of a corner. More...
 
struct  p8est_neighbor_transform_t
 Generic interface for transformations between a tree and any of its neighbors. More...
 

Macros

#define P8EST_DIM   3
 The spatial dimension.
 
#define P8EST_FACES   (2 * P8EST_DIM)
 The number of faces of an octant. More...
 
#define P8EST_CHILDREN   8
 The number of children of an octant. More...
 
#define P8EST_HALF   (P8EST_CHILDREN / 2)
 The number of children/corners touching one face.
 
#define P8EST_EDGES   12
 The number of edges around an octant.
 
#define P8EST_INSUL   27
 The size of insulation layer.
 
#define P8EST_ONLY_P8_LAND(x)   && (x)
 Only use logical AND term in 3D.
 
#define P8EST_ONLY_P8_COMMA(x)   , (x)
 Only use comma and expression in 3D.
 
#define P8EST_DIM_POW(a)   ((a) * (a) * (a))
 Exponentiate with dimension.
 
#define P8EST_FTRANSFORM   9
 size of face transformation encoding
 
#define P8EST_STRING   "p8est"
 p8est identification string
 
#define P8EST_ONDISK_FORMAT   0x3000009
 Increase this number whenever the on-disk format for p8est_connectivity, p8est, or any other 3D data structure changes. More...
 

Typedefs

typedef struct p8est_connectivity p8est_connectivity_t
 This structure holds the 3D inter-tree connectivity information. More...
 

Enumerations

enum  p8est_connect_type_t {
  P8EST_CONNECT_SELF = 30 ,
  P8EST_CONNECT_FACE = 31 ,
  P8EST_CONNECT_EDGE = 32 ,
  P8EST_CONNECT_ALMOST = P8EST_CONNECT_EDGE ,
  P8EST_CONNECT_CORNER = 33 ,
  P8EST_CONNECT_FULL = P8EST_CONNECT_CORNER
}
 Characterize a type of adjacency. More...
 
enum  p8est_connectivity_encode_t {
  P8EST_CONN_ENCODE_NONE = SC_IO_ENCODE_NONE ,
  P8EST_CONN_ENCODE_LAST
}
 Typedef for serialization method. More...
 

Functions

int p8est_connect_type_int (p8est_connect_type_t btype)
 Convert the p8est_connect_type_t into a number. More...
 
const char * p8est_connect_type_string (p8est_connect_type_t btype)
 Convert the p8est_connect_type_t into a const string. More...
 
size_t p8est_connectivity_memory_used (p8est_connectivity_t *conn)
 Calculate memory usage of a connectivity structure. More...
 
void p8est_neighbor_transform_coordinates (const p8est_neighbor_transform_t *nt, const p4est_qcoord_t self_coords[P8EST_DIM], p4est_qcoord_t neigh_coords[P8EST_DIM])
 Transform from self's coordinate system to neighbor's coordinate system. More...
 
void p8est_neighbor_transform_coordinates_reverse (const p8est_neighbor_transform_t *nt, const p4est_qcoord_t neigh_coords[P8EST_DIM], p4est_qcoord_t self_coords[P8EST_DIM])
 Transform from neighbor's coordinate system to self's coordinate system. More...
 
void p8est_connectivity_get_neighbor_transforms (p8est_connectivity_t *conn, p4est_topidx_t tree_id, p8est_connect_type_t boundary_type, int boundary_index, sc_array_t *neighbor_transform_array)
 Fill an array with the neighbor transforms based on a specific boundary type. More...
 
int p8est_connectivity_face_neighbor_corner_set (int c, int f, int nf, int set)
 Transform a corner across one of the adjacent faces into a neighbor tree. More...
 
int p8est_connectivity_face_neighbor_face_corner (int fc, int f, int nf, int o)
 Transform a face corner across one of the adjacent faces into a neighbor tree. More...
 
int p8est_connectivity_face_neighbor_corner (int c, int f, int nf, int o)
 Transform a corner across one of the adjacent faces into a neighbor tree. More...
 
int p8est_connectivity_face_neighbor_face_edge (int fe, int f, int nf, int o)
 Transform a face-edge across one of the adjacent faces into a neighbor tree. More...
 
int p8est_connectivity_face_neighbor_edge (int e, int f, int nf, int o)
 Transform an edge across one of the adjacent faces into a neighbor tree. More...
 
int p8est_connectivity_edge_neighbor_edge_corner (int ec, int o)
 Transform an edge corner across one of the adjacent edges into a neighbor tree. More...
 
int p8est_connectivity_edge_neighbor_corner (int c, int e, int ne, int o)
 Transform a corner across one of the adjacent edges into a neighbor tree. More...
 
p8est_connectivity_tp8est_connectivity_new (p4est_topidx_t num_vertices, p4est_topidx_t num_trees, p4est_topidx_t num_edges, p4est_topidx_t num_ett, p4est_topidx_t num_corners, p4est_topidx_t num_ctt)
 Allocate a connectivity structure. More...
 
p8est_connectivity_tp8est_connectivity_new_copy (p4est_topidx_t num_vertices, p4est_topidx_t num_trees, p4est_topidx_t num_edges, p4est_topidx_t num_corners, const double *vertices, const p4est_topidx_t *ttv, const p4est_topidx_t *ttt, const int8_t *ttf, const p4est_topidx_t *tte, const p4est_topidx_t *eoff, const p4est_topidx_t *ett, const int8_t *ete, const p4est_topidx_t *ttc, const p4est_topidx_t *coff, const p4est_topidx_t *ctt, const int8_t *ctc)
 Allocate a connectivity structure and populate from constants. More...
 
p8est_connectivity_tp8est_connectivity_bcast (p8est_connectivity_t *conn_in, int root, sc_MPI_Comm comm)
 Broadcast a connectivity structure that exists only on one process to all. More...
 
void p8est_connectivity_destroy (p8est_connectivity_t *connectivity)
 Destroy a connectivity structure. More...
 
void p8est_connectivity_set_attr (p8est_connectivity_t *conn, size_t bytes_per_tree)
 Allocate or free the attribute fields in a connectivity. More...
 
int p8est_connectivity_is_valid (p8est_connectivity_t *connectivity)
 Examine a connectivity structure. More...
 
int p8est_connectivity_is_equal (p8est_connectivity_t *conn1, p8est_connectivity_t *conn2)
 Check two connectivity structures for equality. More...
 
int p8est_connectivity_sink (p8est_connectivity_t *conn, sc_io_sink_t *sink)
 Write connectivity to a sink object. More...
 
sc_array_t * p8est_connectivity_deflate (p8est_connectivity_t *conn, p8est_connectivity_encode_t code)
 Allocate memory and store the connectivity information there. More...
 
int p8est_connectivity_save (const char *filename, p8est_connectivity_t *connectivity)
 Save a connectivity structure to disk. More...
 
p8est_connectivity_tp8est_connectivity_source (sc_io_source_t *source)
 Read connectivity from a source object. More...
 
p8est_connectivity_tp8est_connectivity_inflate (sc_array_t *buffer)
 Create new connectivity from a memory buffer. More...
 
p8est_connectivity_tp8est_connectivity_load (const char *filename, size_t *bytes)
 Load a connectivity structure from disk. More...
 
p8est_connectivity_tp8est_connectivity_new_unitcube (void)
 Create a connectivity structure for the unit cube.
 
p8est_connectivity_tp8est_connectivity_new_periodic (void)
 Create a connectivity structure for an all-periodic unit cube.
 
p8est_connectivity_tp8est_connectivity_new_rotwrap (void)
 Create a connectivity structure for a mostly periodic unit cube. More...
 
p8est_connectivity_tp8est_connectivity_new_drop (void)
 Create a connectivity structure for a five-trees geometry with a hole. More...
 
p8est_connectivity_tp8est_connectivity_new_twocubes (void)
 Create a connectivity structure that contains two cubes.
 
p8est_connectivity_tp8est_connectivity_new_twotrees (int l_face, int r_face, int orientation)
 Create a connectivity structure for two trees being rotated w.r.t. More...
 
p8est_connectivity_tp8est_connectivity_new_twowrap (void)
 Create a connectivity structure that contains two cubes where the two far ends are identified periodically.
 
p8est_connectivity_tp8est_connectivity_new_rotcubes (void)
 Create a connectivity structure that contains a few cubes. More...
 
p8est_connectivity_tp8est_connectivity_new_brick (int m, int n, int p, int periodic_a, int periodic_b, int periodic_c)
 An m by n by p array with periodicity in x, y, and z if periodic_a, periodic_b, and periodic_c are true, respectively.
 
p8est_connectivity_tp8est_connectivity_new_shell (void)
 Create a connectivity structure that builds a spherical shell. More...
 
p8est_connectivity_tp8est_connectivity_new_sphere (void)
 Create a connectivity structure that builds a solid sphere. More...
 
p8est_connectivity_tp8est_connectivity_new_torus (int nSegments)
 Create a connectivity structure that builds a revolution torus. More...
 
p8est_connectivity_tp8est_connectivity_new_byname (const char *name)
 Create connectivity structure from predefined catalogue. More...
 
p8est_connectivity_tp8est_connectivity_refine (p8est_connectivity_t *conn, int num_per_dim)
 Uniformly refine a connectivity. More...
 
void p8est_expand_face_transform (int iface, int nface, int ftransform[])
 Fill an array with the axis combination of a face neighbor transform. More...
 
p4est_topidx_t p8est_find_face_transform (p8est_connectivity_t *connectivity, p4est_topidx_t itree, int iface, int ftransform[])
 Fill an array with the axis combination of a face neighbor transform. More...
 
void p8est_find_edge_transform (p8est_connectivity_t *connectivity, p4est_topidx_t itree, int iedge, p8est_edge_info_t *ei)
 Fills an array with information about edge neighbors. More...
 
void p8est_find_corner_transform (p8est_connectivity_t *connectivity, p4est_topidx_t itree, int icorner, p8est_corner_info_t *ci)
 Fills an array with information about corner neighbors. More...
 
void p8est_connectivity_complete (p8est_connectivity_t *conn)
 Internally connect a connectivity based on tree_to_vertex information. More...
 
void p8est_connectivity_reduce (p8est_connectivity_t *conn)
 Removes corner and edge information of a connectivity such that enough information is left to run p8est_connectivity_complete successfully. More...
 
void p8est_connectivity_permute (p8est_connectivity_t *conn, sc_array_t *perm, int is_current_to_new)
 p8est_connectivity_permute Given a permutation perm of the trees in a connectivity conn, permute the trees of conn in place and update conn to match. More...
 
void p8est_connectivity_reorder (sc_MPI_Comm comm, int k, p8est_connectivity_t *conn, p8est_connect_type_t ctype)
 Reorder a connectivity using METIS. More...
 
sc_array_t * p8est_connectivity_reorder_newid (sc_MPI_Comm comm, int k, p8est_connectivity_t *conn, p8est_connect_type_t ctype, sc_array_t *newid)
 Reorder a connectivity using METIS. More...
 
void p8est_connectivity_join_faces (p8est_connectivity_t *conn, p4est_topidx_t tree_left, p4est_topidx_t tree_right, int face_left, int face_right, int orientation)
 p8est_connectivity_join_faces This function takes an existing valid connectivity conn and modifies it by joining two tree faces that are currently boundary faces. More...
 
int p8est_connectivity_is_equivalent (p8est_connectivity_t *conn1, p8est_connectivity_t *conn2)
 p8est_connectivity_is_equivalent This function compares two connectivities for equivalence: it returns true if they are the same connectivity, or if they have the same topology. More...
 
int p8est_connectivity_read_inp_stream (FILE *stream, p4est_topidx_t *num_vertices, p4est_topidx_t *num_trees, double *vertices, p4est_topidx_t *tree_to_vertex)
 Read an ABAQUS input file from a file stream. More...
 
p8est_connectivity_tp8est_connectivity_read_inp (const char *filename)
 Create a p4est connectivity from an ABAQUS input file. More...
 

Variables

const int p8est_face_corners [6][4]
 Store the corner numbers 0..7 for each tree face.
 
const int p8est_face_edges [6][4]
 Store the edge numbers 0..12 for each tree face.
 
const int p8est_face_dual [6]
 Store the face numbers in the face neighbor's system.
 
const int p8est_face_permutations [8][4]
 Store only the 8 out of 24 possible permutations that occur.
 
const int p8est_face_permutation_sets [3][4]
 Store the 3 occurring sets of 4 permutations per face.
 
const int p8est_face_permutation_refs [6][6]
 For each face combination store the permutation set. More...
 
const int p8est_face_edge_permutations [8][4]
 Store only the 8 out of 24 possible permutations that occur.
 
const int p8est_face_edge_permutation_sets [3][4]
 Store the 3 occurring sets of 4 permutations per face.
 
const int p8est_edge_faces [12][2]
 Store the face numbers 0..5 adjacent to each tree edge.
 
const int p8est_edge_corners [12][2]
 Store the corner numbers 0..8 for each tree edge.
 
const int p8est_edge_edge_corners [12][8]
 Store the edge corner numbers 0..1 for the corners touching a tree edge or -1 if combination is invalid.
 
const int p8est_edge_face_corners [12][6][2]
 Store the face corner numbers 0..3 for the faces touching a tree edge. More...
 
const int p8est_edge_face_edges [12][6]
 Store the face edge numbers 0..3 for the faces touching a tree edge. More...
 
const int p8est_corner_faces [8][3]
 Store the face numbers 0..5 for each tree corner.
 
const int p8est_corner_edges [8][3]
 Store the edge numbers 0..11 for each tree corner.
 
const int p8est_corner_face_corners [8][6]
 Store the face corner numbers for the faces touching a tree corner. More...
 
const int p8est_corner_edge_corners [8][12]
 Store the edge corner numbers for the edges touching a tree corner. More...
 
const int p8est_child_edge_faces [8][12]
 Store the faces for each child and edge, can be -1.
 
const int p8est_child_corner_faces [8][8]
 Store the faces for each child and corner, can be -1.
 
const int p8est_child_corner_edges [8][8]
 Store the edges for each child and corner, can be -1.
 

Detailed Description

The connectivity defines the coarse topology of the forest.

A 3D forest consists of one or more octrees, each of which a logical cube. Each tree has a local coordinate system, which defines the origin and the direction of its x-, y-, and z-axes as well as the numbering of its faces, edges, and corners. Each tree may connect to any other tree (including itself) across any of its faces, edges and/or corners, where the neighbor may be arbitrarily rotated and/or flipped. The p8est_connectivity data structure stores these connections.

We impose the following requirement for consistency of p8est_balance :

Note
If a connectivity implies natural connections between trees that are edge neighbors without being face neighbors, these edges shall be encoded explicitly in the connectivity. If a connectivity implies natural connections between trees that are corner neighbors without being edge or face neighbors, these corners shall be encoded explicitly in the connectivity. Please see the documentation of p8est_connectivity_t for the exact encoding convention.

We provide various predefined connectivities by dedicated constructors, such as

Macro Definition Documentation

◆ P8EST_CHILDREN

#define P8EST_CHILDREN   8

The number of children of an octant.

also the nmber of corners

◆ P8EST_FACES

#define P8EST_FACES   (2 * P8EST_DIM)

The number of faces of an octant.

Note
for uniform naming reasons, an octant is represented by the datatype p8est_quadrant_t

◆ P8EST_ONDISK_FORMAT

#define P8EST_ONDISK_FORMAT   0x3000009

Increase this number whenever the on-disk format for p8est_connectivity, p8est, or any other 3D data structure changes.

The format for reading and writing must be the same.

Typedef Documentation

◆ p8est_connectivity_t

This structure holds the 3D inter-tree connectivity information.

Identification of arbitrary faces, edges and corners is possible.

The arrays tree_to_* are stored in z ordering. For corners the order wrt. zyx is 000 001 010 011 100 101 110 111. For faces the order is -x +x -y +y -z +z. They are allocated [0][0]..[0][N-1]..[num_trees-1][0]..[num_trees-1][N-1]. where N is 6 for tree and face, 8 for corner, 12 for edge. If a face is on the physical boundary it must connect to itself.

The values for tree_to_face are in 0..23 where ttf % 6 gives the face number and ttf / 6 the face orientation code. The orientation is determined as follows. Let my_face and other_face be the two face numbers of the connecting trees in 0..5. Then the first face corner of the lower of my_face and other_face connects to a face corner numbered 0..3 in the higher of my_face and other_face. The face orientation is defined as this number. If my_face == other_face, treating either of both faces as the lower one leads to the same result.

It is valid to specify num_vertices as 0. In this case vertices and tree_to_vertex are set to NULL. Otherwise the vertex coordinates are stored in the array vertices as [0][0]..[0][2]..[num_vertices-1][0]..[num_vertices-1][2]. Vertex coordinates are optional and not used for inferring topology.

The edges are stored when they connect trees that are not already face neighbors at that specific edge. In this case tree_to_edge indexes into ett_offset. Otherwise the tree_to_edge entry must be -1 and this edge is ignored. If num_edges == 0, tree_to_edge and edge_to_* arrays are set to NULL.

The arrays edge_to_* store a variable number of entries per edge. For edge e these are at position [ett_offset[e]]..[ett_offset[e+1]-1]. Their number for edge e is ett_offset[e+1] - ett_offset[e]. The entries encode all trees adjacent to edge e. The size of the edge_to_* arrays is num_ett = ett_offset[num_edges]. The edge_to_edge array holds values in 0..23, where the lower 12 indicate one edge orientation and the higher 12 the opposite edge orientation.

The corners are stored when they connect trees that are not already edge or face neighbors at that specific corner. In this case tree_to_corner indexes into ctt_offset. Otherwise the tree_to_corner entry must be -1 and this corner is ignored. If num_corners == 0, tree_to_corner and corner_to_* arrays are set to NULL.

The arrays corner_to_* store a variable number of entries per corner. For corner c these are at position [ctt_offset[c]]..[ctt_offset[c+1]-1]. Their number for corner c is ctt_offset[c+1] - ctt_offset[c]. The entries encode all trees adjacent to corner c. The size of the corner_to_* arrays is num_ctt = ctt_offset[num_corners].

The *_to_attr arrays may have arbitrary contents defined by the user.

Note
If a connectivity implies natural connections between trees that are edge neighbors without being face neighbors, these edges shall be encoded explicitly in the connectivity. If a connectivity implies natural connections between trees that are corner neighbors without being edge or face neighbors, these corners shall be encoded explicitly in the connectivity.

Enumeration Type Documentation

◆ p8est_connect_type_t

Characterize a type of adjacency.

Several functions involve relationships between neighboring trees and/or quadrants, and their behavior depends on how one defines adjacency: 1) entities are adjacent if they share a face, or 2) entities are adjacent if they share a face or corner, or 3) entities are adjacent if they share a face, corner or edge. p8est_connect_type_t is used to choose the desired behavior. This enum must fit into an int8_t.

Enumerator
P8EST_CONNECT_SELF 

No balance whatsoever.

P8EST_CONNECT_FACE 

Balance across faces only.

P8EST_CONNECT_EDGE 

Balance across faces and edges.

P8EST_CONNECT_ALMOST 

= CORNER - 1.

P8EST_CONNECT_CORNER 

Balance faces, edges, corners.

P8EST_CONNECT_FULL 

= CORNER.

◆ p8est_connectivity_encode_t

Typedef for serialization method.

Enumerator
P8EST_CONN_ENCODE_LAST 

Invalid entry to close the list.

Function Documentation

◆ p8est_connect_type_int()

int p8est_connect_type_int ( p8est_connect_type_t  btype)

Convert the p8est_connect_type_t into a number.

Parameters
[in]btypeThe balance type to convert.
Returns
Returns 1, 2 or 3.

◆ p8est_connect_type_string()

const char * p8est_connect_type_string ( p8est_connect_type_t  btype)

Convert the p8est_connect_type_t into a const string.

Parameters
[in]btypeThe balance type to convert.
Returns
Returns a pointer to a constant string.

◆ p8est_connectivity_bcast()

p8est_connectivity_t * p8est_connectivity_bcast ( p8est_connectivity_t conn_in,
int  root,
sc_MPI_Comm  comm 
)

Broadcast a connectivity structure that exists only on one process to all.

On the other processors, it will be allocated using p8est_connectivity_new.

Parameters
[in]conn_inFor the root process the connectivity to be broadcast, for the other processes it must be NULL.
[in]rootThe rank of the process that provides the connectivity.
[in]commThe MPI communicator.
Returns
For the root process this is a pointer to conn_in. Else, a pointer to a newly allocated connectivity structure with the same values as conn_in on the root process.

◆ p8est_connectivity_complete()

void p8est_connectivity_complete ( p8est_connectivity_t conn)

Internally connect a connectivity based on tree_to_vertex information.

Periodicity that is not inherent in the list of vertices will be lost.

Parameters
[in,out]connThe connectivity needs to have proper vertices and tree_to_vertex fields. The tree_to_tree and tree_to_face fields must be allocated and satisfy p8est_connectivity_is_valid (conn) but will be overwritten. The edge and corner fields will be freed and allocated anew.

◆ p8est_connectivity_deflate()

sc_array_t * p8est_connectivity_deflate ( p8est_connectivity_t conn,
p8est_connectivity_encode_t  code 
)

Allocate memory and store the connectivity information there.

Parameters
[in]connThe connectivity structure to be exported to memory.
[in]codeEncoding and compression method for serialization.
Returns
Newly created array that contains the information.

◆ p8est_connectivity_destroy()

void p8est_connectivity_destroy ( p8est_connectivity_t connectivity)

Destroy a connectivity structure.

Also destroy all attributes.

Examples
simple/simple3.c.

◆ p8est_connectivity_edge_neighbor_corner()

int p8est_connectivity_edge_neighbor_corner ( int  c,
int  e,
int  ne,
int  o 
)

Transform a corner across one of the adjacent edges into a neighbor tree.

This version expects the neighbor edge and orientation separately.

Parameters
[in]cA corner number in 0..7.
[in]eAn edge 0..11 that touches the corner c.
[in]neA neighbor edge that is on the other side of e.
[in]oThe orientation between tree boundary edges e and ne.
Returns
Corner number seen from the neighbor.

◆ p8est_connectivity_edge_neighbor_edge_corner()

int p8est_connectivity_edge_neighbor_edge_corner ( int  ec,
int  o 
)

Transform an edge corner across one of the adjacent edges into a neighbor tree.

Parameters
[in]ecAn edge corner number in 0..1.
[in]oThe orientation of a tree boundary edge connection.
Returns
The edge corner number seen from the other tree.

◆ p8est_connectivity_face_neighbor_corner()

int p8est_connectivity_face_neighbor_corner ( int  c,
int  f,
int  nf,
int  o 
)

Transform a corner across one of the adjacent faces into a neighbor tree.

This version expects the neighbor face and orientation separately.

Parameters
[in]cA corner number in 0..7.
[in]fA face number that touches the corner c.
[in]nfA neighbor face that is on the other side of f.
[in]oThe orientation between tree boundary faces f and nf.
Returns
The number of the corner seen from the neighbor tree.

◆ p8est_connectivity_face_neighbor_corner_set()

int p8est_connectivity_face_neighbor_corner_set ( int  c,
int  f,
int  nf,
int  set 
)

Transform a corner across one of the adjacent faces into a neighbor tree.

It expects a face permutation index that has been precomputed.

Parameters
[in]cA corner number in 0..7.
[in]fA face number that touches the corner c.
[in]nfA neighbor face that is on the other side of f.
[in]setA value from p8est_face_permutation_sets that is obtained using f, nf, and a valid orientation: ref = p8est_face_permutation_refs[f][nf]; set = p8est_face_permutation_sets[ref][orientation];
Returns
The corner number in 0..7 seen from the other face.

◆ p8est_connectivity_face_neighbor_edge()

int p8est_connectivity_face_neighbor_edge ( int  e,
int  f,
int  nf,
int  o 
)

Transform an edge across one of the adjacent faces into a neighbor tree.

This version expects the neighbor face and orientation separately.

Parameters
[in]eA edge number in 0..11.
[in]fA face 0..5 that touches the edge e.
[in]nfA neighbor face that is on the other side of f.
[in]oThe orientation between tree boundary faces f and nf.
Returns
The edge's number seen from the neighbor.

◆ p8est_connectivity_face_neighbor_face_corner()

int p8est_connectivity_face_neighbor_face_corner ( int  fc,
int  f,
int  nf,
int  o 
)

Transform a face corner across one of the adjacent faces into a neighbor tree.

This version expects the neighbor face and orientation separately.

Parameters
[in]fcA face corner number in 0..3.
[in]fA face that the face corner fc is relative to.
[in]nfA neighbor face that is on the other side of f.
[in]oThe orientation between tree boundary faces f and nf.
Returns
The face corner number relative to the neighbor's face.

◆ p8est_connectivity_face_neighbor_face_edge()

int p8est_connectivity_face_neighbor_face_edge ( int  fe,
int  f,
int  nf,
int  o 
)

Transform a face-edge across one of the adjacent faces into a neighbor tree.

This version expects the neighbor face and orientation separately.

Parameters
[in]feA face edge number in 0..3.
[in]fA face number that touches the edge e.
[in]nfA neighbor face that is on the other side of f.
[in]oThe orientation between tree boundary faces f and nf.
Returns
The face edge number seen from the neighbor tree.

◆ p8est_connectivity_get_neighbor_transforms()

void p8est_connectivity_get_neighbor_transforms ( p8est_connectivity_t conn,
p4est_topidx_t  tree_id,
p8est_connect_type_t  boundary_type,
int  boundary_index,
sc_array_t *  neighbor_transform_array 
)

Fill an array with the neighbor transforms based on a specific boundary type.

This function generalizes all other inter-tree transformation objects

Parameters
[in]connConnectivity structure.
[in]tree_idThe number of the tree.
[in]boundary_typeType of boundary connection (self, face, edge, corner).
[in]boundary_indexThe index of the boundary.
[in,out]neighbor_transform_arrayArray of the neighbor transforms.

◆ p8est_connectivity_inflate()

p8est_connectivity_t * p8est_connectivity_inflate ( sc_array_t *  buffer)

Create new connectivity from a memory buffer.

This function aborts on malloc errors.

Parameters
[in]bufferThe connectivity is created from this memory buffer.
Returns
The newly created connectivity, or NULL on format error of the buffered connectivity data.

◆ p8est_connectivity_is_equal()

int p8est_connectivity_is_equal ( p8est_connectivity_t conn1,
p8est_connectivity_t conn2 
)

Check two connectivity structures for equality.

Returns
Returns true if structures are equal, false otherwise.

◆ p8est_connectivity_is_equivalent()

int p8est_connectivity_is_equivalent ( p8est_connectivity_t conn1,
p8est_connectivity_t conn2 
)

p8est_connectivity_is_equivalent This function compares two connectivities for equivalence: it returns true if they are the same connectivity, or if they have the same topology.

The definition of topological sameness is strict: there is no attempt made to determine whether permutation and/or rotation of the trees makes the connectivities equivalent.

Parameters
[in]conn1a valid connectivity
[out]conn2a valid connectivity

◆ p8est_connectivity_is_valid()

int p8est_connectivity_is_valid ( p8est_connectivity_t connectivity)

Examine a connectivity structure.

Returns
Returns true if structure is valid, false otherwise.

◆ p8est_connectivity_join_faces()

void p8est_connectivity_join_faces ( p8est_connectivity_t conn,
p4est_topidx_t  tree_left,
p4est_topidx_t  tree_right,
int  face_left,
int  face_right,
int  orientation 
)

p8est_connectivity_join_faces This function takes an existing valid connectivity conn and modifies it by joining two tree faces that are currently boundary faces.

Parameters
[in,out]connconnectivity that will be altered.
[in]tree_lefttree that will be on the left side of the joined faces.
[in]tree_righttree that will be on the right side of the joined faces.
[in]face_leftface of tree_left that will be joined.
[in]face_rightface of tree_right that will be joined.
[in]orientationthe orientation of face_left and face_right once joined (see the description of p8est_connectivity_t to understand orientation).

◆ p8est_connectivity_load()

p8est_connectivity_t * p8est_connectivity_load ( const char *  filename,
size_t *  bytes 
)

Load a connectivity structure from disk.

Parameters
[in]filenameName of the file to read.
[out]bytesSize in bytes of connectivity on disk or NULL.
Returns
Returns valid connectivity, or NULL on file error.

◆ p8est_connectivity_memory_used()

size_t p8est_connectivity_memory_used ( p8est_connectivity_t conn)

Calculate memory usage of a connectivity structure.

Parameters
[in]connConnectivity structure.
Returns
Memory used in bytes.

◆ p8est_connectivity_new()

p8est_connectivity_t * p8est_connectivity_new ( p4est_topidx_t  num_vertices,
p4est_topidx_t  num_trees,
p4est_topidx_t  num_edges,
p4est_topidx_t  num_ett,
p4est_topidx_t  num_corners,
p4est_topidx_t  num_ctt 
)

Allocate a connectivity structure.

The attribute fields are initialized to NULL.

Parameters
[in]num_verticesNumber of total vertices (i.e. geometric points).
[in]num_treesNumber of trees in the forest.
[in]num_edgesNumber of tree-connecting edges.
[in]num_ettNumber of total trees in edge_to_tree array.
[in]num_cornersNumber of tree-connecting corners.
[in]num_cttNumber of total trees in corner_to_tree array.
Returns
A connectivity structure with allocated arrays.

◆ p8est_connectivity_new_byname()

p8est_connectivity_t * p8est_connectivity_new_byname ( const char *  name)

Create connectivity structure from predefined catalogue.

Parameters
[in]nameInvokes connectivity_new_* function. brick235 brick (2, 3, 5, 0, 0, 0) periodic periodic rotcubes rotcubes rotwrap rotwrap shell shell sphere sphere twocubes twocubes twowrap twowrap unit unitcube
Returns
An initialized connectivity if name is defined, NULL else.

◆ p8est_connectivity_new_copy()

p8est_connectivity_t * p8est_connectivity_new_copy ( p4est_topidx_t  num_vertices,
p4est_topidx_t  num_trees,
p4est_topidx_t  num_edges,
p4est_topidx_t  num_corners,
const double *  vertices,
const p4est_topidx_t ttv,
const p4est_topidx_t ttt,
const int8_t *  ttf,
const p4est_topidx_t tte,
const p4est_topidx_t eoff,
const p4est_topidx_t ett,
const int8_t *  ete,
const p4est_topidx_t ttc,
const p4est_topidx_t coff,
const p4est_topidx_t ctt,
const int8_t *  ctc 
)

Allocate a connectivity structure and populate from constants.

The attribute fields are initialized to NULL.

Parameters
[in]num_verticesNumber of total vertices (i.e. geometric points).
[in]num_treesNumber of trees in the forest.
[in]num_edgesNumber of tree-connecting edges.
[in]num_cornersNumber of tree-connecting corners.
[in]verticesCoordinates of the vertices of the trees.
[in]ttvThe tree-to-vertex array.
[in]tttThe tree-to-tree array.
[in]ttfThe tree-to-face array (int8_t).
[in]tteThe tree-to-edge array.
[in]eoffEdge-to-tree offsets (num_edges + 1 values). This must always be non-NULL; in trivial cases it is just a pointer to a p4est_topix value of 0.
[in]ettThe edge-to-tree array.
[in]eteThe edge-to-edge array.
[in]ttcThe tree-to-corner array.
[in]coffCorner-to-tree offsets (num_corners + 1 values). This must always be non-NULL; in trivial cases it is just a pointer to a p4est_topix value of 0.
[in]cttThe corner-to-tree array.
[in]ctcThe corner-to-corner array.
Returns
The connectivity is checked for validity.

◆ p8est_connectivity_new_drop()

p8est_connectivity_t * p8est_connectivity_new_drop ( void  )

Create a connectivity structure for a five-trees geometry with a hole.

The geometry is a 3D extrusion of the two drop example, and covers [0, 3]*[0, 2]*[0, 3]. The additional dimension is Y.

Examples
simple/simple3.c.

◆ p8est_connectivity_new_rotcubes()

p8est_connectivity_t * p8est_connectivity_new_rotcubes ( void  )

Create a connectivity structure that contains a few cubes.

These are rotated against each other to stress the topology routines.

Examples
points/generate_points2.c, and simple/simple3.c.

◆ p8est_connectivity_new_rotwrap()

p8est_connectivity_t * p8est_connectivity_new_rotwrap ( void  )

Create a connectivity structure for a mostly periodic unit cube.

The left and right faces are identified, and bottom and top rotated. Front and back are not identified.

Examples
points/generate_points2.c, and simple/simple3.c.

◆ p8est_connectivity_new_shell()

p8est_connectivity_t * p8est_connectivity_new_shell ( void  )

Create a connectivity structure that builds a spherical shell.

It is made up of six connected parts [-1,1]x[-1,1]x[1,2]. This connectivity reuses vertices and relies on a geometry transformation. It is thus not suitable for p8est_connectivity_complete.

Examples
simple/simple3.c.

◆ p8est_connectivity_new_sphere()

p8est_connectivity_t * p8est_connectivity_new_sphere ( void  )

Create a connectivity structure that builds a solid sphere.

It is made up of two layers and a cube in the center. This connectivity reuses vertices and relies on a geometry transformation. It is thus not suitable for p8est_connectivity_complete.

Examples
simple/simple3.c.

◆ p8est_connectivity_new_torus()

p8est_connectivity_t * p8est_connectivity_new_torus ( int  nSegments)

Create a connectivity structure that builds a revolution torus.

This connectivity reuses vertices and relies on a geometry transformation. It is thus not suitable for p8est_connectivity_complete.

This connectivity reuses ideas from disk2d connectivity. More precisely the torus is divided into segments around the revolution axis, each segments is made of 5 trees (à la disk2d). The total number of trees if 5 times the number of segments.

This connectivity is meant to be used with p8est_geometry_new_torus

Parameters
[in]nSegmentsnumber of trees along the great circle
Examples
simple/simple3.c.

◆ p8est_connectivity_new_twotrees()

p8est_connectivity_t * p8est_connectivity_new_twotrees ( int  l_face,
int  r_face,
int  orientation 
)

Create a connectivity structure for two trees being rotated w.r.t.

each other in a user-defined way.

Parameters
[in]l_faceindex of left face
[in]r_faceindex of right face
[in]orientationorientation of trees w.r.t. each other

◆ p8est_connectivity_permute()

void p8est_connectivity_permute ( p8est_connectivity_t conn,
sc_array_t *  perm,
int  is_current_to_new 
)

p8est_connectivity_permute Given a permutation perm of the trees in a connectivity conn, permute the trees of conn in place and update conn to match.

Parameters
[in,out]connThe connectivity whose trees are permuted.
[in]permA permutation array, whose elements are size_t's.
[in]is_current_to_newif true, the jth entry of perm is the new index for the entry whose current index is j, otherwise the jth entry of perm is the current index of the tree whose index will be j after the permutation.

◆ p8est_connectivity_read_inp()

p8est_connectivity_t * p8est_connectivity_read_inp ( const char *  filename)

Create a p4est connectivity from an ABAQUS input file.

This utility function reads a basic ABAQUS file supporting element type with the prefix C2D4, CPS4, and S4 in 2D and of type C3D8 reading them as bilinear quadrilateral and trilinear hexahedral trees respectively.

A basic 2D mesh is given below. The *Node section gives the vertex number and x, y, and z components for each vertex. The *Element section gives the 4 vertices in 2D (8 vertices in 3D) of each element in counter clockwise order. So in 2D the nodes are given as:

4 3 +----------------—+

+----------------—+ 1 2

and in 3D they are given as:

8 7 +------------------—+ |\ |\ | \ | \ | \ | \ | \ | \ | 5+------------------—+6 | | | | +-—|-------------—+ | 4\ | 3 \ | \ | \ | \ | \ | | | +------------------—+ 1 2

*Heading
box.inp
*Node
1, 5, -5, 5
2, 5, 5, 5
3, 5, 0, 5
4, -5, 5, 5
5, 0, 5, 5
6, -5, -5, 5
7, -5, 0, 5
8, 0, -5, 5
9, 0, 0, 5
10, 5, 5, -5
11, 5, -5, -5
12, 5, 0, -5
13, -5, -5, -5
14, 0, -5, -5
15, -5, 5, -5
16, -5, 0, -5
17, 0, 5, -5
18, 0, 0, -5
19, -5, -5, 0
20, 5, -5, 0
21, 0, -5, 0
22, -5, 5, 0
23, -5, 0, 0
24, 5, 5, 0
25, 0, 5, 0
26, 5, 0, 0
27, 0, 0, 0
*Element, type=C3D8, ELSET=EB1
1, 6, 19, 23, 7, 8, 21, 27, 9
2, 19, 13, 16, 23, 21, 14, 18, 27
3, 7, 23, 22, 4, 9, 27, 25, 5
4, 23, 16, 15, 22, 27, 18, 17, 25
5, 8, 21, 27, 9, 1, 20, 26, 3
6, 21, 14, 18, 27, 20, 11, 12, 26
7, 9, 27, 25, 5, 3, 26, 24, 2
8, 27, 18, 17, 25, 26, 12, 10, 24

This function reads a mesh from filename and returns an associated p4est connectivity.

Parameters
[in]filenamefile to read the connectivity from
Returns
an allocated connectivity associated with the mesh in filename

◆ p8est_connectivity_read_inp_stream()

int p8est_connectivity_read_inp_stream ( FILE *  stream,
p4est_topidx_t num_vertices,
p4est_topidx_t num_trees,
double *  vertices,
p4est_topidx_t tree_to_vertex 
)

Read an ABAQUS input file from a file stream.

This utility function reads a basic ABAQUS file supporting element type with the prefix C2D4, CPS4, and S4 in 2D and of type C3D8 reading them as bilinear quadrilateral and trilinear hexahedral trees respectively.

A basic 2D mesh is given below. The *Node section gives the vertex number and x, y, and z components for each vertex. The *Element section gives the 4 vertices in 2D (8 vertices in 3D) of each element in counter clockwise order. So in 2D the nodes are given as:

4 3 +----------------—+

+----------------—+ 1 2

and in 3D they are given as:

8 7 +------------------—+ |\ |\ | \ | \ | \ | \ | \ | \ | 5+------------------—+6 | | | | +-—|-------------—+ | 4\ | 3 \ | \ | \ | \ | \ | | | +------------------—+ 1 2

*Heading
box.inp
*Node
1, 5, -5, 5
2, 5, 5, 5
3, 5, 0, 5
4, -5, 5, 5
5, 0, 5, 5
6, -5, -5, 5
7, -5, 0, 5
8, 0, -5, 5
9, 0, 0, 5
10, 5, 5, -5
11, 5, -5, -5
12, 5, 0, -5
13, -5, -5, -5
14, 0, -5, -5
15, -5, 5, -5
16, -5, 0, -5
17, 0, 5, -5
18, 0, 0, -5
19, -5, -5, 0
20, 5, -5, 0
21, 0, -5, 0
22, -5, 5, 0
23, -5, 0, 0
24, 5, 5, 0
25, 0, 5, 0
26, 5, 0, 0
27, 0, 0, 0
*Element, type=C3D8, ELSET=EB1
1, 6, 19, 23, 7, 8, 21, 27, 9
2, 19, 13, 16, 23, 21, 14, 18, 27
3, 7, 23, 22, 4, 9, 27, 25, 5
4, 23, 16, 15, 22, 27, 18, 17, 25
5, 8, 21, 27, 9, 1, 20, 26, 3
6, 21, 14, 18, 27, 20, 11, 12, 26
7, 9, 27, 25, 5, 3, 26, 24, 2
8, 27, 18, 17, 25, 26, 12, 10, 24

This code can be called two ways. The first, when vertex==NULL and tree_to_vertex==NULL, is used to count the number of trees and vertices in the connectivity to be generated by the .inp mesh in the stream. The second, when vertices!=NULL and tree_to_vertex!=NULL, fill vertices and tree_to_vertex. In this case num_vertices and num_trees need to be set to the maximum number of entries allocated in vertices and tree_to_vertex.

Parameters
[in,out]streamfile stream to read the connectivity from
[in,out]num_verticesthe number of vertices in the connectivity
[in,out]num_treesthe number of trees in the connectivity
[out]verticesthe list of vertices of the connectivity
[out]tree_to_vertexthe tree_to_vertex map of the connectivity
Returns
0 if successful and nonzero if not

◆ p8est_connectivity_reduce()

void p8est_connectivity_reduce ( p8est_connectivity_t conn)

Removes corner and edge information of a connectivity such that enough information is left to run p8est_connectivity_complete successfully.

The reduced connectivity still passes p8est_connectivity_is_valid.

Parameters
[in,out]connThe connectivity to be reduced.

◆ p8est_connectivity_refine()

p8est_connectivity_t * p8est_connectivity_refine ( p8est_connectivity_t conn,
int  num_per_dim 
)

Uniformly refine a connectivity.

This is useful if you would like to uniformly refine by something other than a power of 2.

Parameters
[in]connA valid connectivity
[in]num_per_dimThe number of new trees in each direction. Must use no more than P8EST_OLD_QMAXLEVEL bits.
Returns
a refined connectivity.

◆ p8est_connectivity_reorder()

void p8est_connectivity_reorder ( sc_MPI_Comm  comm,
int  k,
p8est_connectivity_t conn,
p8est_connect_type_t  ctype 
)

Reorder a connectivity using METIS.

This function takes a connectivity conn and a parameter k, which will typically be the number of processes, and reorders the trees such that if every processes is assigned (num_trees / k) trees, the communication volume will be minimized. This is intended for use with connectivities that contain a large number of trees. This should be done BEFORE a p8est is created using the connectivity. This is done in place: any data structures that use indices to refer to trees before this procedure will be invalid. Note that this routine calls metis and not parmetis because the connectivity is copied on every process. A communicator is required because I'm not positive that metis is deterministic. ctype determines when an edge exist between two trees in the dual graph used by metis in the reordering.

Parameters
[in]commMPI communicator.
[in]kif k > 0, the number of pieces metis will use to guide the reordering; if k = 0, the number of pieces will be determined from the MPI communicator.
[in,out]connconnectivity that will be reordered.
[in]ctypedetermines when an edge exists in the dual graph of the connectivity structure.

◆ p8est_connectivity_reorder_newid()

sc_array_t * p8est_connectivity_reorder_newid ( sc_MPI_Comm  comm,
int  k,
p8est_connectivity_t conn,
p8est_connect_type_t  ctype,
sc_array_t *  newid 
)

Reorder a connectivity using METIS.

This is the same form of p8est_connectivity_reorder but it takes an initialized sc array newid as extra argument. In this way, the users can map old indices to new indices in the case it is necessary (for instance to retrieve high-order nodes previously stored in an array with old indices).

Parameters
[in]commMPI communicator.
[in]kif k > 0, the number of pieces metis will use to guide the reordering; if k = 0, the number of pieces will be determined from the MPI communicator.
[in,out]connconnectivity that will be reordered.
[in]ctypedetermines when an edge exists in the dual graph of the connectivity structure.
[in,out]newidarray that maps old tree indices to new ones. newid has to be an sc_array and it has to be initialized (non-NULL) with element size of size_t (using sc_array_new (sizeof (size_t))). Input length arbitrary, output length modified.

◆ p8est_connectivity_save()

int p8est_connectivity_save ( const char *  filename,
p8est_connectivity_t connectivity 
)

Save a connectivity structure to disk.

Parameters
[in]filenameName of the file to write.
[in]connectivityValid connectivity structure.
Returns
Returns 0 on success, nonzero on file error.

◆ p8est_connectivity_set_attr()

void p8est_connectivity_set_attr ( p8est_connectivity_t conn,
size_t  bytes_per_tree 
)

Allocate or free the attribute fields in a connectivity.

Parameters
[in,out]connThe conn->*_to_attr fields must either be NULL or previously be allocated by this function.
[in]bytes_per_treeIf 0, tree_to_attr is freed (being NULL is ok). If positive, requested space is allocated.

◆ p8est_connectivity_sink()

int p8est_connectivity_sink ( p8est_connectivity_t conn,
sc_io_sink_t *  sink 
)

Write connectivity to a sink object.

Parameters
[in]connThe connectivity to be written.
[in,out]sinkThe connectivity is written into this sink.
Returns
0 on success, nonzero on error.

◆ p8est_connectivity_source()

p8est_connectivity_t * p8est_connectivity_source ( sc_io_source_t *  source)

Read connectivity from a source object.

Parameters
[in,out]sourceThe connectivity is read from this source.
Returns
The newly created connectivity, or NULL on error.

◆ p8est_expand_face_transform()

void p8est_expand_face_transform ( int  iface,
int  nface,
int  ftransform[] 
)

Fill an array with the axis combination of a face neighbor transform.

Parameters
[in]ifaceThe number of the originating face.
[in]nfaceEncoded as nface = r * 6 + nf, where nf = 0..5 is the neigbbor's connecting face number and r = 0..3 is the relative orientation to the neighbor's face. This encoding matches p8est_connectivity_t.
[out]ftransformThis array holds 9 integers. [0]..[2] The coordinate axis sequence of the origin face, the first two referring to the tangentials and the third to the normal. A permutation of (0, 1, 2). [3]..[5] The coordinate axis sequence of the target face. [6]..[8] Edge reversal flags for tangential axes (boolean); face code in [0, 3] for the normal coordinate q: 0: q' = -q 1: q' = q + 1 2: q' = q - 1 3: q' = 2 - q

◆ p8est_find_corner_transform()

void p8est_find_corner_transform ( p8est_connectivity_t connectivity,
p4est_topidx_t  itree,
int  icorner,
p8est_corner_info_t ci 
)

Fills an array with information about corner neighbors.

Parameters
[in]connectivityConnectivity structure.
[in]itreeThe number of the originating tree.
[in]icornerThe number of the originating corner.
[in,out]ciA p8est_corner_info_t structure with initialized array.

◆ p8est_find_edge_transform()

void p8est_find_edge_transform ( p8est_connectivity_t connectivity,
p4est_topidx_t  itree,
int  iedge,
p8est_edge_info_t ei 
)

Fills an array with information about edge neighbors.

Parameters
[in]connectivityConnectivity structure.
[in]itreeThe number of the originating tree.
[in]iedgeThe number of the originating edge.
[in,out]eiA p8est_edge_info_t structure with initialized array.

◆ p8est_find_face_transform()

p4est_topidx_t p8est_find_face_transform ( p8est_connectivity_t connectivity,
p4est_topidx_t  itree,
int  iface,
int  ftransform[] 
)

Fill an array with the axis combination of a face neighbor transform.

Parameters
[in]connectivityConnectivity structure.
[in]itreeThe number of the originating tree.
[in]ifaceThe number of the originating tree's face.
[out]ftransformThis array holds 9 integers. [0]..[2] The coordinate axis sequence of the origin face. [3]..[5] The coordinate axis sequence of the target face. [6]..[8] Edge reversal flag for axes t1, t2; face code for n;
See also
p8est_expand_face_transform.
Returns
The face neighbor tree if it exists, -1 otherwise.

◆ p8est_neighbor_transform_coordinates()

void p8est_neighbor_transform_coordinates ( const p8est_neighbor_transform_t nt,
const p4est_qcoord_t  self_coords[P8EST_DIM],
p4est_qcoord_t  neigh_coords[P8EST_DIM] 
)

Transform from self's coordinate system to neighbor's coordinate system.

Parameters
[in]ntA neighbor transform.
[in]self_coordsInput quadrant coordinates in self coordinates.
[out]neigh_coordsCoordinates transformed into neighbor coordinates.

◆ p8est_neighbor_transform_coordinates_reverse()

void p8est_neighbor_transform_coordinates_reverse ( const p8est_neighbor_transform_t nt,
const p4est_qcoord_t  neigh_coords[P8EST_DIM],
p4est_qcoord_t  self_coords[P8EST_DIM] 
)

Transform from neighbor's coordinate system to self's coordinate system.

Parameters
[in]ntA neighbor transform.
[in]neigh_coordsInput quadrant coordinates in self coordinates.
[out]self_coordsCoordinates transformed into neighbor coordinates.

Variable Documentation

◆ p8est_corner_edge_corners

const int p8est_corner_edge_corners[8][12]
extern

Store the edge corner numbers for the edges touching a tree corner.

Is -1 for invalid combinations.

◆ p8est_corner_face_corners

const int p8est_corner_face_corners[8][6]
extern

Store the face corner numbers for the faces touching a tree corner.

Is -1 for invalid combinations.

◆ p8est_edge_face_corners

const int p8est_edge_face_corners[12][6][2]
extern

Store the face corner numbers 0..3 for the faces touching a tree edge.

Is -1 for invalid combinations of indices

◆ p8est_edge_face_edges

const int p8est_edge_face_edges[12][6]
extern

Store the face edge numbers 0..3 for the faces touching a tree edge.

Is -1 for invalid combinations of indices

◆ p8est_face_permutation_refs

const int p8est_face_permutation_refs[6][6]
extern

For each face combination store the permutation set.

The order is [my_face][neighbor_face]