Today, I worked against numerical fuzz, and branch cuts. In particular I was working to construct Chebfuns for the coupler curve of a four-bar mechanism. According to Charles Wampler, a world expert on mechanism design via synthesis, I should have been able to parameterize the two components of a coupler curve by taking the positive and negative branches of the solution to a quadratic -- straight up quadratic formula. At first, it appeared he was wrong.

On constructing the chebfuns

To create the chebfuns for the coupler curve, we proceed in what is probably the typical way.

1. Write the loop equations for the mechanism.
2. Exploit the loop equations and complex conjugation of isotropic coordinates to get at a quadratic.
3. Solve the quadratic
4. Profit.

No problem, after a bit of mucking about in the weeds, I was able to get some decent pics, like this one:

See the problem? The colors are supposed to correspond to the distinct loops, not strange pieces of it. And, there's all this jitter in the pile of points around transitions between colors. Both of these are problematic. I will discuss the solutions below.

Correcting the colors

The colors come from the positive and negative branches of the square root when solving a quadratic. Remember the quadratic formula?

(-b ± \sqrt(b^2 - 4ac)) / (2a)

Ok, so when we are parameterizing a, b, and c by a variable t, and the variable t doesn't respect the branch cut for sqrt, garbage happens. Let me explain.

We are told as children that you can't take the square root of a negative number, which is a blatant lie. You can. You can also take the square root of a complex number, it works fine. But, the function deals with angles, and angles are periodic, and there are multiple roots of a number. So, we have to be clear about which root we are returning.

The Matlab documentation says that it will return square roots such that the angle is in the interval [-pi, pi]. But, I had set t to be in the interval from [0, 2pi]. So, my intervals were getting chopped up!

The fix was to change my chebfuns to use the interval [-pi pi].

Correcting the fuzz

The second problem to fix was the jitter and garbage points.

If you take the square root of a small real number, plus small error in a random complex direction, you will get unexpected results. This is exactly what was happening in my code. The solution was to wrap the discriminant of the quadratic in a call to real, to discard the imaginary noise.

Note to the scared reader: the discriminant is indeed a real number theoretically, so doing this is sound.

After adding the call, things look very nice.

Conclusion

My coupler curve plots are super nice now. Each piece is its own color, and I have an explicit parameterization of it in terms of the angle of the drive joint. What I do from here is another story.

Here is the plot after fixing it. Thanks for reading!