C 中 LAPACKE dgels_ 的分段错误
Segmentation fault with LAPACKE dgels_ in C
我想用 C 中的 LAPACK 的 dgels_ 解决最小二乘问题 |Ax-b|->min，但我遇到了分段错误（我知道有 a similar problem，但是代码完全不同，答案不适用于我的问题）。
我已经定位到问题了，就是执行dgels_的时候。
代码：
#include <stdio.h>
#include <lapacke.h>
#define COL 3
#define ROW 4
int main()
{
char transa ='N';
lapack_int m, n, nrhs, lda, ldb, info;
m=ROW;
n=COL;
nrhs=1;
lda=ROW;
ldb=ROW;
double work[COL];
double A [COL*ROW] =
{ 1.1, 4.2, 1.7,
2.5, 2.1, 2.8,
3.4, 4.2, 8.5,
4.4, 5.2, 7.8 };
double b[ROW] =
{ 1.5,
2.1,
3.8,
3.4 };
printf("test 1 \n");
dgels_(&transa, &m, &n, &nrhs, A, &lda, b, &ldb, work, 0, &info);
printf("test 2 \n");
return 0;
}
打印第一个"test 1"然后出现段错误。
编辑：
我之前尝试用值编译它，即使用 dgels_(transa, m, n, nrhs, A, lda, b, ldb, work, 0, info)。
但是，然后我收到很多错误消息：
lapack_error.c: In function ‘main’:
lapack_error.c:32:13: warning: passing argument 1 of ‘dgels_’ makes pointer from integer without a cast [-Wint-conversion]
dgels_(transa, m, n, nrhs, A, lda, b, ldb, work, 0, info);
^~~~~~
In file included from /usr/include/lapacke.h:143:0,
from lapack_error.c:2:
/usr/include/lapacke.h:14793:6: note: expected ‘char *’ but argument is of type ‘char’
void LAPACK_dgels( char* trans, lapack_int* m, lapack_int* n, lapack_int* nrhs,
^
lapack_error.c:32:21: warning: passing argument 2 of ‘dgels_’ makes pointer from integer without a cast [-Wint-conversion]
dgels_(transa, m, n, nrhs, A, lda, b, ldb, work, 0, info);
^
In file included from /usr/include/lapacke.h:143:0,
from lapack_error.c:2:
/usr/include/lapacke.h:14793:6: note: expected ‘int *’ but argument is of type ‘int’
void LAPACK_dgels( char* trans, lapack_int* m, lapack_int* n, lapack_int* nrhs,
^
lapack_error.c:32:24: warning: passing argument 3 of ‘dgels_’ makes pointer from integer without a cast [-Wint-conversion]
dgels_(transa, m, n, nrhs, A, lda, b, ldb, work, 0, info);
^
In file included from /usr/include/lapacke.h:143:0,
from lapack_error.c:2:
/usr/include/lapacke.h:14793:6: note: expected ‘int *’ but argument is of type ‘int’
void LAPACK_dgels( char* trans, lapack_int* m, lapack_int* n, lapack_int* nrhs,
^
lapack_error.c:32:27: warning: passing argument 4 of ‘dgels_’ makes pointer from integer without a cast [-Wint-conversion]
dgels_(transa, m, n, nrhs, A, lda, b, ldb, work, 0, info);
^~~~
In file included from /usr/include/lapacke.h:143:0,
from lapack_error.c:2:
/usr/include/lapacke.h:14793:6: note: expected ‘int *’ but argument is of type ‘int’
void LAPACK_dgels( char* trans, lapack_int* m, lapack_int* n, lapack_int* nrhs,
^
lapack_error.c:32:36: warning: passing argument 6 of ‘dgels_’ makes pointer from integer without a cast [-Wint-conversion]
dgels_(transa, m, n, nrhs, A, lda, b, ldb, work, 0, info);
^~~
In file included from /usr/include/lapacke.h:143:0,
from lapack_error.c:2:
/usr/include/lapacke.h:14793:6: note: expected ‘int *’ but argument is of type ‘int’
void LAPACK_dgels( char* trans, lapack_int* m, lapack_int* n, lapack_int* nrhs,
^
lapack_error.c:32:44: warning: passing argument 8 of ‘dgels_’ makes pointer from integer without a cast [-Wint-conversion]
dgels_(transa, m, n, nrhs, A, lda, b, ldb, work, 0, info);
^~~
In file included from /usr/include/lapacke.h:143:0,
from lapack_error.c:2:
/usr/include/lapacke.h:14793:6: note: expected ‘int *’ but argument is of type ‘int’
void LAPACK_dgels( char* trans, lapack_int* m, lapack_int* n, lapack_int* nrhs,
^
lapack_error.c:32:58: warning: passing argument 11 of ‘dgels_’ makes pointer from integer without a cast [-Wint-conversion]
dgels_(transa, m, n, nrhs, A, lda, b, ldb, work, 0, info);
^~~~
In file included from /usr/include/lapacke.h:143:0,
from lapack_error.c:2:
/usr/include/lapacke.h:14793:6: note: expected ‘int *’ but argument is of type ‘int’
void LAPACK_dgels( char* trans, lapack_int* m, lapack_int* n, lapack_int* nrhs,
然后我尝试查找示例程序并找到 this one（使用 sgesv，但我认为它可能与 dgels 类似）。并将 dgels_(transa, m, n, nrhs, A, lda, b, ldb, work, 0, info) 重写为 dgels_(&transa, &m, &n, &nrhs, A, &lda, b, &ldb, work, 0, &info) 没有更多的错误消息，所以我认为这是正确的方法。
0 不是 dgels. Indeed, it must be a pointer to an integer in c, since all arguments are called by reference in fortran, while all arguments are called by value in c. Moreover, the value of lwork specifies the length of the array work, which must be larger than 1 的参数 lwork 的合法值：LWORK >= max( 1, MN + max( MN, NRHS ) )。唯一的例外是 lwork=-1：在这种特殊情况下，函数 returns 为 work[0] 中的数组 work 优化大小。
例如，可以尝试以下几行：
lapack_int lwork=n;
if (m<n){lwork=m;}
lwork=-1;
double query;
dgels_(&transa, &m, &n, &nrhs, A, &lda, b, &ldb, &query, &lwork, &info);
lwork=(int)query;
printf("the optimal size is %d \n",lwork);
work=malloc(lwork*sizeof(double));
if(work==NULL){fprintf(stderr,"maloc failed\n");exit(1);}
dgels_(&transa, &m, &n, &nrhs, A, &lda, b, &ldb, work, &lwork, &info);
printf("test 2 \n");
printf("%g %g %g\n",b[0],b[1],b[2]);
有两个 dgels_ 调用：第一个 returns work 的正确大小。然后分配work，再次调用dgels_进行最小二乘最小化
结果可能与常规 C 开发人员的预期不同：Fortran and C order of dimensions for multidimensional arrays are different.The Lapacke wrapper provide the functions LAPACKE_dgels() and LAPACKE_dgels_work()，参数为 LAPACK_COL_MAJOR/LAPACK_ROW_MAJOR 和 transa。始终通过执行测试确保为这些参数使用正确的值！如果失败请尝试不同的值...
这是一个基于您的示例代码，展示了通过 Lapacke 界面调用 dgels 的不同方式。可以通过gcc main.c -o main -llapacke -llapack -lblas -lm -Wall.
编译
#include <stdio.h>
#include <lapacke.h>
#define COL 3
#define ROW 4
int main()
{
char transa ='N';
lapack_int m, n, nrhs, lda, ldb, info;
m=ROW;
n=COL;
nrhs=1;
lda=ROW;
ldb=ROW;
double* work;
double A [COL*ROW] =
{ 1.1, 4.2, 1.7, 2.5,
2.1, 2.8, 3.4, 4.2,
8.5, 4.4, 5.2, 7.8 };
double b[ROW] =
{ 39.3,
27.4,
29.3,
42.1 };
printf("test 1 \n");
lapack_int lwork=n;
if (m<n){lwork=m;}
lwork=-1;
double query;
dgels_(&transa, &m, &n, &nrhs, A, &lda, b, &ldb, &query, &lwork, &info);
lwork=(int)query;
printf("the optimal size is %d \n",lwork);
work=malloc(lwork*sizeof(double));
if(work==NULL){fprintf(stderr,"maloc failed\n");exit(1);}
dgels_(&transa, &m, &n, &nrhs, A, &lda, b, &ldb, work, &lwork, &info);
printf("test 2 \n");
printf("%g %g %g\n",b[0],b[1],b[2]);
free(work);
double Aa [COL*ROW] =
{ 1.1, 4.2, 1.7, 2.5,
2.1, 2.8, 3.4, 4.2,
8.5, 4.4, 5.2, 7.8 };
double bb[ROW] =
{ 39.3,
27.4,
29.3,
42.1 };
//with lapacke
info= LAPACKE_dgels(LAPACK_COL_MAJOR,transa, m, n, nrhs, Aa, lda, bb, ldb);
printf("%g %g %g\n",bb[0],bb[1],bb[2]);
//
double Aaa [COL*ROW] =
{ 1.1, 4.2, 1.7,
2.5, 2.1, 2.8,
3.4, 4.2, 8.5,
4.4, 5.2, 7.8 };
double bbb[ROW] =
{ 16.3,
17.9,
45.8,
46 };
//with lapacke
info= LAPACKE_dgels(LAPACK_COL_MAJOR,'T', n, m, nrhs, Aaa, n, bbb, ldb);
printf("%g %g %g\n",bbb[0],bbb[1],bbb[2]);
double Aaaa [COL*ROW] =
{ 1.1, 4.2, 1.7,
2.5, 2.1, 2.8,
3.4, 4.2, 8.5,
4.4, 5.2, 7.8 };
double bbbb[ROW] =
{ 16.3,
17.9,
45.8,
46 };
// it is still possible to allocate the buffer yourself once for all if LAPACKE_dgels_work is used.
// it can be useful if dgels is used many times, using the same transa, m,n,nrhs,lda,ldb but different A or b.
info= LAPACKE_dgels_work(LAPACK_COL_MAJOR,'T', n, m, nrhs, Aaaa, n, bbbb, ldb,&query,-1);
lwork=(int)query;
printf("the optimal size is %d \n",lwork);
work=malloc(lwork*sizeof(double));
if(work==NULL){fprintf(stderr,"maloc failed\n");exit(1);}
info= LAPACKE_dgels_work(LAPACK_COL_MAJOR,'T', n, m, nrhs, Aaaa, n, bbbb, ldb,work,lwork);
free(work);
printf("%g %g %g\n",bbbb[0],bbbb[1],bbbb[2]);
return 0;
}
我想用 C 中的 LAPACK 的 dgels_ 解决最小二乘问题 |Ax-b|->min，但我遇到了分段错误（我知道有 a similar problem，但是代码完全不同，答案不适用于我的问题）。 我已经定位到问题了，就是执行dgels_的时候。
代码：
#include <stdio.h>
#include <lapacke.h>
#define COL 3
#define ROW 4
int main()
{
char transa ='N';
lapack_int m, n, nrhs, lda, ldb, info;
m=ROW;
n=COL;
nrhs=1;
lda=ROW;
ldb=ROW;
double work[COL];
double A [COL*ROW] =
{ 1.1, 4.2, 1.7,
2.5, 2.1, 2.8,
3.4, 4.2, 8.5,
4.4, 5.2, 7.8 };
double b[ROW] =
{ 1.5,
2.1,
3.8,
3.4 };
printf("test 1 \n");
dgels_(&transa, &m, &n, &nrhs, A, &lda, b, &ldb, work, 0, &info);
printf("test 2 \n");
return 0;
}
打印第一个"test 1"然后出现段错误。
编辑： 我之前尝试用值编译它，即使用 dgels_(transa, m, n, nrhs, A, lda, b, ldb, work, 0, info)。 但是，然后我收到很多错误消息：
lapack_error.c: In function ‘main’:
lapack_error.c:32:13: warning: passing argument 1 of ‘dgels_’ makes pointer from integer without a cast [-Wint-conversion]
dgels_(transa, m, n, nrhs, A, lda, b, ldb, work, 0, info);
^~~~~~
In file included from /usr/include/lapacke.h:143:0,
from lapack_error.c:2:
/usr/include/lapacke.h:14793:6: note: expected ‘char *’ but argument is of type ‘char’
void LAPACK_dgels( char* trans, lapack_int* m, lapack_int* n, lapack_int* nrhs,
^
lapack_error.c:32:21: warning: passing argument 2 of ‘dgels_’ makes pointer from integer without a cast [-Wint-conversion]
dgels_(transa, m, n, nrhs, A, lda, b, ldb, work, 0, info);
^
In file included from /usr/include/lapacke.h:143:0,
from lapack_error.c:2:
/usr/include/lapacke.h:14793:6: note: expected ‘int *’ but argument is of type ‘int’
void LAPACK_dgels( char* trans, lapack_int* m, lapack_int* n, lapack_int* nrhs,
^
lapack_error.c:32:24: warning: passing argument 3 of ‘dgels_’ makes pointer from integer without a cast [-Wint-conversion]
dgels_(transa, m, n, nrhs, A, lda, b, ldb, work, 0, info);
^
In file included from /usr/include/lapacke.h:143:0,
from lapack_error.c:2:
/usr/include/lapacke.h:14793:6: note: expected ‘int *’ but argument is of type ‘int’
void LAPACK_dgels( char* trans, lapack_int* m, lapack_int* n, lapack_int* nrhs,
^
lapack_error.c:32:27: warning: passing argument 4 of ‘dgels_’ makes pointer from integer without a cast [-Wint-conversion]
dgels_(transa, m, n, nrhs, A, lda, b, ldb, work, 0, info);
^~~~
In file included from /usr/include/lapacke.h:143:0,
from lapack_error.c:2:
/usr/include/lapacke.h:14793:6: note: expected ‘int *’ but argument is of type ‘int’
void LAPACK_dgels( char* trans, lapack_int* m, lapack_int* n, lapack_int* nrhs,
^
lapack_error.c:32:36: warning: passing argument 6 of ‘dgels_’ makes pointer from integer without a cast [-Wint-conversion]
dgels_(transa, m, n, nrhs, A, lda, b, ldb, work, 0, info);
^~~
In file included from /usr/include/lapacke.h:143:0,
from lapack_error.c:2:
/usr/include/lapacke.h:14793:6: note: expected ‘int *’ but argument is of type ‘int’
void LAPACK_dgels( char* trans, lapack_int* m, lapack_int* n, lapack_int* nrhs,
^
lapack_error.c:32:44: warning: passing argument 8 of ‘dgels_’ makes pointer from integer without a cast [-Wint-conversion]
dgels_(transa, m, n, nrhs, A, lda, b, ldb, work, 0, info);
^~~
In file included from /usr/include/lapacke.h:143:0,
from lapack_error.c:2:
/usr/include/lapacke.h:14793:6: note: expected ‘int *’ but argument is of type ‘int’
void LAPACK_dgels( char* trans, lapack_int* m, lapack_int* n, lapack_int* nrhs,
^
lapack_error.c:32:58: warning: passing argument 11 of ‘dgels_’ makes pointer from integer without a cast [-Wint-conversion]
dgels_(transa, m, n, nrhs, A, lda, b, ldb, work, 0, info);
^~~~
In file included from /usr/include/lapacke.h:143:0,
from lapack_error.c:2:
/usr/include/lapacke.h:14793:6: note: expected ‘int *’ but argument is of type ‘int’
void LAPACK_dgels( char* trans, lapack_int* m, lapack_int* n, lapack_int* nrhs,
然后我尝试查找示例程序并找到 this one（使用 sgesv，但我认为它可能与 dgels 类似）。并将 dgels_(transa, m, n, nrhs, A, lda, b, ldb, work, 0, info) 重写为 dgels_(&transa, &m, &n, &nrhs, A, &lda, b, &ldb, work, 0, &info) 没有更多的错误消息，所以我认为这是正确的方法。
0 不是 dgels. Indeed, it must be a pointer to an integer in c, since all arguments are called by reference in fortran, while all arguments are called by value in c. Moreover, the value of lwork specifies the length of the array work, which must be larger than 1 的参数 lwork 的合法值：LWORK >= max( 1, MN + max( MN, NRHS ) )。唯一的例外是 lwork=-1：在这种特殊情况下，函数 returns 为 work[0] 中的数组 work 优化大小。
例如，可以尝试以下几行：
lapack_int lwork=n;
if (m<n){lwork=m;}
lwork=-1;
double query;
dgels_(&transa, &m, &n, &nrhs, A, &lda, b, &ldb, &query, &lwork, &info);
lwork=(int)query;
printf("the optimal size is %d \n",lwork);
work=malloc(lwork*sizeof(double));
if(work==NULL){fprintf(stderr,"maloc failed\n");exit(1);}
dgels_(&transa, &m, &n, &nrhs, A, &lda, b, &ldb, work, &lwork, &info);
printf("test 2 \n");
printf("%g %g %g\n",b[0],b[1],b[2]);
有两个 dgels_ 调用：第一个 returns work 的正确大小。然后分配work，再次调用dgels_进行最小二乘最小化
结果可能与常规 C 开发人员的预期不同：Fortran and C order of dimensions for multidimensional arrays are different.The Lapacke wrapper provide the functions LAPACKE_dgels() and LAPACKE_dgels_work()，参数为 LAPACK_COL_MAJOR/LAPACK_ROW_MAJOR 和 transa。始终通过执行测试确保为这些参数使用正确的值！如果失败请尝试不同的值...
这是一个基于您的示例代码，展示了通过 Lapacke 界面调用 dgels 的不同方式。可以通过gcc main.c -o main -llapacke -llapack -lblas -lm -Wall.
#include <stdio.h>
#include <lapacke.h>
#define COL 3
#define ROW 4
int main()
{
char transa ='N';
lapack_int m, n, nrhs, lda, ldb, info;
m=ROW;
n=COL;
nrhs=1;
lda=ROW;
ldb=ROW;
double* work;
double A [COL*ROW] =
{ 1.1, 4.2, 1.7, 2.5,
2.1, 2.8, 3.4, 4.2,
8.5, 4.4, 5.2, 7.8 };
double b[ROW] =
{ 39.3,
27.4,
29.3,
42.1 };
printf("test 1 \n");
lapack_int lwork=n;
if (m<n){lwork=m;}
lwork=-1;
double query;
dgels_(&transa, &m, &n, &nrhs, A, &lda, b, &ldb, &query, &lwork, &info);
lwork=(int)query;
printf("the optimal size is %d \n",lwork);
work=malloc(lwork*sizeof(double));
if(work==NULL){fprintf(stderr,"maloc failed\n");exit(1);}
dgels_(&transa, &m, &n, &nrhs, A, &lda, b, &ldb, work, &lwork, &info);
printf("test 2 \n");
printf("%g %g %g\n",b[0],b[1],b[2]);
free(work);
double Aa [COL*ROW] =
{ 1.1, 4.2, 1.7, 2.5,
2.1, 2.8, 3.4, 4.2,
8.5, 4.4, 5.2, 7.8 };
double bb[ROW] =
{ 39.3,
27.4,
29.3,
42.1 };
//with lapacke
info= LAPACKE_dgels(LAPACK_COL_MAJOR,transa, m, n, nrhs, Aa, lda, bb, ldb);
printf("%g %g %g\n",bb[0],bb[1],bb[2]);
//
double Aaa [COL*ROW] =
{ 1.1, 4.2, 1.7,
2.5, 2.1, 2.8,
3.4, 4.2, 8.5,
4.4, 5.2, 7.8 };
double bbb[ROW] =
{ 16.3,
17.9,
45.8,
46 };
//with lapacke
info= LAPACKE_dgels(LAPACK_COL_MAJOR,'T', n, m, nrhs, Aaa, n, bbb, ldb);
printf("%g %g %g\n",bbb[0],bbb[1],bbb[2]);
double Aaaa [COL*ROW] =
{ 1.1, 4.2, 1.7,
2.5, 2.1, 2.8,
3.4, 4.2, 8.5,
4.4, 5.2, 7.8 };
double bbbb[ROW] =
{ 16.3,
17.9,
45.8,
46 };
// it is still possible to allocate the buffer yourself once for all if LAPACKE_dgels_work is used.
// it can be useful if dgels is used many times, using the same transa, m,n,nrhs,lda,ldb but different A or b.
info= LAPACKE_dgels_work(LAPACK_COL_MAJOR,'T', n, m, nrhs, Aaaa, n, bbbb, ldb,&query,-1);
lwork=(int)query;
printf("the optimal size is %d \n",lwork);
work=malloc(lwork*sizeof(double));
if(work==NULL){fprintf(stderr,"maloc failed\n");exit(1);}
info= LAPACKE_dgels_work(LAPACK_COL_MAJOR,'T', n, m, nrhs, Aaaa, n, bbbb, ldb,work,lwork);
free(work);
printf("%g %g %g\n",bbbb[0],bbbb[1],bbbb[2]);
return 0;
}