Метод Гаусса с выбором главного элемента
Может оказаться так, что система (1.1) имеет единственное решение, хотя какой либо из миноров матрицы А равен нулю. Заранее неизвестно, что все угловые миноры матрицы А не равны нулю. В этом случае можно использовать метод Гаусса с выбором главного элемента.
1. Выбор главного элемента по столбцу, когда на k-ом шаге исключения в качестве главного элемента выбирают максимальный по модулю коэффициент при неизвестном хk в уравнениях с номерами i=k,k+1,…,m. Затем уравнение, соответствующее выбранному коэффициенту, меняют местами с k-ым уравнением системы, чтобы ведущий элемент занял место коэффициента .После перестановки исключение неизвестного хk выполняют как в схеме единственного деления.
ПРИМЕР. Пусть дана система второго порядка
Предположим, что½а21½>½а11½, тогда переставим уравнения
и применяем первый шаг прямого хода метода Гаусса. В этом случае имеет место перенумерация строк.
2. Выбор главного элемента по строке, т.е. производится перенумерация неизвестных системы.
При ½а12½>½ а11½, на первом шаге вместо неизвестного х1 исключают х2:
.
К этой системе применяем первый шаг прямого хода метода Гаусса.
3. Поиск главного элемента по всей матрице заключается в совместном применении методов 1 и 2. Всё это приводит к уменьшению вычислительной погрешности, но может замедлить процесс решения задачи.
Метод Холецкого (метод квадратных корней)
Пусть дана система
Ах=b, (1.11)
где А − симметричная положительно определенная матрица.
Тогда решение системы (1.11) проводится в два этапа:
1. Симметричная матрица А представляется как произведение двух матриц
А = L∙ LТ.
Рассмотрим метод квадратных корней на примере системы 4-го порядка:
.
Перемножаем матрицы в правой части разложения и сравниваем с элементами в левой части:
, , , , ,
, , .
2. Решаем последовательно две системы
Ly=b,
LTx=у.
Замечания
1) Под квадратным корнем может получиться отрицательное число, следовательно в программе необходимо предусмотреть использование правил действия с комплексными числами.
2) Возможно переполнение, если угловые элементы близки к нулю.
Практическая часть
Задание на контрольную работу
Тема контрольной работы «Решение систем линейных алгебраических уравнений точными методами»
1. Метод Гаусса.
1.1 Решение системы линейных алгебраических уравнений методом Гаусса.
1.2 Решение системы линейных алгебраических уравнений методом Гаусса с помощью программы PROCEDURE SIMQ.
2. Метод квадратных корней (Холецкого)
2.1 Решение системы линейных алгебраических уравнений методом квадратных корней (Холецкого).
2.2 Решение системы линейных алгебраических уравнений методом квадратных корней (Холецкого) с помощью программы Procedure holet.
Программы можно реализовывать в среде Free Pascal или как проект консольного приложения в среде Delphi.
Варианты заданий
Номер варианта задания выбирается по порядковому номеру списка группы.
Таблица 1
Номер варианта | Матрица А коэффициентов системы | Столбец свободных членов b | ||||
1,84 | 2,25 | 2,53 | -6,09 | |||
2,32 | 2,60 | 2,82 | -6,98 | |||
1,83 | 2,06 | 2,24 | -5,52 | |||
2,58 | 2,93 | 3,13 | -6,66 | |||
1,32 | 1,55 | 1,58 | -3,58 | |||
2,09 | 2,25 | 2,34 | -5,01 | |||
2,18 | 2,44 | 2,49 | -4,34 | |||
2,17 | 2,31 | 2,49 | -3,91 | |||
3,15 | 3,22 | 3,17 | -5,27 | |||
1,54 | 1,70 | 1,62 | -1,97 | |||
3,69 | 3,73 | 3,59 | -3,74 | |||
2,45 | 2,43 | 2,25 | -2,26 | |||
Продолжение таблицы 1 | ||||||
1 2 3 | ||||||
1,53 | 1,61 | 1,43 | -5,13 | |||
2,35 | 2,31 | 2,07 | -3,69 | |||
3,83 | 3,73 | 3,45 | -5,98 | |||
2,36 | 2,37 | 2,13 | 1,48 | |||
2,51 | 2,40 | 2,10 | 1,92 | |||
2,59 | 2,41 | 2,06 | 2,16 | |||
3,43 | 3,38 | 3,09 | 5,52 | |||
4,17 | 4,00 | 3,65 | 6,93 | |||
4,30 | 4,10 | 3,67 | 7,29 | |||
3,88 | 3,78 | 3,45 | 10,41 | |||
3,00 | 2,79 | 2,39 | 8,36 | |||
2,67 | 2,37 | 1,96 | 7,62 | |||
3,40 | 3,26 | 2,90 | 13,05 | |||
2,64 | 2,39 | 1,96 | 10,30 | |||
4,64 | 4,32 | 3,85 | 17,89 | |||
2,53 | 2,36 | 1,93 | 12,66 | |||
3,95 | 4,11 | 3,66 | 21,97 | |||
2,78 | 2,43 | 1,94 | 13,93 | |||
2,16 | 1,96 | 1,56 | 13,16 | |||
3,55 | 3,23 | 2,78 | 21,73 | |||
4,85 | 4,47 | 3,97 | 29,75 | |||
2,69 | 2,47 | 2,07 | 19,37 | |||
2,73 | 2,39 | 1,92 | 19,43 | |||
2.93 | 2,52 | 2,02 | 20,80 | |||
3,72 | 3,47 | 3,06 | 30,74 | |||
4,47 | 4,10 | 3,63 | 36,80 | |||
4,96 | 4,53 | 4,01 | 40,79 | |||
4,35 | 4,39 | 3,67 | 40,15 | |||
4,04 | 3,65 | 3,17 | 36,82 | |||
3,14 | 2,69 | 2,17 | 28,10 | |||
4,07 | 3,79 | 3,37 | 40,77 | |||
2,84 | 2,44 | 1,95 | 27,68 | |||
4,99 | 4,50 | 3,97 | 49,37 | |||
3,19 | 2,89 | 2,47 | 33,91 | |||
4,43 | 4,02 | 3,53 | 47,21 | |||
3,40 | 2,92 | 2,40 | 32,92 | |||
2,57 | 2,26 | 1,84 | 28,66 | |||
4,47 | 4,03 | 3,57 | 50,27 | |||
4,89 | 4,40 | 3,87 | 55,03 | |||
Продолжение таблицы 1 | ||||||
2,83 | 2,50 | 2,08 | 33,28 | |||
3,00 | 2,55 | 2,07 | 33,59 | |||
3,72 | 3,21 | 2,68 | 43,43 | |||
3,78 | 3,44 | 3,02 | 46,81 | |||
4,33 | 3,88 | 3,39 | 53,43 | |||
4,76 | 4,24 | 3,71 | 58,73 | |||
4,59 | 4,24 | 3,82 | 59,54 | |||
4,83 | 4,36 | 3,88 | 62,33 | |||
4,06 | 3,53 | 3,01 | 52,11 | |||
4,56 | 4,20 | 3,78 | 61,86 | |||
3,21 | 2,73 | 2,25 | 42,98 | |||
4,58 | 4,04 | 3,52 | 61,67 | |||
3,75 | 3,39 | 2,97 | 53,38 | |||
4,18 | 3,70 | 3,22 | 59,28 | |||
4,43 | 3,88 | 3,36 | 62,62 | |||
2,95 | 2,58 | 2,16 | 44,16 | |||
5,11 | 4,62 | 4,14 | 46,68 | |||
4,38 | 3,82 | 3,30 | 65,34 | |||
2,93 | 2,55 | 2,14 | 46,41 | |||
3,47 | 2,98 | 2,50 | 54,78 | |||
4,78 | 4,22 | 3,70 | 75,81 | |||
3,74 | 3,36 | 2,94 | 63,26 | |||
4,02 | 3,51 | 3,04 | 67,51 | |||
4,18 | 3,61 | 3,09 | 70,03 | |||
4,07 | 4,28 | 3,87 | 84,43 | |||
5,30 | 4,79 | 4,32 | 95,45 | |||
5,11 | 4,54 | 4,03 | 91,69 | |||
4,90 | 4,50 | 4,09 | 94,18 | |||
3,79 | 3,27 | 2,81 | 71,57 | |||
4,01 | 3,43 | 2,91 | 75,45 | |||
4,25 | 3,84 | 3,43 | 86,07 | |||
3,86 | 3,34 | 2,87 | 77,12 | |||
5,40 | 4,82 | 4,30 | 108,97 | |||
3,35 | 2,94 | 2,53 | 70,69 | |||
5,41 | 4,88 | 4,41 | 115,38 | |||
3,88 | 3,30 | 2,78 | 81,07 | |||
3,05 | 2,64 | 2,23 | 67,17 | |||
4,14 | 3,61 | 3,14 | 91,43 | |||
5,63 | 5,03 | 4,52 | 125,40 | |||
Пример решения задачи
2.3.1 Решение системы линейных алгебраических уравнений методом Гаусса.
Решить СЛАУ методом Гаусса.
1. Прямой ход
Разделим все члены первого уравнения на 8, получим:
.
Вычтем из второго и третьего уравнения первое:
.
Разделим все члены второго уравнения на , а третьего на , получим :
.
Вычтем из третьего уравнения второе :
.
Прямой ход завершен.
2.Обратный ход.
Выразим из последнего уравнения системы :
.
Из второго уравнения выразим :
.
Из первого уравнения выразим :
.
Получаем следующее решение системы:
.
2.3.2 Решение системы линейных алгебраических уравнений методом Гаусса с помощью программы PROCEDURE SIMQ
Procedure SIMQ(Nn:Integer;Var Aa:TMatr;
Var Bb:TVector;Var Ks:Integer);
Входные параметры:
Nn − размерность системы;
Aa − матрица коэффициентов системы N *N типа TMatr = array[1..n, 1..n] of real;
Bb− вектор размерности N типа TVector= array[1..n] of real, содержит правые части системы.
Выходные параметры:
Bb − решение системы;
ks – признак правильности решения (код ошибки): если ks=0 то в массиве Bb содержится решение системы.
Пример. Решить систему уравнений
Схема алгоритма приведена на рисунке 1.
Примерный текст процедуры:
PROCEDURE SIMQ(Nn:Integer;Var Aa:TMatr;Var Bb:TVector;Var Ks:Integer);
Label M1;
Const Eps=1e-21;
Var Max,U,V : Real; I,J,K1,L : Integer;
Begin
For I:=1 To Nn Do Aa[i,Nn+1]:=Bb[i];
For I:=1 To Nn Do
Begin
Max:=Abs(Aa[i,i]); K1:=I;
For L:=I+1 To Nn Do If (Abs(Aa[l,i])>Max) Then
Begin
Max:=Abs(Aa[l,i]); K1:=L;End;
If(Max<Eps) Then Begin Ks:=1; Goto M1;
End Else Ks:=0;
If K1<>I Then
For J:=I To Nn+1 Do
Begin U:=Aa[i,j]; Aa[i,j]:=Aa[k1,j]; Aa[k1,j]:=U; End;
V:=Aa[i,i];
For J:=I To Nn+1 Do Aa[i,j]:=Aa[i,j]/V;
For L:=i+1 To Nn Do Begin
V:=Aa[l,i]; For J:=I+1 To Nn+1 Do Aa[l,j]:=Aa[l,j]-Aa[i,j]*V;
End; End;
Bb[nn]:=Aa[Nn,Nn+1];
For I:=Nn-1 Downto 1 Do Begin Bb[i]:=Aa[i,nn+1];
For J:=I+1 To Nn Do Bb[i]:=Bb[i]-Aa[i,j]*Bb[j];
End;
M1:End;
Вычисления по программе привели к следующим результатам:
x1=0.100000E+0001 x2=0.200000Е+0001
x3=0.З00000Е+0001
признак выхода 0
Рисунок 1 – Схема алгоритма метода Гаусса
Варианты заданий для решения систем линейных алгебраических уравнений методом Гаусса приведены в таблице 1.
2.3.3 Решение системы линейных алгебраических уравнений метод квадратных корней (Холецкого)
Решить СЛАУ методом квадратных корней.
1. Пусть
,
Находим элементы матрицы L:
, ,
, , .
Таким образом разложение матрицы А имеет вид:
2. Последовательно решаем две системы Ly=b и LTx=y.
Решением первой системы является вектор ,
а решением второй системы вектор .
2.3.4 Решение системы линейных алгебраических уравнений метод квадратных корней (Холецкого) с помощью программы Procedure holet
Procedure holet (A1: Tmatr; var B1: Tvector;
N: integer);
Входные параметры:
N − размерность системы;
A1− матрица коэффициентов системы N *N +1, типа Tmatr = array[1..n, 1..n+1] of real, N +1 столбец используется для копирования левой части для упрощения вычислений;
B1− вектор свободных членов системы типа Tvector= array[1..n] of real.
Выходные параметры:
B1− решение системы.
Пример. Решить систему уравнений:
Схема алгоритма приведена на рисунке 2.
Здесь произвольная матрица коэффициентов при неизвестных в системе уравнений умножением на транспонированную к ней матрицу приводится к симметричной.
Примерный текст программы:
uses crt;
type
Tmatr=array[1..3,1..4]of real;
Tvector=array[1..3]of real;
const
{матрица размерности N*N , N+1 столбец используется для копирования
левой части СЛАУ для упрощения вычислений}
A:tmatr= ((1,1,1,0),
(1,0,-1,0),
(1,2,1,0));
B:Tvector= (6,-2,8);
var
k,N:integer;
{процедура для вывода матрицы}
Procedure vivod_matr(A1:Tmatr;N :integer);
var
i,j:integer;
begin
For i:=1 to N do begin
For j:=1 to N do write(A1[i,j]:15,' ');
writeln;
End;
writeln;
End;
{процедура для вывода вектора}
Procedure vivod_vectr(B1:Tvector; N :integer);
var
j:integer;
begin
For j:=1 to N do writeln('x',j,'=',B1[j]:0,' ');
writeln;
End;
{метод Холецкого}
{
Входные параметры:
N − размерность системы;
A1− матрица коэффициентов системы N *N +1, N +1 столбец используется для копирования левой части СЛАУ для упрощения вычислений;
B1− вектор свободных членов системы.
Выходные параметры:
B1− решение системы.
B1-при вызове процедуры содержит вектор свободных членов СЛАУ,
Выходные параметры:
B1- содержит решение СЛАУ
}
Procedure holet(A1:Tmatr;var B1:Tvector;N :integer);
var i,j,t:integer;
c,L:Tmatr;
y:Tvector;
buf,summ:real;
Begin
For i:=1 to N do
For j:=1 to N +1 do
begin
c[i,j]:=0;
L[i,j]:=0;
y[i]:=0;
End;
{Умножение матрицы на транспонированную}
For i:=1 to N do
For j:=1 to N do
begin
summ:=0.0;
For t:=1 to N do
summ:=A1[t,j]*A1[t,i]+summ;
c[i,j]:=summ;
End;
{умножение правой части на транспонированную м-цу}
For i:=1 to N do
For j:=1 to N do
y[i]:=A1[j,i]*B1[j]+y[i];
For i:=1 to N do For j:=1 to N do
begin
A1[i,j]:=c[i,j];
B1[i]:=y[i];
End;
For i:=1 to N do
begin
For j:=1 to i do
begin
summ:=0;
For t:=1 to j-1 do
summ:=summ+L[i,t]*L[j,t];
if i<>j then
L[i,j]:=(A1[i,j]-summ)/L[j,j]
else
L[i,i]:=sqrt(A1[i,i]-summ);
End;
End;
For i:=1 to N do L[i,N +1]:=B1[i];
B1[1]:=L[1,N +1]/L[1,1];
For i:=2 to N do
begin
For j:=1 to i-1 do
L[i,N +1]:=L[i,N +1]-L[i,j]*B1[j];
B1[i]:=L[i,N +1]/L[i,i];
End;
For i:=1 to N do
begin
For j:=i+1 to N do
begin
L[i,j]:=L[j,i];
L[j,i]:=0;
End;
L[i,N +1]:=B1[i];
End;
B1[N ]:=L[N ,N +1]/L[N ,N ];
For i:=N -1 downto 1 do
begin
For j:=i+1 to N do
L[i,N +1]:=L[i,N +1]-L[i,j]*B1[j];
B1[i]:=L[i,N +1]/L[i,i];
End;
End;
begin
N:=3;
writeln('dannaya SLAU');
{матрица задается типизированной константой}
vivod_matr(A,N);
writeln('metod Holetckogo');
holet(A,B,N);
vivod_vectr (B, N);
End.
Вычисления по программе привели к следующим результатам:
x1= 0.1000000E+0001
x2= 0.2000000E+0001
x3= 0.3000000E+0001
Рисунок 2 – Схема алгоритма метода Холецкого