bugGNU Octave - Bugs: bug #41454, Inaccurate/incorrect result for...

 
 

bug #41454: Inaccurate/incorrect result for inverse of symmetric matrix

Submitter:  Mattias Linde <linde>
Submitted:  Tue 04 Feb 2014 11:32:42 AM UTC
   
 
Category:  Octave Function Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Incorrect Result
Status:  None Assigned to:  None
Originator Name:  Open/Closed:  * Open
Release:  * dev Operating System:  * GNU/Linux
Fixed Release:  None Planned Release:  None
* Mandatory Fields

Add a New Comment Rich Markup
   

Sat 30 Jan 2016 05:59:31 AM UTC, comment #2: 

I can reproduce the error, but get yet another set of answers


matrixSize =

   239   239

matrixRank =  239
matrixCond =    4.9337e+07
machineEps =    2.2204e-16
withBackslashEye =

   2.5171e-09   3.8788e-08

withInv =

   5.2161e-01   4.8754e-08

symmetryHs = 0
symmetryMyInv =    2.5151e-08
symmetryOctInv =  0.22344


However, note also that I get


norm (inv(octInv) - Hs)
ans =    1.0759e-07


Perhaps the problem is just much less well conditioned than the condition number suggests.  I note that the top left 145x145 block is just the identity.  I assume that makes the condition number very good, despite the complexity in the remainder of the matrix.  (This is where I should learn some linear algebra...  Could either of you tell me the definition of the condition number?)

Lachlan Andrew <lachlan>
Mon 17 Mar 2014 08:50:57 AM UTC, comment #1: 

I can reproduce that error. Anyhow, I'm confused about the small differences that occured.


octave invprecision.m
GNU Octave, version 3.8.0
Copyright (C) 2013 John W. Eaton and others.
This is free software; see the source code for copying conditions.
There is ABSOLUTELY NO WARRANTY; not even for MERCHANTABILITY or
FITNESS FOR A PARTICULAR PURPOSE.  For details, type 'warranty'.

Octave was configured for "i686-pc-linux-gnu".

Additional information about Octave is available at http://www.octave.org.

Please contribute if you find this software useful.
For more information, visit http://www.octave.org/get-involved.html

Read http://www.octave.org/bugs.html to learn how to submit bug reports.
For information about changes from previous versions, type 'news'.

matrixSize =

   239   239

matrixRank =  239
matrixCond =    4.9337e+07
machineEps =    2.2204e-16
withBackslashEye =

   3.7558e-09   3.4613e-08

withInv =

   3.8905e-01   3.5550e-08

symmetryHs = 0
symmetryMyInv =    2.2922e-08
symmetryOctInv =  0.18704


Max Görner <maxg>
Tue 04 Feb 2014 11:32:42 AM UTC, original submission:  


Double precision calculations, machine epsilon at 2.2204e-16.
Input is a symmetric matrix with condition number 4.9337e+07
so at best I cant expect results with more precision than e-09.

AB = BA = Identity, where B = inverse(A)

Since A is symmetric, B should also be symmetric.

Problems I see is that:

  • B is not close to symmetric.
  • norm(AB-I) is 5e-1 where as norm(BA-I) is 5e-8


Problem seen in all the versions I've tested: 3.2.4,
3.6.0 and build from dev sources as of today.

Script:

load('hs.mat')

matrixSize = size(Hs)
matrixRank = rank(Hs)
matrixCond = cond(Hs)

machineEps = eps

I = eye( size(Hs) );

myInv  = Hs \ I;
octInv = inv(Hs);

% AB = BA = I
withBackslashEye = [ norm(full( Hs*myInv-I))  norm(full(myInv*Hs-I)) ]
withInv  = [ norm(full( Hs*octInv-I))  norm(full(octInv*Hs-I)) ]

symmetryHs     = norm(full(Hs-Hs') )
symmetryMyInv  = norm(full(myInv-myInv'))
symmetryOctInv = norm(full(octInv-octInv'))


Results from octave

matrixSize =

   239   239

matrixRank =  239
matrixCond =    4.9337e+07
machineEps =    2.2204e-16
withBackslashEye =

   3.6139e-09   5.4825e-08

withInv =

   5.0890e-01   5.1131e-08

symmetryHs = 0
symmetryMyInv =    2.9754e-08
symmetryOctInv =  0.17683


Results from matlab 8.0:

withBackslashEye =
   1.0e-07 *
    0.0332    0.4951
withInv =
   1.0e-07 *
    0.0332    0.4951
symmetryHs =
     0
symmetryMyInv =
   2.1488e-08
symmetryOctInv =
   2.1488e-08



I also tested changing from using sparse matrices, Hs = full(Hs) and that made things far worse.

Mattias Linde <linde>

 

(Note: upload size limit is set to 16384 kB, after insertion of the required escape characters.)

Attach Files:
   
   
Comment:
   

Attached Files
file #30454:  invprecision.m added by linde (434B - text/x-objcsrc)
file #30455:  hs.mat added by linde (7KiB - application/octet-stream)

 

Depends on the following items: None found

Items that depend on this one: None found

 

Carbon-Copy List
  • -email is unavailable- added by lachlan (Posted a comment)
  • -email is unavailable- added by maxg (Posted a comment)
  • -email is unavailable- added by linde (Submitted the item)
  •  

    There are 0 votes so far. Votes easily highlight which items people would like to see resolved in priority, independently of the priority of the item set by tracker managers.

    Only group members can vote.

     

    Follow 2 latest changes.

    Date Changed by Updated Field Previous Value => Replaced by
    2014-02-04 linde Attached File- Added invprecision.m, #30454
        Attached File- Added hs.mat, #30455

    Back to the top

    Powered by Savane 3.13-758e.
    Corresponding source code