1

Not quite sure what this means. "Warning: Matrix is singular to working precision."

I have a 3x4 matrix called matrix bestM matrix Q is 3x3 of bestM and matrix m is the last column of bestM

I would like to do C = -Inverse matrix of Q * matrix m and I get that warning and C =[Inf Inf Inf] which isn't right because i am calculating for the camera center in the world

bestM = [-0.0031 -0.0002 0.0005 0.9788;
         -0.0003 -0.0006 0.0028 0.2047;
         -0.0000 -0.0000 0.0000 0.0013];   

Q = bestM(1:3,1:3);
m = bestM(:,4);

X = inv(Q);
C = -X*m;
disp(C);
Amro
  • 123,847
  • 25
  • 243
  • 454
ealeon
  • 12,074
  • 24
  • 92
  • 173

2 Answers2

5

A singular matrix can be thought of as the matrix equivalent of zero, when you try to invert 0 it blows up (goes to infinity) which is what you are getting here. user 1281385 is absolutely wrong about using the format command to increase precision; the format command is used to change the format of what is shown to you. In fact the very first line of the help command for format says

format does not affect how MATLAB computations are done.

dvreed77
  • 2,217
  • 2
  • 27
  • 42
1

As found here, a singular matrix is one that does not have an inverse. As dvreed77 already pointed out, you can think of this as 1/0 for matrices.

Why I'm answering, is to tell you that using inv explicitly is almost never a good idea. If you need the same inverse a few hundred times, it might be worth it, however, in most circumstances you're interested in the product C:

C = -inv(Q)*m

which can be computed much more accurately and faster in Matlab using the backslash operator:

C = -Q\m

Type help slash for more information on that. And even if you happen to find yourself in a situation where you really need the inverse explicitly, I'd still advise you to avoid inv:

invQ = Q\eye(size(Q))

Below is a little performance test to demonstrate one of the very few situations where the explicit inverse can be handy:

% This test will demonstrate the one case I ever encountered where
% an explicit inverse proved useful. Unfortunately, I cannot disclose
% the full details without breaking the law, but roughly, it came down
% to this: The (large) design matrix A, a result of a few hundred 
% co-registrated images, needed to be used to solve several thousands 
% of systems, where the result matrices b came from processing the 
% images one-by-one. 
%
% That means the same design matrix was re-used thousands of times, to 
% solve thousands of systems at a time. To add to the fun, the images
% were also complex-valued, but I'll leave that one out of consideration
% for now :)  

clear; clc

% parameters for this demo
its = 1e2;
sz  = 2e3;
Bsz = 2e2;

% initialize design matrix
A = rand(sz);

% initialize cell-array to prevent allocating memory from consuming 
% unfair amounts of time in the first loop. 
% Also, initialize them, NOT copy them (as in D=C,E=D), because Matlab 
% follows a lazy copy-on-write scheme, which would influence the results
C = {cellfun(@(~) zeros(sz,Bsz), cell(its,1), 'uni', false)   zeros(its,1)};
D = {cellfun(@(~) zeros(sz,Bsz), cell(its,1), 'uni', false)   zeros(its,1)};
E = {cellfun(@(~) zeros(sz,Bsz), cell(its,1), 'uni', false)   zeros(its,1)};

% The impact of rand() is the same in both loops, so it has no 
% effect, it just gives a longer total run time. Still, we do the 
% rand explicitly to *include* the indexing operation in the test. 
% Also, caching will most definitely influence the results, because 
% any compiler (JIT), even without optimizations, might recognize the 
% easy performance gain when the code computes the same array over and 
% over again. It probably will, but we have no control over when and 
% wherethat happens. So, we prevent that from happening at all, by 
% re-initializing b at every iteration. 

% The assignment to cell is a necessary part of the demonstration; 
% it is the desired output of the whole calculation. Assigning to cell 
% instead of overwriting 'ans' takes some time, which is to be included 
% in the demonstration, again for cache reasons: the extra time is now
% guaranteed to be equal in both loops, so it really does not matter -- 
% only the total run time will be affected.

% Direct computation
start = tic;
for ii = 1:its    
    b = rand(sz,Bsz);    
    C{ii,1} = A\b;
    C{ii,2} = max(max(abs( A*C{ii,1}-b )));
end
time0 = toc(start);
[max([C{:,2}])   mean([C{:,2}])   std([C{:,2}])]

% LU factorization (everyone's
start = tic;
[L,U,P] = lu(A, 'vector');
for ii = 1:its    
    b = rand(sz,Bsz);    
    D{ii,1} = U\(L\b(P,:));
    D{ii,2} = max(max(abs( A*D{ii,1}-b )));
end
time1 = toc(start);
[max([D{:,2}])   mean([D{:,2}])   std([D{:,2}])]


% explicit inv
start = tic;
invA = A\eye(size(A)); % NOTE: DON'T EVER USE INV()!
for ii = 1:its
    b = rand(sz,Bsz);
    E{ii,1} = invA*b;
    E{ii,2} = max(max(abs( A*E{ii,1}-b )));
end
time2 = toc(start);
[max([E{:,2}])   mean([E{:,2}])   std([E{:,2}])]

speedup0_1 = (time0/time1-1)*100
speedup1_2 = (time1/time2-1)*100
speedup0_2 = (time0/time2-1)*100

Results:

% |Ax-b|
1.0e-12 * % max.   mean     st.dev.
          0.1121   0.0764   0.0159 % A\b   
          0.1167   0.0784   0.0183 % U\(L\b(P,;))
          0.0968   0.0845   0.0078 % invA*b

speedup0_1 = 352.57 % percent
speedup1_2 =  12.86 % percent
speedup0_2 = 410.80 % percent

It should be clear that an explicit inverse has its uses, but just as a goto construct in any language -- use it sparingly and wisely.

Rody Oldenhuis
  • 37,726
  • 7
  • 50
  • 96
  • 2
    Rody, even if you need to solve your system *a few hundred times* you do explicit factorization instead of '\', and not compute inverse. I see no other use of matrix inverse unless your task is to compute a matrix inverse. – angainor Oct 14 '12 at 06:00
  • @angainor: what do you mean with "explicit factorization"? How would you solve a few hundred/thousand systems `Ax=B`, with `A` constant and `B` a few hundred/thousand different m-by-n matrices? – Rody Oldenhuis Oct 14 '12 at 06:23
  • 1
    You run `lu` or `chol` to express `A` as a product of triangular matrices. Solving the system is then equivalent to two trivial triangular solves. – angainor Oct 14 '12 at 06:27
  • See [this thread](http://stackoverflow.com/questions/12831889/efficient-way-to-solve-for-x-in-ax-b-in-matlab-when-both-a-and-b-are-big-matrice/12835140#12835140) – angainor Oct 14 '12 at 06:48
  • Rody, this is simply wrong. First, try to do your tests for matrices larger than 5x5. Say, A and B are 2000x2000. On my PC, LU factorization + backward-forward takes 1s. Computing inverse and invA*B takes 1.5s. Complexity of L/U solves is O(n^2), so for a B-matrix it is the same as matrix-matrix multiply. And it is BLAS3... Also, I have no idea what 'considered too restrictive for most matrices' means. Can you explain? Many problems result in positive-definite matrices, and `chol` is very useful. Do not use it on random matrices, it is *wrong* and does not show anything. Create an SPD matrix. – angainor Oct 15 '12 at 07:55
  • Not to mention the uselessness of this approach for sparse matrices, but it is a different story. Finally, if this does not convince you, computing inverse is less numerically accurate than Gaussian elimination. See help of `inv`, but also a lot of articles lying around. This is more or less basic numerical linear algebra rules... (e.g. [here](http://www.mathworks.com/matlabcentral/newsreader/view_thread/141), [here](http://www.johndcook.com/blog/2010/01/19/dont-invert-that-matrix/)). And yes, I know linking articles from Internet does not prove anything. – angainor Oct 15 '12 at 08:00
  • @angainor: LU is O(n^2), so is a multiply. The thing is, you have to do *two* solves for `chol` or `LU` (plus one multiply when using `P`), and only *one* multiply with `invA`. See my edit. If these results do not convince *you*, well, then you're not doing science anymore :) And my remark that `chol` is too restrictive is based on the fact that `chol` is not a *generally applicable method* -- it has requirements on the matrix' properties, which makes it applicable only to a small subset of all possible matrices, however common they may be; I wasn't considering special or corner-cases. – Rody Oldenhuis Oct 15 '12 at 08:34
  • @angainor: Plus, this test is intended to show that solving systems *repeatedly* for different result vectors, but **with always the same design matrix**, the explicit inverse is useful. It is rare, but this is a situation I encountered IRL, and where I now **know** you need the explicit inverse :) – Rody Oldenhuis Oct 15 '12 at 08:38
  • Rody, I have edited your answer and added a re-written benchmark with the results. About the complexity, you do not do anything twice. Those are triangular matrices - total number of non-zero entries is the same. The performance hit, as you see in my benchmark, is rather small. I do not suspect anyone would want to sacrifice the *provably better* numerical performance of Gaussian elimination for a 5% speedup. – angainor Oct 15 '12 at 09:30
  • @angainor: Have you read my second reply? I think the confusion comes from you misunderstanding what I was trying to show. My test was intended to show how ***repeated solves with changing b*** are sped-up by using explicit inverse. – Rody Oldenhuis Oct 15 '12 at 10:17
  • But they are not. I removed `b=rand` from the loop because it is irrelevant - I do not want to time `rand`. It can equally well be there. I do perform the solves, don't I? Or do I precalculate something? Well, `b(p,:)`, but if you want just put it into the loop - the performance difference will be 10-15%, not 5. And only because MATLAB uses column storage. And anyway you can easily permute the whole problem you solve and avoid explicit permutation of `b` in the loop. That is just *dof* numbering. You should change it if it helps you to solve your system faster. Standard practice. – angainor Oct 15 '12 at 10:22
  • The point is you used mat-mat multiply for permutation. You do not need to do this - this is a trivial row permutation. – angainor Oct 15 '12 at 10:26
  • Rody, did we arrive at any conclusion, or am I missing something? I seriously can not see any advantage. – angainor Oct 15 '12 at 12:35
  • @angainor: Yes, I see -- this whole discussion comes from a misunderstanding. I intended to do something other than what you think. See my latest edit. I went wrong on the permutation `P`, my bad. Also I did a few things a bit less efficient than strictly necessary. I've last done this 2 years ago :) And last: we seem to share a similar irregularity of being at our computers, and a similar tendency for diving in the true depths of a problem, which results in **way too long** comment threads all over the place. Perhaps, e-mail (or SO chat) is a decent alternative for the future? – Rody Oldenhuis Oct 15 '12 at 12:45
  • @RodyOldenhuis: The whole point of my "test" was to give more of a demonstration of a rare corner case (for general matrices) where the explicit inv *can* be useful. The *point* I was trying to make in the answer was that you should avoid it at all times! That seems a bit lost in translation now, so in principle, we are not in disagreement. – Rody Oldenhuis Oct 15 '12 at 12:48
  • Good :) Maybe those comments are indeed too long - since they ask to avoid them ;) next time we can take it elsewhere .. – angainor Oct 15 '12 at 12:56
  • As I was just discussing this, I am on my mission to clarify: [**Don't ever use `A\eye(size(A))`**](http://ideone.com/rAbn9S)! – knedlsepp Feb 16 '15 at 01:01
  • @knedlsepp indeed, since I posted this answer I have learned that there are always better alternatives to it. Thanks :) – Rody Oldenhuis Feb 16 '15 at 08:58
  • @RodyOldenhuis: Would you mind deleting *And even if you happen to find yourself in a situation where you really need the inverse explicitly, I'd still advise you to avoid `inv`.* then? As in fact `inv` is a better alternative than `A\eye(size(A))`. – knedlsepp Feb 16 '15 at 10:39
  • @knedlsepp not always though. The best solution is always highly specific to the problem, and I contend it's *neither* `inv()` *or* backslash-with-`eye`. You have to admit that your test is anything but exhaustive. In any case, this answer is not the best place to put all of the specifics that need to be covered on this discussion, and moreover, I think it's better SO practice to post your own answer :) – Rody Oldenhuis Feb 16 '15 at 11:41
  • @RodyOldenhuis: I could of course post my own answer, but there is a high possibility of it being ignored by the community because of RepOverflow. It is preposterous to say that my tests aren't exhaustive. Would you expect to try *all square matrices*? My tests cover all matrices listed in the MATLAB `gallery` which are listed as *ill-conditioned*, as well as simple well-conditioned matrices generated by `rand`. Certainly there will be better approaches to both `inv` and `mldivide`-`eye`. But this simply gives the false impression that `mldivide`-`eye` is somehow better than `inv`. – knedlsepp Feb 16 '15 at 12:02
  • @RodyOldenhuis: For the sake of scientific reasoning: http://ideone.com/ZQu2Gj This includes **all** `gallery` matrices, each for varying sizes. There is no reason to tell people to use `mldivide`-`eye` instead of `inv`. – knedlsepp Feb 16 '15 at 12:44
  • @knedlsepp ...downvoting because I disagree with you? What good sportsmanship! Look, the answer is not the accepted one, doesn't even answer the question, is 2.5 years old and the question has low views. This is simply not the place for this discussion. We agree on the overarching matter, please don't start a flamewar on the details. – Rody Oldenhuis Feb 16 '15 at 13:29
  • @knedlsepp by the way, *you* can edit my answer you know :) just leave the original in. – Rody Oldenhuis Feb 16 '15 at 13:30
  • @RodyOldenhuis: This is not a flamewar. I think I am acting exactly in the spirit of SO. I find a problem in an answer. I try to contact the OP to explain to him there is a problem with his post. The OP ignores the matter, so I decide this answer is harmful to anybody learning numerics and decide to downvote. I *could* edit your answer, but that is not in the spirit of SO, as one should only change answers to clarify, but **never** to change the OP's original intent! This would only result in me being banned. – knedlsepp Feb 16 '15 at 13:32
  • *doesn't even answer the question* is also a perfectly good reason to downvote BTW! Did you even have a look at the tests? – knedlsepp Feb 16 '15 at 13:34
  • If you don't care to execute it on your machine: All the matrices of `gallery`; each for varying sizes: http://ideone.com/ZQu2Gj 261 times `inv` gives better results than `mldivide`-`eye`. 111 times `mldivide`-`eye` gives better results than `inv`. To me this a compelling case *for* `inv` not *against*! – knedlsepp Feb 16 '15 at 16:43