Tutorial de Álgebra Lineal
Sage Days 15
Mayo 16, 2009

Robert A. Beezer

Universidad de Puget Sound

Gráfico: Gidon Eshel, Univ. de Chicago

Resolviendo Sistema de Ecuaciones

Crear una matriz y un vector.

$A$ es una matriz no singular $3\times 3$

$b$ es un 3-slot vector.

{{{id=0| A = matrix([[2,1,1/3],[-1,6,2],[1/2,1,8]]) A /// [ 2 1 1/3] [ -1 6 2] [1/2 1 8] }}} {{{id=3| b = vector([14/3, 2, -6]) b /// (14/3, 2, -6) }}}

Resolver $Ax=b$.

Solve commands ("right" is location of solution vector).

{{{id=4| A.solve_right(b) /// (2, 1, -1) }}} {{{id=7| print A.det() print A.det() == 0 /// 299/3 False }}} {{{id=5| A.inverse() /// [ 6/13 -1/13 0] [ 27/299 95/598 -1/23] [-12/299 -9/598 3/23] }}}

Vectors are rows or columns as appropriate, compute $A^{-1}b$

{{{id=8| A.inverse()*b /// (2, 1, -1) }}}

Can "divide" by a matrix

{{{id=9| b/A /// (770/299, 44/897, -20/23) }}} {{{id=12| b/A.transpose() /// (2, 1, -1) }}}

Augment and row-reduce

{{{id=16| R = A.augment(matrix(b).transpose()) R /// [ 2 1 1/3 14/3] [ -1 6 2 2] [ 1/2 1 8 -6] }}} {{{id=13| R.echelon_form() /// [ 1 0 0 2] [ 0 1 0 1] [ 0 0 1 -1] }}}

Properties of Matrices

$B$ is a $6\times 5$ matrix of rank 4 over the rationals $\QQ$

{{{id=19| B = matrix(QQ,[ [10,0,3,8,7], [-16,-1,-4,-10,-13], [-6,1,-3,-6,-6], [0,2,-2,-3,-2], [3,0,1,2,3], [-1,-1,1,1,0] ]) B /// [ 10 0 3 8 7] [-16 -1 -4 -10 -13] [ -6 1 -3 -6 -6] [ 0 2 -2 -3 -2] [ 3 0 1 2 3] [ -1 -1 1 1 0] }}} {{{id=22| B.right_kernel() /// Vector space of degree 5 and dimension 1 over Rational Field Basis matrix: [ 1 -3/2 1/2 -1 -1/2] }}} {{{id=23| B.right_kernel().basis() /// [ (1, -3/2, 1/2, -1, -1/2) ] }}} {{{id=25| B.left_kernel().basis() /// [ (1, 0, 3, -1, 3, 1), (0, 1, -2, 1, 1, -1) ] }}} {{{id=26| B.row_space() /// Vector space of degree 5 and dimension 4 over Rational Field Basis matrix: [ 1 0 0 0 2] [ 0 1 0 0 -3] [ 0 0 1 0 1] [ 0 0 0 1 -2] }}} {{{id=27| B.column_space() /// Vector space of degree 6 and dimension 4 over Rational Field Basis matrix: [ 1 0 0 0 -1/4 -1/4] [ 0 1 0 0 -1/4 3/4] [ 0 0 1 0 -1/4 -9/4] [ 0 0 0 1 0 1] }}} {{{id=28| print "Rank", B.rank() print "Right Nullity", B.right_nullity() print "Columns", B.ncols() print print "Rank", B.rank() print "Left Nullity", B.left_nullity() print "Rows", B.nrows() /// Rank 4 Right Nullity 1 Columns 5 Rank 4 Left Nullity 2 Rows 6 }}}

$C$ is a random $50\times 50$ matrix over $\QQ$

{{{id=35| C = random_matrix(QQ, 50, num_bound=10, den_bound=10) C /// 50 x 50 dense matrix over Rational Field (type 'print C.str()' to see all of the entries) }}}

Column 5 of the matrix.

{{{id=39| C[4] /// (3/2, 2/5, 0, 8/5, -2/5, -4/7, 7/4, -3/5, -2/3, 1/9, 7, 2, 0, 5/4, -7/2, -10/9, 5/7, 2, 0, 1, -9/2, 1/5, 7/2, 1, -10/3, -1, -9, -10, -5/3, 0, -1/3, -9/10, -2/5, -7/6, 5/8, 0, 3/4, -2, -3/2, -5/3, -1/3, -2, -1, 1, 1/10, 1/5, 5/4, 9/10, -8, 9/7) }}}

Python slicing for entries, show() for looks

{{{id=29| show(C[20:30, 25:35]) ///

\left(\begin{array}{rrrrrrrrrr} -\frac{7}{6} & 1 & \frac{8}{5} & \frac{5}{3} & -\frac{4}{7} & \frac{5}{3} & -\frac{9}{7} & -\frac{7}{4} & -5 & -\frac{5}{4} \\ -\frac{7}{8} & -\frac{4}{3} & -7 & -\frac{9}{4} & \frac{3}{5} & -\frac{2}{7} & \frac{10}{3} & -1 & 6 & -\frac{7}{9} \\ 0 & -2 & 1 & -\frac{5}{4} & \frac{1}{10} & 0 & 5 & 1 & -\frac{4}{3} & -\frac{7}{6} \\ \frac{5}{4} & 0 & \frac{4}{9} & \frac{5}{4} & 10 & -\frac{3}{8} & \frac{3}{4} & \frac{3}{2} & \frac{5}{9} & -2 \\ \frac{5}{4} & -7 & -\frac{2}{9} & 1 & -\frac{8}{9} & \frac{10}{3} & -5 & -\frac{7}{10} & 1 & -\frac{10}{3} \\ 3 & \frac{2}{5} & \frac{5}{6} & 0 & \frac{4}{9} & -\frac{1}{4} & -\frac{3}{5} & -\frac{1}{3} & 7 & \frac{10}{3} \\ -\frac{7}{9} & 5 & 2 & \frac{6}{7} & -\frac{5}{8} & -\frac{5}{4} & \frac{3}{2} & \frac{1}{2} & \frac{10}{7} & \frac{8}{9} \\ 1 & \frac{3}{4} & \frac{3}{10} & 4 & 0 & -\frac{1}{2} & \frac{3}{4} & 5 & 5 & 0 \\ 1 & 0 & \frac{1}{2} & \frac{7}{10} & 5 & -\frac{3}{2} & \frac{3}{2} & 0 & 0 & -\frac{5}{2} \\ -\frac{7}{5} & \frac{3}{4} & -\frac{7}{4} & 0 & -\frac{2}{3} & -\frac{9}{4} & 1 & \frac{10}{9} & -\frac{3}{5} & -\frac{5}{7} \end{array}\right)

}}} {{{id=32| C.det() /// -76244937800987769506792361955814039582568713211026336805667315205782315518927329287001287987018970228899180767535170733449091911372851595339602228326815772756109043116232618772908486501980400082698339828069631891/896263010485092356716641366987847676083715121677435387928825553754389078504362694654236659984035455684732518400000000000000000000000000000000000000000000000000 }}} {{{id=33| C.trace() /// 9187/2520 }}} {{{id=44| C.norm() /// }}}

EigenStuff

$D$ is a $4\times 4$ matrix, not diagonalizable

{{{id=45| D = matrix(QQ, [ [-2,1,-2,-4], [12,1,4,9], [6,5,-2,-4], [3,-4,5,10] ]) show(D) ///

\left(\begin{array}{rrrr} -2 & 1 & -2 & -4 \\ 12 & 1 & 4 & 9 \\ 6 & 5 & -2 & -4 \\ 3 & -4 & 5 & 10 \end{array}\right)

}}}

Build a characteristic polynomial for $D$ (then simplify)

{{{id=48| var('t') S = D - t*identity_matrix(4) print S print p(t)=S.det() p(t) /// [ -t - 2 1 -2 -4] [ 12 -t + 1 4 9] [ 6 5 -t - 2 -4] [ 3 -4 5 -t + 10] t^4 - 7*t^3 + 18*t^2 - 20*t + 8 }}} {{{id=49| p.find_root(-10, 10) /// 1.9999948114026529 }}} {{{id=52| show(p.factor()) ///

{\left(t - 2\right)}^{3} {\left(t - 1\right)}

}}}

Easier, but less instructive$\dots$

{{{id=53| q(t) = D.charpoly('t') q /// t |--> (((t - 7)*t + 18)*t - 20)*t + 8 }}} {{{id=55| D.eigenvalues() /// [1, 2, 2, 2] }}} {{{id=56| diag, evecs = D.eigenmatrix_right() print "Diagonal matrix with eigenvalues" print diag print print "Matrix with eigenvectors as columns" print evecs /// Diagonal matrix with eigenvalues [1 0 0 0] [0 2 0 0] [0 0 2 0] [0 0 0 2] Matrix with eigenvectors as columns [ 1 1 0 0] [-3 -2 0 0] [-3 1 0 0] [ 0 -2 0 0] }}}

Jordan Canonical Form

{{{id=57| D.jordan_form() /// [1|0 0 0] [-+-----] [0|2 1 0] [0|0 2 1] [0|0 0 2] }}} {{{id=60| J, P = D.jordan_form(transformation = True) print J print print P /// [1|0 0 0] [-+-----] [0|2 1 0] [0|0 2 1] [0|0 0 2] [ 1 5 -8 1] [ -3 -10 21 0] [ -3 5 2 0] [ 0 -10 11 1] }}}

Check the results

{{{id=61| P^(-1)*D*P /// [1 0 0 0] [0 2 1 0] [0 0 2 1] [0 0 0 2] }}}

Manipulate the basis for the Jordan form representation

{{{id=66| Q = P.transpose().gram_schmidt()[0].transpose() Q /// [ 1 75/19 20/29 15/14] [ -3 -130/19 70/87 0] [ -3 155/19 -50/87 5/14] [ 0 -10 -65/87 5/7] }}}

Orthogonal?

{{{id=67| Q.transpose()*Q /// [ 19 0 0 0] [ 0 4350/19 0 0] [ 0 0 175/87 0] [ 0 0 0 25/14] }}}

Resulting representation?

{{{id=68| Q^(-1)*D*Q /// [ 1 20/19 -2915/1653 -535/133] [ 0 2 1 -698/609] [ 0 0 2 1] [ 0 0 0 2] }}}

Decompositions

LU, QR, SVD, Jordan Canonical Form, Smith Normal Form, Cholesky Decomposition, $\dots$

Convert $D$ to a matrix $E$ over the reals (double precision), $\RR$, to obtain QR decomposition

{{{id=69| E = D.change_ring(RDF) E /// [-2.0 1.0 -2.0 -4.0] [12.0 1.0 4.0 9.0] [ 6.0 5.0 -2.0 -4.0] [ 3.0 -4.0 5.0 10.0] }}} {{{id=81| ortho, triangular = E.QR() print "Orthogonal" print ortho print print "Upper Triangular" print triangular print /// Orthogonal [ -0.14396315015 -0.206755085285 0.83647166171 0.486664263392] [ 0.863778900898 0.11873886424 0.366775046298 -0.324442842262] [ 0.431889450449 -0.661782341253 -0.372388950068 0.486664263392] [ 0.215944725225 0.710772502024 -0.164674510583 0.648885684523] Upper Triangular [ 13.8924439894 2.0154841021 3.95898662912 8.78175215913] [ -0.0 -6.24001793541 5.76589282015 11.6505245045] [ -0.0 0.0 -0.284437791007 -0.202100535715] [ -0.0 -0.0 -0.0 -0.324442842262] }}}

Checks

{{{id=85| (ortho.transpose()*ortho - identity_matrix(4)).norm() /// 5.57101754397e-16 }}} {{{id=82| (ortho*triangular - E).norm() /// 5.59659361322e-15 }}}

Vector Spaces

Sage is so much more than numerical computation.

Can work naturally with vector spaces and modules over a variety of fields and rings.

{{{id=84| F. = FiniteField(3^2) /// }}} {{{id=89| F /// Finite Field in a of size 3^2 }}}

$V$ is a 3-dimensional vector space over $F$

"Generator" of all 2-D subspaces of $V$

{{{id=93| subs = V.subspaces(2) /// }}}

Python "list comprehension"

{{{id=95| all_subs = [U for U in subs] /// }}}

Grab one of the subspaces, #42

{{{id=101| all_subs[42] /// Vector space of degree 3 and dimension 2 over Finite Field in a of size 3^2 Basis matrix: [ 1 0 2*a + 2] [ 0 1 2] }}}

How many such subspaces?

How many basis matrices in echelon_form?

{{{id=104| all_subs[86] /// Vector space of degree 3 and dimension 2 over Finite Field in a of size 3^2 Basis matrix: [1 a 0] [0 0 1] }}} {{{id=97| len(all_subs) /// }}}

$=(3^2)^2 + (3^2)^1 + 1$

Accuracy

Octave and Matlab emphasize numerical results - everything is a floating point number.

Their rational form is deceptive. Example by William Stein:

octave:1> format rat;
octave:2> a = [-86/17,40/29,-68/43,-20/11;-24/17,-1/38,-2/25,49/17]
a =
-86/17 40/29 -68/43 -20/11
-24/17 -1/38 -2/25 49/17
octave:3> rref(a)
ans =
1 0 155/2122 -725/384
0 1 -152/173 -6553/795

and in Matlab:

>> format rat;
>> a = [-86/17,40/29,-68/43,-20/11;-24/17,-1/38,-2/25,49/17]

a =
-86/17 40/29 -68/43 -20/11
-24/17 -1/38 -2/25 49/17

>> rref(a)

ans =
1 0 13/178 -725/384
0 1 -152/173 -1426/173

Now in Sage:

{{{id=107| F = matrix(2,[-86/17, 40/29, -68/43, -20/11, -24/17, -1/38, -2/25, 49/17]) show(F.echelon_form()) /// }}}

Entry in lower right corner:

{{{id=130| print N(-6553/795, digits = 9), " Octave" print N(-1426/173, digits=9), " Matlab" print N(-30037214/3644069, digits = 9), " Sage" /// }}}

Speed

Matrices with symbolic entries

{{{id=54| var('x y') n=6 entries = [x^i-y^j+i+j for i in range(1,n+1) for j in range(1,n+1)] G = matrix(SR, n, entries) G /// }}} {{{id=114| G.det() /// }}} {{{id=115| G.det().simplify_full() /// }}} {{{id=116| time G.det() /// }}}