In spite of apparent simplicity fast generation of circles and other curves has always been a subject of numerous theoretical and practical studies, finding applications, among others, in digital plotting, graphical display and numerical machine tool control, compare [12, 15, 21, 25]. Leaving aside the fundamental algoritm of Bresenham  and techniques based on spline functions , in this paper we focus on digital differential analyzers (DDA), see for instance [11, 16, 19]. In section 2 we present brief but exhaustive discussion of digital differential analyzers from a unified point of view, treating them as special one-step difference numerical schemes. Then, in section 3, we introduce and discuss a new DDA scheme, cheap and accurate, based on a two-step numerical scheme.
In terms of the natural arc parameter the equation of the circle of radius is , . The corresponding differential equation is
Considering various discretizations of (1) we can obtain any DDA algorithm.
2 One-step numerical approximations
In this section we present a unified approach to digital differential analysers representing them by one-step numerical schemes. Therefore they can be represented by a general matrix difference equation of the first order:
where depend on the time step . In order to approximate the system (1) these coefficients have to satisfy
In other words, , i.e.,
The qualitative asymptotic behaviour of the discrete solution is determined by the characteristic equation
and its roots (eigenvalues of the matrix) are given by
2.1 Logarithmic spirals
The most popular DDA schemes have the form
It is convenient to denote
i.e., , . Then,
Denoting , , we finally get
From computational point of view the polynomial form of
is preferred. DDA algorithms can be classified by the order of these polynomials. Below we present several schemes of this kind with the corresponding values of and (the exact algorithm should have and ).
First order simultaneous DDA algorithm
This is the conventional DDA method .
Second order simultaneous DDA algorithm
Third order simultaneous DDA algorithm
The best circular interpolator of the third order
We are going to show that this is the best spiral-like circular interpolator of the third order in.
Proof: We assume
Then we compute :
Equating to zero five first coefficients we obtain a system of 5 equations for 5 unknowns, which can be solved easily:
Then, the coefficient by is given by . The last equation of (16) yields two possibilities: either or . Corresponding coefficients by are given by and , respectively. Therefore the best accuracy is attained in the first case, which leads to (14) (the second case has another disadvantage: coefficients and are irrational).
Interestingly enough, considering -expansion of the function we get the same final results (although the starting system of equations is more complicated). The case yields while in the second case, , we get .
2.2 Elliptical deformations
Discrete points generated by (2) lie on an ellipse if and only if and , which is equivalent to
The conditions are known as Barkhausen criteria .
First order sequential DDA algorithm
2.3 Elliptical spirals
Two complex eigenvalues () are obtained when
and, for the generated circle deforms into an elliptical spiral. The following example can be found in .
Second order sequential DDA algorithm
One can easily verify, that . Therefore (19) is satisfied for .
2.4 Exactly circular interpolation
Digital differential analyzer of the form (2) yields an exact circle (up to round-off errors) if and only if , and . Then, the matrix can be parameterized by a single parameter
where depends on time step . In other words,
Implicit midpoint rule
In the case of linear equations the implicit midpoint rule can be represented in an explicit form. Indeed,
can be rewritten in a matrix form as follows
Thus we obtain another algortihm of the form (6), where
see . Note that and , so the implicit midpoint rule generates the circle exactly. However, this algorithm has a relatively high cost (multiplications in every step). Actually the scheme (27) was a motivation for the derivation of the Matsushiro algorithm (13), see . Indeed, (13) is the Taylor truncation of (27) up to the third order.
3 New fast circular interpolator
, applied to a general ordinary differential equation(where and is a given function) reads:
In the case of the circular interpolation and . Therefore we obtain:
compare . This is a two step method and we have to prescribe not only but also . As usual, we can rewrite this scheme in a one-step form increasing the dimension of the matrix problem, namely:
3.1 First integrals
In this section we are going to show that (29) preserves circular trajectories exactly provided that initial conditions are appropriately chosen.
If is a first integral of the system (29), then is also a first integral. What is more, in this case both first integrals are linearly dependent, namely
Proof: Using (29) we compute
If , then, obviously, .
Motivated by above results, we define:
For any satisfying (29) we have
Proof: A simple straightforward calculation.
All points generated by (29) lie exactly (up to round-off errors) on the circle provided that initial conditions satisfy
Proof: We easily verify that
3.2 Low computational cost
The scheme (29) has a very low computation cost, to be compared only with the simplest one-step scheme (10). Indeed, both (29) and (10) can be implemented in terms of two additions and two shifts. However, (29) produces an exact circle while the accuracy of (10) is rather limited.
There is only one more complicated computation, at the start, when initial conditions are are determined. However, taking into account the Taylor expansion
we see that even in that case one can use few shifts to produce . Indeed, let and . Then
where only integer coefficients are left.
The case (i.e., ) corresponds to the circle approximated by regular dodecagon (because then , , compare (36)). Assuming we may confine ourselves to the first three terms.
Characteristic equation reads
hence eigenvalues of are given by:
All the eigenvalues lie on the unit circle: for
. Corresponding eigenvectors are denoted by, , namely
( is given by (35) and the form of other eigenvectors is not important). The initial condition can always be represented as a linear combination of eigenvectors: , where we assume (compare Theorem 3.3) and , are small perturbations. Using (42) and (33) we obtain
Taking into account , we see that the initial perturbations remain unchanged during the evolution.
3.4 Exact circular interpolator
The scheme (29) produces exact circular trajectory but the period is modified. For instance, in the case of regular dodecagon () the total change of after 12 steps is while the exact value is, of course, . In general, the period of the circular interpolator (29) is given by
where is a given function of .
For any all points generated by (45) lie exactly (up to round-off errors) on the circle provided that initial conditions satisfy
Proof: It is enough to repeat the proof of Theorem 3.3 (with replaced by in appropriate places).
The period of the circular interpolator (45) is
If , then scheme
generates the circle exactly (up to round-off errors) preserving the period.
Proof: The scheme (45) will preserve the period exactly if and . Taking into account well known trigonometric identities
we conclude that scheme (45) is exact for .
The scheme (45) can serve as a starting point for deriving cheaper algorithms by taking -polynomials instead of . All these schemes preserve exactly the circular trajectory. The best approximation of the period are obtained for Taylor truncations. In particular, taking
we get the best scheme (45) of order 3. The multiplication by can be replaced by shifts because
and we can approximate by given by
where is chosen to assure the required accuracy.
We discussed DDA circular interpolators based on one-step numerical schemes. We have shown that the best cheap interpolator of the third order is given by (14), see Lemma 2.2. This algorithm can be expressed in terms of shifts and additions only, and it has quite good accuracy. Exact interpolators based on one-step methods are more expensive because multiplications in every step are required.
In section 3 we presented a family of circular interpolators, based on a two-step method, namely the explicit midpoint rule. All of them interpolate the circle exactly (up to round-off errors). The simplest interpolator, given by (29), is very cheap and is very good for graphical applications (when the period of the circular motion is irrelevant). The scheme (48) preserves exactly also the period but is more expensive. Taking an appropriate polynomial we can control the accuracy and computational cost of the scheme (45), obtaining interpolators suitable for different purposes.
Acknowledgement. The first author (J.L.C.) is partly supported by the National Science Centre (NCN) grant no. 2011/01/B/ST1/05137.
-  C.Baumgarten, G.Farin: Approximation of logarithmic spirals, Computer Aided Geometric Design 14 (6) (1997) 515-532.
-  C.A.Bergren: A simple algorithm for circular interpolation, Control Engineering 18 (9) (1971) 57-59.
-  J.E.Bresenham: A linear algorithm for incremental digital display of circular arcs, Commun. ACM 20 (2) (1977) 100-106.
-  I.Chami: A high precision digital differential analyzer for circular interpolation, Tishreen Univ. J. Stud. Sci. Res., Eng. Sci. Ser. 29 (1) (2007).
-  J.L.Cieśliński, B.Ratkiewicz: On simulations of the classical harmonic oscillator equation by difference equations, Adv. Difference Equ. 2006 (2006) 40171 (17pp).
-  J.L.Cieśliński: On the exact discretization of the classical harmonic oscillator equation, J. Difference Equ. Appl. 17 (2011) 1673-1694.
-  P.E.Danielsson: Incremental curve generation, IEEE Trans. Comput. C-19 (9) (1970) 783-793.
-  X.Y.Fan, Y.H.Guo, S.C.Li: New digital differential analyzer interpolation algorithm, Advanced Science Letters 6 (2012) 692-695.
-  E.Hairer, S.P.Nørsett, G.Wanner: Solving ordinary differential equations I. Nonstiff problems, Springer, Berlin 1987.
-  A.Iserles: A first course in the numerical analysis of differential equations, second edition, Cambridge Univ. Press 2009.
-  Y.Koren: Interpolator for a computer numerical control system, IEEE Trans. Comput. C-25 (1) (1976) 32-37.
-  F.S.Lim, Y.S.Wong, M.Rahman: Circular interpolators for numerical control: A comparison of the modified DDA techniques and an LSI interpolator, Computers in Industry 18 (1992) 41-52.
-  N.Matsushiro: A new digital differential analyzer for circle generation, IEICE Trans. Inf. Syst. E81-D (2) (1998) 239-242.
-  N.Matsushiro, I.Oyake: Method of and device for circle generation, US Patent No. 4,999,797 (March 1991).
-  P.G.McCrea, P.W.Baker: On digital differential analyzer (DDA) circle generation for computer graphics, IEEE Trans. Comput. C-24 (11) (1975) 1109-1110.
-  L.V.Moroz: A method for studying characteristics of digital sine-cosine generators, Information Extraction and Processing 36 (112) (2012) 84-90 [in Ukrainian].
-  W.M.Newman, R.F.Sproull: Principles of interactive computer graphics, McGraw-Hill 1979.
-  L.A.Piegl, W.Tiller: Circle approximation using integral B-splines, Computer-Aided Design 35 (2003) 601-607.
-  C.S.Turner: Recursive discrete-time sinusoidal oscillators, IEEE Signal Processing Magazine, May 2003.
-  W.P.Wang, C.Y.Wang: Difference method for generation of circular arcs and ellipses, Computer-Aided Design 21 (1) (1989) 33-37.
-  L.Yong-Kui: Algorithm for circle approximation and generation, Computer-Aided Design 25 (3) (1993) 169-171.