c中的矩阵行列式使用高斯消除核心转储错误
matrix determinant in c using Gauss elimination Core dumped error
我正在尝试用 C 编写一个简单的控制台应用程序,它将使用高斯消元法计算矩阵的行列式。经过大量测试后,我发现我的程序因为核心转储而无法运行 error.After 2 天的编辑和撤消,我找不到问题所在。
我们非常欢迎任何帮助。
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
int recherche_pivot(int k, int n, float *A)
{
int i, j;
if (A[((k - 1) * n + k) - 1] != 0)
{
return k;
}
else
{ //parcours du reste de la colonne
for (i = k + 1; i <= n; i++)
{
if (A[((k - 1) * n + i) - 1] != 0)
{
return i;
}
}
return -1;
}
}
void fois(int n, float p, int i, float * A, float *b, float * x)
{
int a;
for (a = 1; a <= n; a++)
{
x[a - 1] = A[((i - 1) * n + a) - 1] * p;
}
x[n] = b[i - 1] * p;
}
void afficher_system(int n, float * X, float *b)
{
int i, j;
for (i = 1; i <= n; i++)
{
for (j = 1; j <= n; j++)
printf("%f ", X[((i - 1) * n + j) - 1]);
printf(" | %f", b[i - 1]);
printf("nn");
}
printf("nnnn");
}
void saisirmatrice(int n, float *A)
{
int i, j;
for (i = 1; i <= n; i++)
for (j = 1; j <= n; j++)
scanf("%f", &A[((i - 1) * n + j) - 1]);
}
void affichermatrice(int n, float *A)
{
int i, j;
for (i = 1; i <= n; i++)
for (j = 1; j <= n; j++)
printf("A[%d][%d] = %fn", i, j, A[((i - 1) * n + j) - 1]);
}
void elemination(int n, int k, float *b, float *A)
{
int i, l, j;
float * L, piv;
L = (float *) malloc((n) * sizeof(float));
for (i = k + 1; i <= n; i++)
{
piv = -1 * (A[((i - 1) * n + k) - 1] / A[((k - 1) * n + k) - 1]);
fois(n, piv, k, A, b, L);
//afficher_vecteur(n,L);
for (j = 1; j <= n; j++)
{
A[((i - 1) * n + j) - 1] = A[((i - 1) * n + j) - 1] + L[j - 1];
}
b[i - 1] = b[i - 1] + L[n];
afficher_system(n, A, b);
}
}
void permutter(int n, float * A, int i, int j, float * b)
{
int a;
float t[n + 1];
for (a = 1; a <= n; a++)
{
t[a - 1] = A[((i - 1) * n + a) - 1];
A[((i - 1) * n + a) - 1] = A[((j - 1) * n + a) - 1];
A[((j - 1) * n + a) - 1] = t[a - 1];
}
t[n] = b[i - 1];
b[i - 1] = b[j - 1];
b[j - 1] = t[n];
}
void main()
{
float * A, det, *L, *R, *b, s;
int i, j, i0, n, k, stop = 0;
printf("Veuillez donner la taille de la matrice");
scanf("%d", &n);
A = (float *) malloc(sizeof(float) * (n * n));
L = (float*) malloc(n * sizeof(float));
R = (float*) malloc(n * sizeof(float));
b = (float*) malloc(n * sizeof(float));
printf("Veuillez remplir la matrice");
saisirmatrice(n, A);
det = 1;
stop = 0;
k = 1;
do
{
do
{
i0 = recherche_pivot(k, n, A);
if (i0 == k)
{
//Elémination
elemination(n, k, b, A);
k++;
}
else if (i0 == -1)
{
stop = 1;
}
else
{ //cas ou ligne pivot=i0 != k
//permutation
det = -det;
permutter(n, A, k, i0, b);
//elemination
elemination(n, k, b, A);
//afficher_matrice(n,A);
k++;
}
} while ((k <= n) && (stop == 0));
} while (stop == 1 || k == n);
for (i = 1; i < n; i++)
{
det = det * A[((i - 1) * n + i) - 1];
}
printf("Le determinant est :%f", det);
free(A);
free(L);
free(R);
free(b);
}
上面的代码有很多问题。由于数组在 C 中是 zero-indexed,因此您应该从零开始计算矩阵的行和列,而不是从 1 开始计算,然后在 array-indexing 时尝试转换。不需要转换 malloc()
的结果,最好使用标识符而不是显式类型作为 sizeof
运算符的参数:
A = malloc(sizeof(*A) * n * n));
你在main()
中为L
和R
分配了space,然后直到程序结束时才使用这些指针free
]d。然后在 elemination()
函数中为 L
分配;但是你从来没有 free
这个内存,所以你有内存泄漏。您还为 main()
中的 b
分配了 space,但是在将 b
传递给 elemination()
函数之前,您没有在 b
中存储任何值。这势必会出问题。
这里本来就不需要动态分配;我建议使用可变长度数组来存储矩阵的元素。这些从 C99 开始可用,并且可以让您避免所有分配问题。
recherche_pivot()
函数有问题,你比较的地方:
if(A[((k - 1) * n + i) - 1] != 0) {}
这是一个问题,因为数组元素是一个浮点值,它是算术运算的结果;此值不应直接与 0 进行比较。我建议选择一个合适的 DELTA
值来表示零范围,而不是比较:
#define DELTA 0.000001
...
if (fabs(A[((k - 1) * n + i) - 1]) < DELTA) {}
在 permutter()
函数中,您使用数组 float t[n];
来保存临时值。但是这里不需要数组,因为您不需要在交换后保存这些临时值;而只是使用 float t;
。此外,当您交换 b[]
中的值时,您使用 t[n]
来存储临时值,但这是越界的。
elemination()
函数可能应该遍历所有行(第 k
行除外),而不是从第 k
行开始,或者它应该从第 k+1
行。实际上,第 k
行用于消除自身。最后,您在 main()
中用于执行高斯消元的实际算法已损坏。除其他事项外,调用 permutter(n, A, k, i0, b);
将第 k
行与第 i0
行交换,但 i0
是第 k
行的主元列。这没有意义。
实际上看起来您想要做的不仅仅是使用此代码计算行列式,因为您有 b
,它是线性系统的常数向量。问题标题中提到的任务不需要这样做。此外,对于任何 1X1 行列式,您的代码似乎都给出了 1 的结果。这是不正确的;本例应该是单个数字的值。
用于计算行列式的高斯消去法要求您跟踪执行了多少 row-interchanges,并且您要保留每个行相乘的任何因子的 运行 乘积。将一行的倍数添加到另一行以替换该行不会更改行列式的值,这就是下面的 reduce()
函数中使用的操作。最终结果是减少矩阵中对角线项的乘积,乘以 -1,每 row-interchange 操作一次,除以所有用于计算的因子的乘积缩放单个行。在这种情况下,没有这样的因素,所以结果只是简化矩阵的对角线元素与符号校正的乘积。这是原问题贴出的代码使用的方法。
这里有太多问题,我刚刚编写了一个实现该算法的新程序。我认为它与您试图完成的目标很接近,至少在精神上是这样。我确实为矩阵的大小添加了一些输入验证,检查以确保用户输入了正数,并在输入错误时提示 re-entry。填充矩阵的输入循环将受益于类似的输入验证。另请注意,输入大小存储在 signed int
中,以允许检查负输入,成功的输入被转换并存储在 size_t
类型的变量中,这是一个 unsigned
保证包含任何数组索引的整数类型。这是索引数组时使用的正确类型,您会注意到在整个程序中使用 size_t
。
#include <stdio.h>
#include <math.h>
#include <stdbool.h>
#define DELTA 0.000001
void show_matrix(size_t mx_sz, double mx[mx_sz][mx_sz]);
void interchange(size_t r1, size_t r2, size_t mx_sz, double mx[mx_sz][mx_sz]);
void reduce(double factor, size_t r1, size_t r2,
size_t mx_sz, double mx[mx_sz][mx_sz]);
size_t get_pivot(size_t row, size_t mx_sz, double mx[mx_sz][mx_sz]);
double find_det(size_t mx_sz, double mx[mx_sz][mx_sz]);
int main(void)
{
size_t n;
int read_val, c;
printf("Enter size of matrix: ");
while (scanf("%d", &read_val) != 1 || read_val < 1) {
while ((c = getchar()) != '\n' && c != EOF) {
continue; // discard extra characters
}
printf("Enter size of matrix: ");
}
n = (size_t) read_val;
double matrix[n][n];
printf("Enter matrix elements:\n");
for (size_t i = 0; i < n; i++) {
for (size_t j = 0; j < n; j++) {
scanf("%lf", &matrix[i][j]);
}
}
printf("You entered:\n");
show_matrix(n, matrix);
putchar('\n');
double result = find_det(n, matrix);
show_matrix(n, matrix);
putchar('\n');
printf("Determinant: %f\n", result);
return 0;
}
void show_matrix(size_t n, double mx[n][n])
{
for (size_t i = 0; i < n; i++) {
for (size_t j = 0; j < n; j++) {
printf("%7.2f", mx[i][j]);
}
putchar('\n');
}
}
/* interchange rows r1 and r2 */
void interchange(size_t r1, size_t r2, size_t mx_sz, double mx[mx_sz][mx_sz])
{
double temp;
for (size_t j = 0; j < mx_sz; j++) {
temp = mx[r1][j];
mx[r1][j] = mx[r2][j];
mx[r2][j] = temp;
}
}
/* add factor * row r1 to row r2 to replace row r2 */
void reduce(double factor, size_t r1, size_t r2,
size_t mx_sz, double mx[mx_sz][mx_sz])
{
for (size_t j = 0; j < mx_sz; j++) {
mx[r2][j] += (factor * mx[r1][j]);
}
}
/* returns pivot column, or mx_sz if there is no pivot */
size_t get_pivot(size_t row, size_t mx_sz, double mx[mx_sz][mx_sz])
{
size_t j = 0;
while (j < mx_sz && fabs(mx[row][j]) < DELTA) {
++j;
}
return j;
}
double find_det(size_t mx_sz, double mx[mx_sz][mx_sz])
{
size_t pivot1, pivot2;
size_t row;
double factor;
bool finished = false;
double result = 1.0;
while (!finished) {
finished = true;
row = 1;
while (row < mx_sz) {
// determinant is zero if there is a zero row
if ((pivot1 = get_pivot(row - 1, mx_sz, mx)) == mx_sz ||
(pivot2 = get_pivot(row, mx_sz, mx)) == mx_sz) {
return 0.0;
}
if (pivot1 == pivot2) {
factor = -mx[row][pivot1] / mx[row - 1][pivot1];
reduce(factor, row - 1, row, mx_sz, mx);
finished = false;
} else if (pivot2 < pivot1) {
interchange(row - 1, row, mx_sz, mx);
result = -result;
finished = false;
}
++row;
}
}
for (size_t j = 0; j < mx_sz; j++) {
result *= mx[j][j];
}
return result;
}
样本session:
Enter size of matrix: oops
Enter size of matrix: 0
Enter size of matrix: -1
Enter size of matrix: 3
Enter matrix elements:
0 1 3
1 2 0
0 3 4
You entered:
0.00 1.00 3.00
1.00 2.00 0.00
0.00 3.00 4.00
1.00 2.00 0.00
-0.00 -3.00 -9.00
0.00 0.00 -5.00
Determinant: 5.000000
我正在尝试用 C 编写一个简单的控制台应用程序,它将使用高斯消元法计算矩阵的行列式。经过大量测试后,我发现我的程序因为核心转储而无法运行 error.After 2 天的编辑和撤消,我找不到问题所在。 我们非常欢迎任何帮助。
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
int recherche_pivot(int k, int n, float *A)
{
int i, j;
if (A[((k - 1) * n + k) - 1] != 0)
{
return k;
}
else
{ //parcours du reste de la colonne
for (i = k + 1; i <= n; i++)
{
if (A[((k - 1) * n + i) - 1] != 0)
{
return i;
}
}
return -1;
}
}
void fois(int n, float p, int i, float * A, float *b, float * x)
{
int a;
for (a = 1; a <= n; a++)
{
x[a - 1] = A[((i - 1) * n + a) - 1] * p;
}
x[n] = b[i - 1] * p;
}
void afficher_system(int n, float * X, float *b)
{
int i, j;
for (i = 1; i <= n; i++)
{
for (j = 1; j <= n; j++)
printf("%f ", X[((i - 1) * n + j) - 1]);
printf(" | %f", b[i - 1]);
printf("nn");
}
printf("nnnn");
}
void saisirmatrice(int n, float *A)
{
int i, j;
for (i = 1; i <= n; i++)
for (j = 1; j <= n; j++)
scanf("%f", &A[((i - 1) * n + j) - 1]);
}
void affichermatrice(int n, float *A)
{
int i, j;
for (i = 1; i <= n; i++)
for (j = 1; j <= n; j++)
printf("A[%d][%d] = %fn", i, j, A[((i - 1) * n + j) - 1]);
}
void elemination(int n, int k, float *b, float *A)
{
int i, l, j;
float * L, piv;
L = (float *) malloc((n) * sizeof(float));
for (i = k + 1; i <= n; i++)
{
piv = -1 * (A[((i - 1) * n + k) - 1] / A[((k - 1) * n + k) - 1]);
fois(n, piv, k, A, b, L);
//afficher_vecteur(n,L);
for (j = 1; j <= n; j++)
{
A[((i - 1) * n + j) - 1] = A[((i - 1) * n + j) - 1] + L[j - 1];
}
b[i - 1] = b[i - 1] + L[n];
afficher_system(n, A, b);
}
}
void permutter(int n, float * A, int i, int j, float * b)
{
int a;
float t[n + 1];
for (a = 1; a <= n; a++)
{
t[a - 1] = A[((i - 1) * n + a) - 1];
A[((i - 1) * n + a) - 1] = A[((j - 1) * n + a) - 1];
A[((j - 1) * n + a) - 1] = t[a - 1];
}
t[n] = b[i - 1];
b[i - 1] = b[j - 1];
b[j - 1] = t[n];
}
void main()
{
float * A, det, *L, *R, *b, s;
int i, j, i0, n, k, stop = 0;
printf("Veuillez donner la taille de la matrice");
scanf("%d", &n);
A = (float *) malloc(sizeof(float) * (n * n));
L = (float*) malloc(n * sizeof(float));
R = (float*) malloc(n * sizeof(float));
b = (float*) malloc(n * sizeof(float));
printf("Veuillez remplir la matrice");
saisirmatrice(n, A);
det = 1;
stop = 0;
k = 1;
do
{
do
{
i0 = recherche_pivot(k, n, A);
if (i0 == k)
{
//Elémination
elemination(n, k, b, A);
k++;
}
else if (i0 == -1)
{
stop = 1;
}
else
{ //cas ou ligne pivot=i0 != k
//permutation
det = -det;
permutter(n, A, k, i0, b);
//elemination
elemination(n, k, b, A);
//afficher_matrice(n,A);
k++;
}
} while ((k <= n) && (stop == 0));
} while (stop == 1 || k == n);
for (i = 1; i < n; i++)
{
det = det * A[((i - 1) * n + i) - 1];
}
printf("Le determinant est :%f", det);
free(A);
free(L);
free(R);
free(b);
}
上面的代码有很多问题。由于数组在 C 中是 zero-indexed,因此您应该从零开始计算矩阵的行和列,而不是从 1 开始计算,然后在 array-indexing 时尝试转换。不需要转换 malloc()
的结果,最好使用标识符而不是显式类型作为 sizeof
运算符的参数:
A = malloc(sizeof(*A) * n * n));
你在main()
中为L
和R
分配了space,然后直到程序结束时才使用这些指针free
]d。然后在 elemination()
函数中为 L
分配;但是你从来没有 free
这个内存,所以你有内存泄漏。您还为 main()
中的 b
分配了 space,但是在将 b
传递给 elemination()
函数之前,您没有在 b
中存储任何值。这势必会出问题。
这里本来就不需要动态分配;我建议使用可变长度数组来存储矩阵的元素。这些从 C99 开始可用,并且可以让您避免所有分配问题。
recherche_pivot()
函数有问题,你比较的地方:
if(A[((k - 1) * n + i) - 1] != 0) {}
这是一个问题,因为数组元素是一个浮点值,它是算术运算的结果;此值不应直接与 0 进行比较。我建议选择一个合适的 DELTA
值来表示零范围,而不是比较:
#define DELTA 0.000001
...
if (fabs(A[((k - 1) * n + i) - 1]) < DELTA) {}
在 permutter()
函数中,您使用数组 float t[n];
来保存临时值。但是这里不需要数组,因为您不需要在交换后保存这些临时值;而只是使用 float t;
。此外,当您交换 b[]
中的值时,您使用 t[n]
来存储临时值,但这是越界的。
elemination()
函数可能应该遍历所有行(第 k
行除外),而不是从第 k
行开始,或者它应该从第 k+1
行。实际上,第 k
行用于消除自身。最后,您在 main()
中用于执行高斯消元的实际算法已损坏。除其他事项外,调用 permutter(n, A, k, i0, b);
将第 k
行与第 i0
行交换,但 i0
是第 k
行的主元列。这没有意义。
实际上看起来您想要做的不仅仅是使用此代码计算行列式,因为您有 b
,它是线性系统的常数向量。问题标题中提到的任务不需要这样做。此外,对于任何 1X1 行列式,您的代码似乎都给出了 1 的结果。这是不正确的;本例应该是单个数字的值。
用于计算行列式的高斯消去法要求您跟踪执行了多少 row-interchanges,并且您要保留每个行相乘的任何因子的 运行 乘积。将一行的倍数添加到另一行以替换该行不会更改行列式的值,这就是下面的 reduce()
函数中使用的操作。最终结果是减少矩阵中对角线项的乘积,乘以 -1,每 row-interchange 操作一次,除以所有用于计算的因子的乘积缩放单个行。在这种情况下,没有这样的因素,所以结果只是简化矩阵的对角线元素与符号校正的乘积。这是原问题贴出的代码使用的方法。
这里有太多问题,我刚刚编写了一个实现该算法的新程序。我认为它与您试图完成的目标很接近,至少在精神上是这样。我确实为矩阵的大小添加了一些输入验证,检查以确保用户输入了正数,并在输入错误时提示 re-entry。填充矩阵的输入循环将受益于类似的输入验证。另请注意,输入大小存储在 signed int
中,以允许检查负输入,成功的输入被转换并存储在 size_t
类型的变量中,这是一个 unsigned
保证包含任何数组索引的整数类型。这是索引数组时使用的正确类型,您会注意到在整个程序中使用 size_t
。
#include <stdio.h>
#include <math.h>
#include <stdbool.h>
#define DELTA 0.000001
void show_matrix(size_t mx_sz, double mx[mx_sz][mx_sz]);
void interchange(size_t r1, size_t r2, size_t mx_sz, double mx[mx_sz][mx_sz]);
void reduce(double factor, size_t r1, size_t r2,
size_t mx_sz, double mx[mx_sz][mx_sz]);
size_t get_pivot(size_t row, size_t mx_sz, double mx[mx_sz][mx_sz]);
double find_det(size_t mx_sz, double mx[mx_sz][mx_sz]);
int main(void)
{
size_t n;
int read_val, c;
printf("Enter size of matrix: ");
while (scanf("%d", &read_val) != 1 || read_val < 1) {
while ((c = getchar()) != '\n' && c != EOF) {
continue; // discard extra characters
}
printf("Enter size of matrix: ");
}
n = (size_t) read_val;
double matrix[n][n];
printf("Enter matrix elements:\n");
for (size_t i = 0; i < n; i++) {
for (size_t j = 0; j < n; j++) {
scanf("%lf", &matrix[i][j]);
}
}
printf("You entered:\n");
show_matrix(n, matrix);
putchar('\n');
double result = find_det(n, matrix);
show_matrix(n, matrix);
putchar('\n');
printf("Determinant: %f\n", result);
return 0;
}
void show_matrix(size_t n, double mx[n][n])
{
for (size_t i = 0; i < n; i++) {
for (size_t j = 0; j < n; j++) {
printf("%7.2f", mx[i][j]);
}
putchar('\n');
}
}
/* interchange rows r1 and r2 */
void interchange(size_t r1, size_t r2, size_t mx_sz, double mx[mx_sz][mx_sz])
{
double temp;
for (size_t j = 0; j < mx_sz; j++) {
temp = mx[r1][j];
mx[r1][j] = mx[r2][j];
mx[r2][j] = temp;
}
}
/* add factor * row r1 to row r2 to replace row r2 */
void reduce(double factor, size_t r1, size_t r2,
size_t mx_sz, double mx[mx_sz][mx_sz])
{
for (size_t j = 0; j < mx_sz; j++) {
mx[r2][j] += (factor * mx[r1][j]);
}
}
/* returns pivot column, or mx_sz if there is no pivot */
size_t get_pivot(size_t row, size_t mx_sz, double mx[mx_sz][mx_sz])
{
size_t j = 0;
while (j < mx_sz && fabs(mx[row][j]) < DELTA) {
++j;
}
return j;
}
double find_det(size_t mx_sz, double mx[mx_sz][mx_sz])
{
size_t pivot1, pivot2;
size_t row;
double factor;
bool finished = false;
double result = 1.0;
while (!finished) {
finished = true;
row = 1;
while (row < mx_sz) {
// determinant is zero if there is a zero row
if ((pivot1 = get_pivot(row - 1, mx_sz, mx)) == mx_sz ||
(pivot2 = get_pivot(row, mx_sz, mx)) == mx_sz) {
return 0.0;
}
if (pivot1 == pivot2) {
factor = -mx[row][pivot1] / mx[row - 1][pivot1];
reduce(factor, row - 1, row, mx_sz, mx);
finished = false;
} else if (pivot2 < pivot1) {
interchange(row - 1, row, mx_sz, mx);
result = -result;
finished = false;
}
++row;
}
}
for (size_t j = 0; j < mx_sz; j++) {
result *= mx[j][j];
}
return result;
}
样本session:
Enter size of matrix: oops
Enter size of matrix: 0
Enter size of matrix: -1
Enter size of matrix: 3
Enter matrix elements:
0 1 3
1 2 0
0 3 4
You entered:
0.00 1.00 3.00
1.00 2.00 0.00
0.00 3.00 4.00
1.00 2.00 0.00
-0.00 -3.00 -9.00
0.00 0.00 -5.00
Determinant: 5.000000