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()中为LR分配了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