Notes on math, coding, and other stuff

Tag: c

Here’s something I’ve been working on for the past few weeks for one of my courses, CS488 – Intro to Computer Graphics. For the final project, you’re allowed to do any OpenGL or raytracing project, as long as it has 10 reasonable graphics related objectives. Here’s a video of mine:

A screenshot:

It’s a simple game where you control a boat and go around a lake collecting coins. When you collect a coin, there’s a bomb that spawns and follows you around. You die when you hit a bomb. Also if two bombs collide then they both explode (although you can’t see that in the video).

Everything is implemented in bare-metal OpenGL, so none of those modern game engines or physics engines. It’s around 1000-ish lines of C++ (difficult to count because there’s a lot of donated code).

Edit (8/10/2016) – I received an Honorable Mention for this project!

CS488 – Introduction to Computer Graphics

For those that haven’t heard about CS488, it’s one of the “big three” — fourth year CS courses with the heaviest workload and with large projects (the other two being Real-time and Compilers). It’s one of the hardest courses at Waterloo, but also probably the most rewarding and satisfying course I’ve taken.

There are four assignments, each walking you step by step through graphics techniques, like drawing a cube with OpenGL, or building a puppet with hierarchical modelling, or writing a simple ray tracer. Then there’s the final project, where you can choose to make something with OpenGL or extend your ray tracer. The class is split 50/50, about half the class did OpenGL and the other half did a ray tracer. I personally feel that OpenGL gives you more room to be creative and create something unique whereas ray tracing projects end up implementing a mix of different algorithms.

The first two assignments weren’t too bad (I estimate it took me about 10 hours each), but some time during assignment 3 I realized I was spending a lot of time in the lab, so I got an hours tracking app on my phone to track exactly how much time I was spending working on this course. Assignments 3 and 4 each took me 15 hours. I spent 35 hours on my final project, over a period of 3 weeks. I chose relatively easy objectives that I was confident I could do well, which left time to polish the game and do a few extra objectives. I’m not sure what the average is for time spent on the final project, but it’s common to spend 50-100 hours. Bottom line: you can put in potentially unbounded amounts of time to try to get the gold medal, but the effort actually required to get a good grade is quite reasonable.

Now the bad part about this course (obviously not the instructor’s fault) is OpenGL is so incredibly difficult to work with. Even to draw a line on the screen, you have to deal with a lot of low level concepts like vertex array objects, vertex buffer objects, uniform attributes to pass to shaders, stuff like that. It doesn’t help that when something goes wrong in a shader (which runs on the GPU), there’s no way to pass an error message back to the CPU so you can print out variables and debug it. It also doesn’t help that there’s a lot of incompatible OpenGL versions, and code you find in an online tutorial could be subtly broken for the version you’re using. On the other hand, working with OpenGL really makes you appreciate modern game engines like Unity which takes care of all the low level stuff for you.

Position is easy, just represent it with a point in 3D space. But how do you specify its orientation — which direction it’s pointing?

At first glance, it seems a vector will do. After all, a vector points in some direction, right? If the plane is pointing east, represent its orientation by a unit vector pointing east.

Unfortunately, we quickly run into trouble when we try to roll. If we’re facing east, and we roll 90 degrees, we’re still facing east. Clearly we’re missing something.

Euler Angles

When real pilots talk about their orientation, they talk about roll, yaw, pitch. Pitch is going up or down, yaw is going left or right, roll is, well, roll.

Any change in orientation can be described by some combination of roll, yaw, pitch. This is the basis for Euler Angles. We use three angles to represent the airplane’s orientation.

This is all fine and dandy if we want to represent the orientation of a static object in space. But when we try to adjust our orientation, we start to run into problems.

You’re thinking, this should be simple! When we turn left or right, we just increment the yaw variable, right? Yes, it seems to work, at least initially. You can turn left and right, up and down, and roll around.

Implement it in Unity and play around a bit, however, and you begin to notice that things don’t quite behave the way you expect.

In this animation, I’m holding down the right button:

The plane does rotate to the right, but it’s not rotating relative to itself. Instead it’s rotating around some invisible y-axis. If it was rotating relative to itself, the green arrow shouldn’t be moving.

The problem becomes more and more severe when the pitch of the plane becomes higher and higher. The worst case is when the airplane is pointing straight up: then roll and yaw become the same thing! This is called gimbal lock: we have lost a degree of freedom and we can only rotate in 2 dimensions! Definitely not something desirable if we’re controlling a plane or spaceship.

It turns out that no matter what we do, we will suffer from some form of gimbal lock. As long as we use Euler Angles, there is one direction where if we turn too far, everything starts to screw up.

Practical Introduction to Quaternions

All is not lost, however. There is a way to represent orientation that represents all axes equally and does not suffer from gimbal lock. This mythical structure is called the quaternion. Unlike Euler Angles which describe your orientation relative to a fixed set of axes, quaternions do not rely on any fixed axis.

The drawback is that quaternions are unintuitive to understand for humans. There is no way to “look” at a quaternion and be able to visualize what rotation it represents. Fortunately for us, it’s not that difficult to make use of quaternions, even if we can’t visualize quaternions.

There is a lot of theory behind how quaternions work, but in this article, I will gloss over the theory and give a quick primer to quaternions, just the most common facts you need to use them. At the same time, I will implement the operations I describe in C#, so I can integrate them with Unity. If you don’t know C#, you can freely ignore the code.

Definition

A quaternion is an ordered pair of 4 real numbers (w,x,y,z). We write this as

The letters i,j,k are not variables. Rather, they are independent axes. If you like, you can think of the quaternions as a 4 dimensional vector space.

The defining property of quaternions is:

Play around with it a bit and you can derive 6 more identites:

If you’ve worked with complex numbers, this should seem familiar. Instead of 2 parts of a complex number (the real and imaginary parts), we have 4 parts for a quaternion.

The similarity doesn’t end here. Multiplying complex numbers represents a rotation in 2 dimensions. Similarly, multiplying by a quaternion represents a rotation in 3D.

One curious thing to note: we have and . We switched around the terms and the product changed. This means that multiplying quaternions is kind of like multiplying matrices — the order matters. So multiplication is not commutative.

In Unity, the input is not given to us as a single true/false value, but a float between -1 and 1. So holding right increases the LeftRight input gradually until it reaches 1, avoiding a sudden jump in movement.

What’s ToUnityQuaternion? Well, it turns out that Unity already has a Quaternion class that does everything here and much more, so all this could have literally been implemented in one line if we wanted.

Anyways, let’s see the result.

As you can see, holding right turns the plane relative to itself now, and the green arrow stays still. Hooray!

Recently I was invited to compete in the CMOQR — a qualifier contest for the Canadian Math Olympiad. The contest consisted of eight problems, and contestants were allowed about a week’s time to submit written solutions via email.

After a few days, I was able to solve all of the problems except one — the second part of the seventh problem:

Seven people participate in a tournament, in which each pair of players play one game, and one player is declared the winner and the other the loser. A triplet ABC is considered cyclic if A beats B, B beats C, and C beats A.

Can you always separate the seven players into two rooms, so that neither room contains a cyclic triplet?

(Note: the first half of the problem asked the same question for six people — and it’s not too difficult to prove that no matter what, we can put them into two rooms so that neither the first nor the second room contains a cyclic triplet.)

But what happens when we add another person? Can we still put four people in one room, and three people in the other, so that neither rooms contain a cyclic triplet?

There are two possibilities here:

One, it’s always possible. No matter what combinations of wins and losses have occurred, we can always separate them into two rooms in such a way. To prove this, we’ll need to systematically consider all possible combinations, and one by one, verify that the statement is possible for each of the cases.

Two, it’s not always possible. Then there is some counterexample — some combination of wins and losses so that no matter how we separate them, one of the rooms has a cyclic triplet. This is easier to prove: provided that we have the counterexample, we just have to verify that indeed, this case is a counterexample to the statement.

But there’s a problem. Which of the cases does the solution fall into? That is, should we look for a quick solution by counterexample, or look for some mathematical invariant that no counterexample can exist?

Brute Force?

It would be really helpful if we knew the counterexample, or knew for sure what the counterexample was. What if we wrote a computer program to check all the cases? After all, there are only 7 people in the problem, and 7 choose 2 or 21 games played. Then since each game is either won by one player or the other, there are only 2^21 combinations overall (although that does count some duplicates). And 2^21 is slightly over two million cases to check — completely within the bounds of brute force.

So I coded up a possibility-checker. Generate all 2^21 possible arrangements, then for each one, check all possible ways to separate them into two rooms. If it turns out that no matter how we arrange them, a cyclic triplet persists, then display the counterexample. Simple.

I ran the program. It quickly cycled through every possible arrangement, three seconds later exiting without producing a counterexample.

Alright. So there’s no counterexample. I would have to find some nice mathematical invariant, showing that no matter what, there is always some way to group the players so that neither room has a cyclic triplet.

But no such invariant came. I tried several things, but in each attempt couldn’t quite show that the statement held for every case. I knew that there was no counterexample, but I couldn’t prove it. But why? There must be some tricky way to show that no counterexample existed; whatever it was, I couldn’t find it.

Brute Force poorly implemented

Reluctantly, as the deadline came and passed, I submitted my set of solutions without solving the problem. When the solutions came out a week later, the solution to this problem did not contain any tricky way to disprove the counterexample. Instead, what I found was this:

Let be seven players. Let beat when the difference .

Huh? A counterexample, really? Let’s look at it.

Everything is symmetric — we can ‘cycle’ the players around without changing anything. Also, if we take four players, two of them are consecutive. Let them be and .

At this point everything falls into place: in any subset of four players, three of them are cyclic.

But wait … my program had not found any counterexamples! And right here is a counterexample! The culprit was obvious (the reader may have foreseen this by now) — of course, there had to be a problem with my program.

Running my code through a debugger, I found a logic error in the routine converting binary numbers to array configurations, meaning that not all possible configurations were tried. As a result, the counterexample slipped through the hole.

After fixing the code, the program found not one, but a total of 7520 (although not necessarily distinct) counterexamples. Most of them had no elegant structure, but the solution’s configuration was among them.

When to Start Over?

It is true that the program could have been better written, better debugged. But how could you know whether a counterexample existed and your program didn’t find it, or if no counterexample existed at all?

In hindsight, it seems that writing the brute force program made me worse off than if I hadn’t written it at all. After the program ran without finding a single counterexample, I was confident that no counterexample existed, and set out about proving that, instead of looking for counterexamples or symmetry.

When you are stuck on such a math problem — that is, after making a bit of progress you get stuck — it might be profitable to start over. More often than I would like, I prove a series of neat things, without being able to prove the desired result. Then a look at the solutions manual reveals that a very short solution — one or two steps — lay in the opposite direction.

I’ll put an end to my philosophical musings of the day. Fortunately, the cutoff for the CMOQR was low enough that even without solving every single problem, I was still able to meet the cutoff.

Often, as I am minding my daily business, I would suddenly have an interesting thought, or some idea. Being likely to forget it a few minutes later, I would want to write the idea down. But being a lazy and slightly disorganized person, I would grab the nearest pencil and scrap paper, scribble down a few words, and throw the piece of paper aside somewhere.

Usually I would never retrieve that piece of paper again. But occasionally I do want to read the note I wrote. After a few minutes of searching I would usually find the piece of paper I needed. If I didn’t, oh well, it probably wasn’t that important anyway.

But, as you can imagine, this quickly gets out of hand.

Now you would like to do do the note taking on the computer, instead on small pieces of paper. After all, losing files is harder on the computer right?

But how would you do this? In the computer world it’s a bit harder to ‘throw aside’ a file. You have to give the file a name, and put it in some directory, essentially forcing you to be organized. Or if you don’t, files clutter up your computer desktop even faster than sheets of paper clutter up your actual desk.

So being the geek I am, I wrote up my experimental, software solution. Then writing a note would go as follows. You open up command prompt, and enter the note command:

Type up whatever you want, and then hit save:

So your note is saved in some obscure directory of your choosing. Now the Windows 7 search is pretty good (Windows Vista works too) that as long as your directory is indexed, Windows-F and typing pretty much anything in the search bar would quickly bring you to the file (it searches file contents):

Searching for “phone ashley” takes me directly to the note. The program just takes what we’ve written, and saves it in a unique file that also has a datestamp on it. For example the first note written on a day would have _1 appended to it, then _2, and so on.

I’m not yet sure how practical this idea would be, but here’s the C++ code for my program:

Problem 143 of Project Euler is an interesting geometry problem, and one notorious for having one of the lowest solvers to time released ratios of all the problems on the site: under 600 solvers in over 3 years as time of writing.

So what is the problem?

The problem is quite simple. We have a triangle , with all angles smaller than . Point is inside the triangle such that is at its minimum ( can always be determined uniquely).

We define a Torricelli triangle as one where , , , , , are all integers.

Generate all Torricelli triangles with .

The Fermat point

In triangle , with being the point such that is minimized, is called the Fermat point of triangle .

We can use Ptolemy’s theorem to locate the Fermat point in the triangle.

Let be a point such that is equilateral, then construct the circumcircle of :

We now show that must be at the intersection of and the circumcircle of .

Assume that does not lie on the circumcircle of . Then, according to Ptolemy’s theorem:

Equality only occurs when is cyclic. Since is equilateral, we have so we can divide it out to get

Adding to both sides:

Now if did lie on the circumcircle of , we would have , so must lie on the circumcircle of in order to be optimal.

So we know is on the circle. But now we assume that is not on . As we had before, .

Next if is not on , then , or by substitution,

Of course if were actually on , then , and the sum would be optimal.

This proves the optimal place for to be on the intersection of and the circumcircle of .

If we repeat this for the other three sides, we notice that is the intersection of the three circumcircles, and also the intersection of , , and :

Further, we see that quadrilaterals , , and are all cyclic. As , , and similarly and .

Designing an algorithm

We have enough information now to start working on an algorithm. Let us come back to the previous diagram:

So knowing that the central angles are all , we can apply the cosine law ():

A similar formula applies to sides and . We call a pair if is a perfect square.

We have found a Torricelli triangle if for three integers , , , all of , , are all pairs.

This leaves us with an outline of an algorithm:

Generate all pairs under the limit (with and being a square)

Sort the list of pairs and index them (to be explained in step 3)

For each pair in the list, search through the list to check if there exists some where is a pair and is a pair. We index the list to drastically improve searching time.

If an triple is found and , then mark as found. This is easily implemented as a bit array of size 120000 with all bits initially 0 and flipped to 1 when a sum is encountered.

Add up the sums and output

Step 1 takes up most of the time in this algorithm. To do this by ‘brute force’, we need an outer loop from 1 up to 120000, then another inner loop again from 1 to 120000 (actually a bit less), which is essentially operations, which is in the billions. An implementation of this algorithm would take about 3 minutes.

Parameterizing the a^2 + ab + b^2 = c^2 equation

It is well known that all primitive triples of the form can be generated by , , where and are coprime and are of opposite parity.

Can a similar parameterization be found for the equation equation?

Letting , and and dividing by , we get

which is the equation of an ellipse.

Originally we wanted to find integer solutions to the equation . This is equivalent to finding all rational points on the ellipse :

It is easy to see why: if and are rational points such that then and where , , are positive integers, and we arrive back at the original form of the equation.

Also in order for , , to be positive, we will only consider points on the ellipse in the first quadrant, with both and being positive. We do this by choosing a point on the ellipse, and from there combine the equations of the ellipse and the line.

Let us choose the point to be the first point of the line. Then the equation of the line is , where is valid only if it is positive and greater than 1.

Substituting into the equation of the ellipse:

Simplifying:

Now evaluating for :

Notice now that for and to be rational, we just need the slope to be rational. So we can write as where and are positive integers and .

Simplifying the expressions for and :

Since we defined and , we have a valid triple if , , and .

If and are not coprime, then neither will , , and be coprime as they will share a common divisor. So in order for to be primitive we will need and to be coprime.

Still this isn’t quite sufficient. Consider what would happen if . We have which can be written as which is divisible by 3. Next which is divisible by 3. Finally, which can be written as , again divisible by 3! So if , then the resulting triple is not primitive!

But it turns out that we can ignore all cases where . Given a triple , with , it is possible to find an identical parameterization giving .

As , we have and also . Thus we can have and where and are positive integers:

Multiplying by 3, , and . Combining the two and substituting we get

If we substitute and for and in the triple , we get:

So this proves that when , the case is nonprimitive and has already been covered by a smaller (possibly primitive) case. We are done.

Implementation in C++

Now using the new parameterization scheme, we can generate all pairs in about operations instead of operations. So we can loop up to , which cuts down the time down to about 2 seconds:

An important, yet basic problem in computational geometry is the following task: given a polygon (a sequence of points) and a point, determine if the point is inside, or outside the polygon.

To a human eye, this problem is trivial. After all, it seems obvious that point is outside the polygon, whereas point is inside it.

Computers, on the other hand, lack this simple intuition. It is not entirely obvious how to go about creating an algorithm that correctly determines this.

Triangles containing the origin

Let us first focus on an easier version of the problem. Instead of any arbitrary polygon, let our polygon be limited to three points, that is, a triangle. Instead of any point in space, let us determine if the triangle contains the origin.

We observe that we can draw a ray from our point outwards. Any direction will do. If the point is inside the triangle, then the ray will eventually hit a side of it.

On the other hand, if the point is outside the triangle, the ray may hit two sides of the triangle, or it may hit none. However, it will never hit exactly one side:

Aha; this is the approach. Now we can transpose our diagram to a Cartesian plane. For simplicity, let our ray be a line directly upwards from the origin, that is, the entire y-axis above the origin.

Our algorithm now takes the three line segments of the triangle, and counts how many of them intersect the y-axis above the origin.

There are three cases. In the first case, the triangle doesn’t intersect the y-intercept at all above the origin, thus the triangle doesn’t contain the origin:

In the second scenario, the triangle intersects the y-axis twice above the origin, in which case the triangle still does not contain the origin:

In the last case, the triangle intersects the y-axis once above the origin, in which case the origin is in the triangle:

So the algorithm is, for every line segment (pair of coordinates) in the triangle:

Calculate the y-intercept, say,

If this line segment doesn’t cross the y-axis above the origin so throw it away. Otherwise, continue.

Make sure that the line segment actually crosses the y-axis at all. Or, check that the two points are on opposite sides of the y-axis.

You have a board of 7 squares in a single row, three black, and three red pieces. At the initialization of the game, the three black pieces are placed on the first (leftmost) three squares, and the three red pieces on the last (rightmost) three squares. This leaves one square in the middle.

The aim is to get all the red pieces onto the left side, and the black pieces on the right side, or in other words to reverse the position of all the pieces. This is done by executing a sequence of two possible moves:

One, a piece may move onto an empty square adjacent to it.

Two, a piece may hop over a piece of the opposite color adjacent to it, provided that the square on the other side of the second piece is empty.

Furthermore a move may only bring a piece closer to its destination; that is, it is not allowed to backtrack.

This particular instance of the problem (of three pieces each and seven squares) can be completed in 15 moves:

While this sequence completes the puzzle in 15 moves, it is not so obvious that this is the shortest sequence, that there is no shorter sequence that completes the puzzle in fewer than 15 moves. It is possible to prove a generalized result mathematically.

A formula for the number of moves

Let us generalize the game. We define an n-game to be an instance of the puzzle, with n pieces of each color separated by an empty square. For an n-game, the board contains squares.

Proposition: to complete an n-game, exactly moves are required.

In order for the n black pieces to completely trade places with the n red pieces, each black piece must hop over, or be hopped over by n red pieces.

Consider the pieces to be paired up as follows:

In the sequence of the game, each piece trades places with its partner. (Note that this is correct because the order of pieces of the same color cannot change).

It can be seen that the relative displacement between two members of a pair is . This is because each piece travels squares to get to its destination. In other words a piece travels squares relative to its partner.

A move by either member of a pair adds a relative displacement of 1, while a hop adds a relative displacement of 2.

We have already shown that each pair must have, between them, n hops. This is equivalent to a relative displacement of 2n. Two additional moves are required to bring the relative displacement up to . Therefore, moves are required for a pair to exchange places.

As there are n pairs in the game, the total number of moves required to complete the game is , as desired.

Implementation of a solver

Using the technique of a recursive depth-first search, it is a fairly easy programming exercise to write a solver for this puzzle. The solver produces the actual sequence of moves to complete an n-game.