13

I'm having trouble solving a system of the form Ax=B

The solution to the system should be 

x = inv(A)*B

However, this doesn't work.

I get the following error message when I try the above line of code:

Warning: Matrix is close to singular or badly scaled.
     Results may be inaccurate. RCOND = 1.156482e-018.

It seems that matlab is having trouble inverting the matrix that I've specified. I tried to verify that the inverse function was working properly by typing in inv(A)*A

This should give the identity matrix, however I got the same error and some garbage numbers.

This is the A matrix that I'm using:

A = [5/2   1/2  -1     0     0    -1/2  -1/2   0     0
     1/2   1/2   0     0     0    -1/2  -1/2   0     0 
    -1     0     5/2  -1/2  -1     0     0    -1/2   1/2
     0     0    -1/2   1/2   0     0     0     1/2  -1/2
     0     0    -1     0     3/2  -1/2   1/2   0     0
    -1/2  -1/2   0     0    -1/2   2     0    -1     0  
    -1/2  -1/2   0     0     1/2   0     1     0     0 
     0     0    -1/2   1/2   0    -1     0     2     0 
     0     0     1/2  -1/2   0     0     0     0     1]

Any ideas as to why this isn't working? I also tried to convert A to a sparse matrix (sparse(A)), and then run the inverse command. No dice.

Stewie Griffin
  • 14,889
  • 11
  • 39
  • 70
Tushar Garg
  • 769
  • 3
  • 8
  • 16
  • 3
    This is fundamentally a conceptual, math problem, not a programming one. Your syntax is correct. Warnings are not errors, they're warnings. – Karl Knechtel Nov 02 '11 at 02:43
  • Thanks for your response Karl. I know that the inverse of a matrix multiplied by itself yields the identity matrix. For some reason this isn't working in matlab. This leads me to believe that there is an issue with the way I am inputting the data and not just a mathematics problem. – Tushar Garg Nov 02 '11 at 02:43
  • 6
    You don't seem to understand that not all matrices have an inverse. That is a mathematics problem. – Karl Knechtel Nov 02 '11 at 02:46
  • 2
    Even if a matrix is invertable, you shouldn't use inv(A) to solve a linear problem. Use mldivide instead. –  Feb 03 '15 at 02:54

4 Answers4

21

The problem is indeed in your mathematics. The matrix you provided isn't of full rank, so it isn't invertible. You could verify that manually (haven't taken the time to do so), but MATLAB already points this out by showing that warning.

Since you are working with floating point numbers, this sometimes causes other subtle problems, one of which you can see in the result of det(A), which is in the order of 1e-16, i.e. machine precision or 0 in practice.

You can see that this Matrix is not of full rank by executing the rank function: rank(A) = 8. For a 9x9 matrix, this indeed means that the matrix is not invertible for doubles (as the rank function accounts for machine precision).

If you want to use MATLAB to get a result that corresponds to a manual calculation, you can use the Symbolic Toolbox and its vpa (variable precision arithmetic) to work around possible numerical problems at the cost of a slower calculation.

B = [5  1 -2  0  0 -1 -1  0  0;
     1  1  0  0  0 -1 -1  0  0;
    -2  0  5 -1 -2  0  0 -1  1;
     0  0 -1  1  0  0  0  1 -1;
     0  0 -2  0  3 -1  1  0  0;
    -1 -1  0  0 -1  4  0 -2  0;
    -1 -1  0  0  1  0  2  0  0;
     0  0 -1  1  0 -2  0  4  0;
     0  0  1 -1  0  0  0  0  2];
A = B/2;
size(A)    % = [9 9]
det(A)     % = -1.38777878078145e-17
rank(A)    % = 8
C = vpa(A);
det(C)     % = 0.0
rank(C)    % = 8

Both with VPA and floating points you will get that the rank is 8, the size is [9 9] and the determinant is practically 0, i.e. singular or not invertible. Changing a few entries might make your matrix regular (non-singular), but it is not guaranteed to work and it will solve a different problem.

To solve your actual problem A*x=b for x, you can try to use mldivide (a.k.a. the backslash operator) or a Moore-Penrose pseudo-inverse:

x1 = A\b;
x2 = pinv(A)*b;

But do remember that such a system does not have a unique solution, so both the pseudo-inverse and the backslash operator may (and in this case will) return very different solutions, whether any of them is acceptable really depends on your application.

Egon
  • 4,757
  • 1
  • 23
  • 38
  • Anyone reading this and not knowing the mathematics should read the Matlab docs to `mldivide` and `pinv`. If the matrix is rank-deficient by one, as in this case, the problem might lend itself to specifying a further constraint such as norm(x) == 1. – Unapiedra Jun 30 '14 at 13:14
  • @Unapiedra I'm working on a similar case where my matrix A is rank-deficient by one. I'm using the approach pinv(A)*b. I'd greatly appreciate it if you could elaborate a bit on what you mean by specifying the constraint norm(x) == 1. What is x? Where is this expression to be implemented. What does this do? – André Foote Feb 01 '16 at 19:17
11

It means exactly what it says. The matrix is singular, which means it can't really be inverted. Not all matrices can.

In geometrical terms, you have a matrix that transforms one 9-dimensional object into another, but flattens one dimension out completely. That can't be undone; there's no way to tell how far to pull things out in that direction.

Karl Knechtel
  • 62,466
  • 11
  • 102
  • 153
  • Okay, lets pretend it is indeed singular. What if I change the first 3 entries in the matrix to some random integers. lets say, 5 1 and 3. If the matrix was indeed singular before (determinant = 0), changing a few of the values should change its determinant as well. After I compute the inverse on the new matrix, I get the same error. So now whats going on?? – Tushar Garg Nov 02 '11 at 03:02
2

The matrix is singular, consider B=2*A below:

B = [5  1 -2  0  0 -1 -1  0  0;
     1  1  0  0  0 -1 -1  0  0;
    -2  0  5 -1 -2  0  0 -1  1;
     0  0 -1  1  0  0  0  1 -1;
     0  0 -2  0  3 -1  1  0  0;
    -1 -1  0  0 -1  4  0 -2  0;
    -1 -1  0  0  1  0  2  0  0;
     0  0 -1  1  0 -2  0  4  0;
     0  0  1 -1  0  0  0  0  2]

det(B)

0
Alex
  • 5,863
  • 2
  • 29
  • 46
  • Okay, thanks for your response. Interestingly enough when I plug in your matrix into matlab, I get the following for det(B): -1.7764e-015. I'm guessing this has to do with the variable type, but how did you get yours to display 0 – Tushar Garg Nov 02 '11 at 03:22
  • 1
    You are seeing floating point errors creeping in. Anything that close to 0 is 0, especially if you start off with values close to 1. – Alex Nov 02 '11 at 03:26
  • How do I change this? (just out of curiosity) – Tushar Garg Nov 02 '11 at 03:30
1

bicgstab(A,b,tol,maxit), an iterative solver, was able to solve a singular linear system A*x=b for a singular matrix A: 

size(A)=[162, 162] 
rank(A)=14 
cond(A)=4.1813e+132

I used: 

tol=1e-10; 
maxit=100;

None of the above-mentioned (including svd, \, inv, pinv, gmres) worked for me but bicgstab did a good job. bicgstab converged at iteration 4 to a solution with relative residual 1.1e-11. It works fast for sparse matrices.

See documentation here: https://uk.mathworks.com/help/matlab/ref/bicgstab.html

Wolfie
  • 27,562
  • 7
  • 28
  • 55
Feruza
  • 11
  • 1