Note that there are some explanatory texts on larger screens.

plurals
  1. PO
    primarykey
    data
    text
    <p><strong>Original answer:</strong> I'm not sure if you will like how mathematical courses typically introduce matrices. As a programmer you might be happier with grabbing any decent 3D graphics book. It should certainly have very concrete 3x3 matrices. Also find out the ones that will teach you <a href="http://en.wikipedia.org/wiki/Projective_transformations" rel="nofollow noreferrer">projective transformations</a> (projective geometry is a very beautiful area of low-dimensional geometry and easy to program).</p> <h1>Mini-course in matrix math with Python 3</h1> <p><strong>Contents:</strong></p> <ol> <li><strong>Matrices</strong> <code>[Vector, __add__, reflect_y, rotate, dilate, transform]</code></li> <li><strong>Matrices: Overloaded</strong> <code>[Matrix, __add__, __str__, __mul__, zero, det, inv, __pow__]</code></li> <li><strong>Bonus: Complex numbers</strong></li> <li><strong>Matrices: The (R)evolution</strong>. It's already in the making (there's a summary at the end)</li> </ol> <p><strong>Preface:</strong> Based on my teaching experience, I think that the courses referenced by others are very good <em>courses</em>. That means if your goal is understanding matrices as mathematicians do, than you should by all means get the whole course. But if your goals are more modest, here's my try at something more tailored to your needs (but still written with the goal to convey many theoretical concepts, kind of contradicting my original advice.)</p> <p><strong>How to use:</strong></p> <ul> <li>This post is long. You might consider printing this and going slow, like one part a day.</li> <li>Code is essential. This is a course for programmers. Exercises are essential too. </li> <li>You should <strong><a href="http://mit.edu/~unknot/www/matrices.py" rel="nofollow noreferrer">take a look at the code companion</a></strong> which contains all this code and much more</li> <li>It's "2 for the price of 1" special: you can also learn Python 3 here. And complex numbers.</li> <li>I'll highly value any attempt to read this (do I officially qualify for the longest post ever?), so feel free to comment if you don't understand something (and also if you do).</li> </ul> <h1>1. Matrices</h1> <h2>Vectors</h2> <p>Before matrices come vectors. You sure know how to handle the 2- and 3- dimensional vectors:</p> <pre><code>class Vector: """This will be a simple 2-dimensional vector. In case you never encountered Python before, this string is a comment I can put on the definition of the class or any function. It's just one of many cool features of Python, so learn it here! """ def __init__(self, x, y): self.x = x self.y = y </code></pre> <p>now you can write</p> <pre><code>v = Vector(5, 3) w = Vector(7, -1) </code></pre> <p>but it's not much fun by itself. Let's add more useful methods:</p> <pre><code> def __str__(self: 'vector') -&gt; 'readable form of vector': return '({0}, {1})'.format(self.x, self.y) def __add__(self:'vector', v: 'another vector') -&gt; 'their sum': return Vector(self.x + v.x, self.y + v.y) def __mul__(self:'vector', number: 'a real number') -&gt; 'vector': '''Multiplies the vector by a number''' return Vector(self.x * number, self.y * number) </code></pre> <p>That makes things more interesting as we can now write:</p> <pre><code>print(v + w * 2) </code></pre> <p>and get the answer <code>(19, 1)</code> nicely printed as a vector (if the examples look unfamiliar, think how this code would look in C++). </p> <h2>Transformations</h2> <p>Now it's all cool to be able to write <code>1274 * w</code> but you need more vector operations for the graphics. Here are some of them: you can flip the vector around <code>(0,0)</code> point, you can reflect it around <code>x</code> or <code>y</code> axis, you can rotate it clockwise or counterclockwise (it's a good idea to draw a picture here). </p> <p>Let's do some simple operations:</p> <pre><code> ... def flip(self:'vector') -&gt; 'vector flipped around 0': return Vector(-self.x, -self.y) def reflect_x(self:'vector') -&gt; 'vector reflected around x axis': return Vector(self.x, -self.y) print(v.flip(), v.reflect_x()) </code></pre> <ul> <li><strong>Question:</strong> is it possible to express <code>flip(...)</code> using the operations I had below? What about <code>reflect_x</code>?</li> </ul> <p>Now you may wonder why I omitted <code>reflect_y</code>. Well, it's because I want you to stop for a moment and write your own version of it. Ok, here's mine:</p> <pre><code> def reflect_y(self:'vector') -&gt; 'vector reflected around y axis': return self.flip().reflect_x() </code></pre> <p>See, if you look how this function computes, it's actually quite trivial. But suddenly an amazing thing happened: I was able to write a transformation using only the existing transformations <code>flip</code> and <code>reflect_x</code>. For all I care, <code>reflect_y</code> could be defined in a derived class without access to <code>x</code> and <code>y</code> and it would still work!</p> <p>Mathematicians would call these functions <em>operators</em>. They would say that <code>reflect_y</code> is an operator obtained by <em>composition</em> of operators <code>flip</code> and <code>reflect_x</code> which is denoted by <code>reflect_y = flip ○ reflect_x</code> (you should see the small circle, a Unicode symbol <code>25CB</code>).</p> <ul> <li><strong>Note:</strong> I will quite freely use the <code>=</code> symbol from now to denote that two operations produce the same result, like in the paragraph above. This is a "mathematical <code>=</code>", which <a href="http://en.wikipedia.org/wiki/List_of_undecidable_problems" rel="nofollow noreferrer">cannot be expressed as a program</a>.</li> </ul> <p>So if I do </p> <pre><code>print(v.reflect_y()) </code></pre> <p>I get the result <code>(-5, 3)</code>. Go and picture it!</p> <ul> <li><strong>Question:</strong> Consider a composition <code>reflect_y ◦ reflect_y</code>. How would you name it? </li> </ul> <h2>Rotations</h2> <p>Those operations were nice and useful, but you are probably wondering why am so slow to introduce rotations. Ok, here I go:</p> <pre><code> def rotate(self:'vector', angle:'rotation angle') -&gt; 'vector': ?????? </code></pre> <p>At this point, if you know how to rotate vectors, you should go on and fill in the question marks. Otherwise please bear with me for one more simple case: counterclockwise rotation by <code>90</code> degrees. This one is not hard to draw on a piece of paper:</p> <pre><code> def rotate_90(self:'vector') -&gt; 'rotated vector': new_x = - self.y new_y = self.x return Vector(new_x, new_y) </code></pre> <p>Trying </p> <pre><code>x_axis = Vector(1, 0) y_axis = Vector(0, 1) print(x_axis.rotate_90(), y_axis.rotate_90()) </code></pre> <p>now gives <code>(0, 1) (-1, 0)</code>. Run it yourself!</p> <ul> <li><strong>Question:</strong> Prove that <code>flip = rotate_90 ◦ rotate_90</code>.</li> </ul> <p>Anyway, I won't hide the secret ingredient for longer:</p> <pre><code>import math # we'll need math from now on ... class Vector: ... def rotate(self:'vector', angle:'rotation angle') -&gt; 'rotated vector': cos = math.cos(angle) sin = math.sin(angle) new_x = cos * self.x - sin * self.y new_y = sin * self.x + cos * self.y return Vector(new_x, new_y) </code></pre> <p>Now let's try something along the lines:</p> <pre><code>print(x_axis.rotate(90), y_axis.rotate(90)) </code></pre> <p>If you expect the same result as before, <code>(0, 1) (-1, 0)</code>, you're bound to be disappointed. That code prints:</p> <pre><code>(-0.448073616129, 0.893996663601) (-0.893996663601, -0.448073616129) </code></pre> <p>and boy, is it ugly!</p> <ul> <li><p><strong>Notation:</strong> I will say that we applied operation <code>rotate(90)</code> to <code>x</code> in the example above. The knowledge we gained is that <code>rotate(90) != rotate_90</code>.</p></li> <li><p><strong>Question:</strong> What happened here? How to express <code>rotate_90</code> in terms of <code>rotate</code>? How to express <code>flip</code> in terms of <code>rotate</code>?</p></li> </ul> <h2>Dilations</h2> <p>Those rotations are certainly useful, but they are not everything you need to do even the 2D graphics. Consider the following transformations: </p> <pre><code> def dilate(self:'vector', axe_x:'x dilation', axe_y:'y dilation'): '''Dilates a vector along the x and y axes''' new_x = axe_x * self.x new_y = axe_y * self.y return Vector(new_x, new_y) </code></pre> <p>This <code>dilate</code> thing dilates the <code>x</code> and <code>y</code> axes in a possibly different way. </p> <ul> <li><strong>Exercise:</strong> Fill in the question marks in <code>dilate(?, ?) = flip</code>, <code>dilate(?, ?) = reflect_x</code>.</li> </ul> <p>I will use this <code>dilate</code> function to demonstrate a thing mathematicians call <em>commutativity</em>: that is, for every value of parameters <code>a</code>, <code>b</code>, <code>c</code>, <code>d</code> you can be sure that </p> <pre><code>dilate(a, b) ◦ dilate(c, d) = dilate(c, d) ◦ dilate(a, b) </code></pre> <ul> <li><p><strong>Exercise:</strong> Prove it. Also, is it true that for all possible values of parameters those below would hold?</p> <ul> <li><code>rotate(a) ◦ rotate(b) = rotate(b) ◦ rotate(a)</code></li> <li><code>dilate(a, b) ◦ rotate(c) = rotate(c) ◦ dilate(a, b)</code></li> <li><code>rotate(a) ◦ __mul__(b) = __mul__(b) ◦ rotate(a)</code></li> </ul></li> </ul> <h2>Matrices</h2> <p>Let's summarize all the stuff we had around here, our <em>operators on vector <code>x</code></em></p> <ul> <li><code>flip</code>, <code>reflect_x</code>, <code>*</code>, <code>rotate(angle)</code>, <code>dilate(x, y)</code></li> </ul> <p>from which one could make some really crazy stuff like </p> <ul> <li><code>flip ◦ rotate(angle) ◦ dilate(x, y) ◦ rotate(angle_2) ◦ reflect_y + reflect_x = ???</code></li> </ul> <p>As you create more and more complicated expressions, one would hope for some kind of order that would suddenly reduce all possible expressions to a useful form. Fear not! Magically, every expression of the form above can be simplified to</p> <pre><code> def ???(self:'vector', parameters): '''A magical representation of a crazy function''' new_x = ? * self.x + ? * self.y new_y = ? * self.x + ? * self.y return Vector(new_x, new_y) </code></pre> <p>with some numbers and/or parameters instead of <code>?</code>s.</p> <ul> <li><strong>Example:</strong> Work out what the values of '?' are for <code>__mul__(2) ◦ rotate(pi/4)</code></li> <li><strong>Another example:</strong> Same question for <code>dilate(x, y) ◦ rotate(pi/4)</code></li> </ul> <p>This allows us to write a universal function</p> <pre><code> def transform(self:'vector', m:'matrix') -&gt; 'new vector': new_x = m[0] * self.x + m[1] * self.y new_y = m[2] * self.x + m[3] * self.y return Vector(new_x, new_y) </code></pre> <p>which would take any 4-tuple of numbers, called <strong>matrix</strong>, and <em>apply</em> it to vector <code>x</code>. Here's an example:</p> <pre><code>rotation_90_matrix = (0, -1, 1, 0) print(v, v.rotate_90(), v.transform(rotation_90_matrix)) </code></pre> <p>which prints <code>(5, 3) (-3, 5) (-3, 5)</code>. Note that if you apply <code>transform</code> with any matrix to origin, you still get origin:</p> <pre><code>origin = Vector(0, 0) print(origin.transform(rotation_90_matrix)) </code></pre> <ul> <li><strong>Exercise:</strong> what are the tuples <code>m</code> that describe <code>flip</code>, <code>dilate(x, y)</code>, <code>rotate(angle)</code>?</li> </ul> <p>As we part with the <code>Vector</code> class, here's an exercise for those who want to test both their vector math knowledge and Pythonic skills:</p> <ul> <li><strong>The final battle:</strong> Add to the <code>Vector</code> class all vector operations that you can come up with (how many of standard operators can you overload for vectors? Check out my answer).</li> </ul> <h1>2. Matrices: Overloaded</h1> <p>As we found out in the previous section, a matrix can be thought of a shorthand that allows us to encode a vector operation in a simple way. For example, <code>rotation_90_matrix</code> encodes the rotation by 90 degrees.</p> <h2>Matrix objects</h2> <p>Now as we shift our attention from vectors to matrices, we should by all means have a class for matrix as well. Moreover, in that function <code>Vector.transform(...)</code> above the role of the matrix was somewhat misrepresented. It's more usual for <code>m</code> to be fixed while vector changes, so from now on our transformations will be methods of matrix class:</p> <pre><code>class Matrix: def __init__(self:'new matrix', m:'matrix data'): '''Create a new matrix. So far a matrix for us is just a 4-tuple, but the action will get hotter once The (R)evolution happens! ''' self.m = m def __call__(self:'matrix', v:'vector'): new_x = self.m[0] * v.x + self.m[1] * v.y new_y = self.m[2] * v.x + self.m[3] * v.y return Vector(new_x, new_y) </code></pre> <p>If you don't know Python, <code>__call__</code> overloads the meaning of <code>(...)</code> for matrices so I can use the standard notation for a matrix <em>acting</em> on a vector. Also, the matrices are usually written using a single uppercase letter:</p> <pre><code>J = Matrix(rotation_90_matrix) print(w, 'rotated is', J(w)) </code></pre> <ul> <li><strong>Exercise:</strong> repeat this example with matrices from the previous exercise. </li> </ul> <h2>Addition</h2> <p>Now, let's find out what else we can do with matrices. Remember that matrix <code>m</code> is really just a way to encode an operaton on vectors. Note that for two functions <code>m1(x)</code> and <code>m2(x)</code> I can create a new function (using <a href="http://en.wikipedia.org/wiki/Lambda_calculus" rel="nofollow noreferrer">lambda notation</a>) <code>m = lambda x: m1(x) + m2(x)</code>. It turns out if <code>m1</code> and <code>m2</code> were enconded by matrices, <em>you can also encode this <code>m</code> using matrices</em>!</p> <ul> <li><strong>Exercise:</strong> Think through any difficulties you might have with this statement.</li> </ul> <p>You just have to add its data, like <code>(0, 1, -1, 0) + (0, 1, -1, 0) = (0, 2, -2, 0)</code>. Here's how to add two tuples in Python, with some very useful and highly Pythonic techniques:</p> <pre><code> def __add__(self:'matrix', snd:'another matrix'): """This will add two matrix arguments. snd is a standard notation for second argument. (i for i in array) is Python's powerful list comprehension. zip(a, b) is used to iterate over two sequences together """ new_m = tuple(i + j for i, j in zip(self.m, snd.m)) return Matrix(new_m) </code></pre> <p>Now we can write expressions like <code>J + J</code> or even <code>J + J + J</code>, but to see the results we have to figure out how to print a Matrix. A possible way would be to print a 4-tuple of numbers, but let's take a hint from the <code>Matrix.__call__</code> function that the numbers should be organized into a <code>2x2</code> block:</p> <pre><code> def as_block(self:'matrix') -&gt; '2-line string': """Prints the matrix as a 2x2 block. This function is a simple one without any advanced formatting. Writing a better one is an exercise. """ return ('| {0} {1} |\n' .format(self.m[0], self.m[1]) + '| {0} {1} |\n' .format(self.m[2], self.m[3]) ) </code></pre> <p>If you look at this function in action you'll notice there is some room for improvement:</p> <pre><code>print((J + J + J).as_block()) </code></pre> <ul> <li><strong>Exercise:</strong> write a nicer function <code>Matrix.__str__</code> that will round the numbers and print them in the fields of fixed length.</li> </ul> <p>Now you should be able to write the matrix for rotation:</p> <pre><code>def R(a: 'angle') -&gt; 'matrix of rotation by a': cos = math.cos(a) sin = math.sin(a) m = ( ????? ) return Matrix(m) </code></pre> <ul> <li><p><strong>Exercise:</strong> Examine the code for <code>Vector.rotate(self, angle)</code> and fill in the question marks. Test with</p> <pre><code>from math import pi print(R(pi/4) + R(-pi/4)) </code></pre></li> </ul> <h2>Multiplication</h2> <p>The most important thing we can do with one-parameter functions is compose them: <code>f = lambda v: f1(f2(v))</code>. How to mirror that with matrices? This requires us to examine how <code>Matrix(m1) ( Matrix(m2) (v))</code> works. If you expand it, you'll notice that </p> <pre><code>m(v).x = m1[0] * (m2[0]*v.x + m2[1]*v.y) + m1[1] * (m2[2]*v.x + m2[3]*v.y) </code></pre> <p>and similarly for <code>m(v).y</code>, which, if you open the parentheses, looks suspiciously similar to <code>Matrix.__call__</code> using a new tuple <code>m</code>, such that <code>m[0] = m1[0] * m2[0] + m1[2] * m2[2]</code>. So let's take this as a hint for a new definiton:</p> <pre><code> def compose(self:'matrix', snd:'another matrix'): """Returns a matrix that corresponds to composition of operators""" new_m = (self.m[0] * snd.m[0] + self.m[1] * snd.m[2], self.m[0] * snd.m[1] + self.m[1] * snd.m[3], ???, ???) return Matrix(new_m) </code></pre> <ul> <li><p><strong>Exercise:</strong> Fill in the question marks here. Test it with</p> <pre><code>print(R(1).compose(R(2))) print(R(3)) </code></pre></li> <li><p><strong>Math exercise:</strong> Prove that <code>R(a).compose(R(b))</code> is always the same as <code>R(a + b)</code>.</p></li> </ul> <p>Now let me tell the truth: this <code>compose</code> function is actually how mathematicians decided to <em>multiply</em> matrices. This makes sense as a notation: <code>A * B</code> is a matrix that decribes operator <code>A ○ B</code>, and as we'll see next there are deeper reasons to call this 'multiplication' as well.</p> <p>To start using multiplication in Python all we have to do is to order it so in a <code>Matrix</code> class:</p> <pre><code> class Matrix: ... __mul__ = compose </code></pre> <ul> <li><strong>Exercise:</strong> Compute <code>(R(pi/2) + R(pi)) * (R(-pi/2) + R(pi))</code>. Try to find the answer yourself first on a piece of paper.</li> </ul> <h2>Rules for <code>+</code> and <code>*</code></h2> <p>Let's make some good name for the matrix that corresponds to the <code>dilate(a, b)</code> operator. Now there's nothing wrong with <code>D(a, b)</code>, but I'll use a chance to introduce a standard notation:</p> <pre><code>def diag(a: 'number', b: 'number') -&gt; 'diagonal 2x2 matrix': m = (a, 0, 0, b) return Matrix(m) </code></pre> <p>Try <code>print(diag(2, 12345))</code> to see why it's called a <em>diagonal</em> matrix. </p> <p>As the composition of operations was found before to be not always commutative, <code>*</code> operator won't be always commutative for matrices either. </p> <ul> <li><strong>Exercise:</strong> go back and refresh the <em>commutativity</em> thing if necessary. Then give examples of matrices <code>A</code>, <code>B</code>, made from <code>R</code> and <code>diag</code>, such that <code>A * B</code> is not equal to <code>B * A</code>.</li> </ul> <p>This is somewhat strange, since multiplication for numbers is always commutative, and raises the question whether <code>compose</code> really deserves to be called <code>__mul__</code>. Here's quite a lot of rules that <code>+</code> and <code>*</code> <em>do</em> satisfy:</p> <ol> <li><code>A + B = B + A</code></li> <li><code>A * (B + C) = A * B + A * C</code></li> <li><code>(A + B) * C = A * C + B * C</code></li> <li><code>(A * B) * C = A * (B * C)</code></li> <li><p>There is an operation called <code>A - B</code> and <code>(A - B) + B = A</code></p> <ul> <li><strong>Exercise:</strong> Prove these statements. How to define <code>A - B</code> in terms of <code>+</code>, <code>*</code> and <code>diag</code>? What does <code>A - A</code> equal to? Add the method <code>__sub__</code> to the class <code>Matrix</code>. What happens if you compute <code>R(2) - R(1)*R(1)</code>? What should it be equal to?</li> </ul></li> </ol> <p>The <code>(A * B) * C = A * (B * C)</code> equality is called <em>associativity</em> and is especially nice since it means that we don't have to worry about putting parentheses in an expression of the form <code>A * B * C</code>:</p> <pre><code>print(R(1) * (diag(2,3) * R(2))) print((R(1) * diag(2,3)) * R(2)) </code></pre> <p>Let's find analogues to regular numbers <code>0</code> and <code>1</code> and subtraction:</p> <pre><code>zero = diag(0, 0) one = diag(1, 1) </code></pre> <p>With the following easily verifiable additions:</p> <ol> <li><code>A + zero = A</code></li> <li><code>A * zero = zero</code></li> <li><code>A * one = one * A = A</code></li> </ol> <p>the rules become complete, in the sense that there is a short name for them: <em>ring axioms</em>. Mathematicians thus would say that matrices form a <em>ring</em>, and they indeed always use symbols <code>+</code> and <code>*</code> when talking about rings, and so shall we.</p> <p>Using the rules it's possible to easily compute the expression from the previous section:</p> <pre><code>(R(pi/2) + R(pi)) * (R(-pi/2) + R(pi)) = R(pi/2) * R(-pi/2) + ... = one + ... </code></pre> <ul> <li><strong>Exercise:</strong> Finish this. Prove that <code>(R(a) + R(b)) * (R(a) - R(b)) = R(2a) - R(2b)</code>. </li> </ul> <h2>Affine Transformations</h2> <p>Time to return to how we defined matrices: they are a shortcut to some operations you can do with vectors, so it's something you can actually draw. You might want to take a pen or look at the materials that others suggested to see examples of different plane transformations.</p> <p>Among the transformations we'll be looking for the <em>affine</em> ones, those who look 'the same' everywhere (no bending). For example, a rotation around some point <code>(x, y)</code> qualifies. Now this one cannot be expressed as <code>lambda v: A(v)</code>, but in can be written in the form <code>lambda v: A(v) + b</code> for some matrix <code>A</code> and vector <code>b</code>. </p> <ul> <li><strong>Exercise:</strong> find the <code>A</code> and <code>b</code> such that a rotation by <code>pi/2</code> around the point <code>(1, 0)</code> has the form above. Are they unique?</li> </ul> <p>Note that for every vector there is an affine transformation which is a <em>shift</em> by the vector. </p> <p>An affine transformation may stretch or dilate shapes, but it should do in the same way everywhere. Now I hope you believe that the area of any figure changes by a constant number under the transformation. For a transformation given by matrix <code>A</code> this coeffiecient is called the <em>determinant</em> of <code>A</code> and can be computed applying the formula for an area to two vectors <code>A(x_axis)</code> and <code>A(y_axis)</code>:</p> <pre><code> def det(self: 'matrix') -&gt; 'determinant of a matrix': return self.m[0]*self.m[3] - self.m[1] * self.m[2] </code></pre> <p>As a sanity check, <code>diag(a, b).det()</code> is equal to <code>a * b</code>.</p> <ul> <li><strong>Exercise:</strong> Check this. What happens when one of arguments is 0? When it's negative?</li> </ul> <p>As you can see, the determinant of rotation matrix is always the same:</p> <pre><code>from random import random r = R(random()) print (r, 'det =', r.det()) </code></pre> <p>One interesting thing about <code>det</code> is that it is multiplicative (it kind of follows from the definition if you meditate long enough):</p> <pre><code>A = Matrix((1, 2, -3, 0)) B = Matrix((4, 1, 1, 2)) print(A.det(), '*', B.det(), 'should be', (A * B).det()) </code></pre> <h2>Inverse</h2> <p>A useful thing you can do with matrices is write a system of two linear equations</p> <pre><code>A.m[0]*v.x + A.m[1]*v.y = b.x A.m[2]*v.x + A.m[3]*v.y = b.y </code></pre> <p>in a simpler way: <code>A(v) = b</code>. Let's solve the system as they teach in (some) high schools: multiply first equation by <code>A.m[3]</code>, second by -A.m<a href="http://en.wikipedia.org/wiki/Projective_transformations" rel="nofollow noreferrer">1</a> and add (if in doubt, do this on a piece of paper) to solve for <code>v.x</code>.</p> <p>If you really tried it, you should have got <code>A.det() * v.x = (A.m[3]) * b.x + (-A.m[1]) * b.y</code>, which suggests that you can always get <code>v</code> by multiplying <code>b</code> by some other matrix. This matrix is called <em>inverse</em> of <code>A</code>:</p> <pre><code> def inv(self: 'matrix') -&gt; 'inverse matrix': '''This function returns an inverse matrix when it exists, or raises ZeroDivisionError when it doesn't. ''' new_m = ( self.m[3] / self.det(), -self.m[1] / self.det(), ????? ) return Matrix(new_m) </code></pre> <p>As you see, this method fails loudly when determinant of matrix is zero. If you really want you can catch this expection with:</p> <pre><code>try: print(zero.inv()) except ZeroDivisionError as e: ... </code></pre> <ul> <li><strong>Exercise:</strong> Finish the method. Prove that inverse matrix doesn't exist when <code>self.det() == 0</code>. Write the method to divide matrices and test it. Use the inverse matrix to solve an equation <code>A(v) = x_axis</code> (<code>A</code> was defined above).</li> </ul> <h2>Powers</h2> <p>The main property of inverse matrix is that <code>A * A.inv()</code> always equals to <code>one</code></p> <ul> <li><strong>Exercise:</strong> check that yourself. Explain why that should be so from the definition of inverse matrix.</li> </ul> <p>That's why mathematicians denote <code>A.inv()</code> by <code>A</code><sup>-1</sup>. How about we write a nice function to use <code>A ** n</code> notation for <code>A</code><sup>n</sup>? Note that the naive <code>for i in range(n): answer *= self</code> cycle is O(|n|) which is certainly too slow, because this can be done with a complexity of <code>log |n|</code>:</p> <pre><code> def __pow__(self: 'matrix', n:'integer') -&gt; 'n-th power': '''This function returns n-th power of the matrix. It does it more efficiently than a simple for cycle. A while loop goes over all bits of n, multiplying answer by self ** (2 ** k) whenever it encounters a set bit. ... </code></pre> <ul> <li><p><strong>Exercise:</strong> Fill in the details in this function. Test it with</p> <p><code>X, Y = A ** 5, A ** -5</code> <code>print (X, Y, X * Y, sep = '\n')</code></p></li> </ul> <p>This function only works for integer values of <code>n</code>, even though for some matrices we can also define a fractional power, such as square root (in other words, a matrix <code>B</code> such that <code>B * B = A</code>).</p> <ul> <li><strong>Exercise:</strong> Find a square root of <code>diag(-1, -1)</code>. Is this the only possible answer? Find an example of matrix that <em>doesn't</em> have a square root.</li> </ul> <h1>Bonus: Complex numbers</h1> <p>Here I'm going to introduce you to the subject in exactly one section! Since it's a complex subject, I'm likely to fail, so please forgive me in advance.</p> <p>First, similarly to how we have matrices <code>zero</code> and <code>one</code>, we can make a matrix out of any real number by doing <code>diag(number, number)</code>. Matrices of that form can be added, subtracted, multiplied, inverted and the results would mimic what happens with the numbers themselves. So for all practical purposes, one can say that, e.g., <code>diag(5, 5)</code> <em>is</em> 5.</p> <p>However, Python doesn't know yet how to handle expressions of the form <code>A + 1</code> or <code>5 * B</code> where <code>A</code> and <code>B</code> are matrices. If you're interested, you should by all means go and do the following exercise or look at my implementation (which uses a cool Python feature called <em>decorator</em>); otherwise, just know that it's been implemented.</p> <ul> <li><strong>Exercise for gurus:</strong> Change the operators in a <code>Matrix</code> class so that in all standard operations where one of operands is a matrix and another a number, the number is automatically converted to the <code>diag</code> matrix. Also add comparison for equality.</li> </ul> <p>Here's an example test:</p> <pre><code>print( 3 * A - B / 2 + 5 ) </code></pre> <p>Now here's the first interesting <strong>complex number</strong>: the matrix <code>J</code>, introduced in the beginning and equal to <code>Matrix((0, 1, -1, 0))</code>, has a funny property that <code>J * J == -1</code> (try it!). That means <code>J</code> is certainly not a normal number, but, as I just said, matrices and numbers easily mix together. For example,</p> <pre><code>(1 + J) * (2 + J) == 2 + 2 * J + 1 * J + J * J = 1 + 3 * J </code></pre> <p>using the rules listed some time before. What happens if we test this in Python?</p> <pre><code>(1 + J) * (2 + J) == 1 + 3*J </code></pre> <p>That should happily say <code>True</code>. Another example:</p> <pre><code>(3 + 4*J) / (1 - 2*J) == -1 + 2*J </code></pre> <p>As you might have guessed, the mathematicians don't call those 'crazy numbers', but they do something similar - they call expressions of the form <code>a + b*J</code> <em>complex numbers</em>. Because those are still instances of our <code>Matrix</code> class, we can do quite a lot of operations with those: addition, subtraction, multiplication, division, power - it's all already implemented! Aren't matrices amazing?</p> <p>I have overlooked the question of how to print the result of operation like <code>E = (1 + 2*J) * (1 + 3*J)</code> so that it looks like an expression with <code>J</code> rather than a <code>2x2</code> matrix. If you examine it carefully, you'll see that you need to print the left column of that matrix in the format <code>... + ...J</code> (just one more nice thing: it's exactly <code>E(x_axis)</code>!) Those who know the difference between <code>str()</code> and <code>repr()</code> should see it's natural to name a function that would produce expression of such form as <code>repr()</code>.</p> <ul> <li><p><strong>Exercise:</strong> Write the function <code>Matrix.__repr__</code> that would do exactly that and try some tests with it, like <code>(1 + J) ** 3</code>, first computing the result on paper and then trying it with Python.</p></li> <li><p><strong>Math question:</strong> What is the determinant of <code>a + b*J</code>? If you know what the <em>absolute value</em> of complex number is: how they are connected? What is the absolute value of <code>a</code>? of <code>a*J</code>?</p></li> </ul> <h1>3. Matrices: The (R)evolution</h1> <p>In the final part of this trilogy we will see that everything is a matrix. We'll start with general <code>M x N</code> matrices, and find out how vectors can be thought of as <code>1 x N</code> matrices and why numbers are the same as diagonal matrices. As a side note we'll explore the complex numbers as <code>2 x 2</code> matrices.</p> <p>Finally, we will learn to write affine and projective transformations using matrices.</p> <p>So the classes planned are <code>[MNMatrix, NVector, Affine, Projective]</code>.</p> <p>I guess if you was able to bear with me until here, you could be interested in this sequel, so I'd like to hear if I should continue with this (and where, since I'm pretty much sure I'm beyond what considered reasonable length of a single document).</p>
    singulars
    1. This table or related slice is empty.
    plurals
    1. This table or related slice is empty.
    1. This table or related slice is empty.
    1. This table or related slice is empty.
    1. VO
      singulars
      1. This table or related slice is empty.
    2. VO
      singulars
      1. This table or related slice is empty.
    3. VO
      singulars
      1. This table or related slice is empty.
 

Querying!

 
Guidance

SQuiL has stopped working due to an internal error.

If you are curious you may find further information in the browser console, which is accessible through the devtools (F12).

Reload