## 1 Introduction

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 [7] and techniques based on spline functions [22], 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

(1) |

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:

(2) |

where depend on the time step . In order to approximate the system (1) these coefficients have to satisfy

(3) |

In other words, , i.e.,

The qualitative asymptotic behaviour of the discrete solution is determined by the characteristic equation

(4) |

and its roots (eigenvalues of the matrix

) are given by(5) |

where .

### 2.1 Logarithmic spirals

The most popular DDA schemes have the form

(6) |

It is convenient to denote

(7) |

i.e., , . Then,

Denoting , , we finally get

(8) |

###### Corollary 2.1.

From computational point of view the polynomial form of

is preferred. DDA algorithms can be classified by the order of these polynomials

[16]. 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

(10) |

This is the conventional DDA method [16].

#### Second order simultaneous DDA algorithm

#### Third order simultaneous DDA algorithm

(12) |

compare [6].

#### Matsushiro’s analyzer

#### The best circular interpolator of the third order

Interestingly enough, the Matsushiro algoritm can be easily improved by changing the last term, see [8, 20]:

(14) |

We are going to show that this is the best spiral-like circular interpolator of the third order in

.###### Lemma 2.2.

Proof: We assume

(15) |

Then we compute :

Equating to zero five first coefficients we obtain a system of 5 equations for 5 unknowns, which can be solved easily:

(16) |

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

(17) |

The conditions are known as Barkhausen criteria [23].

#### First order sequential DDA algorithm

### 2.3 Elliptical spirals

Two complex eigenvalues () are obtained when

(19) |

and, for the generated circle deforms into an elliptical spiral. The following example can be found in [16].

#### Second order sequential DDA algorithm

(20) |

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

(22) |

where depends on time step . In other words,

(23) |

where . Taking we obtain the exact numerical scheme (compare [21, 24])

(24) |

Exact numerical schemes can be constructed for any matrix differential equation of the first order, see [9, 10].

The Taylor expansion of (24) up to the third order yields the scheme (12), compare [6]. We recall that this is not the best scheme of the third order, see Lemma 2.2.

#### Implicit midpoint rule

In the case of linear equations the implicit midpoint rule can be represented in an explicit form. Indeed,

(25) |

can be rewritten in a matrix form as follows

(26) |

Thus we obtain another algortihm of the form (6), where

(27) |

see [15]. 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 [17]. Indeed, (13) is the Taylor truncation of (27) up to the third order.

## 3 New fast circular interpolator

The main result of our paper consists in applying a two step method to the circular interpolation. The explicit midpoint rule (see [14]), known also as 2-step Nyström method [13]

, applied to a general ordinary differential equation

(where and is a given function) reads:(28) |

In the case of the circular interpolation and . Therefore we obtain:

(29) |

compare [20]. 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:

(30) |

### 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.

###### Lemma 3.1.

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

(31) |

Motivated by above results, we define:

(32) |

###### Lemma 3.2.

For any satisfying (29) we have

(33) |

Proof: A simple straightforward calculation.

###### Theorem 3.3.

All points generated by (29) lie exactly (up to round-off errors) on the circle provided that initial conditions satisfy

(34) |

### 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

(37) |

we see that even in that case one can use few shifts to produce . Indeed, let and . Then

(38) |

and

(39) |

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.

### 3.3 Stability

Characteristic equation reads

(40) |

hence eigenvalues of are given by:

(41) |

All the eigenvalues lie on the unit circle: for

. Corresponding eigenvectors are denoted by

, , namely(42) |

( 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

(43) |

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

(44) |

In order to preserve exactly not only the trajectory but also the period of the circular motion, we can use a nonstandard modification of (29), following the approach of [9, 10]:

(45) |

where is a given function of .

###### Theorem 3.5.

For any all points generated by (45) lie exactly (up to round-off errors) on the circle provided that initial conditions satisfy

(46) |

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

(47) |

###### Theorem 3.6.

If , then scheme

(48) |

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

(49) |

we conclude that scheme (45) is exact for .

The exact interpolator (48) uses two multiplication by at every step, hence it is relatively expensive (but two times cheaper than its one-step conterpart (24)).

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

(50) |

we get the best scheme (45) of order 3. The multiplication by can be replaced by shifts because

(51) |

and we can approximate by given by

(52) |

where is chosen to assure the required accuracy.

## 4 Conclusions

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.

## References

- [1]
- [2]
- [3]
- [4]
- [5] C.Baumgarten, G.Farin: Approximation of logarithmic spirals, Computer Aided Geometric Design 14 (6) (1997) 515-532.
- [6] C.A.Bergren: A simple algorithm for circular interpolation, Control Engineering 18 (9) (1971) 57-59.
- [7] J.E.Bresenham: A linear algorithm for incremental digital display of circular arcs, Commun. ACM 20 (2) (1977) 100-106.
- [8] I.Chami: A high precision digital differential analyzer for circular interpolation, Tishreen Univ. J. Stud. Sci. Res., Eng. Sci. Ser. 29 (1) (2007).
- [9] 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).
- [10] J.L.Cieśliński: On the exact discretization of the classical harmonic oscillator equation, J. Difference Equ. Appl. 17 (2011) 1673-1694.
- [11] P.E.Danielsson: Incremental curve generation, IEEE Trans. Comput. C-19 (9) (1970) 783-793.
- [12] X.Y.Fan, Y.H.Guo, S.C.Li: New digital differential analyzer interpolation algorithm, Advanced Science Letters 6 (2012) 692-695.
- [13] E.Hairer, S.P.Nørsett, G.Wanner: Solving ordinary differential equations I. Nonstiff problems, Springer, Berlin 1987.
- [14] A.Iserles: A first course in the numerical analysis of differential equations, second edition, Cambridge Univ. Press 2009.
- [15] Y.Koren: Interpolator for a computer numerical control system, IEEE Trans. Comput. C-25 (1) (1976) 32-37.
- [16] 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.
- [17] N.Matsushiro: A new digital differential analyzer for circle generation, IEICE Trans. Inf. Syst. E81-D (2) (1998) 239-242.
- [18] N.Matsushiro, I.Oyake: Method of and device for circle generation, US Patent No. 4,999,797 (March 1991).
- [19] 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.
- [20] L.V.Moroz: A method for studying characteristics of digital sine-cosine generators, Information Extraction and Processing 36 (112) (2012) 84-90 [in Ukrainian].
- [21] W.M.Newman, R.F.Sproull: Principles of interactive computer graphics, McGraw-Hill 1979.
- [22] L.A.Piegl, W.Tiller: Circle approximation using integral B-splines, Computer-Aided Design 35 (2003) 601-607.
- [23] C.S.Turner: Recursive discrete-time sinusoidal oscillators, IEEE Signal Processing Magazine, May 2003.
- [24] W.P.Wang, C.Y.Wang: Difference method for generation of circular arcs and ellipses, Computer-Aided Design 21 (1) (1989) 33-37.
- [25] L.Yong-Kui: Algorithm for circle approximation and generation, Computer-Aided Design 25 (3) (1993) 169-171.

Comments

There are no comments yet.