From 70eaf6fc1e533b54928bf8aaf3002db535994fec Mon Sep 17 00:00:00 2001 From: schillic Date: Wed, 2 Aug 2023 19:31:14 +0200 Subject: [PATCH 1/6] add radius_ball function --- docs/src/lib/sets/Ball1.md | 1 + docs/src/lib/sets/Ball2.md | 1 + docs/src/lib/sets/BallInf.md | 1 + docs/src/lib/sets/Ballp.md | 1 + src/Sets/Ball1.jl | 17 +++++++++++++++++ src/Sets/Ball2.jl | 17 +++++++++++++++++ src/Sets/BallInf.jl | 17 +++++++++++++++++ src/Sets/Ballp.jl | 17 +++++++++++++++++ test/Sets/Ball1.jl | 3 +++ test/Sets/Ball2.jl | 6 ++++++ test/Sets/BallInf.jl | 3 +++ test/Sets/Ballp.jl | 7 ++++++- 12 files changed, 90 insertions(+), 1 deletion(-) diff --git a/docs/src/lib/sets/Ball1.md b/docs/src/lib/sets/Ball1.md index e5d411be91..eb1fa4d9d7 100644 --- a/docs/src/lib/sets/Ball1.md +++ b/docs/src/lib/sets/Ball1.md @@ -11,6 +11,7 @@ Ball1 ∈(::AbstractVector, ::Ball1, ::Bool=false) vertices_list(::Ball1) center(::Ball1) +radius_ball(::Ball1) rand(::Type{Ball1}) constraints_list(::Ball1) translate(::Ball1, ::AbstractVector) diff --git a/docs/src/lib/sets/Ball2.md b/docs/src/lib/sets/Ball2.md index 7ccc558b96..a4dc6f2c1e 100644 --- a/docs/src/lib/sets/Ball2.md +++ b/docs/src/lib/sets/Ball2.md @@ -10,6 +10,7 @@ Ball2 σ(::AbstractVector, ::Ball2) ∈(::AbstractVector, ::Ball2) center(::Ball2) +radius_ball(::Ball2) rand(::Type{Ball2}) sample(::Ball2{N}, ::Int) where {N} translate(::Ball2, ::AbstractVector) diff --git a/docs/src/lib/sets/BallInf.md b/docs/src/lib/sets/BallInf.md index 09cfc07416..5143875ba9 100644 --- a/docs/src/lib/sets/BallInf.md +++ b/docs/src/lib/sets/BallInf.md @@ -7,6 +7,7 @@ CurrentModule = LazySets ```@docs BallInf center(::BallInf) +radius_ball(::BallInf) radius(::BallInf, ::Real=Inf) radius_hyperrectangle(::BallInf) radius_hyperrectangle(::BallInf, ::Int) diff --git a/docs/src/lib/sets/Ballp.md b/docs/src/lib/sets/Ballp.md index fa5af55a31..655c295631 100644 --- a/docs/src/lib/sets/Ballp.md +++ b/docs/src/lib/sets/Ballp.md @@ -10,6 +10,7 @@ Ballp ρ(::AbstractVector, ::Ballp) ∈(::AbstractVector, ::Ballp) center(::Ballp) +radius_ball(::Ballp) rand(::Type{Ballp}) translate(::Ballp, ::AbstractVector) translate!(::Ballp, ::AbstractVector) diff --git a/src/Sets/Ball1.jl b/src/Sets/Ball1.jl index 5cd306bb10..34bd72d449 100644 --- a/src/Sets/Ball1.jl +++ b/src/Sets/Ball1.jl @@ -74,6 +74,23 @@ function center(B::Ball1) return B.center end +""" + radius_ball(B::Ball1) + +Return the ball radius of a ball in the 1-norm. + +### Input + +- `B` -- ball in the 1-norm + +### Output + +The ball radius. +""" +function radius_ball(B::Ball1) + return B.radius +end + """ vertices_list(B::Ball1) diff --git a/src/Sets/Ball2.jl b/src/Sets/Ball2.jl index a2d74139e8..2d39d0e812 100644 --- a/src/Sets/Ball2.jl +++ b/src/Sets/Ball2.jl @@ -86,6 +86,23 @@ function center(B::Ball2) return B.center end +""" + radius_ball(B::Ball2) + +Return the ball radius of a ball in the 2-norm. + +### Input + +- `B` -- ball in the 2-norm + +### Output + +The ball radius. +""" +function radius_ball(B::Ball2) + return B.radius +end + """ ρ(d::AbstractVector, B::Ball2) diff --git a/src/Sets/BallInf.jl b/src/Sets/BallInf.jl index 9e68d56bae..0297a9b216 100644 --- a/src/Sets/BallInf.jl +++ b/src/Sets/BallInf.jl @@ -156,6 +156,23 @@ function center(B::BallInf) return B.center end +""" + radius_ball(B::BallInf) + +Return the ball radius of a ball in the infinity norm. + +### Input + +- `B` -- ball in the infinity norm + +### Output + +The ball radius. +""" +function radius_ball(B::BallInf) + return B.radius +end + """ σ(d::AbstractVector, B::BallInf) diff --git a/src/Sets/Ballp.jl b/src/Sets/Ballp.jl index 0fc25d354c..0e8caff21e 100644 --- a/src/Sets/Ballp.jl +++ b/src/Sets/Ballp.jl @@ -94,6 +94,23 @@ function center(B::Ballp) return B.center end +""" + radius_ball(B::Ballp) + +Return the ball radius of a ball in the p-norm. + +### Input + +- `B` -- ball in the p-norm + +### Output + +The ball radius. +""" +function radius_ball(B::Ballp) + return B.radius +end + """ σ(d::AbstractVector, B::Ballp) diff --git a/test/Sets/Ball1.jl b/test/Sets/Ball1.jl index d9a0bb792c..c0f6091d41 100644 --- a/test/Sets/Ball1.jl +++ b/test/Sets/Ball1.jl @@ -59,6 +59,9 @@ for N in [Float64, Rational{Int}, Float32] # center @test center(b) == N[0, 0] + # radius_ball + @test LazySets.radius_ball(b) == N(2) + # boundedness @test isbounded(b) diff --git a/test/Sets/Ball2.jl b/test/Sets/Ball2.jl index eb3783527d..795ef636f1 100644 --- a/test/Sets/Ball2.jl +++ b/test/Sets/Ball2.jl @@ -61,6 +61,12 @@ for N in [Float64, Float32] # unicode constructor @test ○(center(b), radius(b)) == b + # center + @test center(b) == N[0, 0] + + # radius_ball + @test LazySets.radius_ball(b) == N(2) + # boundedness @test isbounded(b) && isboundedtype(typeof(b)) diff --git a/test/Sets/BallInf.jl b/test/Sets/BallInf.jl index 18e130c265..2809bd9cfa 100644 --- a/test/Sets/BallInf.jl +++ b/test/Sets/BallInf.jl @@ -62,6 +62,9 @@ for N in [Float64, Rational{Int}, Float32] b = BallInf(c, N(2)) @test center(b) == c && center(b, 1) == N(1) && center(b, 2) == N(2) + # radius_ball + @test LazySets.radius_ball(b) == N(2) + # support vector for single entry vector svec = σ(SingleEntryVector(2, 3, N(2)), BallInf(zeros(N, 3), N(2))) @test svec[1] ∈ Interval(N[-2, 2]) && svec[2] == N(2) && svec[3] ∈ Interval(N[-2, 2]) diff --git a/test/Sets/Ballp.jl b/test/Sets/Ballp.jl index ba4365f226..a170c3356d 100644 --- a/test/Sets/Ballp.jl +++ b/test/Sets/Ballp.jl @@ -34,7 +34,12 @@ for N in [Float64, Float32] @test Ballp(N(Inf), N[3], N(4)) isa BallInf # center - @test center(b) == b.center + @test center(b) == N[0, 0] + + # radius_ball + @test LazySets.radius_ball(b) == N(2) + + # an_element @test an_element(b) isa AbstractVector{N} # boundedness From ea91911d96fc8a0f57d68b75fa32b67b9cc04139 Mon Sep 17 00:00:00 2001 From: schillic Date: Thu, 3 Aug 2023 06:30:59 +0200 Subject: [PATCH 2/6] add ball_norm function --- docs/src/lib/sets/Ball1.md | 5 +++-- docs/src/lib/sets/Ball2.md | 5 +++-- docs/src/lib/sets/BallInf.md | 1 + docs/src/lib/sets/Ballp.md | 5 +++-- src/Sets/Ball1.jl | 18 ++++++++++++++++++ src/Sets/Ball2.jl | 18 ++++++++++++++++++ src/Sets/BallInf.jl | 18 ++++++++++++++++++ src/Sets/Ballp.jl | 17 +++++++++++++++++ test/Sets/Ball1.jl | 3 +++ test/Sets/Ball2.jl | 3 +++ test/Sets/BallInf.jl | 3 +++ test/Sets/Ballp.jl | 3 +++ 12 files changed, 93 insertions(+), 6 deletions(-) diff --git a/docs/src/lib/sets/Ball1.md b/docs/src/lib/sets/Ball1.md index eb1fa4d9d7..ef9d389587 100644 --- a/docs/src/lib/sets/Ball1.md +++ b/docs/src/lib/sets/Ball1.md @@ -6,12 +6,13 @@ CurrentModule = LazySets ```@docs Ball1 +center(::Ball1) +radius_ball(::Ball1) +ball_norm(::Ball1) σ(::AbstractVector, ::Ball1) ρ(::AbstractVector, ::Ball1) ∈(::AbstractVector, ::Ball1, ::Bool=false) vertices_list(::Ball1) -center(::Ball1) -radius_ball(::Ball1) rand(::Type{Ball1}) constraints_list(::Ball1) translate(::Ball1, ::AbstractVector) diff --git a/docs/src/lib/sets/Ball2.md b/docs/src/lib/sets/Ball2.md index a4dc6f2c1e..004a004ceb 100644 --- a/docs/src/lib/sets/Ball2.md +++ b/docs/src/lib/sets/Ball2.md @@ -6,11 +6,12 @@ CurrentModule = LazySets ```@docs Ball2 +center(::Ball2) +radius_ball(::Ball2) +ball_norm(::Ball2) ρ(::AbstractVector, ::Ball2) σ(::AbstractVector, ::Ball2) ∈(::AbstractVector, ::Ball2) -center(::Ball2) -radius_ball(::Ball2) rand(::Type{Ball2}) sample(::Ball2{N}, ::Int) where {N} translate(::Ball2, ::AbstractVector) diff --git a/docs/src/lib/sets/BallInf.md b/docs/src/lib/sets/BallInf.md index 5143875ba9..accf06bdf8 100644 --- a/docs/src/lib/sets/BallInf.md +++ b/docs/src/lib/sets/BallInf.md @@ -8,6 +8,7 @@ CurrentModule = LazySets BallInf center(::BallInf) radius_ball(::BallInf) +ball_norm(::BallInf) radius(::BallInf, ::Real=Inf) radius_hyperrectangle(::BallInf) radius_hyperrectangle(::BallInf, ::Int) diff --git a/docs/src/lib/sets/Ballp.md b/docs/src/lib/sets/Ballp.md index 655c295631..ec2b097060 100644 --- a/docs/src/lib/sets/Ballp.md +++ b/docs/src/lib/sets/Ballp.md @@ -6,11 +6,12 @@ CurrentModule = LazySets ```@docs Ballp +center(::Ballp) +radius_ball(::Ballp) +ball_norm(::Ballp) σ(::AbstractVector, ::Ballp) ρ(::AbstractVector, ::Ballp) ∈(::AbstractVector, ::Ballp) -center(::Ballp) -radius_ball(::Ballp) rand(::Type{Ballp}) translate(::Ballp, ::AbstractVector) translate!(::Ballp, ::AbstractVector) diff --git a/src/Sets/Ball1.jl b/src/Sets/Ball1.jl index 34bd72d449..154b1b24c4 100644 --- a/src/Sets/Ball1.jl +++ b/src/Sets/Ball1.jl @@ -91,6 +91,24 @@ function radius_ball(B::Ball1) return B.radius end +""" + ball_norm(B::Ball1) + +Return the characteristic norm of a ball in the 1-norm. + +### Input + +- `B` -- ball in the 1-norm + +### Output + +The characteristic norm, which is `1`. +""" +function ball_norm(B::Ball1) + N = eltype(B) + return one(N) +end + """ vertices_list(B::Ball1) diff --git a/src/Sets/Ball2.jl b/src/Sets/Ball2.jl index 2d39d0e812..a4962a0363 100644 --- a/src/Sets/Ball2.jl +++ b/src/Sets/Ball2.jl @@ -103,6 +103,24 @@ function radius_ball(B::Ball2) return B.radius end +""" + ball_norm(B::Ball2) + +Return the characteristic norm of a ball in the 2-norm. + +### Input + +- `B` -- ball in the 2-norm + +### Output + +The characteristic norm, which is `2`. +""" +function ball_norm(B::Ball2) + N = eltype(B) + return N(2) +end + """ ρ(d::AbstractVector, B::Ball2) diff --git a/src/Sets/BallInf.jl b/src/Sets/BallInf.jl index 0297a9b216..67be1b27d6 100644 --- a/src/Sets/BallInf.jl +++ b/src/Sets/BallInf.jl @@ -173,6 +173,24 @@ function radius_ball(B::BallInf) return B.radius end +""" + ball_norm(B::BallInf) + +Return the characteristic norm of a ball in the infinity norm. + +### Input + +- `B` -- ball in the infinity norm + +### Output + +The characteristic norm, which is `Inf`. +""" +function ball_norm(B::BallInf) + N = eltype(B) + return N(Inf) +end + """ σ(d::AbstractVector, B::BallInf) diff --git a/src/Sets/Ballp.jl b/src/Sets/Ballp.jl index 0e8caff21e..b6cd2b726e 100644 --- a/src/Sets/Ballp.jl +++ b/src/Sets/Ballp.jl @@ -111,6 +111,23 @@ function radius_ball(B::Ballp) return B.radius end +""" + ball_norm(B::Ballp) + +Return the characteristic norm of a ball in the p-norm. + +### Input + +- `B` -- ball in the p-norm + +### Output + +The characteristic norm, which is `p`. +""" +function ball_norm(B::Ballp) + return B.p +end + """ σ(d::AbstractVector, B::Ballp) diff --git a/test/Sets/Ball1.jl b/test/Sets/Ball1.jl index c0f6091d41..1bc876f0f2 100644 --- a/test/Sets/Ball1.jl +++ b/test/Sets/Ball1.jl @@ -62,6 +62,9 @@ for N in [Float64, Rational{Int}, Float32] # radius_ball @test LazySets.radius_ball(b) == N(2) + # ball_norm + @test LazySets.ball_norm(b) == N(1) + # boundedness @test isbounded(b) diff --git a/test/Sets/Ball2.jl b/test/Sets/Ball2.jl index 795ef636f1..3851057ecf 100644 --- a/test/Sets/Ball2.jl +++ b/test/Sets/Ball2.jl @@ -67,6 +67,9 @@ for N in [Float64, Float32] # radius_ball @test LazySets.radius_ball(b) == N(2) + # ball_norm + @test LazySets.ball_norm(b) == N(2) + # boundedness @test isbounded(b) && isboundedtype(typeof(b)) diff --git a/test/Sets/BallInf.jl b/test/Sets/BallInf.jl index 2809bd9cfa..72aa99c697 100644 --- a/test/Sets/BallInf.jl +++ b/test/Sets/BallInf.jl @@ -65,6 +65,9 @@ for N in [Float64, Rational{Int}, Float32] # radius_ball @test LazySets.radius_ball(b) == N(2) + # ball_norm + @test LazySets.ball_norm(b) == N(Inf) + # support vector for single entry vector svec = σ(SingleEntryVector(2, 3, N(2)), BallInf(zeros(N, 3), N(2))) @test svec[1] ∈ Interval(N[-2, 2]) && svec[2] == N(2) && svec[3] ∈ Interval(N[-2, 2]) diff --git a/test/Sets/Ballp.jl b/test/Sets/Ballp.jl index a170c3356d..c4c0bc971f 100644 --- a/test/Sets/Ballp.jl +++ b/test/Sets/Ballp.jl @@ -39,6 +39,9 @@ for N in [Float64, Float32] # radius_ball @test LazySets.radius_ball(b) == N(2) + # ball_norm + @test LazySets.ball_norm(b) == N(3) + # an_element @test an_element(b) isa AbstractVector{N} From 6e6e05723f5cd62f586f19c40862e9a80b771c02 Mon Sep 17 00:00:00 2001 From: schillic Date: Thu, 3 Aug 2023 09:03:12 +0200 Subject: [PATCH 3/6] add AbstractBallp interface --- docs/src/lib/interfaces.md | 15 ++++++++++ src/Interfaces/AbstractBallp.jl | 31 ++++++++++++++++++++ src/Interfaces/AbstractCentrallySymmetric.jl | 5 ++-- src/LazySets.jl | 1 + src/Sets/Ball2.jl | 4 +-- src/Sets/Ballp.jl | 4 +-- 6 files changed, 53 insertions(+), 7 deletions(-) create mode 100644 src/Interfaces/AbstractBallp.jl diff --git a/docs/src/lib/interfaces.md b/docs/src/lib/interfaces.md index b34ab8c478..8707db84e7 100644 --- a/docs/src/lib/interfaces.md +++ b/docs/src/lib/interfaces.md @@ -553,3 +553,18 @@ dim(::AbstractPolynomialZonotope) * [Polynomial zonotope (DensePolynomialZonotope)](@ref def_DensePolynomialZonotope) * [Simplified sparse polynomial zonotope (SimpleSparsePolynomialZonotope)](@ref def_SimpleSparsePolynomialZonotope) * [Sparse polynomial zonotope (SparsePolynomialZonotope)](@ref def_SparsePolynomialZonotope) + +## [Balls in the p-norm (AbstractBallp)](@id def_AbstractBallp) + +A ball is a centrally-symmetric set with a characteristic p-norm. + +```@docs +AbstractBallp +``` + +### Implementations + +* [Ball1](@ref def_Ball1) +* [Ball2](@ref def_Ball2) +* [BallInf](@ref def_BallInf) +* [Ballp](@ref def_Ballp) diff --git a/src/Interfaces/AbstractBallp.jl b/src/Interfaces/AbstractBallp.jl new file mode 100644 index 0000000000..7079c8d6f8 --- /dev/null +++ b/src/Interfaces/AbstractBallp.jl @@ -0,0 +1,31 @@ +export AbstractBallp + +""" + AbstractBallp{N} <: AbstractCentrallySymmetric{N} + +Abstract type for p-norm balls. + +### Notes + +See [`Ballp`](@ref) for a standard implementation of this interface. + +Every concrete `AbstractBallp` must define the following methods: + +- `center(::AbstractBallp)` -- return the center +- `radius_ball(::AbstractBallp)` -- return the ball radius +- `ball_norm(::AbstractBallp)` -- return the characteristic norm + +The subtypes of `AbstractBallp`: + +```jldoctest; setup = :(using LazySets: subtypes) +julia> subtypes(AbstractBallp) +2-element Vector{Any}: + Ball2 + Ballp +``` + +There are two further set types implementing the `AbstractBallp` interface, but +they also implement other interfaces and hence cannot be subtypes: `Ball1` and +`BallInf`. +""" +abstract type AbstractBallp{N} <: AbstractCentrallySymmetric{N} end diff --git a/src/Interfaces/AbstractCentrallySymmetric.jl b/src/Interfaces/AbstractCentrallySymmetric.jl index 56f09a957d..8dd04b6df9 100644 --- a/src/Interfaces/AbstractCentrallySymmetric.jl +++ b/src/Interfaces/AbstractCentrallySymmetric.jl @@ -21,9 +21,8 @@ The subtypes of `AbstractCentrallySymmetric`: ```jldoctest; setup = :(using LazySets: subtypes) julia> subtypes(AbstractCentrallySymmetric) -3-element Vector{Any}: - Ball2 - Ballp +2-element Vector{Any}: + AbstractBallp Ellipsoid ``` """ diff --git a/src/LazySets.jl b/src/LazySets.jl index 158bc5cfad..10bf2eba52 100644 --- a/src/LazySets.jl +++ b/src/LazySets.jl @@ -75,6 +75,7 @@ include("Interfaces/AbstractPolygon.jl") include("Interfaces/AbstractSingleton.jl") include("Interfaces/AbstractHPolygon.jl") include("Interfaces/AbstractAffineMap.jl") +include("Interfaces/AbstractBallp.jl") # ============================= # Types representing basic sets diff --git a/src/Sets/Ball2.jl b/src/Sets/Ball2.jl index a4962a0363..43ce4e7084 100644 --- a/src/Sets/Ball2.jl +++ b/src/Sets/Ball2.jl @@ -6,7 +6,7 @@ export Ball2, volume """ - Ball2{N<:AbstractFloat, VN<:AbstractVector{N}} <: AbstractCentrallySymmetric{N} + Ball2{N<:AbstractFloat, VN<:AbstractVector{N}} <: AbstractBallp{N} Type that represents a ball in the 2-norm. @@ -52,7 +52,7 @@ julia> σ([1.0, 2, 3, 4, 5], B) 0.3370999312316211 ``` """ -struct Ball2{N<:AbstractFloat,VN<:AbstractVector{N}} <: AbstractCentrallySymmetric{N} +struct Ball2{N<:AbstractFloat,VN<:AbstractVector{N}} <: AbstractBallp{N} center::VN radius::N diff --git a/src/Sets/Ballp.jl b/src/Sets/Ballp.jl index b6cd2b726e..7ca9ecd9d4 100644 --- a/src/Sets/Ballp.jl +++ b/src/Sets/Ballp.jl @@ -4,7 +4,7 @@ import Base: rand, export Ballp """ - Ballp{N<:AbstractFloat, VN<:AbstractVector{N}} <: AbstractCentrallySymmetric{N} + Ballp{N<:AbstractFloat, VN<:AbstractVector{N}} <: AbstractBallp{N} Type that represents a ball in the p-norm, for ``1 ≤ p ≤ ∞``. @@ -54,7 +54,7 @@ julia> σ([1.0, 2, 3, 4, 5], B) 0.33790011086518895 ``` """ -struct Ballp{N<:AbstractFloat,VN<:AbstractVector{N}} <: AbstractCentrallySymmetric{N} +struct Ballp{N<:AbstractFloat,VN<:AbstractVector{N}} <: AbstractBallp{N} p::N center::VN radius::N From ffc9690437434c1c86ef19e29970c2c9855454ad Mon Sep 17 00:00:00 2001 From: schillic Date: Thu, 3 Aug 2023 09:03:30 +0200 Subject: [PATCH 4/6] low/high for AbstractBallp --- src/Interfaces/AbstractBallp.jl | 32 ++++++++++++++++++++++++++++++++ src/Sets/Ball1.jl | 16 ++++++++++++++++ src/Sets/BallInf.jl | 16 ++++++++++++++++ test/Sets/Ball1.jl | 4 ++++ test/Sets/Ball2.jl | 4 +++- test/Sets/BallInf.jl | 6 +++--- test/Sets/Ballp.jl | 4 ++++ 7 files changed, 78 insertions(+), 4 deletions(-) diff --git a/src/Interfaces/AbstractBallp.jl b/src/Interfaces/AbstractBallp.jl index 7079c8d6f8..68ccaea804 100644 --- a/src/Interfaces/AbstractBallp.jl +++ b/src/Interfaces/AbstractBallp.jl @@ -29,3 +29,35 @@ they also implement other interfaces and hence cannot be subtypes: `Ball1` and `BallInf`. """ abstract type AbstractBallp{N} <: AbstractCentrallySymmetric{N} end + +function low(B::AbstractBallp) + return _low_AbstractBallp(B) +end + +function _low_AbstractBallp(B::LazySet) + return center(B) .- radius_ball(B) +end + +function low(B::AbstractBallp, i::Int) + return _low_AbstractBallp(B, i) +end + +function _low_AbstractBallp(B::LazySet, i::Int) + return center(B, i) - radius_ball(B) +end + +function high(B::AbstractBallp) + return _high_AbstractBallp(B) +end + +function _high_AbstractBallp(B::LazySet) + return center(B) .+ radius_ball(B) +end + +function high(B::AbstractBallp, i::Int) + return _high_AbstractBallp(B, i) +end + +function _high_AbstractBallp(B::LazySet, i::Int) + return center(B, i) + radius_ball(B) +end diff --git a/src/Sets/Ball1.jl b/src/Sets/Ball1.jl index 154b1b24c4..76a6b00562 100644 --- a/src/Sets/Ball1.jl +++ b/src/Sets/Ball1.jl @@ -109,6 +109,22 @@ function ball_norm(B::Ball1) return one(N) end +function low(B::Ball1) + return _low_AbstractBallp(B) +end + +function low(B::Ball1, i::Int) + return _low_AbstractBallp(B, i) +end + +function high(B::Ball1) + return _high_AbstractBallp(B) +end + +function high(B::Ball1, i::Int) + return _high_AbstractBallp(B, i) +end + """ vertices_list(B::Ball1) diff --git a/src/Sets/BallInf.jl b/src/Sets/BallInf.jl index 67be1b27d6..941cabc608 100644 --- a/src/Sets/BallInf.jl +++ b/src/Sets/BallInf.jl @@ -191,6 +191,22 @@ function ball_norm(B::BallInf) return N(Inf) end +function low(B::BallInf) + return _low_AbstractBallp(B) +end + +function low(B::BallInf, i::Int) + return _low_AbstractBallp(B, i) +end + +function high(B::BallInf) + return _high_AbstractBallp(B) +end + +function high(B::BallInf, i::Int) + return _high_AbstractBallp(B, i) +end + """ σ(d::AbstractVector, B::BallInf) diff --git a/test/Sets/Ball1.jl b/test/Sets/Ball1.jl index 1bc876f0f2..b8d5e40877 100644 --- a/test/Sets/Ball1.jl +++ b/test/Sets/Ball1.jl @@ -65,6 +65,10 @@ for N in [Float64, Rational{Int}, Float32] # ball_norm @test LazySets.ball_norm(b) == N(1) + # low/high/extrema + @test extrema(b) == (low(b), high(b)) == (N[-2, -2], N[2, 2]) + @test extrema(b, 1) == (low(b, 1), high(b, 1)) == (N(-2), N(2)) + # boundedness @test isbounded(b) diff --git a/test/Sets/Ball2.jl b/test/Sets/Ball2.jl index 3851057ecf..d4fc1248be 100644 --- a/test/Sets/Ball2.jl +++ b/test/Sets/Ball2.jl @@ -153,16 +153,18 @@ for N in [Float64, Float32] b4 = Ball2(N[4, 3, 2, 1], N(2)) @test project(b4, [2, 4]) == Ball2(N[3, 1], N(2)) - # low/high + # low/high/extrema B = Ball2(N[1, 2], N(2)) l1 = low(B, 1) l2 = low(B, 2) h1 = high(B, 1) h2 = high(B, 2) + @test extrema(B, 1) == (l1, h1) @test l1 == N(-1) && l2 == N(0) && h1 == N(3) && h2 == N(4) @test box_approximation(B) ≈ Hyperrectangle(; low=[l1, l2], high=[h1, h2]) @test low(B) == [N(-1), N(0)] @test high(B) == [N(3), N(4)] + @test extrema(B) == (low(B), high(B)) # reflect @test reflect(b4) == Ball2(N[-4, -3, -2, -1], N(2)) diff --git a/test/Sets/BallInf.jl b/test/Sets/BallInf.jl index 72aa99c697..45cfe13dac 100644 --- a/test/Sets/BallInf.jl +++ b/test/Sets/BallInf.jl @@ -108,10 +108,10 @@ for N in [Float64, Rational{Int}, Float32] vl = vertices_list(b) @test vl == [center(b)] - # high and low + # low/high/extrema b = BallInf(N[1, 2], N(1)) - @test high(b) == N[2, 3] - @test low(b) == N[0, 1] + @test extrema(b) == (low(b), high(b)) == (N[0, 1], N[2, 3]) + @test extrema(b, 1) == (low(b, 1), high(b, 1)) == (N(0), N(2)) # isflat @test !isflat(BallInf(N[1, 1], N(1))) && isflat(BallInf(N[1, 1], N(0))) diff --git a/test/Sets/Ballp.jl b/test/Sets/Ballp.jl index c4c0bc971f..af54b198c6 100644 --- a/test/Sets/Ballp.jl +++ b/test/Sets/Ballp.jl @@ -42,6 +42,10 @@ for N in [Float64, Float32] # ball_norm @test LazySets.ball_norm(b) == N(3) + # low/high/extrema + @test extrema(b) == (low(b), high(b)) == (N[-2, -2], N[2, 2]) + @test extrema(b, 1) == (low(b, 1), high(b, 1)) == (N(-2), N(2)) + # an_element @test an_element(b) isa AbstractVector{N} From 0f4f8cdc1181f6c6982e6b5857ac3e23ef20f87e Mon Sep 17 00:00:00 2001 From: schillic Date: Thu, 3 Aug 2023 09:03:33 +0200 Subject: [PATCH 5/6] minkowski_sum of AbstractBallp --- src/ConcreteOperations/minkowski_sum.jl | 33 ++++++++++++++++++++++++ test/ConcreteOperations/minkowski_sum.jl | 16 ++++++++++++ 2 files changed, 49 insertions(+) diff --git a/src/ConcreteOperations/minkowski_sum.jl b/src/ConcreteOperations/minkowski_sum.jl index 833b733a80..9b65a816c6 100644 --- a/src/ConcreteOperations/minkowski_sum.jl +++ b/src/ConcreteOperations/minkowski_sum.jl @@ -507,3 +507,36 @@ function minkowski_sum(P1::SparsePolynomialZonotope, E = cat(expmat(P1), expmat(P2); dims=(1, 2)) return SparsePolynomialZonotope(c, G, GI, E) end + +# Given two balls of the same p-norm, their Minkowski sum is again a p-norm ball, +# which can be seen as follows: +# +# σ_Bp(d) = c + r \\frac{d}{‖d‖_p} +# +# ρ_Bp(d) +# = ⟨c + r \\frac{d}{‖d‖_p}, d⟩ +# = ⟨c, d⟩ + \\frac{r}{‖d‖_p} * ⟨d, d⟩ +# +# ρ_Bp¹⊕Bp²(d) +# = ρ_Bp¹(d) + ρ_Bp²(d) +# = ⟨c¹, d⟩ + \\frac{r¹}{‖d‖_p} * ⟨d, d⟩ + ⟨c², d⟩ + \\frac{r²}{‖d‖_p} * ⟨d, d⟩ +# = ⟨c¹ + c², d⟩ + \\frac{r¹ + r²}{‖d‖_p} * ⟨d, d⟩ +# = ρ_Bp³(d) +# +# where Bp³ = Ball(center(Bp¹) + center(Bp²), radius(Bp¹) + radius(Bp²)) +function minkowski_sum(B1::BT1, B2::BT2) where {BT<:AbstractBallp,BT1<:BT,BT2<:BT} + p = ball_norm(B1) + if ball_norm(B2) != p + throw(ArgumentError("this method only applies to balls of the same norm")) + end + return Ballp(p, center(B1) + center(B2), radius_ball(B1) + radius_ball(B2)) +end + +for B in (:Ball1, :BallInf) + @eval begin + # see AbstractBallp method + function minkowski_sum(B1::$B, B2::$B) + return $B(center(B1) + center(B2), radius_ball(B1) + radius_ball(B2)) + end + end +end diff --git a/test/ConcreteOperations/minkowski_sum.jl b/test/ConcreteOperations/minkowski_sum.jl index 8c1b5bb420..d83874b5e8 100644 --- a/test/ConcreteOperations/minkowski_sum.jl +++ b/test/ConcreteOperations/minkowski_sum.jl @@ -5,4 +5,20 @@ for N in [Float64, Float32, Rational{Int}] minkowski_sum(X, Y) minkowski_sum(X, Y; algorithm=Polyhedra.FourierMotzkin()) end + + for B in (Ball1, BallInf) + B1 = B(N[1, 2], N(3)) + B2 = B(N[4, 5], N(6)) + @test minkowski_sum(B1, B2) == B(N[5, 7], N(9)) + end +end + +for N in [Float64, Float32] + B1 = Ball2(N[1, 2], N(3)) + B2 = Ball2(N[4, 5], N(6)) + @test minkowski_sum(B1, B2) == Ball2(N[5, 7], N(9)) + + B1 = Ballp(N(3), N[1, 2], N(3)) + B2 = Ballp(N(3), N[4, 5], N(6)) + @test minkowski_sum(B1, B2) == Ballp(N(3), N[5, 7], N(9)) end From 883829a185bd2a6b1ec120bc854ac84652d713c2 Mon Sep 17 00:00:00 2001 From: schillic Date: Thu, 3 Aug 2023 13:22:07 +0200 Subject: [PATCH 6/6] generalize code from Ballp to AbstractBallp --- docs/src/lib/interfaces.md | 8 ++ docs/src/lib/sets/Ballp.md | 8 +- src/Interfaces/AbstractBallp.jl | 134 ++++++++++++++++++++++++++++++++ src/Sets/Ballp.jl | 131 ------------------------------- 4 files changed, 147 insertions(+), 134 deletions(-) diff --git a/docs/src/lib/interfaces.md b/docs/src/lib/interfaces.md index 8707db84e7..762624da91 100644 --- a/docs/src/lib/interfaces.md +++ b/docs/src/lib/interfaces.md @@ -562,6 +562,14 @@ A ball is a centrally-symmetric set with a characteristic p-norm. AbstractBallp ``` +This interface defines the following functions: + +```@docs +σ(::AbstractVector, ::AbstractBallp) +ρ(::AbstractVector, ::AbstractBallp) +∈(::AbstractVector, ::AbstractBallp) +``` + ### Implementations * [Ball1](@ref def_Ball1) diff --git a/docs/src/lib/sets/Ballp.md b/docs/src/lib/sets/Ballp.md index ec2b097060..4ea6d218ea 100644 --- a/docs/src/lib/sets/Ballp.md +++ b/docs/src/lib/sets/Ballp.md @@ -9,9 +9,6 @@ Ballp center(::Ballp) radius_ball(::Ballp) ball_norm(::Ballp) -σ(::AbstractVector, ::Ballp) -ρ(::AbstractVector, ::Ballp) -∈(::AbstractVector, ::Ballp) rand(::Type{Ballp}) translate(::Ballp, ::AbstractVector) translate!(::Ballp, ::AbstractVector) @@ -33,3 +30,8 @@ Inherited from [`AbstractCentrallySymmetric`](@ref): * [`an_element`](@ref an_element(::AbstractCentrallySymmetric)) * [`extrema`](@ref extrema(::AbstractCentrallySymmetric)) * [`extrema`](@ref extrema(::AbstractCentrallySymmetric, ::Int)) + +Inherited from [`AbstractBallp`](@ref): +* [`σ`](@ref σ(::AbstractVector, ::AbstractBallp)) +* [`ρ`](@ref ρ(::AbstractVector, ::AbstractBallp)) +* [`∈`](@ref ∈(::AbstractVector, ::AbstractBallp)) diff --git a/src/Interfaces/AbstractBallp.jl b/src/Interfaces/AbstractBallp.jl index 68ccaea804..0b747976f0 100644 --- a/src/Interfaces/AbstractBallp.jl +++ b/src/Interfaces/AbstractBallp.jl @@ -61,3 +61,137 @@ end function _high_AbstractBallp(B::LazySet, i::Int) return center(B, i) + radius_ball(B) end + + +""" + σ(d::AbstractVector, B::AbstractBallp) + +Return the support vector of a ball in the p-norm in a given direction. + +### Input + +- `d` -- direction +- `B` -- ball in the p-norm + +### Output + +The support vector in the given direction. +If the direction has norm zero, the center of the ball is returned. + +### Algorithm + +The support vector of the unit ball in the ``p``-norm along direction ``d`` is: + +```math +σ(d, \\mathcal{B}_p^n(0, 1)) = \\dfrac{\\tilde{v}}{‖\\tilde{v}‖_q}, +``` +where ``\\tilde{v}_i = \\frac{|d_i|^q}{d_i}`` if ``d_i ≠ 0`` and +``\\tilde{v}_i = 0`` otherwise, for all ``i=1,…,n``, and ``q`` is the conjugate +number of ``p``. +By the affine transformation ``x = r\\tilde{x} + c``, one obtains that +the support vector of ``\\mathcal{B}_p^n(c, r)`` is + +```math +σ(d, \\mathcal{B}_p^n(c, r)) = \\dfrac{v}{‖v‖_q}, +``` +where ``v_i = c_i + r\\frac{|d_i|^q}{d_i}`` if ``d_i ≠ 0`` and ``v_i = 0`` +otherwise, for all ``i = 1, …, n``. +""" +function σ(d::AbstractVector, B::AbstractBallp) + p = ball_norm(B) + q = p / (p - 1) + v = similar(d) + N = promote_type(eltype(d), eltype(B)) + @inbounds for (i, di) in enumerate(d) + v[i] = di == zero(N) ? di : abs.(di) .^ q / di + end + vnorm = norm(v, p) + if isapproxzero(vnorm) + svec = copy(center(B)) + else + svec = center(B) .+ radius_ball(B) .* (v ./ vnorm) + end + return svec +end + +""" + ρ(d::AbstractVector, B::AbstractBallp) + +Evaluate the support function of a ball in the p-norm in the given direction. + +### Input + +- `d` -- direction +- `B` -- ball in the p-norm + +### Output + +Evaluation of the support function in the given direction. + +### Algorithm + +Let ``c`` and ``r`` be the center and radius of the ball ``B`` in the p-norm, +respectively, and let ``q = \\frac{p}{p-1}``. Then: + +```math +ρ(d, B) = ⟨d, c⟩ + r ‖d‖_q. +``` +""" +function ρ(d::AbstractVector, B::AbstractBallp) + p = ball_norm(B) + q = p / (p - 1) + return dot(d, center(B)) + radius_ball(B) * norm(d, q) +end + +""" + ∈(x::AbstractVector, B::AbstractBallp) + +Check whether a given point is contained in a ball in the p-norm. + +### Input + +- `x` -- point/vector +- `B` -- ball in the p-norm + +### Output + +`true` iff ``x ∈ B``. + +### Notes + +This implementation is worst-case optimized, i.e., it is optimistic and first +computes (see below) the whole sum before comparing to the radius. +In applications where the point is typically far away from the ball, a fail-fast +implementation with interleaved comparisons could be more efficient. + +### Algorithm + +Let ``B`` be an ``n``-dimensional ball in the p-norm with radius ``r`` and let +``c_i`` and ``x_i`` be the ball's center and the vector ``x`` in dimension +``i``, respectively. +Then ``x ∈ B`` iff ``\\left( ∑_{i=1}^n |c_i - x_i|^p \\right)^{1/p} ≤ r``. + +### Examples + +```jldoctest +julia> B = Ballp(1.5, [1.0, 1.0], 1.) +Ballp{Float64, Vector{Float64}}(1.5, [1.0, 1.0], 1.0) + +julia> [0.5, -0.5] ∈ B +false + +julia> [0.5, 1.5] ∈ B +true +``` +""" +function ∈(x::AbstractVector, B::AbstractBallp) + @assert length(x) == dim(B) "a vector of length $(length(x)) is " * + "incompatible with a set of dimension $(dim(B))" + N = promote_type(eltype(x), eltype(B)) + p = ball_norm(B) + sum = zero(N) + @inbounds for i in eachindex(x) + sum += abs(center(B, i) - x[i])^p + end + return _leq(sum^(one(N) / p), radius_ball(B)) +end diff --git a/src/Sets/Ballp.jl b/src/Sets/Ballp.jl index 7ca9ecd9d4..5032445975 100644 --- a/src/Sets/Ballp.jl +++ b/src/Sets/Ballp.jl @@ -128,137 +128,6 @@ function ball_norm(B::Ballp) return B.p end -""" - σ(d::AbstractVector, B::Ballp) - -Return the support vector of a ball in the p-norm in a given direction. - -### Input - -- `d` -- direction -- `B` -- ball in the p-norm - -### Output - -The support vector in the given direction. -If the direction has norm zero, the center of the ball is returned. - -### Algorithm - -The support vector of the unit ball in the ``p``-norm along direction ``d`` is: - -```math -σ(d, \\mathcal{B}_p^n(0, 1)) = \\dfrac{\\tilde{v}}{‖\\tilde{v}‖_q}, -``` -where ``\\tilde{v}_i = \\frac{|d_i|^q}{d_i}`` if ``d_i ≠ 0`` and -``\\tilde{v}_i = 0`` otherwise, for all ``i=1,…,n``, and ``q`` is the conjugate -number of ``p``. -By the affine transformation ``x = r\\tilde{x} + c``, one obtains that -the support vector of ``\\mathcal{B}_p^n(c, r)`` is - -```math -σ(d, \\mathcal{B}_p^n(c, r)) = \\dfrac{v}{‖v‖_q}, -``` -where ``v_i = c_i + r\\frac{|d_i|^q}{d_i}`` if ``d_i ≠ 0`` and ``v_i = 0`` -otherwise, for all ``i = 1, …, n``. -""" -function σ(d::AbstractVector, B::Ballp) - p = B.p - q = p / (p - 1) - v = similar(d) - N = promote_type(eltype(d), eltype(B)) - @inbounds for (i, di) in enumerate(d) - v[i] = di == zero(N) ? di : abs.(di) .^ q / di - end - vnorm = norm(v, p) - if isapproxzero(vnorm) - svec = B.center - else - svec = @.(B.center + B.radius * (v / vnorm)) - end - return svec -end - -""" - ρ(d::AbstractVector, B::Ballp) - -Evaluate the support function of a ball in the p-norm in the given direction. - -### Input - -- `d` -- direction -- `B` -- ball in the p-norm - -### Output - -Evaluation of the support function in the given direction. - -### Algorithm - -Let ``c`` and ``r`` be the center and radius of the ball ``B`` in the p-norm, -respectively, and let ``q = \\frac{p}{p-1}``. Then: - -```math -ρ(d, B) = ⟨d, c⟩ + r ‖d‖_q. -``` -""" -function ρ(d::AbstractVector, B::Ballp) - q = B.p / (B.p - 1) - return dot(d, B.center) + B.radius * norm(d, q) -end - -""" - ∈(x::AbstractVector, B::Ballp) - -Check whether a given point is contained in a ball in the p-norm. - -### Input - -- `x` -- point/vector -- `B` -- ball in the p-norm - -### Output - -`true` iff ``x ∈ B``. - -### Notes - -This implementation is worst-case optimized, i.e., it is optimistic and first -computes (see below) the whole sum before comparing to the radius. -In applications where the point is typically far away from the ball, a fail-fast -implementation with interleaved comparisons could be more efficient. - -### Algorithm - -Let ``B`` be an ``n``-dimensional ball in the p-norm with radius ``r`` and let -``c_i`` and ``x_i`` be the ball's center and the vector ``x`` in dimension -``i``, respectively. -Then ``x ∈ B`` iff ``\\left( ∑_{i=1}^n |c_i - x_i|^p \\right)^{1/p} ≤ r``. - -### Examples - -```jldoctest -julia> B = Ballp(1.5, [1.0, 1.0], 1.) -Ballp{Float64, Vector{Float64}}(1.5, [1.0, 1.0], 1.0) - -julia> [0.5, -0.5] ∈ B -false - -julia> [0.5, 1.5] ∈ B -true -``` -""" -function ∈(x::AbstractVector, B::Ballp) - @assert length(x) == dim(B) "a vector of length $(length(x)) is " * - "incompatible with a set of dimension $(dim(B))" - N = promote_type(eltype(x), eltype(B)) - sum = zero(N) - for i in eachindex(x) - sum += abs(B.center[i] - x[i])^B.p - end - return _leq(sum^(one(N) / B.p), B.radius) -end - """ rand(::Type{Ballp}; [N]::Type{<:Real}=Float64, [dim]::Int=2, [rng]::AbstractRNG=GLOBAL_RNG, [seed]::Union{Int, Nothing}=nothing)