## Posts Tagged ‘algorithm’

### Avoiding loops and routing failures with A*

#### The problem

With GraphPaper, generating path with the A* algorithm I ran into some interesting cases where the algorithm would fail to generate a path between the start and end points:

Note that the path should go around, not directly across the big rectangle in the center.

There’s a few things being masked here, as GraphPaper will always try to produce a result, even if it’s non-optimal. So we need to disable:

• The mechanism that makes a path directly from start to end if there’s a failure in generating a path
• Path optimization

With these disabled, and some annotations showing the order in which the path segments were generated, we can get a better idea of what’s happening:

Here we see that the algorithm generates a path that loops around its adjacent object (1-12), the path moves away from the object (13-14), then moves right back towards the object (15-18), then we’re out of accessible routing points, so we have a failure and we end the path by going directly to the endpoint. It’s important to note that routing points are not grid points here, they’re points designated around objects (the rectangles with the dots in the center) and anchors (the rectangles where paths start and end).

#### Why is this happening?

This might seem like a bug in the A* implementation but, as far as I can tell, the algorithm is performing as it should and the issue here comes down to the path cost computation done in the main loop. A* doesn’t actually define how to compute cost, but GraphPaper uses run-of-the-mill Euclidean distance as a metric.

In the A* main loop, the cost of going from our current point to a given point, n, is given by:

g(n) is the total length of the path from the starting point to n (sum of the length of all pieces of the path). In practice, we keep track of the current path length, so we only need to compute and add the straight-line length from current to n:

h(n) is the straight-line length from n to the goal:

The issues we’re seeing surface as a result of f(n) being based solely on distance computations, note that:

• We may have a good point, n, it’s close to our goal but requires a large jump from our current location, the large increase in g(n) leads to another, closer, point being prioritized
• We may have a bad point, n, it’s further from our goal but only requires a small jump from our current location, h(n) increases a bit but the relatively small increase in g(n), leads to this point being prioritized over a better one

The implication here is also that shorter jumps become preferable, even if we end up going in the opposite direction of the goal, as such jumps are safer (i.e. less likely to lead to an blocked/invalid path).

#### Smarter routing

To improve the routing, I wanted to adjust the cost computation such that moving to points that and in the direction of (and thus closer) to the goal are prioritized, so I introduced a new term, t(n), in the cost computation:

To compute t(n), we first compute the dot product of 2 normalized vectors:

• vecCurrentTo(endpt): from current to the endpoint, representing the ideal direction we should head in
• vecCurrentTo(n): from current to n, representing the direction to n

The dot product ranges from [-1, 1]; tending to -1 when n is in the opposite direction of the ideal, and tending to 1 when n is in the same direction as the ideal.

We then scale the dot product result as we can’t just bias our overall cost, f(n), by a range this small, so we scale by the straight-line length from current to n. Using this length as a scaling factor also serves to bias towards longer paths in good/positive directions directions.

In code, the functions g(n), h(n) and t(n) components are computed like this:

```// g(n) = length/cost of _startPoint to _vp + _currentRouteLength const currentToVisibleLength = (new Line(currentPoint, visiblePt)).getLength(); let gn = currentToVisibleLength + _currentRouteLength; // h(n) = length/cost of _vp to _endPoint let hn = (new Line(visiblePt, _endPoint)).getLength(); // t(n) = // a. get the relationship between the 2 vectors (dot product) // b. scale to give influence (scale by currentToVisibleLength) // .. using currentToVisibleLength as scaling factor is influence towards longer paths in good directions vs shorter paths in bad directions const vecToTargetIdeal = (new Vec2(_endPoint.getX() - currentPoint.getX(), _endPoint.getY() - currentPoint.getY())).normalize(); const vecVisibleToEndpt = (new Vec2(visiblePt.getX() - currentPoint.getX(), visiblePt.getY() - currentPoint.getY())).normalize(); const tn = vecToTargetIdeal.dot(vecVisibleToEndpt) * currentToVisibleLength;```

With the updated path cost computation, we successful get a path to the endpoint and, even with path optimization off, we get a much better path overall with less unnecessary points in the path and no looping around objects:

It’s also noting how this performs as the position of the object changes (again, path optimization is disabled and the mechanism that makes a path directly from start to end on routing failure is also disabled):

Without t(n)

With t(n)

In pretty much every case, applying the t(n) factor to the path computation produces a better path with no failures in routing to the endpoint.

### A* path optimization

The A* algorithm, using a straight-line distance heuristic function, is great in terms of performance, but yields a number of cases where the paths produced are not optimal.

For example, here is a path created in ScratchGraph using GraphPaper (the underlying library used for creating the connectors):

I suspect the non-optimal result is even more pronounced here given that the path is determined based on specific routing points that exist around the objects (you can spot these as the inflection points in the path), instead of a uniform grid. That said, I don’t want to use a uniform grid; while it’s likely a non-issue with A*, there is also the cost of computing which routing points are accessible vs blocked, and that cost grows quickly as the possible number of routing points grow.

An approach that works well, and doesn’t drastically increase the cost of path computation, is simplifying the generated path by checking if corresponding points in the path are visible to each other. If so, we can remove the intermediate points, and simply connect those 2 points together. Here’s a snipped from the GraphPaper codebase:

```/** * * @param {Point[]} _pointsInRoute * @param {Function} _arePointsVisibleToEachOther */ optimize: function(_pointsInRoute, _arePointsVisibleToEachOther) { let start = 0; let end = _pointsInRoute.length - 1; while(true) { if((end-start) <= 1) { start++; end = _pointsInRoute.length - 1; if(start >= _pointsInRoute.length-2) { break; } } if(_arePointsVisibleToEachOther(_pointsInRoute[start], _pointsInRoute[end])) { _pointsInRoute.splice(start + 1, (end-start) - 1); end = _pointsInRoute.length - 1; } else { end--; } } }```

The function works as follows:

• Begin with the points `start` (first point in the path) and `end` (last point in the path).
• Relative to `start`, check if the other points in the path are visible to it (in the code above, we iterate backwards from the endpoint). If we find that `start` is visible to another point, we eliminate the intermediate points from the path.
• Once we get to `end` being the point after `start`, update `start` to the next point in the path and reset `end` to whatever the last point in the path is.
• Repeat the latter 2 steps until we’ve checked all corresponding points in the path (`start` is the point directly preceding `end`).

Optimizing the path shown above with this method yields an optimal path:

In terms of performance, the cost is based on the number of points in the path and the cost of whatever computation the `_arePointsVisibleToEachOther()` function does. For my use-cases, paths have relatively few points and `_arePointsVisibleToEachOther()` consists of a number of fairly fast line-line intersection checks, so the method is fairly cheap and there’s no significant decrease in performance when generating paths.

### Brute-force convex hull construction

I’ve been experimenting a bit with convex hull constructions and below I’ll explain how to do a brute-force construction of a hull.

It’s worth noting up-front that the brute-force method is slow, O(n3) worst case complexity. So why bother? I think there are a few compelling reasons:

• The brute-force method expresses the fundamental solution, which gives you the basic building blocks and understanding to approach more complex solutions
• It’s faster to implement
• It’s still a viable solution when n is small, and n is usually small.

#### What is a convex hull?

You can find a formal definition on Wikipedia. Informally, and specific to computational geometry, the convex hull is a convex polygon in which all points are either vertices of said polygon or enclosed within the polygon.

#### Brute-force construction

• Iterate over every pair of points (p,q)
• If all the other points are to the right (or left, depending on implementation) of the line formed by (p,q), the segment (p,q) is part of our result set (i.e. it’s part of the convex hull)

Here’s the top-level code that handles the iteration and construction of resulting line segments:

`/** * Compute convex hull */var computeConvexHull = function() { console.log("--- "); for(var i=0; i<points.length; i++) { for(var j=0; j<points.length; j++) { if(i === j) { continue; } var ptI = points[i]; var ptJ = points[j]; // Do all other points lie within the half-plane to the right var allPointsOnTheRight = true; for(var k=0; k<points.length; k++) { if(k === i || k === j) { continue; } var d = whichSideOfLine(ptI, ptJ, points[k]); if(d < 0) { allPointsOnTheRight = false; break; } } if(allPointsOnTheRight) { console.log("segment " + i + " to " + j); var pointAScreen = cartToScreen(ptI, getDocumentWidth(), getDocumentHeight()); var pointBScreen = cartToScreen(ptJ, getDocumentWidth(), getDocumentHeight()); drawLineSegment(pointAScreen, pointBScreen); } } }};`

The “secret sauce” is the whichSideOfLine() method:

`/** * Determine which side of a line a given point is on */var whichSideOfLine = function(lineEndptA, lineEndptB, ptSubject) { return (ptSubject.x - lineEndptA.x) * (lineEndptB.y - lineEndptA.y) - (ptSubject.y - lineEndptA.y) * (lineEndptB.x - lineEndptA.x);};`

This is a bit of linear algebra derived from the general equation for a line.

The result represents the side of a line a point is one, based on the sign of the result. We can check if the point is on the left or on the right, it doesn’t matter as long as there is consistency and the same check is done for all points.

#### How it looks

I made a few diagrams to show the first few steps in the algorithm, as segments constituting the convex hull are found. The shaded area represents our success case, where all other points are to the right of the line formed by the points under consideration. Not shown are the failure cases (i.e. one or more points are on the left of the line formed by the points under consideration).

#### Code and Demo

You can play around with constructing a hull below by double-clicking to add vertices.

You can find the code on GitHub.

### Implementing a “did you mean…?” function

A while ago I became interested in how one would go about implementing something akin to Google Search’s “did you mean…” function. This answer from StackOverflow provides a good overview of what Google does, which involves looking at a user’s incorrect entry as well as a subsequent correction provided by the user. Data mining both the incorrect entry and correction, for millions of users and billions (trillions?) of entries, Google Search can thus make an intelligent guess as to what a user really meant when an incorrect entry is submitted. While most applications can’t do this at Google Search’s scale or generality (which is a span of terms from across the entire web), I can certainly see this model working for smaller applications, which only need to deal with a smaller subset of terms (this blog for example only needs to handle terms that are in my posts). Remove the data capture aspect, providing a fixed database of terms from which to search, and implementation becomes even simpler!

Digging deeper into implementation details, there’s the problem of figuring out how closely an input string (X) matches each string in a database of terms (Ti). The closeness or distance here is the edit distance between the 2 strings, or the minimum number of edits it takes to turn one string into another. There are a number of edit distance algorithms, but the Levenshtein distance seems to be a popular choice.

For a simple did-you-mean suggester, computing the edit distance is the crux of the method.

To test things out, I did a simple project in C++. Using a straightforward implementation of the Levenshtein distance and a vector of 112 chemical elements (up to copernicium), I wrote a program that would prompt the user for the name of an element, if the element was found it would output “ELEMENT FOUND”, if not it would suggest the name of an element based on the user’s input.

### Includes

`#include <iostream>#include <string>#include <vector>`

### Levenshtein distance implementation

`int ld(const std::string& strA, const std::string& strB){    int lenA = strA.length() + 1;    int lenB = strB.length() + 1;    int** mat;    mat = new int*[lenA];    for(int i=0; i<lenA; i++)    {        mat[i] = new int[lenB];    }    for(int i=0; i<lenA; i++)    {        mat[i][0] = i;    }        for(int j=0; j<lenB; j++)    {        mat[0][j] = j;    }    for(int i=1; i<=lenA-1; i++)    {        for(int j=1; j<=lenB-1; j++)        {            if(strA[i-1] == strB[j-1])            {                mat[i][j] = mat[i-1][j-1];            }            else            {                mat[i][j] = std::min(mat[i-1][j-1]+1, std::min(mat[i-1][j] + 1, mat[i][j-1] + 1) );            }        }    }    int ret = mat[lenA-1][lenB-1];    // memory cleanup    for(int i=0; i<lenA; i++)    {        delete [] mat[i];    }    delete [] mat;    return ret;}`

### Function to construct std::vector of chemical elements

`std::vector<std::string> make_elements_vector(){    std::vector<std::string> elements;    elements.push_back("hydrogen");    elements.push_back("helium");    elements.push_back("lithium");    elements.push_back("beryllium");    elements.push_back("boron");    elements.push_back("carbon");    elements.push_back("nitrogen");    elements.push_back("oxygen");    elements.push_back("fluorine");    elements.push_back("neon");    elements.push_back("sodium");    elements.push_back("magnesium");    elements.push_back("aluminium");    elements.push_back("silicon");    elements.push_back("phosphorus");    elements.push_back("sulfur");    elements.push_back("chlorine");    elements.push_back("argon");    elements.push_back("potassium");    elements.push_back("calcium");    elements.push_back("scandium");    elements.push_back("titanium");    elements.push_back("vanadium");    elements.push_back("chromium");    elements.push_back("manganese");    elements.push_back("iron");    elements.push_back("cobalt");    elements.push_back("nickel");    elements.push_back("copper");    elements.push_back("zinc");    elements.push_back("gallium");    elements.push_back("germanium");    elements.push_back("arsenic");    elements.push_back("selenium");    elements.push_back("bromine");    elements.push_back("krypton");    elements.push_back("rubidium");    elements.push_back("strontium");    elements.push_back("yttrium");    elements.push_back("zirconium");    elements.push_back("niobium");    elements.push_back("molybdenum");    elements.push_back("technetium");    elements.push_back("ruthenium");    elements.push_back("rhodium");    elements.push_back("palladium");    elements.push_back("silver");    elements.push_back("cadmium");    elements.push_back("indium");    elements.push_back("tin");    elements.push_back("antimony");    elements.push_back("tellurium");    elements.push_back("iodine");    elements.push_back("xenon");    elements.push_back("caesium");    elements.push_back("barium");    elements.push_back("lanthanum");    elements.push_back("cerium");    elements.push_back("praseodymium");    elements.push_back("neodymium");    elements.push_back("promethium");    elements.push_back("samarium");    elements.push_back("europium");    elements.push_back("gadolinium");    elements.push_back("terbium");    elements.push_back("dysprosium");    elements.push_back("holmium");    elements.push_back("erbium");    elements.push_back("thulium");    elements.push_back("ytterbium");    elements.push_back("lutetium");    elements.push_back("hafnium");    elements.push_back("tantalum");    elements.push_back("tungsten");    elements.push_back("rhenium");    elements.push_back("osmium");    elements.push_back("iridium");    elements.push_back("platinum");    elements.push_back("gold");    elements.push_back("mercury");    elements.push_back("thallium");    elements.push_back("lead");    elements.push_back("bismuth");    elements.push_back("polonium");    elements.push_back("astatine");    elements.push_back("radon");    elements.push_back("francium");    elements.push_back("radium");    elements.push_back("actinium");    elements.push_back("thorium");    elements.push_back("protactinium");    elements.push_back("uranium");    elements.push_back("neptunium");    elements.push_back("plutonium");    elements.push_back("americium");    elements.push_back("curium");    elements.push_back("berkelium");    elements.push_back("californium");    elements.push_back("einsteinium");    elements.push_back("fermium");    elements.push_back("mendelevium");    elements.push_back("nobelium");    elements.push_back("lawrencium");    elements.push_back("rutherfordium");    elements.push_back("dubnium");    elements.push_back("seaborgium");    elements.push_back("bohrium");    elements.push_back("hassium");    elements.push_back("meitnerium");    elements.push_back("darmstadtium");    elements.push_back("roentgenium");    elements.push_back("copernicium");    return elements;}`

### Application Logic

`int main(int argc, char* argv[]){    std::vector<std::string> elements = make_elements_vector();    std::cout << "What element are you attempting to find? ";    std::string inputStr;    std::cin >> inputStr;    size_t minDistIndex = 0;    int minDist = INT_MAX;    for(size_t i=0; i<elements.size(); i++)    {        int dist = ld(elements[i], inputStr);        if(dist < minDist)        {            minDist = dist;            minDistIndex = i;        }    }    if(minDist == 0)    {        std::cout << "ELEMENT FOUND!" << std::endl;    }    else    {        std::string dym = "Did you mean " + elements[minDistIndex] + "?";        std::cout << dym << std::endl;    }    return 0;}`

### Search for “tillium”

This little demo works surprisingly well and suggestions are more-or-less inline with what you’d expect. There are obviously limitations as you think about applying this to other domains as language, context, etc. are not taken into consideration, but as a simple suggester it holds up pretty well and is perhaps a nice addition to a number of search methods in a variety of applications (the majority of which don’t seem to implement anything of the sort).

### firesync and the copy-delete-rename problem

An interesting problem with firesync popped up a while ago. I was syncing files on my laptop and a file didn’t get updated. Thinking this was odd, I tried to sync again, and got the same problem. So I looked at the file modification times of the files on the 2 computers and noticed the problem. What happened was…

1. I had a file (we’ll call it fileA) and made a copy of it (we’ll call the copy fileB).
2. I deleted a file (we’ll call it fileC)
3. I renamed fileB to the name of fileC (thus replacing fileC with fileB)

Unfortunately, when fileB was made, Windows set the file modification time of fileB to that of fileA and fileA had a modification time <= the modification time of fileC. So when firesync saw the file it looked like fileC didn’t need to be updated.

It’s a weird and complex little problem, but the good news it that when the file copy was done, Windows gave fileB a newer creation time. So it’s a somewhat easy fix that’ll be implemented in the next version of firesync.