bugGNU Octave - Bugs: bug #32966, spdiags is slow and buggy

 
 

bug #32966: spdiags is slow and buggy

Submitter:  Marco Caliari <caliari>
Submitted:  Thu 31 Mar 2011 01:21:38 PM UTC
   
 
Category:  Octave Function Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Incorrect Result
Status:  Fixed Assigned to:  None
Originator Name:  Open/Closed:  * Closed
Release:  * dev Operating System:  * GNU/Linux
Fixed Release:  None Planned Release:  None
* Mandatory Fields

Add a New Comment Rich Markup
   

Jump to the original submission

Thu 19 Dec 2013 12:16:20 AM UTC, comment #6: 

Marco,

I rewrote spdiags extraction routine to use diag() in C++.  It is about 20% faster than the previous code.  I also overhauled the function (http://hg.savannah.gnu.org/hgweb/octave/rev/190ef1764d30).  Closing report.

Rik <rik5>
Group administrator
Wed 18 Dec 2013 04:21:23 PM UTC, comment #5: 

Marco,

I applied one of your patches that changes diag() to return an empty matrix rather than issue a warning when the requested diagonal is out of range (http://hg.savannah.gnu.org/hgweb/octave/rev/a86d608c413c).  I didn't apply the second patch which speeds up spdiags by using diag() internally because it fails on the first example given in the Matlab documentation (http://www.mathworks.com/help/matlab/ref/spdiags.html).

--Rik

Rik <rik5>
Group administrator
Thu 14 Mar 2013 07:08:21 AM UTC, comment #4: 

I'm sorry Rik, but I won't be able to contribute to Octave for the next few months. I'm looking forward to start again.

Thanks,

Marco

Marco Caliari <caliari>
Group Member
Tue 12 Mar 2013 11:19:16 PM UTC, comment #3: 

Marco, this bug report is 1.5 years old now.  Do you want to rebase your patches onto a recent Mercurial changeset and close this report?

Rik <rik5>
Group administrator
Fri 30 Sep 2011 04:06:19 PM UTC, comment #2: 

Marco,

This wasn't the case when you filed this bug report, but you are part of the Octave team now and have commit privileges to Mercurial.  Do you want to just apply the patches in this bug report yourself and close it?  The performance gains seem impressive.

Rik <rik5>
Group administrator
Fri 08 Apr 2011 01:54:26 PM UTC, comment #1: 

I'm submitting two patches: diag.patch fixes the compatibility issue 3). spdiags.m.patch rewrites the function calls spdiags (A) and spdiags(A, C) and uses diag to extract diagonals. This fixes issues 1) and 2). Infact, using the testfile testspdiags.m (which creates sparse matrices of random sizes and extract random diagonals), I see

+verbatim-

> testspdiags(10) % random sparse matrices up to size 10

warning: spdiagsorig failed
...
warning: spdiagsorig failed
cpuorigtot =  10.013
cputot =  0.85205
-verbatim-


> testspdiags(100) % random sparse matrices up to size 100
warning: spdiagsorig failed
...
warning: spdiagsorig failed
cpuorigtot =  11.421
cputot =  5.7804



> testspdiags(1000) % random sparse matrices up to size 1000
cpuorigtot =  326.95
cputot =  32.586


Because of issue 2) in the original spdiags, sometimes you get the warning. Two compatibility issues remain: number 4) (which I think is a bug in ml). Then, if A is a vector, ml gives an error, whereas the original spdiags in Octave (and mine) works. When you call spdiags with 1 or 2 input args you want to extract diagonals (no ambiguity as in diag): so a vector is just a rectangular matrix.

Marco

(file #23138, file #23139, file #23140)

Marco Caliari <caliari>
Group Member
Thu 31 Mar 2011 01:21:38 PM UTC, original submission:  

Dear maintainers,

stimulated by the comment

## FIXME: we can do this without a loop if we are clever

in spdiags.m, I found the following:

1) diag is much faster then spdiags in extracting diagonals. I get

> v=sprand(1000,1000,0.5);
> tic,spdiags(v,0);,toc
Elapsed time is 0.072 seconds.
> tic,diag(v,0);,toc
Elapsed time is 0.0021 seconds.


2) spdiags fails in the following case:

> spdiags([0,1;0,0],2)
error: spdiags: A(I,J,...) = X: dimensions mismatch

whereas it works with [1,1;1,1].

3) diag is not compatible:

> diag(ones(2),2)
error: diag: requested diagonal out of range

whereas Matlab gives empty matrix 0-by-1.

If 3) is fixed, I can easily rewrite spdiags.m to avoid problems 2) and 1) (without testing if the diagonal is out of range). However, I'm not clever and I cannot do it without a loop :-(

4) there is another incompatible behaviour:

> spdiags(ones(2),-1)
ans =

   1
   0
> spdiags(ones(2,3),-1)
ans =

   1
   0

whereas in Matlab

>> spdiags(ones(2),-1)
ans =
     1
     0
>> spdiags(ones(2,3),-1)
ans =
     0
     1

I do not understand this Matlab (7.6.0) behaviour.

If we take a decision about 4) and someone fixes 3), I can work on 1) and 2).

Cheers,

Marco

Marco Caliari <caliari>
Group Member

 

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

Attach Files:
   
   
Comment:
   

Attached Files
file #23138:  diag.patch added by caliari (1KiB - text/x-patch)
file #23139:  spdiags.m.patch added by caliari (2KiB - text/x-patch)
file #23140:  testspdiags.m added by caliari (612B - text/x-objcsrc)

 

Depends on the following items: None found

Items that depend on this one: None found

 

Carbon-Copy List
  • -email is unavailable- added by mtmiller (Updated the item)
  • -email is unavailable- added by rik5 (Updated the item)
  • -email is unavailable- added by caliari (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 9 latest changes.

    Date Changed by Updated Field Previous Value => Replaced by
    2013-12-19 rik5 StatusPatch Submitted Fixed
        Open/ClosedOpen Closed
    2013-12-18 rik5 Summaryspdiags is slow and buggy, diag is not compatible spdiags is slow and buggy
    2013-10-28 mtmiller CategoryNone Octave Function
        Release3.4.0 dev
    2011-08-04 rik5 StatusNone Patch Submitted
    2011-04-08 caliari Attached File- Added diag.patch, #23138
        Attached File- Added spdiags.m.patch, #23139
        Attached File- Added testspdiags.m, #23140

    Back to the top

    Powered by Savane 3.13-d3ae.
    Corresponding source code