My heart sank because I can remember a hot autumn day (or did it last a week?) decades ago, when my brother and I were sitting in a business math class at Merrimack College being led through this exact topic by a teacher-in-training. It ranked among the most boring episodes in my educational experience. And pointless, too. Why were we going through all this? Because we had matrices to solve, I guess.
More recently I've posted about the more interesting, visual, interactive applications of matrices in 3D graphics:
So I worked through a few problems in Gaussian elimination to refresh my memory and predictably enough, I thought about automating this boring stuff with Python. Sure, you can just use Numpy's Linear Algebra module to solve them. Let's solve this example (from Wikipedia):
2x + y - z = 8
-3x - y + 2z = -11
-2x + y + 2z = -3
You have to put the coefficients and constants into matrices like this and type a few words to get the solution:
That means x = 2, y = 3 and z = -1. Check it by plugging those values into the original equations.
But maybe using modern-day solvers is cheating. This is math class, after all. Best to learn to make a solver yourself. I started by dividing the whole first row by the element in the first column. To use a list's values and their place in the list (their index), Python has a handy function called "enumerate." For example, declare a list like this and print the index and value of each item:
Using enumerate, it was easy to divide all the values in the first row (index 0) by the value of the top left entry.
#get top left entry to be 1
#by dividing row by A[0][0]
divisor = A[0][0]
for i,v in enumerate(A[0]):
A[0][i] = v / divisor
The next step was making the other elements in that column 0 by adding the negative of those values. Like in our example, the first item in the second row is -3. So to everything in that row we add 3 times the corresponding entry in the first row. Do the same thing for the third row. So having the index and the value of each entry in the matrix is crucial.
#add row 2 and 3 to additive inverse to make value 0
for i in [1,2]:
addinv = -1*A[i][0]
for ind,value in enumerate(A[i]):
A[i][ind] = value + addinv*A[0][ind]
Then you do the same thing to the second row, second column and make the other entries 0, and the same thing for the third row, third column, etc. I could smell a loop. You just need to have a list of rows and take out the one that's on the diagonal:
rowList = [x for x in range(len(A))]
rowList.remove(j)
I also needed to make sure we weren't dividing by zero. The whole function ended up being less than 20 lines:
But does it work? Let's use our Wikipedia example. I entered the matrix in augmented form, meaning it was a 3x4 matrix:
My function returned the entire matrix, but if you look at the last values in each row, those are the solutions for x, y and z.
The great thing is this function even works for bigger matrices. If you had a system of 4 unknowns, you'd need 4 equations like this:
4w + x + 2y + 3z = 6
3w - 6x + 5y + 4z = 1
9w + 7x - 7y + 8z = 4
w + x + y - z = 1
That would be cruel to make a student do by hand, but using our function, the worst part is just typing it all in:
The last numbers in each row are the solutions for w, x, y and z. Checked it with numpy:
>>>
[[-1.01875]
[ 1.8375 ]
[ 1.75625]
[ 1.575 ]]
One of the questions the students asked was to find the equation of a parabola given 3 points like this:
We know it's a parabola, so the solution is in the form
y = ax2 + bx + c
and we know 3 pairs of x and y coordinates, so it's a system of equations to solve for a, b and c:
(2) = a(-2)2 + b(-2) + c
(-1) = a(1)2 + b(1) + c
(-6) = a(2)2 + b(2) + c
or
4a - 2b + c = 2
a + b + c = -1
4a + 2b + c = -6
No problem solving this with gauss():
The solutions are a = -1, b = -2 and c = 2. So the equation of the parabola is
y = -x2 - 2x + 2
There are plenty of computer-based and online tools to solve linear systems using matrices. We should be teaching our math students how to make the tools themselves.
No comments:
Post a Comment