Define a full polytope structure

[serenity4]

So far we never define a full polytope structure, but maybe it would be nice to have one. So that if you're interested in an object with all its vertices, edges, faces and so on (essentially, a mesh in 3D) without relying on conventions, there's a type available and you don't have to pass around all k-faces individually. We don't even need to encode the value of the vertices if we stick to an abstract polytope (which only describes connectivity, but not data).

[juliohm]

It would be really nice to see how a prototype implementation of this. I wonder how it scales and what the public API would look like. Looking forward to learning more about it.

[serenity4]

For simple shapes, when it's possible to use a mapping to find all k-faces from vertices, it's probably going to just be inefficient memory-wise. Although, if you don't really care about efficiency at that level, it would be easier for querying faces should you the indices you want. One good thing though is that the resulting polytope would be portable, particularly when the polytope is complex: in computer graphics, you have control over all faces: vertices, edges, "faces" so if you want to keep the geometry you need to save everything.

Implementing it for any number of dimensions might be tricky. In 3D, it might look something like

struct DensePolyhedron{F<:AbstractVector{<: Polytope{2, Dim, T}}, E<:AbstractVector{<: Polytope{1, Dim, T}}, V<:AbstractVector{<: Polytope{0, Dim, T}}, Dim, T} <: Polyhedron{Dim, T}
    faces::F
    edges::E
    vertices::V
end

which, if you only have the same polytope type inside each k-face (only Triangles, Segments and Points), may be OK, but may quickly become a mess if you mix several k-polytopes types for k-faces. Taking simplices, though, would look reasonable for some applications. Maybe, instead of using N fields for N-polytopes, you could use a N-tuple where each item is a list:

struct DensePolyhedron{F<:AbstractVector{<: Polytope{2, Dim, T}}, E<:AbstractVector{<: Polytope{1, Dim, T}}, V<:AbstractVector{<: Polytope{0, Dim, T}}, Dim, T} <: Polyhedron{Dim, T}
    faces::Tuple{F, E, V} # in the mathematical sense (0-, 1-, 2-faces)
end

so that the number of fields stays fixed. But still, we would probably need N+2 type parameters at least to keep it efficient. I don't know, maybe that can be a problem. Particularly if you consider high dimensions, but then you also just have the problem that the number of k-faces (to describe useful geometry) would be insane.

[serenity4]

I think it might be worth it to see if we could simply use graphs as an internal representation. For example, if you have a polyhedron with V vertices, E edges and F faces:

• you begin from a list of vertex indices 1:V, instancing a first graph G1
• you form E edges between the graph vertices of G1, that represent the geometric edges
• you create a new graph G2 from edge indices 1:E
• you form F edges between the graph vertices of G2, that represent the geometric "faces"

and so on for higher polytopes. This would be good since it would allow reuse of k-face indices to construct (k+1)-faces, which is convenient when you have a 3D geometry that shares a lot of edges/faces for example.

[serenity4]

Just to give an idea of the number of vertices/edges/faces of 3D Voronoi triangulations, here are some plots (produced from p.46 of this source; when only bounds are available, I took the average of min and max bounds):

Gist: https://gist.github.com/serenity4/7783951c517ecee5e40fb90768434c41

[juliohm]

Nice plots! Sorry for the hiatus today. I was finishing a tutorial and I am finally free again to jump back into Meshes.jl development.

[serenity4]

So, after looking up UnstructuredMesh in detail, this one in theory allows a full polytope structure. It's just that people need to fill the connec field with connectivity information regarding all faces, and that vertex data is required.
I believe though that with a hierarchical graph-based approach, it would be closer to the mathematical theory and more efficient for parametric dimensions > 3.

[juliohm]

I was considering the graph-like structure as well. It is very common in other mesh implementations out there (e.g. GMSH) as it allows efficient traversal of ranks. From what I understand we only need to store the connectivities of the highest-dimensional polytopes (e.g. tetrahedron) and provide adjacency and incidence lists to navigate the ranks.

After I finish the triangulation algorithm, I will have more time to think about it. Feel free to share proposals and we can brainstorm optimal data structures.

[serenity4]

Sorry for the delay, I've been working on some geometric algebra package in parallel (something that would ressemble Grassmann.jl, MIT-licensed, more user-friendly, but probably- hopefully not by much- less performant).

As for proposals, I've quickly drafted an implementation

using LightGraphs

abstract type AbstractPolytope{N} end

struct FullPolytope{N} <: AbstractPolytope{N}
    graph::NTuple{N,SimpleDiGraph}
end

"""
Build a `FullPolytope` with `n_verts` from the structure provided in `indices`.

The nth argument passed as `indices` describes connectivity from the (n-1)th face to the nth face.

## Example

`` `julia
julia> FullPolytope(4,
    [(1, 2), (2, 3), (3, 1), (2, 3), (3, 4), (4, 2)], # edges (segments)
    [(1, 2, 3), (4, 5, 6)] # 2 triangle faces from edges 1, 2, 3 and 4, 5, 6 respectively
)

julia> FullPolytope(4,
    [(1, 2, 3, 1), (2, 3, 4, 2)], # edges (closed chains)
    [(1,), (2,)]) # faces (made from a closed chain -> PolySurface with no inner chains)
`` `
"""
function FullPolytope(n_verts, indices::Vararg{<:AbstractVector{<:Tuple}})
    graphs = SimpleDiGraph[]
    n_nodes = n_verts
    g = SimpleDiGraph(n_verts)
    for vec in indices
        for (edge_src, edge_dests) in enumerate(vec)
            for dest in edge_dests
                add_edge!(g, edge_src, dest)
            end
        end
        push!(graphs, g)
        n_nodes = length(vec)
        g = SimpleDiGraph(n_nodes)
    end
    FullPolytope{length(indices)}(tuple(graphs...))
end

paramdim(::Type{FullPolytope{N}}) where {N} = N

### Benchmarks
using BenchmarkTools

@benchmark FullPolytope(4, [(1, 2), (2, 3), (3, 1), (2, 3), (3, 4), (4, 2)], [(1, 2, 3), (4, 5, 6)]) # ~7 µs
@benchmark FullPolytope(4, [(1, 2, 3, 1), (2, 3, 4, 2)], [(1,), (2,)]) # ~12 µs

It introduces LightGraphs as a dependency. The reason for that is that it would be useful for querying graph properties (e.g. determine all cycles of the graph, which indicate chain structures- chain of vertices, a.k.a. Chain-like, or chains of edges, faces or any polytope). We can remove it if we don't want to exploit the graph structure, though.
This would suggest a rather different approach than the one taken so far, using only combinatorial information (which is why the AbstractPolytope type comes back!). I find it interesting from a mathematical viewpoint, but the current structure of the package is not adapted to that.

Right now it's quite slow, indices tuples can be of any length so type inference struggles a bit. If we restrict the length of each tuple, or provide constructors for some predefined polytopes, it should speed it up. However since these "predefined polytopes" (e.g. Triangle, Chain) currently use vertex information, we'd need abstract polytope implementations instead (except if we go for the approach given below).

A probably better approach, in that it would follow the same logic as the rest of the package, would be to declare FullPolytope as a Polytope with additional Dim and T parameters, and an additional vertices field.

By the way, I don't know about the name FullPolytope, I was also thinking of DensePolytope, but ultimately it's just a polytope that explicits all its structure.

[juliohm]

Thanks for sharing! I think we want to strike a full polytope structure as a mesh with incidence and adjacency lists as I mentioned, without explicitly creating graph objects. I have a guess that the code will look much simpler that way, but we need to evolve these ideas over time.

What package you are working on for geometric algebra? It would be nice to contribute and learn there as well. :)

[serenity4]

I think we want to strike a full polytope structure as a mesh with incidence and adjacency lists

Are those incidence and adjacency lists the same as implemented in LightGraphs ? In that case, if we don't want to exploit the structure of the graph we can use them yes. Else, I think reusing what LightGraphs provides could reduce the amount of code we have to write (things like checking for cycles, shortest path between two nodes - edges, vertices, faces...). All that because a polytope can be represented as a graph (in fact it simplifies a lot of things when we think about it as a graph). IIRC LightGraphs is quite fast so I'd love to see this graph structure exploited in its purest form and still retain performance. It would be nice to see how both approaches perform (using lists or using a graph).

Sure! It's at https://github.com/serenity4/GeometricAlgebra.jl. It's still young and far from polished though. I'll put links to what I found were good introductions to geometric algebra tomorrow, so check out the README if you're interested ;)

[juliohm]

I am looking forward to playing more with geometric algebra. :) Are you on Zulip as well? We could exchange messages there too.

[juliohm]

I think this issue overlaps with #87 . The ultimate goal is to be able to represent the connectivities efficiently without copying vertices. There is a lot of literature on this, and I will try to address things as I progress with my research.

Closing this one for now.

[serenity4]

Yes and no. This issue was not about defining an optimal polytope structure. It was about defining a full structure without implicit assumptions so that the poset it represents can be exploited. There is a potential for doing arbitrary queries on the structure of the polytope without having to switch between (efficient) data structures.

[juliohm]

Can you compare this full structure you have in mind to what we currently have in SimpleMesh?

…

On Tue, Mar 30, 2021, 18:25 Cédric Belmant ***@***.***> wrote: Yes and no. This issue was not about defining an optimal polytope structure. It was about defining a full structure without implicit assumptions so that the poset it represents can be exploited. There is a potential for doing arbitrary queries on the structure of the polytope without having to switch between (efficient) data structures. — You are receiving this because you modified the open/close state. Reply to this email directly, view it on GitHub <#12 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/AAZQW3JKS4LPHJGAXU6LSBTTGI6VNANCNFSM4TPESHMQ> .

[serenity4]

SimpleMesh allows any structure but it is highly inconvenient to work with in a general setting. For example, in a polytope a 2-face is defined by its relation to 1-faces, but the current connectivity logic only associates k-faces to vertices. If I want to find face loops (an extension of Chain but for general faces, not just 0-faces aka vertices), I don't see any convenient way to do that in SimpleMesh. If I know k-faces from their relations to their (k-1)-faces, it's a lot easier.

[juliohm]

Makes total sense and I think some of the data structures I mentioned will help with that too.

…

On Tue, Mar 30, 2021, 18:36 Cédric Belmant ***@***.***> wrote: SimpleMesh allows any structure but it is highly inconvenient to work with in a general setting. A 2-face for example is defined by its relation to 1-faces, but the current connectivity logic only associates k-faces to vertices. If I want to find face loops (an extension of Chain but for general faces, not just 0-faces aka vertices), I don't see any convenient way to do that in SimpleMesh. If I know k-faces from their relations to their (k-1)-faces, it's a lot easier. — You are receiving this because you modified the open/close state. Reply to this email directly, view it on GitHub <#12 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/AAZQW3MNDRMZYEYL4RALHV3TGI745ANCNFSM4TPESHMQ> .