Decomposing an Affine Map
In this section we present an illustrative example of the decomposed image of a linear map.
Preliminaries: Polygon, Linear Map, and Plotting
Consider the matrix-valued function $Φ(θ) := \begin{pmatrix} \cos (θ) & -\sin (θ) \\ \sin (θ) & \cos (θ) \end{pmatrix}$, $θ ∈ [π, π]$.
using LazySets, LazySets.Approximations, Plots
Φ(theta) = [cos(theta) -sin(theta); sin(theta) cos(theta)]
Φ (generic function with 1 method)
Now define an arbitrary convex polygon with five vertices with $\operatorname{CH}$ denoting the convex hull operation,
\[\mathcal{X} := \operatorname{CH}\big( (1, 0.5), (1.1, 0.2), (1.4, 0.3), (1.7, 0.5), (1.4, 0.8) \big)\]
This set can be defined as:
X = VPolygon([[1.0, 0.5], [1.1, 0.2], [1.4, 0.3], [1.7, 0.5], [1.4, 0.8]])
VPolygon{Float64, Vector{Float64}}([[1.0, 0.5], [1.1, 0.2], [1.4, 0.3], [1.7, 0.5], [1.4, 0.8]])
You can as well define the convex hull of the one element sets (singletons) via
C = CHArray([Singleton([1.0, 0.5]), Singleton([1.1, 0.2]), Singleton([1.4, 0.3]), Singleton([1.7, 0.5]), Singleton([1.4, 0.8])])
Observe that C
is just a lazy convex hull, whereas X
is a polygon in vertex representation.
Applying the linear map $Φ(π/4) · \mathcal{X}$, we get a new polygon $\mathcal{X}'$ which is the counter-clockwise turn of $\mathcal{X}$ by $θ \triangleq 45°$. In this package the linear map is not computed explicitly but only wrapped in a LinearMap
instance.
Xp = Φ(pi/4) * X
typeof(Xp)
LinearMap{Float64, VPolygon{Float64, Vector{Float64}}, Float64, Matrix{Float64}}
Let us plot the two polygons, $\mathcal{X}$ in green and $\mathcal{X}'$ in blue.
example = plot(X, color="green")
plot!(example, Xp, 1e-2, color="blue")
Note that we have passed 1e-2
as additional argument for the LinearMap
set ($\mathcal{X}'$) because by default such a set is just plotted as its box (or hyperrectangle) approximation. The value 1e-2
is the precision up to which the set is (over-)approximated with a polygon, which in this case is sufficient to obtain the actual set again. If such value is not specified, it defaults to 1e-3
.
Cartesian Decomposition
Next we want to decompose $\mathcal{X}$ into a Cartesian product of intervals. That is, we project it to the $x$-axis and $y$-axis and then compose these intervals again: $\hat{\mathcal{X}} := \hat{\mathcal{X}}_1 × \hat{\mathcal{X}}_2$.
Xhat = overapproximate(X) # approximation of X with an axis-aligned polygon
plot!(example, Xhat, color="gray", alpha=0.3)
Decomposed Image of a Linear Map
Now let us compute the linear map for the box approximation, and let us call it $\mathcal{Y} := Φ(π/4) · \hat{\mathcal{X}}$. This will be a diamond-like shape (the box turned by 45°
).
Y = Φ(pi/4) * Xhat
plot!(example, Y, color="yellow", alpha=0.3)
However, we want our approximation be again a Cartesian product of intervals, so we have to overapproximate this diamond-like shape again: $\hat{\mathcal{Y}} := \hat{\mathcal{X}}' = \hat{\mathcal{X}}_1' × \hat{\mathcal{X}}_2'$
Xhatp = overapproximate(Y)
plot!(example, Xhatp, color="gray", alpha=0.3)
As we can see, the resulting box $\hat{\mathcal{X}}'$ is not a tight overapproximation of $\mathcal{X}'$. We can, however, gain precision by reducing the angle by which we turn the set, e.g., making two smaller turns. Why not try it out?