Wyznacznik macierzy kwadratowej

Niech A będzie macierzą stopnia n (czyli macierzą kwadratową n x n), o elementach (Subscript[A, ik]). Utwórzmy iloczyn  (Subscript[A, 1(Subscript[α, 1])]) (Subscript[A, 2(Subscript[α, 2])]) ... (Subscript[A, (Subscript[nα, n])]) , gdzie ((Subscript[α, 1 ]),(Subscript[α, 2 ]), ... (Subscript[α, n ])) jest permutacją liczb {1, 2,  ... , n}.
Ten iloczyn zawiera dokładnie jeden element każdego wiersza i każdej kolumny macierzy A. Mnożymy go przez znak permutacji sgn((Subscript[α, 1 ]) (Subscript[α, 2 ])... (Subscript[α, n ])) i tworzymy sumę po wszystkich permutacjach

(∑_(
wszystkie permutacje 
(Subscript[α, 1 ]) (Subscript[α, 2 ]) ...(Subscript[α, n ])   ))  sgn((Subscript[α, 1 ]) (Subscript[α, 2 ])...(Subscript[α, n ]))(Subscript[A, 1(Subscript[α, 1])]) (Subscript[A, 2(Subscript[α, 2])])...(Subscript[A, (Subscript[nα, n])]) ,
 
którą nazywamy wyznacznikiem macierzy A i oznaczamy Det(A).

Uwaga: zwykle liczymy wyznaczniki, korzystając z twierdzeń, a nie wprost z definicji.

Kluczowe w tej definicji jest zrozumienie, czym jest znak permutacji ! 
 sgn((Subscript[α, 1 ]) (Subscript[α, 2 ])...(Subscript[α, n ]))to (-1)^(liczba przestawień w ciągu  ((Subscript[α, 1 ]) (Subscript[α, 2 ])...(Subscript[α, n ])), która daje w wyniku ciąg (1 2 ... n))

P4=Permutations[{1,2,3,4}]

Przykład permutacji nieparzystej:
sgn(4231) =(-1)^(liczba przestawień w ciągu (4,2,3,1), która daje ciąg (1,2,3,4). Widać, że wystarczy jedno przestawienie (1⟷4), bo od razu daje (1234).
Można rozważać także inne przestawienia, ale ich liczba będzie zawsze nieparzysta ! 
Na przykład:
najpierw (1⟷3) daje (4213),
potem (2⟷1) daje (4123),
potem (4⟷1) daje (1423), 
potem (4⟷2) daje (1243),
potem (4⟷3) daje (1234)

sgn4231=Signature[{4,2,3,1}]

Przykład permutacji parzystej
sgn(1342) = (-1)^(liczba przestawień w ciągu (1,3,4,2), która daje ciąg (1,2,3,4). Widać, że wystarczą dwa przestawienia:
najpierw (4⟷2) daje (1,3,2,4), 
potem (3⟷2) daje (1,2,3,4).
Można rozważać także inne przestawienia, ale ich liczba będzie zawsze parzysta ! 
Na przykład:
najpierw (2⟷1) daje (2341),
potem (4⟷1) daje (2314),
potem (3⟷1) daje (2134), 
potem (1⟷2) daje (1234)

sgn1342=Signature[{1,3,4,2}]

Macierz 1 x 1

A1={{a11}};MatrixForm[A1]

(* Dla jednego elementu mamy tylko permutację tożsamościową z zerową, a więc parzystą liczbą przestawień ! *)

Det[A1]

Macierz 2 x 2

A2={{a11,a12},{a21,a22}};MatrixForm[A2]

(* Teraz mamy dwie permutacje: parzystą (12) i nieparzystą (21)*)

Det[A2]

Macierz 3 x 3

A3={{a11,a12,a13},{a21,a22,a23},{a31,a32,a33}};MatrixForm[A3]

(* W tym przypadku spośród sześciu permutacji trzy są parzyste i trzy nieparzyste *)

Det[A3]

(obrazek z Wikipedii dla tzw. schematu Sarrusa)

Macierz 4 x 4

A4={{a11,a12,a13,a14},{a21,a22,a23,a24},{a31,a32,a33,a34},{a41,a42,a43,a44}};MatrixForm[A4]

Det[A4]

(* Zwłaszcza ten ostatni wzór niełatwo zapamiętać i są prostsze sposoby liczenia wyznaczników dla większych macierzy. *)

Ostrzeżenie : nie istnieje schemat Sarrusa dla macierzy 4 x 4 !!!

(* Spróbujmy teraz powtórzyć otrzymany powyżej wynik na Det[A4], korzystając wprost z definicji - jest to bardzo dobre ćwiczenie z operacji na listach *)

P4

P4[[1]]

P4[[1]][[1]]

P4[[1]][[2]]

Signature[P4[[1]]]

Length[P4]

MyDetA4=Sum[Signature[P4[[i]]]*A4[[1,P4[[i]][[1]]]]*A4[[2,P4[[i]][[2]]]]*A4[[3,P4[[i]][[3]]]]*A4[[4,P4[[i]][[4]]]],{i,1,Length[P4]}]

Det[A4]==MyDetA4

W przypadku konkretnych kwadratowych macierzy A i B mamy

A={{1,0,3,4},{2,4,-2,8},{1,0,-4,6},{3,-7,9,1}}

MatrixForm[A]

B={{-2,1,8,0},{7,0,-2,2},{0,1,2,3},{2,1,-9,3}}

MatrixForm[B]

Det[A]

Det[B]

(* Zachodzi ogólnie dla macierzy kwadratowych: Det[A.B]=Det[B.A]=Det[A]*Det[B]. Na przykład:  *)

Det[A.B]

Det[B.A]

Det[A]*Det[B]

Ślad (trace) macierzy kwadratowej: suma elementów na przekątnej macierzy

Tr[A2]

Tr[A3]

Tr[A4]

W przypadku naszych kwadratowych macierzy A i B

Tr[A]

Tr[B]

(* Zachodzi Tr[A.B]=Tr[B.A]. Na przykład: *)

Tr[A.B]

Tr[B.A]

H=A+I*B;

Tr[H]

(* Przypomnienie: sprzężenie hermitowskie *)

G=ConjugateTranspose[H]

Tr[G]

Szczególne przypadki macierzy

(* macierz o rozmiarach m x n i wszystkich elementach identycznych  *)

Cst=ConstantArray[3,{3,5}];MatrixForm[Cst]

(* Przy pomocy ConstantArray można zbudować zerową macierz kwadratową o rozmiarze n *)

NullMatrix[n_]:=ConstantArray[0,{n,n}]

Nm4=NullMatrix[4]

MatrixForm[Nm4]

(* macierz jednostkowa (identity matrix) *)

Id4=IdentityMatrix[4]

MatrixForm[IdentityMatrix[4]]

(* Zachodzi *)

A4.Id4==A4

Id4.A4==A4

(* Macierz diagonalna zbudowana ze wskazanych elementów *)

Diag4=DiagonalMatrix[{d11,d22,d33,d44}]

MatrixForm[Diag4]

Macierz odwrotna (inverse matrix) istnieje dla macierzy kwadratowej A  o niezerowym wyznaczniku. Oznaczamy ją symbolem (A^(-1)). Z definicji zachodzi (A^(-1)).A=A.(A^(-1))=Id (macierz jednostkowa).

Inverse[A2]

Inverse[A3]

InvA=Inverse[A]

InvB=Inverse[B]

MatrixForm[InvA.A]

MatrixForm[A.InvA]

MatrixForm[InvB.B]

MatrixForm[B.InvB]

(* Ilustracja twierdzenia o wyznaczniku macierzy odwrotnej: Jeśli istnieje macierz odwrotna A^(-1), to zachodzi 
Det[A^(-1)]=1/Det[A] 
 *)

Det[Inverse[A]]==1/Det[A]

Det[Inverse[B]]==1/Det[B]

(* Ilustracja twierdzenia: Jeśli istnieją macierze odwrotne A^(-1) i B^(-1) (oczywiście o tych samych rozmiarach), to zachodzi 
(A.B)^(-1)= B^(-1).A^(-1)
 *)

Simplify[Inverse[A.B]==Inverse[B].Inverse[A]]

MatrixForm[Inverse[A.B]]

MatrixForm[Inverse[B].Inverse[A]]

(* Poniżej pokazuję na przykładach macierzy A i B, że dwie kolejne operacje (1) liczenia macierzy odwrotnej i (2) transponowania macierzy można wykonywać w dowolnej kolejności *)

Inverse[Transpose[A]]==Transpose[Inverse[A]]

MatrixForm[Inverse[Transpose[A]]]

MatrixForm[Transpose[Inverse[A]]]

Inverse[Transpose[B]]==Transpose[Inverse[B]]

MatrixForm[Inverse[Transpose[B]]]

MatrixForm[Transpose[Inverse[B]]]

Układy równań liniowych Cramera. Twierdzenie Cramera

Rozważmy układ n równań liniowych na n niewiadomych x1, x2, ..., xn. Układ taki można zapisać w postaci:
A11*x1 + A12*x2 + ... + A1n*xn = B1,
A21*x1 + A22*x2 + ... + A2n*xn = B2,
...
An1*x1 + An2*x2 + ... + Ann*xn = Bn.

Współczynniki Aik tworzą macierz kwadratową - tzw. macierz współczynników. Mówimy też o kolumnie wyrazów wolnych
B=(({{B1}, {B2}, {...}, {Bn}}))  oraz o kolumnie niewiadomych X = (({{X1}, {X2}, {...}, {Xn}})) .
W postaci macierzowej układ równań można zapisać w postaci A.X = B. 
Ten zapis jest  absolutnie równoważny zapisowi stosowanemu powyżej (z jawnym układem równań, bo równość dwóch macierzy oznacza równość wszystkich odpowiadających sobie elementów macierzy.
  
Jeśli Det[A] jest różny od zera, to mamy do czynienia z układem Cramera.

Twierdzenie Cramera:
Układ równań z Det[A] różnym od zera ma jednoznaczne rozwiązanie, które możemy zapisać w postaci    X= A^(-1).B

Równoważne sformułowanie (nieoczywiste bez znajomości liczenia macierzy odwrotnych):
Niewiadomą Xs otrzymujemy jako iloraz dwóch wyznaczników:
Xs=Det[As]/Det[A], gdzie 
As jest macierzą A, w której s-ta kolumna została zastąpiona przez kolumnę wyrazów wolnych, czyli B.

(* Przykład *)

A={{2,3,1,-2,1},{3,-2,4,-6,2},{7,2,-1,-1,3},{6,0,1,-9,4},{5,-1,2,-8,5}}

MatrixForm[A]

Det[A]

(* Mamy układ Cramera ! *)

MatrixForm[A]

B={{1},{-2},{-1},{3},{0}}

MatrixForm[B]

X={{X1},{X2},{X3},{X4},{X5}}

MatrixForm[X]

ukladrownan=A.Flatten[X]==Flatten[B]

Flatten[X]

Flatten[B]

Solve[ukladrownan,Flatten[X]]

(* Flatten[B] i Flatten[X] "spłaszcza" wielowymiarowe listy do postaci prostych list jednowymiarowych *)

(* Nie ma potrzeby samemu tworzyć układ równań i rozwiązywać go przy pomocy ogólnej instrukci Solve, bo mamy narzędzie bardziej specjalistyczne; LinearSolve, gdzie podajemy jako argumenty macierz współczynników oraz kolumnę wyrazów wolnych. *)

LinearSolve[A,B]

(* Oczywiście możemy też użyć macierzy odwrotnej do A *)

Inverse[A].B

Rząd macierzy

Rozważmy macierz A o m wierszach i n kolumnach. Jeśli wykreślić m-k wierszy i n-k kolumn to pozostanie z niej macierz kwadratowa stopnia k ( k x k). Wyznacznik z tej macierzy nazywamy minorem stopnia k wyjętym z macierzy A. Jeśli istnieje różny od zera minor stopnia r wyjęty z macierzy A, a wszystkie minory stopnia r+1 są równe zero, to mówimy, że macierz A jest rzędu r i piszemy Rank[A]= r.

Twierdzenia, które zilustrujemy na ćwiczeniach.
Rząd macierzy nie ulega zmianie, jeśli:
(1) pomnożyć dowolną kolumnę przez liczbę rożną od zera
(2) do elementów jednej kolumny dodać odpowiednie elementy innej kolumny pomnożone przez jakąś stałą
(3) przestawić kolumny
(4) transponować macierz

A={{1,2,3},{4,5,6},{7,8,9},{0,0,0}}

MatrixForm[A]

AA={{1,2,3},{4,5,6},{7,8,9}}

MatrixForm[AA]

Det[AA]

AAA={{1,2},{4,5}}

MatrixForm[AAA]

Det[AAA]

(* Dlatego ... *)

MatrixRank[A]

Mamy do dyspozycji użyteczną instrukcję, 
     RowReduce[A],
która daje uproszczoną (tzw. zredukowaną) macierz o tym samym rzędzie co macierz A.

MatrixForm[RowReduce[A]]

(* Po wykonaniu tej instrukcji widać, że z macierzy A nie da się wykroić macierzy kwadratowej o niezerowym wyznaczniku dla rozmiaru macierzy większego od 2 *)

Ale uwaga !!!

B={{1,1,b},{1,1,1},{1,2,1}}

MatrixForm[B]

Det[B]

MatrixRank[B]

To jest niepoprawny wynik !!!
Wynik poprawny to:
 rz(B)=3, jeśli b ≠ 1,
 rz(B)=2, jeśli b = 1.

Twierdzenie Kroneckera-Capelli’ego.
Dotyczy ono najogólniejszego układu m równań na n niewiadomych x1, x2, ..., xn. Układ taki można zapisać w postaci:
A11*x1 + A12*x2 + ... + A1n*xn = B1,
A21*x1 + A22*x2 + ... + A2n*xn = B2,
...
Am1*x1 + Am2*x2 + ... + Amn*xn = Bm.

Współczynniki Aik tworzą macierz (niekoniecznie kwadratową) - tzw. macierz współczynników. Mówimy też o kolumnie wyrazów wolnych
B=(({{B1}, {B2}, {...}, {Bm}}))  oraz o kolumnie niewiadomych X = (({{X1}, {X2}, {...}, {Xn}})) .
W postaci macierzowej układ równań można zapisać w postaci:
A.X = B.

Dalej tworzymy macierz DAB przez dołączenie do macierzy A kolumny wyrazów wolnych, B. 

Twierdzenie mówi, że rozwiązanie układu A.X=B istnieje tylko wtedy, gdy Rank[A]=Rank[DAB]. Ponadto, gdy Rank[A]=Rank[DAB]=n (liczba niewiadomych), to istnieje dokładnie jedno rozwiązanie; gdy Rank[A]=Rank[DAB] =r < n, to istnieje nieskończenie wiele rozwiązań zależnych od n-r parametrów.

(* Przykład 1: zadanie 3.2a z książki Henryka Arodzia 
i Krzysztofa Rościszewskiego, 
Algebra i geometria analityczna w zadaniach, Wydawnictwo Znak, Kraków, 2005
({{\[Piecewise], ({{2x1+x2+x3+x4+x5, =1,}, {-x1-x2+x3-x4+x5, =3,}, {0x1-x2+3x3-x4+3x5, =7.}})}})                *)

A={{2,1,1,1,1},{-1,-1,1,-1,1},{0,-1,3,-1,3}}

B={1,3,7}

MatrixForm[A]

DAB={{2,1,1,1,1,1},{-1,-1,1,-1,1,3},{0,-1,3,-1,3,7}}

MatrixForm[DAB]

rA=MatrixRank[A]

rDAB=MatrixRank[DAB]

(* rA=rDAB, więc układ ma rozwiązanie zależne od (5-2=3) parametrów. *)

MatrixForm[RowReduce[A]]

(* Jako nieosobliwy minor możemy wybrać macierz 2 x 2 z lewego górnego rogu A i sprowadzić nasz układ do następującej postaci:

({{\[Piecewise], ({{2x1+x2, = 1-x3-x4-x5,}, {-x1-x2, = 3-x3+x4-x5.}})}})

Stąd ostatecznie dostaniemy

x1=4-2*x3-2*x5
x2=-7+3*x3-x4+3*x5
*)

x1=4-2*x3-2*x5 (* to dostałem, dodając dwa równania stronai *)

x2=Simplify[1-x3-x4-x5-2*x1] (* liczę x2 z pierwszego równania *)

(* Sprawdzamy, czy nasze rozwiązanie na x1 i x2 (zalezne od x3, x4 i x5) spełnia układ równań *)

Simplify[2x1+x2+x3+x4+x5==1]

Simplify[-x1-x2+x3-x4+x5==3]

Simplify[0*x1-x2+3*x3-x4+3*x5==7]

(* Teraz sformułujemy układ równań i zobaczymy, jak komendy programu Mathematica sobie z nim poradzą *)

Clear[x1,x2,x3,x4,x5];

ukladKC=A.{x1,x2,x3,x4,x5}==B

(* Co daje komenda Solve dla naszego układu ? *)

Solve[{2 x1+x2+x3+x4+x5,-x1-x2+x3-x4+x5,-x2+3 x3-x4+3 x5}=={1,3,7},{x1,x2,x3,x4,x5}]

(* Dostaliśmy rozwiązanie na x4 i x5 zależne od x1, x2 i x3, co wygląda rozsądnie. *)

x4=-1-((3 x1)/2)-x2;x5=2-(x1/2)-x3;

(* Sprawdzamy, czy to rozwiązanie na x4 i x5, uzyskane przy pomocy Solve i zalezne od x1, x2 i x3, spełnia układ równań *)

Simplify[2x1+x2+x3+x4+x5==1]

Simplify[-x1-x2+x3-x4+x5==3]

Simplify[0*x1-x2+3*x3-x4+3*x5==7]

Solve[{2 x1+x2+x3+x4+x5,-x1-x2+x3-x4+x5,-x2+3 x3-x4+3 x5}=={1,3,7},{x1,x2}]

(* A co daje LinearSolve ? *)

LinearSolve[A,B]

(* Na pewno nie jest to pełne rozwiązanie ! *)

Zobaczmy jeszcze na zakończenie, jak działa twierdzenie, w przypadku układu, który nie posiada rozwiązania, czyli układu sprzecznego

(* Przykład 2: "z głowy" *)

Clear[DAB,X,x1,x2,x3,x4,x5,ukladKC,rA,rDAB];

A={{1,2,3},{1,2,3},{3,1,-1}}

MatrixForm[A]

B={1,2,0}

MatrixForm[B]

DAB={{1,2,3,1},{1,2,3,2},{3,1,-1,0}}

MatrixForm[DAB]

MatrixRank[A]

MatrixRank[DAB]

ukladKC=A.{x1,x2,x3}==B

Solve[ukladKC,{x1,x2,x3}]

(* czyli Solve nie znalazł żadnego rozwiązania *)

LinearSolve[A,B]

(* To jest poprawny wynik ! *)

Problem własny

Jeśli zachodzi
                        A.X = λ X, 
 to znaczy

(({{A11, A12, ..., A1n}, {A21, A22, ..., A2n}, {..., ..., ..., ...}, {An1, An2, ..., Ann}}))(({{X1}, {X2}, {...}, {Xn}}))=λ (({{X1}, {X2}, {...}, {Xn}}))

dla pewnej liczby λ i wektora X= (({{X1}, {X2}, {...}, {Xn}})) ≠(({{0}, {0}, {...}, {0}})) ,

to liczbę λ nazywamy wartością własną macierzy A,
a wektor X nazywamy wektorem własnym macierzy A związanym z wartością własną  λ. 

Składowe wektora własnego można uważać za rozwiązanie liniowego i jednorodnego układu n równań na n niewiadomych

                    (A - λ (Subscript[I, d]) ). X = 0,
                    
gdzie (Subscript[I, d]) jest macierzą jednostkową (identity matrix) n x n.
                    
Z twierdzenia Kroneckera-Capelli’ego wynika, że taki układ może mieć niezerowe rozwiązanie wtedy i tylko wtedy, gdy

                            Det[(A - λ (Subscript[I, d]) )]    = 0.

Równanie to nazywamy charakterystycznym (lub wiekowym), a wielomian charakterystyczny dla zmiennej λ stanowi lewą stronę tego równania.

(* Przykład 1: macierz 2 x 2  *)

Clear[X,n];

A={{0,1},{1,0}}

MatrixForm[A]

Id2=IdentityMatrix[2]

Det[A-λ*Id2]

(* To jest wielomian charakterystyczny (wiekowy) *)

rch=Det[A-λ*Id2]==0

(* To jest równanie charakterystyczne na λ *)

Solve[rch,λ]

(* Mamy dwie wartości własne i możemy poszukać "ręcznie" wektorów własnych z nimi związanych *)

(* λ=1 *)

Solve[A.{x1,x2}=={x1,x2},{x1,x2}]

(* Rozwiązanie nie może być i nie jest jednoznaczne; oznacza ono, że każdy wektor postaci {x1,x1} jest wektorem własnym A do wartości własnej 1 *)

A.{x1,x1}=={x1,x1}

(* λ= -1 *)

Solve[A.{x1,x2}==-{x1,x2},{x1,x2}]

(* Rozwiązanie oznacza, że każdy wektor postaci {x1,-x1} jest wektorem własnym A do wartości własnej -1 *)

A.{x1,-x1}==-{x1,-x1}

(* Oczywiście Mathematica może wszystko to za nas zrobić sama. Oto kilka użytecznych komend: *)

(* Wielomian charakterystyczny A w zmiennej λ *)

CharacteristicPolynomial[A,λ]

(* Wartości własne macierzy A *)

Eigenvalues[A]

(* Wektory własne macierzy A *)

Eigenvectors[A]

(* Mathematica podaje mozliwie najprostszą postać wektora własnego. Jest jasne, ze jeśli X jest wektorem własnym, to także c*X jest wektorem własnym, gdzie c jest dowolną liczbą różną od zera. *)

(* Dwie powyższe komendy można zebrać w jedną:  Eigensystem *)

pw=Eigensystem[A]

pw[[1]] (* lista wartości własnych *)

nww=Length[pw[[1]]]

pw[[2]] (* lista wektorów własnych *)

Do[Print[A.pw[[2]][[i]]==pw[[1]][[i]]*pw[[2]][[i]]],{i,1,nww}]

(* Czyli równanie własne jest spełnione dla obu wartości własnych macierzy A*)

Diagonalizacja macierzy

Dla niediagonalnej macierzy A szukamy takiej nieosobliwej macierzy U, by zachodziło
  (U^(-1 ))A U = (({{(Subscript[λ, 1]), 0, ..., 0}, {0, (Subscript[λ, 2]), ..., 0}, {0, 0, ..., 0}, {0, 0, 0, (Subscript[λ, n])}})) ≡ D. Znalezienie takiej macierzy U ogromnie ułatwia policzenie (A^n) dla dużych wartości n.  Zachodzi bowiem
  A = U D (U^(-1 )), więc (A^n) = (((U D (U^(-1)) ))^n) = U D (U^(-1)) U D (U^(-1)) ... U D (U^(-1)) U D (U^(-1)) =
  U (D^n) (U^(-1)), a liczenie (D^n) jest bardzo proste. (Dlaczego ?)

Przypadek ogólny diagonalizacji

Twierdzenie: Macierz stopnia n można zdiagonalizować wtedy i tylko wtedy, gdy ma ona n liniowo niezależnych wektorów własnych. Nie oznacza to, że wszystkie wartości własne macierzy muszą być różne, bo nawet w tym przypadku można niekiedy znaleźć niezależne liniowo wektory! Z takich wektorów własnych budujemy macierz U (kolumna obok kolumny). Zachodzi wtedy


                    (U^(-1 ))A U = (({{(Subscript[λ, 1]), 0, ..., 0}, {0, (Subscript[λ, 2]), ..., 0}, {0, 0, ..., 0}, {0, 0, 0, (Subscript[λ, n])}})) ,

gdzie (Subscript[λ, 1]), (Subscript[λ, 2]), ...,  (Subscript[λ, n]) są wartościami własnymi (niekoniecznie różnymi) macierzy A. 

Uwaga: wektory własne, z których tworzymy macierz U nie muszą być unormowane.

(* ClearAll["Global`*"]  *)

(* Przykład 2: macierz 3 x 3, która ma trzy różne wartości własne.

H. Arodź, K. Rościszewski, Algebra i geometria analityczna w zadaniach, Wydawnictwo Znak, Kraków 2006 *)

A={{1,1,1 },{0,2,2 },{0,0,3}}

MatrixForm[A]

pw=Eigensystem[A]

Um=Array[U,{3,3}];

Um[[All,1]]=pw[[2]][[1]]

(* Nie musimy normować wektorów własnych *)

Um[[All,2]]=pw[[2]][[2]]

Um[[All,3]]=pw[[2]][[3]]

MatrixForm[Simplify[Inverse[Um].A.Um]]

(* ClearAll["Global`*"] *)

(* Przykład 3: macierz 3 x 3, która ma podwójną wartość własną

http://en.wikipedia.org/wiki/Diagonalizable_matrix#Diagonalization
*)

A={{-1,3,-1 },{-3,5,-1 },{-3,3,1}}
{{-1,3,-1},{-3,5,-1},{-3,3,1}}

MatrixForm[A]

pw=Eigensystem[A]

Um=Array[U,{3,3}];

Um[[All,1]]=pw[[2]][[1]]

(* Nie musimy normować wektorów własnych *)

Um[[All,2]]=pw[[2]][[2]]

Um[[All,3]]=pw[[2]][[3]]

MatrixForm[Simplify[Inverse[Um].A.Um]]

Czy rzeczywistą macierz diagonalizuje zawsze macierz rzeczywista?

(* ClearAll["Global`*"]  *)

(* Przykład 4: macierz 2 x 2, której nie da się zdiagonalizować przy pomocy rzeczywistej macierzy. Trzeba użyć macierzy zespolonej !

http://en.wikipedia.org/wiki/Diagonalizable_matrix#Diagonalization
*)

A={{0,1},{-1,0}}

MatrixForm[A]

pw=Eigensystem[A]

Um=Array[U,{2,2}];

Um[[All,1]]=pw[[2]][[1]]

(* Nie musimy normować wektorów własnych *)

Um[[All,2]]=pw[[2]][[2]]

MatrixForm[Um]

(* Macierz diagonalizująca macierz rzeczywistą jest zespolona *)

MatrixForm[Simplify[Inverse[Um].A.Um]]

(* Macierz A ma dwa liniowo niezależne wektory własne, więc dała się jej zdiagonalizować !  *)

Czy wszystkie macierze można zdiagonalizować ?

(* ClearAll["Global`*"]  *)

(* Przykład 5: macierz, której nie da się zdiagonalizować. 

http://en.wikipedia.org/wiki/Diagonalizable_matrix#Diagonalization
*)

A={{0,1},{0,0}}

MatrixForm[A]

pw=Eigensystem[A]

Um=Array[U,{2,2}];

Um[[All,1]]=pw[[2]][[1]]

Um[[All,2]]=pw[[2]][[2]]

MatrixForm[Simplify[Inverse[Um].A.Um]]

(* Ta macierz nie ma dwóch liniowo niezależnych wektorów własnych, więc nie da się jej zdiagonalizować !  *)

Macierze hermitowskie (symetryczne)

Przypomnienie: macierze hermitowskie to takie macierze, które jednoczesna operacja transpozycji i sprzężenia zespolonego (czyli sprzężenie hermitowskie) pozostawia bez zmiany. Rzeczywista macierz hermitowska jest macierzą symetryczną.

Macierze hermitowskie grają bardzo ważną role w fizyce kwantowej!

Przytoczę najprostszą postać twierdzenia, którego ogólna postać przedstawiana jest na wykładzie z algebry.

Zakładamy, że hermitowska macierz A  (A=(SuperDagger[A])) o wymiarach n x n ma n różnych wartości własnych (Subscript[λ, 1]) < (Subscript[λ, 2]) < ... < (Subscript[λ, n]) . Niech wektory (Subscript[v, 1]), (Subscript[v, 2]), ..., (Subscript[v, n. ])będą odpowiednimi wektorami własnymi. Można pokazać, że wszystkie wartości własne są rzeczywiste oraz że zachodzi

                                     (SuperDagger[(Subscript[v, k])]) (Subscript[v, i]) = 0 dla i ≠ k, 
                                     
czyli iloczyn skalarny dwóch różnych wektorów własnych jest równy zero. Załóżmy jeszcze, że każdy wektor jest unormowany do jedności:       
                                       (SuperDagger[(Subscript[v, i])]) (Subscript[v, i]) = 1 .

Z takich wektorów własnych budujemy macierz U (kolumna obok kolumny). Zachodzi wtedy 

                    (U^(† ))A U = (({{(Subscript[λ, 1]), 0, ..., 0}, {0, (Subscript[λ, 2]), ..., 0}, {0, 0, ..., 0}, {0, 0, 0, (Subscript[λ, n])}})) .

(* ClearAll["Global`*"]  *)

(* Przykład 6: macierz hermitowska 2 x 2 *)

A={{0,I},{-I,0}}

MatrixForm[A]

pw=Eigensystem[A]

Um=Array[U,{2,2}]

v1=Normalize[pw[[2]][[1]]]

v2=Normalize[pw[[2]][[2]]]

MatrixForm[v1]

MatrixForm[v2]

Um[[All,1]]=v1

Um[[All,2]]=v2

MatrixForm[Um]

MatrixForm[ConjugateTranspose[Um].A.Um]

(* ClearAll["Global`*"]  *)

(* Przykład 7: macierz hermitowska 3 x 3 *)

A={{1,0,1-I},{0,1,3},{1+I,3,2}}

MatrixForm[A]

pw=Eigensystem[A]

nww=Length[pw[[1]]]

Do[Print[A.pw[[2]][[i]]==pw[[1]][[i]]*pw[[2]][[i]]],{i,1,nww}]

Um=Array[U,{3,3}]

v1=Normalize[pw[[2]][[1]]]

v2=Normalize[pw[[2]][[2]]]

v3=Normalize[pw[[2]][[3]]]

Simplify[Conjugate[v1].v1]

Simplify[Conjugate[v2].v1]

Simplify[Conjugate[v3].v1]

Simplify[Conjugate[v1].v2]

Simplify[Conjugate[v2].v2]

Simplify[Conjugate[v3].v2]

Simplify[Conjugate[v1].v3]

Simplify[Conjugate[v2].v3]

Simplify[Conjugate[v3].v3]

MatrixForm[v1]
MatrixForm[v2]
MatrixForm[v3]

Um[[All,1]]=v1

Um[[All,2]]=v2

Um[[All,3]]=v3

MatrixForm[Simplify[Um.ConjugateTranspose[Um]]]

(* To oznacza, że macierz U jest unitarna ! *)

pw[[1]]

MatrixForm[Simplify[ConjugateTranspose[Um].A.Um]]