Метод исключения Гаусса

Алгоритм

Метод Гаусса - классический метод решения системы линейных алгебраических уравнений (СЛАУ). Рассмотрим систему линейных уравнений с действительными постоянными коэффициентами

\begin{align*} \left\{ \begin{array}{ccc} a_{11} \cdot x_{1} + a_{12} \cdot x_{2} + & \cdots & + a_{1 n} \cdot x_{n}=y_{1} \\ a_{21} \cdot x_{1} + a_{22} \cdot x_{2} + & \cdots & + a_{2 n} \cdot x_{n}=y_{2} \\ & \cdots & \\ a_{n 1} \cdot x_{1} + a_{n 2} \cdot x_{2} + & \cdots & + a_{n n} \cdot x_{n}=y_{n} \end{array} \right. \end{align*}

или в матричной форме

\begin{align*} A x = y, \end{align*}

где

\begin{align*} A=\left( \begin{array}{ccc} a_{11} & \dots & a_{1 n} \\ & \dots & \\ a_{n 1} & \dots & a_{n n} \end{array} \right), \quad x=\left( \begin{array}{c} x_{1} \\ \vdots \\ x_{n} \end{array} \right), \quad y=\left( \begin{array}{c} y_{1} \\ \vdots \\ y_{m} \end{array} \right). \end{align*}

Метод Гаусса решения системы линейных уравнений включает в себя 2 стадии:

  • последовательное (прямое) исключение;
  • обратная подстановка.

Последовательное исключение

Исключения Гаусса основаны на идее последовательного исключения переменных по одной до тех пор, пока не останется только одно уравнение с одной переменной в левой части. Затем это уравнение решается относительно единственной переменной. Таким образом, систему уравнений приводят к треугольной (ступенчатой) форме. Для этого среди элементов первого столбца матрицы выбирают ненулевой (а чаще максимальный) элемент и перемещают его на крайнее верхнее положение перестановкой строк. Затем нормируют все уравнения, разделив его на коэффициент ai1, где i– номер столбца.

\begin{align*} \left\{ \begin{array}{ccc} x_{1}+\frac{a_{12}}{a_{11}} \cdot x_{2} + & \cdots & + \frac{a_{1 n}}{a_{11}} \cdot x_{n}=\frac{y_{1}}{a_{11}} \\ x_{1}+\frac{a_{22}}{a_{21}} \cdot x_{2} + & \cdots & + \frac{a_{2 n}}{a_{21}} \cdot x_{n}=\frac{y_{2}}{a_{21}} \\ & \cdots & \\ x_{1}+\frac{a_{n 2}}{a_{n 1}} \cdot x_{2} + & \cdots & + \frac{a_{n n}}{a_{n 1}} \cdot x_{n}=\frac{y_{n}}{a_{n 1}} \end{array} \right. \end{align*}

Затем вычитают получившуюся после перестановки первую строку из остальных строк:

\begin{align*} \left\{ \begin{array}{r c r c r c l} x_{1} & + & \frac{a_{12}}{a_{11}} \cdot x_{2} & + \cdots + & \frac{a_{1 n}}{a_{11}} \cdot x_{n} & = & \frac{y_{1}}{a_{11}} \\ 0 & + & \left(\frac{a_{22}}{a_{21}}-\frac{a_{12}}{a_{11}}\right) \cdot x_{2} & + \cdots + & \left(\frac{a_{2 n}}{a_{21}}-\frac{a_{1 n}}{a_{11}}\right) \cdot x_{n} & = & \left(\frac{y_{2}}{a_{21}}-\frac{y_{1}}{a_{11}}\right) \\ && & \cdots & \\ 0 & + & \left(\frac{a_{n 2}}{a_{n 1}}-\frac{a_{12}}{a_{11}}\right) \cdot x_{2} & + \cdots + & \left(\frac{a_{n n}}{a_{n 1}}-\frac{a_{1 n}}{a_{11}}\right) \cdot x_{n} & = & \left(\frac{y_{n}}{a_{n 1}}-\frac{y_{1}}{a_{11}}\right) \end{array} \right. \end{align*}

Получают новую систему уравнений, в которой заменены соответствующие коэффициенты.

\begin{align*} \left\{ \begin{array}{r c r c r c r} x_{1} & + & a_{12}^{\prime} \cdot x_{2} & + \cdots + & a_{1 n}^{\prime} \cdot x_{n} & = & y_{1}^{\prime} \\ 0 & + & a_{22}^{\prime} \cdot x_{2} & + \cdots + & a_{2 n}^{\prime} \cdot x_{n} & = & y_{2}^{\prime} \\ && & \cdots & \\ 0 & + & a_{n 2}^{\prime} \cdot x_{2} & + \cdots + & a_{n n}^{\prime} \cdot x_{n} & = & y_{n}^{\prime} \end{array} \right. \end{align*}

После того, как указанные преобразования были совершены, первую строку и первый столбец мысленно вычёркивают и продолжают указанный процесс для всех последующих уравнений пока не останется уравнение с одной неизвестной:

\begin{align*} \left\{ \begin{array}{r c r c r c r c r} x_{1} & + & a_{12}^{\prime} \cdot x_{2} & + & a_{13}^{\prime} \cdot x_{3} & + \cdots + & a_{1 n}^{\prime} \cdot x_{n} & = & y_{1}^{\prime} \\ 0 & + & x_{2} & + & a_{23}^{\prime \prime} \cdot x_{3} & + \cdots + & a_{2 n}^{\prime \prime} \cdot x_{n} & = & y_{2}^{\prime \prime} \\ 0 & + & 0 & + & x_{3} & + \cdots + & a_{3 n}^{\prime \prime \prime} \cdot x_{n} & = & y_{3}^{\prime \prime \prime} \\ && && & \cdots & \\ 0 & + & 0 & + & 0 & + \cdots + & x_{n} & = & y_{n}^{n \prime} \end{array} \right. \end{align*}

Обратная подстановка

Обратная подстановка предполагает подстановку полученного на предыдущем шаге значения переменной \(xn\) в предыдущие уравнения:

\begin{align*} \begin{array}{rcl} x_{n-1} & = & y_{n-1}^{(n-1)^{\prime}}-a_{(n-1) n}^{(n-1) \prime} \cdot x_{n} \\ x_{n-2}+a_{(n-2)(n-1)}^{(n-2) \prime} \cdot x_{n-1} & = & y_{n-2}^{(n-2)^{\prime}}-a_{(n-2) n}^{(n-2) \prime} \cdot x_{n} \\ \cdots \\ x_{2}+a_{23}^{\prime \prime} \cdot x_{3}+\cdots+a_{2(n-1)}^{\prime \prime} \cdot x_{n-1} & = & y_{2}^{\prime \prime}-a_{2 n}^{\prime \prime} \cdot x_{n} \\ x_{1}+a_{12}^{\prime} \cdot x_{2}+a_{13}^{\prime} \cdot x_{3}+\cdots+a_{1(n-1)}^{\prime} \cdot x_{n-1} & = & y_{1}^{\prime}-a_{1 n}^{\prime} \cdot x_{n} \end{array} \end{align*}

Эта процедура повторяется для всех оставшихся решений:

\begin{align*} \begin{array}{rcl} x_{n-2} & = & \left(y_{n-2}^{(n-2)^{\prime}}-a_{(n-2) n}^{(n-2)^{\prime}} \cdot x_{n}\right)-a_{(n-2)(n-1)}^{(n-2) \prime} \cdot x_{n-1} \\ & \cdots & \\ x_{2}+a_{23}^{\prime \prime} \cdot x_{3}+\cdots & = & \left(y_{2}^{\prime \prime}-a_{2 n}^{\prime \prime} \cdot x_{n}\right)-a_{2(n-1)}^{\prime \prime} \cdot x_{n-1} \\ x_{1}+a_{12}^{\prime} \cdot x_{2}+a_{13}^{\prime} \cdot x_{3}+\cdots & = & \left(y_{1}^{\prime}-a_{1 n}^{\prime} \cdot x_{n}\right)-a_{1(n-1)}^{\prime} \cdot x_{n-1} \end{array} \end{align*}

Использование

Представьте, что мы хотим решить следующую линейную систему уравнений:

\begin{align*} \left( \begin{array}{cc} 2 & 1 \\ -1 & 1 \end{array} \right) \left( \begin{array}{c} x_{1} \\ x_{2} \end{array} \right) = \left( \begin{array}{c} 5 \\ 2 \end{array} \right). \end{align*}

Тогда код будет выглядеть так:

// example_gauss_elimination.cpp

#include <iostream>
#include "../src/numerary.hpp" // Numerary library

using namespace std;
using namespace numerary;

/* The main function */
int main() {

    double **a = new double*[2];
    double *y = new double[2];
    double *x = new double[2];
    short int is_solved;

    for (int i = 0; i < 2; i ++)
        a[i] = new double[2];


    // Initialize matrix A
    a[0][0] = 2;
    a[0][1] = 1;

    a[1][0] = -1;
    a[1][1] = 1;

    // Initialize vector y
    y[0] = 5;
    y[1] = 2;

    is_solved = Numerary::linear_systems_of_equations(a, y, x, 2, "gauss");

    if (is_solved == 1) {
        cout << "System solved!" << endl;
        cout << "x = (" << x[0] << ", " << x[1] << ")" << endl;
    } else {
        cout << "Method is not allowed!";
    }

    for (int i = 0; i < 2; i++) delete[] a[i];

    delete[] a;
    delete[] x;
    delete[] y;

    return 0;
}