diff --git a/src/elements.jl b/src/elements.jl index 31f2bef6..d14de3c7 100644 --- a/src/elements.jl +++ b/src/elements.jl @@ -341,6 +341,20 @@ end lift(v::VStruct{T}) where {T} = constructor(v)(pushbefore(v.a, zero(T))) lift(v::AbstractVector{T}) where {T} = pushbefore(v, one(T)) +function MA.promote_operation( + ::typeof(translate), + ::Type{V}, + ::Type{<:AbstractVector{T}}, +) where {S,V<:AbstractVector{S},T} + return similar_type(V, MA.promote_operation(+, S, T)) +end +function MA.mutable_operate!( + ::typeof(translate), + a::AbstractVector, + b::AbstractVector, +) + return MA.mutable_operate!(+, a, b) +end translate(p::AbstractVector, v) = p + v translate(r::VStruct, v) = r diff --git a/src/repelemop.jl b/src/repelemop.jl index 8049d900..35ad3381 100644 --- a/src/repelemop.jl +++ b/src/repelemop.jl @@ -18,14 +18,14 @@ export translate function htranslate(p::HRep{T}, v::Union{AbstractVector{S}}) where {S, T} f = (i, h) -> translate(h, v) d = FullDim(p) - Tout = Base.promote_op(+, T, S) + Tout = MA.promote_operation(+, T, S) similar(p, d, Tout, hmap(f, d, Tout, p)...) end function vtranslate(p::VRep{T}, v::Union{AbstractVector{S}}) where {S, T} f = (i, u) -> translate(u, v) d = FullDim(p) - Tout = Base.promote_op(+, T, S) + Tout = MA.promote_operation(+, T, S) similar(p, d, Tout, vmap(f, d, Tout, p)...) end @@ -54,6 +54,32 @@ function translate(p::Polyhedron, v) end end +function MA.promote_operation( + ::typeof(translate), + P::Type{<:Rep{S}}, + ::Type{<:AbstractVector{T}}, +) where {S,T} + U = MA.promote_operation(+, S, T) + return similar_type(P, U) +end +function MA.mutable_operate!(::typeof(translate), rep::VRepresentation, v) + return mutable_operate_elements!(translate, rep, v) +end +function MA.mutable_operate!(::typeof(translate), rep::HRepresentation, v) + return mutable_operate_elements!(translate, rep, v) +end +function MA.mutable_operate!(::typeof(translate), p::Polyhedron, v) + # The redundancy should not change with translation so we give + # the current redundancy information. + if hrepiscomputed(p) + sethrep!(MA.operate!(translate, hrep(p), v), hredundancy(p)) + end + if vrepiscomputed(p) + # The redundancy should not change with translation + setvrep!(MA.operate!(translate, vrep(p), v), vredundancy(p)) + end +end + ###### # IN # ###### diff --git a/src/vecrep.jl b/src/vecrep.jl index 499f9aae..3b43896e 100644 --- a/src/vecrep.jl +++ b/src/vecrep.jl @@ -180,6 +180,18 @@ convexhull!(v::PointsHull, p::AbstractVector) = push!(v.points, p) @vecrepelem PointsHull Point points +MA.mutability(::Type{<:PointsHull}) = MA.IsMutable() +function mutable_operate_elements!( + op::F, + v::PointsHull, + args::Vararg{Any,N}, +) where {F<:Function,N} + for i in eachindex(v.points) + v.points[i] = MA.operate!(op, v.points[i], args...) + end + return v +end + """ vrep(lines::LineIt, rays::RayIt; d::FullDim) @@ -233,6 +245,25 @@ convexhull!(v::RaysHull, r::Ray) = push!(v.rays, r) @vecrepelem RaysHull Ray rays @subrepelem RaysHull Line lines +MA.mutability(::Type{<:RaysHull}) = MA.IsMutable() +function mutable_operate_elements!( + ::typeof(translate), + ::RaysHull, + ::AbstractVector, +) +end +function mutable_operate_elements!( + op::F, + v::RaysHull, + args::Vararg{Any,N}, +) where {F<:Function,N} + mutable_operate_elements!(op, v.lines, args...) + for i in eachindex(v.rays) + v.rays[i] = MA.operate!(op, v.rays[i], args...) + end + return v +end + vreptype(::Type{RaysHull{T, AT, D}}) where {T, AT, D} = Hull{T, AT, D} """ @@ -279,6 +310,17 @@ convexhull!(v::Hull, r::VStruct) = convexhull!(v.rays, r) @subrepelem Hull Line rays @subrepelem Hull Ray rays +MA.mutability(::Type{<:Hull}) = MA.IsMutable() +function mutable_operate_elements!( + op::F, + v::Hull, + args::Vararg{Any,N}, +) where {F<:Function,N} + mutable_operate_elements!(op, v.points, args...) + mutable_operate_elements!(op, v.rays, args...) + return v +end + fulltype(::Type{<:Union{Hull{T, AT, D}, PointsHull{T, AT, D}, LinesHull{T, AT, D}, RaysHull{T, AT, D}}}) where {T, AT, D} = Hull{T, AT, D} dualtype(::Type{<:Intersection{T, BT, D}}, ::Type{AT}) where {T, AT, BT, D} = Hull{T, AT, D} diff --git a/test/repelemop.jl b/test/repelemop.jl new file mode 100644 index 00000000..f52f7434 --- /dev/null +++ b/test/repelemop.jl @@ -0,0 +1,43 @@ +using Test + +import MutableArithmetics +const MA = MutableArithmetics + +using Polyhedra + +function translate_test(T::Type) + o = one(T) + v = [o, -o] + # We need to create independent copies of `o` in case these are `BigInt`. + pr = convexhull([1o, 2o]) + @test MA.promote_operation(translate, typeof(pr), typeof(v)) == typeof(pr) + @test MA.mutability(typeof(pr)) isa MA.IsMutable + @test MA.mutability(typeof(pr), translate, typeof(pr), typeof(v)) isa MA.IsMutable + rr = conichull([3o, -o]) + @test MA.promote_operation(translate, typeof(rr), typeof(v)) == typeof(rr) + @test MA.mutability(typeof(rr)) isa MA.IsMutable + @test MA.mutability(typeof(rr), translate, typeof(rr), typeof(v)) isa MA.IsMutable + vr = pr + rr + @test MA.promote_operation(translate, typeof(vr), typeof(v)) == typeof(vr) + alloc_test(() -> MA.promote_operation(translate, typeof(vr), typeof(v)), 0) + @test MA.mutability(typeof(vr)) isa MA.IsMutable + @test MA.mutability(typeof(vr), translate, typeof(vr), typeof(v)) isa MA.IsMutable + alloc_test(() -> MA.mutability(typeof(vr), translate, typeof(vr), typeof(v)), 0) + vr2 = MA.operate!(translate, vr, v) + rec_identical(vr, vr2) + alloc_test(() -> MA.mutable_operate!(+, a, b), 0) + alloc_test(() -> MA.mutable_operate!(translate, a, b), 0) + alloc_test(() -> MA.mutable_operate!(translate, vr.points.points[1], v), 0) + alloc_test(() -> MA.mutable_operate!(translate, vr.points, v), 0) + alloc_test(() -> MA.mutable_operate!(translate, vr.rays, v), 0) + alloc_test(() -> MA.mutable_operate!(translate, vr, v), 0) + alloc_test(() -> MA.operate!(translate, vr, v), 0) + +# hr = HalfSpace([1o, 1], 2o) ∩ HyperPlane([-o, 2o, 0o]) +# @test MA.promote_operation(translate, typeof(hr), typeof(v)) == typeof(vr) +# @test MA.mutability(translate, typeof(hr), typeof(v)) isa MA.IsMutable +end + +@testset "Translate" begin + translate_test(Float64) +end diff --git a/test/utils.jl b/test/utils.jl index 5667d213..b06061c0 100644 --- a/test/utils.jl +++ b/test/utils.jl @@ -1,6 +1,35 @@ using Test using Polyhedra +# Allocating size for allocating a `BigInt`. +# Half size on 32-bit. +const BIGINT_ALLOC = Sys.WORD_SIZE == 64 ? 48 : 24 + +function alloc_test(f, n) + f() # compile + @test n == @allocated f() +end + +function rec_identical(vr::Polyhedra.LinesHull, vr2::Polyhedra.LinesHull) + @test vr === vr2 + @test vr.lines === vr2.lines +end +function rec_identical(vr::Polyhedra.RaysHull, vr2::Polyhedra.RaysHull) + @test vr === vr2 + @test vr2.rays === vr.rays + rec_identical(vr.lines, vr2.lines) +end +function rec_identical(vr::Polyhedra.PointsHull, vr2::Polyhedra.PointsHull) + @test vr === vr2 + @test vr2.points === vr.points +end + +function rec_identical(vr::Polyhedra.Hull, vr2::Polyhedra.Hull) + @test vr2 === vr + rec_identical(vr.points, vr2.points) + rec_identical(vr.rays, vr2.rays) +end + _isapprox(x::Real, y::Real) = _isapprox(promote(x, y)...) _isapprox(x::T, y::T) where {T<:Real} = x == y _isapprox(x::T, y::T) where {T<:AbstractFloat} = y < x+1024*eps(T) && x < y+1024*eps(T)