1 // This file is part of Eigen, a lightweight C++ template library 2 // for linear algebra. 3 // 4 // Copyright (C) 2012 Desire Nuentsa Wakam <desire.nuentsa_wakam (at) inria.fr> 5 // 6 // This Source Code Form is subject to the terms of the Mozilla 7 // Public License v. 2.0. If a copy of the MPL was not distributed 8 9 #define EIGEN_NO_DEBUG_SMALL_PRODUCT_BLOCKS 10 #include "sparse.h" 11 #include <Eigen/SPQRSupport> 12 13 14 template<typename MatrixType,typename DenseMat> 15 int generate_sparse_rectangular_problem(MatrixType& A, DenseMat& dA, int maxRows = 300, int maxCols = 300) 16 { 17 eigen_assert(maxRows >= maxCols); 18 typedef typename MatrixType::Scalar Scalar; 19 int rows = internal::random<int>(1,maxRows); 20 int cols = internal::random<int>(1,rows); 21 double density = (std::max)(8./(rows*cols), 0.01); 22 23 A.resize(rows,cols); 24 dA.resize(rows,cols); 25 initSparse<Scalar>(density, dA, A,ForceNonZeroDiag); 26 A.makeCompressed(); 27 return rows; 28 } 29 30 template<typename Scalar> void test_spqr_scalar() 31 { 32 typedef SparseMatrix<Scalar,ColMajor> MatrixType; 33 MatrixType A; 34 Matrix<Scalar,Dynamic,Dynamic> dA; 35 typedef Matrix<Scalar,Dynamic,1> DenseVector; 36 DenseVector refX,x,b; 37 SPQR<MatrixType> solver; 38 generate_sparse_rectangular_problem(A,dA); 39 40 Index m = A.rows(); 41 b = DenseVector::Random(m); 42 solver.compute(A); 43 if (solver.info() != Success) 44 { 45 std::cerr << "sparse QR factorization failed\n"; 46 exit(0); 47 return; 48 } 49 x = solver.solve(b); 50 if (solver.info() != Success) 51 { 52 std::cerr << "sparse QR factorization failed\n"; 53 exit(0); 54 return; 55 } 56 //Compare with a dense solver 57 refX = dA.colPivHouseholderQr().solve(b); 58 VERIFY(x.isApprox(refX,test_precision<Scalar>())); 59 } 60 void test_spqr_support() 61 { 62 CALL_SUBTEST_1(test_spqr_scalar<double>()); 63 CALL_SUBTEST_2(test_spqr_scalar<std::complex<double> >()); 64 } 65
