Problem doing a left division on two matrices

Hi everyone,

Me and a collaegue of mine have different results when doing a left division of 2 matrices.

Here are the 2 matrices :

C =

   2.000000000000000  -0.005154639175258
  -0.005154639175258   0.000013285152514

SecondMembre_DIFF =

 -17.239054990531251
   0.044430554099307

The operation is :

AB_DIFF=C\SecondMembre_DIFF.

Results are :

- For me (Matlab R2014b 64 bits) : AB_DIFF =

   1.0e+04 *

   0.005266876946754
   2.377985922156901

- For my colleague (Matlab R2015b 64 bits) : AB_DIFF =

   NaN
   NaN

Also when I do the same operation in a different way:

AB_DIFF=inv(C)*SecondMembre_DIFF,

I've the same result as my colleague (NaN matrix), which seems normal because det(C)=0, so inv(C)=inf.

Anyone knows the reason of this difference? Is it a bug on my version of Matlab?

Thanks in advance for your help!

Answer 1

[John D'Errico]

It is not a bug in MATLAB. It is a bug in your incomplete understanding of floating point arithmetic, as well as your incomplete understanding of linear algebra. In no special order...

1. NEVER use det to test a matrix as being singular or not.

C = [2.000000000000000  -0.005154639175258
    -0.005154639175258   0.000013285152514];

det(C)
ans =
      8.95527277550504e-16

det(1e12*C)
ans =
          895533710.718155

Is C any less singular than is 1e12*C? Of course not!!!!!!

Instead, try rank or cond on C to learn about its singularity. (You can use tools like svd to learn even more.)

rank(C)
ans =
     1

rank(1e12*C)
ans =
     1

cond(C)
ans =
      4.46670368349431e+15

cond(1e12*C)
ans =
      4.46667327796135e+15

See that both rank and cond are able to show that C is essentially numerically singular.

I repeat myself. NEVER USE DET FOR THESE THINGS. In fact, det is rarely the right tool for any computation. (There are a few special cases where it can be useful, IF used with the proper care.)

2. Be careful how you pass around numbers to a colleague. For example, copy and paste from the command window is NOT sufficient to pass around the exact numbers in the array C. In fact, tiny errors in how you create those arrays is sufficient to change the results by a huge amount. So as I pasted in the matrix C, I probably used a slightly different matrix than did you or your colleague. This probably explains much of the differences you have seen.

d = [ -17.239054990531251
       0.044430554099307];

C\d
ans =
         -8.63478457507908
         -5.91974696762184

As you can see, I got different results than did you. Why? Again, look at the condition number. LEARN WHAT THE CONDITION NUMER TELLS YOU. Trying to solve a linear system that has a large condition number will magnify any TINY floating point errors in the inputs by the condition number. What was the condition number? It was on the order of 1e15. So what should you expect? Random GARBAGE for results.

3. Finally, use backslash to solve linear systems, not inv. In some cases, pinv is a good idea, especially when the system is singular, as it is here.

pinv(C)*d
ans =
       -8.61947023977726
      0.0222151294839632

It is important to understand the difference between the various solutions, why they are different, but that is potentially the subject of a class in linear algebra, and I've already written a lot.