用LAPACKE求解线性方程组

Question

我试图解一些线性方程(对称的，三对角的，正的)。我得用LAPACKE。我的代码如下：

#include <lapacke.h>
#include <stdio.h>



void print_mtrx(double * mtrx, int n, int m)
{
    int i, j;

    for(i = 0; i < n; i++)
    {
        for(j = 0; j < m; j++)
        {
            printf("%f ", mtrx[i*m+j]);
        }
        printf("\n");
    }
    printf("\n");
}

int main()
{
    double matrix[5*5] = {
        2,  0,  0,  0,  0,
        0,  2,  0,  0,  0,
        0,  0,  2,  0,  0,
        0,  0,  0,  2,  0,
        0,  0,  0,  0,  2
    };

    double rozw[5] = {1,2,3,4,5};

    double matrix2[5*5] = {
        7,  0,  0,  0,  0,
        0,  7,  0,  0,  0,
        0,  0,  7,  0,  0,
        0,  0,  0,  7,  0,
        0,  0,  0,  0,  7
    };



    LAPACKE_dptsv(LAPACK_COL_MAJOR, 5, 5, matrix, matrix2, rozw, 5);

    print_mtrx(matrix, 5, 5);
    print_mtrx(matrix2, 5, 5);
    print_mtrx(rozw, 5, 1);

}

LAPACKE的函数似乎什么也不做，没有任何错误。主要的问题是，我不知道函数参数代表什么。我搜索了很久，但没有真正的文档。下面是我设法找到或猜到的：

• int matrix_order -- LAPACK_COL_MAJOR或LAPACK_ROW_MAJOR，如何在内存中表示矩阵
• lapack_int n-矩阵的大小(即列数)
• lapack_int nrhs --不确定，可能是矢量b的大小
• 双*d-方程矩阵
• 双份*电子--不知道。
• 双*b-d方程解的向量
• lapack_int ldb -- b的领先方向(为什么？它是否与nrhs相同，而nrhs本身与n相同？

我怎样才能找到这些论点的真正意义呢？我如何使我的代码工作？

Answer 1

当谈到BLAS和/或LAPACK的文档时，英特尔可能是最全面的。您可以查找?ptsv的文档，这说明了每个参数的作用。

(提示:在谷歌搜索BLAS或LAPACK时，一定要删除s/d/c/z前缀。)

下面是相关的片段：

该程序对X求解真实或复杂的线性方程组A*X = B，其中A是一个nbyn对称/Hermitian正定三对角矩阵，矩阵B的列是单独的右手边，X的列是相应的解。 将A分解为A = L*D*LT (实香精)或A = L*D*LH (复合风味)，然后用A的因式来求解方程组A*X = B。 输入参数 n：矩阵A；n ≥ 0的阶数。 nrhs：B中的右侧数，列数；nrhs ≥ 0。 d：数组，维度至少max(1, n)。包含三对角矩阵A的对角线元素。 e，b：数组：e(n - 1)，b(ldb,*)。数组e包含A的(n - 1)次对角线元素。数组b包含矩阵B，其列是方程组的右侧。b的第二个维度至少必须是max(1,nrhs)。 ldb：b的领先维度；ldb ≥ max(1, n)。 输出参数 d：由对角线矩阵D的n对角线元素从A的L*D*LT (实)/ L*D*LH (复)分解中覆盖的。 e：由单元双对角线因子L的(n - 1)次对角线元素从A分解而来的覆盖。 b：由解矩阵X覆盖。 info：如果是info = 0，执行是成功的。如果是info = -i，则i-th参数具有一个非法值。如果是info = i，则i阶的前导子(以及矩阵A本身)不是正定的，而且解还没有计算出来。除非i = n，否则分解尚未完成。

Answer 2

因此，我们必须查看纯LAPACK (8f.html#af1bd4c731915bd8755a4da8086fd79a8)的文档，并忽略不正确(如果是LAPACKE)的注释，即LDB大于或等于max(1，N)。

正确的程序如下：

#include <lapacke.h>
#include <stdio.h>



void print_mtrx(double * mtrx, int n, int m)
{
    int i, j;

    for(i = 0; i < n; i++)
    {
        for(j = 0; j < m; j++)
        {
            printf("%f ", mtrx[i*m+j]);
        }
        printf("\n");
    }
    printf("\n");
}

int main()
{
    double diagonal[5] = {5,1,5,1,5};
    double subdiagonal[4] = {0,0,0,0};

    double solution[5] = {1,2,3,4,5};


    LAPACKE_dptsv(LAPACK_ROW_MAJOR, 5 /*size of matrix*/, 1 /*number of columns in solution*/,
                  diagonal, subdiagonal, solution, 1 /*leading dimension of solution vector*/);

    print_mtrx(solution, 5, 1);
}