Студопедия КАТЕГОРИИ: АвтоАвтоматизацияАрхитектураАстрономияАудитБиологияБухгалтерияВоенное делоГенетикаГеографияГеологияГосударствоДомЖурналистика и СМИИзобретательствоИностранные языкиИнформатикаИскусствоИсторияКомпьютерыКулинарияКультураЛексикологияЛитератураЛогикаМаркетингМатематикаМашиностроениеМедицинаМенеджментМеталлы и СваркаМеханикаМузыкаНаселениеОбразованиеОхрана безопасности жизниОхрана ТрудаПедагогикаПолитикаПравоПриборостроениеПрограммированиеПроизводствоПромышленностьПсихологияРадиоРегилияСвязьСоциологияСпортСтандартизацияСтроительствоТехнологииТорговляТуризмФизикаФизиологияФилософияФинансыХимияХозяйствоЦеннообразованиеЧерчениеЭкологияЭконометрикаЭкономикаЭлектроникаЮриспунденкция |
Метод LU разложения матрицы
LU – разложение матрицы – это представление квадратной матрицы А в виде произведения лево-треугольной матрицы L и верхнетреугольной матрицы U. Такое разложение обосновывается следующей теоремой: Теорема: Любая квадратная матрицы А, главные миноры которой отличны от нуля, представима в виде произведения двух матриц A=LU, где , Причем разложение единственно с точностью до диагональных элементов.
На практике для удобства элементы u11, u22, …,unn принимают равными единице. После того, как разложение A=LU получено, для нахождения решения СЛАУ можно использовать два шага: 1. Найти набор Y из решения системы L*Y=B (прямой ход МГ) 2. Найти значения переменных Х их решения U*X=Y (обратный ход МГ) Отметим, что это разложение A=LU помогает тогда, когда требуется найти решение СЛАУ при разных правых частях уравнений. В случае составления программы алгоритм следующий: 1. L – нулевая матрица, U – единичная. 2. Для i=1..n 3. Для j=1..n 4. Если i>=j, то иначе 5. Для i=1..n
6. Для i=n..1
Программа, написанная на языке Pascal, имеет вид: Uses crt; Const n1=10; Var k,i,j,n:integer; l,a,u:array[1..n1,1..n1] of real; b,y,x:array[1..n1] of real; sum:real; {определяем типы используемых массивов} BEGIN Write(‘n=’);read(n);
{БЛОК ВВОДА МАТРИЦЫ А} For i:=1 to n do begin For j:=1 to n do begin write('введите a[',i,',',j,'] '); readln(a[I,j]); end; write('введите b[',i,'] ');read(b[i]); end;
{БЛОК ПОКАЗА МАТРИЦЫ А НА ЭКРАНЕ} Writeln(‘Исходная матрица А’); For i:=1 to n do begin For j:=1 to n do write(a[i,j],’ ‘); Writeln; end;
{БЛОК НАХОЖДЕНИЯ МАТРИЦ L и U} For i:=1 to n do u[i,i]:=1; {матрица U - единичная} For i:=1 to n do For j:=1 to n do If i>=j then Begin For k:=1 to j-1 do sum:=sum+l[i,k]*u[k,j]; l[I,j]:=a[I,j]-sum; sum:=0; end else begin For k:=1 to i-1 do sum:=sum+l[i,k]*u[k,j]; u[i,j]:=(a[i,j]-sum)/l[i,i]; sum:=0; end;
{БЛОК ПОКАЗА МАТРИЦЫ L НА ЭКРАНЕ} Writeln(‘Преобразованная матрица L’); For i:=1 to n do begin For j:=1 to n do write(l[i,j],’ ‘); Writeln; End;
{БЛОК ПОКАЗА МАТРИЦЫ U НА ЭКРАНЕ} Writeln(‘Преобразованная матрица U’) For i:=1 to n do begin For j:=1 to n do write(u[i,j],’ ‘); Writeln; End; {БЛОК ПОДСЧЕТА И ПОКАЗА ПЕРЕМЕННЫХ Y} For i:=1 to n do Begin For j:=1 to i-1 do sum:= sum+l[i,j]*y[j]; Y[i]:=(b[i]-sum)/l[I,i]; writeln(‘y[’,I,’]=’,y[i]); Sum:=0; End;
{БЛОК ПОДСЧЕТА И ПОКАЗА ПЕРЕМЕННЫХ X} For i:=n downto 1 do begin For j:=i+1 to n do sum:=sum+u[i,j]*x[j]; X[i]:=y[i]-sum; writeln(‘x[’,i,’]=’,x[i]); sum:=0; End; Readln; END. Пример 2.2: Дана матрица А. Найти элементы матриц L и U. 3 2 1 A= 1 1 -1 4 -1 5 Решение: 0 0 0 1 0 0 1.Пусть L= 0 0 0, U= 0 1 0 0 0 0 0 0 1 2.i=1 3. j=1 4. 1>=1 then l11=a11=3 3.j=2 4. 1>=2 else u12=(a12-0)/l11 = 2/3
3.j=3 4. 1>=3 else u13=(a13-0)/l11 = 1/3 2.i=2 3. j=1 4. 2>=1 then l21=a21=1 3.j=2 4. 2>=2 then l22=a22-l21*u12 = 1/3 3.j=3 4. 2>=3 else u23=(a23-l21*u13)/l22 = -4 2.i=3 3. j=1 4. 3>=1 then l31=a31=4 3.j=2 4. 3>=2 then l32=a32-l31*u12 = -11/3 3.j=3 4. 3>=3 then l33=(a33-l31*u13-l32*u23) = -11
Таким образом, получены матрицы 3 0 0 1 2/3 1/3 L= 1 1/3 0, U= 0 1 -4 4 -11/3 -11 0 0 1
5.i=1 y1=b1/l11=5/3 5.i=2 y2=(b2-l21*y1)/l22=(0-1*5/3)/ 1/3 = -5 5.i=3 y3=(b3-l31*y1-l32*y2)/l33 = (3-4*5/3-(-11/3)*(-5))/(-11) = 2 Таким образом, получен столбец Y=(5/3, -5, 2)
6.i=3 x3=y3=2 6.i=2 x2=(y2-u23*x3)=(-5-(-4)*2)=3 6.i=1 x1=(y1-u12*x2-u13*x3)=5/3-2/3*3-1/3*2=-1 Таким образом, получены неизвестные х1=-1, х2 = 3, х3 = 2. 3. Нахождение обратной матрицы методом Гаусса Ненулевая матрица А называется обратимой, если существует матрица А-1, называемая обратной, и выполнено условие: А-1∙А = А∙А-1 = Е Пусть Х – обратная матрица, тогда А∙Х = Е. Представляя матрицы Х и Е в виде хi=(x1i x2i … xni)T и еi=(0 … 1 … 0)Т, получим n систем A∙xi = ei, где 1<=i<=n, где в качестве неизвестных векторов выступают столбцы искомой матрицы Х, а в качестве известных векторов правой части системы – поочередно столбцы единичной матрицы. Разложение матрицы А достаточно сделать один раз. Пример 2.3: Дана матрица А. Найти обратную к ней. 3 2 1 A= 1 1 -1 4 -1 5 Решение: Приведем аналитическое решение. Итак, необходимо решить три системы (количество решаемых систем равно порядку исходной матрицы, n=3) любым удобным способом.
1) 3x11+2x21+x31=1 x11+2/3x21+1/3x31=1/3 x11+ x21 – x31=0 ó x11+ x21 - x31=0 ó 4x11 – x21+5x31=0 4x11 - x31 + 5x31=0
x11+2/3x21+1/3x31=1/3 x11+2/3x21+1/3x31=1/3 1/3x21 - 4/3x31=-1/3 ó x21 - 4x31=-1 ó -11/3x21+11/3x31=-4/3 -11/3x21+11/3x31=-4/3
x11+2/3x21+1/3x31=1/3 x11+2/3x21+1/3x31=1/3 x11=-4/11 x21 - 4x31=-1 ó x21 - 4x31=-1 x21=9/11 -11x31= -5 x31=5/11 x31=5/11
2) 3x12+2x22+x32=0 x12+2/3x22+1/3x32=0 x12+ x22 – x32=1 ó x12+ x22 - x32=1 ó 4x12 – x22+5x32=0 4x12 - x32 + 5x32=0
x12+2/3x22+1/3x32=0 x12+2/3x22+1/3x32=0 1/3x22 - 4/3x32=1 ó x22 - 4x32=3 ó -11/3x22+11/3x32=0 -11/3x22+11/3x32=0
x12+2/3x22+1/3x32=0 x12+2/3x22+1/3x32=0 x12=1 x22 - 4x32=1 ó x22 - 4x32=3 x22=-1 -11x32= 11 x32=-1 x32=-1
3) 3x13+2x23+x33=0 x13+2/3x23+1/3x33=0 x13+ x23 – x33=0 ó x13+ x23 - x33=0 ó 4x13 – x23+5x33=1 4x13 - x33 + 5x33=1
x13+2/3x23+1/3x33=0 x13+2/3x23+1/3x33=0 1/3x23 - 4/3x33=0 ó x23 - 4x33=0 ó -11/3x23+11/3x33=1 -11/3x23+11/3x33=1
x13+2/3x23+1/3x33=0 x13+2/3x23+1/3x33=0 x13=3/11 x23 - 4x33=0 ó x23 - 4x33=0 x23=-4/11 -11x33= 1 x33=-1/11 x33=-1/11
Таким образом, обратная матрица к матрице A имеет вид: -4/11 1 3/11 X= 9/11 -1 -4/11 5/11 -1 -1/11
Проверка: А∙Х = Е
3 2 1 -4/11 1 3/11 1 0 0 1 1 -1 × 9/11 -1 -4/11 = 0 1 0 4 -1 5 5/11 -1 -1/11 0 0 1 В случае составления программы алгоритм следующий:
9. Программа, написанная на языке Pascal, имеет вид: Uses crt; Const n1=10; Var z,k,i,j,n:integer; l,a,x,b,a1,b1:array[1..n1,1..n1] of real; sum:real; BEGIN clrscr; write('n=');read(n); For i:=1 to n do For j:=1 to n do begin Write('vvedite a[',i,',',j,']='); readln(a[i,j]); end; Writeln(' matrix À'); For i:=1 to n do Begin For j:=1 to n do write(a[i,j],' '); Writeln; End; For z:=1 to n do b[z,z]:=1; For k:=1 to n-1 do For i:=k+1 to n do Begin l[i,k]:=a[i,k]/a[k,k]; for z:=1 to n do b[i,z]:=b[i,z]-l[i,k]*b[k,z]; For j:=k+1 to n do a[i,j]:=a[i,j]-l[i,k]*a[k,j]; End; For i:=1 to n do For j:=1 to n do Begin A1[i,j]:=a[i,j]/a[i,i]; for z:=1 to n do B1[i,z]:=b[i,z]/a[I,i]; end; for z:=1 to n do For k:=n downto 1 do begin For j:=k+1 to n do sum:=sum+a1[k,j]*x[j,z]; x[k,z]:=(b1[k,z]-sum)/a1[k,k]; sum:=0; end; writeln; Writeln('matrix À^(-1)'); For i:=1 to n do begin For j:=1 to n do write(x[i,j],' '); Writeln; end; Readln; END.
|
||
Последнее изменение этой страницы: 2018-04-12; просмотров: 912. stydopedya.ru не претендует на авторское право материалов, которые вылажены, но предоставляет бесплатный доступ к ним. В случае нарушения авторского права или персональных данных напишите сюда... |