For anybody viewing this on github, there is a rendering issue with some of the latex (one equation is missing and the column matrices are not displaying properly). It should display properly if you download the notebook and open it with jupyter.
Any other issues please contact j.blood14@imperial.ac.uk
For most day to day applications, or at GCSE, we can solve simultaneous equations by substitution. For most of us, when we only have a couple of values to find this method is not only adequate, but it is easier and less time consuming to carry out with pen and paper.
When we start to encounter simultaneous equations with multiple values to solve, substitution can be laborious and inelegant. Luckily, there is a quicker and more pleasing method to solve such equations: matrices!
For this example we will just use two unknown variables, $x$ and $y$, but this method can be extended to as many as you like.
In order to find these values using matrices we need to build... matrices. We will separate our coefficients into square matrix $A$, our $x$ and $y$ values into column matrix $X$, and our results into column matrix $B$. Multiplying matrix $A$ by matrix $X$ will give us matrix $B$. In general terms it looks like this: (If you want more information about multiplying matrices, please visit Wikipedia.)
$$AX=B$$Our matrices will look like this:
$$\begin{pmatrix}3 & -1\\5 & 2\\ \end{pmatrix} \begin{pmatrix}x\\y\\ \end{pmatrix}=\begin{pmatrix}3\\38\\ \end{pmatrix}$$Take a few minutes to make sure you can understand where these values are coming from.
Any number divided by itself is equal to 1. This is the same as the number multiplied by 1 divided by the number. For example:
$$7 \times \frac{1}{7} = 1 $$This is all well and good for integers, but what about matrices?
1 divided by a number is called it's inverse. And for a number, $n$, can be written as $n^{-1}$. Matrices also have this inversion property whereby a matrix multiplied by its inverse is equal to one. This will come in very useful for our general equation above. If we multiply both sides by the inverse of matrix A we get:
$$AXA^{-1} = A^{-1}B$$And as we now know that $AA^{-1} = 1$ we get:
$$ X = A^{-1}B $$i.e. Our values for x and y is the column matrix we get from multiplying B by our inverse of A.
To invert a $\begin{pmatrix} a&b\\c&d\\ \end{pmatrix}$ a and d must swap places, and the sign is changed on b and c. Each of these values is then divided by the determinant of the matrix ($ad - bc$). We are only very briefly going over these concepts of linear algebra here, and if you wish to know more about invertible matrices then there are a variety of books and online video lessons that can help you. For now, let us invert our matrix, A:
$$A^{-1}= \frac{1}{11} \begin{pmatrix}2&1\\-5&3\\ \end{pmatrix} $$Which is equal to:
$$A^{-1}=\begin{pmatrix}\frac{2}{11}&\frac{1}{11}\\ \frac{-5}{11}&\frac{3}{11}\\ \end{pmatrix}$$Now we have our inverted matrix, we just need to multiply it by B:
$$X = \begin{pmatrix}\frac{2}{11}&\frac{1}{11}\\ \frac{-5}{11}&\frac{3}{11}\\ \end{pmatrix} \begin{pmatrix}3\\ 38\\ \end{pmatrix} = \begin{pmatrix}\frac{2}{11}3 + \frac{1}{11}38\\ \frac{-5}{11}3+\frac{3}{11}38\\ \end{pmatrix} = \begin{pmatrix}3\\9\\\end{pmatrix}$$$$X=\begin{pmatrix}3\\9\\\end{pmatrix}$$This means: $$x=3$$ $$y=9$$
Thankfully we don't need to do all of this by hand, because it can start to take a long time. To carry out this task in python, we are using NumPy. NumPy is an extremely powerful library in python to work with n-dimensional arrays (the object type is known as ndarray). Using the NumPy functionality we should be able to write a function to solve simultaneous equations in much less than one second. In the example we are using a 2x2 array (you can think of array and matrix being interchangeable terms here).
I will walk you through the functions that we need and then let you write a general function. Try not to cheat, but I will upload a separate sim_eq.py script to github to show you how I did it.
First things first, we need to import numpy as is shown below.
In [ ]:
import numpy as np #import the NumPy library
Using the examples from above: $$ 3x - y = 3 $$ $$ 5x + 2y = 38 $$
Let's make our $A$ (2x2) and $B$ (1x2) arrays (we are trying to find $X$!)
In [ ]:
A=np.array([[3,-1],
[5,2]])
B=np.array([3,38])
To invert an array we use:
np.linalg.inv()
In [ ]:
A_inverted=np.linalg.inv(A)
A_inverted
NumPy will store the values as floats, whereas a human might write them as fractions (as I did above).
The next step is to multiply $A^{-1}$ by $B$ to get $X$, which we do in NumPy by using the dot (scalar) product.
np.dot(A_inverted,B)
In [ ]:
X=np.dot(A_inverted,B)
X
You might now by wondering why this looks like a row vector, when you were probably expecting a column vector. NumPy treats each of these as identical (the eagle eyed among you might notice that I used a row vector when defining B, try using a column vector and see what happens).
If needed we can now assign our $x$ and $y$ values in the normal way and our simultaneous equations are solved.
In [ ]:
x,y=X[0],X[1]
print('x =',int(x),'and y =',int(y))
Now try to write a general function that takes in a list of NumPy array (or lists, which you can convert to NumPy arrays if you're feeling adventurous) where the first $k$ number of elements in each list/array are the coefficients and the last element of each list is the result. This should return an array of values for the variables and it should be generalised so that you can solve for any number of unknown variables. Remember you need the same number of equations as variables, and that 0 is a valid coefficient. Good luck!