Решение систем линейных уравнений методом Гаусса

Алгоритмизация / Решение систем линейных уравнений методом Гаусса

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

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

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

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

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

Иллюстрирующий пример

Пусть дана система уравнений
Система линейных уравнений (пример)
или в матричной форме
Система линейных уравнений в матричной форме (пример)
Выбираем строку с максимальным коэффициентом ai1 и меняем ее с первой.
Перемещение строк в матрице коэффициентов
Нормируем уравнения относительно коэффициента при x1:
Нормирование коэффициентов (пример)
Нормирование уравнений
Вычитаем 1 уравнение из 2 и 3:
Вычитание уравнений
Выбираем строку с наибольшим коэффициентом при ai2 (уравнение 1 не рассматривается) и перемещаем ее на место 2.
Перемещение строк коэффициентов (пример)
Нормируем 2 и 3 уравнения относительно коэффициента при x2

gauss15
Вычитаем уравнение 2 из 3
gauss16
Нормируем уравнение 3 относительно коэффициента при x3
gauss17
Откуда получаем x3=2. Подставляем полученное значение в уравнения 2 и 1 получаем
gauss18
Подставляя полученное значение x2=5 в уравнение 1, найдем
gauss19
Таким образом, решением системы уравнений будет вектор
Вектор решения системы линейных уравнений
Реализация

#include <iostream>
using namespace std;
// Вывод системы уравнений
void sysout(double **a, double *y, int n) {
  for (int i = 0; i < n; i++) {
    for (int j = 0; j < n; j++){
      cout << a[i][j] << "*x" << j;
      if (j < n - 1) {
        cout << " + ";
      }
    }
  cout << " = " << y[i] << endl;
  }
  return;
}
double * gauss(double **a, double *y, int n) {
  double *x, max;
  int k, index;
  const double eps = 0.00001;  // точность
  x = new double[n];
  k = 0;
  while (k < n) {
  // Поиск строки с максимальным a[i][k]
    max = abs(a[k][k]);
    index = k;
    for (int i = k + 1; i < n; i++) {
      if (abs(a[i][k]) > max) {
        max = abs(a[i][k]);
        index = i;
      }
    }
  // Перестановка строк
    if (max < eps) {
    // нет ненулевых диагональных элементов
      cout << "Решение получить невозможно из-за нулевого столбца " ;
      cout << index << " матрицы A" << endl;
      return 0;
    }
    for (int j = 0; j < n; j++) {
      double temp = a[k][j];
      a[k][j] = a[index][j];
      a[index][j] = temp;
    }
    double temp = y[k];
    y[k] = y[index];
    y[index] = temp;
  // Нормализация уравнений
    for (int i = k; i < n; i++) {
      double temp = a[i][k];
      if (abs(temp) < eps) continue; // для нулевого коэффициента пропустить
      for (int j = 0; j < n; j++) {
        a[i][j] = a[i][j] / temp;
      }
      y[i] = y[i] / temp;
      if (i == k)  continue; // уравнение не вычитать само из себя
      for (int j = 0; j < n; j++) {
        a[i][j] = a[i][j] - a[k][j];
      }
      y[i] = y[i] - y[k];
    }
    k++;
  }
  // обратная подстановка
  for (k = n - 1; k >= 0; k--) {
    x[k] = y[k];
    for (int i = 0; i < k; i++) {
      y[i] = y[i] - a[i][k] * x[k];
    }
  }
  return x;
}
int main() {
  double **a, *y, *x;
  int n;
  system("chcp 1251");
  system("cls");
  cout << "Введите количество уравнений: ";
  cin >> n;
  a = new double*[n];
  y = new double[n];
  for (int i = 0; i < n; i++) {
    a[i] = new double[n];
    for (int j = 0; j < n; j++) {
      cout << "a[" << i << "][" << j << "]= ";
      cin >> a[i][j];
    }
  }
  for (int i = 0; i < n; i++) {
    cout << "y[" << i << "]= ";
    cin >> y[i];
  }
  sysout(a, y, n);
  x = gauss(a, y, n);
  for (int i = 0; i < n; i++){
    cout << "x[" << i << "]=" << x[i] << endl;
  }
  cin.get(); cin.get();
  return 0;
}

Результат выполнения
Пример
Назад

Комментариев к записи: 33


    • Елена Вставская

      Да, есть такая недоработка. Проверка совместности системы в примере отсутствует.


  • Вместо
    const double eps = 0.00001;  // точность
    можно прописать
    const double eps = DBL_EPSILON;  // точность
    // DBL_EPSILON определено во float.h


    • Елена Вставская

      Да, но для одной константы придется вставлять в проект целую библиотеку float.h.


      • Это не библиотека, нет идущих вместе никаких исполняемых модулей, это просто определения констант, точнее макросов препроцессора, т.е. вместо вашего 0.00001 просто подставится число, в итоге даже размер исполняемого файла не изменится.



  • Алексей

    пишет [Error] ‘abs’ was not declared in this scope, объявляю и появляется куча ошибок, можете сказать как поправить?


    • Елена Вставская

      Попробуйте подключить <cmath> или воспользоваться функцией fabs().


  • Дмитрий

    А есть ли какой-то способ решить систему из 2х уровнений, но с 3мя неизвестными ? К примеру, если нам дана сумма их a+b+c=const(s) и const(1)*a+const(2)*b+const(3)*c = const(s)*const(4) ? Ну или, что-то подобное, Спасибо.


    • Дмитрий, для того чтобы система имела единственное решение необходимым условием является одинаковое количество уравнений и неизвестных. В случае если количество неизвестных больше, чем количество уравнений, система будет иметь бесконечное множество решений.


      • Добрый день,

        а можно ли перепрограммировать текущий код, так чтобы программа могла также решать системы с неодинаковым количеством переменных и количеством строк, как например сделали на этом сайте:

        http://math.semestr.ru/gauss/gauss.php

        ?

        С Уважением,
        Игорь.


        • Елена Вставская

          Если количество переменных больше количества уравнений, то часть переменных, превышающих количество уравнений, задается в качестве свободных.
          Если количество уравнений превышает количество переменных, то такая система вообще может не иметь решения в случае если она не является вырожденной.

          Пример вырожденной системы уравнений:

          x1 + 2*x2 = 3

          3*x1 + 6*x2 = 9


  • Здравствуйте! Спасибо за код. У меня вылазят ошибки sysout: идентификатор не найден

    и gauss: идентификатор не найден

    хотя все сделал, как у вас. В чем может быть хитрость?



      • Ответить без Вашего кода довольно сложно. У меня приведенный пример компилируется и работает.


  • Анастасия

    Можно еще вопрос? в моей системе программа не выдает 0.
    т.е. у меня система:

    10×1-7×2=7;
    -3×1+2×2+6×3=4;
    5×1-x2+5×3=6.

    ответ, посчитанный вручную :

    x1=0;
    x2=-1;
    x3=1.

    программа выдает:

    x1=1.11022e-016;
    x2=-1;
    x3=1.

    почему так может быть?


  • Добрый день,

    после x = gauss(a, y, n); в главной функции main, нужно добавить проверку исключения, а именно:

    sysout(a, y, n);
    x = gauss(a, y, n);
    if (x != 0)
    {
    for (int i = 0; i < n; i++) {
    cout << "x[" << i << "]=" << x[i] << endl;
    }

    еще меня вот что смутило, например есть такая система уравнений: (в скобочках я имею в виду нижний регистр)

     

    x(2)+3x(3)-x(4)=2

    2*x(1)-3*x(2) +2*x(3)=-2

    2*x(1)-4*x(2)+x(4)=3

    -2*x(1)+5*x(2)-3*x(3)+3*x(4)=1

    результат вычисления на бумаге у меня получился такой:
    x1=23,3

    x2=-10,2

    x3=7

    x4=8.8

    в программе такой:

    x[0]=-2.25833
    x[1]=-1.25
    x[2]=1.38333
    x[3]=2.8

    где может быть ошибка, давайте вместе разберем?


    • Спасибо за замечание, Игорь!

      Ошибку исправила, сейчас Ваша система решается верно.


  • Добрый день,

    а если бы, например, было 4 уравнения с 5-ю неизвестными (x1,x2,x3,x4,x5), то это можно запрограммировать?

    (так как вроде по формулировке — бесконечное множество решений)

    С Уважением,
    Игорь.


    • Да, так и есть. Количество уравнений должно соответствовать количеству неизвестных, причем система должна быть невырожденной, чтобы она имела единственное решение.


  • Добрый день,

    почему после нормирования 2 и 3 уравнения относительно коэффициента при x2 в первой строке вместо [1 0,4 0,2] * [x1] = [9,4] стало [1 0,5 0,25]*[x1] = [10]?

    С Уважением,
    Игорь.







        • Для однозначного решения количество уравнений должно соответствовать количеству неизвестных. Поэтому для прямоугольных матриц однозначного решения не существует.



          • Все равно не поняла вопроса. Какую задачу требуется решить?


          • Z2   
            Решение СЛАУ по модулю, по модулю 2
            Z2 

            Z2 


          • Тут метод Гаусса вряд ли поможет. Это — тема для отдельной статьи.


Добавить комментарий

Ваш e-mail не будет опубликован. Обязательные поля помечены *