I am writing a program about geometry. I don’t want to get into details yet, I will do it when I have more to show.
While coding it was necessary to get the coordinates of the intersection of two circles. It is a deceptively simple problem, and without thinking much about it you can find several solutions on the Internet [SOF1] [SEX1] [WMW1].
But these solutions, sometimes start from geometrical reasoning for several cases, sometimes they are not general enough (considering a coordinate system that facilitates the solution for instance, what though ingenious do not solve the general problem without additional steps, in general not described ) or both [WMW1].
I decided to make my own solution and started as everyone analyzing the cases, apparently obvious: when circles are too far to cross each other, when they touch in one point (this is called osculating circles), when you cross at two points. Simple, right? But when specifying the conditions to enable writing the program other cases emerged: when circles have the same center but different radii, meaning that they don’t cross, when the center of one circle is inside another but they touch at one or two points.
I start wondering that this was turning into a mess of ifs, even before this part of the program was ready to compile [1]. Wouldn’t exist another case, or would my ifs had conditions that were more generic and encompassed other more specific conditions contained in them? In that case the program would never reach the specific conditions if the general ones were before them in code. As always when programming, something banal is not so banal if you look closely.
The Art of UNIX programming
I decided to use a thing I learned many time ago: that in programming frequently the formalization of a problem in solid mathematical terms makes it much simpler.
In the book The Art of UNIX Programming [RAY1] [2], Eric Raymond writes:
One subtle but powerful way to promote compactness in a design is to organize it around a strong core algorithm addressing a clear formal definition of the problem, avoiding heuristics and fudging.
“Formalization often clarifies a task spectacularly. It is not enough for a programmer to recognize that bits of his task fall within standard computer-science categories — a little depth-first search here and a quicksort there. The best results occur when the nub of the task can be formalized, and a clear model of the job at hand can be constructed. It is not necessary that ultimate users comprehend the model. The very existence of a unifying core will provide a comfortable feel, unencumbered with the why-in-hell-did-they-do-that moments that are so prevalent in using Swiss-army-knife programs.”
-- Doug McIlroy
This is an often-overlooked strength of the Unix tradition. Many of its most effective tools are thin wrappers around a direct translation of some single powerful algorithm.
Perhaps the clearest example of this is diff(1), the Unix tool for reporting differences between related files. This tool and its dual, patch(1), have become central to the network-distributed development style of modern Unix. A valuable property of diff is that it seldom surprises anyone. It doesn't have special cases or painful edge conditions, because it uses a simple, mathematically sound method of sequence comparison. This has consequences:
“By virtue of a mathematical model and a solid algorithm, Unix diff contrasts markedly with its imitators. First, the central engine is solid, small, and has never needed one line of maintenance. Second, the results are clear and consistent, unmarred by surprises where heuristics fail.”
-- Doug McIlroy
Thus, people who use diff can develop an intuitive feel for what it will do in any given situation without necessarily understanding the central algorithm perfectly. Other well-known examples of this special kind of clarity achieved through a strong central algorithm abound in Unix:
All three of these programs are so bug-free that their correct functioning is taken utterly for granted, and compact enough to fit easily in a programmer's hand. Only a part of these good qualities are due to the polishing that comes with a long service life and frequent use; most of it is that, having been constructed around a strong and provably correct algorithmic core, they never needed much polishing in the first place.
I decided to do exactly that and treat the simple problem of two circles intersection formal and exactly instead of thinking about each special case.
Geometry (according Descartes)
I wanted a comprehensive solution that assured me that I had all the cases. I preferred a Cartesian geometry approach from the beginning. For that is enough to write the equations of two circles in an arbitrary position, but if you do that you will get lost in a sea of letters when you try to do the math, that because you will have two coordinates of each center and will carry them all the time, and the relevant information is only the circles radii and the distance between the centers, everything else depends on coordinate choice.
You can use an algebraic computation program [3]. If you use Maxima and enter the general equations and tell it to solve them you won’t get anything [4], just an empty list. On the other side, if an adequate coordinate system is chosen the problem is greatly simplified [5]. It is enough to choose the origin as one of the centers and the $x$ axis on the straight line that connect the centers. If the original case is not that we can always convert one coordinate system into another.
With this coordinate system we have the following equations:
\[{x}^{2}+{y}^{2}={r}_{1}^{2}\]
\[{\left( x-d\right) }^{2}+{y}^{2}={r}_{2}^{2}\]
Maxima can solve this system easily [6] (although it take a few minutes) finding the intersection points coordinates:
\[x=\frac{{d}^{2}-{r}_{2}^{2}+{r}_{1}^{2}}{2\,d},y=-\frac{\sqrt{-{d}^{4}+2\,{r}_{2}^{2}\,{d}^{2}+2\,{r}_{1}^{2}\,{d}^{2}-{r}_{2}^{4}+2\,{r}_{1}^{2}\,{r}_{2}^{2}-{r}_{1}^{4}}}{2\,d}\]
\[x=\frac{{d}^{2}-{r}_{2}^{2}+{r}_{1}^{2}}{2\,d},y=\frac{\sqrt{-{d}^{4}+2\,{r}_{2}^{2}\,{d}^{2}+2\,{r}_{1}^{2}\,{d}^{2}-{r}_{2}^{4}+2\,{r}_{1}^{2}\,{r}_{2}^{2}-{r}_{1}^{4}}}{2\,d}\]
note that the solution is given only in terms of the radii and the distance between the centers, no coordinates of circles centers. But this still seems complicated. Lets call the term under radical the determinant and factor it, so we get:
\[det=-\left( d-{r}_{2}-{r}_{1}\right) \,\left( d-{r}_{2}+{r}_{1}\right) \,\left( d+{r}_{2}-{r}_{1}\right) \,\left( d+{r}_{2}+{r}_{1}\right) \]
This is interesting because it reduces the number of multiplications and sums [7], and because it gives us the conditions in a natural way. Simplifying and using the notation ${s}_{r}={r}_{2}+{r}_{1}$ and ${d}_{r}={r}_{2}-{r}_{1}$, we get:
\[det=-\left( {d}_{r}-d\right) \,\left( {d}_{r}+d\right) \,\left( {s}_{r}-d\right) \,\left( {s}_{r}+d\right) \]
Expressing the solution this way, we have:
\[x=\frac{dr\,sr-{d}^{2}}{2\,d},y=-\frac{\sqrt{det}}{2\,d}\]
\[x=\frac{dr\,sr-{d}^{2}}{2\,d},y=\frac{\sqrt{det}}{2\,d}\]
This will gives us the conditions from the determinant, because as it is under the radical must be positive or there is no solution. As it is a four term product the signal depends on the terms, and as ${s}_{r}$ is always positive (as ${r}_{1}$ and ${r}_{2}$ are always positive) and ${d}$ is also always positive because it is a distance.
As for ${d}_{r}$ it is ${r}_{1}-{r}_{2}$ but this seem arbitrary as ${r}_{2}$ can be either the radius that is bigger or smaller than ${r}_{1}$, it is only a question of what name I give each radius. So it could be positive or negative, but ${d}_{r}$ appears in only two terms $\left( {d}_{r}-d\right) \,\left( {d}_{r}+d\right) $ and expanding we have ${d}_{r}^{2}-{d}^{2}$, that is, the ${d}_{r}$ signal do not change the result (because it appears in two terms and gets squared in the expansion), so we can take ${d}_{r}$ without sign, simplifying our reasoning.
Thus the sign of ${det}$ depends on the terms that are subtractions that can be positive or negative. $\left( {d}_{r}-d\right) $ and $\left( {s}_{r}-d\right)$. Then we have 4 cases:
| ${s}_{r}>d$ | ${s}_{r}<d$ | |
| ${d}_{r}>d$ | $det<0$ | $det>0$ |
| ${d}_{r}>d$ | $det>0$ | $det<0$ |
It is now much simpler to specify the cases being sure they are comprehensive. Lets see the meaning of each case:
- ${s}_{r}<d,{d}_{r}>d$: this case has a solution, but given the limitations on the values of ${r}_{1}$, ${r}_{2}$ and ${d}$ it never happens, if we have ${r}_{2}+{r}_{1}<d$ and ${r}_{2}-{r}_{1}>d$ we would have ${r}_{2}+{r}_{1}<{d}<{r}_{2}-{r}_{1}$, or ${r}_{2}+{r}_{1}<{r}_{2}-{r}_{1}$ what is impossible for ${r}_{1}$ and ${r}_{2}$ positive;
- ${s}_{r}>d,{d}_{r}>d$: this case is possible, and has no solution, because one circle in inside the other and the smaller do not touch the bigger because it is too small for they to cross;
- ${s}_{r}<d,{d}_{r}<d$: this case is possible, and has no solution because the circles have their centers too far away for them to touch;
- ${s}_{r}>d,{d}_{r}<d$: this case is possible, and it has solution, the circles touch in two points.
We still have one case to analyze:
- $det=0$: we have only one solution, the circles osculate, touching in only one point. For that it is enough that one of the terms be 0 (as they are being multiplied), meaning ${s}_{r}=d$ or ${d}_{r}=d$ as the terms with sums cannot be 0. If ${s}_{r}=d$ we have circles osculating from outside,and if ${d}_{r}=d$ we have circles osculating from inside;
Besides the values of $det$ we still have two cases to analyze:
- $d=0, {r}_{1}\not={r}_{2}$: as $d$ is on the solution divisor, we have a division by 0 and there is no solution. This represents two concentric circles one having its radius smaller than the other;
- $d=0, {r}_{1}={r}_{2}$: now in the solution expression besides 0 on denominator we have 0 on the numerator. Lets understand what this means, the equation of the two circles is the same, meaning we have only one equation and the solution for the intersection is the set of all points that satisfy the equations, that is the entire circle. It is not the case that it has no solution, it has an infinite number of them. In this case in mathematics we say that the system is compatible, but undetermined. This represents two concentric circles with the same radius.
But it is not over yet, we need to reconvert the coordinates for the system that we had before simplifying the coordinate system. This will give us the general response.
We can consider the transformation between two orthonormal coordinate systems as a translation followed by a rotation. The translation shifts the origin position and rotation turns the axes until they are superimposed. Translation is only a sum (or subtraction), rotation is a little more complex and generally is expressed by a rotation matrix applied to coordinates [8].
\[\begin{bmatrix}x\prime\cr y\prime\end{bmatrix}=\begin{bmatrix}\mathrm{cos}\left( \theta\right) & \mathrm{sin}\left( \theta\right) \cr -\mathrm{sin}\left( \theta\right) & \mathrm{cos}\left( \theta\right) \end{bmatrix} . \left( \begin{bmatrix}x\cr y\end{bmatrix}-\begin{bmatrix}{x}_{0}\cr {y}_{0}\end{bmatrix}\right) \]
Where $(x\prime, y\prime)$ are the coordinates in the simplified coordinate system and $(x,y)$ are the coordinates of the original more general system. As we have the result in terms of $(x\prime, y\prime)$ we need invert this expression [9]. This can be done in several ways, what gives us the following result:
\[\begin{bmatrix}\mathrm{cos}\left( \theta\right) & -\mathrm{sin}\left( \theta\right) \cr \mathrm{sin}\left( \theta\right) & \mathrm{cos}\left( \theta\right) \end{bmatrix} . \begin{bmatrix}x\prime\cr y\prime\end{bmatrix}+\begin{bmatrix}{x}_{0}\cr {y}_{0}\end{bmatrix}=\begin{bmatrix}x\cr y\end{bmatrix}\]
A lifestyle
Now we have the generalized solution. Not just that, with the analysis we made the problem was elegantly solved, we know the solution is comprehensive and this has suggested to us an efficient and expressive way to implement it.
What I want to convey here is less the solution of the two circle intersection than the way of approaching the problem, that is to create a formal and solid model and not to create a patchwork of particular of extreme cases.
But it is incredibly complicated to convince other people of that. With time I realized that the reluctance of people in accepting more formal processes was not the result of a healthy skepticism, but it was the result of unfamiliarity of these structures and algorithms. They saw them as theoretical abstractions without contact with reality and preferred a concrete approach, through a stumbling heuristics.
I have a belief, very different from the current fads. The belief that software must be built like a cathedral, aiming at something bigger. Some will argument that my analogy is flawed and that cathedrals were not made for money, to those I remember that many fundamental software either.
Most will say I am crazy saying that time is money, that we need to have quick results so we can make money. I will say that a little less speed and a little more consideration will avoid money waste. Avoiding a vicious circle of corrections of errors or inadequacies that need to be made even quicker than the original product, and this urgency creates new problems in increasingly faster and decreasingly effective cycles, where corrections introduce problems that frequently are even worse. And all could be avoided. If you bother with invasive update systems like the ones from Microsoft and Adobe, maybe you understand what I mean.
What do I propose? It is simple, if software is not pleasant to make, if it is not seem elegant, beautiful and robust at the end then I don’t want to have anything to do with it. This works for me and kept me employed and productive for decades. If you think that to make a career in programming you cannot think this way, you are failing to do something that could be very pleasant. Basically I propose this as a lifestyle.
The art of the fugue
Decades ago I saw a Brazilian singer/comedian interview where he talked about a technical issue said that that was not art. I believe that a piece of technology or science is art, an processor project, a FPGA, the proof of a mathematical theorem, or a computer program can be art.
But there are some differences, one of them is that these things do not have the original goal of being appreciated as art. But neither a cathedrals, they are built for practical end, but that does not stop them of being made as a work of art for reaching that. In this particular point a computer program can be built this way.
As artists uses art to express themselves, programmers also can use programs to express themselves and when this is done we have real art, intentionally.
Another difference is that a piece of technology or science is a form of art that cannot be properly appreciated easily, it requires effort and knowledge.
Let me draw a parallel: Bach in the end of his life was seen by many, like his son Carl Phillip Emmanuel as a relic, his polyphonic music full of canons and fugues was considered out of fashion. It was felt that it was too sophisticated for the public to appreciate. His employer the priest where he was kapellmeister thought the same, that the sophistication of his music instead of praising God, distracted the congregation from him.
Carl Phillip Emanuel called him to the court of Frederick the Great, of Prussia, where he was treated with a mix of reverence and derision, as they wanted to prove that his greatness though acknowledged was obsolete [GAI1].
Today however, when someone mentions Bach, everyone thinks first in Johann Sebastian and not in Carl Philipp Emanuel or Johann Christian, to the extent that I can refer to him as Bach without fear of being ambiguous.
It was true that Bach music was sophisticated, it was true that the public do not appreciated it in its entirety, but this was a problem with the public and not with Bach. When I was a teenager I liked Bach, but I even liked of the compositions of Carl Phillip Emmanuel a little more, so one day I read a book: Gödel, Escher, Bach: an Eternal Golden Braid [HOF1], where it was taught me among other things the details of polyphony, and the multiple voices separated and united again when I heard a fugue and Bach never sounded the same again.
With time this kind of experience repeated countless times and I realized that each time I learned about something that I didn’t like (or didn’t like so much at first), this made me take a fresh look at it that led me to appreciate.
I came to the conclusion that the problem with something that I didn’t like was generally in my blindness and ignorance and not in the thing itself, that is how I came to appreciate opera for instance, but this was a felling that was not different of the feelings I had when I learned to read numbers of more than two digits, or when I could design a circuit with a relay that memorized its state, in an a mechanical analogue to a flip-flop, or when I understood binary numbers, or when I understood why the integral that gives Fourier series coefficients works, or when I discovered that the Taylor series expansion of the exponential gives a sum of sines and cosines if it is applied to an imaginary angle.
These are precious moments of my life, so precious as the smell of churros at the zoo or the memory of my grandmother. These are the moments that show me that there is no shame in ignorance, there is only shame in believing that the ignorance of something is in any way a virtue, losing an opportunity to learn. It is wonderful that there is always something to learn, or life would be sad.
LeoNerd
The feeling of rapture, that I discovered with these things was no different at all from the feeling of rapture with Leonardo pictures (I loved in particular the two versions of the Virgin of the Rocks), and it is interesting that he didn’t differentiate in importance his art and his technique not only in painting or sculpture, but in the project of devices like the ornithopter and the parachute. When Bill Gates spent a fortune acquiring some Leonardo manuscripts [WIKI1] I understood well why, we can say that he was a nerd and to him everything was art, unlike that Brazilian singer/comedian thought.
Notes
[1]. This is a code smell called complex conditionals. This is always an indication of an heuristic approach and future problems.
[2]. The electronic version of this book use the Creative Commons license, reserving the right to publish on paper to the publisher of the book. So it is easy to find online versions, like those in [RAY1] and [RAY2].
[3]. This is generally called CAS (computer algebra system). There are several of them, like: Maxima (that is the GPL version of MIT Macsyma), Maple, Mathematica, Reduce (that was proprietary but today is open with a BSD license), some electronic calculators like the TI Nspire CAS (that used a version of the old Derive that was a program small and light enough so that it could be implemented in a electronic calculator), and the HP50g.
[4]. On Maxima we would do that to solve:
circle1:(x-x[1])^2+(y-y[1])^2=r[1]^2;
circle2:(x-x[2])^2+(y-y[2])^2=r[2]^2;
solve([circle1,circle2],[x,y]);
[5]. Yes, I know, it was exactly what I criticized above, but I will show how to extend to the general case later.
[6]. To solve the simplified system on Maxima:
normalized_circle1:x^2+y^2=r[1]^2;
normalized_circle2:(x-d)^2+y^2=r[2]^2;
solve([normalized_circle1,normalized_circle2],[x,y]);
[7]. This similar to Horner’s method of polynomial evaluation, where you can calculate the value of a polynomial of degree $n$ with only $n$ sums and $n$ multiplications.
[8]. To understand why a matrix you need to know a little of linear algebra. Basically we apply a transformation to each point where each coordinate on the new system depends on its original coordinates multiplied by two different factors, one for x and another for y. The coordinate value is transformed like a straight line equation, where if I walk in the direction of one of the original coordinate axes I move a little in each axis of the new coordinate system, the more inclined is the line I move in relation to the new axis the less I move in the new axis. This proportionality factors are the sines and cosines of the inclination between the axes. There are a countless ways to understand that, by projection, by scalar product, etc. You have much to gain in insight exploring these ways, but this is out of the scope of this post.
[9]. You can multiply both sides by the inverse of the rotation matrix. It is trivial to invert a 2 by 2 matrix (yes, linear algebra again, as my teacher used to say “1 over the determinant plus the transposed cofactor matrix”). But we can gain in insight and insight and admiration for the consistency beauty exploring other ways. To invert the expression also means a rotation in the other direction. As sine is an odd function (meaning $sin(\theta)=-sin(-\theta)$) and cosine is an even function (meaning $cos(\theta)=cos(-\theta)$) simply substituting $\theta$ by $=-\theta$ in the rotation matrix we have:
\[\begin{bmatrix}\mathrm{cos}\left( -\theta\right) & \mathrm{sin}\left( -\theta\right) \cr -\mathrm{sin}\left( -\theta\right) & \mathrm{cos}\left( -\theta\right) \end{bmatrix}=\begin{bmatrix}\mathrm{cos}\left( \theta\right) & -\mathrm{sin}\left( \theta\right) \cr \mathrm{sin}\left( \theta\right) & \mathrm{cos}\left( \theta\right) \end{bmatrix}\]
that is also the inverse. Note that the inverse is equal of the to the transposed (a matrix like that is called orthogonal).
References
[SOF1]: python - Circle-circle intersection points - Stack Overflow. Available at: <http://stackoverflow.com/questions/3349125/circle-circle-intersection-points>. Accessed at: 21 MAI 2014.
[SEX1]: intersection - How can I find the points at which two circles intersect? - Mathematics Stack Exchange. Available at: <http://math.stackexchange.com/questions/256100/how-can-i-find-the-points-at-which-two-circles-intersect>. Accessed at: 21 MAI 2014.
[WMW1]: Circle-Circle Intersection -- from Wolfram MathWorld. Available at: <http://mathworld.wolfram.com/Circle-CircleIntersection.html>. Accessed at: 21 MAI 2014.
[RAY1]: RAYMOND, Eric S. The Art of UNIX Programming. Boston: Addison-Wesley, 2003. p. 111, cap. 4.
[RAY2]: The Art of UNIX Programming: Compactness and Orthogonality: Compactness and the Strong Single Center. Available at: <http://www.faqs.org/docs/artu/ch04s02.html#id2894555>. Accessed at: 26 MAI 2014.
[RAY3]: The Art of Unix Programming : Eric Steven Raymond : Free Download & Streaming : Internet Archive. Available at: <https://archive.org/details/ost-computer-science-the_art_of_unix_programming-1>. Accessed at: 26 MAI 2014.
[GAI1]: GAINES, James R. Evening in the Palace of Reason: Bach Meets Frederick the Great in the Age of Enlightenment. New York: Harper Perennial, 2006.
[HOF1]: HOFSTADTER, Douglas R. Gödel, Escher, Bach: an eternal Golden Braid. London: Penguin Books, 1980.
[WIKI1]: Codex Leicester. In: Wikipedia. Available at: <http://en.wikipedia.org/wiki/Codex_Leicester>. Accessed at: 02 JUN 2014.
Comments