bugGNU Octave - Bugs: bug #51246, pinv(0) different from Matlab

 
 

bug #51246: pinv(0) different from Matlab

Submitted by:  Guillaume <gyom>
Submitted on:  Fri 16 Jun 2017 12:26:04 PM UTC  
 
Category:  Octave Function Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Regression
Status:  Fixed Assigned to:  None
Originator Name:  Guillaume Open/Closed:  Closed
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.

 

( Jump to the original submission)

Fri 06 Apr 2018 05:47:45 PM UTC, comment #19:

It looks like this has been resolved, but nevertheless I'd like to chime in.

In reference to comment 4:

"I would still open a discussion on Octave Maintainers list just to see if there is a good mathematical reason for returning either 0 or Inf. "

There is a good mathematical reason for returning pinv(0) = 0. Although there are different mathematical notions of generalized inverse, the most widely used is the Moore-Penrose inverse (https://en.wikipedia.org/wiki/Moore%E2%80%93Penrose_inverse). For any matrix A, there exists a unique matrix B with the following four properties:

  1. A*B*A = A
  2. B*A*B = B
  3. A*B is Hermitian
  4. B*A is Hermitian

B is the Moore-Penrose inverse of A. It is easy to verify that the Moore-Penrose inverse of a zero matrix is a zero matrix.

If you look at the documentation of pinv in Matlab, this is precisely what is returned there. This is what Octave should return as well.

Travis Arnold <teerav42>
Fri 03 Nov 2017 05:10:35 PM UTC, comment #18:

I checked in a change along with BIST tests here (http://hg.savannah.gnu.org/hgweb/octave/rev/d3dc76efb38b).

I also found that we were not calculating the tolerance correctly for single matrices. Apparently there was a copy and paste error from the code for double matrices and we were using

in calculations instead of the "float" version of epsilon.

Rik <rik5>
Project Administrator
Fri 03 Nov 2017 04:15:23 PM UTC, comment #17:

>> pinv ([0 1e-15; 0 0], 1e-14)


ans =

0 0
0 0

>> pinv (zeros (2,2))


ans =

0 0
0 0

>> version


ans =

'9.3.0.713579 (R2017b)'

>>

Looks OK.

Michael Godfrey <godfrey>
Project Member
Fri 03 Nov 2017 04:09:16 PM UTC, comment #16:

Hmm, okay. What about the following:

I think numbers below the tolerance should be rounded to zero and this will be equivalent to

which should be a 2x2 zero matrix according to the results so far.

Rik <rik5>
Project Administrator
Fri 03 Nov 2017 03:55:32 PM UTC, comment #15:

Last pinv (D) needs fixing:
D = diag ([0 0 0])

D =

0 0 0
0 0 0
0 0 0

pinv (D)

ans =

0 0 0
0 0 0
0 0 0

>> version


ans =

'9.3.0.713579 (R2017b)'

>>
>> D = full (diag ([0 0 0]))


D =

0 0 0
0 0 0
0 0 0

>> pinv (D)


ans =

0 0 0
0 0 0
0 0 0

>>

Michael Godfrey <godfrey>
Project Member
Fri 03 Nov 2017 03:29:27 PM UTC, comment #14:

@Michael: Thanks. I've got a cset that addresses this issue. I know it is a corner case, but I'm actually not pleased that Octave has been returning the wrong answer for zero vectors and scalars for many years.

A final test to run through Matlab for diagonal matrices

In Octave, this creates a special Diagonal matrix and the result is

But if I use full matrices, I get NaN

So what does Matlab do?

Rik <rik5>
Project Administrator
Fri 03 Nov 2017 10:36:29 AM UTC, comment #13:

>> pinv (complex (0,0))

pinv ([complex(0,0), complex(0,0), complex(0,0)])

ans =

0

ans =

0
0
0

>> version


ans =

'9.3.0.713579 (R2017b)'

>>

Michael Godfrey <godfrey>
Project Member
Fri 03 Nov 2017 02:46:19 AM UTC, comment #12:

I've been looking at resolving the issues with inv(0) and am now re-visiting this bug. What does Matlab return for

Rik <rik5>
Project Administrator
Tue 20 Jun 2017 09:37:18 AM UTC, comment #11:

So Matlab and Octave agree on the 2x2 matrices. Is the patch of comment #9 ready to be pushed or are there other issues to take care of?

Guillaume <gyom>
Sat 17 Jun 2017 10:11:53 AM UTC, comment #10:
Marco Caliari <caliari>
Project Member
Fri 16 Jun 2017 06:23:08 PM UTC, comment #9:

Attached is a quick hack (inverse.diff) which applies only to double matrices. It explicitly checks for the case of a scalar which is either 0 or Inf and returns the correct output value.

What does Matlab return for these two cases

Rik <rik5>
Project Administrator
Fri 16 Jun 2017 05:36:24 PM UTC, comment #8:

Any fix also needs to respect bug #49690 where the inverse of Inf should be 0 for Matlab compatibility. My guess is that is where the problem in inverse for the input 0 has crept in. See cset 22154.

Rik <rik5>
Project Administrator
Fri 16 Jun 2017 05:00:48 PM UTC, comment #7:

The code for the pseudo_inverse of a double matrix is in liboctave/array/dMatrix.cc. It looks like it calculates the svd, and then also calculates an inverse. My guess is the inverse code is where the problem creeps in.

Rik <rik5>
Project Administrator
Fri 16 Jun 2017 04:47:32 PM UTC, comment #6:

The pseudoinverse of 0 is 0 (wikipedia, for instance). And the pseudoinverse of a null vector is the transposed null vector (this answers comment #5).
inv(0) was Inf on 4.2.0 and this is correct.

Marco Caliari <caliari>
Project Member
Fri 16 Jun 2017 04:39:21 PM UTC, comment #5:

pinv is a weird function.

What does Matlab do for this test case?

Rik <rik5>
Project Administrator
Fri 16 Jun 2017 04:25:30 PM UTC, comment #4:

I would still open a discussion on Octave Maintainers list just to see if there is a good mathematical reason for returning either 0 or Inf.

And from the results of your test it looks like there are two Matlab incompatibilities.

Rik <rik5>
Project Administrator
Fri 16 Jun 2017 04:09:37 PM UTC, comment #3:

In Matlab:

I agree this is a design choice but I would say that 100% compatible with Matlab is required for functions as fundamental as inv & pinv.

Guillaume <gyom>
Fri 16 Jun 2017 03:50:23 PM UTC, comment #2:

Leaving aside what Matlab does for the moment, what would you like to see for the inverse of a scalar?

In general, the inverse would just be 1/x for a scalar. So,

If I follow this analogy, pinv(0) should return 1/0

which is what the function seems to be doing.

Of course, maybe it is better to be Matlab compatible for this case. Or maybe pinv should follow inv.

There is definitely room for interpretation here, rather than being absolutely the wrong answer. Maybe you should post on the Octave Maintainer's list to gather some other opinions.

Rik <rik5>
Project Administrator
Fri 16 Jun 2017 02:20:48 PM UTC, comment #1:

I confirm and changed to regression, since it was 0 on 4.0.2. Moreover

is correct on dev.

Marco Caliari <caliari>
Project Member
Fri 16 Jun 2017 12:26:04 PM UTC, original submission:

pinv(0) returns Inf in Octave and 0 in Matlab.

Guillaume <gyom>

 

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

Attach Files:
   
   
Comment:
   

Attached Files
file #40928:  inverse.diff added by rik5 (1KiB - text/x-patch)

 

Depends on the following items: None found

Items that depend on this one: None found

 

Carbon-Copy List
  • -email is unavailable- added by teerav42 (Posted a comment)
  • -email is unavailable- added by godfrey (Posted a comment)
  • -email is unavailable- added by rik5 (Posted a comment)
  • -email is unavailable- added by caliari (Posted a comment)
  • -email is unavailable- added by gyom (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 6 latest changes.

    Date Changed by Updated Field Previous Value => Replaced by
    2017-11-03 rik5 StatusIn Progress => Fixed
        Open/ClosedOpen => Closed
    2017-06-16 rik5 Attached File- => Added inverse.diff, #40928
    2017-06-16 rik5 StatusNone => In Progress
        Summarypinv(0) => pinv(0) different from Matlab
    2017-06-16 caliari Item GroupMatlab Compatibility => Regression

    Back to the top


    Powered by Savane 3.4