матричный определитель в c с использованием исключения Гаусса Ошибка сбрасывания ядра

Question

Я пытаюсь создать простое консольное приложение на C, которое будет вычислять детерминант матрицы с использованием исключения Гаусса. после многих тестов я узнал, что моя программа не работает из-за ошибки ядра. После двух дней редактирования и отмены я не смог найти проблему. Любая помощь более чем приветствуется.

#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); 
} 

Answer 1

В приведенном выше кодексе есть много проблем. Поскольку массивы с нулевой индексацией на C, вы должны считать строки и столбцы ваших матриц, начиная с нуля, вместо подсчета с 1, а затем попытаться преобразовать при индексировании массива. Там нет необходимости подавать результат в malloc(), и лучше использовать идентификатор вместо явного типа в качестве аргумента для оператора sizeof:

A = malloc(sizeof(*A) * n * n)); 

Вы выделить пространство для L и R в main(), и то никогда не используйте эти указатели до конца программы, когда они free d. Затем вы назначаете для L функцию elemination(); но вы никогда не free этой памяти, поэтому у вас есть утечка памяти. Вы также выделяете место для b в main(), но вы не храните никаких значений в b перед передачей его функции elemination(). Это неизбежно вызовет проблемы.

Здесь нет необходимости в динамическом размещении; Я предлагаю использовать массив переменной длины для хранения элементов матрицы. Они доступны с C99 и позволят избежать всех проблем с распределением.

Существует проблема в функции recherche_pivot(), где вы сравните:

if(A[((k - 1) * n + i) - 1] != 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 детерминант. Это неверно; в этом случае это должно быть значение единственного числа.

Метод исключения Гаусса для расчета детерминанты требует, чтобы вы отслеживали, сколько выполняется ряд развязок строк, и что вы сохраняете запущенный продукт любых факторов, по которым умножаются отдельные строки. Добавление нескольких строк из одной строки в другую для замены этой строки не изменяет значение детерминанта, и это операция, используемая в функции reduce() ниже. Конечным результатом является произведение диагональных элементов в приведенной матрице, умноженное на -1 один раз для каждой операции обмена строк, деленный на произведение всех факторов, используемых для масштабирования отдельных строк. В этом случае таких факторов нет, поэтому результат является просто произведением диагональных элементов приведенной матрицы со знаковой поправкой. Это метод, используемый кодом, опубликованным в исходном вопросе.

Здесь было так много проблем, что я только что написал новую программу, которая реализует этот алгоритм. Я думаю, что это близко, по крайней мере по духу, к тому, что вы пытались выполнить. Я добавил некоторую проверку ввода для размера матрицы, убедившись, что пользователь вводит положительное число и запрашивает повторную запись, если вход плохой. Входной цикл, который заполняет матрицу, выиграет от аналогичной проверки ввода. Также обратите внимание, что размер ввода сохраняется в 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; 
} 

Пример сессии:

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