## An error with area models

Wednesday, November 28th, 2007 at 2:10 pm

The subset of the 1-dimensional real line is an important model in computational geometry because it can be defined accurately in terms of an even-sized array of floats, and because it can be used as a building block for more sophisticated geometric structures.

The diagram below shows illustrates such a subset, A which happens to have 3 components and could be expressed by the array:

```   size(A) = 6;
A[0] = 0.5;  A[1] = 1.7;  A[2] = 2.9;
A[3] = 4.1;  A[4] = 5.8;  A[5] = 6.0;
```

which means it lies inside the interval [0.0, 6.0].

For convenience, subsets are assumed to be closed, which means they contain the endpoints of their connected segments. A contains(x) function could be written as:

```bool contains(A, x)
{
for (i = 0; i < size(A); i += 2)
{
if (A[i] <= x) and (x <= A[i + 1])
return True;
}
return False;
}
```

It's easy to calculate the compliment A' of this subset within an interval using the following routine:

```complement(A, lo, hi)
{
if (size(A) == 0) or (A[size(A) - 1] != hi)
A.append(hi);
else
A.eraseAt(size[A] - 1);

if (size(A) == 0) or (A[0] != lo)
A.insertAt(0, lo);
else
A.eraseAt(0);
}
```

The (size(A) == 0) condition is superfluous, because after the first paragraph size(A) will not be even, so cannot be 0.

(There’s no generally accepted pseudocode for array handling, but knowing what needs to be done should make it possible to understand what it means. The statements len(A), count(A), A.size(), and size(A) all mean the same thing.)

This function is not strictly the complement set, because if A' (drawn on the lower half of the above diagram) is the value of A after this function has been executed, then its intersection with A is non-empty because it contains five boundary points.

There are a lot of other problems with this representation, which I once wrote up in a sixty page chapter of a book I began on computational geometry before I got the current work. But one of the good things about it is it’s a self-inverting function. That is, if you apply compliment() twice, you come back to the value you started out with.

I have constructed a representation of a 2-dimensional subset of the plane using two arrays of these 1-dimensional subsets by arranging them in a weave, like so:

The area is broken down into these square cells between the 1-dimensional lines with the endpoints of the subsets connected up in the obvious manner. As you can see, this approximates the area bounded by the polygonal model (the triangle drawn with a green dotted line), which in this case has been placed in the most inconvenient way possible — one of its edges has been perfectly aligned with a horizontal 1-dimensional subset.

As you would expect, it’s pretty rare for this to happen; if the spacing between the lines is a random-looking value, like 1.928374746, then the odds of hitting one of the edges of the polygon are pretty low. However, when it does, things shouldn’t go particularly wrong…

Unless you try making the complement of the area.

There is an important consistency requirement with regards to these area models, which in my code are called “weaves” (the 1-dimensional subsets are known as “fibres”). It is that the crossing points are consistent. That is, for any vertical fibre A aligned with value x = a, and any horizontal fibre B aligned with the value y = b, then:

```    contains(A, b) == contains(B, a)
```

This condition is satisfied by the diagram above. Up till now I had implicitly found the complement of the area by applying complement() to each fibre in turn. This worked nearly always, and I had never suspected that there was a problem. Unfortunately, when you do it for such an instance where endpoints of subset segments lie exactly on a cross-fibre, the consistency requirement fails, and the rest of the algorithm which depends on this requirement being true breaks down in catastrophic and unpredictable ways.

Some of these failures were recovered from by additional code later on. However, once we began running our fancy algorithms on multi-processor machines, these fixes stopped working. I am therefore in a bit of a pickle.

Of course, the quick fix would be to ensure that there always is a good fudge factor to prevent the polygon (or any of its offsets) becoming perfectly aligned with one of the fibres. Multiplying the weave spacing by 0.999928376 would almost certainly make any given bug instantly disappear, whilst leaving the underlying design flaw in place.

For most people this would the best course of action, because they didn’t write the algorithm, they don’t understand it, more invasive fixes are likely to make things worse, and it’s not their code anyway. If they’re a very good programmer, they will have moved on up the company before it reappears. Workload tends to be pretty constant all the time, so there is nothing to be gained from fixing a very marginal problem absolutely given all the costs it entails.

That’s assuming you’re sure it is a marginal problem, and not just hoping it is. However, I do note that it is unusual for a programmer, unless he is young and naive, to believe he has a stake in the software to the extent that he expects that if it goes wrong in 20 years time it will still be his problem. It might not be true, but if he believes it, his results will be different.

So I better stop blogging and continue fixing this now-known issue now, having discussed it with the world (even if no one in the world understands it; it helps to talk it through). It’s not going very easily; it’s going to get slower before it gets faster because have to leave some duplicate calculations floating around before I work out how to avoid them. I hope everyone understands.

### 1 Comment

• 1. Freesteel » Blog Ar&hellip replies at 5th December 2007, 10:03 am :

[…] Freesteel Two CAM programmers on the loose

« An error with area models

That partial book
December 5t […]