diff --git a/src/toporelations/iterable.jl b/src/toporelations/iterable.jl new file mode 100644 index 000000000..1793d31c0 --- /dev/null +++ b/src/toporelations/iterable.jl @@ -0,0 +1,277 @@ +""" + <:AbstractRelation{Q,P} + +Read as "\$typename Q's of P" + +See also: [`Bounding`](@ref), [`Shared`](@ref), [`Adjacent`](@ref) +""" +abstract type AbstractRelation{Q,P} end + +(T::Type{<:AbstractRelation{Q,P} where {Q,P}})(msh::HalfEdgeMesh) = T(topology(msh)) + +struct TopologyIterator{T<:AbstractRelation} + relation::T + e::HalfEdge +end + +Base.IteratorSize(::Type{<:TopologyIterator}) = Base.SizeUnknown() +Base.eltype(::Type{<:TopologyIterator}) = Int + +struct ParametricDimension{N} end + +const VERTEX = ParametricDimension{0} +const SEGMENT = ParametricDimension{1} +const FACE = ParametricDimension{2} +const CELL = ParametricDimension{3} + +_dim(::Type{ParametricDimension{N}}) where N = N + +struct Bounding{Q<:ParametricDimension,P<:ParametricDimension,T<:HalfEdgeTopology} <: AbstractRelation{Q,P} + t::T + + function Bounding{Q,P}(t::T) where {Q<:ParametricDimension,P<:ParametricDimension,T<:HalfEdgeTopology} + _dim(P) ≤ paramdim(t) || throw(DomainError((P,T), "topology with rank $(ParametricDimension{paramdim(t)}) has no $P. `P` must be ≤ 2")) + _dim(Q) < _dim(P) || throw(DomainError((Q,P), "cannot calculate bounding $Q's of $P")) + P === VERTEX && throw(DomainError(P, "the boundary of a VERTEX is not well-defined")) + + return new{Q,P,T}(t) + end +end + +function (B::Bounding{Q,FACE})(elem::Int) where {Q} + return TopologyIterator(B, half4elem(B.t, elem)) +end + +function (B::Bounding{Q,SEGMENT})(edge::Int) where {Q} + return TopologyIterator(B, half4edge(B.t, edge)) +end + +function Base.iterate(itr::TopologyIterator{<:Bounding{SEGMENT}}) + if isnothing(itr.e.elem) + return nothing + else + return edge4pair(itr.relation.t, (itr.e.head, itr.e.half.head)), itr.e.next + end +end +function Base.iterate(itr::TopologyIterator{<:Bounding{SEGMENT}}, state::HalfEdge) + if state.next === itr.e.next + return nothing + else + return edge4pair(itr.relation.t, (state.head, state.half.head)), state.next + end +end + +function Base.iterate(itr::TopologyIterator{<:Bounding{VERTEX,FACE}}) + if isnothing(itr.e.elem) + return nothing + else + return itr.e.head, itr.e + end +end + +function Base.iterate(itr::TopologyIterator{<:Bounding{VERTEX,FACE}}, state::HalfEdge) + n = state.next + if n === itr.e + return nothing + else + return n.head, n + end +end + +Base.IteratorSize(::Type{<:TopologyIterator{<:Bounding{VERTEX,SEGMENT}}}) = Base.HasLength() +Base.length(::TopologyIterator{<:Bounding{VERTEX,SEGMENT}}) = 2 + +Base.iterate(itr::TopologyIterator{<:Bounding{VERTEX,SEGMENT}}) = itr.e.head, itr.e.half +function Base.iterate(itr::TopologyIterator{<:Bounding{VERTEX,SEGMENT}}, state::HalfEdge) + if state === itr.e + return nothing + else + return state.head, state.half + end +end + +struct Shared{Q<:ParametricDimension,P<:ParametricDimension,T<:HalfEdgeTopology} <: AbstractRelation{Q,P} + t::T + + function Shared{Q,P}(t::T) where {Q<:ParametricDimension,P<:ParametricDimension,T<:HalfEdgeTopology} + _dim(Q) ≤ paramdim(t) || throw(DomainError((Q,T), "topology with rank $(ParametricDimension{paramdim(t)}) has no $Q. `Q` must be ≤ 2")) + _dim(Q) > _dim(P) || throw(DomainError((Q,P), "cannot calculate shared $Q's of $P")) + + return new{Q,P,T}(t) + end +end + +function (S::Shared{Q,FACE})(elem::Int) where {Q} + return TopologyIterator(S, half4elem(S.t, elem)) +end + +function (S::Shared{Q,SEGMENT})(edge::Int) where {Q} + return TopologyIterator(S, half4edge(S.t, edge)) +end + +function (S::Shared{Q,VERTEX})(vert::Int) where {Q} + return TopologyIterator(S, half4vert(S.t, vert)) +end + +function Base.iterate(itr::TopologyIterator{<:Shared{SEGMENT,VERTEX}}) + e = itr.e + si = edge4pair(itr.relation.t, (e.head, e.half.head)) + return si, (e, :prev) +end + +function Base.iterate(itr::TopologyIterator{<:Shared{SEGMENT,VERTEX}}, (e, dir)::Tuple{HalfEdge,Symbol}) + if dir === :prev + if !isnothing(e.elem) + h = e.prev.half + if h.elem != itr.e.elem + si = edge4pair(itr.relation.t, (h.head, h.half.head)) + return si, (h, :prev) + else + return nothing + end + elseif !isnothing(itr.e.half.elem) + h = itr.e.half.next + si = edge4pair(itr.relation.t, (h.head, h.half.head)) + return si, (h, :next) + else + return nothing + end + else # dir === :next + h = e.half + if !isnothing(h.elem) + n = h.next + si = edge4pair(itr.relation.t, (n.head, n.half.head)) + return si, (n, :next) + else + return nothing + end + end +end + +function Base.iterate(itr::TopologyIterator{<:Shared{FACE,VERTEX}}) + if isnothing(itr.e.elem) + h = itr.e.half + return h.elem, (h, :next) + else + h = itr.e + return h.elem, (h, :prev) + end +end + +function Base.iterate(itr::TopologyIterator{<:Shared{FACE,VERTEX}}, (e, dir)::Tuple{HalfEdge,Symbol}) + if dir === :prev + h = e.prev.half + if !isnothing(h.elem) + if h.elem != itr.e.elem + return h.elem, (h, :prev) + else + return nothing + end + elseif !isnothing(itr.e.half.elem) + return itr.e.half.elem, (itr.e.half, :next) + else + return nothing + end + else # dir === :next + h = e.next.half + if !isnothing(h.elem) + return h.elem, (h, :next) + else + return nothing + end + end +end + +function Base.iterate(itr::TopologyIterator{<:Shared{FACE,SEGMENT}}) + e = itr.e + if isnothing(e.elem) + return e.half.elem, e.half + else + return e.elem, e + end +end + +function Base.iterate(itr::TopologyIterator{<:Shared{FACE,SEGMENT}}, e::HalfEdge) + if e === itr.e.half || isnothing(e.half.elem) + return nothing + else + return e.half.elem, e.half + end +end + +struct Adjacent{Q<:ParametricDimension,P<:ParametricDimension,T<:HalfEdgeTopology} <: AbstractRelation{Q,P} + t::T + + function Adjacent{Q,P}(t::T) where {Q<:ParametricDimension,P<:ParametricDimension,T<:HalfEdgeTopology} + _dim(Q) ≤ paramdim(t) || throw(DomainError((Q,T), "topology with rank $(ParametricDimension{paramdim(t)}) has no $Q. `Q` and `P` must be ≤ 2")) + _dim(Q) == _dim(P) || throw(DomainError((Q,P), "cannot calculate adjacent $Q's of $P")) + + return new{Q,P,T}(t) + end +end + +# Convenience constructor to avoid repeating `Q` +Adjacent{Q}(t::HalfEdgeTopology) where Q = Adjacent{Q,Q}(t) +Adjacent{Q}(msh::HalfEdgeMesh) where Q = Adjacent{Q,Q}(topology(msh)) + +function (A::Adjacent{VERTEX})(vert::Int) + return TopologyIterator(A, half4vert(A.t, vert)) +end + +function (A::Adjacent{FACE})(elem::Int) + return TopologyIterator(A, half4elem(A.t, elem)) +end + +function Base.iterate(itr::TopologyIterator{<:Adjacent{VERTEX}}) + e = itr.e + return e.half.head, (e, :prev) +end + +function Base.iterate(itr::TopologyIterator{<:Adjacent{VERTEX}}, (e, dir)::Tuple{HalfEdge,Symbol}) + if dir === :prev + if !isnothing(e.elem) + h = e.prev.half + if h.elem != itr.e.elem + return h.half.head, (h, :prev) + else + return nothing + end + elseif !isnothing(itr.e.half.elem) + h = itr.e.half.next + return h.half.head, (h, :next) + else + return nothing + end + else # dir === :next + h = e.half + if !isnothing(h.elem) + n = h.next + return n.half.head, (n, :next) + else + return nothing + end + end +end + +function Base.iterate(itr::TopologyIterator{<:Adjacent{FACE}}) + if isnothing(itr.e.elem) + return nothing + else + n = itr.e + while isnothing(n.half.elem) + n = n.next + end + return n.half.elem, n.next + end +end +function Base.iterate(itr::TopologyIterator{<:Adjacent{FACE}}, state::HalfEdge) + while isnothing(state.half.elem) && state !== itr.e + state = state.next + end + if state === itr.e + return nothing + else + return state.half.elem, state.next + end +end + diff --git a/test/toporelations_iterators.jl b/test/toporelations_iterators.jl new file mode 100644 index 000000000..d09e2fcd4 --- /dev/null +++ b/test/toporelations_iterators.jl @@ -0,0 +1,302 @@ +@testset "HalfEdgeTopology topological iteration" begin +# +# 3 ------------ 4 +# / \ / +# / \ 2 / +# / \ / +# / \ / +# / 1 \ / +# / \ / +# 1 ------------ 2 +# This is a simple 2 triangle topology + simple = HalfEdgeTopology(connect.([(1, 2, 3), (4, 3, 2)])) + +# 6 +# / \ +# / \ +# / \ +# / 3 \ +# / \ +# / \ +# 1 ------------ 3 +# / \ / \ +# / \ 1 / \ +# / \ / \ +# / \ / \ +# / 4 \ / 2 \ +# / \ / \ +# 5 ------------ 2 ------------ 4 +# This topology has one triangle (1) which is fully surrounded (neighbors on every edge) + triforce = HalfEdgeTopology(connect.([(1, 2, 3), (4, 3, 2), (1, 3, 6), (1, 2, 5)])) + +# 3 +# / | \ +# / | \ +# / | \ +# / | \ +# / 1 | 2 \ +# / | \ +# 1 ------2------ 4 +# \ | / +# \ 4 | 3 / +# \ | / +# \ | / +# \ | / +# \ | / +# 5 +# This topology has one vertex (2) which is fully surrounded (the connected edges are a cycle) + diamond = HalfEdgeTopology(connect.([(1, 2, 3), (4, 3, 2), (2, 4, 5), (1, 2, 5)])) + + nonnothing_elem = x -> !isnothing(x.elem) + @testset "Bounding" begin + @test_throws DomainError Bounding{FACE,CELL}(simple) + @test_throws DomainError Bounding{FACE,FACE}(simple) + @test_throws DomainError Bounding{ParametricDimension{-1},VERTEX}(simple) + + @testset "Bounding{SEGMENT,FACE}" begin + ∂ = Boundary{2,1}(simple) + B = Bounding{SEGMENT,FACE}(simple) + for i in 1:nfaces(simple,2) + @test issetequal(B(i), ∂(i)) + end + @test isnothing(iterate(TopologyIterator(B, half4pair(simple, (1,2)).half))) + for e in filter(nonnothing_elem, simple.halfedges) + @test issetequal(TopologyIterator(B, e), ∂(e.elem)) + end + + ∂ = Boundary{2,1}(triforce) + B = Bounding{SEGMENT,FACE}(triforce) + for i in 1:nfaces(triforce,2) + @test issetequal(B(i), ∂(i)) + end + for e in filter(nonnothing_elem, triforce.halfedges) + @test issetequal(TopologyIterator(B, e), ∂(e.elem)) + end + + ∂ = Boundary{2,1}(diamond) + B = Bounding{SEGMENT,FACE}(diamond) + for i in 1:nfaces(diamond,2) + @test issetequal(B(i), ∂(i)) + end + for e in filter(nonnothing_elem, diamond.halfedges) + @test issetequal(TopologyIterator(B, e), ∂(e.elem)) + end + end + + @testset "Bounding{VERTEX,FACE}" begin + ∂ = Boundary{2,0}(simple) + B = Bounding{VERTEX,FACE}(simple) + for i in 1:nfaces(simple,2) + @test issetequal(B(i), ∂(i)) + end + @test isnothing(iterate(TopologyIterator(B, half4pair(simple, (1,2)).half))) + for e in filter(nonnothing_elem, simple.halfedges) + @test issetequal(TopologyIterator(B, e), ∂(e.elem)) + end + + ∂ = Boundary{2,0}(triforce) + B = Bounding{VERTEX,FACE}(triforce) + for i in 1:nfaces(triforce,2) + @test issetequal(B(i), ∂(i)) + end + for e in filter(nonnothing_elem, triforce.halfedges) + @test issetequal(TopologyIterator(B, e), ∂(e.elem)) + end + + ∂ = Boundary{2,0}(diamond) + B = Bounding{VERTEX,FACE}(diamond) + for i in 1:nfaces(diamond,2) + @test issetequal(B(i), ∂(i)) + end + for e in filter(nonnothing_elem, diamond.halfedges) + @test issetequal(TopologyIterator(B, e), ∂(e.elem)) + end + end + + @testset "Bounding{VERTEX,SEGMENT}" begin + ∂ = Boundary{1,0}(simple) + B = Bounding{VERTEX,SEGMENT}(simple) + for i in 1:nfacets(simple) + @test issetequal(B(i), ∂(i)) + end + for e in simple.halfedges + @test issetequal(TopologyIterator(B, e), ∂(edge4pair(simple, (e.head, e.half.head)))) + end + + ∂ = Boundary{1,0}(triforce) + B = Bounding{VERTEX,SEGMENT}(triforce) + for i in 1:nfacets(triforce) + @test issetequal(B(i), ∂(i)) + end + for e in triforce.halfedges + @test issetequal(TopologyIterator(B, e), ∂(edge4pair(triforce, (e.head, e.half.head)))) + end + + ∂ = Boundary{1,0}(diamond) + B = Bounding{VERTEX,SEGMENT}(diamond) + for i in 1:nfacets(diamond) + @test issetequal(B(i), ∂(i)) + end + for e in diamond.halfedges + @test issetequal(TopologyIterator(B, e), ∂(edge4pair(diamond, (e.head, e.half.head)))) + end + end + end + + @testset "Shared" begin + @test_throws DomainError Bounding{CELL,FACE}(simple) + @test_throws DomainError Bounding{FACE,FACE}(simple) + + @testset "Shared{SEGMENT,VERTEX}" begin + 𝒞 = Coboundary{0,1}(simple) + S = Shared{SEGMENT,VERTEX}(simple) + for i in 1:nvertices(simple) + @test issetequal(S(i), 𝒞(i)) + end + for e in simple.halfedges + @test issetequal(TopologyIterator(S, e), 𝒞(e.head)) + end + + 𝒞 = Coboundary{0,1}(triforce) + S = Shared{SEGMENT,VERTEX}(triforce) + for i in 1:nvertices(triforce) + @test issetequal(S(i), 𝒞(i)) + end + for e in triforce.halfedges + @test issetequal(TopologyIterator(S, e), 𝒞(e.head)) + end + + 𝒞 = Coboundary{0,1}(diamond) + S = Shared{SEGMENT,VERTEX}(diamond) + for i in 1:nvertices(diamond) + @test issetequal(S(i), 𝒞(i)) + end + for e in diamond.halfedges + @test issetequal(TopologyIterator(S, e), 𝒞(e.head)) + end + end + + @testset "Shared{FACE,VERTEX}" begin + 𝒞 = Coboundary{0,2}(simple) + S = Shared{FACE,VERTEX}(simple) + for i in 1:nvertices(simple) + @test issetequal(S(i), 𝒞(i)) + end + for e in simple.halfedges + @test issetequal(TopologyIterator(S, e), 𝒞(e.head)) + end + + 𝒞 = Coboundary{0,2}(triforce) + S = Shared{FACE,VERTEX}(triforce) + for i in 1:nvertices(triforce) + @test issetequal(S(i), 𝒞(i)) + end + for e in triforce.halfedges + @test issetequal(TopologyIterator(S, e), 𝒞(e.head)) + end + + 𝒞 = Coboundary{0,2}(diamond) + S = Shared{FACE,VERTEX}(diamond) + for i in 1:nvertices(diamond) + @test issetequal(S(i), 𝒞(i)) + end + for e in diamond.halfedges + @test issetequal(TopologyIterator(S, e), 𝒞(e.head)) + end + end + + @testset "Shared{FACE,SEGMENT}" begin + 𝒞 = Coboundary{1,2}(simple) + S = Shared{FACE,SEGMENT}(simple) + for i in 1:nfacets(simple) + @test issetequal(S(i), 𝒞(i)) + end + for (e,i) in zip(simple.halfedges, repeat(1:nfacets(simple), inner=2)) + @test issetequal(TopologyIterator(S, e), 𝒞(i)) + end + + 𝒞 = Coboundary{1,2}(triforce) + S = Shared{FACE,SEGMENT}(triforce) + for i in 1:nfacets(triforce) + @test issetequal(S(i), 𝒞(i)) + end + for (e,i) in zip(triforce.halfedges, repeat(1:nfacets(triforce), inner=2)) + @test issetequal(TopologyIterator(S, e), 𝒞(i)) + end + + 𝒞 = Coboundary{1,2}(diamond) + S = Shared{FACE,SEGMENT}(diamond) + for i in 1:nfacets(diamond) + @test issetequal(S(i), 𝒞(i)) + end + for (e,i) in zip(diamond.halfedges, repeat(1:nfacets(diamond), inner=2)) + @test issetequal(TopologyIterator(S, e), 𝒞(i)) + end + end + end + + @testset "Adjacent" begin + @test_throws DomainError Adjacent{CELL,FACE}(simple) + @test_throws DomainError Adjacent{VERTEX,SEGMENT}(simple) + + @testset "Adjacent{VERTEX,VERTEX}" begin + 𝒜 = Adjacency{0}(simple) + A = Adjacent{VERTEX}(simple) + for i in 1:nvertices(simple) + @test issetequal(A(i), 𝒜(i)) + end + for e in simple.halfedges + @test issetequal(TopologyIterator(A, e), 𝒜(e.head)) + end + + 𝒜 = Adjacency{0}(triforce) + A = Adjacent{VERTEX}(triforce) + for i in 1:nvertices(triforce) + @test issetequal(A(i), 𝒜(i)) + end + for e in triforce.halfedges + @test issetequal(TopologyIterator(A, e), 𝒜(e.head)) + end + + 𝒜 = Adjacency{0}(diamond) + A = Adjacent{VERTEX}(diamond) + for i in 1:nvertices(diamond) + @test issetequal(A(i), 𝒜(i)) + end + for e in diamond.halfedges + @test issetequal(TopologyIterator(A, e), 𝒜(e.head)) + end + end + + @testset "Adjacent{FACE,FACE}" begin + 𝒜 = Adjacency{2}(simple) + A = Adjacent{FACE}(simple) + for i in 1:nfaces(simple,2) + @test issetequal(A(i), 𝒜(i)) + end + @test isnothing(iterate(TopologyIterator(A, half4pair(simple, (1,2)).half))) + for e in filter(nonnothing_elem, simple.halfedges) + @test issetequal(TopologyIterator(A, e), 𝒜(e.elem)) + end + + 𝒜 = Adjacency{2}(triforce) + A = Adjacent{FACE}(triforce) + for i in 1:nfaces(triforce,2) + @test issetequal(A(i), 𝒜(i)) + end + for e in filter(nonnothing_elem, triforce.halfedges) + @test issetequal(TopologyIterator(A, e), 𝒜(e.elem)) + end + + 𝒜 = Adjacency{2}(diamond) + A = Adjacent{FACE}(diamond) + for i in 1:nfaces(diamond,2) + @test issetequal(A(i), 𝒜(i)) + end + for e in filter(nonnothing_elem, diamond.halfedges) + @test issetequal(TopologyIterator(A, e), 𝒜(e.elem)) + end + end + end +end +