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

 
 

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

Submitted by:  Mattias Linde <linde>
Submitted on:  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

Add a New Comment (Rich Markup)
   

You are not logged in

Please log in, so followups can be emailed to you.

 

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)
  •  

    Do you think this task is very important?
    If so, you can add your encouragement to it.
    This task has 0 encouragements so far.

    Only project 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.5