bugGNU Octave - Bugs: bug #52593, qr(A,'matrix') and qr(A,'vector')...

 
 

bug #52593: qr(A,'matrix') and qr(A,'vector') return the wrong result with a tall skinny matrix

Submitter:  Joseph Young <josyoun>
Submitted:  Wed 06 Dec 2017 12:01:44 AM UTC
   
 
Category:  Octave Function Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Matlab Compatibility
Status:  Fixed Assigned to:  None
Originator Name:  josyoun Open/Closed:  * Closed
Release:  * 4.2.1 Operating System:  * GNU/Linux
Fixed Release:  None Planned Release:  None
* Mandatory Fields

Add a New Comment Rich Markup
   

Jump to the original submission

Thu 07 Dec 2017 08:42:52 AM UTC, comment #7: 

@Rik: it works for me.
@Joseph: in case you want to stay with the stable branch, you can do someething like


A = [-3.6 0 9.72;3.4 1.9 8.67;3.8 1.7 0;-1.6 4 0;-1.6 4 0];
[Q, R, P] = qr(A);
p = find (P == 1)' - (0:2) * 3


Maybe you already found even a better way.

Marco Caliari <caliari>
Group Member
Thu 07 Dec 2017 12:54:52 AM UTC, comment #6: 

I checked in this change (http://hg.savannah.gnu.org/hgweb/octave/rev/c9d229f1db04) on the development branch of Octave.  I took the time to rewrite the documentation for qr() as well since it was a bit of a confusing mess.

In order to see the changes you will need to build Octave from development sources available through Mercurial.  See the Octave wiki for instructions on getting the source code and building it.

Marking bug as fixed and closing report.

Rik <rik5>
Group administrator
Wed 06 Dec 2017 11:14:56 PM UTC, comment #5: 

One other thing, I think Marco found the right piece of code where the bug is.  Thanks Marco!

Joseph Young <josyoun>
Wed 06 Dec 2017 10:49:44 PM UTC, comment #4: 

Hi, Rik.  Thanks for taking a look at this.  I would contend that it's partially a MATLAB compatibility issue, but the result for  [Q R p]=qr(A,'matrix') is actually incorrect.  Octave returns a permutation vector p when a matrix is explicitly requested.  I would call returning an economy sized factorization for [Q R p]=qr(A,'vector') when MATLAB returns the full factorization a compatibility issue.

Anyway, maybe it's not a big deal in how it's categorized, but I wanted to bring this up in case there's a difference in how and when issues are prioritized internally.  At the end of the day, I think the same bit of code will end up fixing both bugs.

Joseph Young <josyoun>
Wed 06 Dec 2017 06:04:51 PM UTC, comment #3: 

I have the start of a cset for this.  I realized the same thing that Marco did which is that the original coder misunderstood when an economy versus standard QR factorization should be made. 

I've marked this as a Matlab Compatibility issue.

Rik <rik5>
Group administrator
Wed 06 Dec 2017 09:01:46 AM UTC, comment #2: 

The problem starts at line 23812 of libinterp/dldfcn/qr.cc


template <typename T>
static typename octave::math::qr<T>::type
qr_type (int nargin, int nargout)
{
  return ((nargout == 0 || nargout == 1)
          ? octave::math::qr<T>::raw
          : (nargin == 2) ? octave::math::qr<T>::economy : octave::math::qr<T>::std);
}


It is not possible to select the type by just considering the number of input and output arguments. qr(A,0) is type economy, while qr(A,'matrix') and qr(A,'vector') are std. And then there are qr(A,B,0), qr(A,B,'matrix'), qr(A,B,'vector').

Marco Caliari <caliari>
Group Member
Wed 06 Dec 2017 07:34:48 AM UTC, comment #1: 

Confirmed in dev as well.

Marco Caliari <caliari>
Group Member
Wed 06 Dec 2017 12:01:44 AM UTC, original submission:  

It looks like qr(A,'matrix') and qr(A,'vector') return the wrong result when the matrix is tall and skinny.  For example,

> A=[-3.6 0 9.72;3.4 1.9 8.67;3.8 1.7 0;-1.6 4 0;-1.6 4 0];
> [Q R p] = qr(A,'vector')

Q =

   0.74626  -0.49565  -0.14732
   0.66565   0.55568   0.16516
   0.00000   0.57352   0.27055
   0.00000  -0.24148   0.66250
   0.00000  -0.24148   0.66250

R =

   13.02487   -0.42334    1.26473
    0.00000    6.62577    0.09892
    0.00000    0.00000    6.07377

p =

   3   1   2

and

> A=[-3.6 0 9.72;3.4 1.9 8.67;3.8 1.7 0;-1.6 4 0;-1.6 4 0];
> [Q R p] = qr(A,'matrix')

Q =

   0.74626  -0.49565  -0.14732
   0.66565   0.55568   0.16516
   0.00000   0.57352   0.27055
   0.00000  -0.24148   0.66250
   0.00000  -0.24148   0.66250

R =

   13.02487   -0.42334    1.26473
    0.00000    6.62577    0.09892
    0.00000    0.00000    6.07377

p =

   3   1   2

In both cases, the skinny factorization is given, which normally requires the argument 0 to be passed.  In addition, when matrix a matrix permutation is requested, a vector permutation is given.  Conversely, in MATLAB, we have:

>> A = [-3.6 0 9.72;3.4 1.9 8.67;3.8 1.7 0;-1.6 4 0;-1.6 4 0];
>> [Q R p] = qr(A,'vector')                                  


Q =

   -0.7463    0.4957    0.1473   -0.2964   -0.2964
   -0.6656   -0.5557   -0.1652    0.3323    0.3323
         0   -0.5735   -0.2706   -0.5468   -0.5468
         0    0.2415   -0.6625    0.5373   -0.4627
         0    0.2415   -0.6625   -0.4627    0.5373


R =

  -13.0249    0.4233   -1.2647
         0   -6.6258   -0.0989
         0         0   -6.0738
         0         0         0
         0         0         0


p =

     3     1     2

>> [Q R p] = qr(A,'matrix')


Q =

   -0.7463    0.4957    0.1473   -0.2964   -0.2964
   -0.6656   -0.5557   -0.1652    0.3323    0.3323
         0   -0.5735   -0.2706   -0.5468   -0.5468
         0    0.2415   -0.6625    0.5373   -0.4627
         0    0.2415   -0.6625   -0.4627    0.5373


R =

  -13.0249    0.4233   -1.2647
         0   -6.6258   -0.0989
         0         0   -6.0738
         0         0         0
         0         0         0


p =

     0     1     0
     0     0     1
     1     0     0

which I contend should be expected result.  Really, I'm looking for a way to get a vector permutation matrix as well as the entire Q matrix, not the skinny factorization and this doesn't appear to be currently possible due to a bug.

For reference:

GNU Octave Version: 4.2.1
MATLAB Version: 9.3.0.713579 (R2017b)

both where run on Linux.  Thanks for the hard work!

Joseph Young <josyoun>

 

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

Attach Files:
   
   
Comment:
   

No files currently attached

 

Depends on the following items: None found

Items that depend on this one: None found

 

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

    Date Changed by Updated Field Previous Value => Replaced by
    2017-12-07 rik5 StatusConfirmed Fixed
        Open/ClosedOpen Closed
    2017-12-06 rik5 Item GroupNone Matlab Compatibility
    2017-12-06 caliari StatusNone Confirmed

    Back to the top

    Powered by Savane 3.13-02a9.
    Corresponding source code