Geometry Algorithm Difficulties

In implementing the geometry algorithms the largest problem I’ve had is overcoming the fact that floating point numbers are not numbers. This issue crops up in most of the lower level algorithms such as testing if a point is inside a polygon or doing the boolean difference of two polygons etc.

The solution is to use a fuzzy representation of polygons. That is if a point is within \epsilon of a polygon then it is inside the polygon. This method works fairly well and can be implemented but:

1. I didn’t write my functions using this method from the beginning so I’ve had to gradually retrofit my functions which is a pain and is fairly difficult because testing becomes a nightmare.

2. The implementation itself can be a bit tricky. In point of fact the selection of \epsilon depends on the scale of a polygon (e.g. If the edges of a polygon have length 1000 then \epsilon should be larger than when a polygon has edges of length 0.001)

For the above two reasons I’ve been contemplating using CGAL for implementation of low level functions. The only troubles are I’d need to rewrite my high level algorithms in C++ and debugging/testing algorithms in Python is much much better than C++.

That’s my current quandary. As of now I’m still hacking on the Python code and fixing bugs but I’ve also created a Subversion branch for future C++ development using CGAL.

Youtube Video Ratings

When viewing most youtube videos you can rate the video as in the following example.

The overall video rating, I assume, is the mean of the ratings, that is
\textrm{video rating} = \frac{\Sigma rating_i}{\# \textrm{number of ratings}}.

A problem with this rating technique is sensitivity to outliers e.g. people voting 0 or 5. One solution is for youtube to use the median instead of the mean. The median can be computed in linear time but the computation is significant, see selection algorithms.

But would using the median actually change the overall video rating? We need to analyze how people rate videos. Is someone’s rating of a video independent of the previous ratings?

It is possible that people will give a high rating if they believe the video has been rated too low and vice versa. In this scenario the mean is the same as the median!

It is also possible that people rate a video independent of the previous ratings in which case using the median is desirable.

To distinguish between which of these scenarios is more correct we can give the following two tests.

1. Use a common video in both tests and a large random set of users. This is so we can compare the tests against each other.

Test1: The video rating is given by the mean of the ratings and a user can see the video rating before they rate it.

Test2: The user cannot see the video rating before they rate the video and again the video rating is given by the mean of the ratings.

If the overall video rating is the same in both tests then it strongly suggests that users rate videos independent of the previous ratings. Of course if the distribution of user ratings is normal than the mean and median will be the same the same in Test1 and Test2. One way to modify the experiment if this happens is to add a few outlier ratings to Test1 and then see if the users correct the rating to what it would be without the outliers.

Now if the video rating in Test1 agrees with the median of the video ratings in Test2, then we know that using the mean is a good substitute for the median.

Protected: Vert-Solve Running Time

I’m interested in the running time of Poly-Solve and Vert-Solve. First note that because Poly-Solve and Vert-Solve use the same combinatorial approach they have the same big O running time.

Because I gave a more detailed listing for Vert-Solve, I’ll analyze Vert-Solve’s running time. For clarity I’ve reposted the Vert-Solve algorithm below.

function vert_solve(P = {set of puzzle pieces}, s = silhouette)
if P == ∅: return true

for (&p, v), w ∈ Pxs:
#By x I mean the Cartesian product, so (&p, v), w varies over all possible vertex pairs from P and s.

  g=position_vertex_pair((&p, v), w) #create a generator to position the vertex, v, so that it coincides with w, that is translate p so that v and w have the same coordinates. Note that there are several possible orientations of v‘s line segments (that’s why a generator is used).
    P’ = Pp #remove puzzle piece p from the set of puzzle pieces.

   while g.next_orientation() == true: #continue loop, until all valid (by valid I mean p is a subset of s) possible orientations of v‘s lines segments are exhausted.
          s’ = s\p #subtract puzzle piece p from shape s and set s’ equal to the shape resulting from the subtraction.

           if vert_solve(P’, s’) == true: return true
end for
return false

Let n = |P| + |s|, where |x| means the number of line segments of polygon x.
Note that the number of line segments is at most 2 * number of vertices.
For the set P, |P| means the total number of line segments for all the polygons in P.

Part A.

1. Let R(n) denote the running time of vert-solve for an input with n line segments i.e. n = |P| + |s|.

2. |Pxs| ≤ n2.

3. position_vertex_pair((&p, v), w) takes O(1) time, since it only creates a generator.

4. P’ = P – p, takes O(n) time (simply scan through the list P looking for vertices of p).

5. g.next_orientation() is the same algorithm as match_vertices which I gave in my post on 2D Geometric Algorithms. The the running time for match_vertices is O(m4), where m = |p| + |s| ≤ n. Because we’re doing big O running time I’ll let m = n, so g.next_orientation() takes O(n4) time.

6. s’ = s\p takes O(n3), see my post on 2D Geometric Algorithms for the running time analysis. Note that |s’| ≤ |s| – 3.

7. The running time for vert_solve(P’, s’) is R(n – 1) since |P’| + |s’| must be at least one smaller than |P| + |s|. This is because we’ve removed |p| line segments from P and added at most |p| – 1 line segments to s.

Part B. The while loop given by 5 runs at most n times and each loop takes n4 + n3 + R(n – 1) time. So the total time is
n^5 + n^4 + n\cdot R(n - 1).

The for loop given by 2 runs at most n2 times, so the total time is R(n) = n^2 \cdot (O(1) + n + n^5 + n^4 + n\cdot R(n - 1). Which simplifies to R(n) = n^7 + n^3R(n-1) (we only care about big O running time so the lower terms get absorbed into the n^7 term). Expanding yields n^7 + n^3(n-1)^7 + n^3(n-1)^3(n-2)^7 + \ldots + (n!)^3 \cdot 1^7 which is \le n^7 + n^3n^7+(n(n-1))^3n^7+\ldots + (n!)^3n^7 which equals n^7(1+n^3+(n(n-1))^3 + \ldots + (n!)^3), which is \le n^7(1 + n \cdot (n!)^3) \le n^7 + n^8 \cdot (n!)^3.

Thus the total running times is O(n^8{n!}^3). Hopefully I can improve this to exponential because n^8{n!}^3 is terrible!


Today I spent time learning how to juggle. I’ve wanted to juggle for at least three or four years but I never read any articles on how to juggle. As I was browsing around the internet watching various videos I found this video of Chris Bliss juggling. Then I searched around for a how to on juggling and found an excellent article.

And thus I was off and running. It took about four hours to get the basic routine of juggling three tennis balls. At this point I can do 10 catches.

I’ve uploaded a video of me doing a 4 toss juggle.


If it isn’t clear what I said in the video the transcript is “Oh, wait you’ve got to replace the battery”. By battery I meant the camera battery.

And that was how my Saturday was spent.

Blog Categories

Currently WordPress has a very simply category mechanism where there exists a simple list of categories and the blog author can select any number of them for a post.

For example in assigning categories for blog posts I’ve often used the category Computational Geometry, this category is a strict subset of Computer Science, so I’d like readers to somehow be able to access these posts when they’re browsing the computer science category.

From the author’s point of view this category style is very easy to use. But for a reader it doesn’t allow enough interconnection between categories. For example the category computational geometry is a subcategory of computer science but when browsing the computer science posts there is no way to access the computational geometry posts.

What I want is to browse the computer science category and have the subcategories and computer science posts listed in the same way as Computer Science on Wikipedia.

Yet another project I need to do 🙂

2D Geometric Algorithms

For implementing Vert-Solve or Poly-Solve several lower level algorithms are needed. In particular I want to describe the most important ones which are:

1. subset(s1, s2) #Test if s1 is a subset of s2.

2. match_vertices(v1, v2) #Orient vertex v1 and its respective shape so that v1 and one of v1‘s line segments coincide with v2 and one of v2‘s line segments.

3. subtract_shape(s1, s2) #Subtract s1 from s2 and return the new shape.

All three algorithms make use of more basic algorithms such as line intersection.

1. The following python code implements subset. I note that shape.nodes is a list of points defining shape (i.e. the list of vertices). The function segment_contained_in tests if the line segment (node1, node2) is a subset of shape2), it will return true if (node1, node2) are part of shape2’s boundary. The function inside_shape test if (r, c) is strickly inside shape2 that is if (r,c) is on the boundary of shape2, then inside_shape((r, c), shape2) returns false.

def subset(shape1, shape2):
    """returns true if shape1 is a subset of shape2, returns false otherwise"""
    for (r,c) in shape1.nodes:
        if not inside_shape((r,c),shape2) and not on_boundary((r,c),shape2):
            return False
    for (node1,node2) in shape1.line_segments:
        if not segment_contained_in((node1,node2), shape2):
            return False
    return True

Let n = max(|shape1.line_segments|, |shape2.line_segments|), where || denotes set size. The number of vertices in shape1 or shape2 is strickly less than n. The running time of inside_shape and on_boundary is linear in n; segment_contained_in also has linear running time in n. Thus subset has O(2n2) running time.

2. I use the following code for matching vertices. It makes use of python generators.

def match_vertices(i, j, shape1, shape2):
    vertex1 = shape1.vertices[i]
    vertex2 = shape2.vertices[j]
    (point1,segments1) = vertex1
    (point2,segments2) = vertex2
    #translate vertex1 ontop of vertex2
    shape1.translate((point2[0] - point1[0], point2[1] - point1[1]))
    #for each segment1 in vertex1 and for each segment2, make segment1 parallel to segment2. This way we 
    #try all orientations that can make shape1 contained inside of p_shape_to_fit.
    for segment1 in segments1:
        for segment2 in segments2:
            #translate so point2 is the origin
            tsegment1 = translate(segment1, (-point2[0], -point2[1]))
            tsegment2 = translate(segment2, (-point2[0], -point2[1]))
            vec1 = tsegment1[0]
            if(points_equal(tsegment1[0], [0, 0])): vec1 = tsegment1[1]
            vec2 = tsegment2[0]
            if(points_equal(tsegment2[0], [0, 0])): vec2 = tsegment2[1]
            #theta1 is the angle that segment1 makes with the positive x-axis.
            #0 <= theta1 < 2pi
            theta1 = angle(vec1)
            #theta2 is the angle that segment2 makes with the positive x-axis
            #0 <= theta2 < 2pi
            theta2 = angle(vec2)
            rotation_angle = theta2 - theta1
            p1 = rotate2d(vec1, rotation_angle)
            p1 = (p1[0] + point2[0], p1[1] + point2[1])
            if not inside_shape(p1, shape2) and not on_boundary(p1, shape2):
            nodes = shape1.nodes[:]
            line_segments = shape1.line_segments[:]
            verts = shape1.vertices[:]
            #rotate segment1 into segment2
            shape1.rotate(point2, rotation_angle)
            if subset(shape1, shape2):
                yield True
            shape1.nodes = nodes
            shape1.line_segments = line_segments
            shape1.vertices = verts
    yield False

Let n = max(|shape1.line_segments|, |shape2.line_segments|).
The number of line segments attached to vertex1 is at most n. It appears that there should only be two line segments attached to vertex1 but consider the case of a disconnected shape. In the shape below one of the vertices has 4 line segments attached to it.

Disconnected Polygon

The only important call for running time analysis inside the inner loop is subset(shape1, shape2), which as I’ve shown in 1. has n2 running time. Thus the total running time is O(n4<). In non pathological cases where |segments1| = 2 and |segments2| = 2, the running time is O(n2).

The above gives the total running time of the generator but I’m also interested in the running time for one iteration of match_vertices. That is what is running time for one iteration of match_vertices. For this we need to find the running time between yield statements. The worst case scenario is that there are no matches, that is the the first call to match_vertices returns false. In which case the call takes O(n4) time.

The running time for getting the next iteration given by generator

3. The following is the python code for subtract_shape

def subtract_shape(p_shape1, p_shape2):
    new_line_segments = []
    if not subset(p_shape1, p_shape2): return None
    for line_segment1 in p_shape1.line_segments:
        split_seg = split_segment(line_segment1, p_shape2)
        #the part of the line segment that is outside p_shape2
        outside_part = split_seg['outside']
        #None of the line segment should be outside of p_shape2 since p_shape1 is a subset of p_shape2
        assert(outside_part == [])
        #the boundary part of a line segment does not form part of the boundary for the new shape, so we can ignore this part of the line segment
        boundary_part = split_seg['boundary']
        #this is the part of the line segment that will be included in the boundary of the new shape
        inside_part = split_seg['inside']
    for line_segment2 in p_shape2.line_segments:
        split_seg = split_segment(line_segment2, p_shape1)
        #the part of the line segment that is outside p_shape1
        outside_part = split_seg['outside']
        #the boundary part of a line segment does not form part of the boundary for the new shape, so we can ignore this part of the line segment
        boundary_part = split_seg['boundary']
        #this is the part of the line segment that will be included in the boundary of the new shape
        inside_part = split_seg['inside']
        #none of the line segments in p_shape2 should be in p_shape1 since p_shape1 is a subset of p_shape2
    #remove zero length line segments
    #the zero length line segments are a result of the floating point errors.
    for seg in new_line_segments:
        if(nodes_equal(seg[0], seg[1])):
    return shape('', '', new_line_segments)

I believe split_segment runs in quadratic time with respect to |p_shape2.line_segments|. In which case subtract_shape has O(n3) running time.

Protected: Algorithm Heuristics

In my past posts on solving tangram puzzles I’ve given two slightly different algorithms. Neither algorithm as I gave it had any heuristics. I aim to show you a few simple heuristics that should work in non-pathological cases.

The most obvious heuristic for the Poly-Solve algorithm is to sort the puzzle pieces by size (or bounding box) before you pass them into position_puzzle_pieces. When I implemented Poly-Solve this did substantially speed up the running time. The reason this heuristic works is because large puzzle pieces have fewer possible orientations inside the silhouette which implies that fewer combinations need to be tried. The pathological case for this algorithm is when the largest puzzle piece must be fitted last, which when you think about it, is truly pathological.

For the Vert-Solve heuristics I define P to be the set of puzzle pieces, where each element is a vertex v together with a pointer to its puzzle piece.

1. One heuristic for the Vert-Solve algorithm is to try using vertices from the puzzle pieces that are perfect matches. By perfect match for a vertex, v, I mean v is identical modulo translations, rotations, and reflections to a vertex in the silhouette. With this heuristic the complexity of the new silhouette should be reduced each time vert_solve recurses.

2. Another possible heuristic for Vert-Solve is to sort the vertices in P from largest to smallest, where by size I mean the length of the shortest line segment of the vertex. And also sort the vertices in s. Thus when we select the first tuple ((&p,v), w) from Pxs, the v, w are the largest possible vertices. The logic behind this heuristic is similar to the above heuristic for Poly-Solve.

3. We can combine 1 and 2 to create a new sorting algorithm for P, using weights for the importance of 1 and 2. In fact there are lots of ways to sort P but the running time of the sorting algorithm can be important.

That is all for today. It’s easy to think of more heuristics for both algorithms but I need to test them before posting because otherwise I’m only speculating.