c - 使用 Lapack 在 C 中求解矩阵

Question

所以我有两个矩阵；一个是 2 x 2，另一个是 2 x 1。我想使用外部库 lapack 来解决线性系统，似乎我需要调用函数 dgesv_() 这是

   Ax = B

而且我必须求解 x。

所以我真的很困惑他们所说的调用函数的符号以及如何将它转换成我现在拥有的两个数组。

 #include <stdlib.h>
 #include <stdio.h>


static double Angstroms[2];
static double Energy[2];
static double ax[2][2];

void file_input ();
void polynomial ();


int main () {

     file_input ();
     polynomial ();
     return 0;
}

 void file_input () {

     float a, b;
     int i;

     FILE * in_file = fopen("H2Mini.txt", "r");
     if (outfile == NULL) {
         printf ("Error file does not exist");
         exit (-1);
     }
     for (i = 0; i <= 1; i++) {
          fscanf(in_file, "%f %f\n", &a, &b);
          Angstroms[i] = a;
          Energy [i] = b;
     }
     fclose(in_file);

}

 void polynomial () {

      int i;
      FILE * outfile = fopen("PolyTest2.txt", "w");
      if (outfile == NULL) {
         printf ("Error file does not exist");
         exit (-1);
      }
     for (i = 0; i <= 1; i++) {
         ax[i][0] = 1;
         fprintf (outfile, "%.8f ", ax[i][0]);

     }
     fprintf (outfile, "\n");
     for (i = 0; i <= 1; i ++) {
         ax[i][1] = Angstroms[i];
         fprintf (outfile, "%.8f ", ax[i][1]);
     }
  }

所以我的文件是这个

    [2.00000000 3.00000000
     6.00000000 5.00000000]

ax 数组看起来像这样

    [1.00000000 1.00000000
     2.00000000 6.00000000]

能量阵列是这样的

    [3.00000000 5.00000000]

我的 ax 数组是 Ax = b 方程中的 Ax 项，我的 b 项是能量数组。

我确实看过他们的功能文档，实现起来有点困惑。在 dgesv (n, nrhs, a, lda, ipiv, b, ldb, info)

此代码的任何真正清晰的示例都会非常有帮助!

Answer 1

我不确定你所说的“斧头阵列看起来像这样”是什么意思？ Ax 应该是一个 vector ，不是吗？

但是一般来说，使用 C 包装器 LAPACKE 调用 dgesv 应该如下所示:

LAPACKE_dgesv(
            LAPACK_ROW_MAJOR, // The storage format you use, usually row major in c.
            n,                // Dimensions of A (nxn)
            1,                // Number of b's, here always one because you want to solve for a single right hand side.
            A,                // The input matrix A, it will be destroyed in the call so best make a copy.
            n,                // LDA basically the length of each column of A in memory. Normally this is also n if A is nxn.
            p,                // An additional array where lapacke saves pivot ordering. Probably not of any interest to you. Just give it an array of size n. 
            bToX,             // On input input the rhs b, on output this will be the solution vector x
            1 );