bugGNU Octave - Bugs: bug #58944, [z,r,p]=qr(a,b) unexpectedly...

 
 

bug #58944: [z,r,p]=qr(a,b) unexpectedly returns z as q instead of q' * b

Submitter:  A.R. Burgers <arb>
Submitted:  Thu 13 Aug 2020 02:35:13 PM UTC
   
 
Category:  Octave Function Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Incorrect Result
Status:  In Progress Assigned to:  None
Originator Name:  Open/Closed:  * Open
Release:  * dev Operating System:  * Any
Fixed Release:  None Planned Release:  None
* Mandatory Fields

Add a New Comment Rich Markup
   

Jump to the original submission

Wed 27 Jan 2021 04:57:17 PM UTC, comment #14: 

@Rik: Documenting the case in comment #2 is still missing on the default branch.
It would be good if someone who actually knows what this function does (i.e., someone other than me) could update the documentation.

Markus Mützel <mmuetzel>
Group administrator
Wed 27 Jan 2021 04:11:54 PM UTC, comment #13: 

Are we done with this bug report?  Okay to close?

Rik <rik5>
Group administrator
Tue 26 Jan 2021 03:47:15 PM UTC, comment #12: 

@Rik: Thanks for clearing this up for me.
The input check was missing in Octave 6, too. So, I pushed this change to stable:
https://hg.savannah.gnu.org/hgweb/octave/rev/4d3cfe193f65

The other part - [C, R, E] = qr (A, B) with sparse input - only works on default.
It would be nice if someone who knows what they are writing about could update the documentation.

Markus Mützel <mmuetzel>
Group administrator
Tue 26 Jan 2021 12:28:58 AM UTC, comment #11: 

@Markus: Sorry, I think I misled you.  The 3-output argument syntax is valid for dense matrices as long as there is just a single matrix A.  This is what polyfit.m is relying on and, as I said, it is valid.  The extra argument '0' in


[q,r,p] = qr (a, 0);


is not a matrix B, but a flag asking for an "economy"-sized return value.

So, I think the input validation needs to look like


if (have_b && nargout > 2)
  error (...)


Rik <rik5>
Group administrator
Mon 25 Jan 2021 12:33:11 PM UTC, comment #10: 

The attached patch tightens the input validation and removes the BISTs that would fail now.

However, `test polyfit` fails with that change. It seems to rely on that (undocumented) extension...

(file #50782)

Markus Mützel <mmuetzel>
Group administrator
Mon 25 Jan 2021 02:04:04 AM UTC, comment #9: 

@Markus: Those tests can be removed when the input validation is made tighter.  It will be very obvious because those tests will fail during 'make check' immediately.

Rik <rik5>
Group administrator
Sun 24 Jan 2021 04:44:36 PM UTC, comment #8: 

Typo in the markup. Repeating the previous post:

There are currently several BISTs like the following ones in qr.cc where `a`
is dense:

%! [q,r,p] = qr (a, 0); t(j++) = __testqr (q, r, a, p);
%! [q,r,p] = qr (a',0); t(j++) = __testqr (q, r, a', p);


So is this an un-documented extension in Octave? Or should we also remove
these tests?

Markus Mützel <mmuetzel>
Group administrator
Sun 24 Jan 2021 04:43:38 PM UTC, comment #7: 

There are currently several BISTs like the following ones in qr.cc where `a` is dense:

%! [q,r,p] = qr (a, 0);  t(j++) = __testqr (q, r, a, p);
%! [q,r,p] = qr (a',0);  t(j++) = __testqr (q, r, a', p);
-vebatim-

So is this an un-documented extension in Octave? Or should we also remove these tests?

Markus Mützel <mmuetzel>
Group administrator
Thu 14 Jan 2021 05:08:17 PM UTC, comment #6: 

Yes, let's keep the bug open.  We can close it when the input validation is corrected to reject


[C, R, E] = qr (A, B)


when A, B are not sparse.  And yes, the documentation should be updated to mention this form.

Rik <rik5>
Group administrator
Thu 14 Jan 2021 07:30:25 AM UTC, comment #5: 

Ah. Thank you. I misunderstood (and didn't find that explained in the help text).

>> A = sparse(rand(3,2));
>> B = sparse(rand(3,1));
>> [C, R, E] = qr (A, B);
>> X = E*(R\C)
X =

  -0.1954
   0.4075

>> A\B
ans =

Compressed Column Sparse (rows = 2, cols = 1, nnz = 2 [100%])

  (1, 1) -> -0.1954
  (2, 1) -> 0.4075


That looks good to me.

Should we close this as fixed? Or keep it open to document that we now support this syntax for `qr`?

What about the syntax with two full matrices as input and three output arguments (which is also undocumented)? What does it return? Should that be rejected?

Markus Mützel <mmuetzel>
Group administrator
Thu 14 Jan 2021 07:19:13 AM UTC, comment #4: 

There is no "solution" of A*X=B, only a "least square solution". Can you please compare A\B directly with X?

Marco Caliari <caliari>
Group Member
Wed 13 Jan 2021 02:40:36 PM UTC, comment #3: 

I tried the following with Octave at hg id 512591ccb174 (i.e., after the changes for bug #57033 ):

>> A = rand(3,2);
>> B = rand(3,1);
>> [C, R, E] = qr (sparse(A), sparse(B));
>> X = E*(R\C);
>> A*X
ans =

   0.1617
   0.2306
   0.2725

>> B
B =

   0.802331
   0.021141
   0.069525


It looks like the identity described in comment #2 doesn't hold with our current implementation.

Is that a problem?

Markus Mützel <mmuetzel>
Group administrator
Fri 21 Aug 2020 10:27:40 AM UTC, comment #2: 

I think that


[C, R, E] = qr (A, B)


should work only for sparse matrices and returns the fill-reducing ordering matrix E, such that A*X=B iff X=E*(R\C).
It is related to bug #57033, especially comments #6 and #9.

Marco Caliari <caliari>
Group Member
Thu 13 Aug 2020 09:38:32 PM UTC, comment #1: 

This may not be implemented in Octave.  Looking at the documentation, the calling forms are:


 -- [Q, R] = qr (A)
 -- [Q, R, P] = qr (A) # non-sparse A
 -- X = qr (A) # non-sparse A
 -- R = qr (A) # sparse A
 -- [C, R] = qr (A, B)
 -- [...] = qr (..., 0)
 -- [...] = qr (..., "vector")
 -- [...] = qr (..., "matrix")


There is a three-output form for a single input matrix A, but there is no listing for a three-output form with 2 inputs.  In the main body of the documentation there is also no discussion of a form such as


[C, R, P] = qr (A, B)


It gets weirder.  Matlab won't accept


[z,r,p] = qr (a,b);


and complains that the matrices must be sparse.  If I do that with


[z,r,p] = qr (sparse(a),sparse(b));


then it works and returns a 3x1 vector.

However, Octave goes the other way.  it will accept full matrices with a permutation output, but not sparse ones:


[z,r,p] = qr(sparse(a),sparse(b));
error: qr: Permutation output is not supported for sparse input


I think we need to decide what should happen and then go implement that.  Should the [z,r,p] form be allowed with full matrices?  Should the [z,r,p] = qr (a,b) form work with sparse matrices?

Rik <rik5>
Group administrator
Thu 13 Aug 2020 02:35:13 PM UTC, original submission:  

according to the docs [z,r] = qr(a,b) returns the vector z=q'*b, where q is the matrix returned by [q,r]=qr(a). This works as advertized, however [z,r,p] = qr(a,b) unexpectedly (for me) returns z as the matrix q. This applies to all of octave 5/6/7


ver = version
a = rand(3,2);
b = rand(3,1);
[z,r] = qr(a,b);
size_z_no_p_requested = size(z)
[z,r,p] = qr(a,b);
size_z_p_requested = size(z)



ver = 5.2.1
size_z_no_p_requested =
   3   1

size_z_p_requested =
   3   3


A.R. Burgers <arb>

 

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

Attach Files:
   
   
Comment:
   

Attached Files

 

Depends on the following items: None found

Items that depend on this one: None found

 

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

    Date Changed by Updated Field Previous Value => Replaced by
    2021-01-25 mmuetzel Attached File- Added bug58944_qr_input_output_validation.patch, #50782
    2021-01-14 rik5 StatusNeed Info In Progress
    2020-08-13 rik5 StatusNone Need Info

    Back to the top

    Powered by Savane 3.13-3230.
    Corresponding source code