The NumPy or SciPy packages are very helpful for solving linear systems in Python. For instructions on how to install them look here.

We will use mostly NumPy when working with linear algebra. We could just as well have used SciPy, which contains all the functions that are in NumPy, as well as other more advanced functions that are not included in NumPy.

Want to get better at programming, check out our page Tips and tricks for coding.

Page content

Matrices in Python

Like many other programming languages, the indices in Python arrays starts at 0. This is commonly known, but may catch you by surprise if you're used to the indices in e.g. Matlab.


In Python, matrices are not it's own thing, but rather a list of list/nested lists. Let's look at an example:

Our matrix of choice will be To represent this matrix in Python, we will consider the matrix as two independent lists,  and , put together in one list. Written in code, it looks like this:

A = [[1, 2, 3], [4, 5, 6]]
# To make the matrices easier to read while coding, we can place the different lists vertical to eachother instead of horizontally.
B = [[1, 2, 3],
	   [4, 5, 6]]


# Printing both of these matrices will result in the same output:
[[1, 2, 3], [4, 5, 6]]

# To alter the output to a more readable format, for-loops are very helpful:
[[1, 2, 3], 
[4, 5, 6]]

Now, let's see how to obtain different values from our matrix:

A = [[1, 2, 3], 
	[4, 5, 6]]

print("A =", A) 				# Entire matrix
print("A[1] =", A[1])      		# 2nd row
print("A[1][2] =", A[1][2])   	# 3rd element of 2nd row
print("A[0][-1] =", A[0][-1])   # Last element of 1st Row

column = [];        			# Empty list
for row in A:
  column.append(row[2])   		# Adding the element in the 3rd column of every row.

print("3rd column =", column)

The script above will result in the following output:

A = [[1, 2, 3], [4, 5, 6]]
A[1] = [4, 5, 6]
A[1][2] = 6
A[0][-1] = 3
3rd column = [3, 6]

Matrices in Python using NumPy

First of all, we have to import the NumPy library. For the sake of this tutorial, we will only do it once, but you can assume it is done everywhere NumPy is being used.

import numpy as np

A lot of great, more detailed information about the different functions in NumPy can be found by searching in their user manual. This wiki will only describe a limited number of functions and their functionalities.

To make an array (or matrix) using NumPy, we will use the function numpy.array, and simply use the same syntax as before, but now as a function parameter.

There is another NumPy function for making matrices, numpy.matrix, but this is no longer recommended to use. The class may be removed in the future.


A = np.array([[1, 2, 3], [4, 5, 6]])		# Array of ints

B = np.array([[1.1, 2, 3], [4, 5, 6]])		# Array of floats

C = np.array([[1, 2, 3], [4, 5, 6]], dtype=complex)		# Array of complex numbers

When printing the matrices above using print, the output will look like this:

[[1 2 3]					# A
 [4 5 6]]

[[1.1 2.  3. ]				# B
 [4.  5.  6. ]]

[[1.+0.j 2.+0.j 3.+0.j]		# C
 [4.+0.j 5.+0.j 6.+0.j]]

Notice how the matrices are printed in an easily readable format without the use of for-loops. This is another benefit of using NumPy arrays.

There are several other functions in NumPy that can create specific matrices. Here are some of them:

  • numpy.zeros returns a new array with given shape and type, filled with zeros.

    zeros
    A = np.zeros((2, 3))
    
    # When printed:
    
    [[0. 0. 0.]
     [0. 0. 0.]]
  • numpy.ones returns a new array with given shape and type, filled with ones.

    ones
    A = np.ones((2, 3))
    
    # When printed:
    
    [[1. 1. 1.]
     [1. 1. 1.]]
  • numpy.empty returns a new array with given shape and type, without initalizing entries.

    empty
    A = np.empty((2, 3))
    
    # When printed:
    
    [[ 1.33360294e+241  1.00200641e-310  3.42196482e-210]		# Random numbers
     [ 2.42869050e-313  1.59166660e-308 -1.95565444e-310]]
  • numpy.full returns a new array with given shape and type, filled with fill_value.

    full
    A = np.full((2, 3), 2019)
    
    B = np.full((2, 3), np.inf)
    
    # When printed:
    
    [[2019 2019 2019]		# A
     [2019 2019 2019]]	
    
    [[inf inf inf]			# B
     [inf inf inf]]

Linear algebra using NumPy

NumPy has several useful functions for built-in linear algebra. For a full list, check out the NumPy documentation on linear algebra. Optionally, and as mentioned at the top of this page, see SciPys documentation on linear algebra.

Dot product/Matrix multiplication

To compute the dot product of two arrays, A and B, we can use the function numpy.dotHowever, as recommended by NumPy, if both arrays are 2-D arrays it is preffered to use numpy.matmul or "@ B". As an example, we will now compute the dot product of two matrices, A and B, in all three ways. Let 

A = np.array([[1, 2], [3, 4]])
B = np.array([[5, 6], [7, 8]])

# Using "np.dot"
print(np.dot(A, B))

# Using "np.matmul"
print(np.matmul(A, B))

# Using "@"
print(A @ B)


# All three version will produce the same output:
[[19 22]
 [43 50]]

The @-operator and numpy.matmul uses the the same computing approach, where as numpy.dot is different. For more detailed information about the differences, read the function documentations provided by NumPy.

Solving a system of linear equations

Let the set of equations we would like to solve be and . Let's first write it as matrices:

The function we will use to solve the equation is numpy.linalg.solveIt takes two parameters, a and b, where both are arrays (remember that our matrix is an array type). numpy.linalg.solve will then computes the “exact” solution, x, of the well-determined, i.e., full rank, linear matrix equation ax = b .

a = np.array([[5, 2],[1,4]])    # Left-hand side
b = np.array([33, 21])          # Right-hand side

print(np.linalg.solve(a,b))     # Printing the solution


# Output:
[5. 4.]		# In our example this evaluates to x = 5, y = 4, which is the correct answer!


# You may also save the values directly into variables:
x, y = np.linalg.solve(a,b)

Eigenvalues and eigenvectors

numpy.linalg.ein is a function that computes the eigenvalues and eigenvectors of a square array, returning two arrays consisting of the eigenvalues and eigenvectors respectivly. Lets look at an example. We want to find the eigenvalues and eigenvectors of matrix A, where.

A = np.array([[1, 3], [4, -3]])

w, v = np.linalg.eig(A)		# "w" are the eigenvalues, "v" the eigenvectors

print(w)
print(v)

 
# Output:
[ 3. -5.]						# w
[[ 0.83205029 -0.4472136 ]		# v
 [ 0.5547002   0.89442719]]

The eigenvectors are arranged so that column v[:,i] is the eigenvector corresponding to the eigenvalue w[i], i.e. is the eigenvector to the eigenvalue 3, and is the eigenvector to the eigenvalue -5.

Note that the function returns the eigenvectors with a unit length. An easier way to present the eigenvectors above is as and , respectively.

Determinant

To compute the determinant of a square array, we will use the function numpy.linalg.det. It takes an array (or several arrays combined) as input, and returns the determinant. As an example, we will compute the determinant of matrix A, where.

A = np.array([[1, 2], [3, 4]])

print(np.linalg.det(A))


# Output:
-2.0000000000000004
Example: Finding the determinant of multiple arrays at once
A = np.array([[1, 2], [3, 4]])
B = np.array([[1, 2], [2, 1]])
C = np.array([[4, 4], [2, 3]])

 # To compute the determinant of multiple arrays at once, we have to put our arrays into a single array as input.
 # This will result in an array of length 3 as output, with the determinants respectivly.
print(np.linalg.det([A, B, C]))    


# Output:
[-2. -3.  4.] 		# Evaluates to: det(A) = -2, det(B) = -3 and det(C) = 4.

Inverting a matrix

To invert a matrix, we can use the function numpy.linalg.inv. It takes a square array as input, and returns the invers of the given array. As an example, we will compute the invers of matrix A, where 

A = np.array([[1, 2], [3, 4]])

print(np.linalg.inv(A))


# Output:
[[-2.   1. ]
 [ 1.5 -0.5]]
Example: Check the answer
# To check if the invers we got in the example above is correct, we simply compute the dot product of A and its invers.
A = np.array([[1, 2], [3, 4]])
B = np.array([[-2., 1. ], [1.5, -0.5]])		# This is the presumed invers of A.

print(A @ B)	# Matrix product/ dot product


# Output:
[[1. 0.]		# This is the identity matrix, which implies that B really is the invers of A.
 [0. 1.]]		# Remember how the invers works, where "B @ A" would have given us the same answer.

Exercises and solutions

If you want have a look at some exercises and solutions regarding NumPy, SciPy, matrices and linear algebra, visit Exercises and solutions, linear algebra in Python.

The exercises provided will have a varying degree of difficulty, thus including more advanced and detailed methods than described here.