user7829436
user7829436

Reputation:

How to check a matrix is not singular with a computer

I would like to know a robust way to check if a matrix is not singular using a computer. I know using the determinant (requiring it to be non zero) can be misleading as it cannot distinguish between a case where a matrix is indeed singular and due to numerical errors you get a very small value (e.g.~10^-12) and a case like 10^-12*I which again gives a very small determinant whereas the matrix is definitely not singular (it is orthogonal).

There is a good link (How to find out if a matrix is singular?) arguing one can use the condition of the matrix, or in other words the ratio of biggest singular value over smallest singular value.

However isn't this problematic again? The 2x2 matrix [10^8 0; 0 10^-8] is orthogonal and therefore definitely not singular, however its singular values are 10^8 and 10^-8 (condition number 10^16), therefore according to the link above it would be classified as singular.

Is a correct way to normalise the rows before singular value decomposition and then simply check the smallest singular value to be small (e.g. smaller than 10^-7)?

Upvotes: 4

Views: 977

Answers (2)

Giorgos Altanis
Giorgos Altanis

Reputation: 2760

The condition number does not decide whether the matrix is singular or not, it shows whether the solutions obtained are robust w.r.t. RHS of the linear system or not. A nonsingular matrix might have very bad condition.

A matrix A is singular if any of its columns can be expressed as a linear combination of the remaining columns. This is equivalent to saying that A is nonsingular if and only if it is full rank. So a rank-revealing factorization should be used. The recommended way (at least this is what I was taught back then!) to do this is the Rank-revealing QR factorization.

Upvotes: 0

Luis Mendo
Luis Mendo

Reputation: 112699

The condition number of a matrix measures how sensitive the linear system A x = b is with respect to perturbations in b. A large condition number means that a relative perturbation in b can be greatly amplified in the solution x.

The term relative perturbation here means how much the original and perturbed vectors differ, compared to the size of the original vector. Specifically, let b1 denote the perturbed version of b and x1 the corresponding perturbed solution. The relative perturbation in b (or in x) is defined as norm(b-b1)/norm(b) (or norm(x-x1)/norm(x)).

With this definition, the significance of the condition number can be phrased as follows: a large condition number means that norm(x1-x)/norm(x) can be much larger than norm(b1-b)/norm(x). (For a proof of this result see Applied Linear Algebra (3rd ed.) by B. Noble and J.W. Daniels, page 271; or this Q&A at Mathematics Stack Exchange.)

Your example matrix is

A = [10^8 0; 0 10^-8];

with

>> cond(A)
ans =
     1.0000e+16

Consider the following system using this matrix:

b = [1; 0]; % original b
x = A\b; % original solution
b1 = b + 0.01; % perturbed b
x1 = A\b1; % perturbed solution

This gives the solutions x = [1e-8; 0] (original) and x1 = [1.01e-08; 1e6] (perturbed). The relative perturbation in the solution is

>> norm(x-x1)/norm(x)
ans =
     1.0000e+14

As you can see, it is much larger than the relative perturbation that was introduced in b,

>> norm(b-b1)/norm(b)
ans =
   0.0141

Note that for other choices of b the relative perturbation may not get so dramatically amplified. The condition number characterizes the worst-case behaviour over all possible choices of b.

On the other hand, consider the row-normalized version of A:

B = A;
B(1,:) = B(1,:)/norm(B(1,:));
B(2,:) = B(2,:)/norm(B(2,:));

This is just the identity matrix:

>> B
B =
     1     0
     0     1

which, of course, is as well-conditioned as it gets. So now the relative perturbation in x is not amplified with respect to that in b:

y = B\b;
y1 = B\b1;

gives

>> norm(y-y1)/norm(y)
ans =
   0.0141

which is the same as the relative perturbation in b.

Upvotes: 2

Related Questions