tulip-control / polytope

Geometric operations on polytopes of any dimension
https://pypi.org/project/polytope
Other
73 stars 19 forks source link

Polytope.reduce and removal of possibly overlapping polytopes #82

Open victoraalves opened 1 year ago

victoraalves commented 1 year ago

Hi everyone,

I have a quick question: Does the reduce function, responsible for removing redundant inequalities from the H-rep, indirectly removes overlapping polytopes (if there are any?)

The docstring for this function gives some information but I wanted to double check. I am working in a package that uses polytope as a requirement and knowing that reduce would do this for me seems faster than a workaround I currently have implemented.

https://github.com/tulip-control/polytope/blob/931698601fd35f18e12ed5e20e10ee1dfc3721e6/polytope/polytope.py#L1051

"""Remove redundant inequalities from the hyperplane representation.

Uses the algorithm described at [1],
by solving one LP for each facet.

[1] https://www.cs.mcgill.ca/~fukuda/soft/polyfaq/node24.html

Warning:
  - nonEmptyBounded == 0 case is not tested much.

@type poly: L{Polytope} or L{Region}

@return: Reduced L{Polytope} or L{Region} object
"""

Thanks!

victoraalves commented 1 year ago

If this is a port of MPT's v3.0 reduce.m, then it has the same functionality I would assume? Per MPT's documentation:

image

necozay commented 1 year ago

reduce function does not eliminate redundant polytopes in a region. When we designed the package, the intention was that regions are not explicitly defined but only generated by taking unions of polytopes. Hence, there is a flag inside union that will remove redundant polytopes (it also tries to merge two polytopes if their union is convex but I am not sure if this flag is extensively tested).

For instance in the following piece of code reduce(reg) will have redundant polytopes but reg2 won't.

p0 = pc.Polytope(A0,b0)
p1 = pc.Polytope(A0,b0)
lst = []
lst.append(p0)
lst.append(p1)
reg = pc.Region(lst)
print(pc.reduce(reg))
reg2 = pc.union(p0,p1,True)
print(reg2)
victoraalves commented 1 year ago

Thank you for the detailed explanation. I noticed that indeed, union with the flag=True will try to merge polytopes if their union is convex, as well as removing redundant polytopes. However, I have noticed that this can be quite slow if I have a polytope.Region with several polytopes.

One idea that is also slow, but faster than using union with flag=True is to take the Region's bounding box, remove the original region from the bounding box, generating a "hollow" polytopic region, and then removing this hollow region from the original bounding box (this is done removing one polytope from the region at a time). The result of this is a non-overlapping polytopic region. The process animated:

overlapping_polytopes

This elegant "trick" is credited to David Vinson, in his Ph.D. dissertation "A new measure of process operability for the improved steady-state design of chemical processes. Ph.D. Thesis, Lehigh University, USA, 2000."

I have implemented this trick but it's the current bottleneck of my application in terms of computational time. This is why I opened the discussion so I could gage more information about Polytope's algorithms.

Nevertheless, I wanted to thank you for developing Polytope, it has been a very useful package!