LAPACKE矩阵求逆分段故障C
LAPACKE matrix inversion segmentation fault C
我正在尝试使用 lapacke 库编写代码来反转 C 中的复杂矩阵。但是我遇到了一个分段错误,它似乎取决于矩阵的大小 N。更重要的是,每次编译程序或触摸任何东西时,发生分段错误的大小都会有所不同。这让我觉得代码在某个地方试图访问分配不当或被禁止的内存。不幸的是,我不明白这是怎么发生的,因为它似乎与 LAPACKE 函数主题本身有关。事实上,当函数 /*MatrixComplexInv(invA,A,N);*/
(其中调用 LAPACKE 函数进行反转)被注释时,分段错误 不会 发生。
下面是可以编译的工作代码和运行自己的代码。
#include <stdio.h>
#include <lapacke.h>
#include <complex.h>
#include <stdlib.h>
#include <math.h>
void Ctranspose( double complex *, double complex * ,int );
void MatrixComplexInv(double complex *, double complex *, int );
int main(int argc, const char * * argv) {
int i,j,k,N = 4;/*if N> bigger than a small number 4,6,7.. it gives segmentation fault*/
double complex *A = calloc(N*N,sizeof(double complex)),
*b = calloc(N*N,sizeof(double complex)),
*Ap =calloc(N*N,sizeof(double complex));
double complex *invA =calloc(N*N,sizeof(double complex));
for(i=0;i<N;i++){
for(j=0;j<N;j++){
A[i*N+j] = 1+sin(i*j)*i+I*j;
Ap[i*N+j] = 1+sin(i*j)*i+I*j;
}
}
/*Segmentation fault in this function, due to
*
LAPACKE_zgetrf(LAPACK_ROW_MAJOR, n, n, tempA , n,&n);
LAPACKE_zgetri(LAPACK_ROW_MAJOR, n, tempA , n, &n );
*
* both.
*/
MatrixComplexInv(invA,A,N);
for(i=0;i<N;i++){
for(j=0;j<N;j++){
for(k = 0;k<N;k++){
b[i*N+j]+=invA[i*N + k]*Ap[k*N + j];
}
printf("(%lf,%lf)\t", creal(b[i*N + j]),cimag(b[i*N + j]));/*tests that the result produces the inverse matrix A^{-1}A = 1*/
}
printf("\n");
}
return 0;
}
void Ctranspose( double complex *Transposed, double complex *M ,int n)
{
int i,j;
for(i=0;i<n;i++)
for(j=0;j<n;j++) Transposed[i+n*j] = M[i*n+j];
}
void MatrixComplexInv(double complex *invA, double complex *A, int n)
{
double complex *tempA = (double complex*) malloc( n*n*sizeof(double complex) );
Ctranspose(tempA,A,n);
/*SEGMENTATION HAPPEN IN THESE TWO FUNCTIONS*/
LAPACKE_zgetrf(LAPACK_ROW_MAJOR, n, n, tempA , n,&n);
LAPACKE_zgetri(LAPACK_ROW_MAJOR, n, tempA , n, &n );
Ctranspose(invA,tempA,n);
free(tempA);
}
在 LAPACKE_zgetrf(LAPACK_ROW_MAJOR, n, n, tempA , n,&n);
中,LAPACKE_zgetrf points to n, a single integer. On the contrary, the argument ipiv
should be a pointer to an array of integer of dimension max(m,n), to store pivot indices 的最后一个参数这可以解释分段错误。
LAPACKE_zgetrf()
计算的 ipiv
也必须作为输入提供给 LAPACKE_zgetri()
,以获得正确的逆矩阵。
我正在尝试使用 lapacke 库编写代码来反转 C 中的复杂矩阵。但是我遇到了一个分段错误,它似乎取决于矩阵的大小 N。更重要的是,每次编译程序或触摸任何东西时,发生分段错误的大小都会有所不同。这让我觉得代码在某个地方试图访问分配不当或被禁止的内存。不幸的是,我不明白这是怎么发生的,因为它似乎与 LAPACKE 函数主题本身有关。事实上,当函数 /*MatrixComplexInv(invA,A,N);*/
(其中调用 LAPACKE 函数进行反转)被注释时,分段错误 不会 发生。
下面是可以编译的工作代码和运行自己的代码。
#include <stdio.h>
#include <lapacke.h>
#include <complex.h>
#include <stdlib.h>
#include <math.h>
void Ctranspose( double complex *, double complex * ,int );
void MatrixComplexInv(double complex *, double complex *, int );
int main(int argc, const char * * argv) {
int i,j,k,N = 4;/*if N> bigger than a small number 4,6,7.. it gives segmentation fault*/
double complex *A = calloc(N*N,sizeof(double complex)),
*b = calloc(N*N,sizeof(double complex)),
*Ap =calloc(N*N,sizeof(double complex));
double complex *invA =calloc(N*N,sizeof(double complex));
for(i=0;i<N;i++){
for(j=0;j<N;j++){
A[i*N+j] = 1+sin(i*j)*i+I*j;
Ap[i*N+j] = 1+sin(i*j)*i+I*j;
}
}
/*Segmentation fault in this function, due to
*
LAPACKE_zgetrf(LAPACK_ROW_MAJOR, n, n, tempA , n,&n);
LAPACKE_zgetri(LAPACK_ROW_MAJOR, n, tempA , n, &n );
*
* both.
*/
MatrixComplexInv(invA,A,N);
for(i=0;i<N;i++){
for(j=0;j<N;j++){
for(k = 0;k<N;k++){
b[i*N+j]+=invA[i*N + k]*Ap[k*N + j];
}
printf("(%lf,%lf)\t", creal(b[i*N + j]),cimag(b[i*N + j]));/*tests that the result produces the inverse matrix A^{-1}A = 1*/
}
printf("\n");
}
return 0;
}
void Ctranspose( double complex *Transposed, double complex *M ,int n)
{
int i,j;
for(i=0;i<n;i++)
for(j=0;j<n;j++) Transposed[i+n*j] = M[i*n+j];
}
void MatrixComplexInv(double complex *invA, double complex *A, int n)
{
double complex *tempA = (double complex*) malloc( n*n*sizeof(double complex) );
Ctranspose(tempA,A,n);
/*SEGMENTATION HAPPEN IN THESE TWO FUNCTIONS*/
LAPACKE_zgetrf(LAPACK_ROW_MAJOR, n, n, tempA , n,&n);
LAPACKE_zgetri(LAPACK_ROW_MAJOR, n, tempA , n, &n );
Ctranspose(invA,tempA,n);
free(tempA);
}
在 LAPACKE_zgetrf(LAPACK_ROW_MAJOR, n, n, tempA , n,&n);
中,LAPACKE_zgetrf points to n, a single integer. On the contrary, the argument ipiv
should be a pointer to an array of integer of dimension max(m,n), to store pivot indices 的最后一个参数这可以解释分段错误。
LAPACKE_zgetrf()
计算的 ipiv
也必须作为输入提供给 LAPACKE_zgetri()
,以获得正确的逆矩阵。