From f564afb233862b25fa874a8b6a4ad62e214d70ca Mon Sep 17 00:00:00 2001 From: Jan Weidner Date: Fri, 5 May 2017 22:06:51 +0200 Subject: [PATCH 1/4] convexhulls work again --- src/convexhulls.jl | 51 ++++++++++++++++++++++------------------------ test/REQUIRE | 1 - test/runtests.jl | 2 +- 3 files changed, 25 insertions(+), 29 deletions(-) delete mode 100644 test/REQUIRE diff --git a/src/convexhulls.jl b/src/convexhulls.jl index 42ac2b9..4ee1165 100644 --- a/src/convexhulls.jl +++ b/src/convexhulls.jl @@ -1,5 +1,3 @@ - - Base.eltype(fg::AFG) = eltype(typeof(fg)) Base.length(fg::AFG) = length(vertices(fg)) nvertices(fg::AFG) = length(fg) @@ -19,7 +17,7 @@ Base.deleteat!(c::AbstractFlexibleGeometry, i) = (deleteat!(vertices(c), i); c) Base.copy{FG <: AFG}(fl::FG) = FG(copy(vertices(fl))) push(fl::AFG, pt) = push!(copy(fl), pt) -vertices{T}(x::AbstractFlexibleGeometry{T}) = x._ +vertices(x::AbstractFlexibleGeometry) = x._ vertices(s::Simplex) = Tuple(s) standard_cube_vertices(::Type{Val{1}}) = [Vec(0), Vec(1)] @@ -41,37 +39,36 @@ end _combine_vcat(vert_1, vert_last) end end +@generated function vertices_rettype(r::HyperRectangle) + N = ndims(r) + T = eltype(r) + RET = NTuple{(2^N), Vec{N,T}} + :($RET) +end -@generated function vertices{N,T}(r::HyperRectangle{N,T}) - ret_type = NTuple{(2^N), Vec{N,T}} - quote - o = origin(r) - v = widths(r) - f(sv) = o + sv .* v - tuple(map(f, standard_cube_vertices(Val{N}))...)::$ret_type - end +function vertices(r::HyperRectangle) + N = ndims(r) + ret_type = vertices_rettype(r) + o = origin(r) + v = widths(r) + f(sv) = o + sv .* v + tuple(map(f, standard_cube_vertices(Val{N}))...)::ret_type end vertices(c::HyperCube) = vertices(convert(HyperRectangle, c)) -#vertices(s::AbstractConvexHull) = Tuple(s) - -vertexmat{M, T}(s::Simplex{M, T}) = Mat{length(T)}(vcat(vertices(s)...)) -vertexmat{M}(s::AbstractGeometry{M}) = Mat{M}(vcat(vertices(s)...)) -function vertexmat(s::AbstractFlexibleGeometry) - verts = vcat(vertices(s)...) - M, N = spacedim(s), nvertices(s) - Mat{M, N}(verts) -end +vertexmat(s) = hcat(vertices(s)...) vertexmatrix(s::AbstractConvexHull) = Matrix(vertexmat(s))::Matrix{numtype(s)} -convert{S <: Simplex}(::Type{S}, fs::FlexibleSimplex) = S(tuple(vertices(fs)...)) -convert{F <: AFG}(::Type{F}, s::Simplex) = F(collect(vertices(s))) -convert{FS <: FlexibleSimplex}(::Type{FS}, f::FS) = f -convert{FG <: AFG, FS <: FlexibleSimplex}(::Type{FG}, f::FS) = FG(vertices(f)) -convert{R <: HyperRectangle}(::Type{R}, c::HyperCube) = R(origin(c), widths(c)) -convert{F <: FlexibleConvexHull}(::Type{F}, s::Simplex) = F(collect(vertices(s))) -convert{F <: FlexibleConvexHull}(::Type{F}, c) = F(collect(vertices(c))) +convert(S::Type{<: Simplex}, fs::FlexibleSimplex) = S(tuple(vertices(fs)...)) +convert(F::Type{<: AFG}, s::Simplex) = F(collect(vertices(s))) +convert(FG::Type{<: AFG}, f::FlexibleSimplex) = FG(vertices(f)) +convert(R::Type{<: HyperRectangle}, c::HyperCube) = R(origin(c), widths(c)) +convert(F::Type{<:FlexibleConvexHull}, s::Simplex) = F(collect(vertices(s))) +convert(F::Type{<:FlexibleConvexHull}, s::FlexibleSimplex) = F(vertices(s)) +convert(F::Type{FlexibleConvexHull}, c::FlexibleConvexHull) = c +convert(F::Type{<:FlexibleConvexHull}, c::FlexibleConvexHull) = F(vertices(c)) +convert(F::Type{<:FlexibleConvexHull}, c) = F(collect(vertices(c))) function Base.isapprox(s1::AbstractConvexHull, s2::AbstractConvexHull;kw...) isapprox(vertexmat(s1), vertexmat(s2); kw...) diff --git a/test/REQUIRE b/test/REQUIRE deleted file mode 100644 index 94e516f..0000000 --- a/test/REQUIRE +++ /dev/null @@ -1 +0,0 @@ -BaseTestNext diff --git a/test/runtests.jl b/test/runtests.jl index 949e218..01f1104 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -4,6 +4,7 @@ import Base.Test.@inferred @testset "GeometryTypes" begin + include("convexhulls.jl") include("polygons.jl") include("hyperrectangles.jl") include("faces.jl") @@ -17,7 +18,6 @@ import Base.Test.@inferred include("simplices.jl") include("lines.jl") include("polygons.jl") - #include("convexhulls.jl") #include("gjk.jl") include("cylinder.jl") end From c3318e59bce82de418e5d31606f27e2d2fbef148 Mon Sep 17 00:00:00 2001 From: Jan Weidner Date: Fri, 5 May 2017 23:53:13 +0200 Subject: [PATCH 2/4] gjk passes again --- src/convexhulls.jl | 1 + src/gjk.jl | 30 +++++++++++++++++++------- src/simplices.jl | 14 ++++++------- test/gjk.jl | 52 ++++++++++++++++++++++------------------------ test/runtests.jl | 2 +- 5 files changed, 56 insertions(+), 43 deletions(-) diff --git a/src/convexhulls.jl b/src/convexhulls.jl index 4ee1165..fa1d5ae 100644 --- a/src/convexhulls.jl +++ b/src/convexhulls.jl @@ -68,6 +68,7 @@ convert(F::Type{<:FlexibleConvexHull}, s::Simplex) = F(collect(vertices(s))) convert(F::Type{<:FlexibleConvexHull}, s::FlexibleSimplex) = F(vertices(s)) convert(F::Type{FlexibleConvexHull}, c::FlexibleConvexHull) = c convert(F::Type{<:FlexibleConvexHull}, c::FlexibleConvexHull) = F(vertices(c)) +convert(::Type{F}, x::F) where {F <: FlexibleConvexHull} = x convert(F::Type{<:FlexibleConvexHull}, c) = F(collect(vertices(c))) function Base.isapprox(s1::AbstractConvexHull, s2::AbstractConvexHull;kw...) diff --git a/src/gjk.jl b/src/gjk.jl index cf69cd1..feb5fb3 100644 --- a/src/gjk.jl +++ b/src/gjk.jl @@ -3,19 +3,34 @@ using Base.Cartesian -type_immutable{T, n}(::Type{FlexibleSimplex{T}}, ::Type{Val{n}}) = Simplex{n,T} +type_immutable(::Type{FlexibleSimplex{T}}, ::Type{Val{n}}) where {T,n}= Simplex{n,T} type_immutable(s, v=Val{length(s)}) = type_immutable(typeof(s), v) make_immutable(s, v=Val{length(s)}) = convert(type_immutable(s, v), s) +#@generated function with_immutable{n}(f, s, max_depth::Type{Val{n}}=Val{10}) +# quote +# l = length(s) +# @nif $n d->(d == l) d->(f(make_immutable(s, Val{d}))) d->error("length(s) = $s < max_depth = {n}") +# end +#end + """ with_immutable{n}(f, s, max_depth::Type{Val{n}}=Val{10}) Apply function f to immutable variant of s without introducing a type instability. """ -@generated function with_immutable{n}(f, s, max_depth::Type{Val{n}}=Val{10}) - quote - l = length(s) - @nif $n d->(d == l) d->(f(make_immutable(s, Val{d}))) d->error("length(s) = $s < max_depth = {n}") +function with_immutable(f, s) + l = length(s) + if l == 1 + f(make_immutable(s, Val{1})) + elseif l == 2 + f(make_immutable(s, Val{2})) + elseif l == 3 + f(make_immutable(s, Val{3})) + elseif l == 4 + f(make_immutable(s, Val{4})) + else + f(make_immutable(s, Val{l})) end end @@ -80,7 +95,7 @@ end support_vector(ch, v) = support_vector_max(ch,v)[1] """ -s = shrink!(s::FlexibleSimplex, pt, atol=0.) + s = shrink!(s::FlexibleSimplex, pt, atol=0.) Drop as many vertices from s as possible, while keeping pt inside. """ @@ -96,7 +111,7 @@ end """ -pt_best, dist = gjk0(c) + pt_best, dist = gjk0(c) Compute the distance of a convex set to the origin using gjk algorithm. If c is not convex, the algorithm may or may not yield an optimal solution. @@ -132,7 +147,6 @@ function gjk0(c, push!(sim, wk) pt_best::T, sqd = proj_sqdist(zero_point, sim) sqd == 0 && return pt_best, norm(pt_best) - #@show weights(pt_best, sim) shrink!(sim, pt_best, atol) end end diff --git a/src/simplices.jl b/src/simplices.jl index 69f7551..ef677c8 100644 --- a/src/simplices.jl +++ b/src/simplices.jl @@ -94,7 +94,7 @@ translation(s::Simplex) = first(vertices(s)) v1 = vert[1] end tupl = Expr(:tuple) - for i = 2:m + @inbounds for i = 2:m sym = Symbol("diff$i") push!(init.args, :($sym = vert[$i]-v1)) append!(tupl.args, (:($sym[$j]) for j = 1:N)) @@ -115,7 +115,7 @@ simplex_face(s::Simplex, i::Int) = Simplex(deleteat(vertices(s), i)) """ proj_sqdist{T}(pt::T, s::LineSegment{T}, best_sqd=eltype(T)(Inf)) """ -function proj_sqdist{T}(pt::T, s::LineSegment{T}, best_sqd=eltype(T)(Inf)) +function proj_sqdist{T}(pt::T, s::LineSegment{T}, best_sqd=Inf) v0, v1 = vertices(s) pt = pt - v0 v = v1 - v0 @@ -143,7 +143,7 @@ distance is greater then best_sqd, the behaviour of pt_proj is not defined. distance is greater then best_sqd is returned instead. """ -function proj_sqdist{nv,T}(pt::T, s::Simplex{nv,T}, best_sqd = eltype(T)(Inf)) +function proj_sqdist{nv,T}(pt::T, s::Simplex{nv,T}, best_sqd=Inf) w = weights(pt, s) best_proj = Vec(vertexmat(s) * w) # at this point best_proj lies in the subspace spanned by s, @@ -152,7 +152,7 @@ function proj_sqdist{nv,T}(pt::T, s::Simplex{nv,T}, best_sqd = eltype(T)(Inf)) if sqd >= best_sqd # pt is far away even from the subspace spanned by s return best_proj, best_sqd elseif any(w .< 0) # pt is closest to point inside a face of s - @inbounds for i in 1:length(w) + for i in 1:length(w) if w[i] < 0 proj, sqd = proj_sqdist(pt, simplex_face(s, i), best_sqd) if sqd < best_sqd @@ -167,17 +167,17 @@ function proj_sqdist{nv,T}(pt::T, s::Simplex{nv,T}, best_sqd = eltype(T)(Inf)) end end -sqdist(pt, s, best = Inf) = proj_sqdist(pt, s, best)[2] +sqdist(pt, s, best=Inf) = proj_sqdist(pt, s, best)[2] """ proj_sqdist(p::Vec, q::Vec, best_sqd=eltype(p)(Inf)) """ -@inline function proj_sqdist(p::Vec, q::Vec, best_sqd = eltype(p)(Inf)) +@inline function proj_sqdist(p::Vec, q::Vec, best_sqd=Inf) q, min(best_sqd, sqnorm(p-q)) end """ proj_sqdist{T}(pt::T, s::Simplex{1, T}, best_sqd=eltype(T)(Inf)) """ -@inline function proj_sqdist{T}(pt::T, s::Simplex{1, T}, best_sqd=eltype(T)(Inf)) +@inline function proj_sqdist{T}(pt::T, s::Simplex{1, T}, best_sqd=Inf) proj_sqdist(pt, translation(s), best_sqd) end diff --git a/test/gjk.jl b/test/gjk.jl index 3b60bd7..ab18f9b 100644 --- a/test/gjk.jl +++ b/test/gjk.jl @@ -1,6 +1,5 @@ -import GeometryTypes: ⊖, support_vector_max -import GeometryTypes: gjk0 import GeometryTypes: type_immutable, with_immutable, make_immutable +import GeometryTypes: gjk0, ⊖, support_vector_max @testset "gjk" begin @@ -11,35 +10,35 @@ import GeometryTypes: type_immutable, with_immutable, make_immutable @test gjk(c1,c2) ≈5 @test min_euclidean(c1,c2) ≈5 - c1 = Simplex(Vec(-1.,0,0)) - c2 = Simplex(Vec(4.,0,0)) + c1 = Simplex(Vec(-1.,0.,0.)) + c2 = Simplex(Vec(4.,0.,0.)) @test gjk(c1,c2) ≈5 @test min_euclidean(c1,c2) ≈5 - c1 = FlexibleConvexHull([Vec(0.,0), Vec(0.,1), Vec(1.,0),Vec(1.,1)]) + c1 = FlexibleConvexHull([Vec(0.,0.), Vec(0.,1.), Vec(1.,0.),Vec(1.,1.)]) c2 = Simplex(Vec(4.,0.5)) - @test gjk(c1, c2) ≈3 - @test min_euclidean(c1,c2) ≈3 + @test gjk(c1, c2) ≈ 3 + @test min_euclidean(c1,c2) ≈ 3 - pt1 = Vec(1,2,3.) - pt2 = Vec(3,4,5.) + pt1 = Vec(1.,2.,3.) + pt2 = Vec(3.,4.,5.) @test gjk(pt1, pt2) ≈norm(pt1-pt2) @test min_euclidean(pt1, pt2) ≈norm(pt1-pt2) end @testset "gjk intersecting lines" begin - c1 = Simplex(Vec(1,1.), Vec(1, 2.)) + c1 = Simplex(Vec(1.,1.), Vec(1., 2.)) @test gjk(c1, c1) == 0. @test min_euclidean(c1,c1) == 0. - c2 = Simplex(Vec(1,1.), Vec(10, 2.)) + c2 = Simplex(Vec(1.,1.), Vec(10., 2.)) @test gjk(c1, c2) == 0. @test min_euclidean(c1,c2) == 0. - c3 = Simplex(Vec(0, 1.), Vec(2,2.)) + c3 = Simplex(Vec(0., 1.), Vec(2.,2.)) md = vertices(c1 ⊖ c3) @test md == [Vec(1.0,0.0),Vec(-1.0,-1.0),Vec(1.0,1.0),Vec(-1.0,0.0)] - @test gjk0(FlexibleConvexHull(md)) == (Vec(0, 0.), 0.) + @test gjk0(FlexibleConvexHull(md)) == (Vec(0., 0.), 0.) @test gjk(c1, c3) == 0. @test min_euclidean(c1,c3) == 0. @@ -47,13 +46,13 @@ import GeometryTypes: type_immutable, with_immutable, make_immutable @testset "Cube" begin c = HyperCube(Vec(0.5,0.5,0.5), 1.) - @test min_euclidean(Vec(2,2,2.), c) ≈ gjk(Vec(2,2,2.), c) ≈ √(3/4) + @test min_euclidean(Vec(2.,2.,2.), c) ≈ gjk(Vec(2.,2.,2.), c) ≈ √(3/4) - s = Simplex(Vec(1, 0.5, 0.5), Vec(1,2,3.)) + s = Simplex(Vec(1., 0.5, 0.5), Vec(1.,2.,3.)) @test 0 <= min_euclidean(s, c) <= 1e-14 @test 0 <= gjk(s, c) <= 1e-14 - s = Simplex(Vec(2,2,2.), Vec(2,3,2.), Vec(2,2,7.), Vec(3,4,5.)) + s = Simplex(Vec(2.,2.,2.), Vec(2.,3.,2.), Vec(2.,2.,7.), Vec(3.,4.,5.)) @test min_euclidean(s,c) ≈ gjk(s, c) ≈ √(3/4) end @@ -61,18 +60,18 @@ import GeometryTypes: type_immutable, with_immutable, make_immutable @testset "support_vector_max" begin r = HyperRectangle(Vec(-0.5, -1.), Vec(1., 2.)) - @test support_vector_max(r, Vec(1,0.)) == (Vec(0.5,-1.), 0.5) - @test support_vector_max(r, Vec(2,0.)) == (Vec(0.5,-1.), 1.) - @test support_vector_max(r, Vec(-1,0.)) == (Vec(-0.5,-1.), 0.5) - @test support_vector_max(r, Vec(0, 1.)) == (Vec(-0.5,1.), 1.) - @test support_vector_max(r, Vec(1, 1.)) == (Vec(0.5,1.), 1.5) - @test support_vector_max(FlexibleConvexHull(r), Vec(1, 1.)) == (Vec(0.5,1.), 1.5) - - c1 = Simplex(Vec(1,1.), Vec(1, 2.)) - c3 = Simplex(Vec(0, 1.), Vec(2,2.)) + @test support_vector_max(r, Vec(1.,0.)) == (Vec(0.5,-1.), 0.5) + @test support_vector_max(r, Vec(2.,0.)) == (Vec(0.5,-1.), 1.) + @test support_vector_max(r, Vec(-1.,0.)) == (Vec(-0.5,-1.), 0.5) + @test support_vector_max(r, Vec(0., 1.)) == (Vec(-0.5,1.), 1.) + @test support_vector_max(r, Vec(1., 1.)) == (Vec(0.5,1.), 1.5) + @test support_vector_max(FlexibleConvexHull(r), Vec(1., 1.)) == (Vec(0.5,1.), 1.5) + + c1 = Simplex(Vec(1.,1.), Vec(1., 2.)) + c3 = Simplex(Vec(0., 1.), Vec(2.,2.)) md = c1 ⊖ c3 fh = FlexibleConvexHull(md) - for v in [Vec(1,0.), Vec(12,-10.), Vec(0,-1.), Vec(1,1.)] + for v in [Vec(1.,0.), Vec(12.,-10.), Vec(0.,-1.), Vec(1.,1.)] v_md, s_md = support_vector_max(md, v) v_fh, s_fh = support_vector_max(fh, v) @test s_md ≈ s_fh @@ -80,7 +79,6 @@ import GeometryTypes: type_immutable, with_immutable, make_immutable end end - @testset "make immutable" begin T = Vec{2, Float64} n = 3 diff --git a/test/runtests.jl b/test/runtests.jl index 01f1104..d2be602 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -16,8 +16,8 @@ import Base.Test.@inferred include("hypersphere.jl") include("typeutils.jl") include("simplices.jl") + include("gjk.jl") include("lines.jl") include("polygons.jl") - #include("gjk.jl") include("cylinder.jl") end From 45599c744bd53b1835cec678525bdd6c53d834a0 Mon Sep 17 00:00:00 2001 From: Jan Weidner Date: Fri, 5 May 2017 23:59:46 +0200 Subject: [PATCH 3/4] fix --- src/simplices.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/simplices.jl b/src/simplices.jl index ef677c8..3776a0a 100644 --- a/src/simplices.jl +++ b/src/simplices.jl @@ -94,7 +94,7 @@ translation(s::Simplex) = first(vertices(s)) v1 = vert[1] end tupl = Expr(:tuple) - @inbounds for i = 2:m + for i = 2:m sym = Symbol("diff$i") push!(init.args, :($sym = vert[$i]-v1)) append!(tupl.args, (:($sym[$j]) for j = 1:N)) @@ -152,7 +152,7 @@ function proj_sqdist{nv,T}(pt::T, s::Simplex{nv,T}, best_sqd=Inf) if sqd >= best_sqd # pt is far away even from the subspace spanned by s return best_proj, best_sqd elseif any(w .< 0) # pt is closest to point inside a face of s - for i in 1:length(w) + @inbounds for i in 1:length(w) if w[i] < 0 proj, sqd = proj_sqdist(pt, simplex_face(s, i), best_sqd) if sqd < best_sqd From 950c794fbf4030de13831ba13037dada378c9e37 Mon Sep 17 00:00:00 2001 From: Jan Weidner Date: Sat, 6 May 2017 00:03:45 +0200 Subject: [PATCH 4/4] replace _ fieldnames --- src/convexhulls.jl | 2 +- src/types.jl | 4 ++-- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/src/convexhulls.jl b/src/convexhulls.jl index fa1d5ae..6f585e0 100644 --- a/src/convexhulls.jl +++ b/src/convexhulls.jl @@ -17,7 +17,7 @@ Base.deleteat!(c::AbstractFlexibleGeometry, i) = (deleteat!(vertices(c), i); c) Base.copy{FG <: AFG}(fl::FG) = FG(copy(vertices(fl))) push(fl::AFG, pt) = push!(copy(fl), pt) -vertices(x::AbstractFlexibleGeometry) = x._ +vertices(x::AbstractFlexibleGeometry) = x.vertices vertices(s::Simplex) = Tuple(s) standard_cube_vertices(::Type{Val{1}}) = [Vec(0), Vec(1)] diff --git a/src/types.jl b/src/types.jl index 8d4d483..dcbca13 100644 --- a/src/types.jl +++ b/src/types.jl @@ -191,7 +191,7 @@ FlexibleConvexHull{T} Represents the convex hull of finitely many points. The number of points is not fixed. """ immutable FlexibleConvexHull{T} <: AFG{T} - _::Vector{T} + vertices::Vector{T} end """ @@ -200,7 +200,7 @@ FlexibleSimplex{T} Represents a Simplex whos number of vertices is not fixed. """ immutable FlexibleSimplex{T} <: AFG{T} - _::Vector{T} + vertices::Vector{T} end """