使用Eigen解Ax=b線性方程組
阿新 • • 發佈:2018-12-18
#include <iostream> using namespace std; #include <ctime> // Eigen 部分 #include <Eigen/Core> // 稠密矩陣的代數運算(逆,特徵值等) #include <Eigen/Dense> #define MATRIX_SIZE 100 /**************************** * 本程式演示了 Eigen 基本型別的使用 ****************************/ int main( int argc, char** argv ) { // 解方程 // 我們求解 A * x = b 這個方程 // 直接求逆自然是最直接的,但是求逆運算量大 Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic > A1; A1 = Eigen::MatrixXd::Random( MATRIX_SIZE, MATRIX_SIZE ); Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic > b1; b1 = Eigen::MatrixXd::Random( MATRIX_SIZE, 1 ); clock_t time_stt = clock(); // 計時 // 直接求逆 Eigen::Matrix<double,MATRIX_SIZE,1> x = A1.inverse()*b1; cout <<"time use in normal inverse is " << 1000* (clock() - time_stt)/(double)CLOCKS_PER_SEC << "ms"<< endl; cout<<x<<endl; // QR分解colPivHouseholderQr() time_stt = clock(); x = A1.colPivHouseholderQr().solve(b1); cout <<"time use in Qr decomposition is " <<1000*(clock() - time_stt)/(double)CLOCKS_PER_SEC <<"ms" << endl; cout <<x<<endl; //QR分解fullPivHouseholderQr() time_stt = clock(); x = A1.fullPivHouseholderQr().solve(b1); cout <<"time use in Qr decomposition is " <<1000*(clock() - time_stt)/(double)CLOCKS_PER_SEC <<"ms" << endl; cout <<x<<endl; /* //llt分解 要求矩陣A正定 time_stt = clock(); x = A1.llt().solve(b1); cout <<"time use in llt decomposition is " <<1000*(clock() - time_stt)/(double)CLOCKS_PER_SEC <<"ms" << endl; cout <<x<<endl;*/ /*//ldlt分解 要求矩陣A正或負半定 time_stt = clock(); x = A1.ldlt().solve(b1); cout <<"time use in ldlt decomposition is " <<1000*(clock() - time_stt)/(double)CLOCKS_PER_SEC <<"ms" << endl; cout <<x<<endl;*/ //lu分解 partialPivLu() time_stt = clock(); x = A1.partialPivLu().solve(b1); cout <<"time use in lu decomposition is " <<1000*(clock() - time_stt)/(double)CLOCKS_PER_SEC <<"ms" << endl; cout <<x<<endl; //lu分解(fullPivLu() time_stt = clock(); x = A1.fullPivLu().solve(b1); cout <<"time use in lu decomposition is " <<1000*(clock() - time_stt)/(double)CLOCKS_PER_SEC <<"ms" << endl; cout <<x<<endl; return 0; }