This story actually happened in a few hours, but the ride was the must fun I have had in the last few months(professionally), so I thought in sharing it.

## The task

I had to implement the Graham’s Scan in Scala. Given that the algorithm has been around since the 70’s, you’d think it would be a walk in the park. It was more like a game of rapid chess.

1. Find

p_{0}in Q with minimum y-coordinate (and minimum x-coordinate if there are ties).

2. Sorted the remaining points of Q (that is, Q − {p_{0}}) by polar angle in counterclockwise order with respect top_{0}.

3. TOP [S] = 0 ▷ Lines 3-6 initialize the stack to contain, from bottom to top, first three points.

4. PUSH (p_{0}, S)

5. PUSH (p_{1}, S)

6. PUSH (p_{2}, S)

7.fori= 3tom▷ Perform test for each point p_{3}, …, p_{m}.

8.do whilethe angle between NEXT_TO_TOP[S], TOP[S], andp_{i}makes a nonleft turn ▷ remove if not a vertex

9.doPOP(S)

10. PUSH (S,p_{i})

11.returnS

As you can see, the “good stuff” happens in a stack. All the versions describing the algorithm that I could found, were written this way. So instead of starting from the abstraction of the algorithm, I had to abstract an opinionated interpretation, to come up with a functional version.

First things first, I created a data set of points for plotting in Octave. Octave already has a convex hull function, so I could **see** the solution before hand.

X = [0.410036859415169 0.16990953347486892 0.7438936112344358 0.9879763068249436 0.09850362854534644 0.27745454018107196 0.2463477117635604 0.0369350463381225 0.7356804213480419 0.979020964716884 0.4104133424281955 0.9566561367393451 0.36405116056182674 0.6743515475246237 0.2014896354108069 0.5937635359142525 0.9990876238009762 0.4753053860797861 0.15127539785912514 0.3109492347362207 0.3710641858259751 0.4981026670484303 0.472172885898987 0.9804204746053471 0.20208956184140103 0.7386385686849415 0.3182368631449125 0.4529941453549946 0.3035123394853474 0.3785220082038615 0.6028255171683854 0.4754941518809386 0.02885621173133124 0.619495881160238 0.03310900644639769 0.2791084654454382 0.8674411544744653 0.393894922057599 0.271592427317065 0.31136184917695364]; Y = [0.0964872038572564 0.8409000674315462 0.5513623119954284 0.9591944778976642 0.3152016460947692 0.16362690913082834 0.31570455065654346 0.5538361485253371 0.7670638114116081 0.015668444677785942 0.3467907029082512 0.07399794456597841 0.12236562106445759 0.9981954727371612 0.9308695268083278 0.8560770018856956 0.6579161973082933 0.5959520081450651 0.13661844095029363 0.6514314696732396 0.19298645312277207 0.00994355691826554 0.7010697212377054 0.1634450128412981 0.3793185478977382 0.868604478585853 0.4590548622219143 0.18377356484707785 0.2866066023269005 0.8213505729786262 0.38257319264118383 0.4566803134679752 0.4748548161642837 0.26842322605193014 0.5464910173382869 0.6503654422352256 0.9950193150372821 0.730931140118583 0.5869777402729242 0.06419667711655552]; k = convhull(X,Y); plot (X(k), Y(k), "r-", X, Y, "b+");

## Starting

First, I found a cool way of sorting the points by polar angle:

val coordList = origin :: coords.filterNot(_ == origin). sortBy(point => atan2(point._2 - origin._2, point._1 - origin._1))

Just to be sure, I plotted the sorted points without any further processing:

In another article, I found this amazing formula for deciding if a point makes a left turn in relation with a line created by the other two:

// isLeft(): tests if a point is Left|On|Right of an infinite line. // Input: three points P0, P1, and P2 // Return: > 0 for P2 left of the line through P0 and P1 // = 0 for P2 on the line // < 0 for P2 right of the line // See: Algorithm 1 on Area of Triangles inline float isLeft( Point P0, Point P1, Point P2 ) { return (P1.x - P0.x)*(P2.y - P0.y) - (P2.x - P0.x)*(P1.y - P0.y); }

So, P2 is the new point, P1 the last point added to the hull, and P0 the second to last.

## Misinterpretation.

In my first version I was doing the following:

- Add the first three points to the hull.
- After every new point, check if it makes a left turn in relation with the two previous
- If it does, add it to the hull.
- If it doesn’t, remove the last point from the solution and add the new one to the hull.

Sounds reasonable right?

I wrote the code and plotted the result in Octave.

Yes! No! It found all of the points of the hull, but it kept also the wrong ones. Is hard when a perfectly reasonable code doesn’t do you think it should do. Anyway, I knew that the algorithm was fine, I knew that the points were sorted in a proper way. So it had to be my interpretation of the algorithm. After a careful reading of the pseudo code above, I realized that before adding a new point, I needed to scan the already found vertices in relation with the new point, and removing P1 every time a non-left angle were found.

## foldRight

The coolest moment for me was realizing that instead of the stack, I could just use foldRight. That gave away the solution.

case class Point(x: Double, y: Double) { //improves readability def goesLeft(p0: Point, p1: Point): Boolean = (p1.x - p0.x) * (this.y - p0.y) - (this.x - p0.x) * (p1.y - p0.y) > 0 } def scan(coords: List[Point]): List[Point] = { //instead of using PUSH and POP, foldRight def addToHull(hull: List[Point], new_point: Point): List[Point] = //the new point is always added but filtering all the non-left turn vertices from the current hull new_point :: hull.foldRight(List.empty[Point]) { case (p1, currentHull@(p0 :: _)) => if (new_point.goesLeft(p0, p1)) p1 :: currentHull else currentHull case (p, currentHull) => p :: currentHull } val min = coords.minBy(_.y) //adding the min to close the polygon on octave min :: coords.sortBy(point => atan2(point.y - min.y, point.x - min.x)) .foldLeft(List.empty[Point])(addToHull) }

I created a project in GitHub, that also has a java 8 solution, not very functional. What the test suites actually do is to generate random samples and a Octave file with the solution for plotting.