2.4. Стандартные операторы мatlab для решения систем линейных алгебраических уравнений.
Наиболее простой способ решения системы линейных алгебраических уравнений можно получить с помощью оператора x=A\B . Для этого необходимо ввести матрицы A и B, а затем поделить их с помощью обратной косой черты.
За простотой использования этого оператора скрыт сложный алгоритм, который практически использует все современные методы расчетной работы с матрицами. Поэтому анализ задачи и поиск методов решения выполняет МATLAB.
Оператор solve позволяет решить систему уравнений, если ее записать в ' ' во входных параметрах оператора: [x,y,z,s]=solve('x+2*y+z+4*s=13','2*x+4*z+3*s=28','4*x+2*y+2*z+s=20','-…
3*x+y+3*z+2*s=6'). Его чаще используют для решения системы символьными методами.
При вводе матрицы B будут решаться столько систем уравнений сколько столбцов в матрице B.
Если ввести матрицу A, то можно получить матрицы L,U и P командой МATLAB: [L,U,P] = lu(A).
Полезно знать основные операции с матрицами:
+ - * ^ (сложение, вычитание, умножение и возведение в целую степень. В показателе степени может быть матрица).
' ( транспонирование матрицы).
При всех операциях необходимо согласовывать размеры матриц.
При умножении матриц не следует нарушать правила умножения матриц.
Примеры, описанных операций C = A + B, C = A - B, C = A*B, C = A^2, C=A’.
При умножении матриц можно перемножать строки, столбцы и матрицы в таких допустимых комбинациях:
Матрица А на матрицу В – получим матрицу С.
Матрица А на столбец В - получим столбец С.
Строка А на матрицу В – получим строку С.
Строка А на столбец В - получим скаляр С.
Столбец А на строку В – получим матрицу С.
Наиболее трудоемкую операцию обращения матрицы можно выполнить так:
С = A^-1, поделить единичную матрицу на А( E/A), использовать оператор inv(A).
Единичную
матрицу
можно
создать так: E=eye(3,3).
Создание нулевой матрицы zeros(3,3)
и матрицы состоящей из единиц ones(3,3).
К матрице можно обращаться с одним входным аргументом: A(5). Счет элементов в матрице идет последовательно по столбцам.
Удалить элементы матрицы, строку или столбец можно с помощью присваивания пустой матрицы [ ].
Знак двоеточие (:) означает необходимость использовать элементы строки или столбца от заданного до конца строки или столбца. Например, для матрицы A:
A=[4 -1 1;4 -8 1;-2 1 5];
A(1,:)=[ ] -удалили первую строку и C=A(1:3,2:3) –выделили блок, состоящий из второй и третей строк.
Приведем алгоритм умножения матриц в виде процедуры записанной на языке Турбо Паскаль. Он позволяет детально проследить все возможные комбинации при умножении матриц.
Procedure UmnMat (Var a,b,c:mat; Var M,N,L:integer);
Var k,I,j:integer;
Var s:rial;
Begin
For k:=1 to M do
For j:=1 to L do
Begin s:=0;
For i:=1 to N do
S:=s+A[K,I]*B[I,J];
C[K,J]:=S
End;
End.
Обращение к процедуре UmnMat(A,B,C,M,N,L);
Пусть
размерность матрицы А
,
а матрицы В
.
Произведение: матрица С размером
.
Реализован общий случай 1.
Если L=1,
то реализован случай 2 и результат имеет
вид столбца
.
Если
,
то имеем 3 случай и результат будет
строка
.
Если
и
,
то имеем 4 случай и результат будет
скаляр в
.
Если
,
то имеем 5 случай А (
)
и В(
),
а С(
),
т.е. матрица.
Порядок написания m, n, L, как фактических параметров, при обращении к процедуре
UmnMat должен соблюдаться, поскольку параметры циклов жёстко связанны с параметрами m,n,L.
Задание для самостоятельной работы:
Решить систему линейных алгебраических уравнений методами Гаусса, Гаусса-Зейделя и LU- разложения.
а)
Программа метода Гаусса-Зейделя.
function [X]=gausseid(A,B,P,delta,max1)
%A матрица [N,N] системы уравнений
%B матрица [N,1]
%P начальное приближение решения [N,1]
%err текущая ошибка
%delta допустимая ошибка решения
%iter полученное число итераций
%max1 заданное число итераций
% A=[4 -1 1;4 -8 1;-2 1 5];
% B=[7 -21 15];
%Обращение [X,iter]=gausseid(A,B,[1;2;2],0.00001,20)
N=length(B);
for k=1:max1
for j=1:N
if j==1
X(1)=(B(1)-A(1,2:N)*P(2:N))/A(1,1);
elseif j==N
X(N)=(B(N)-A(N,1:N-1)*(X(1:N-1))')/A(N,N);
else
X(j)=(B(j)-A(j,1:j-1)*X(1:j-1)-A(j,j+1:N)*P(j+1:N))/A(j,j);
end
end
err=abs(norm(X'-P));
relerr=err/(norm(X)+eps);
P=X';
if (err<delta)|(relerr<delta)
break
end
end
X=X';
k