-
Notifications
You must be signed in to change notification settings - Fork 87
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Introducing iteration on vertices #1133
base: master
Are you sure you want to change the base?
Conversation
Thank you @ghyatzo for the contribution. To facilitate the usage of this iterator, we just need to implement a new For the |
Hi, yes absolutely can do. Regarding flatten, I've done some tries, but the
flatten requires an iterator of iterators, which means that we have to build an iterator for each subelement in the multigeometry. |
@ghyatzo I meant that we simply need to define the function eachvertex(p::Polytope) = (v for v in vertices(p)) The geometry that will benefit the most is the |
The reason why I didn't use simply a generator It should still be generic enough since underneat we are calling Please let me know if there is anything I overlooked or wrong with my reasoning. |
I've improved the flatten approach using a generator of generators and not an array as I was naively doing before. julia> eachvertex_custom(e::Union{Mesh,Polytope,Multi}) = VertexItr(e)
julia> eachvertex_gen(p::Polytope) = (v for v in vertices(p))
julia> eachvertex_flat(e::Multi) = Iterators.flatten(eachvertex_gen(g) for g in e.geoms)
eachvertex_flat (generic function with 1 method)
julia> m
MultiNgon
├─ Triangle((x: 0.0119639 m, y: 0.139913 m), (x: 0.347587 m, y: 0.958666 m), (x: 0.145339 m, y: 0.109886 m))
└─ Hexagon((x: 0.25018 m, y: 0.85677 m), ..., (x: 0.57332 m, y: 0.3497 m))
julia> @btime begin
s = 0
for _ in eachvertex_flat($m)
s+=1
end
s
end
1.104 μs (39 allocations: 2.95 KiB)
9
julia> @btime begin
s = 0
for _ in eachvertex_custom($m)
s+=1
end
s
end
497.851 ns (9 allocations: 288 bytes)
9 those 9 allocations are from an Anyway, I completely understand prefering the simpler solution. Also: is it enough to just have the |
Isn't ok to define eachvertex(m::Multi) = (v for g in m.geoms for v in vertices(g)) in the |
is shorthand for
from lowering the code:
and has similar performances
|
Yes, but please restrict to
Interesting. Thanks for clarifying this. |
Appreciate if you can keep it simple for now without an auxiliary iterator struct. Please define the We can reconsider the auxiliary iterator struct in the future. I understand it allocates a fixed amount no matter the size of the |
Alright, will do. shall I make a separate PR and keep this as draft with the struct approach? |
Feel free to update the PR here directly. You can copy the VertexIter struct in the comments to save it for future reference. |
For future reference: struct VertexItr{T}
el::T
end
_v_iterate(el::Polytope, i) =
(@inline; (i - 1) % UInt < nvertices(el) % UInt ? (@inbounds vertex(el, i), i + 1) : nothing)
_v_iterate(el::Mesh, i) =
(@inline; (i - 1) % UInt < nvertices(el) % UInt ? (@inbounds vertex(el, i), i + 1) : nothing)
Base.iterate(itr::VertexItr{<:Mesh}, i=1) = _v_iterate(itr.el, i)
Base.iterate(itr::VertexItr{<:Polytope}, i=1) = _v_iterate(itr.el, i)
Base.iterate(itr::VertexItr{<:MultiPolytope}, state=(1, 1)) = begin
ig, ivg = state
ig > length(itr.el.geoms) && return nothing
is = _v_iterate(itr.el.geoms[ig], ivg)
is === nothing && return Base.iterate(itr, (ig + 1, 1))
v, ivg = is
return (v, (ig, ivg))
end
Base.length(itr::VertexItr{<:Mesh}) = nvertices(T)
Base.length(itr::VertexItr{<:Polytope}) = nvertices(T)
Base.length(itr::VertexItr{<:MultiPolytope}) = sum(nvertices, itr.el.geoms)
Base.IteratorSize(itr::VertexItr) = Base.HasLength()
Base.eltype(::VertexItr) = Point
eachvertex(e::T) where {T} = error("Vertex iterator not implemented for type $T")
eachvertex(e::Union{Mesh,Polytope,MultiPolytope}) = VertexItr(e) |
I would like some pointers WRT the new tests, I'm a bit unsure where best to put them I've also took the liberty of rewriting some |
You can add tests in the corresponding test files. For instance,
I shared a comment above about that. Just a test that uses the
Does it lead to the same type inference? I know that list comprehensions are very efficient, and that the compiler can do magic with them. I am not sure if this is the case with the |
I'll do some tests, but I expect it does, since I have the sneaking suspicion that doing Edit:
|
I am not sure how to exclude the allocation of the object being iterated over, like one can do with
but
another example:
|
The vertices of the |
Seems like all the allocations arrive from here: function iterate(g::Generator, s...)
@inline
y = iterate(g.iter, s...)
y === nothing && return nothing
y = y::Tuple{Any, Any} # try to give inference some idea of what to expect about the behavior of the next line
return (g.f(y[1]), y[2]) # <--- this line here in particular
end which is the code for iterating generators in Base. in particular, I've seen these allocations using julia> try ttt() catch err; err.errors[1] end
Allocating runtime call to "jl_f_getfield" in ./Base.jl:49
| getproperty(x, f::Symbol) = (@inline; getfield(x, f))
Stacktrace:
[1] getproperty
@ ./Base.jl:49 [inlined]
[2] iterate
@ ./generator.jl:48 [inlined]
[3] var"##ttt#446"()
@ Main ./REPL[190]:2
julia> try ttt() catch err; err.errors[2] end
Allocating runtime call to "jl_get_nth_field_checked" in ./tuple.jl:31
| getindex(@nospecialize(t::Tuple), i::Int) = getfield(t, i, @_boundscheck)
Stacktrace:
[1] getindex
@ ./tuple.jl:31 [inlined]
[2] iterate
@ ./generator.jl:48 [inlined]
[3] var"##ttt#446"()
@ Main ./REPL[190]:2
julia> try ttt() catch err; err.errors[3] end
Dynamic dispatch in ./generator.jl:48
| return (g.f(y[1]), y[2])
Stacktrace:
[1] iterate
@ ./generator.jl:48 [inlined]
[2] var"##ttt#446"()
@ Main ./REPL[190]:2
julia> try ttt() catch err; err.errors[4] end
Allocating runtime call to "jl_f_getfield" in ./Base.jl:49
| getproperty(x, f::Symbol) = (@inline; getfield(x, f))
Stacktrace:
[1] getproperty
@ ./Base.jl:49 [inlined]
[2] iterate
@ ./generator.jl:45 [inlined]
[3] var"##ttt#446"()
@ Main ./REPL[190]:3
julia> try ttt() catch err; err.errors[5] end
Dynamic dispatch to function iterate in ./generator.jl:45
| y = iterate(g.iter, s...)
Stacktrace:
[1] iterate
@ ./generator.jl:45 [inlined]
[2] var"##ttt#446"()
@ Main ./REPL[190]:3
julia> try ttt() catch err; err.errors[6] end
Allocating runtime call to "jl_f_getfield" in ./Base.jl:49
| getproperty(x, f::Symbol) = (@inline; getfield(x, f))
Stacktrace:
[1] getproperty
@ ./Base.jl:49 [inlined]
[2] iterate
@ ./generator.jl:48 [inlined]
[3] var"##ttt#446"()
@ Main ./REPL[190]:3
julia> try ttt() catch err; err.errors[7] end
Allocating runtime call to "jl_get_nth_field_checked" in ./tuple.jl:31
| getindex(@nospecialize(t::Tuple), i::Int) = getfield(t, i, @_boundscheck)
Stacktrace:
[1] getindex
@ ./tuple.jl:31 [inlined]
[2] iterate
@ ./generator.jl:48 [inlined]
[3] var"##ttt#446"()
@ Main ./REPL[190]:3
julia> try ttt() catch err; err.errors[8] end
Dynamic dispatch in ./generator.jl:48
| return (g.f(y[1]), y[2])
Stacktrace:
[1] iterate
@ ./generator.jl:48 [inlined]
[2] var"##ttt#446"()
@ Main ./REPL[190]:3 not completely sure what is going on. Seems like the compiler is giving up in trying to infer stuff and fallbacks to runtime allocation. In this case |
That is very interesting. The If that is the case, we have a good motive to introduce the custom iterator in a new Thank you for the deep analysis of memory usage, it really helps with these design choices. 💯 |
I am a doofus, I was checking the allocations using non const global variables, that is why it was struggling. julia> @allocated begin
for v in eachvertex(T) # const Triangle
end
end
0
julia> @allocated begin
for v in eachvertex(H) # const Hexagon
end
end
0
julia> @allocated begin
for v in eachvertex(M) # const Multi([T, H])
end
end
3312 # <--- this still allocates.
julia> @allocated begin
for v in eachvertex(RG) # const RegularGrid
end
end
0
julia> @allocated begin
for v in eachvertex(SG) # const StructuredGrid
end
end
0 The issues still seems to be the compiler choking on the nested generator at some point. |
But if the nested generators work well with non global variables, we are set, right? It is simpler code with the desired behavior. |
Yeah, but the code I shared above is with constants, not global variables. julia> forloops(g) = for v in eachvertex(g) end
julia> @allocated forloops(m)
3024 |
Got it. Thanks for checking!
Em sex., 15 de nov. de 2024, 12:01, cschen ***@***.***>
escreveu:
… Yeah, but the code I shared above is with constants, not global variables.
And it still allocates. Also inside a function barrier:
julia> forloops(g) = for v in eachvertex(g) end
julia> @allocated forloops(m)3024
—
Reply to this email directly, view it on GitHub
<#1133 (comment)>,
or unsubscribe
<https://github.com/notifications/unsubscribe-auth/AAZQW3NX7H7AXDKNBQVYMSD2AYEGLAVCNFSM6AAAAABRVEG6CCVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDINZZGA4DOMZSHE>
.
You are receiving this because you commented.Message ID:
***@***.***>
|
I have been hacking away at it, but I can't seem to get rid of the type instability in the the compiler correctly infers that Which incurs in dynamic dispatch calls to find the number of vertices through I've tried multiple approaches but I fear that without adding some extra information in the structure itself to help the compiler Couldn't find a way around it yet, let me know if you have any insights. |
This should be fine for most practical purposes. The case where Multi
becomes relevant is with PolyArea, a different type of Ngon that is not
aware of the number N of vertices at compile time. If you convert both the
Triangle and the Hexagon to PolyArea you get zero allocations?
Although we support Multi of Ngon, and other types of Multi not common
elsewhere, the main use case is with PolyArea.
Em sáb., 16 de nov. de 2024, 19:23, cschen ***@***.***>
escreveu:
… I have been hacking away at it, but I can't seem to get rid of the type
instability in the Multi iteration.
for example, using as an example m = Multi([T, H]) where T is a Triangle
and H is a Hexagon.
the compiler correctly infers that m is of type Multi{M, C, <:Ngon}, but
when iterating over the geometries it gives up on infering the type of each
geometry. It knows that m.geoms is of type Vector{Ngon{N}} where N but it
stops there.
Which incurs in dynamic dispatch calls to find the number of vertices
through nvertices and therefore allocations (although much less than the
generators).
I've tried multiple approaches but I fear that without adding some extra
information in the structure itself to help the compiler
having complete inference of all elements inside the vector might be hard.
Couldn't find a way around it yet, let me know if you have any insights.
—
Reply to this email directly, view it on GitHub
<#1133 (comment)>,
or unsubscribe
<https://github.com/notifications/unsubscribe-auth/AAZQW3N2JXY5OYRMNQNUMZL2A7AVBAVCNFSM6AAAAABRVEG6CCVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDIOBQHAZTANRZGA>
.
You are receiving this because you commented.Message ID:
***@***.***>
|
Hi,
following #1127 and #1131 I'm starting this draft PR.
I've added the main iteration logic for single geometry and multi geometry types.
This si incomplete as I would require some input on the following points
in particular, I am not very familiar with all the nomenclature and type system structure used in this package, so I would appreciate an indication of which types we are interested in and/or for which onese we want to be able to call
vertices
.WRT to point 2, I've gone for stuffing all iteration related code into a single file for now, let me know if this has to change.
cheers!