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

Размещено на http://www.allbest.ru

Размещено на http://www.allbest.ru

Белорусский государственный университет

Факультет прикладной математики и информатики

Кафедра вычислительной математики

Лабораторная работа №1

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

Подготовил: студент 2 курса 2 группы

Зинькович И.А.

Преподаватель: Самусенко А.В.

Минск

2014

Теоретический материал

Пусть есть система:

Прямой ход состоит из n 1 шагов исключения.

1-й шаг. Целью этого шага является исключение неизвестного x₁ из уравнений с номерами i = 2, 3, …, n. Предположим, что коэффициент a₁₁ 0. Будем называть его главным элементом 1-го шага.

Найдем величины q_i1 = a_i1/a₁₁ (i = 2, 3, …, n), называемые множителями 1-го шага. Вычтем последовательно из второго, третьего, …, n-го уравнений системы первое уравнение, умноженное соответственно на q₂₁, q₃₁, …, q_n1. Это позволит обратить в нуль коэффициенты при x₁ во всех уравнениях, кроме первого. В результате получим эквивалентную систему

a₁₁x₁ + a₁₂x₂ + a₁₃x₃ + … + a_1nx_n = b_1,

a₂₂⁽¹⁾x₂ + a₂₃⁽¹⁾x₃ + … + a_2n⁽¹⁾x_n = b₂^(1),

a₃₂⁽¹⁾x₂ + a₃₃⁽¹⁾x₃ + … + a_3n⁽¹⁾x_n = b₃⁽¹⁾,

a_n2⁽¹⁾x₂ + a_n3⁽¹⁾x₃ + … + a_nn⁽¹⁾x_n = b_n⁽¹⁾.

в которой a_ij⁽¹⁾ и b_ij⁽¹⁾ вычисляются по формулам

a_ij⁽¹⁾ = a_ij ? q_i1a_1j, b_i⁽¹⁾ = b_i ? q_i1b₁.

2-й шаг. Целью этого шага является ислючение неизвестного x₂ из уравнений с номерами i = 3, 4, …, n. Пусть a₂₂⁽¹⁾ ? 0, где a₂₂⁽¹⁾ - коэффициент, называемый главным (или ведущим) элементом 2-го шага. Вычислим множители 2-го шага

q_i2 = a_i2⁽¹⁾ / a₂₂⁽¹⁾ (i = 3, 4, …, n)

и вычтем последовательно из третьего, четвертого, …, n-го уравнения системы второе уравнение, умноженное соответственно на q₃₂, q₄₂, …, q_m2. В результате получим систему

a₁₁x₁ + a₁₂x₂ + a₁₃x₃ + … + a_1nx_n = b₁ ,

a₂₂⁽¹⁾x₂ + a₂₃⁽¹⁾x₃ + … + a_2n⁽¹⁾ x_n = b₂⁽¹⁾ ,

a₃₃⁽²⁾x₃ + … + a_3n⁽²⁾x_n = b₃⁽²⁾ ,

a_n3⁽²⁾x₃ + … + a_nn⁽²⁾x_n = b_n⁽²⁾ .

Здесь коэффициенты a_ij⁽²⁾ и b_ij⁽²⁾ вычисляются по формулам: a_ij⁽²⁾ = a_ij⁽¹⁾ - q_i2a_2j⁽¹⁾, b_i⁽²⁾ = b_i⁽¹⁾ - q_i2b₂⁽¹⁾.

Аналогично проводятся остальные шаги. Опишем очередной k-й шаг.

k-й шаг. В предположении, что главный (ведущий) элемент k-го шага a_kk^(k-1) отличен от нуля, вычислим множители k-го шага q_ik = a_ik^(k-1) / a_kk^(k-1) (i = k + 1, …, n) и вычтем последовательно из (k + 1)-го, …, n-го уравнений полученной на предыдущем шаге системы k-e уравнение, умноженное соответственно на q_k+1,k, q_k+2,k, …, q_nk.

После (n - 1)-го шага исключения получим систему уравнений

a₁₁x₁ + a₁₂x₂ + a₁₃x₃ + … + a_1nx_n = b₁ ,

a₂₂⁽¹⁾x₂ + a₂₃⁽¹⁾x₃ + … + a_2n⁽¹⁾x_n = b₂⁽¹⁾ ,

a₃₃⁽²⁾x₃ + … + a_3n⁽²⁾x_n = b₃⁽²⁾ ,

a_nn^(n-1)x_n = b_n^(n-1) .

матрица A^(n-1) которой является верхней треугольной. На этом вычисления прямого хода заканчиваются.

Обратный ход. Из последнего уравнения системы находим xn. Подставляя найденное значение x_n в предпоследнее уравнение, получим x_n-1. Осуществляя обратную подстановку, далее последовательно находим x_n-1, x_n-2, …, x₁. Вычисления неизвестных здесь проводятся по формулам

x_n = b_n^(n-1) / a_nn^(n-1),

x_k = (b_n^(k-1) - a_k,k+1^(k-1)x_k+1 - … - a_kn^(k-1)x_n) / a_kk^(k-1), (k = n - 1, …, 1).

Необходимость выбора главных элементов. Заметим, что вычисление множителей, а также обратная подстановка требуют деления на главные элементы a_kk^(k-1). Поэтому если один из главных элементов оказывыется равным нулю, то схема единственного деления не может быть реализована. Здравый смысл подсказывает, что и в ситуации, когда все главные элементы отличны от нуля, но среди них есть близкие к нулю, возможен неконтролируемый рост погрешности.

линейный алгебраический уравнение гаусс

Условие задания

Решить методом Гаусса систему линейных уравнений Ax = b, где A = D + C*k. Матрицы D, C, b взяты из №41, k = 0,2.

Код программы (main.cpp)

#include <fstream>

#include <cmath>

using namespace std;

double round(double x) {

return (x < 0) ? ceil(x-0.5): floor(x+0.5);

}

int main(){

int SIZE;

double** A;

double** copy;

double* b;

double* result;

double* delta;

double k = 0.2;

ifstream in("in.txt");

in >> SIZE;

A = new double*[SIZE];

copy = new double*[SIZE];

for (int i = 0; i < SIZE; i++){

A[i] = new double[SIZE];

copy[i] = new double[SIZE];

}

double** D = new double*[SIZE];

for (int i = 0; i < SIZE; i++)

D[i] = new double[SIZE];

for (int i = 0; i < SIZE; i++)

for (int j = 0; j < SIZE; j++)

in >> D[i][j];

double** C = new double*[SIZE];

for (int i = 0; i < SIZE; i++)

C[i] = new double[SIZE];

for (int i = 0; i < SIZE; i++)

for (int j = 0; j < SIZE; j++)

in >> C[i][j];

for (int i = 0; i < SIZE; i++)

for (int j = 0; j < SIZE; j++){

A[i][j] = D[i][j] + C[i][j] * k;

copy[i][j] = A[i][j];

}

b = new double[SIZE];

result = new double[SIZE];

delta = new double[SIZE];

for (int i = 0; i < SIZE; i++){

in >> b[i];

result[i] = b[i];

}

for (int i = 0; i < SIZE; ++i) {

if (A[i][i]!= 0) {

result[i] /= A[i][i];

for (int j = SIZE - 1; j >= i; --j) {

A[i][j] /= A[i][i];

}

}

else {

return 0;

}

for (int j = i + 1; j < SIZE; ++j) {

for (int k = i; k < SIZE; k++) {

if (i == k)

continue;

A[j][k] -= A[j][i] * A[i][k];

}

result[j] -= result[i] * A[j][i];

A[j][i] = 0;

}

}

for (int j = SIZE - 1; j >= 0; j--){

for (int i = 0; i < j; i++){

result[i] -= A[i][j] * result[j];

A[i][j] = 0;

}

}

ofstream out("out.txt");

FILE* f = fopen("out.txt", "w");

for (int i = 0; i < SIZE; i++){

result[i] = round(result[i] * 100000) / 100000;

}

for (int i = 0; i < SIZE; i++)

fprintf(f, "%.5lf\n", result[i]);

for (int i = 0; i < SIZE; i++)

for (int j = 0; j < SIZE; j++) {

double sum = 0;

for (int k = 0; k < SIZE; k++)

sum += copy[i][k] * result[k];

delta[i] = abs(b[i] - sum);

}

for (int i = 0; i < SIZE; i++)

fprintf(f, "\n%.10lf", delta[i]);

fclose(f);

}

В результате получили X:

0.9169573287

0.2508817999

-1.1237128647

0.0277923094

1.0872656987

-0.2322239583

Округляя до пяти знаков после запятой:

0.91696

0.25088

-1.12371

0.02779

1.08727

-0.23222

При вычислении AX - B, где X - округлённые данные, получим:

0.0000061000

0.0000061000

0.0000057000

0.0000022000

0.0000035000

0.0000031000

Размещено на Allbest.ru