## How to get all the solutions to an underdetermined system of linear equations using lsqnonneg()?

### MUKTIKANTA (view profile)

on 13 Oct 2018
Latest activity Commented on by MUKTIKANTA

on 14 Oct 2018

### Bruno Luong (view profile)

I am trying find out all the positive solutions to an underdetermined system of linear equations and my code is as follows:
format rat
a=[3 22 119 620 3 1 3 3 2 20 20 22 22 119 119 3 3 3 3 3 3 1 1 1 20 20 20 20 2 1 2 20 20 3 3 3 3 3 3 1 1 1 2 2 1;
3 0 0 0 3 1 3 3 2 0 0 0 0 0 0 3 3 3 3 3 3 1 1 1 0 0 0 0 2 1 2 0 0 3 3 3 3 3 3 1 1 1 2 2 1;
0 22 0 0 3 1 0 0 0 20 20 22 22 0 0 3 3 3 3 0 0 1 1 1 20 20 20 20 2 1 0 20 20 3 3 3 3 3 3 1 1 1 2 2 1;
0 0 119 0 0 0 3 0 0 20 20 0 0 119 119 3 3 0 0 3 3 1 0 0 20 20 20 20 0 0 2 20 20 3 3 3 3 3 3 1 1 1 2 2 1;
0 0 0 620 0 0 0 3 2 0 0 22 22 119 119 0 0 3 3 3 3 0 1 1 20 20 20 20 2 1 2 20 20 3 3 3 3 3 3 1 1 1 2 2 1];
b=[620;3;22;119;620];
w=lsqnonneg(a,b);
Each time I run this code, it gives me only one solution. How can I get all the solution? Is there any method available in MATLAB to achieve this?

Walter Roberson

### Walter Roberson (view profile)

on 13 Oct 2018
That is rather the wrong function to use.
An underdetermined system of linear equations describes hyperplane (or possibly a plane or a line.)
You have 45 variables and 5 equations. Given any arbitrary values for the last 40 variables, exact values of the first 5 can be calculated.
The set of all positive solutions could only be described by boundary equations. I would have to work out how many of those you would need. It is late at night here and one avenue of thought is suggesting that you might need 45!/5! so I would want to think about that when I am more awake. Maybe it is 45!/40! ?

### Bruno Luong (view profile)

on 13 Oct 2018
Edited by Bruno Luong

### Bruno Luong (view profile)

on 13 Oct 2018

The number of solutions is infinity. Explanation:
>> size(a)
ans =
5 45
In your case you have 5 eqt and 45 unknowns, and
>> rank(a)
ans =
5
confirms your matrix is full row-rank. meaning all the 5 rows of A are independent vectors in R^45.
One particular solution can be obtained with backshlash (or lsqnonneg if you prefer)
>> x0=a\b;
>> fprintf('x0 = %s\n', mat2str(x0));
x0 = [0;0;5.93708057585062e-17;0.77741935483871;0;0;0;0;0;0;0;0.863636363636363;0;0.974789915966387;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0.999999999999998;0;0;0;0;0;0;0;0;0;0;0]
The kernel of a, meaning a basis of equation A*X = 0 can be obtained with NULL command.
>> X=null(a);
>> size(X)
ans =
45 40
It tells that the kernel has dimension 40 and a maths theorem tells it's equals to size(a,2)-rank(a).
To confirm X is the kernel, you can check with frobenus norm of max values (they are 0 at the numerical precision)
>> y=a*X;
>> norm(y,'fro')
ans =
1.5079e-13
>> max(max(abs(y)))
ans =
7.5533e-14
Now all solutions of a*x=b can be written as
x = x0 + X*r
where r is any arbitrary vector of size 40 x 1. The set of solutions is an affine space of dimension 40 and there are of course an infinity vectors in this huge space. To check it try:
r = rand(40,1); % any random vector
x = x0 + X*r; % general form of the solution a*x=b
check if a*x = b
>> max(abs(a*x-b))
ans =
1.1369e-13
If you add the constraints (x >= 0), meaning
x(i) >= 0 for all i = 1,...45
Then the general solution is the intersection of 45 half-hyperspace and the affine space of dimension 40.
It will be a polytope of dimension 40 living in R^45. It is still huge. This polytope is actually a an interior of the convex combination of n vertices, n can be as high as:
>> nchoosek(45,40)
ans =
1221759
each of them is one particular solution of
a*x = b
x >= 0
There is no way you can enumerate them all and make on top of that an arbitrary sub-convex combination if that is your intention.

Walter Roberson

### Walter Roberson (view profile)

on 13 Oct 2018
Thanks, Bruno. I would have found it tough to pull all those steps together.
Bruno Luong

### Bruno Luong (view profile)

on 13 Oct 2018
Note: what I wrote is not entirely correct, sometime the polytope is unbounded if no other constraint on x is introduced and it can not express as the convex combination of finite set of vertexes.
However I hop it discourage people trying to beat the market with such brute force method.
MUKTIKANTA

### MUKTIKANTA (view profile)

on 14 Oct 2018
Thanks, Bruno for such a well-explained answer.