Convex set: Finding the A and b matrices from Ax=b when the points of the convex set is known.

Background:

I want to create convex sets from points given in a map. E.g. if I have some points distinguishing land and water in a harbor and create a convex set from this. Then I would use the convex set to solve a optimaztion problem.

Specific:

Is there a generic way of finding the A and b matrices from the equation Ax=b where the points of the convex set is known?

Say I have the convex set below with points:

P =

0 0

1.0000 -1.0000

1.5000 -0.5000

1.5000 0.5000

1.0000 1.0000

0 0

How could I extract the matrices A and b such that I could make ineqalities for an opimization problem?

I've already created a convex set with A = [-1,0; 1,0; -1,1; -1, -1/2]; and b = [300,200,300,600]'; which gave me the convex set below:

but is it possible to solve it the other way?

Appriciate any help!

Thanks

Answer 1

[John D'Errico]

@Matt J has a set of tools for working with such objects.

https://www.mathworks.com/matlabcentral/fileexchange/30892-analyze-n-dimensional-convex-polyhedra?s_tid=ta_fx_results

In there, he has the utility vert2lcon, which can convert a set of vertices from the polytope into a set of linear constraints.

Or, you can go further back, and use the @Michael Kleder code, vert2con, which is less sophisticated, but should still work fine on a simple convex hull in 2-d.

https://www.mathworks.com/matlabcentral/fileexchange/7895-vert2con-vertices-to-constraints?s_tid=srchtitle

[Matt J]

Note that Michael Kleder's code gave me problems in some cases, which is partly what motivated me to write my toolset. I can't remember if vert2con was ever inaccurate, but this definitely fails:

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

[Henrik Rohde Nordhus]

Fantastic! Thank you so much for the quick and helpful response!

[Henrik Rohde Nordhus]

Hi, again! I'm writing my code in Python now and having some trouble with converting the vert2lcon function. Do you know if a sufficient code exists for Python?

Thanks!

[Matt J]

I don't know if there's any ready-made Python equivalent to vert2lcon, but it should be possible to make one for someone willing to put in the effort. Everything used internally by vert2lcon are just standard linear algebra tools and convhulln, which does appear to have a scipy equivalent,

https://docs.scipy.org/doc/scipy/reference/generated/scipy.spatial.ConvexHull.html

You also have the option of calling specific Matlab functions from Python,

https://www.mathworks.com/help/matlab/matlab_external/call-matlab-functions-from-python.html

[Henrik Rohde Nordhus]

Thanks! I tried to use a amatlab.engine but got some errors and couldn't find a solution online that had encountered the same as me... I ended up making my own function, however, it doesn't work quite as well and sometime it return numerically unstable matrices. Hopefully it will work soon!

Answer 2

[Matt J]

Make a change of variables in the problem by rewriting x in terms of an unknown vector of coefficients

where are the given vertices of the convex set.

Answer 3

[Matt J]

If you only need to address the 2D case, then you can use this:

function [A,b]=vert2ineq2D(a)
%Finds the expression of a 2D convex polygon as a set of linear inequalities from
%a set of vertices
%
%  [A,b]=vert2ineq2D(a)
%
%in:
%
% a: Nx2 matrix whos rows are polygon vertex coordinates. The rows must
%    must be ordered so that adjacent rows correspond to adjacent vertices
%    (which will trivially be the case for triangles).
%
%out:
%
% [A,b]: Nx1 vector and Nx2 matrix such that A*x<=b if and only if x is a point
%        inside the polygon
%
%
%Example: Detect whether a point is in a triangle
%
%  %%%data
%  Vertices=[0 0; 1 0; 0 1]; %vertices of a triangle
%  p1=[.5;.25]; %a point inside the triangle
%  p2=[.5;-.25];%a point outside the triangle
%
%  [A,b]=vert2ineq2D(Vertices);
%
%  >>all(A*p1<=b)  %Test if p1 is in the triangle. 
% 
%    ans =
% 
%       1
% 
%  >>all(A*p2<=b)   %Test if  p2 in the triangle.
% 
%   ans =
% 
%       0



centroid=mean(a).';
R=[0 1; -1 0];

A=diff([a;a(1,:)])*R;

b=sum(A.*a,2);

ii=(A*centroid>=b);

b(ii)=-b(ii);
A(ii,:)=-A(ii,:);