Solving system of equations with constraints
27 views (last 30 days)
Show older comments
Deepa Maheshvare
on 24 Jan 2019
Edited: Deepa Maheshvare
on 24 Jan 2019
I am solving a system of linear equations with constraints using linprog.
Aeq = [-1 0 1 0 0 0 0;1 -1 0 -1 0 0 0;0 1 0 0 1 0 0;-1 0 0 0 0 1 0;1 0 0 0 0 0 1];
rank(Aeq)
size(A)
beq = [0 0 0 0 0];
size(b)
A=[];
b=[];
f = [0 0 0 -1 0 0 0];
lb = [0 -10 0 0 -10 -10 -10];
ub = [10 10 10 10 10 10 10];
x = linprog(f,A,b,Aeq,beq,lb,ub)
The optimal solution that I obtain is,
x =
0
-10
0
10
10
0
0
From what I observe, the values in the solution vector are mostly the values of the bounds.
Could someone suggest how to obtain other solutions?
For instance,
x = [ 4.857 ; 5.143; 4.857;10;5.143;4.857;4.857] is also a solution.
I am not sure how other solutions can be obtained.
0 Comments
Accepted Answer
John D'Errico
on 24 Jan 2019
Linear programming returns the optimal solution. You want another optimal solution? How many optimal solutions are there? :) Ok, lets get serious.
You have a very simple LP problem. Maximiize variable 4, subject to a set of equality constraints, and bound constraints.
Linprog finds this solution:
xlinprog =
0
-10
0
10
10
0
0
with x(4) at the upper bound. So, what flexibility do you have?
Aeq
Aeq =
-1 0 1 0 0 0 0
1 -1 0 -1 0 0 0
0 1 0 0 1 0 0
-1 0 0 0 0 1 0
1 0 0 0 0 0 1
Aeq*x == 0 is a system of 5 equations in 7 unknowns. So there will be two degrees of freedom in there.
format long g
Anull = null(Aeq)
Anull =
0.297698713811669 -0.354487221322
-0.342890338788254 -0.48945793762951
0.297698713811669 -0.354487221322
0.640589052599922 0.134970716307511
0.342890338788254 0.48945793762951
0.297698713811669 -0.354487221322
-0.297698713811669 0.354487221322
That is, we can add any linear combination of the two columns of null(Aeq) to the solution, and get another solution that lies on the equality constraint. Efectively, we can create a family of quasi-solutions that look like
xlinprog + s*Anull(:,1) + t*Anull(:,2)
The problem is to do so, while x(4) remains unchanged. What does that mean? We can now formulate one other simple equation:
Anull(4,1)*s + Anull(4,2)*t == 0
Thus, given two parameters s and t, we can add any amount of these two vectors to the linprog solution, as long as that additional equality holds true. In turn, that lets us eliminate t, as
t = -s*Anull(4,1)/Anull(4,2)
So the family of all "solutions" to your problem are of the form
xlinprog + s*Anull(:,1) - s*Anull(:,2)*Anull(4,1)/Anull(4,2)
That should make it clear there is a one dimensional family of probably infinitely many solutions, but those other solutions will not lie at a vertex of the solution space, but be interior to that domain. Since linprog will return a solution that always lies at a boundary vertex, it gave you what it did.
xfamily = @(s) xlinprog + s*Anull(:,1) - s*Anull(:,2)*Anull(4,1)/Anull(4,2)
lb = [0 -10 0 0 -10 -10 -10]';
ub = [10 10 10 10 10 10 10]';
Now, we can quickly observe that xfamily(0) is the linprog solution, and that the entire set of validly equivalent solutions lies in the interval on s of [0,5.05014177670167].
[lb,xfamily(0), xfamily(5.05014177670167),ub]
ans =
0 0 10 10
-10 -10 7.1054e-15 10
0 0 10 10
0 10 10 10
-10 10 -7.1054e-15 10
-10 0 10 10
-10 0 -10 10
Anything between those two values will result in a validly equivalent solution, that remains within the bounds.
Easy peasy.
0 Comments
More Answers (1)
Torsten
on 24 Jan 2019
Linear optimization codes that use the simplex method always return corner points of the LP problem as solution.
Maybe
options = optimoptions('linprog','Algorithm','interior-point')
can help.
1 Comment
See Also
Categories
Find more on Linear Least Squares in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!