improve numerical accuracy of polytope.reduce (polytope.intersect)

[johnyf]

polytope.intersect calls polytope.reduce to simplify the H-representation after intersection, which basically concatenates constraints. When the numbers in given matrices have different orders, reduce does not work well. Below is an example where a1 intersection a2 and a2 intersection a1 lead to two different results (with the cvxopt conelp solver). a1 and a2 are two continuous propositions for a problem where the proposition preserving partition and the abstraction results tulip returns do not make sense. Though root cause seems to be reduce.

A short-term solution would be to include some warnings in docstrings/tutorial so that users can transform their system to have continuous state variables of similar magnitude. A mid-term solution might be to introduce some conditioning within the polytope package based on the size of input matrices. Long-term solution might be to interface with better LP solvers, which should partially alleviate this issue.

a1 = pc.Polytope.from_box([[0., 30000.],[0., 25.]])

a2 = pc.Polytope.from_box([[0., 40000.],[0., 25.]])

pc.reduce(pc.Polytope(vstack((a2.A,a1.A)),vstack((a2.b,a1.b))))
Out[86]: 
Single polytope 
[[ 1.  0.] |    [[ 30000.]
 [ 0.  1.] x <=  [    25.]
 [-1. -0.] |     [     0.]
 [-0. -1.]]|     [     0.]]

pc.reduce(pc.Polytope(vstack((a1.A,a2.A)),vstack((a1.b,a2.b))))
Out[87]: 
Single polytope 
[[ 1.  0.] |    [[ 40000.]
 [ 0.  1.] x <=  [    25.]
 [-1. -0.] |     [     0.]
 [-0. -1.]]|     [     0.]]

[pejaab]

I would greatly appreciate it, if you could elaborate a little more what is meant by "so that users can transform their system to have continuous state variables of similar magnitude", as I have been struggling with exactly this instability for quite a while and just now realized that you already have an open issue... any pointers would be more than welcome.

[johnyf]

I only faintly remember my thoughts at the time. I think that what I had in mind was to start by applying some geometric transformation that stretches / squeezes dimensions to avoid ill-conditioned shapes (as the polytopes in the OP are--slivers).

However, in the presence of dynamics, a transformation that avoids such extreme geometry will also transform the dynamics, and could in principle lead to ill-conditioning of the dynamics. So there is a trade-off. Perhaps geometric preprocessing of this sort would be more applicable to kinematic analyses.

Regarding the phrase "state variables of similar magnitude", I may have been referring to the polytope that defines the entire domain of interest (bounded). If this polytope is a box with some ratio of sides too large or small, then this means that one variable (coordinate) can take much larger values than some other variable.