Python Matrix Module

Matrices are an important part of linear algebra. By arranging scalars in a rectangular manner, one can elegantly encode vector transformations like scaling, rotating, shearing and squashing, solve systems of linear equations and represent general vector space homomorphisms.
However, as powerful as matrices are, when actually applying the theoretical tools one has to calculate specific values. Doing so by hand can be done, yet gets cumbersome quite quickly when dealing with any matrices which contain more than a few rows and columns.

So, even though there are a lot of other implementations already present, I set out to write a Python matrix module containing a matrix class capable of basic matrix arithmetic (matrix multiplication, transposition, …) together with a set of functions which perform some higher-level matrix manipulation like applying Gaussian elimination, calculating the reduced row echelon form, determinant, inversion and rank of a given matrix.

Module source code can be seen below and downloaded. When saved as matrix.py in the current working directory, one can import the module as follows.

>>> import matrix
>>> A = matrix.Matrix([[13,  1, 20, 18],
...                    [ 9, 24,  0,  9],
...                    [14, 22,  5, 18],
...                    [19,  9, 15, 14]])
>>> print A**-1
        -149/1268  -67/634   83/1268 171/1268
          51/1268 239/1902 -105/1268 -33/1268
Matrix(    73/634 803/4755  -113/634 -87/3170 )
          13/1268  -75/634  197/1268 -83/1268

Matrices are defined over a field, typically \mathbb{F} = \mathbb{R} in theoretical use, though for my implementation I chose not to use a double data structure, as it lacked the conceptual precision in numbers like a third. As one cannot truly represent a large portion of the reals anyways, I chose to use \mathbb{F} = \mathbb{Q}, which also is a field though can be — to a certain scalar size and precision — accurately represented using fractional data types (Python’s built-in Fraction is used here).

To simplify working with matrices, the implemented matrix class supports operator overloading such that the following expressions — A[i,j], A[i,j]=l, A*B, A*l, A+B, -A, A/l, A+B, A-B, A**-1, A**"t", ~A, A==B, A!=B — all behave in a coherent and intuitive way for matrices A, B, scalars l and indices i, j.

When working with matrices, there are certain rules that must be obeyed, like proper size when adding or multiplying, invertibility when inverting and others. To minimize potential bug track down problems, I tried to include a variety of detailed exceptions (ValueErrors) explaining the program’s failure at that point.

Apart from basic matrix arithmetic, a large part of my module centers around Gaussian elimination and the functions that follow from it. At their heart lies the implementation of GaussianElimination, a function which calculates the reduced row echelon form rref(A) of a matrix together with the transformation matrix T such that T*A = rref(A), a list of all matrix pivot coordinates, the number of row transpositions used to achieve row echelon form and a product of all scalars used to achieve reduced row echelon form.
From this function, rref(A) simply returns the first, rrefT(A) the second parameter. Functions rcef(A) (reduced column echelon form) and rcefS(A) (A*S=rcef(A)) follow from repeated transposition.
Determinant calculation uses characteristic determinant properties (multilinear, alternating and the unit hypercube has hypervolume one).

  • \det \begin{pmatrix}  a_{1 1}&a_{1 2}&\dots&a_{1 n}\\  \vdots&\vdots&\ddots&\vdots&\\  \lambda \cdot a_{i 1}&\lambda \cdot a_{i 2}&\dots&\lambda \cdot a_{i n}\\  \vdots&\vdots&\ddots&\vdots&\\  a_{n 1}&a_{n 2}&\dots&a_{n n}  \end{pmatrix} = \lambda \cdot \det \begin{pmatrix}  a_{1 1}&a_{1 2}&\dots&a_{1 n}\\  \vdots&\vdots&\ddots&\vdots&\\  a_{i 1}&a_{i 2}&\dots&a_{i n}\\  \vdots&\vdots&\ddots&\vdots&\\  a_{n 1}&a_{n 2}&\dots&a_{n n}  \end{pmatrix}
  • \det \begin{pmatrix}  a_{1 1}&a_{1 2}&\dots&a_{1 n}\\  \vdots&\vdots&\ddots&\vdots&\\  a_{i 1}&a_{i 2}&\dots&a_{i n}\\  \vdots&\vdots&\ddots&\vdots&\\  a_{j 1}&a_{j 2}&\dots&a_{j n}\\  \vdots&\vdots&\ddots&\vdots&\\  a_{n 1}&a_{n 2}&\dots&a_{n n}  \end{pmatrix} = -\det \begin{pmatrix}  a_{1 1}&a_{1 2}&\dots&a_{1 n}\\  \vdots&\vdots&\ddots&\vdots&\\  a_{j 1}&a_{j 2}&\dots&a_{j n}\\  \vdots&\vdots&\ddots&\vdots&\\  a_{i 1}&a_{i 2}&\dots&a_{i n}\\  \vdots&\vdots&\ddots&\vdots&\\  a_{n 1}&a_{n 2}&\dots&a_{n n}  \end{pmatrix}
  • \det \begin{pmatrix}1&a_{1 2}&a_{1 3}&\dots&a_{1 n}\\0&1&a_{2 3}&\dots&a_{2 n}\\0&0&1&\dots&a_{3 n}\\\vdots&\vdots&\vdots&\ddots&\vdots\\0&0&0&\dots&1\end{pmatrix} = \det \bold{1}_n = 1

Using these properties, the determinant is equal to the product of the total product of all factors used in transforming the matrix into reduced row echelon form and the permutation parity (minus one to the power of the number of transpositions used).
Questions regarding invertibility and rank can also conveniently be implemented using the above described information.

All in all, this Python module implements basic matrix functionality paired with a bit more involved matrix transformations to build a usable environment for manipulating matrices in Python.


# Python 2.7 code; Jonathan Frech; 9th, 10th, 11th, 12th, 14th of December 2017
# 4th of January 2018: added support for matrix inequality, implemented __ne__

Continue reading

Advertisements

Generic C / Python Polyglot

A polyglot — coming from the Greek words πολύς (many) and γλώττα (tongue) — is a piece of source code which can run in multiple languages, often performing language-dependent tasks.
Even though such a source code’s feature may not seem particularily useful or possible with certain language combinations, trying to bend not one but multiple language’s syntactic rules and overall behavior is still interesting.

An example of a rather simple polyglot would be a Python 2 / Python 3 polyglot — if one counts those as two separate languages. Because of their enormous similarities, one can pick out differences and use those to determine which language runs the source code.

if 1/2>0:print("Python 3")
else:print("Python 2")

Utilizing Python’s version difference regarding integer division and real division is one of the most common ways to tell them apart, as it can also be used to only control a specific part of a program instead of having to write two nearly identical programs (increases the program’s style and cuts on bytes — an important consideration if one golfs a polyglot).

However, polyglots combining languages that are not as similar as Python 2 and Python 3 require more effort. The following is a general Python 2 / C polyglot, meaning that nearly all C and Python 2 programs can be mixed using this template (there are a few rules both languages need to abide which will come apparent later).

#define _\
"""
main(){printf("C");}
#define _"""
#define/*
print"Python 2"
#*/_f

In the above code, main(){printf("C");} can be nearly any C code and print"Python 2" can be nearly any Python 2 code.
Language determination is exclusively done via syntax. A C compiler sees a #define statement and a line continuation \, another two #define statements with a block comment in between and actual compilable C source code (view first emphasis).
Python, on the other hand, treats all octothorps, #, as comments, ignoring the line continuation, and triple-quoted strings, """...""", as strings rather than statements and thus only sees the Python code (view second emphasis).

My first ever polyglot used this exact syntactical language differentiation and solved a task titled Life and Death of Trees (link to my answer).

Brainfuck X

While browsing StackExchange PCG questions and answers, I came across a challenge regarding drawing the swiss flag. In particular, I was interested in benzene’s answer, in which they showcased a Brainfuck dialect capable of creating two-dimensional 24-bit color images. In this post I present this dialect with slight changes of my own, as well as an interpreter I wrote in Python 2.7 (source code is listed below and can also be downloaded).

Brainfuck X generated Swiss flagUrban Müller’s original Brainfuck (my vanilla Brainfuck post can be found here) works similar to a Turing machine, in that the memory consists of a theoretically infinitely large tape with individual cells which can be modified. What allows Brainfuck X (or Braindraw, as benzene called their dialect) to create color images is, that instead of a one-dimensional tape, a three-dimensional tape is used. This tape extends infinitely in two spacial dimensions and has three color planes. Each cell’s value is limited to a byte (an integer value from 0 to 255) which results in a 24-bit color depth.

Adding to Brainfucks eight commands (+-<>[].,), there are two characters to move up and down the tape (^v) and one character to move forwards in the color dimension (*). Starting on the red color plane, continuing with the green and ending in the blue. After the blue color plane, the color planes cycle and the red color plane is selected. benzene’s original language design which I altered slightly had three characters (rgb) to directly select a color plane. Whilst this version is supported by my interpreter (the flag --colorletters is necessary for that functionality), I find my color star more Brainfucky — directly calling color planes by their name seems nearly readable.
Brainfuck’s vanilla eight characters still work in the same way, Brainfuck X can thereby execute any vanilla Brainfuck program. Also, there still is a plaintext output — the tape’s image is a program’s secondary output.

Having executed the final Brainfuck instruction, the interpreter prints out the tape to the terminal — using ANSI escape codes. Because of this, the color depth is truncated in the terminal view, as there are only 216 colors supported.
For the full 24-bit color depth output, I use the highly inefficient Portable Pixmap Format (.ppm) as an output image file format. To open .ppm files, I recommend using the GNU Image Manipulation Program; specifying the output file name is done via the --output flag.

The Swiss flag image above was generated by benzene’s Braindraw code (see their StackExchange answer linked to above); the resulting .ppm file was then scaled and converted using GIMP.
Interpreter command: python brainfuckx.py swiss.bfx -l -o swiss.ppm

Usage

  • Being written in pure Python, the interpreter is completely controlled via the command line. The basic usage is python brainfuckx.py <source code file>; by using certain flags the functionality can be altered.
  • --input <input string>-i <input string> specifies Brainfuck’s input and is given as a byte stream (string).
  • --simplify, -s outputs the source code’s simplified version; the source code with all unnecessary characters removed.
  • --colorstar selects the color star color plane change model which is the default.
  • --colorletters, -l selects the color letter color plane change model.
  • --silent stops the interpreter from outputting warnings, infos and the final tape.
  • --maxcycles <cycles>, -m <cycles> defines the maximum number of cycles the Brainfuck program can run; the default is one million.
  • --watch, -w allows the user to watch the program’s execution.
  • --watchdelay <delay> defines the time in seconds the interpreter sleeps between each watch frame.
  • --watchskip <N> tells the interpreter to only show every Nth cycle of the execution.
  • --output <output file name>, -o <output file name> saves the final tape as a .ppm image file.

# Python 2.7 Code; Jonathan Frech, 24th, 25th of August 2017

Continue reading