Code covered by the BSD License  

Highlights from
Representing Polyhedral Convex Hulls by Vertices or (In)Equalities

5.0

5.0 | 9 ratings Rate this file 19 Downloads (last 30 days) File Size: 7.89 KB File ID: #30892

Representing Polyhedral Convex Hulls by Vertices or (In)Equalities

by Matt J

 

30 Mar 2011 (Updated 07 Mar 2012)

Express bounded polyhedron via equalities/inequalities or vertices.

| Watch this File

File Information
Description

This submission contains VERT2LCON and LCON2VERT, which will find the linear constraints defining a bounded polyhedron in R^n, given its vertices, or vice versa. They are extensions of Michael Kleder's VERT2CON and CON2VERT functions that can handle cases where the polyhedron is not solid in R^n, i.e., where the polyhedron is defined by both equality and inequality constraints. Some other refinements to the original VERT2CON and CON2VERT subroutines, meant to improve speed and reliability, have also been made.
 

SYNTAX:
 
   [A,b,Aeq,beq]=VERT2LCON(V,TOL)

   [V,nr,nre]=LCON2VERT(A,b,Aeq,beq,TOL)
 
The rows of the N x n matrix V are a series of N vertices of the polyhedron
in R^n, defined by the linear constraints
   
    A*x <= b
    Aeq*x = beq
 
For LCON2VERT, Aeq=beq=[] are the default inputs, implying no equality constraints. The output "nr" lists non-redundant inequality constraints, and "nre" lists non-redundant equality constraints. Note further that A,b should be used to specify a solid region, while Aeq and beq should be used to further constrain the region to a sub-manifold of R^n, if needed. You should not use the inequality constraints to express a region in R^n that can be described by equality constraints.
 
The optional TOL argument is a tolerance used for rank-estimation purposes and to handle finite precision issues. Default=1e-10.
 
 
EXAMPLE:
 
Consider V=eye(3) corresponding to the 3D region defined by x+y+z=1, x>=0, y>=0, z>=0.
 
 
 >> [A,b,Aeq,beq]=vert2lcon(eye(3))
  
     A =
 
         0.4082 -0.8165 0.4082
         0.4082 0.4082 -0.8165
        -0.8165 0.4082 0.4082
 
 
     b =
 
         0.4082
         0.4082
         0.4082
 
 
     Aeq =
 
         0.5774 0.5774 0.5774
 
 
     beq =
 
         0.5774

And the reverse conversion is,

   >> V=lcon2vert(A,b,Aeq,beq)
  
          V =
  
              1.0000 0.0000 0.0000
             -0.0000 1.0000 0
                   0 0.0000 1.0000

Acknowledgements

The author wishes to acknowledge the following in the creation of this submission:
VERT2CON - vertices to constraints, CON2VERT - constraints to vertices

MATLAB release MATLAB 7.11 (2010b)
Tags for This File  
Everyone's Tags
Tags I've Applied
Add New Tags Please login to tag files.
Comments and Ratings (22)
14 Apr 2011 Sandro  
14 Apr 2011 Sandro

It is definitely a good job!!! I may mess up the rate submission above.

14 Apr 2011 Matt J

Ah, ok.

14 Apr 2011 Sandro

Sorry Matt, I am not sure how to remove the first rating. Two comments: "unique" function may be better to be replaced by "consolidator"; It may be good to normalize A,b,etc before the final output(I am not quite sure about its necessity.)

14 Apr 2011 Matt J

Hi Sandro. Thanks for the feedback. I'll think about your recommendations. Could you explain the advantage of consolidator over unique? I chose a normalization for Aeq and beq that led to nice whole integers where possible. A and b are produced by Michael Kleder's functions which seem to normalize it this way already.

Regarding your rating. It's fine. The FEX only scores according to your most recent rating. This is to allow people to update their ratings if the author makes improvements to the submission.

16 Apr 2011 Sandro

Hi Matt,

my concern of using consolidator is because I think function unique has pretty high requirement of two vectors to be the same. When run vert2lcon on a large data set, i have seen a large A with many rows are the same up to 6-8 digital number. This enlarges the size of A and it may not be necessary for the user.

"A and b are produced by Michael Kleder's functions which seem to normalize it this way already."

I doubt about it.

Best,
Sandro

14 Jul 2011 Timothy Stratton

Nice very powerful!

05 Sep 2011 Harshavardhan Sundar  
05 Sep 2011 Harshavardhan Sundar

Hi Matt. First of all thanks a lot for the wonderful piece of code. The code works fine for most of the cases. However there are a few occasions where I get erroneous vertices. One such example is below. My A matrix is :
A=[1 0 0
0 1 0
0 0 1
-1 0 0
0 -1 0
0 0 -1
-0.387696185212551 -0.295276274520099 0.703237014010810
5.73326513862449 -0.422448628998633 -0.113292931747631
-3.79255233000727 -0.360645958414035 -0.673892379691217
-1.27188203250127 0.544750863219514 -0.791324376598499
-2.30540199426269 -2.50522755054822 -0.291771858420232
0.616578245027929 -0.654271120193835 0.574406987898971
-1.49508806389885 -2.68877782397954 -0.989225453798092
5.34556895341194 -0.717724903518731 0.589944082263178
-3.40485614479472 -0.0653696838939366 -1.37712939370203
-0.884185847288714 0.840027137739613 -1.49456139060931
-1.91770580905014 -2.20995127602812 -0.995008872431042
0.228882059815378 -0.949547394713933 1.27764400190978
-1.10739187868630 -2.39350154945944 -1.69246246780890
1.94071280861722 -0.783094587412668 -0.787185311438849
4.46138310612323 0.122302234220881 -0.904617308346131
3.42786314436180 -2.92767617954685 -0.405064790167863
5.11668689359657 0.231822491195202 -0.687699919646602
4.23817707472564 -3.11122645297817 -1.10251838554572
-2.52067029750601 -0.905396821633549 0.117431996907282
1.48715033574458 -2.14458159213419 0.382120521270985
-3.17597408497934 -1.01491707860787 -0.0994853917922465
2.29746426610842 -2.32813186556550 -0.315333074106875
-1.03351996176143 -3.04997841376774 0.499552518178267
-0.655303787473336 -0.109520256974320 -0.216917388699529
-0.223206031397585 -3.23352868719905 -0.197901077199593
-1.68882374923476 -3.15949867074206 0.282635129478739
-0.810313930363841 0.183550273431316 0.697453595377860
-0.878509818870922 -3.34304894417337 -0.414818465899121]
b=[6
4
2
0
0
0
-0.154147860211922
16.6671436561829
-8.83179250583593
-1.48693603688105
-8.92116339978305
0.423928365940654
-8.95704617639833
16.5129957959710
-8.67764464562401
-1.33278817666913
-8.76701553957113
0.269780505728732
-8.80289831618640
7.83535115034699
15.1802076193019
7.74598025639987
16.2432152902423
7.71009747978459
-7.34485646895489
-0.0893708939471196
-8.40786413989528
-0.125253670562397
-7.43422736290201
-1.06300767094039
-7.47011013951728
-8.49723503384240
0.0358827766152754
-8.53311781045768]

I used the plotregion command available in MATLAB central to plot the region. I have reasons to believe that the plotregion command is plotting the region accurately. However for this particular example the vertices I get (graphically) from the plotregion command are not coinciding with the vertices obtained from lcon2vert. I also tried con2vert by Michael Kleder with little success.

06 Sep 2011 Matt J

@Harshavardhan - one additional note. In the meantime, you should also modify this line of code

 if ef ~= 1

to this

  if ~all(A*c < b)

This will at least warn you properly that the function is about to fail.

06 Sep 2011 Matt J

@Harshavardhan - This is a resend of a previous comment which seems to have gotten lost by the FEX.

The data set you've provided is quite intriguing! I've traced the problem you're seeing to Michael Kleder's CON2VERT code (my lcon2vert is mainly a wrapper for this).

The problem is that, in order to work, CON2VERT needs to find an initial point c inside the polytope. The simple method that it uses seems to have had a good rate of success historical success, but lo and behold, your particular A,b data set has defeated it!

I know this is the problem however, because

(i) I've checked that the initial point c that it generates does not satisfy A*c<b

(ii) When I manually change CON2VERT, hardcoding the initial c to a point that I know is in your polytope, I get the correct vertices.

I've alerted Michael Kleder to this via the CON2VERT webpage, but he doesn't appear to have been monitoring the page for several years. I doubt he'll be responding with any input.

I have some ideas of my own to improve the initialization scheme, but they will take some time to develop.

In the meantime, if you're in a situation where you already know an initial point, c, inside your polytope, a workaround for you might be to modify CON2VERT so that it accepts this c is an additional input argument.

Beyond that, I don't know what to say. You have difficult data...

08 Sep 2011 Matt J

@Harshavardhan - I've just updated the package with what I believe is a more robust version of LCON2VERT. At the very least, it successfully handles your example data set. Also, it will more reliably report cases where it cannot find a solution.

18 Sep 2011 Harshavardhan Sundar

Hi Matt, Thanks again (Officially)....the latest of your code lcon2vert and vert2lcon are working absolutely fine....

22 Sep 2011 Sandro

hi Matt,

Do you have any exposure of how high the dim could be? I am testing something like x1+x2+...+xn=1, xi>=0, i=1,2,...n. I randomly generate some points within this region and also give the vertices like (0,0,0..,0), (1,0,0,..0), (0,1,0,...0)...to help checking.

It looks like n=10 works, but n=20 is stuck due to qullmx function called by convexhulln.

Is there any possible we can have function like qullmx but more powerful in the sense of dim and speed?

Sandro

23 Sep 2011 Matt J

Hi Sandro,

I'm not having any such difficult in R2011a, so perhaps you need to upgrade. Here is the test that I ran, and no errors resulted

>> [a,b,ae,be]=vert2lcon(eye(20));
>> V=lcon2vert(a,b,ae,be);
>> norm(V'*V-eye(20))

ans =

  5.7516e-015

31 Oct 2011 david

Hi Matt,
I'll be very thankful if you helps me with this problem: When I applied your function i get the following error :
  "??? Maximum recursion limit of 500 reached. Use set(0,'RecursionLimit',N)
to change the limit. Be aware that exceeding your available stack space can
crash MATLAB and/or your computer."

What i have to do in this case ?

31 Oct 2011 Matt J

@david: You need to tell me how I can reproduce the error.

13 Feb 2012 Deepanshu Vasal

Hi Matt,

I am trying to use your code with my dataset, but I am getting singular matrix warning with the vertices generated are out of bounds (Nan, inf).
My data set is as follows :-

A = - [
    
    0 1 0 -2 0 1 0 0 ;
    0 0 0 0 0 1 -2 1 ;
    0 0 1 -2 1 0 0 0 ;
    
    -1 1 1 -1 0 0 0 0 ;
    0 0 0 -1 1 1 -1 0 ;
    
    
    0 0 0 1 -1 0 -1 1 ;
    1 -1 0 -1 1 0 0 0 ;
    0 0 1 -1 0 -1 1 0 ;
    eye(8);
    -eye(8)];

b = - [zeros(8,1); -10*ones(16,1)];

is there something that cannot be handled by this code? I see that the initial point taken is fine (satisfies the constraints.) but cannot figure out where is the problem in the code.
If you come across any problem, please let me know. Also do you know the reference for the algorithm used in this code?

Thanks so much.
Deepanshu
Thanks so much
Deepanshu

14 Feb 2012 Matt J

@Deepanshu - Are you sure the first 8 constraints are not meant to be equality constraints?

In any case, I answered you in an earlier email, but here's the recap: the region described by your inequality constraints appears to be non-solid. Such inequality constraints are numerically-unstable, which is why all MATLAB Optimization Toolbox functions ask you to represent your linear constraints in the form:

Aineq*x<=bineq
Aeq*x=beq

where Aineq*x=bineq describes a SOLID region.

LCON2VERT works similarly. You must rewrite your non-solid region as above and invoke lcon2vert in 4-argument form:

lcon2vert(Aineq,bineq,Aeq,beq)

15 Feb 2012 Deepanshu Vasal

Matt, Thanks for your reply. Having equality constraint was not the problem. If I take
A = - [
    0 0 0 0 0 1 -2 1 ;
    0 0 1 -2 1 0 0 0 ;
    eye(8);
    -eye(8)];
and
b = - [zeros(2,1); -10*ones(16,1)];

it does not work but if I reduce redundant dimensions and put

A = - [
     0 0 0 1 -2 1 ;
     1 -2 1 0 0 0 ;
    eye(6);
    -eye(6)];
and
b = - [zeros(2,1); -10*ones(12,1)];

it works.

but I used lrs package at http://cgm.cs.mcgill.ca/~avis/C/lrs.html

and I got vertices for my problem.

Anyway I appreciate your effort and your such timely response. good luck

Deepanshu

16 Feb 2012 Matt J

@Deepanshu - actually the 2nd case doesn't appear to work for me, either. I think I've found the problem, but it will take me a few days until I can fix it.

09 Apr 2012 Deepanshu Vasal  
Please login to add a comment or rating.
Updates
06 Apr 2011

added LCON2VERT which does the inverse of VERT2CON

08 Sep 2011

Improved the reliability of CON2VERT subroutine, both in terms of performance and of error reporting.

14 Sep 2011

*Further robustness of lcon2vert
*Improved weeding of non-unique constraints in vert2lcon output
*Outputs A,Aeq of vert2lcon will now have normalized rows.

07 Mar 2012

Various bug fixes and improvements in robustness. These address failure cases discovered recently by users.

Tag Activity for this File
Tag Applied By Date/Time
vert2con Matt J 30 Mar 2011 15:13:53
polyhedron Matt J 30 Mar 2011 15:13:53
polytope Matt J 31 Mar 2011 01:51:06
vertices Matt J 31 Mar 2011 01:51:06
convex hull Matt J 31 Mar 2011 08:54:31
convhull Matt J 31 Mar 2011 08:54:59
convhulln Matt J 31 Mar 2011 08:54:59
constraint Matt J 31 Mar 2011 09:09:13
constraints Matt J 31 Mar 2011 09:09:13
hull Matt J 31 Mar 2011 09:09:14
vertex Matt J 31 Mar 2011 09:09:14
convex Matt J 31 Mar 2011 09:09:14
polygon Matt J 31 Mar 2011 09:09:14
optimization Matt J 31 Mar 2011 09:09:14
solid Matt J 31 Mar 2011 09:10:01
vert2lcon Matt J 06 Apr 2011 10:47:39
lcon2vert Matt J 06 Apr 2011 10:47:39
convex analysis Matt J 17 Apr 2011 23:17:32

Contact us at files@mathworks.com