Программирование: Метод Гаусса

[Fenix]

Метод Гаусса
Данный метод применяется для решения системы линейных алгебраических уравнений (далее СЛАУ).
Суть метода:
Прямой ход:
Применяя элементарные преобразования, такие как:
-перестановка строк
-умножение строки на константу
-прибавление к одной строке матрицы другой строки
(аналогично определяются преобразования столбцов)
свести матрицу к верхнетреугольному виду:
а[1][1] a[1][2] … a[1][n] | b[1]
0 a[2][2] … a[2][n] | b[2]
,,,
0 0 … a[m][n] | b[m]

Обратных ход:
Последовательно находим все корни.

Пример:
2 -3 1 | -2
-3 1 1 | 2
2 4 -1 | 3
Сводим к верхнетреугольному виду:
Занулим 1-й столбец, начиная с 2-го элемента.
2 -3 1 | -2
-3 1 1 | 2 +1-я строка *(-3/2)
2 4 -1 | 3 +1-я строка *(-1)
Получили:
2 -3 1 | -2
0 -7/2 5/2 | -1
0 7 -2 | 5
Зануляем 2-й столбец начиная с 3-го элемента:
2 -3 1 | -2
0 -7/2 5/2 | -1
0 7 -2 | 3 +2я-я строка *2
Получили:
2 -3 1 | -2
0 -7/2 5/2 | -1
0 0 3 | 3
Прямой ход закончен.
Обратных ход:
СЛАУ имеет вид:
A[1][1]*X[1] + A[1][2]*X[2] + A[1][3]*X[3]=B[1]
A[2][1]*X[1] + A[2][2]*X[2] + A[2][n]*X[3]=B[2]
A[3][1]*X[1] + A[3][2]*X[2] + A[3][n]*X[3]=B[3]
Подставим коэффициенты:
2*X[1] - 3*X[2] + 1*X[3]=-2
(-7/2)*X[2] + (5/2)*X[3]=-1
3*X[3]=3
Обратных ход заключается в последовательном нахождении переменных Х.
X[3]= 3/3=1
X[2]=(-1- (5/2)*3 )/(-7/2)=1
X[1]=(-2 +3*1 - 1*1)/2=0

Кому не понятно: Википедия: Метод Гаусса

Ниже реализация метода на С++.
Недостатки:
- Не переставляет строки. поэтому на главной диагонали не могут стоять нули.

Входные данные в файле slau.in
Формат входного файла:
1-я строка - размер матрицы, далее сама матрица.

//gauss.cpp

#include <iostream>
#include <fstream>
using namespace std;

int main(){

    fstream ifile;
    ifile.open("slau.in",ios::in);
    if(!ifile.is_open()) return 1;

    // intput
    int size;
    ifile>>size;
    float slau[size][size+1];
    for(int i=0;i<size;i++)
        for(int j=0;j<size+1;j++)
            ifile>>slau[i][j];
    ifile.close();

    // Gauss' method
    // forward
    float factor;

    for(int i=0;i<size-1;i++){
        for(int j=i+1;j<size;j++){
            if(slau[i][i]!=0)
                factor=slau[j][i]/slau[i][i];
            else
                return 2;
            for(int k=i;k<size+1;k++){
                slau[j][k]+=(-factor)*slau[i][k];
            }
        }
    }
    // backward
    float answer[size];
    slau[size-1][size-1]==0?answer[size-1]=0:answer[size-1]=slau[size-1][size]/slau[size-1][size-1];

    float temp;
    for(int i=size-2;i>-1;i--){
        temp=slau[i][size];
        if(slau[i][size]!=0){
            for(int j=size-1;j>i;j--){
                temp-=slau[i][j]*answer[j];
            }
        cout<<endl;
        temp==0?answer[i]=0:answer[i]=slau[i][i]/temp;
        }
    }


    // output
    fstream ofile;
    ofile.open("slau.out",ios::out);
    for(int i=0;i<size;i++){
        for(int j=0;j<size+1;j++){
            j==3?ofile<<'|':ofile;
            ofile<<' '<<slau[i][j]<<' ';
        }
        ofile<<endl;
    }
    for(int i=0;i<size;i++)
        ofile<<"answer["<<i<<"]:"<<answer[i]<<endl;
    ofile<<"size: "<<size<<endl;
    ofile.close();

    return 0;
}

 

[Ваня]

Кошмар, как это все можно понимать?

 

[Fenix]

Вань, нормально

 

[Vladimir_S]

Вот здесь http://www.tehnari.ru/f41/t43929/#post435904 я выкладывал решение той же задачи на Паскале. Между прочим, вариант с нулевым диагональным членом там учтен и обходится путем перестановки строк.

 

[Fenix]

Владимир, на счет перестановки я в курсе, допишу позже, делов-то там. Это больше шпаргалки для меня.