Plot feasible region of a high-dimensional linear programming along some dimensions

8 visualizaciones (últimos 30 días)
Let be a vector of unknowns. Consider the linear system
where are matrices and vectors of known parameters with appropriate dimensions. They are attached in .mat format in the last comment to the question. Let . is bounded.
My objective: I would like to plot the region of values of for which there exists such that satisfies the linear system above. Let us call such a region .
I initially thought about the following strategy to reach my objective.
1) Find the extreme points of .
2) Transform the extreme points of using E to get a description of the extreme points of .
3) Use fill to colour the area traced by the extreme point of .
I'm stuck at 1). I tried to use this without success, by running
[V,nr,nre]=lcon2vert(C,d,A,b,[])
which stops with a long error message that I'm struggling to intepret. On top of that, even if I was able to make the function lcon2vert running, I doubt that 32 unknowns is a manageable size.
I would like to know your opinion on this. Could you think of other strategies? For example an algorithm that find random points in across each iteration may be an idea: as the number of iteration goes to infinity, the scatter of the points gets closer and closer to .
Clarification: The linear system above has at least one solution (see possible_solution_complete.mat)
  7 comentarios
Matt J
Matt J el 14 de Ag. de 2019
For the convenience of the community, I have reposted CT's data in a single .mat file

Iniciar sesión para comentar.

Respuesta aceptada

Bruno Luong
Bruno Luong el 14 de Ag. de 2019
Editada: Bruno Luong el 14 de Ag. de 2019
Something is wrong in my implementation of understanding of your problem, but there is no feasible point found by my code.
A=load('A.txt');
b=load('b.txt');
C=load('C.txt');
d=load('d.txt');
E=load('E.txt');
Aeq = [A, zeros(size(A,1),2); ...
E, -eye(2)
];
beq = [b; ...
zeros(size(E,1),1)
];
N = 360;
Rb = nan(2,N);
for i=0:N-1
theta = 2*pi/N*i;
fi = [zeros(32,1); cos(theta); sin(theta)];
[xi,~, exitflag] = linprog(fi, ...
[C, zeros(size(C,1),2)], d, ...
Aeq, beq, [], []);
if exitflag>=0
Rb(:,i+1) = xi([33 34]);
end
end
Edit: this code no longer applicable, see updated code using mat files below
  18 comentarios
CT
CT el 15 de Ag. de 2019
With 3 dimensions and the matrices attached, the feasible region should not necessarily be in the 3-D unit simplex. This is because the matrices attached have been designed for the 2-D case.
Matt J
Matt J el 16 de Ag. de 2019
It indicates the boundary/domain is a line. But again there might be just an artifact of discretization (or not). Who knows. I increases N to 3600, I still get a line.
This can also be seen by eliminating the equality constraints A*y=b. The solutions have the form y=Na*c+y0 where Na=null(A) and y0 is the given feasible point. The solutions z=[x33,x34] then have the form z=E*(Na*c+y0). However,
>> Na=null(A); rank(E*Na,1e-12)
ans =
1
So, z lies in some 1-dimensional space.

Iniciar sesión para comentar.

Más respuestas (2)

Bruno Luong
Bruno Luong el 14 de Ag. de 2019
Editada: Bruno Luong el 14 de Ag. de 2019
What is your motivation to do that?
You clearly underestimate such task, simply because the number of vertices of simplex grow very fast with the dimension.
I see your inequality constrainst are 72 (C matrix) assuming they are independent and number of variable is 34. Then the number of vertices can go as big as
>> nchoosek(72,34)
Warning: Result may not be exact. Coefficient is greater than 9.007199e+15 and is only accurate to 15 digits
> In nchoosek (line 92)
ans =
3.9656e+20
If you project on subspace defined by A ane E (20+2 rows) then you get a subspace of dimension 12, then the number of vertices is still very large
>> nchoosek(72,34-22)
ans =
1.5363e+13
To store all of them you need
nchoosek(72,34-22)*8/1e12
~123 Tb.
Each of the vertices require to solve a linear system of 34 x 34.
There is very little chance that you could find a method to do it in reasonably time and with enough accuracy.
Generate a random points to cover such shape? Good luck!!!
  2 comentarios
CT
CT el 14 de Ag. de 2019
Thanks. I suspected that, but wanted to double check with experts if there is a simpler way to obtain and plot , which is only two-dimensional. At the end I just want . Is obtaining necessary to get ?
Bruno Luong
Bruno Luong el 14 de Ag. de 2019
Editada: Bruno Luong el 14 de Ag. de 2019
R is projection of P in 2-dimensional subspace, so I think there is no short-cut,
Now if you want an approximation of R you can discretize:
f_i := [0,0,...,0, cos(theta_i), sin(theta_i)].'
where
theta_i is 2*pi*i/N; i=0,1,2....,N-1 and N large enough.
Finding the N solutions of N LP
x_i = argmin f_i'*x such that
equality/inequalities constraints with matrices A, C, E are meet.
After taking coordinated 33/34 of {x_i} you will getsome discretization of the boundary of R, but there is no warranty what so ever that descretization is the true shape.

Iniciar sesión para comentar.


Matt J
Matt J el 17 de Ag. de 2019
I tried to use this without success ... which stops with a long error message that I'm struggling to intepret.
The error message indicates that the region described by C*y<=d is not solid, i.e., it lies in some strict subspace of R^32, to within numerical precision. This is a numerically unstable way of specifying your constraints. It means that some of the inequalities in the system C*x<=d should really be expressed as equalities.
If you think this assessment is inaccurate, try to present an interior point y, one that strictly satisifes the inequalities C*y<d.

Categorías

Más información sobre Solver Outputs and Iterative Display en Help Center y File Exchange.

Etiquetas

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by