使用Lapack求解线性代数方程组（C/C++语言）

#include <iostream>
using namespace std;

#include "mkl_cblas.h"
#include "mkl_lapacke.h"
#include "tnt.h"
using namespace TNT;

#ifdef _DLL //动态连接

#pragma comment(lib,"mkl_intel_c_dll.lib")
#pragma comment(lib,"mkl_intel_thread_dll.lib")
#pragma comment(lib,"mkl_core_dll.lib")
#pragma comment(lib,"libiomp5md.lib")

#else //静态连接

#pragma comment(lib,"mkl_intel_c.lib")
#pragma comment(lib,"mkl_intel_thread.lib")
#pragma comment(lib,"mkl_core.lib")
#pragma comment(lib,"libiomp5md.lib")

#endif // _DLL

int main()
{
 Array2D<double> A(4, 4, 0.0);
 Array1D<double> B(4, 0.0);
 Array1D<int> ipiv(4, 0);
 double CondNum(1.0);
 double *rcond = &CondNum;

 A[0][0] = 1.0; A[0][1] = 4.0; A[0][2] = -2.0; A[0][3] = 3.0;
 A[1][0] = 2.0; A[1][1] = 2.0; A[1][2] = 0.0; A[1][3] = 4.0;
 A[2][0] = 3.0; A[2][1] = 0.0; A[2][2] = -1.0; A[2][3] = 2.0;
 A[3][0] = 1.0; A[3][1] = 2.0; A[3][2] = 2.0; A[3][3] = -3.0;

 B[0] = 6.0;
 B[1] = 2.0;
 B[2] = 1.0;
 B[3] = 8.0;

 

 //计算矩阵A条件数：
 LAPACKE_dgecon(LAPACK_ROW_MAJOR,'1',4,A[0],4,1,rcond);
 cout<<"A的条件数（1-范数）："<<CondNum<<endl;

 

 //求解线性方程组(注意要保证A是可逆的)
 int info = LAPACKE_dgesv(LAPACK_ROW_MAJOR, 4, 1, A[0], 4, &ipiv[0], &B[0], 1);
 cout<<"计算结果："<<endl<<B<<" "<<info<<endl;

 

return 0;

}

 

注 1 ： LAPACKE 为LAPACK的C接口；

注 2 ： LAPACK_dgesv参数含义：

LAPACK_ROW_MAJOR： C、C++矩阵存储方式

4：矩阵行数，列数

1：右端向量列数

A[0]：矩阵指针

4：Leading Deminsion of A: The first dimension of A （注意！）

&ipiv[0]:向量指针

1: Leading Deminsion of B: The first dimension of B （注意！）

注 3 ： LAPACK_dgesv 函数返回值为0，表示求解过程正常，否则返回错误类型。