1 // This file is part of Eigen, a lightweight C++ template library 2 // for linear algebra. 3 // 4 // Copyright (C) 2010 Benoit Jacob <jacob.benoit.1 (at) gmail.com> 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 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 9 10 #include "main.h" 11 #include <Eigen/SVD> 12 13 template<typename MatrixType> void upperbidiag(const MatrixType& m) 14 { 15 const typename MatrixType::Index rows = m.rows(); 16 const typename MatrixType::Index cols = m.cols(); 17 18 typedef Matrix<typename MatrixType::RealScalar, MatrixType::RowsAtCompileTime, MatrixType::ColsAtCompileTime> RealMatrixType; 19 typedef Matrix<typename MatrixType::Scalar, MatrixType::ColsAtCompileTime, MatrixType::RowsAtCompileTime> TransposeMatrixType; 20 21 MatrixType a = MatrixType::Random(rows,cols); 22 internal::UpperBidiagonalization<MatrixType> ubd(a); 23 RealMatrixType b(rows, cols); 24 b.setZero(); 25 b.block(0,0,cols,cols) = ubd.bidiagonal(); 26 MatrixType c = ubd.householderU() * b * ubd.householderV().adjoint(); 27 VERIFY_IS_APPROX(a,c); 28 TransposeMatrixType d = ubd.householderV() * b.adjoint() * ubd.householderU().adjoint(); 29 VERIFY_IS_APPROX(a.adjoint(),d); 30 } 31 32 void test_upperbidiagonalization() 33 { 34 for(int i = 0; i < g_repeat; i++) { 35 CALL_SUBTEST_1( upperbidiag(MatrixXf(3,3)) ); 36 CALL_SUBTEST_2( upperbidiag(MatrixXd(17,12)) ); 37 CALL_SUBTEST_3( upperbidiag(MatrixXcf(20,20)) ); 38 CALL_SUBTEST_4( upperbidiag(Matrix<std::complex<double>,Dynamic,Dynamic,RowMajor>(16,15)) ); 39 CALL_SUBTEST_5( upperbidiag(Matrix<float,6,4>()) ); 40 CALL_SUBTEST_6( upperbidiag(Matrix<float,5,5>()) ); 41 CALL_SUBTEST_7( upperbidiag(Matrix<double,4,3>()) ); 42 } 43 } 44
