Н.И. Яцкин - Линейная алгебра (Теоремы и алгоритмы) (1109879), страница 95
Текст из файла (страница 95)
1Прил. 1Коды Maple-процедурU:=ColumnOperation(E,[k,i],inplace=false);# Перестановка ненулевого элемента# в нужную позицию (второй этап).print(`U=`,U);T:=T.U;# Накопление элементарных преобразований# в матрице перехода.B:=Transpose(U).B.U;# Пересчет матрицы квадратичной формы.print(`B=`,B);fi;fi;fi;U:=Matrix(n,E);U[k..k+1,k..k+1]:=V;# Матрица преобразования,# реализующего второй прием Лагранжа.print(`U=`,U);T:=T.U;# Накопление элементарных преобразований# в матрице перехода.B:=Transpose(U).B.U;# Пересчет матрицы квадратичной формы.print(`B=`,B);elif B[k,k]=0 then# Если диагональ юго-восточного блока ненулевая,# но начальный ее элемент - нулевой,# то проводим подготовку к применению# первого приема Лагранжа.for i from k+1 to n doif B[i,i]<>0 then# Обнаружен ненулевой диагональный элемент.break;fi;od;U:=ColumnOperation(E,[k,i],inplace=false);# Перестановка ненулевого диагонального# элемента в нужную позицию.print(`U=`,U);T:=T.U;# Накопление элементарных преобразований# в матрице перехода.B:=Transpose(U).B.U;# Пересчет матрицы квадратичной формы.print(`B=`,B);fi;579580Коды Maple-процедурU:=Matrix(n,E);# Заготовка для матрицы преобразования,# реализующего первый прием Лагранжа;# далее она заполняется.for j from k+1 to n doU[k,j]:=-B[k,j]/B[k,k];od;if not Equal(U,E) then# Если преобразование не тождественное,# то применяем его к матрице квадратичной формы,# производим "накопление"# в результирующей матрице перехода# и выдаем промежуточные результаты на печать.print(`U=`,U);T:=T.U;# Накопление элементарных преобразований# в матрице перехода.B:=Transpose(U).B.U;# Пересчет матрицы квадратичной формы.print(`B=`,B);fi;fi;od;RETURN(B,T);# Возвращаются: диагональная матрица,# задающая диагональный вид данной квадратичной формы,# и матрица перехода к диагональному виду.end proc;> Quadro[Jacob]:=proc(A::'Matrix'(symmetric))local B,n,i,j,k,Delta,mu,T,AT,sys,sol;#####Процедура Jacob приведенияк конгруэнтному диагональному видусимметрической квадратной матрицы(отвечающей квадратичной форме),в предположении выполнения условия Якоби.n:=RowDimension(A);# Размер матрицы.Delta[0]:=1;# Угловой минор нулевого порядка.Прил.
1Прил. 1Коды Maple-процедурfor i from 1 to n doDelta[i]:=Determinant(SubMatrix(A,1..i,1..i));# Вычисление угловых миноров.print(evaln(Delta[i])=Delta[i]);if i<n and Delta[i]=0 then# Если очередной (не последний) угловой минор# оказывается нулевым, то# выдается сообщение об ошибке.ERROR("Условия Якоби не выполнены!");else# В противном случае вычисляются# диагональные элементы диагональной матрицы,# задающей диагональный вид# данной квадратичной формы.mu[i]:=Delta[i]/Delta[i-1];fi;od;for i from 1 to n doprint(evaln(mu[i])=mu[i]);od;B:=Matrix(n):# Заготовка для диагональной матрицы.# Далее следует ее заполнение.for i from 1 to n doB[i,i]:=mu[i]:od:T:=Matrix(n,symbol=`t`):# Заготовка для унитреугольной# матрицы перехода к диагональному виду.# Далее следует ее заполнение# (элементы выше главной диагонали# остаются неопределенными).for i from 1 to n doT[i,i]:=1;od:for i from 2 to n dofor j from 1 to i-1 doT[i,j]:=0;od;od;581582Коды Maple-процедурprint(`T=`,T);AT:=A.T;# Вычисление произведения# матрицы квадратичной формы# на матрицу перехода.print(map(simplify,AT));sys:={}:# Начало формирования# системы (множества) линейных уравнений# для отыскания неопределенных элементов# матрицы перехода.for i from 1 to n dofor j from i+1 to n dosys:=sys union {AT[i,j]=0};# К (изначально пустому)# множеству уравнений добавляется# очередное, отвечающее# наддиагональному элементу# ранее вычисленного произведения матриц.od;od;print(sys);sol:=solve(sys);# Вычисляется решение указанной системы.print(sol);# Далее в матрицу перехода# заносятся вычисленные значения# наддиагональных элементов.for i from 1 to n dofor j from i+1 to n dofor k from 1 to nops(sol) doif T[i,j]=lhs(sol[k]) thenT[i,j]:=rhs(sol[k]);break;fi;od;od;od;RETURN(B,T);# Возвращаются: диагональная матрица,# задающая диагональный вид данной# квадратичной формы,# и матрица перехода к диагональному виду.end proc;Прил.
1Прил. 1Коды Maple-процедур> Quadro[Signature]:=proc(A::'Matrix'(diagonal))local n,B1,B2,E,T1,T2,T,i,j,s,t,r,pos,neg,zer,str,new_ord;#####Процедура Signature (сценарного типа)приведения диагональной матрицы(отвечающей квадратичной форме)к конгруэнтному нормальному виду(над полем действительных чисел).n:=RowDimension(A);# Размер матрицы.E:=IdentityMatrix(n);# Единичная матрица.s:=0;t:=0;# Заготовки для накопления# индексов инерции.pos:=[];neg:=[];zer:=[];# Заготовки списков номеров переменных,# которым отвечают# положительные (отрицательные, нулевые)# диагональные коэффициенты.# Далее следует просмотр диагонали# и заполнение указанных списков.for i from 1 to n doif A[i,i]>0 thens:=s+1;pos:=[pos[],i];elif A[i,i]<0 thent:=t+1;neg:=[neg[],i];elsezer:=[zer[],i];fi;od:r:=s+t;# Вычисление ранга матрицы.print(`s=`||s,`t=`||t,`r=`||r);print(`pos=`,pos,`neg=`,neg,`zer=`,zer);# Далее определяется# тип квадратичной формы# (в плане знакоопределенности).583584Коды Maple-процедурif s=n thenstr:=`Форма положительно определена.`;elif s<n and s>0 and t=0 thenstr:=`Форма положительно полуопределена.`;elif s<n and s>0 and t>0 thenif r=n thenstr:=`Форма невырожденна и знакопеременна.`;elsestr:=`Форма вырождена и знакопеременна.`;fi;elif s=0 and t=0 thenstr:=`Форма нулевая.`;elif s=o and t<n thenstr:=`Форма отрицательно полуопределена.`;elsestr:=`Форма отрицательно определена.`;fi;print(str);new_ord:=[pos[],neg[],zer[]];# Переупорядочивание переменных:# сначала должны идти переменные,# которым отвечают положительные# диагональные элементы,# потом - те, которым отвечают отрицательные,# и, наконец, - те, которым отвечают нулевые# диагональные элементы.T1:=Matrix(n,E);# Заготовка для матрицы# перестановочного перехода.####Далее следует ее заполнение:столбцы единичной матрицырасполагаются в соответствиис новым порядком.for j from 1 to n doT1[1..n,j..j]:=E[1..n,new_ord[j]..new_ord[j]];od;B1:=Transpose(T1).A.T1;# Вычисление матрицы,# отвечающей диагональному виду,# с правильно упорядочненными по знаку# диагональными элементами.print(`T1=`,T1,`B1=`,B1);T2:=Matrix(n,E):# Заготовка для матрицы# перехода к нормальному виду.Прил.
1Прил. 1Коды Maple-процедур# Далее следует ее заполнение.for i from 1 to n doif i<=s thenT2[i,i]:=1/sqrt(B1[i,i]);elif i>s and i<=r thenT2[i,i]:=1/sqrt(-B1[i,i]);elseT2[i,i]:=1;fi;od:B2:=Transpose(T2).B1.T2;# Вычисление матрицы,# отвечающей нормальному виду# квадратичной формы.print(`T2=`,T2,`B2=`,B2);RETURN(B2,T1.T2,[s,t]);# Возвращаются: матрица нормального вида,# матрица перехода# от диагонального вида к нормальному,# сигнатура (индексы инерции)# данной квадратичной формы.end proc;> save Quadro,"F:/MaplePackages/Quadro.m";# Сохранение пакета.585586Коды Maple-процедурПрил. 13а. Пример применения процедур пакета Quadro(к ТР3 "Диагонализациясимметрических билинейных (квадратичных)форм "; п. 39.1)> restart;with(LinearAlgebra):> read "F:/MaplePackages/Quadro.m"; with(Quadro);[ Jacob , Lagr , Signature ]>A:=Matrix([[0,1,0,1,0,1,0],[1,0,1,1,1,1,1],[0,1,0,0,1,0,1],[1,1,0,0,1,1,1],[0,1,1,1,0,1,0],[1,1,0,1,1,0,1],[0,1,1,1,0,1,0]]);⎡0⎢⎢⎢⎢1⎢⎢0A := ⎢⎢1⎢⎢⎢⎢0⎢⎢1⎢⎢⎣0101111101001011100111011101011011010⎤⎥1 ⎥⎥⎥1 ⎥⎥1 ⎥⎥⎥0 ⎥⎥⎥1 ⎥⎥⎥0 ⎥⎦> (DL,TL):=Lagr(A);# Запуск программы диагонализации по Лагранжу.# Переменным DL и TL будут присвоены значения# диагональной матрицы, конгруэнтной введенной,# и матрицы перехода к диагонализирующему базису.# Следите за заменами переменных# (они регистрируются в матрицах U)# и за преобразованиями данной матрицы# (регистрируются в матрицах B).# Так, на первом шаге применяется# второй прием Лагранжа,# в результате чего на диагонали,# которая была в исходной матрице нулевой,# появляются ненулевые элементы.k=1⎡⎢1⎢⎢1⎢⎢⎢0U= , ⎢⎢0⎢⎢⎢⎢0⎢0⎢⎢⎢⎣0-110000000100000001000000010000000100⎤⎥0 ⎥⎥⎥0 ⎥⎥0 ⎥⎥⎥0 ⎥⎥⎥0 ⎥⎥⎥1 ⎥⎦Прил.
1Коды Maple-процедур⎡2⎢⎢⎢⎢0⎢⎢1B= , ⎢⎢2⎢⎢⎢⎢1⎢⎢2⎢⎢⎣1⎡⎢1⎢⎢⎢⎢0⎢⎢⎢0U= , ⎢⎢⎢⎢0⎢⎢0⎢⎢0⎢⎢⎣0⎡⎢ 2⎢0⎢⎢⎢⎢⎢⎢ 0⎢⎢0B= , ⎢⎢⎢⎢⎢0⎢⎢⎢⎢ 0⎢⎢⎢⎢0⎣0-2101011100101-1201000001000000-2-100100001-12-112-11210101200011100-1-20-101111010-1200010001120-120-1220011015871⎤⎥1 ⎥⎥⎥1 ⎥⎥1 ⎥⎥⎥0 ⎥⎥⎥1 ⎥⎥⎥0 ⎥⎦-100001000-1-10-20-1 ⎤⎥2 ⎥⎥0 ⎥⎥⎥0 ⎥⎥⎥0 ⎥⎥⎥0 ⎥⎥0 ⎥⎥⎥1 ⎥⎦0⎤⎥1 ⎥⎥⎥1 ⎥⎥2 ⎥⎥⎥0 ⎥⎥⎥-1 ⎥⎥2 ⎥⎥⎥0 ⎥⎥⎥-1 ⎥⎥2 ⎥⎦k=2⎡⎢ 1⎢⎢⎢⎢ 0⎢⎢⎢0U= , ⎢⎢⎢⎢ 0⎢⎢ 0⎢⎢ 0⎢⎢⎣0⎡2⎢⎢⎢⎢0⎢⎢0B= , ⎢⎢0⎢⎢⎢⎢0⎢⎢0⎢⎢⎣0000000012100000-200000000-11-1110010000120010000-1-20-1000100000000001000-1-10-200⎤⎥1 ⎥⎥2 ⎥⎥⎥0 ⎥⎥⎥0 ⎥⎥⎥0 ⎥⎥0 ⎥⎥⎥1 ⎥⎦0⎤⎥0 ⎥⎥⎥1 ⎥⎥0 ⎥⎥⎥0 ⎥⎥⎥0 ⎥⎥⎥0 ⎥⎦588Коды Maple-процедурПрил. 1k=3⎡1⎢⎢⎢⎢0⎢⎢0U= , ⎢⎢0⎢⎢⎢⎢0⎢⎢0⎢⎢⎣0⎡⎢2⎢⎢0⎢⎢⎢0B= , ⎢⎢0⎢⎢⎢⎢0⎢0⎢⎢⎢⎣0⎡⎢1⎢⎢0⎢⎢⎢⎢0⎢U= , ⎢⎢⎢⎢0⎢⎢0⎢⎢0⎢⎢⎣0⎡2⎢⎢⎢⎢ 0⎢⎢ 0⎢⎢⎢0B= , ⎢⎢⎢⎢⎢⎢ 0⎢⎢⎢⎢ 0⎢⎢⎣0010000000010000010000000010000000100-20000000-2-10-1000-101-11000100000-1-10-200⎤⎥0 ⎥⎥⎥0 ⎥⎥1 ⎥⎥⎥0 ⎥⎥⎥0 ⎥⎥⎥0 ⎥⎦010000010000000000-12100000-1200100⎤⎥0 ⎥⎥⎥⎥0 ⎥⎥⎥⎥0 ⎥⎥⎥0 ⎥⎥0 ⎥⎥⎥1 ⎥⎦0-2000-200000000000121-121001000001000000-120-3200⎤⎥0 ⎥⎥⎥0 ⎥⎥0 ⎥⎥⎥0 ⎥⎥⎥0 ⎥⎥⎥1 ⎥⎦0⎤⎥0 ⎥⎥⎥0 ⎥⎥⎥⎥1 ⎥⎥⎥⎥0 ⎥⎥⎥⎥0 ⎥⎥⎥⎥0 ⎥⎦k=4⎡1⎢⎢⎢⎢0⎢⎢0U= , ⎢⎢0⎢⎢⎢⎢0⎢⎢0⎢⎢⎣0010000000100000001000000-210000010100⎤⎥0 ⎥⎥⎥0 ⎥⎥-2 ⎥⎥⎥0 ⎥⎥⎥0 ⎥⎥⎥1 ⎥⎦Прил.
1Коды Maple-процедур⎡⎢2⎢⎢0⎢⎢⎢⎢0⎢B= , ⎢⎢0⎢⎢⎢⎢0⎢⎢0⎢⎢⎣00-2000-2000000000001200000000000-21-21-21000012105890⎤⎥0 ⎥⎥⎥0 ⎥⎥⎥⎥0 ⎥⎥⎥⎥-2 ⎥⎥1 ⎥⎥⎥-2 ⎥⎦k=5⎡⎢ 1⎢⎢ 0⎢⎢⎢⎢ 0⎢0U= , ⎢⎢⎢⎢ 0⎢⎢⎢⎢ 0⎢⎢⎣0⎡2⎢⎢⎢⎢ 0⎢⎢ 0⎢⎢⎢0B= , ⎢⎢⎢⎢⎢⎢ 0⎢⎢⎢⎢ 0⎢⎢⎣001000010000100000001000000000-2000-20000000120000000000⎤⎥0 ⎥⎥⎥0 ⎥⎥0 ⎥⎥⎥⎥-1 ⎥⎥⎥⎥0 ⎥⎥⎥1 ⎥⎦00000000-20-3200⎤⎥0 ⎥⎥⎥0 ⎥⎥⎥⎥0 ⎥⎥⎥⎥0 ⎥⎥⎥⎥0 ⎥⎥⎥⎥0 ⎥⎦k=6⎡⎢ 2⎢0⎢⎢⎢⎢ 0⎢⎢⎢0DL , TL := ⎢⎢⎢⎢⎢⎢ 0⎢⎢⎢⎢ 0⎢⎢⎣00-2000-2000012000000000000-20-32000000000⎡⎢0 ⎤ ⎢⎢⎥⎥ ⎢0⎥ ⎢⎥ ⎢0 ⎥⎥ ⎢⎢⎥⎥ ⎢⎢0 ⎥⎥ ⎢⎢⎥, ⎢⎢⎥0 ⎥⎥ ⎢⎢⎥ ⎢⎢⎥0 ⎥⎥ ⎢⎢⎥⎥ ⎢⎢0 ⎥⎦ ⎢⎢⎢⎢⎣-12121-121-1-111-10000010000100000000000-1-21-1-120-1212100 ⎤⎥⎥⎥⎥⎥0 ⎥⎥⎥⎥0 ⎥⎥⎥⎥0 ⎥⎥⎥⎥⎥-1 ⎥⎥⎥⎥0 ⎥⎥⎥1 ⎥⎦590Коды Maple-процедурПрил.
1> (DJ,TJ):=Jacob(A);Δ1 = 0Error, (in Jacob) Условия Якоби не выполнены!> (DS,TS,sg):=Signature(DL);# Нормализация и вычисление сигнатуры.s=2 , t=4 , r=6pos= , [ 1 , 4 ], neg= , [ 2 , 3 , 5 , 6 ], zer=, [ 7 ]Форма вырождена и знакопеременна⎡⎢1⎢⎢0⎢⎢⎢0T1= , ⎢⎢0⎢⎢⎢⎢0⎢⎢0⎢⎢⎣0⎡⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢T2= , ⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎢⎣00010000100000001000000001000000010⎡2⎢⎢0⎤⎢⎢ 0⎥0 ⎥⎥⎢⎢⎥⎢⎢ 0⎥0⎥⎢⎥0 ⎥, B1= , ⎢⎢ 0⎥⎥⎢⎢0⎥⎢⎢ 0⎥⎥0⎥⎢⎢⎥⎥⎢⎢ 01⎦⎢⎢⎣02200000020000002200000022000000220000000000063001200000000000-2000-2000-200000000000-320⎤0 ⎥⎥⎥⎥0 ⎥⎥⎥⎥⎡⎢1⎥⎢00⎥⎢⎢⎥⎥⎢0⎥⎥⎥, B2= , ⎢⎢0⎢⎢0⎥⎥⎥⎢⎢0⎥⎥⎢⎢⎢00⎥⎥⎥⎢⎢⎣0⎥⎥0 ⎥⎥⎥⎥1 ⎥⎦010000000-10000000-10000⎤⎥⎥0 ⎥⎥⎥⎥0 ⎥⎥⎥0 ⎥⎥⎥0 ⎥⎥⎥⎥0 ⎥⎥⎥⎥0 ⎥⎦0000-10000000-100⎤⎥0 ⎥⎥⎥0 ⎥⎥0 ⎥⎥⎥0 ⎥⎥⎥0 ⎥⎥⎥0 ⎥⎦Прил.