# How to use the qr decomposition correctly?

3 messages
Open this post in threaded view
|

## How to use the qr decomposition correctly?

 Hi,everybody.    I'm quite new to numerical computing, especially using boost's lapack-bindings. So I hope to find some generous numerical computing guru among you to set me on the right path.         -- Begin code sample -- //qrf.h#ifndef QRF_H #define QRF_H #if defined(_MSC_VER) && (_MSC_VER >= 1020) # pragma once #endif #include void qr_factorize( const boost::numeric::ublas::matrix& A, boost::numeric::ublas::matrix& Q, boost::numeric::ublas::triangular_matrix& R); #endif//qrf.cpp#include "stdafx.h" #include "qrf.h" #include #include #include #include #define BOOST_NUMERIC_BINDINGS_USE_CLAPACK #include #include #include #include #undef BOOST_NUMERIC_BINDINGS_USE_CLAPACK void qr_factorize( const boost::numeric::ublas::matrix& A, boost::numeric::ublas::matrix& Q, boost::numeric::ublas::triangular_matrix& R) { namespace ublas = boost::numeric::ublas; BOOST_UBLAS_CHECK(A.size1() >= A.size2(), ublas::external_logic()); std::vector tau(A.size2()); ublas::matrix CQ(A); int info; info = boost::numeric::bindings::lapack::geqrf(CQ, tau); BOOST_UBLAS_CHECK(info == 0, ublas::internal_logic()); ublas::triangular_matrix CR = ublas::triangular_adaptor >, ublas::upper>( ublas::project(CQ, ublas::range(0, A.size2()), ublas::range(0, A.size2())) ); ublas::identity_matrix I(A.size1()); ublas::matrix tempQ(I); info = boost::numeric::bindings::lapack::ormqr('L', 'N', CQ, tau, tempQ, boost::numeric::bindings::lapack::optimal_workspace()); BOOST_UBLAS_CHECK(info == 0, ublas::internal_logic()); #if BOOST_UBLAS_TYPE_CHECK BOOST_UBLAS_CHECK(ublas::detail::expression_type_check(ublas::prod(CQ, CR), A), ublas::internal_logic()); #endif ublas::matrix CCQ(CQ); Q.assign_temporary(CCQ); R.assign_temporary(CR); }// test_LAPACK_QRF.cpp : Defines the entry point for the console application. // #include "stdafx.h" #include "test_LAPACK_QRF.h" #include "qrf.h" #include #include #include #ifdef _DEBUG #ifdef _MSC_VER # pragma comment(lib, "libf2cd.lib") # pragma comment(lib, "BLASd.lib") # pragma comment(lib, "clapackd.lib") #endif #else #ifdef _MSC_VER # pragma comment(lib, "libf2c.lib") # pragma comment(lib, "BLAS.lib") # pragma comment(lib, "clapack.lib") #endif #endif #ifdef _DEBUG #define new DEBUG_NEW #endif // The one and only application object CWinApp theApp; using namespace std; namespace ublas = boost::numeric::ublas; int _tmain(int argc, TCHAR* argv[], TCHAR* envp[]) { int nRetCode = 0; HMODULE hModule = ::GetModuleHandle(NULL); if (hModule != NULL) { // initialize MFC and print and error on failure if (!AfxWinInit(hModule, NULL, ::GetCommandLine(), 0)) { // TODO: change error code to suit your needs _tprintf(_T("Fatal Error: MFC initialization failed\n")); nRetCode = 1; } else { // TODO: code your application's behavior here. //ublas::matrix m(3, 3); ublas::matrix m(3, 3); ublas::matrix Q; ublas::triangular_matrix R; m(0, 0) = 12; m(0, 1) = -51; m(0, 2) = 4; m(1, 0) = 6; m(1, 1) = 167; m(1, 2) = -68; m(2, 0) = -4; m(2, 1) = 24; m(2, 2) = -41; qr_factorize(m, Q, R); std::cout << Q << std::endl; std::cout << R << std::endl; std::cout << ublas::prod(Q, R) << std::endl; return 0; } } else { // TODO: change error code to suit your needs _tprintf(_T("Fatal Error: GetModuleHandle failed\n")); nRetCode = 1; } return nRetCode; } -- End code sample -- Compiling this code gives an error from lapack: Error 1 error C2589: '(' : illegal token on right side of '::' d:\3rd\include\boost\numeric\bindings\lapack\ormqr.hpp 189  traits::detail::array work( std::max(1,n_w*32) ); raymond Chen Shanghai[hidden email] _______________________________________________ ublas mailing list [hidden email] http://lists.boost.org/mailman/listinfo.cgi/ublasSent to: [hidden email]