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@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 #include "sparse.h"
9 #include <Eigen/SPQRSupport>
10
11
12 template<typename MatrixType,typename DenseMat>
generate_sparse_rectangular_problem(MatrixType & A,DenseMat & dA,int maxRows=300,int maxCols=300)13 int generate_sparse_rectangular_problem(MatrixType& A, DenseMat& dA, int maxRows = 300, int maxCols = 300)
14 {
15 eigen_assert(maxRows >= maxCols);
16 typedef typename MatrixType::Scalar Scalar;
17 int rows = internal::random<int>(1,maxRows);
18 int cols = internal::random<int>(1,rows);
19 double density = (std::max)(8./(rows*cols), 0.01);
20
21 A.resize(rows,rows);
22 dA.resize(rows,rows);
23 initSparse<Scalar>(density, dA, A,ForceNonZeroDiag);
24 A.makeCompressed();
25 return rows;
26 }
27
test_spqr_scalar()28 template<typename Scalar> void test_spqr_scalar()
29 {
30 typedef SparseMatrix<Scalar,ColMajor> MatrixType;
31 MatrixType A;
32 Matrix<Scalar,Dynamic,Dynamic> dA;
33 typedef Matrix<Scalar,Dynamic,1> DenseVector;
34 DenseVector refX,x,b;
35 SPQR<MatrixType> solver;
36 generate_sparse_rectangular_problem(A,dA);
37
38 int m = A.rows();
39 b = DenseVector::Random(m);
40 solver.compute(A);
41 if (solver.info() != Success)
42 {
43 std::cerr << "sparse QR factorization failed\n";
44 exit(0);
45 return;
46 }
47 x = solver.solve(b);
48 if (solver.info() != Success)
49 {
50 std::cerr << "sparse QR factorization failed\n";
51 exit(0);
52 return;
53 }
54 //Compare with a dense solver
55 refX = dA.colPivHouseholderQr().solve(b);
56 VERIFY(x.isApprox(refX,test_precision<Scalar>()));
57 }
test_spqr_support()58 void test_spqr_support()
59 {
60 CALL_SUBTEST_1(test_spqr_scalar<double>());
61 CALL_SUBTEST_2(test_spqr_scalar<std::complex<double> >());
62 }
63