linear combination and convex hull on sspz

docs/src/lib/sets/SimpleSparsePolynomialZonotope.md

@@ -16,4 +16,7 @@ order(::SimpleSparsePolynomialZonotope)
 linear_map(::AbstractMatrix, ::SimpleSparsePolynomialZonotope)
 minkowski_sum(::SimpleSparsePolynomialZonotope, ::SimpleSparsePolynomialZonotope)
 cartesian_product(::SimpleSparsePolynomialZonotope, ::SimpleSparsePolynomialZonotope)
+linear_combination(::SimpleSparsePolynomialZonotope, ::SimpleSparsePolynomialZonotope)
+convex_hull(::SimpleSparsePolynomialZonotope, ::SimpleSparsePolynomialZonotope)
+convex_hull(::SimpleSparsePolynomialZonotope)
 ```

src/ConcreteOperations/convex_hull.jl

@@ -595,3 +595,38 @@ end
 convex_hull(X::LazySet, ::EmptySet) = X
 convex_hull(::EmptySet, X::LazySet) = X
 convex_hull(∅::EmptySet, ::EmptySet) = ∅
+
+"""
+    convex_hull(P1::SimpleSparsePolynomialZonotope, P2::SimpleSparsePolynomialZonotope)
+
+Computes the convex hull of the simple sparse polynomial zonotopes `P1` and `P2`.
+
+### Input
+
+- `P1` : simple sparse polynomial zonotopes
+- `P2` : simple sparse polynomial zonotopes
+
+### Output
+
+Tightest convex simple sparse polynomial zonotope containing `P1` and `P2`.
+"""
+function convex_hull(P1::SimpleSparsePolynomialZonotope, P2::SimpleSparsePolynomialZonotope)
+    return linear_combination(linear_combination(P1, P1), linear_combination(P2, P2))
+end
+
+"""
+    convex_hull(P::SimpleSparsePolynomialZonotope)
+
+Computes the convex hull of the simple sparse polynomial zonotope `P`.
+
+### Input
+
+- `P` -- simple sparse polynomial zonotopes
+
+### Output
+
+Tightest convex simple sparse polynomial zonotope containing `P`.
+"""
+function convex_hull(P::SimpleSparsePolynomialZonotope)
+    return linear_combination(P, P)
+end

src/ConcreteOperations/linear_combination.jl

@@ -0,0 +1,41 @@
+export linear_combination
+
+"""
+    linear_combination(P1::SimpleSparsePolynomialZonotope, P2::SimpleSparsePolynomialZonotope)
+
+compute the linear combination of simple sparse polynomial zonotopes `P1` and `P2`.
+
+### Input
+
+- `P1` -- simple sparse polynomial zonotope
+- `P2` -- simple sparse polynomial zonotope
+
+### Output
+
+linear combination of `P1` and `P2`.
+
+### Notes
+
+The linear combination of two sets ``P₁`` and ``P₂`` is defined as
+
+```math
+\\{1/2(1+λ)p₁ + 1/2(1-λ)p₂ | p₁ ∈ P₁, p₂ ∈ P₂, λ ∈ [-1, 1]\\}
+```
+"""
+function linear_combination(P1::SimpleSparsePolynomialZonotope, P2::SimpleSparsePolynomialZonotope)
+    c1, c2 = center(P1), center(P2)
+    G1, G2 = genmat(P1), genmat(P2)
+    E1, E2 = expmat(P1), expmat(P2)
+
+    c = 0.5 * (c1 + c2)
+    G = 0.5 * hcat(c1 - c2, G1, G1, G2, -G2)
+
+    N = promote_type(eltype(E1), eltype(E2))
+    E = hcat(vcat(zeros(N, nparams(P1) + nparams(P2)), one(N)),
+             vcat(E1, zeros(N, nparams(P2) + 1, ngens(P1))),
+             vcat(E2, zeros(N, nparams(P2), ngens(P1)), ones(N, 1, ngens(P1))),
+             vcat(zeros(N, nparams(P1), ngens(P2)), E2, zeros(N, 1, ngens(P2))),
+             vcat(zeros(N, nparams(P1), ngens(P2)), E2, ones(N, 1, ngens(P2))))
+
+    return SimpleSparsePolynomialZonotope(c, G, E)
+end

src/LazySets.jl

@@ -159,6 +159,7 @@ include("ConcreteOperations/intersection.jl")
 include("ConcreteOperations/isdisjoint.jl")
 include("ConcreteOperations/issubset.jl")
 include("ConcreteOperations/isstrictsubset.jl")
+include("ConcreteOperations/linear_combination.jl")
 include("ConcreteOperations/minkowski_difference.jl")
 include("ConcreteOperations/minkowski_sum.jl")
 include("Utils/samples.jl")

test/Sets/SimpleSparsePolynomialZonotope.jl

@@ -32,4 +32,25 @@ for N in [Float64, Float32, Rational{Int}]
     @test genmat(CPS) == [1 2 0 0;2 2 0 0;0 0 1 2;0 0 2 2.]
     @test expmat(CPS) == [1 4 0 0;1 2 0 0;0 0 1 4;0 0 1 2]
 
+    _c = N[2.0, 0.0]
+    _g = N[ 0.0  0.5  1.0  0.5  1.0  0.5  1.0  -0.5  -1.0
+           0.0  1.0  1.0  1.0  1.0  1.0  1.0  -1.0  -1.0]
+
+    _e = [0  1  4  1  4  0  0  0  0
+          0  1  2  1  2  0  0  0  0
+          0  0  0  0  0  1  4  1  4
+          0  0  0  0  0  1  2  1  2
+          1  0  0  1  1  0  0  1  1]
+
+    LCS = linear_combination(S, S)
+
+    @test center(LCS) == _c
+    @test genmat(LCS) == _g
+    @test expmat(LCS) == _e
+
+    CH1 = convex_hull(S)
+
+    @test center(CH1) == _c
+    @test genmat(CH1) == _g
+    @test expmat(CH1) == _e
 end