bugGNU Octave - Bugs: bug #57033, Replace CXSPARSE with SPQR

 
 

bug #57033: Replace CXSPARSE with SPQR

Submitted by:  Rik <rik5>
Submitted on:  Thu 10 Oct 2019 04:57:19 PM UTC  
 
Category:  Libraries Severity:  3 - Normal
Priority:  5 - Normal Item Group:  Feature Request
Status:  Patch Submitted Assigned to:  None
Originator Name:  Open/Closed:  Open
Release:  dev Operating System:  Any

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 18 Oct 2019 09:16:30 PM UTC, comment #6: 

Hello,
here are some additional notes on the test script.

Let A be a sparse m x n matrix. The test script is divided into four sections.
In the first section the QR-factorization [Q,R,E] = qr(A) is calculated. E is an array of size n that represents a (fill-reducing) column permutation of A. The script creates a permutation matrix P from E (with Q*R = A*P) and checks whether Q*R*P' - A is small. 
In the second section some overdetermined systems A*x = b (m > n) are solved by calling x = qr(A,b). The script prints norm (A*x - b,2). Usually the values of norm(A*x - b,2) are between 1.0 and 6.0. The n columns of A span a subspace of the m-dimensional vector space. If b lies outside this subspace, there is no x with A * x = b. Only an approximation A*x of b can be calculated.
In the third section some underdetermined systems A*x = b (m < n) are solved by calling x = qr(A,b). The script prints norm(A*x - b,2). Usually the values of norm(A*x - b,2) have the order of magnitude -16. If m < n, there is a high chance that a column-subset of A spans the complete m-dimensional vector space. In this case, there exists x with A*x = b and norm(A*x - b,2) is very small. If no such subset exists, the value of norm (A*x - b,2) is large (> 0.1).
The fourth section checks whether the QR-factorization for a dense A still works.

Simon

Simon <sim79>
Sun 13 Oct 2019 03:14:43 PM UTC, comment #5: 

Thanks, Rik.

The attached patch adds a test for -lspqr to configure. However, I don't know how to check for a function from a template. Instead I checked for the C function "SuiteSparseQR_C" which happens to live in the same library.
Also some more logic could probably be added to make sure that we don't link to CXSparse in case SPQR can be used. (If this is what we want to do.)
There is a comment in configure.ac that the order of the suitesparse flags matters on some systems. The current order works on my system (Ubuntu). I don't know what the "correct" order is though.

The patch applies on top of the patch from comment #0. Feel free to merge the patches if you want to.

During compilation, I see a few shadowed declaration warnings:

../liboctave/numeric/sparse-qr.cc: In member function ‘SPARSE_T octave::math::sparse_qr<SPARSE_T>::sparse_qr_rep::R(bool) const [with SPARSE_T = SparseMatrix]’:
../liboctave/numeric/sparse-qr.cc:371:59: warning: declaration of ‘econ’ shadows a member of ‘octave::math::sparse_qr<SparseMatrix>::sparse_qr_rep’ [-Wshadow]
     sparse_qr<SparseMatrix>::sparse_qr_rep::R (bool econ) const
                                                           ^~~~~
../liboctave/numeric/sparse-qr.cc:130:23: note: shadowed declaration is here
       octave_idx_type econ;
                       ^~~~
../liboctave/numeric/sparse-qr.cc: In member function ‘RET_T octave::math::sparse_qr<SPARSE_T>::sparse_qr_rep::wide_solve(const RHS_T&, octave_idx_type&) const [with RHS_T = SparseMatrix; RET_T = SparseMatrix; SPARSE_T = SparseMatrix; octave_idx_type = long int]’:
../liboctave/numeric/sparse-qr.cc:1059:23: warning: declaration of ‘S’ shadows a member of ‘octave::math::sparse_qr<SparseMatrix>::sparse_qr_rep’ [-Wshadow]
       cholmod_sparse *S;
                       ^
../liboctave/numeric/sparse-qr.cc:132:57: note: shadowed declaration is here
       typename cxsparse_types<SPARSE_T>::symbolic_type *S;
                                                         ^
../liboctave/numeric/sparse-qr.cc: In member function ‘SPARSE_T octave::math::sparse_qr<SPARSE_T>::sparse_qr_rep::R(bool) const [with SPARSE_T = SparseComplexMatrix]’:
../liboctave/numeric/sparse-qr.cc:1542:66: warning: declaration of ‘econ’ shadows a member of ‘octave::math::sparse_qr<SparseComplexMatrix>::sparse_qr_rep’ [-Wshadow]
     sparse_qr<SparseComplexMatrix>::sparse_qr_rep::R (bool econ) const
                                                                  ^~~~~
../liboctave/numeric/sparse-qr.cc:130:23: note: shadowed declaration is here
       octave_idx_type econ;
                       ^~~~

The tests from the .zip file also produce some results and don't crash for me. I'm also not sure how to interpret the results though.

I'll try to cross-build for Windows next.

(file #47661)

Markus Mützel <mmuetzel>
Project Member
Sat 12 Oct 2019 11:48:21 PM UTC, comment #4: 

So I can get it to compile by adding "-lspqr" to the environment variable LDFLAGS and then re-configuring and compiling.

This suggests only that some more work has to be done with the configuration and build system to check for the existence of libspqr during configuration, and then to add "-lspqr" to the LIBS variable.

With the compiled version I was able to run test_spqr.m in the zip file.  I don't know what output it is supposed to produce, but it doesn't segfault or do anything bad.

Rik <rik5>
Project Administrator
Sat 12 Oct 2019 11:18:47 PM UTC, comment #3: 

I also fail to compile on Linux.  This is with suitesparse-dev-5.1.2.  It is failing at the link stage.  It may be that the SQPR routines are in a different package such as libspqr2 (but I have that installed), or that some different link flags need to be set.

Rik <rik5>
Project Administrator
Sat 12 Oct 2019 10:47:42 AM UTC, comment #2: 

I tried to compile with this patch (to eventually try with a cross-build on Windows). However, the linker fails to resolve some references:

  CXXLD    src/octave-cli
/usr/bin/ld: liboctave/.libs/liboctave.so: undefined reference to `cholmod_dense_struct* SuiteSparseQR_qmult<double>(int, cholmod_sparse_struct*, cholmod_dense_struct*, long*, cholmod_dense_struct*, cholmod_common_struct*)'
/usr/bin/ld: liboctave/.libs/liboctave.so: undefined reference to `cholmod_sparse_struct* SuiteSparseQR_min2norm<double>(int, double, cholmod_sparse_struct*, cholmod_sparse_struct*, cholmod_common_struct*)'
/usr/bin/ld: liboctave/.libs/liboctave.so: undefined reference to `long SuiteSparseQR<double>(int, double, long, cholmod_sparse_struct*, cholmod_sparse_struct**, long**, cholmod_sparse_struct**, long**, cholmod_dense_struct**, cholmod_common_struct*)'
/usr/bin/ld: liboctave/.libs/liboctave.so: undefined reference to `long SuiteSparseQR<std::complex<double> >(int, double, long, cholmod_sparse_struct*, cholmod_sparse_struct**, long**, cholmod_sparse_struct**, long**, cholmod_dense_struct**, cholmod_common_struct*)'
/usr/bin/ld: liboctave/.libs/liboctave.so: undefined reference to `cholmod_sparse_struct* SuiteSparseQR_qmult<double>(int, cholmod_sparse_struct*, cholmod_dense_struct*, long*, cholmod_sparse_struct*, cholmod_common_struct*)'
/usr/bin/ld: liboctave/.libs/liboctave.so: undefined reference to `cholmod_dense_struct* SuiteSparseQR_qmult<std::complex<double> >(int, cholmod_sparse_struct*, cholmod_dense_struct*, long*, cholmod_dense_struct*, cholmod_common_struct*)'
/usr/bin/ld: liboctave/.libs/liboctave.so: undefined reference to `cholmod_dense_struct* SuiteSparseQR_min2norm<std::complex<double> >(int, double, cholmod_sparse_struct*, cholmod_dense_struct*, cholmod_common_struct*)'
/usr/bin/ld: liboctave/.libs/liboctave.so: undefined reference to `cholmod_sparse_struct* SuiteSparseQR_min2norm<std::complex<double> >(int, double, cholmod_sparse_struct*, cholmod_sparse_struct*, cholmod_common_struct*)'
/usr/bin/ld: liboctave/.libs/liboctave.so: undefined reference to `cholmod_dense_struct* SuiteSparseQR_min2norm<double>(int, double, cholmod_sparse_struct*, cholmod_dense_struct*, cholmod_common_struct*)'
collect2: error: ld returned 1 exit status
make[2]: *** [Makefile:15663: src/octave-cli] Error 1
make[2]: Leaving directory '/home/osboxes/Repositories/Octave/octave-1/.build'
make[1]: *** [Makefile:27400: all-recursive] Error 1
make[1]: Leaving directory '/home/osboxes/Repositories/Octave/octave-1/.build'
make: *** [Makefile:11032: all] Error 2

My steps were:
1) Update octave repo to hg id 3fe26656e73c
2) Apply the patch from comment #0.
3) Run "make maintainer-clean"
4) Complete "bootstrap"-"configure"-"make all" cycle.

This was on Ubuntu 19.04 with libsuitesparse-dev 5.4.0.
I understand that this is an incomplete patch. Please, let me know if it is too early to try it on a cross-build for Windows.
Otherwise: What do I have to do to get this to compile and link (on Ubuntu-Linux as a first step)?

Markus Mützel <mmuetzel>
Project Member
Thu 10 Oct 2019 04:58:44 PM UTC, comment #1: 

I filed this issue report to keep track of a promising patch to implement SPQR for sparse operations in Octave.  I've added the author of the patch (Simon) to the CC list for this bug, as well as Marco and Kai who I believe are interested in this topic.

Rik <rik5>
Project Administrator
Thu 10 Oct 2019 04:57:19 PM UTC, original submission:  

I appended a test patch and a zip file. The zip file contains the script test_spqr.m with test cases.  Two test cases use SuiteSparse test matrices in Matrix Market Format. So I added the test matrices and MatrixMarketUtils to the zip file. I hope the patch will work on your system. 

The current status of the patch is
1) QR factorization can be calculated for sparse (real|complex) matrices.
2) least square and min2norm solutions can be calculated for sparse (real|complex) matrices.

Remarks
1) I will improve code formatting. I didn't apply the GNU coding standards consistently.
2) Probably the error handling and validation of parameters can be improved.
3) The solve methods for least square and min2norm are moved from class sparse_qr_rep to class sparse_qr. I can change this if necessary.
4) Macros HAVE_SUITESPARSE_SPQR_H and HAVE_SPQR are defined in oct-sparse.h. I don't know how to put them into config.h.
5) I didn't do any tests on Windows or Mac.
6) SPQR stores the results of the QR factorization either in a collection of data structures (H, R, HPinv, HTau) or in one "big" structure of type SuiteSparseQR_factorization. Unfortunately, there is no user-callable SPQR function to extract H, R, HPinv and HTau from a SuiteSparseQR_factorization-object. At least I didn't find one.

The implementation of the R-,C-,Q- and V-methods in class sparse_qr_rep is based on H, R, HPinv and HTau (calculated in the constructor). For the implementation of the solve methods the function SuiteSparseQR_min2norm is very convenient. This function uses the SuiteSparseQR_factorization and does not reuse H, R, HPinv and HTau. So for a matrix A I do one QR factorization in order to get R, Q etc. and a further factorization to solve Ax = b.

In the old implementation of sparse_qr_rep there was only one QR factorization in the constructor. The solve methods and the calculation of R, C, Q and V based on this single factorization. Therefore I think that the old approach is more efficient than my implementation.

I am thankful for tips or comments how to proceed.

Simon

Rik <rik5>
Project Administrator

 

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

Attach Files:
   
   
Comment:
   

Attached Files
file #47661:  bug57033_spqr_configure.patch added by mmuetzel (4KiB - application/octet-stream)
file #47640:  test_spqr_20191006.patch added by rik5 (61KiB - text/x-patch)
file #47641:  test_spqr_20191006.zip added by rik5 (1MiB - application/zip)

 

Depends on the following items: None found

Digest:
   bug dependencies.

 

Carbon-Copy List
  • -email is unavailable- added by sim79 (Posted a comment)
  • -email is unavailable- added by mmuetzel (Posted a comment)
  • -email is unavailable- added by rik5 (Submitted the item)
  • -email is unavailable- added by rik5 (New patch to reimplement sparse operations with SPQR)
  • -email is unavailable- added by rik5 (New patch to reimplement sparse operations with SPQR)
  • -email is unavailable- added by rik5 (New patch to reimplement sparse operations with SPQR)
  •  

    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 7 latest changes.

    Date Changed by Updated Field Previous Value => Replaced by
    2019-10-13 mmuetzel Dependencies- => bugs #51950 is dependent
    2019-10-13 mmuetzel Attached File- => Added bug57033_spqr_configure.patch, #47661
    2019-10-10 rik5 Attached File- => Added test_spqr_20191006.patch, #47640
        Attached File- => Added test_spqr_20191006.zip, #47641
        Carbon-Copy- => Added marco.caliari@univr.it
        Carbon-Copy- => Added k.ohlhus@gmail.com
        Carbon-Copy- => Added simon.hau79@gmail.com

    Back to the top


    Powered by Savane 3.5