my way or the Conway

Another day, another implementation. Today I have implemented Conway’s game of life. It is a very simple model that can show very complicated results. It is also a lot of fun. What are the rule to this game? Well, I will tell you. First, let us start with the components of the system, we have

  • a 2-dimensional grid
  • a set of rules
    • alive and don’t have two or three alive neighbours -> dead
    • dead and three alive neighbours -> alive

Now we only have to define neighbours. In a 2D grid these are the 8 squares around a particular square. I have drawn a little ASCII graphic.

|---+---+---|
| O | O | O |
|---+---+---|
| O | X | O |
|---+---+---|
| O | O | O |
|---+---+---|

The O’s mark the neighbours and the X the particular square in question. Now all we have to do is implement it.

it’s complicated

I have been busy with complex numbers in JavaScript. I wanted to implement the Durand Kerner root finding algorithm. This algorithm requires complex arithmetic in order to work. So, in order to implement something, you need to implement something else. You just have to hope that you have arrived at the lowest level.

What is a complex number? It is a number of the type: \[ z = a + bi \] this is called rectangular form, and it is the form most people are familiar with. There is another form, called polar form, which looks like this: \[ z = re^{i\theta} \]

To translate from polar form to rectangular form we use Euler’s formula: \[ re^{i\theta} = r(\cos \theta + i \sin \theta). \] To translate from rectangular to polar coordinates, we use the special function \(\arctan\): \[ a + bi = \sqrt{a^2 + b^2}e^{i \arctan(b, a)}. \]

Why do we need these two forms? Well, some of the arithmetic operations make more sense in one form than in the other and, moreover, are easier to compute. To implement most of the functions we are familiar with we need to start with the three basic operators:

  • addition
  • multiplication
  • exponentiation

My implementation assumes the rectangular form as canonical. All operators are written for this form. For exponentiation we translate to polar, perform the operations, and translate back to rectangular.

Now for the mathematical definitions: \[ (a + bi) + (c + di) = (a + c) + (b + d)i \] \[ (a + bi) \cdot (c + di) = (a \cdot c – b \cdot d) + (a \cdot d + b \cdot c)i \] To define exponentiation we will first define the logarithm and the exponent. Notice how these operations are very similar to the translation operations defined before. \[ \log(a + bi) = \log(\sqrt{a^2 + b^2}e^{i \arctan(b, a)}) = \frac{1}{2}\log(a^2 + b^2) + \arctan(b, a)i \] \[ e^{a + bi} = e^a e^{bi} = e^a(\cos b + i\sin b) \] Now we can establish the general complex power: \[ (a + bi) ^{c + di} = e^{(c + di)\log(a + bi)} \] Which is defined in terms of our previous operations \(\cdot,\log,\exp\).

And that is most of it. To approximate trigonometric functions we use the identity (also found by using Euler’s formula): \[ \cos z = \frac{e^{iz} + e^{-iz}}{2} \] Calculating the inverse is left as an exercise to the reader. My code can be found here.

brainfark

I wrote my own little brainfuck interpreter. You can read all about brainfuck on the esoteric programming website. They have a whole list of weird programming languages. So what is brainfuck? They tell you all about it on the Esolang website. But I’m gonna tell you anyway.

Brainfuck is a programming language designed to be as small as possible. This does not mean that programs written in brainfuck are neccesarily small. No, it means that the interpreter/compiler can be written in a small number of bytes.

There are only 8 commands in the language:

  • > move pointer right
  • < move pointer left
  • + increase value at pointer
  • - decrease value at pointer
  • . output character at pointer
  • , input character at pointer
  • [ jump past ] if value at pointer is 0
  • ] jump to [ if value at pointer is not 0.

This language is Turing complete, which means that any program may be written in it. This is not recommended however.

I implemented my interpreter in the C programming language, because of its native pointer support. Let’s look at the code:

void eval(char *s, char *p) {
  int i;

  for (; *s; s++)
         if (*s == '>') ++p;
    else if (*s == '<') --p;
    else if (*s == '+') ++*p;
    else if (*s == '-') --*p;
    else if (*s == '.') printf("%c", *p);
    else if (*s == ',') scanf("%c", p);
    else if (*s == '[' && !*p) {
      for (i = 1, s++; i; s++)
        i += (*s == '[') - (*s == ']');
      s--;
    }
    else if (*s == ']' && *p)
      for (i = 1, s--; i; s--)
        i += (*s == ']') - (*s == '[');
}

Doesn’t that look nice and neat? I was very pleased with myself. Until I found out that a lot of people have written brainfuck interpreters in brainfuck itself, that are shorter and of course more elegant. The *s pointer points to the program text. The *p pointer points to the actual data array.

The structure of the program should be familiar if you have ever written any semantic sucking code. We try to match one of the 8 commands and execute logic accordingly. The s-- at the end of the ] is there to ensure the string pointer does not go to far past the matching ].

To actually make this code do something, some additional I/O code needs to be written. Try making an interpreter for yourself and see what happens!

Double trouble

Following up on the post about Lorenz systems I have made a visualization of the double pendulum. Both of these systems are prime examples of systems studied in chaos theory. That is: systems that are highly sensitive to initial conditions.

Many of these systems can also be seen in nature. For example the weather. One of the problems is that numbers in nature are non computable. They simpy can not be estimated to within arbitrary precision by any program. This is kind of a scary thought.

Back to the system at hand. I used the wiki, the explanation there is probably better. But here is my take anyway. The double pendulum can be described by 7 variables:

  • gravity constant
  • length both rods
  • mass both rods
  • angle first rod
  • angle second rod
  • momentum first rod
  • momentum second rod.

The total energy (Langrangian) in the system consists of:

  • linear kinetic energy \(E_k = \frac{1}{2} m v^2\)
  • rotational kinetic energy \(E_r = \frac{1}{2} I \omega^2\)
  • potential energy \(E_p = mgh\)

We use some other identities to fill in these formula’s:

  • moment of inertia thin rod \(I = \frac{1}{12} m l^2\)
  • x-center of mass rod \(x = \frac{l}{2} \sin \theta\)
  • y-center of mass rod \(y = \frac{l}{2} \cos \theta\)
  • 2d velocity \(v = \dot{x} + \dot{y}\)

We throw all this together to get our pendulum applet. Have a happy new year!

Lorenz system

I have been doing some stuff with differential equations. They provide some unexpected results and show random behavior. I made a little visualization with the html canvas of a Lorenz (not Lorentz) systemhttps://en.wikipedia.org/wiki/Lorenz_system

<canvas>
<script>
  c = document.getElementsByTagName("canvas")[0];
  (ctx = c.getContext("2d")), (c.width = 500), (c.height = 500);
  (x = y = z = 1), (a = 10), (b = 28), (c = 8 / 3), (dt = 0.005);
  setInterval(() => {
      [x, y, z] = [
	  x + a * (y - x) * dt,
	  y + (x * (b - z) - y) * dt,
	  z + (x * y - c * z) * dt
      ];
      ctx.fillRect(200 + 5 * x, 200 + 5 * y, 1, 1);
  });
</script>

You can see it live here:https://05dd8515-5617-4479-92db-df9de14327bc.htmlpasta.com/

My mandelbrot

I have tried implementing a zoomable mandelbrot set for a while. Today I created a satisfactory implementation. There is no native support for complex numbers in JavaScript, so we keep both real and imaginary parts of the complex number separate. The mandelbrot recurrence relation is:

z_{n+1} = z_{n} + c

The value of z_ 0 is c. We iterate over these values. The mandelbrot set is the set of all complex numbers for which this recurrence relation converges. We will write a function for this relation and then we are basically done. The rest of the code will simply draw the set and handle zoom events.

function f(a0, b0) {
      [a, b] = [a0, b0]
      for (i = 0; i < N && a*a+b*b<4; i++)
	  [a, b] = [a*a - b*b + a0, 2*a*b + b0]
      return i === N
  }

This is our function. The implementation can be seen live at https://e835fe86-83d5-4143-83c7-6e3fef8c6704.htmlpasta.com/. Use the spacebar to zoom and the arrow keys to move.

Deep Learning and CNN’s, which optimizer to use.

When designing CNN’s, or any neural network for that matter there is always a debate over what you want to use, and if methods for improving neural networks actually work. If methods actually work and if they work for the reason stated or a different reason is always a difficult thing to say. In the following few posts I will go over commonly used methods and what their effect is on the Neural networks. To start off this series I will begin with optimizers in a CNN (convolutional neural network). Using a control network I will go over the different optimizers that exist and go over the effects they have. In order to keep things fair we will use the same network but only change the optimizer.

Continue reading “Deep Learning and CNN’s, which optimizer to use.”

Random number cycles

Random number generation is often done by applying a function recursively to a seed number. The function needs to be bounded to a certain domain. Obtaining a function that is bounded is easily done by applying the modulo operation on the base function. The function must be deterministic. For example the following function

Continue reading “Random number cycles”

Fibonacci for beginners

The Fibonacci sequence is one of the most famous sequences. It appears everywhere in nature and mathematics. It is defined as:

\[ F(1) = 1 \] \[ F(2) = 1 \] \[ F(n) = F(n – 1) + F(n – 2) \]

We don’t want to calculate all these numbers by hand. For n = 100 the number has more than twenty digits. Addition is easy for humans and even easier for computers. Let’s write it with some Python, because everybody loves Python: Continue reading “Fibonacci for beginners”

Ulam’s primes

Searching for prime numbers has been a topic that many mathematicians have dedicated their time to over the years. Being able to generate prime numbers, and test if a number is prime quickly is an important problem in the field of computer science, as this would lead to major problems in encryption.

One of these mathematicians that enjoyed spending their time finding primes was the polish physicist/mathematician Stanislaw Ulam. Although he may be more commonly known for his works on the Monte Carlo method, or the Manhattan project. He also came up with an interesting idea concerning primes.

Continue reading “Ulam’s primes”