API reference
P4estTypes.Connectivity — TypeConnectivity{X,P}Connectivity for a Pxest which holds the mesh information for the roots of Pxest quadtrees or octrees. The parameter X is 4 if the roots are quads (2D aka p4est) and 8 if they are hexes (3D aka p8est).
Fields
pointer: The pointer (of typeP) can be a pointer to either aP4estTypes.P4est.p4est_connectivityor aP4estTypes.P4est.p8est_connectivity. See the help documentation for these types for more information about the underlying p4est structures.
Usage
Connectivity{4}(name::Symbol)Construct a connectivity mesh for the roots of a forest-of-quadtrees using p4est's built-in mesh connectivities. Valid values for name are
:unitsquare: the unit square.:periodic: all-periodic unit square.:rotwrap: periodic unit square (the left and right faces are identified, and bottom and top opposite).:corner: three-tree mesh around a corner.:pillow: two trees on top of each other.:moebius: a five-tree moebius band.:star: six-tree star.:cubed: six sides of a unit cube.:disk_nonperiodic: five-tree flat spherical disk.:icosahedron: for mapping the sphere using an icosahedron (see@doc P4estTypes.P4est.p4est_connectivity_new_icosahedronfor more info).:shell2d: 2D spherical shell.:disk2d: maps a 2D disk.
Connectivity{4}(:disk, periodic_x::Bool, periodic_y::Bool)Create a connectivity structure for a five-tree flat spherical disk. The arguments periodic_x and periodic_y determine if the disk is periodic in the x and y directions, respectively.
See @doc P4estTypes.P4est.p4est_connectivity_new_disk for detailed information.
Connectivity{8}(name::Symbol)Construct a connectivity mesh for the roots of a forest-of-octrees using p8est's built-in mesh connectivities. Valid values for name are
:unitcube: the unit cube.:periodic: an all-periodic unit cube.:rotcubes: contains a few cubes (these are rotated against each other to stress the topology routines).:rotwrap: a mostly periodic unit cube (see@doc P4estTypes.P4est.p8est_connectivity_new_rotwrap).:shell: a spherical shell (see@doc P4estTypes.P4est.p8est_connectivity_new_shell).:sphere: a solid sphere (see@doc P4estTypes.P4est.p8est_connectivity_new_sphere).:twocubes: two cubes.:twowrap: two cubes where the two far ends are identified periodically.
Connectivity{8}(:torus, nsegments)Create a connectivity structure that builds a revolution torus. Here nsegments are the number of trees along the great circle.
See @doc P4estTypes.P4est.p8est_connectivity_new_torus for detailed information.
Connectivity{8}(:torus, nsegments)Create a connectivity structure that builds a revolution torus. Here nsegments are the number of trees along the great circle.
See @doc P4estTypes.P4est.p8est_connectivity_new_torus for detailed information.
Connectivity{X}(:twotrees, l_face, r_face, orientation) where {X}Create a connectivity structure (X=4 for quadtrees and X=8 for octrees) for two trees being rotated with respect to each other in a user-defined way. Here l_face and r_face are the 0-based indices of left and right faces, respectively. The argument orientation gives the orientation code of the trees with respect to each other.
Connectivity{X}(vertices, elements) where {X}Creates a connectivity from the given list of vertices and element-to-vertex connectivity. The parameter set to X=4 is for quads and X=8 for hexes.
vertices: should be a number-of-vertices by 3 matrix where the columns correspond to x, y, and z coordinates (typically thezcoordinate will be zero for a 2D forest).elements: should be a number-of-vertices by 4 or 8 matrix where the columns vertex indices used to define each element. Note that z-ordering should be used, and it should use zero-indexing.
Connectivity{X}(filename::String) where {X}Create a connectivity from an ABAQUS input at filename. The parameter set to X=4 is for quads and X=8 for hexes.
See @doc P4estTypes.P4est.p4est_connectivity_read_inp and @doc P4estTypes.P4est.p8est_connectivity_read_inp for example ABAQUS input files.
See also
brick: a function to create a rectangularConnectivity.connectivity: a function to get the connectivity of aPxest.refine: a function to create a refinedConnectivity.
P4estTypes.Connectivity — Methodfunction Connectivity{X}(vertices::AbstractVector, cells::AbstractVector) where {X}
Example:
'''julia
julia> vertices = [
(0.0, 0.0, 0.0),
(0.0, 0.0, 1.0),
(0.0, 1.0, 0.0),
(0.0, 1.0, 1.0),
(1.0, 0.0, 0.0),
(1.0, 0.0, 1.0),
(1.0, 1.0, 0.0),
(1.0, 1.0, 1.0),
]
julia> cells = [
Int32.((1, 3, 2, 4)),
Int32.((4, 3, 8, 7)),
Int32.((7, 5, 8, 6)),
Int32.((6, 2, 5, 1)),
Int32.((1, 3, 5, 7)),
Int32.((2, 4, 6, 8)),
]
julia> conn = Connectivity{4}(vertices, cells)
'''
6--------8 cell 1 connects vertices 1, 3, 2, 4
/ /| cell 2 connects vertices 4, 3, 8, 7
/ / | .
2--------4 | .
| | | .
| 5 | 7 cell 6 connects vertices 2, 4, 6, 8
| | /
| |/
1--------3P4estTypes.GhostLayer — TypeGhostLayer{X,P}Stores a ghost layer of quadrants that neighbor the domain local to the rank for a Pxest{X}. Also stores the corresponding local domain quadrants, mirrors, that are in other rank's ghost layers.
Fields
pointer: The pointer (of typeP) can be a pointer to either aP4estTypes.P4est.p4est_ghost_tor aP4estTypes.P4est.p8est_ghost_t. See the help documentation for these types for more information about the underlying p4est structures.
See also
ghostlayer: a function used to construct aGhostLayer
P4estTypes.LNodes — TypeLNodes{X,P}Stores a parallel numbering of Lobatto nodes for a Pxest{X}.
Fields
pointer: The pointer (of typeP) can be a pointer to either aP4estTypes.P4est.p4est_lnodes_tor aP4estTypes.P4est.p8est_lnodes_t. See the help documentation for these types for more information about the underlying p4est structures.comm: The MPI Communicator that includes the ranks participating in the lnods.
See also
lnodes: a function used to constructLNodes
P4estTypes.Pxest — TypePxest{X,P,C} <: AbstractArray{P4estTypes.Tree,1}Stores the forest of quadtrees (when X=4) or octrees (when X=8).
This forest of octrees can be accessed in two ways. First, as an array-of-arrays. Each rank holds an array of quadrants for each tree of the Connectivity associated with the forest. (Note, the quadrants are distributed among the ranks. So, each rank will only have access to the quadrants it owns.) Second, using iterateforest to iterate over the volumes, faces, edges, and corners of the forest via callback functions.
Fields
pointer: The pointer (of typeP) can be a pointer to either aP4estTypes.P4est.LibP4est.p4estor aP4estTypes.P4est.LibP4est.p8est. See the help documentation for these types for more information about the underlying p4est structures.connectivity: The connectivity (of typeC) the forest is associated with. This is stored so the connectivity will not be reclaimed by the garbage collector too early.comm: The MPI Communicator that includes the ranks participating in the forest.
See also
pxest: a function that constructs aPxestfrom aConnectivity.iterateforest: a function to iterate over the volumes, faces, edges, and corners of the forest.refine!: refine the quadrants of the forest.coarsen!: coarsen the quadrants of the forest.balance!: two-to-one balance the quadrants of the forest.partition!: partition the quadrants of the forest.ghostlayer: get the ghost layer of quadrants for the forest.lnodes: get a global node numbering.P4estTypes.savevtk: save a VTK representation of the forest.
P4estTypes.QuadrantWrapper — TypeQuadrantWrapper{X,P}Stores a Pxest{X} quadrant (where X=4 indicates a quadrant and X=8 indicates an octant; quadrant is used both as the general term and the term for the 2D object).
Fields
pointer: The pointer (of typeP) can be a pointer to either aP4estTypes.P4est.p4est_quadrantor aP4estTypes.P4est.p8est_quadrant. See the help documentation for these types for more information about the underlying p4est structures.
P4estTypes.Tree — TypeTree{X,P,Q} <: AbstractArray{QuadrantWrapper,1}Stores the quadrants in a tree of a Pxest{X}.
Fields
pointer: The pointer (of typeP) can be a pointer to either aP4estTypes.P4est.p4est_treeor aP4estTypes.P4est.p8est_tree. See the help documentation for these types for more information about the underlying p4est structures.forest: The forest (of typeQ) the tree is associated with. This is stored so the forest will not be reclaimed by the garbage collector too early.
P4estTypes.CONNECT_CORNER — MethodP4estTypes.CONNECT_CORNER(::Val{4})Returns an integer indicating connecting quadrants across faces and corners.
P4estTypes.CONNECT_CORNER — MethodP4estTypes.CONNECT_CORNER(::Val{8})Returns an integer indicating connecting octants across faces, edges, and corners.
P4estTypes.CONNECT_EDGE — MethodP4estTypes.CONNECT_EDGE(::Val{8})Returns an integer indicating connecting octants across faces and edges.
P4estTypes.CONNECT_FACE — MethodP4estTypes.CONNECT_FACE(::Val{4})Returns an integer indicating connecting quadrants across faces.
P4estTypes.CONNECT_FACE — MethodP4estTypes.CONNECT_FACE(::Val{8})Returns an integer indicating connecting octants across faces.
P4estTypes.CONNECT_FULL — MethodP4estTypes.CONNECT_FULL(::Val{4})Returns an integer indicating connecting quadrants across faces and corners.
P4estTypes.CONNECT_FULL — MethodP4estTypes.CONNECT_FULL(::Val{8})Returns an integer indicating connecting octants across faces, edges, and corners.
P4estTypes.balance! — Methodbalance!(forest; kw...)Enforce the two-to-one quadrant size constraint across the forest. By default, this constraint is enforced across faces, edges, and corners.
The keyword arguments (kw...) for the balancing are:
connect: type of constraint enforced which can take the values:P4estTypes.CONNECT_FULL(Val(4)): enforce across face, and corner.P4estTypes.CONNECT_FULL(Val(8)): enforce across face, edge, and corner.P4estTypes.CONNECT_FACE(Val(4)): enforce across face.P4estTypes.CONNECT_FACE(Val(8)): enforce across face.P4estTypes.CONNECT_EDGE(Val(8)): enforce across face and edge.P4estTypes.CONNECT_CORNER(Val(4)): enforce across face and corner.P4estTypes.CONNECT_CORNER(Val(8)): enforce across face, edge, and corner.
init = nothing: callback function with prototypeinit(forest, treeid, quadrant)called for each quadrant created to initialized the user data.replace = nothing: callback function with prototypereplace(forest, treeid, outgoing, incoming)called for eachoutgoingquadrant with their associatedincomingquadrants. Note bothoutgoingandincomingare arrays witheltypeQuadrantWrapper.
See @doc P4estTypes.P4est.p4est_balance_ext and @doc P4estTypes.P4est.p8est_balance_ext for more information about the underlying p4est balance functions.
P4estTypes.brick — Functionbrick(n::NTuple{2, Integer}, p::NTuple{2, Bool}=(false, false))Returns a new Connectivity that is a rectangular n[1]-by-n[2] quadtree connectivity. The brick is periodic in x and y if p[1] and p[2] are true, respectively.
P4estTypes.brick — Functionbrick(l, m, n, p=false, q=false, r=false)Returns a new Connectivity that is a rectangular l-by-m-by-n octree mesh. The brick is periodic in x, y, and z if p, q, and r are true, respectively.
P4estTypes.brick — Functionbrick(n::NTuple{3, Integer}, p::NTuple{3, Bool}=(false, false, false))Returns a new Connectivity that is a rectangular n[1]-by-n[2]-by-n[3] octree mesh. The brick is periodic in x, y, and z if p[1], p[2], and p[3] are true, respectively.
P4estTypes.brick — Functionbrick(l, m, p=false, q=false)Returns a new Connectivity that is a rectangular l-by-m quadtree mesh. The brick is periodic in x and y if p and q are true, respectively.
P4estTypes.coarsen! — Methodcoarsen!(forest; coarsen = (_...) -> false, kw...)Coarsen the quadrants of the forest determined by the coarsen callback. The coarsen(forest, treeid, siblings) callback is called for each set of sibling quadrants local to the rank that are eligible for coarsening. If the callback returns true the siblings will coarsen into one quadrant otherwise they will be untouched.
The other keyword arguments (kw...) for the coarsening are:
recursive = false: iftruecoarsening will be recursive otherwise each set of rank-local siblings will only be visited once.init = nothing: callback function with prototypeinit(forest, treeid, quadrant)called for each quadrant created to initialized the user data.replace = nothing: callback function with prototypereplace(forest, treeid, outgoing, incoming)called for each set ofoutgoingquadrants with their associatedincomingquadrant. Note bothoutgoingandincomingare arrays witheltypeQuadrantWrapper.
See @doc P4estTypes.P4est.p4est_coarsen_ext and @doc P4estTypes.P4est.p8est_coarsen_ext for more information about the underlying p4est coarsening functions.
P4estTypes.comm — Methodcomm(forest)Return the MPI Communicator used by the forest.
P4estTypes.connectivity — Methodconnectivity(forest)Return the Connectivity used by the forest.
P4estTypes.coordinates — Methodcoordinates(quadrant::QuadrantWrapper{4})Returns a tuple of the quadrant's integer coordinates inside its tree.
P4estTypes.coordinates — Methodcoordinates(quadrant::QuadrantWrapper{8})Returns a tuple of the quadrant's integer coordinates inside its tree.
P4estTypes.expand! — Methodexpand!(ghost, forest, nodes)Expand the ghost layer to include all elements with nodes shared with neighboring ranks.
Consider the following forest-of-quadtrees
+-------+---+---+
| | 5 | 6 |
| 2 +---+---+
| | 3 | 4 |
+-------+---+---+
| | |
| 0 | 1 |
| | |
+-------+-------+that is partitioned so rank 0 owns quadrants {0,1} and rank 1 owns quadrants {2, 3, 4, 5, 6}. A fully connected ghost layer on rank 0 would include quadrants {2, 3, 4}. QuadrantWrapper 5 shares a global node with rank 0 but is not in the fully connected ghost layer. This function expands the ghost layer to include quadrants like this.
See sharers to get a list of the global nodes shared with neighboring ranks.
P4estTypes.ghostlayer — Methodghostlayer(forest::Pxest{X}; connection=P4estTypes.CONNECT_FULL(Val(X)))Construct a ghost layer of quadrants that neighbor the local to the rank for the given forest. Here connection determines what neighboring quadrants to include (across face, edge, corner, or full) and can take the values:
P4estTypes.CONNECT_FULL(Val(4)): get face and corner neighbors.P4estTypes.CONNECT_FULL(Val(8)): get face, edge, and corner neighbors.P4estTypes.CONNECT_FACE(Val(4)): get face neighbors.P4estTypes.CONNECT_FACE(Val(8)): get face neighbors.P4estTypes.CONNECT_EDGE(Val(8)): get face and edge neighbors.P4estTypes.CONNECT_CORNER(Val(4)): get face and corner neighbors.P4estTypes.CONNECT_CORNER(Val(8)):get face, edge, and corner neighbors.
P4estTypes.ghosts — Methodghosts(gl::GhostLayer)Returns an array-like structure with the QuadrantWrappers that neighbor the domain of the local rank.
P4estTypes.globalid — Methodglobalid(nodes::LNodes, localid::p4est_locidx_t)Returns the global id associated with the node in LNodes with local id localid.
P4estTypes.iterateforest — Methoditerateforest(forest; kw...)Execute the callbacks passed as keyword arguments for every volume, face, edge, and corner of the rank-local forest.
The keyword arguments (kw...) for the iteration are:
ghost = nothing: the ghost layer associated for the mesh. Used inface,edge, andcornercallbacks for neighboring elements not rank local.volume = nothing: Callback used for every volume (aka quadrant) of the local forest with the prototypevolume(forest, ghost, quadrant, quadid, treeid, userdata).face = nothing: Not implemented yet.edge = nothing: Not implemented yet.corner = nothing: Not implemented yet.userdata = nothing: User data passed to the callbacks.
See @doc P4estTypes.P4est.p4est_iterate and @doc P4estTypes.P4est.p8est_iterate for more information about the iteration.
P4estTypes.lengthofglobalquadrants — Methodlengthofglobalquadrants(forest)Return the number of quadrants distributed across the whole forest.
P4estTypes.lengthoflocalquadrants — Methodlengthoflocalquadrants(forest)Return the number of quadrants local to the rank.
P4estTypes.level — Methodlevel(quadrant::QuadrantWrapper)Returns the level of refinement for the quadrant. Level 0 is the coarsest level and P4estTypes.P4est.P4EST_QMAXLEVEL is the maximum refinement level.
P4estTypes.lnodes — Methodlnodes(forest::Pxest; ghost = nothing, degree = 1)Construct a parallel node numbering for the forest. If the ghost layer ghost is not provided it will be constructed.
A degree > 0 indicates that degree degree Lobotto nodes will be constructed. A degree < 0 indicates that the boundary objects (faces, edges, and corners) will be numbered.
See @doc P4estTypes.P4est.p4est_lnodes_t and @doc P4estTypes.P4est.p8est_lnodes_t for a more detailed discussion of the numbering based on degree.
P4estTypes.mirrors — Methodmirrors(gl::GhostLayer)Returns an array-like structure with the QuadrantWrappers in the local domain that are in neighboring rank's ghost layers.
P4estTypes.offset — Methodoffset(tree::Tree)The cumulative sum of the quadrants over earlier trees on this rank (locals only).
P4estTypes.partition! — Methodpartition!(forest; kw...)Partition the quadrants of the forest. By default this will partition the quadrants equally across the ranks of the forest.
By default sibling elements are split among the ranks. This means they cannot be coarsened with coarsen! and can cause MPI dependent coarsening. If allow_for_coarsening==true then this is avoided by keeping sibling quadrants on the same rank.
A weight(forest, treeid, quadrant) callback may provided (which gives the Float64 weight of each quadrant) for a weighted partitioning of the forest.
Alternatively, the forest may be partitioned to equally distribute the globally numbered nodes via LNodes. This is done by setting lnodes_degree to the node degree. This requires the GhostLayer which if not passed in ghost will be created.
The keyword arguments (kw...) for the partitioning are:
ghost = nothing:GhostLayerused when partitioning byLNodes.lnodes_degree = nothing: partition based onLNodesif this is set to the degree.allow_for_coarsening = false: iftruesibling groups that may be coarsened will be collect on the same rank.weight = nothing: callback that give theFloat64weight of each quadrant to perform a weighted partitioning.
See @doc P4estTypes.P4est.p4est_partition_ext, @doc P4estTypes.P4est.p8est_partition_ext, @doc P4estTypes.P4est.p4est_partition_lnodes, and @doc P4estTypes.P4est.p8est_partition_lnodes, for more information about the underlying p4est partition functions.
P4estTypes.pxest — Methodpxest(connectivity::Connectivity{X}; kw...) where {X}Generate a distributed forest of quadtrees (if X=4) or octrees (if X=8) based on connectivity. Each element of connectivity becomes a tree root.
The connectivity is duplicated on all ranks but the leaves of the forest are split (based on a space-filling curve order) among the ranks.
The keyword arguments (kw...) that control the construction of the forest are:
comm = MPI.COMM_WORLD: the MPI Communicator object of the ranks sharing the forest.min_quadrants = 0: the minimum number of quadrants per rank. (This makes the initial refinement patternMPI.Comm_sizespecific.)min_level = 0: the minimum level of quadrant refinement for the forest.fill_uniform = true: iftruethe forest will be filled with a uniform mesh otherwise it is the coarsest possible mesh.data_type = Nothing: anisbitstypeof the user data stored for each quadrant.init_function = nothing: callback function with prototypeinit_function(forest, treeid, quadrant)called for each quadrant to initialized the user data.
P4estTypes.quadrantstyle — Functionquadrantndims(forest)Returns 2 if the forest quadrants are 2D and 3 if they are 3D.
P4estTypes.quadrantstyle — Methodquadrantstyle(forest)Returns 4 if the forest quadrants are 2D and 8 if they are 3D.
P4estTypes.refine! — Methodrefine!(forest; refine = (_...) -> false, kw...)Refine the quadrants of the forest determined by the refine callback. The refine(forest, treeid, quadrant) callback is called for each quadrant local to the rank. If the callback returns true the quadrant will refine into multiple quadrants otherwise it will be untouched.
The other keyword arguments (kw...) for the refining are:
recursive: iftruerefining will be recursive otherwise each rank-local quadrant will only be visited once.maxlevel = -1: the maximum level of refinement possible during this call.init = nothing: callback function with prototypeinit(forest, treeid, quadrant)called for each quadrant created to initialized the user data.replace = nothing: callback function with prototypereplace(forest, treeid, outgoing, incoming)called for eachoutgoingquadrant with their associatedincomingquadrants. Note bothoutgoingandincomingare arrays witheltypeQuadrantWrapper.
See @doc P4estTypes.P4est.p4est_refine_ext and @doc P4estTypes.P4est.p8est_refine_ext for more information about the underlying p4est refinement functions.
P4estTypes.refine — Methodrefine(c::Connectivity{4}, nedge)Returns a new Connectivity that is c uniformly refined with nedge new trees in each direction.
P4estTypes.refine — Methodrefine(c::Connectivity{8}, nedge)Returns a new Connectivity that is c uniformly refined with nedge new trees in each direction.
P4estTypes.savevtk — Methodsavevtk(prefix, forest; kw...)Save the distributed forest-of-octrees (or quadtrees) forest to a set of VTK files.
A .vtu file with the file name prefix is created per rank storing the rank-local quadrants. Further, .pvtu and .visit collection files are created for ease of importing the mesh into Paraview and Visit, respectively.
The keyword arguments (kw...) are:
scale = 1.0: ascale < 1.0places a visual gap between adjacent quadrants.writetree = true: iftrueinclude the zero-based tree id in VTK cell data.writelevel = true: iftrueinclude the quadrant level in VTK cell data.writerank = true: iftrueinclude the MPI rank in VTK cell data.wraprank = 0: ifwraprank > 0the MPI rank is stored modulowraprank.
P4estTypes.setverbosity — MethodP4estTypes.setverbosity(logpriority::P4estTypes.SC.LP.LogPriority)Sets the verbosity of p4est. See P4estTypes.SC.LP.LogPriority for a list of valid priorities.
P4estTypes.sharers — Methodsharers(nodes::LNodes)Returns a Dict mapping the neighboring rank with the global ids of the nodes shared between it and the local rank.
Note, this Dict does not include an entry for the local rank.
P4estTypes.unsafe_element_nodes — Methodunsafe_element_nodes(nodes::LNodes)Return an array containing the unique continuous node number for each local degree-of-freedom.
See @doc P4estTypes.P4est.p4est_lnodes_t and @doc P4estTypes.P4est.p8est_lnodes_t for a more details.
Note, this unsafely wraps a C array. So, you must ensure that the nodes structure is preserved while using the return value.
P4estTypes.unsafe_face_code — Methodunsafe_face_code(nodes::LNodes)Return an array containing the face code for each quadrant of the mesh. The face code indicates which faces and edges of the quadrant are hanging.
See the p4est functions p4est_lnodes_decode and p8est_lnodes_decode to determine how to decode the face code.
Note, this unsafely wraps a C array. So, you must ensure that the nodes structure is preserved while using the return value.
P4estTypes.unsafe_global_first_quadrant — Methodunsafe_global_first_quadrant(forest::Pxest)Returns 0-based indices into global quadrants. This includes an extra entry at the end of the array so that 1-based range into the global quadrants for rank r can be built with
(global_first_quadrant[r]+1):global_first_quadrant[r+1]Note, this unsafely wraps a C array. So, you must ensure that the forest structure is preserved while using the return value.
P4estTypes.unsafe_global_owned_count — Methodunsafe_global_owned_count(nodes::LNodes)Return an array containing the number of independent nodes owned by each rank.
See @doc P4estTypes.P4est.p4est_lnodes_t and @doc P4estTypes.P4est.p8est_lnodes_t for a more details.
Note, this unsafely wraps a C array. So, you must ensure that the nodes structure is preserved while using the return value.
P4estTypes.unsafe_loaduserdata — Methodunsafe_loaduserdata(quadrant::QuadrantWrapper, type::Type)Return the user data of type type associated with the quadrant.
P4estTypes.unsafe_local_num — Methodunsafe_local_num(quadrant::QuadrantWrapper)Returns the local_num field of the underlying quadrant. This value is only sometimes set so the function is marked unsafe. For example, it is set for quadrants returned by ghosts and mirrors.
P4estTypes.unsafe_mirror_proc_mirrors — Methodunsafe_mirror_proc_mirrors(ghost::GhostLayer)Returns 0-based indices into mirrors. This is used in conjunction with mirror_proc_offsets to get the mirror quadrants associated with each rank. For example
rrange = (mirror_proc_offsets[r]+1):mirror_proc_offsets[r+1]
mirrors(ghost)[mirror_proc_mirrors(rrange)]selects the mirror quadrants associated with rank r.
See @doc P4estTypes.P4est.p4est_ghost_t and @doc P4estTypes.P4est.p8est_ghost_t for a more details.
Note, this unsafely wraps a C array. So, you must ensure that the ghost structure is preserved while using the return value.
P4estTypes.unsafe_mirror_proc_offsets — Methodunsafe_mirror_proc_offsets(ghost::GhostLayer)Returns 0-based indices into mirror_proc_mirrors for each rank. This includes an extra entry at the end of the array so that 1-based range into mirror_proc_mirrors for rank r can be built with
(mirror_proc_offsets[r]+1):mirror_proc_offsets[r+1]See @doc P4estTypes.P4est.p4est_ghost_t and @doc P4estTypes.P4est.p8est_ghost_t for a more details.
Note, this unsafely wraps a C array. So, you must ensure that the ghost structure is preserved while using the return value.
P4estTypes.unsafe_proc_offsets — Methodunsafe_proc_offsets(ghost::GhostLayer)Returns 0-based indices into ghosts for each rank. This includes an extra entry at the end of the array so that 1-based range into ghosts for rank r can be built with
(proc_offsets[r]+1):proc_offsets[r+1]Thus the ghost quadrants associated with rank r can be obtained with
ghosts(ghost)[(proc_offsets[r]+1):proc_offsets[r+1]]See @doc P4estTypes.P4est.p4est_ghost_t and @doc P4estTypes.P4est.p8est_ghost_t for a more details.
Note, this unsafely wraps a C array. So, you must ensure that the ghost structure is preserved while using the return value.
P4estTypes.unsafe_storeuserdata! — Methodunsafe_storeuserdata!(quadrant::QuadrantWrapper, data)Store the user data data associated with the quadrant.
P4estTypes.unsafe_which_tree — Methodunsafe_which_tree(quadrant::QuadrantWrapper)Returns the which_tree field of the underlying quadrant. This value is only sometimes set so the function is marked unsafe.
P4estTypes.uses_mpi — MethodP4estTypes.uses_mpi()Returns true if the p4est C library uses MPI.
P4estTypes.SC.LP.Always — TypeP4estTypes.SC.LP.AlwaysLog priority indicating to log absolutely everything.
P4estTypes.SC.LP.Debug — TypeP4estTypes.SC.LP.DebugLog priority indicating to log any information on the internal state.
P4estTypes.SC.LP.Default — TypeP4estTypes.SC.LP.DefaultThe libsc default log priority.
P4estTypes.SC.LP.Error — TypeP4estTypes.SC.LP.ErrorLog priority indicating to log errors only.
P4estTypes.SC.LP.Essential — TypeP4estTypes.SC.LP.EssentialLog priority indicating to log a few lines at most for a major api function.
P4estTypes.SC.LP.Info — TypeP4estTypes.SC.LP.InfoLog priority indicating to log most relevant things a function is doing.
P4estTypes.SC.LP.LogPriority — TypeP4estTypes.SC.LP.LogPriorityAn abstract type for p4est and libsc log priorities. The follow priorities are available:
P4estTypes.SC.LP.Default: the libsc default.P4estTypes.SC.LP.Always: log absolutely everything.P4estTypes.SC.LP.Trace: prefix file and line number.P4estTypes.SC.LP.Debug: any information on the internal state.P4estTypes.SC.LP.Verbose: information on conditions, decisions.P4estTypes.SC.LP.Info: most relevant things a function is doing.P4estTypes.SC.LP.Statistics: important for consistency/performance.P4estTypes.SC.LP.Essential: a few lines at most for a major api function.P4estTypes.SC.LP.Production: log a few lines max per program.P4estTypes.SC.LP.Error: log errors only.P4estTypes.SC.LP.Silent: never log anything.
P4estTypes.SC.LP.Production — TypeP4estTypes.SC.LP.ProductionLog priority indicating to log a few lines max per program.
P4estTypes.SC.LP.Silent — TypeP4estTypes.SC.LP.SilentLog priority indicating to never log anything.
P4estTypes.SC.LP.Statistics — TypeP4estTypes.SC.LP.StatisticsLog priority indicating to log important for consistency/performance.
P4estTypes.SC.LP.Trace — TypeP4estTypes.SC.LP.TraceLog priority indicating to log trace information with prefix file and line number.
P4estTypes.SC.LP.Verbose — TypeP4estTypes.SC.LP.VerboseLog priority indicating to log information on conditions, decisions.