3.2. Некоторые точные методы решения слау

Метод Гаусса. Если число неизвестных менее 200, то для решения СЛАУ используют метод Гаусса, иначе называемый методом последовательного исключения неизвестных. Суть его заключается в приведении матрицы коэффициентов к треугольному виду. Вычисление состоит из двух этапов – прямого хода и обратного хода (обратной подстановки). Прямой ход выражается в последовательном исключении неизвестных из системы (3.1) для преобразования ее к эквивалентной системе с верхней треугольной матрицей. Вычисление значений неизвестных производят на этапе обратного хода.

Прямой ход состоит из -го шага.

1-й шаг. Целью этого шага является исключение неизвестного из уравнений с номерами. Предположим, что коэффициент(это ведущий элемент 1-го шага), вычислим величины

.

Вычтем последовательно из второго, третьего, ..., -го уравнений системы (3.1) первое уравнение, умноженное соответственно на,,...,. Это позволит обратить в нуль коэффициенты приво всех уравнениях, кроме первого. В результате получим эквивалентную систему

(3.4)

где и() находятся по формулам

, .

Фактически проведенные вычисления означают умножение матрицы на некоторую матрицуследующего вида:

.

Таким образом, исходное уравнение можно записать при помощи матрицы и векторовв виде

.

2-й шаг. Целью этого шага является исключение неизвестного из уравнений с номерами. Предположим, что коэффициент(это ведущий элемент 2-го шага), вычислим величины

.

Вычтем последовательно из третьего, четвертого, ..., -го уравнений системы (3.4) второе уравнение, умноженное соответственно на,,...,. В результате получим эквивалентную систему

где и() находятся по формулам

, .

Аналогично 1-му шагу получается уравнение вида

,

где , а. Матрица преобразованияимеет вид

.

k-й шаг. Целью этого шага является исключение неизвестного из уравнений с номерами. Если ведущий элемент-го шага отличен от нуля (), то

.

Вычтем из -го, ...,-го уравнений полученное на предыдущем шаге системы-е уравнение, умноженное соответственно на,,...,. Получим уравнение относительно матрицыи вектора.

Окончательно после -го шага исключения получим систему уравнений

(3.5)

Матрица системы уравнений (3.5) является верхней треугольной и обозначается. Вычисления прямого хода заканчиваются. В результате имеем

или ,

где матрица является нижней треугольной и имеет вид

.

Таким образом, получено -разложение матрицы. Исходное уравнение (3.1) можно представить как матричное уравнение

или в виде системы

Данная система решается просто, так как и– треугольные матрицы.

Обратный ход. Из последнего уравнения системы (3.5) находим . Подставляя найденное значениев предпоследнее уравнение (3.5), получаем. Осуществляя обратную подстановку, далее последовательно находим,,...,. Вычисления неизвестных здесь производятся по формулам

, (3.6)

(3.7)

при .

Число арифметических операций, необходимых для реализации данного метода, составляет приблизительно . При этом необходимо все время отслеживать неравенство нулю коэффициентов. Более того, даже если все ведущие элементы отличны от нуля, но среди них есть близкие к нулю, то в этом случае возможен неконтролируемый рост погрешности. Избежать этого недостатка метода можно используя масштабирование, т.е. подбор множителей, которые приведут систему уравнений к эквивалентной, но с коэффициентами близкими к единице. Другой способ избежать указанного недостатка состоит в перенумерации на каждом шаге неизвестных так, чтобы ведущий элемент был максимальным по модулю среди оставшихся непреобразованных уравнений.

Отметим, что при реализации метода Гаусса в программной среде MatLab для получения -разложения матрицыдостаточно одной команды –[L,U]=lu(A). При этом матрица представляет собой нижнюю треугольную матрицу, а– верхнюю треугольную.

Приведем программу для решения СЛАУ методом Гаусса. Входные данные: A – матрица размерности (nn) и b – вектор свободных членов.

Программа 3.1.

n=input(‘ n=? ‘)

L=eye(n);

C=A;

C(:,n+1)=b;

for j=1:(n-1),

for i=(j+1):n,

k=C(i,j)/C(j,j); L(i,j)=k;

C(i,:)=C(i,:)-k*C(j,:),

pause,

end,

end

U=C;

U(:,(n+1))=[];

disp(‘ L U L*U ‘)

L, U, L*U,

pause

x=zeros(n,1);

x(n)=C(n,(n+1))/C(n,n),

pause

for i=1:(n-1),

v=C((n-1),:); v(n+1)=[]; s=v*x;

x(n-i)=(C((n-i),(n+1))-s)/C((n-i),(n-i)),

pause,

end

disp(‘ невязка A*x-b ‘)

e=A*x-b

end.

Достаточно часто приходится иметь дело со специальными матрицами, а именно, с положительно-определенными симметричными ленточными матрицами. В этом случае методы вычисления точного решения СЛАУ упрощаются. Если матрица коэффициентов является симметричной и положительно-определенной, то ее разложение в произведение треугольных матриц приводит к равенству

.

При этом требуется определить только одну матрицу, что упрощает вычисления. Эта модификация метода Гаусса называется методом Холецкого или методом квадратного корня. Коэффициенты матрицы определяются из произведениянаи приравнивания соответствующих коэффициентов коэффициентам матрицы.

QR-разложение СЛАУ. Другим методом точного решения СЛАУ, пригодного для численной реализации, является представление матрицы в виде разложения в произведение ортогональной матрицыи верхней треугольной, так называемоеQR-разложение. Для ортогональной матрицы выполняется равенство , а, следовательно, решение уравнения (3.1) сводится к решению уравнения

, (3.8)

что аналогично обратному ходу метода Гаусса.

Для определения соответствующих матриц, так же как в методе Гаусса, применяется последовательное исключение неизвестных. Для исключения неизвестного из второго уравнения системы (3.1), используется матрица

,

где числа ивычисляются по формулам

, .

При этом для данных чисел выполняются следующие соотношения:

, .

В результате умножения слева матрицы на матрицупервый элемент второй строки становится равным нулю. Далее, умножением матрицына матрицу, убираетсяиз третьего уравнения системы (3.1). Для того чтобы убрать-й элемент в-й строке необходимо умножить матрицуна матрицу, в которой отличны от единичной матрицы только четыре элемента: элементы с индексами (,) и (,) равны, а элементы с индексами (,) равны (), с индексами (,) равны.

На -м шаге такого исключения получим верхнюю треугольную матрицу, которая будет равна

,

где

.

Отметим, что матрица является ортогональной, так как представляет собой произведение ортогональных матриц. Таким образом, искомая матрицаравнаи также является ортогональной. Найденное разложение используется для решения уравнения (3.1), а в программной средеMatLab это разложение выполняет оператор [Q,R]=qr(A).

Метод прогонки. Очень часто, особенно при решении дифференциальных уравнений конечно-разностными методами, получаются ленточные матрицы, в простейшем случае трехдиагональные.

Пусть СЛАУ (3.1) имеет вид

, (3.9)

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

.

Аналогично , обозначив

, ,

получим

.

На -м этапе находим равенство

,

подставив правую часть которого в соотношение для , получим

.

Откуда следует, что

(3.10)

с коэффициентами

(3.11)

где . Для равенство имеет вид

,

то есть

.

Найдя из последнего равенства , получим

. (3.12)

По формулам (3.11) найдем идля. Затем, вычислив по формуле (3.12), по формуле (3.10) последовательно определим при .

Решение данной задачи с помощью программы MatLab. В программе 3.2 рассматриваются векторы и размерностью , и размерностью , при этом, при обращении кm-файлу, содержащему приведенную ниже программу, эти вектора и матрицы должны быть заданы. Они представляют собой соответствующие диагонали матрицы системы и ее правую часть.

Программа 3.2

n=…;

m=n+1;

al=zeros(n,1);

be=zeros(n,1);

alfa=-C(1)/B(1);

beta=F(1)/B(1);

al(1)=alfa;

be(1)=beta;

for j=2:n,

d=A(j)*alfa+B(j);

alfa=-C(j)/d;

beta=(F(j)-A(j)*beta)/d;

al(j)=alfa; be(j)=beta;

end;

y(m)=(F(m)-A(n)*be(n))/(A(n)*al(n)+be(n));

for j=1:n,

y(m-j)=y(m+1-j)*al(m-j)+be(m-j);

end;

end.