 #### Matrices and Vectors in Sage

Linear algebra is a fundamental task for mathematical software. Linear algebra, i.e., dealing with vectors and matrices, is easily automated because it involves computations that must be performed according to well-defined formulas and algorithms. To do so, Sage uses two (non-Python) data types: matrix and vector. However, there are many ways to do linear algebra in Sage. We will concentrate on some of them here. Let's get started.

### Vectors and Vector Spaces

We start by defining a vector space that consists of three-element vectors with real numbers as elements: the Euclidean space $\mathbb{R}^3$. We then print the canonical basis for that vector space.
The command VectorSpace takes as first argument the base field (e.g., the reals RR, the rationals QQ, or the complex numbers CC, while the second argument is the dimension of the vector space.
︡13855390-665e-4fea-9553-106dc0a9a307︡{"done":true,"html":"\n

\n \n

#### \n Matrices and Vectors in Sage\n

\nLinear algebra is a fundamental task for mathematical software. Linear algebra, i.e., dealing with vectors and matrices, is easily automated because it involves computations that must be performed according to well-defined formulas and algorithms. To do so, Sage uses two (non-Python) data types: matrix and vector. However, there are many ways to do linear algebra in Sage. We will concentrate on some of them here. Let's get started.
\n

### \n Vectors and Vector Spaces\n

\n We start by defining a vector space that consists of three-element vectors with real numbers as elements: the Euclidean space $\\mathbb{R}^3$. We then print the canonical basis for that vector space.
The command VectorSpace takes as first argument the base field (e.g., the reals RR, the rationals QQ, or the complex numbers CC, while the second argument is the dimension of the vector space.
"} ︠23922888-0df2-482c-91a2-1d3772503366s︠ R3 = VectorSpace(RR,3) print("Basis for space:") (b1,b2,b3)=R3.basis() print b1 print b2 print b3 ︡fcb9a9f4-4b40-41fc-9ead-85b8e608a53d︡{"stdout":"Basis for space:\n"}︡{"stdout":"(1.00000000000000, 0.000000000000000, 0.000000000000000)\n"}︡{"stdout":"(0.000000000000000, 1.00000000000000, 0.000000000000000)\n"}︡{"stdout":"(0.000000000000000, 0.000000000000000, 1.00000000000000)\n"}︡{"done":true}︡ ︠0ee14ec4-87f8-4315-b9dd-b9edd87884a1i︠ %html Let use define some vectors in this vector space, and then look at their linear combination, their norm, scalar multiple, scalar product (aka "dot product") and their cross product: ︡1ed30e38-0c33-4f45-9fc0-592d532881c2︡{"done":true,"html":"\nLet use define some vectors in this vector space, and then look at their linear combination, their norm, scalar multiple, scalar product (aka \"dot product\") and their cross product: \n"} ︠da2162b8-89b5-4e6a-a84f-392970c24501s︠ vector1 = R3([-1,2,7]) vector2 = R3([4,-9,2]) print vector1 print vector2 ︡68a071c2-e831-4fb2-8554-8053ff42a58d︡{"stdout":"(-1.00000000000000, 2.00000000000000, 7.00000000000000)\n"}︡{"stdout":"(4.00000000000000, -9.00000000000000, 2.00000000000000)\n"}︡{"done":true}︡ ︠bd4b32a6-f126-4262-a709-284a8ba3d67fs︠ print("Linear combination:") a,b = var("a,b") print(a*vector1+b*vector2) print((a*vector1+b*vector2).simplify()) ︡e8d761ad-86db-4db2-854f-8d07b8a50159︡{"stdout":"Linear combination:\n"}︡{"stdout":"(-1.00000000000000*a + 4.00000000000000*b, 2.00000000000000*a - 9.00000000000000*b, 7.00000000000000*a + 2.00000000000000*b)\n"}︡{"stdout":"(-1.0*a + 4.0*b, 2.0*a - 9.0*b, 7.0*a + 2.0*b)\n"}︡{"done":true}︡ ︠6d76a2b2-4750-489a-bbcf-43849dd0e51bs︠ print("Norm of vector1:") print vector1.norm() print sqrt(vector1 * vector1) # using scalar product ︡11f525be-df9f-44bd-a30a-a51e469461cc︡{"stdout":"Norm of vector1:\n"}︡{"stdout":"7.34846922834953\n"}︡{"stdout":"7.34846922834953\n"}︡{"done":true}︡ ︠72d1f642-d127-4539-809b-92330ca37ea7s︠ print("Scalar multiplication:") print 2 * vector1 ︡55f4ed49-020d-4160-b816-40c4e236b69b︡{"stdout":"Scalar multiplication:\n"}︡{"stdout":"(-2.00000000000000, 4.00000000000000, 14.0000000000000)\n"}︡{"done":true}︡ ︠6322eb6d-8d96-46c2-9ddf-1e2d91596a62s︠ print("Scalar (dot) products:") print vector1 * vector2 # using operator * print vector1.inner_product(vector2) # using methods print b1 * b2 # scalar product of basis vectors is 0 ︡0b9dcbd2-32d2-4bac-94de-acb1a269cd4f︡{"stdout":"Scalar (dot) products:\n"}︡{"stdout":"-8.00000000000000\n"}︡{"stdout":"-8.00000000000000\n"}︡{"stdout":"0.000000000000000\n"}︡{"done":true}︡ ︠eb7bebee-9755-4f45-a75b-fd3e31045925s︠ print("Vector (cross) product:") print vector1.cross_product(vector2) ︡2ed4412a-0384-468b-b692-8c0ee50c70a5︡{"stdout":"Vector (cross) product:\n"}︡{"stdout":"(67.0000000000000, 30.0000000000000, 1.00000000000000)\n"}︡{"done":true}︡ ︠7e801c93-cfbd-4951-a096-1ad78894b964i︠ %html You can also create a vector by specifying the base field directly (so without the need of the VectorSpace command first) as follows: ︡5a88ba1c-a521-4894-b661-0bf8dc54c66f︡{"done":true,"html":"\n You can also create a vector by specifying the base field directly (so without the need of the VectorSpace command first) as follows:\n"} ︠1f73483c-cb6c-497b-a005-23f58c9113cds︠ u1 = vector(QQ, [1,2/7,10/3]) print u1 u2 = vector(RR, [1,2/7,10/3]) print u2 ︡c1608117-ac8c-4a73-beaf-2b12d86f43fd︡{"stdout":"(1, 2/7, 10/3)\n"}︡{"stdout":"(1.00000000000000, 0.285714285714286, 3.33333333333333)\n"}︡{"done":true}︡ ︠e8391698-154b-46bf-93e7-3a7056d6a3a9i︠ %html Problem: Since vectors are essentially lists, you can "slice" them, get single elements, and change elements as you are used to from lists. So, define a four-dimensional vector with components $2$, $7/12$, $4/5$ and $3$, respectively. Print the vector consisting of only the first three components, assign a different value to the second and third component and print the resulting new vector. ︡fc5ad637-047a-42d5-93df-ed1b5bf0a5d4︡{"done":true,"html":"\nProblem: Since vectors are essentially lists, you can \"slice\" them, get single elements, and change elements as you are used to from lists. So, define a four-dimensional vector with components $2$, $7/12$, $4/5$ and $3$, respectively. Print the vector consisting of only the first three components, assign a different value to the second and third component and print the resulting new vector. \n"} ︠bcd1db7e-c11b-4b38-ba46-1a806f4d0b00︠ ︡63540328-06a1-401c-ba6b-6544a6e27f8a︡ ︠633a864f-b1a2-4d87-9594-1a27ef3b727e︠ ︡1088ee04-22b6-4e00-b83b-212834dd3c0c︡ ︠37f3bde7-c03f-41a1-8870-e892cbb52a02i︠ %html

### Matrices

Matrices work similar to vectors in Sage, only that instead of a list we now have a list of list. Instead of VectorSpace for vectors, we now use MatrixSpace.
︡382bd540-6f47-4354-94d7-4fe043be6a77︡{"done":true,"html":"\n

### \n Matrices\n

\nMatrices work similar to vectors in Sage, only that instead of a list we now have a list of list. Instead of VectorSpace for vectors, we now use MatrixSpace. \n
"} ︠ae07a6b2-5319-4016-ae73-ffd151827469s︠ M4 = MatrixSpace(RR,4) print("Identity matrix:") show(M4.identity_matrix()) ︡0c9928d4-752e-4089-a39c-589732623a4a︡{"stdout":"Identity matrix:\n"}︡{"html":"
$\\displaystyle \\left(\\begin{array}{rrrr}\n1.00000000000000 & 0.000000000000000 & 0.000000000000000 & 0.000000000000000 \\\\\n0.000000000000000 & 1.00000000000000 & 0.000000000000000 & 0.000000000000000 \\\\\n0.000000000000000 & 0.000000000000000 & 1.00000000000000 & 0.000000000000000 \\\\\n0.000000000000000 & 0.000000000000000 & 0.000000000000000 & 1.00000000000000\n\\end{array}\\right)$
"}︡{"done":true}︡ ︠84f72225-7f85-4287-8067-190c2fe899b1s︠ A = M4.matrix([[0,-1,-1,1],[1,1,1,1],[2,4,1,-2],[3,1,-2,2]]) print("Matrix A:") show(A) ︡077a028c-fde8-47d6-971e-2ef4b73a7454︡{"stdout":"Matrix A:\n"}︡{"html":"
$\\displaystyle \\left(\\begin{array}{rrrr}\n0.000000000000000 & -1.00000000000000 & -1.00000000000000 & 1.00000000000000 \\\\\n1.00000000000000 & 1.00000000000000 & 1.00000000000000 & 1.00000000000000 \\\\\n2.00000000000000 & 4.00000000000000 & 1.00000000000000 & -2.00000000000000 \\\\\n3.00000000000000 & 1.00000000000000 & -2.00000000000000 & 2.00000000000000\n\\end{array}\\right)$
"}︡{"done":true}︡ ︠ee190978-c525-4fce-b7dc-db704e93e40ci︠ %html Matrices are great to solve simultaneous linear equations like $$\left( \begin{array}{cccc} 0 & -1 & -1 & 1 \\ 1 & 1 & 1 & 1 \\ 2 & 4 & 1 & -2 \\ 3 & 1 & -2 & 2 \end{array} \right) \left(\begin{array}{c} x_1 \\ x_2 \\ x_3 \\ x_4 \end{array}\right) = \left(\begin{array}{c} 0 \\ 6 \\ -1 \\ 3 \end{array}\right)$$ Here is how you can do it in Sage: ︡7318d4fb-fe65-41be-b6b7-b7132bca10a2︡{"done":true,"html":"\nMatrices are great to solve simultaneous linear equations like\n $$\\left( \\begin{array}{cccc} 0 & -1 & -1 & 1 \\\\ 1 & 1 & 1 & 1 \\\\ 2 & 4 & 1 & -2 \\\\ 3 & 1 & -2 & 2 \\end{array} \\right) \\left(\\begin{array}{c} x_1 \\\\ x_2 \\\\ x_3 \\\\ x_4 \\end{array}\\right) = \\left(\\begin{array}{c} 0 \\\\ 6 \\\\ -1 \\\\ 3 \\end{array}\\right)$$\nHere is how you can do it in Sage:\n"} ︠579ada50-1488-4895-ad32-9e0d58de0150s︠ A = M4.matrix([[0,-1,-1,1],[1,1,1,1],[2,4,1,-2],[3,1,-2,2]]) print("Matrix A:") show(A) b = vector(RR,[0,6,-1,3]) print("Vector b:") show(b) solution = A.solve_right(b) print("Solution to Ax=b:") show(solution) ︡caeb7a07-d087-4a51-88a8-d90b6f4a3b4b︡{"stdout":"Matrix A:\n"}︡{"html":"
$\\displaystyle \\left(\\begin{array}{rrrr}\n0.000000000000000 & -1.00000000000000 & -1.00000000000000 & 1.00000000000000 \\\\\n1.00000000000000 & 1.00000000000000 & 1.00000000000000 & 1.00000000000000 \\\\\n2.00000000000000 & 4.00000000000000 & 1.00000000000000 & -2.00000000000000 \\\\\n3.00000000000000 & 1.00000000000000 & -2.00000000000000 & 2.00000000000000\n\\end{array}\\right)$
"}︡{"stdout":"Vector b:\n"}︡{"html":"
$\\displaystyle \\left(0.000000000000000,\\,6.00000000000000,\\,-1.00000000000000,\\,3.00000000000000\\right)$
"}︡{"stdout":"Solution to Ax=b:\n"}︡{"html":"
$\\displaystyle \\left(2.00000000000000,\\,-1.00000000000000,\\,3.00000000000000,\\,2.00000000000000\\right)$
"}︡{"done":true}︡ ︠b20a5c51-89d3-4bc3-8cdd-038c2d66f1d0i︠ %html As with vectors, there are several ways to create matrices. Have a look at the following possibilities. ︡a4f0a15b-fb5c-44d4-a5c7-1b7109974f22︡{"done":true,"html":"\n As with vectors, there are several ways to create matrices. Have a look at the following possibilities.\n"} ︠7c3936ab-9e67-4b29-8d6a-bebca8b35bc7s︠ B = matrix(QQ,[[1,2,3,4],[2,3,1,2],[-1,-1,-1,-1]]) show(B) C = matrix(RR, 3, 4, [1,2,3,4,5,6,7,8,9,10,11,12]) show(C) ︡faa7e2f7-02f1-4ddb-bade-bbd0bc1471e3︡{"html":"
$\\displaystyle \\left(\\begin{array}{rrrr}\n1 & 2 & 3 & 4 \\\\\n2 & 3 & 1 & 2 \\\\\n-1 & -1 & -1 & -1\n\\end{array}\\right)$
"}︡{"html":"
$\\displaystyle \\left(\\begin{array}{rrrr}\n1.00000000000000 & 2.00000000000000 & 3.00000000000000 & 4.00000000000000 \\\\\n5.00000000000000 & 6.00000000000000 & 7.00000000000000 & 8.00000000000000 \\\\\n9.00000000000000 & 10.0000000000000 & 11.0000000000000 & 12.0000000000000\n\\end{array}\\right)$
"}︡{"done":true}︡ ︠82fab370-9190-4e57-a8af-84ae3cd2340fi︠ %html Problem: Solve the linear equation $$\left( \begin{array}{ccc} 0 & -1 & 1 \\ 1 & 1 & 1 \\ 2 & 4 & 8 \end{array} \right) \left(\begin{array}{c} x_1 \\ x_2 \\ x_3 \end{array}\right) = \left(\begin{array}{c} 0 \\ 6 \\ 1 \end{array}\right)$$ ︡1efd78fa-c539-43cf-a711-e04d320a01a7︡{"done":true,"html":"\n Problem: Solve the linear equation \n $$\\left( \\begin{array}{ccc} 0 & -1 & 1 \\\\ 1 & 1 & 1 \\\\ 2 & 4 & 8 \\end{array} \\right) \\left(\\begin{array}{c} x_1 \\\\ x_2 \\\\ x_3 \\end{array}\\right) = \\left(\\begin{array}{c} 0 \\\\ 6 \\\\ 1 \\end{array}\\right)$$\n"} ︠7b1b6a37-b306-4133-bdcf-6a8de43172b4︠ ︡8dfa0bc9-b3e0-4fda-b2b0-99b6c98d94c4︡ ︠2ee8c79c-027f-4473-9615-824721d561e4︠ ︡f940fec3-43bb-4e2e-bf2e-c665417d4e61︡ ︠094c05f0-a124-492c-965a-b5db806c4a77i︠ %html We can also use the reduced row echelon form, look at the following example. ︡d9294e8e-c980-4ca2-a638-9786b735e813︡{"done":true,"html":"\n We can also use the reduced row echelon form, look at the following example.\n"} ︠8a035787-e4f2-42ec-81cd-d6096cd98277s︠ A = matrix(RR, 3,3,[3,-4,5,1,1,-8,2,1,1]) show(A) ︡ec6c916d-4ae5-4e2b-8577-846abd536fe9︡{"html":"
$\\displaystyle \\left(\\begin{array}{rrr}\n3.00000000000000 & -4.00000000000000 & 5.00000000000000 \\\\\n1.00000000000000 & 1.00000000000000 & -8.00000000000000 \\\\\n2.00000000000000 & 1.00000000000000 & 1.00000000000000\n\\end{array}\\right)$
"}︡{"done":true}︡ ︠d64e6be8-9c89-4528-992d-7b218a8c323ds︠ b = vector(RR,[14,-5,7]) show(b) ︡0af40e43-1da9-41ef-9fd7-0a75363a1626︡{"html":"
$\\displaystyle \\left(14.0000000000000,\\,-5.00000000000000,\\,7.00000000000000\\right)$
"}︡{"done":true}︡ ︠c728a9b8-d0d0-49e8-9183-a3d504626183s︠ B = A.augment(b) show(B) ︡6c9c82f0-1dab-4a72-903a-44d5ec99d0e0︡{"html":"
$\\displaystyle \\left(\\begin{array}{rrrr}\n3.00000000000000 & -4.00000000000000 & 5.00000000000000 & 14.0000000000000 \\\\\n1.00000000000000 & 1.00000000000000 & -8.00000000000000 & -5.00000000000000 \\\\\n2.00000000000000 & 1.00000000000000 & 1.00000000000000 & 7.00000000000000\n\\end{array}\\right)$
"}︡{"done":true}︡ ︠8766bd75-8626-4de3-a0f2-f40ffcb19e8bs︠ print("B in echelon form:") show(B.rref()) ︡be1169c6-f143-442c-8d8a-45abbea2ac9b︡{"stdout":"B in echelon form:\n"}︡{"html":"
$\\displaystyle \\left(\\begin{array}{rrrr}\n1.00000000000000 & 0.000000000000000 & 0.000000000000000 & 3.00000000000000 \\\\\n0.000000000000000 & 1.00000000000000 & 0.000000000000000 & 8.88178419700125 \\times 10^{-16} \\\\\n0.000000000000000 & 0.000000000000000 & 1.00000000000000 & 1.00000000000000\n\\end{array}\\right)$
"}︡{"done":true}︡ ︠33fa08cd-fddf-4aa2-9457-28752b891293i︠ %html How do you interpret this last result?

Problem: Experiment with linear equations using the two methods above. What happens in the case that the system has no solution or infinitely many solution?
︡d80f5f89-f0ee-4651-9e07-ce094d964b7f︡{"done":true,"html":"\n How do you interpret this last result?

\n Problem: Experiment with linear equations using the two methods above. What happens in the case that the system has no solution or infinitely many solution?\n
"} ︠c37da5df-9400-4a67-bcae-845e1e36e3c6︠ ︡1c7e621d-af57-4270-9118-30cf2ed99495︡ ︠19ac8ecc-49c5-46cf-8fa1-a68bde068ddd︠ ︡32821f81-c606-48fb-ae49-11f7fc4ebabb︡ ︠13abb1f9-5782-411f-8efa-3e1780f82cf7︠ ︡4eb41326-f430-4ac1-b3a9-437f224a4af3︡ ︠eab67bc8-810c-44b0-bc54-39e53dd8717b︠ %html We can also use the standard operations on matrices, adding and multiplying them. ︡d937488f-8712-4c6d-954d-fd1c381e1564︡{"done":true,"html":"\n We can also use the standard operations on matrices, adding and multiplying them.\n"} ︠c960c9b0-63fd-4e8e-ba85-280717084541s︠ A = matrix(QQ,[[3,2,1],[4,5,6]]) B = matrix(QQ,[[2,2,2],[1,2,3]]) print("Matrix addition:") show(A+B) ︡622d1f42-0673-4652-8282-2b8a99820bd7︡{"stdout":"Matrix addition:\n"}︡{"html":"
$\\displaystyle \\left(\\begin{array}{rrr}\n5 & 4 & 3 \\\\\n5 & 7 & 9\n\\end{array}\\right)$
"}︡{"done":true}︡ ︠c6936831-34bf-4528-8dfe-753b358ade2cs︠ print("Scalar multiplication:") show(1/2*A) ︡b6842fa2-dd92-4be0-8ff9-66642dccd55b︡{"stdout":"Scalar multiplication:\n"}︡{"html":"
$\\displaystyle \\left(\\begin{array}{rrr}\n\\frac{3}{2} & 1 & \\frac{1}{2} \\\\\n2 & \\frac{5}{2} & 3\n\\end{array}\\right)$
"}︡{"done":true}︡ ︠f39a9c00-d64f-42bf-a109-00daa7ff1461s︠ a,b,c,d,e,f = var("a,b,c,d,e,f") C = matrix(QQ, [[4,2,1],[5,3,7]]) D = matrix(SR, [[a,b],[c,d],[e,f]]) print("Matrix multiplication:") show(C * D) ︡43a6a60b-6b5b-4185-bebf-89bb4de32367︡{"stdout":"Matrix multiplication:\n"}︡{"html":"
$\\displaystyle \\left(\\begin{array}{rr}\n4 \\, a + 2 \\, c + e & 4 \\, b + 2 \\, d + f \\\\\n5 \\, a + 3 \\, c + 7 \\, e & 5 \\, b + 3 \\, d + 7 \\, f\n\\end{array}\\right)$
"}︡{"done":true}︡ ︠c7f4d67b-3682-41b5-86ab-2667db0a05ac︠ x1, x2, x3 = var("x1,x2,x3") X = vector([x1,x2,x3]) print("Multiplying a matrix and a vector:") show(C * X) ︡19c8e2ef-d691-499c-8d0c-f65ef0cd2e65︡{"stdout":"Multiplying a matrix and a vector:\n"}︡{"html":"
$\\displaystyle \\left(4 \\, x_{1} + 2 \\, x_{2} + x_{3},\\,5 \\, x_{1} + 3 \\, x_{2} + 7 \\, x_{3}\\right)$
"}︡{"done":true}︡ ︠83567047-4275-446b-a49f-f2e47ad9f627i︠ %html Sage can find the transpose, the inverse or the determinant of a matrix. ︡2c873ee4-6411-4e0b-af3b-7b0651c9ca60︡{"done":true,"html":"\n Sage can find the transpose, the inverse or the determinant of a matrix.\n"} ︠ccbad577-abdc-48c9-b249-8626fe113563s︠ A = matrix(QQ,[[2,5,4],[3,1,2],[5,4,6]]) print("Matrix A:") show(A) ︡87f700c9-4d83-44a7-a2de-8eb9e8498383︡{"stdout":"Matrix A:\n"}︡{"html":"
$\\displaystyle \\left(\\begin{array}{rrr}\n2 & 5 & 4 \\\\\n3 & 1 & 2 \\\\\n5 & 4 & 6\n\\end{array}\\right)$
"}︡{"done":true}︡ ︠8c94fcb7-f523-4267-8354-c8e688d41b3es︠ print("Determinant of A:") print A.det() print det(A) print("Rank of A:") print A.rank() ︡42d49eea-89e0-498f-8354-bf1bdce4d2ec︡{"stdout":"Determinant of A:\n"}︡{"stdout":"-16\n"}︡{"stdout":"-16\n"}︡{"stdout":"Rank of A:\n"}︡{"stdout":"3\n"}︡{"done":true}︡ ︠951ae80e-66c8-4b51-9eb6-338b46f564b5s︠ print("Transpose of A:") show(A.transpose()) print("Inverse of A:") show(A.inverse()) print("Check:") show(A.inverse() * A) print("Adjoint of A:") show(A.adjoint()) ︡d7d8cc86-5bfd-444a-915e-9d0cf8567ee0︡{"stdout":"Transpose of A:\n"}︡{"html":"
$\\displaystyle \\left(\\begin{array}{rrr}\n2 & 3 & 5 \\\\\n5 & 1 & 4 \\\\\n4 & 2 & 6\n\\end{array}\\right)$
"}︡{"stdout":"Inverse of A:\n"}︡{"html":"
$\\displaystyle \\left(\\begin{array}{rrr}\n\\frac{1}{8} & \\frac{7}{8} & -\\frac{3}{8} \\\\\n\\frac{1}{2} & \\frac{1}{2} & -\\frac{1}{2} \\\\\n-\\frac{7}{16} & -\\frac{17}{16} & \\frac{13}{16}\n\\end{array}\\right)$
"}︡{"stdout":"Check:\n"}︡{"html":"
$\\displaystyle \\left(\\begin{array}{rrr}\n1 & 0 & 0 \\\\\n0 & 1 & 0 \\\\\n0 & 0 & 1\n\\end{array}\\right)$
$\\displaystyle \\left(\\begin{array}{rrr}\n-2 & -14 & 6 \\\\\n-8 & -8 & 8 \\\\\n7 & 17 & -13\n\\end{array}\\right)$