?gesv
求解线性方程组,其中是的矩阵,和为的矩阵。求解过程中,对使用具有部分旋转和行交换的LU分解,将其分解为的形式,其中是置换矩阵,是单位下三角矩阵,是上三角矩阵。并通过LU分解的结果来求解线性方程组。
接口定义
C Interface:
void sgesv_(const int *n,const int *nrhs,float *a,const int *lda,int *ipiv, float *b, const int *ldb,int *info);
void dgesv_(const int *n,const int *nrhs,double *a,const int *lda,int *ipiv, double *b, const int *ldb,int *info);
void cgesv_(const int *n,const int *nrhs,float_Complex *a,const int *lda,int *ipiv, float_Complex *b, const int *ldb,int *info);
void zgesv_(const int *n,const int *nrhs, double_Complex *a,const int *lda,int *ipiv, double_Complex *b, const int *ldb,int *info);
Fortran Interface:
SGESV(n, nrhs, a, lda, ipiv, b, ldb, info);
DGESV(n, nrhs, a, lda, ipiv, b, ldb, info);
CGESV(n, nrhs, a, lda, ipiv, b, ldb, info);
ZGESV(n, nrhs, a, lda, ipiv, b, ldb, info);
参数
参数名 |
类型 |
描述 |
输入/输出 |
---|---|---|---|
n |
整数型 |
矩阵A的阶数,要求n≥0。 |
输入 |
nrhs |
整数型 |
右端项的数量,即矩阵B的列数,要求nrhs≥0。 |
输入 |
a |
|
矩阵维度为(lda,n)。
|
输入/输出 |
lda |
整数型 |
A的leading dimension大小,要求lda≥max(1, n)。 |
输入 |
ipiv |
整数型数组 |
数组维度为n。 存放了置换矩阵P的主元索引,即矩阵第i行与第ipiv(i)行交换。 |
输出 |
b |
|
矩阵维度为(ldb,nrhs)。
|
输入/输出 |
ldb |
整数型 |
B的leading dimension大小,要求ldb≥max(1, n)。 |
输入 |
info |
整数型 |
执行结果:
|
输出 |
依赖
#include "klapack.h"
示例
C Interface:
int n = 5; int nrhs = 2; int lda = 5; int ldb = 5; int ipiv[5]; int info = 0; double a[]={72.1673 , 66.1857 , 64.7644 , 28.0199 , 91.4151, 6.5180 , 62.8483 , 72.4323 , 46.5760 , 8.6928, 28.9821 , 42.1828 , 18.6437 , 99.8612 , 35.6972, 67.9812 , 5.0880 , 85.5035 ,79.2945 , 54.5920, 28.6869 , 49.7512 , 7.5186 ,28.6929 , 84.6041}; double b[]={ 9.4532 , 1.5204 , 2.2127 , 0.9891, 7.1778, 6.8955 , 7.2465, 3.5019 , 8.2268, 3.5287}; dgesv_(&n,&nrhs,a,&lda,ipiv,b,&ldb,&info); /* * Output: * 0.17538061067669766 0.16637572403098155 * -0.11183914210321674 -3.7758100714851153E-002 * 5.5415093265516101E-002 0.15550950667724139 * -1.4901096204948673E-002 -7.3593964815566168E-002 * -0.10693481055391466 -0.15230871441871899 */
Fortran Interface:
PARAMETER(n=5) PARAMETER(nrhs=2) PARAMETER(lda=5) PARAMETER(ldb=5) INTEGER::info =0 INTEGER :: ipiv(n) REAL(8) :: a(n*n) DATA a / 72.1673 , 66.1857 , 64.7644 , 28.0199 , 91.4151, $ 6.5180 , 62.8483 , 72.4323 , 46.5760 , 8.6928, $ 28.9821 , 42.1828 , 18.6437 , 99.8612 , 35.6972, $ 67.9812 , 5.0880 , 85.5035 ,79.2945 , 54.5920, $ 28.6869 , 49.7512 , 7.5186 ,28.6929 , 84.6041 / REAL(8) :: b(n, nrhs) DATA b /9.4532 , 1.5204 , 2.2127 , 0.9891, 7.1778, $ 6.8955 , 7.2465, 3.5019 , 8.2268, 3.5287 / EXTERNAL DGESV CALL DGESV(n, nrhs, a, lda, ipiv, b, ldb, info); * Output: * 0.17538061067669766 0.16637572403098155 * -0.11183914210321674 -3.7758100714851153E-002 * 5.5415093265516101E-002 0.15550950667724139 * -1.4901096204948673E-002 -7.3593964815566168E-002 * -0.10693481055391466 -0.15230871441871899