Cycle detection

Last updated

In computer science, cycle detection or cycle finding is the algorithmic problem of finding a cycle in a sequence of iterated function values.

Contents

For any function f that maps a finite set S to itself, and any initial value x0 in S, the sequence of iterated function values

must eventually use the same value twice: there must be some pair of distinct indices i and j such that xi = xj. Once this happens, the sequence must continue periodically, by repeating the same sequence of values from xi to xj 1. Cycle detection is the problem of finding i and j, given f and x0.

Several algorithms are known for finding cycles quickly and with little memory. Robert W. Floyd's tortoise and hare algorithm moves two pointers at different speeds through the sequence of values until they both point to equal values. Alternatively, Brent's algorithm is based on the idea of exponential search. Both Floyd's and Brent's algorithms use only a constant number of memory cells, and take a number of function evaluations that is proportional to the distance from the start of the sequence to the first repetition. Several other algorithms trade off larger amounts of memory for fewer function evaluations.

The applications of cycle detection include testing the quality of pseudorandom number generators and cryptographic hash functions, computational number theory algorithms, detection of infinite loops in computer programs and periodic configurations in cellular automata, automated shape analysis of linked list data structures, and detection of deadlocks for transactions management in DBMS.

Example

This function defines the cycles {4} and {1, 6, 3}. Functional graph.svg
This function defines the cycles {4} and {1, 6, 3}.

The figure shows a function f that maps the set S = {0,1,2,3,4,5,6,7,8} to itself. If one starts from x0 = 2 and repeatedly applies f, one sees the sequence of values

2, 0, 6, 3, 1, 6, 3, 1, 6, 3, 1, ....

The cycle in this value sequence is 6, 3, 1.

Definitions

Let S be any finite set, f be any function from S to itself, and x0 be any element of S. For any i > 0, let xi = f(xi 1). Let μ be the smallest index such that the value xμ reappears infinitely often within the sequence of values xi, and let λ (the loop length) be the smallest positive integer such that xμ = xλ + μ. The cycle detection problem is the task of finding λ and μ. [1]

One can view the same problem graph-theoretically, by constructing a functional graph (that is, a directed graph in which each vertex has a single outgoing edge) the vertices of which are the elements of S and the edges of which map an element to the corresponding function value, as shown in the figure. The set of vertices reachable from starting vertex x0 form a subgraph with a shape resembling the Greek letter rho (ρ): a path of length μ from x0 to a cycle of λ vertices. [2]

Computer representation

Generally, f will not be specified as a table of values, the way it is shown in the figure above. Rather, a cycle detection algorithm may be given access either to the sequence of values xi, or to a subroutine for calculating f. The task is to find λ and μ while examining as few values from the sequence or performing as few subroutine calls as possible. Typically, also, the space complexity of an algorithm for the cycle detection problem is of importance: we wish to solve the problem while using an amount of memory significantly smaller than it would take to store the entire sequence.

In some applications, and in particular in Pollard's rho algorithm for integer factorization, the algorithm has much more limited access to S and to f. In Pollard's rho algorithm, for instance, S is the set of integers modulo an unknown prime factor of the number to be factorized, so even the size of S is unknown to the algorithm. To allow cycle detection algorithms to be used with such limited knowledge, they may be designed based on the following capabilities. Initially, the algorithm is assumed to have in its memory an object representing a pointer to the starting value x0. At any step, it may perform one of three actions: it may copy any pointer it has to another object in memory, it may apply f and replace any of its pointers by a pointer to the next object in the sequence, or it may apply a subroutine for determining whether two of its pointers represent equal values in the sequence. The equality test action may involve some nontrivial computation: for instance, in Pollard's rho algorithm, it is implemented by testing whether the difference between two stored values has a nontrivial greatest common divisor with the number to be factored. [2] In this context, by analogy to the pointer machine model of computation, an algorithm that only uses pointer copying, advancement within the sequence, and equality tests may be called a pointer algorithm.

Algorithms

If the input is given as a subroutine for calculating f, the cycle detection problem may be trivially solved using only λ + μ function applications, simply by computing the sequence of values xi and using a data structure such as a hash table to store these values and test whether each subsequent value has already been stored. However, the space complexity of this algorithm is proportional to λ + μ, unnecessarily large. Additionally, to implement this method as a pointer algorithm would require applying the equality test to each pair of values, resulting in quadratic time overall. Thus, research in this area has concentrated on two goals: using less space than this naive algorithm, and finding pointer algorithms that use fewer equality tests.

Floyd's tortoise and hare

Floyd's "tortoise and hare" cycle detection algorithm, applied to the sequence 2, 0, 6, 3, 1, 6, 3, 1, ... Tortoise and hare algorithm.svg
Floyd's "tortoise and hare" cycle detection algorithm, applied to the sequence 2, 0, 6, 3, 1, 6, 3, 1, ...

Floyd's cycle-finding algorithm is a pointer algorithm that uses only two pointers, which move through the sequence at different speeds. It is also called the "tortoise and the hare algorithm", alluding to Aesop's fable of The Tortoise and the Hare.

The algorithm is named after Robert W. Floyd, who was credited with its invention by Donald Knuth. [3] [4] However, the algorithm does not appear in Floyd's published work, and this may be a misattribution: Floyd describes algorithms for listing all simple cycles in a directed graph in a 1967 paper, [5] but this paper does not describe the cycle-finding problem in functional graphs that is the subject of this article. In fact, Knuth's statement (in 1969), attributing it to Floyd, without citation, is the first known appearance in print, and it thus may be a folk theorem, not attributable to a single individual. [6]

The key insight in the algorithm is as follows. If there is a cycle, then, for any integers iμ and k ≥ 0, xi = xi + , where λ is the length of the loop to be found, μ is the index of the first element of the cycle, and k is a whole integer representing the number of loops. Based on this, it can then be shown that i = μ for some k if and only if xi = x2i (if xi = x2i in the cycle, then there exists some k such that 2i = i + , which implies that i = ; and if there are some i and k such that i = , then 2i = i + and x2i = xi + ). Thus, the algorithm only needs to check for repeated values of this special form, one twice as far from the start of the sequence as the other, to find a period ν of a repetition that is a multiple of λ. Once ν is found, the algorithm retraces the sequence from its start to find the first repeated value xμ in the sequence, using the fact that λ divides ν and therefore that xμ = xμ + v. Finally, once the value of μ is known it is trivial to find the length λ of the shortest repeating cycle, by searching for the first position μ + λ for which xμ + λ = xμ.

The algorithm thus maintains two pointers into the given sequence, one (the tortoise) at xi, and the other (the hare) at x2i. At each step of the algorithm, it increases i by one, moving the tortoise one step forward and the hare two steps forward in the sequence, and then compares the sequence values at these two pointers. The smallest value of i> 0 for which the tortoise and hare point to equal values is the desired value ν.

The following Python code shows how this idea may be implemented as an algorithm.

deffloyd(f,x0)->(int,int):"""Floyd's cycle detection algorithm."""# Main phase of algorithm: finding a repetition x_i = x_2i.# The hare moves twice as quickly as the tortoise and# the distance between them increases by 1 at each step.# Eventually they will both be inside the cycle and then,# at some point, the distance between them will be# divisible by the period λ.tortoise=f(x0)# f(x0) is the element/node next to x0.hare=f(f(x0))whiletortoise!=hare:tortoise=f(tortoise)hare=f(f(hare))# At this point the tortoise position, ν, which is also equal# to the distance between hare and tortoise, is divisible by# the period λ. So hare moving in cycle one step at a time, # and tortoise (reset to x0) moving towards the cycle, will # intersect at the beginning of the cycle. Because the # distance between them is constant at 2ν, a multiple of λ,# they will agree as soon as the tortoise reaches index μ.# Find the position μ of first repetition.    mu=0tortoise=x0whiletortoise!=hare:tortoise=f(tortoise)hare=f(hare)# Hare and tortoise move at same speedmu+=1# Find the length of the shortest cycle starting from x_μ# The hare moves one step at a time while tortoise is still.# lam is incremented until λ is found.lam=1hare=f(tortoise)whiletortoise!=hare:hare=f(hare)lam+=1returnlam,mu

This code only accesses the sequence by storing and copying pointers, function evaluations, and equality tests; therefore, it qualifies as a pointer algorithm. The algorithm uses O(λ + μ) operations of these types, and O(1) storage space. [7]

Brent's algorithm

Richard P. Brent described an alternative cycle detection algorithm that, like the tortoise and hare algorithm, requires only two pointers into the sequence. [8] However, it is based on a different principle: searching for the smallest power of two 2i that is larger than both λ and μ. For i = 0, 1, 2, ..., the algorithm compares x2i1 with each subsequent sequence value up to the next power of two, stopping when it finds a match. It has two advantages compared to the tortoise and hare algorithm: it finds the correct length λ of the cycle directly, rather than needing to search for it in a subsequent stage, and its steps involve only one evaluation of the function f rather than three. [9]

The following Python code shows how this technique works in more detail.

defbrent(f,x0)->(int,int):"""Brent's cycle detection algorithm."""# main phase: search successive powers of twopower=lam=1tortoise=x0hare=f(x0)# f(x0) is the element/node next to x0.whiletortoise!=hare:ifpower==lam:# time to start a new power of two?tortoise=harepower*=2lam=0hare=f(hare)lam+=1# Find the position of the first repetition of length λtortoise=hare=x0foriinrange(lam):# range(lam) produces a list with the values 0, 1, ... , lam-1hare=f(hare)# The distance between the hare and tortoise is now λ.# Next, the hare and tortoise move at same speed until they agreemu=0whiletortoise!=hare:tortoise=f(tortoise)hare=f(hare)mu+=1returnlam,mu

Like the tortoise and hare algorithm, this is a pointer algorithm that uses O(λ + μ) tests and function evaluations and O(1) storage space. It is not difficult to show that the number of function evaluations can never be higher than for Floyd's algorithm. Brent claims that, on average, his cycle finding algorithm runs around 36% more quickly than Floyd's and that it speeds up the Pollard rho algorithm by around 24%. He also performs an average case analysis for a randomized version of the algorithm in which the sequence of indices traced by the slower of the two pointers is not the powers of two themselves, but rather a randomized multiple of the powers of two. Although his main intended application was in integer factorization algorithms, Brent also discusses applications in testing pseudorandom number generators. [8]

Gosper's algorithm

R. W. Gosper's algorithm [10] [11] finds the period , and the lower and upper bound of the starting point, and , of the first cycle. The difference between the lower and upper bound is of the same order as the period, i.e. .

Advantages

The main feature of Gosper's algorithm is that it never backs up to reevaluate the generator function, and is economical in both space and time. It could be roughly described as a concurrent version of Brent's algorithm. While Brent's algorithm gradually increases the gap between the tortoise and hare, Gosper's algorithm uses several tortoises (several previous values are saved), which are roughly exponentially spaced. According to the note in HAKMEM item 132, [11] this algorithm will detect repetition before the third occurrence of any value, i.e. the cycle will be iterated at most twice. This note also states that it is sufficient to store previous values; however, the provided implementation [10] stores values. For example, assume the function values are 32-bit integers and the second iteration of the cycle ends after at most 232 function evaluations since the beginning (viz. ). Then Gosper's algorithm will find the cycle after at most 232 function evaluations, while consuming the space of 33 values (each value being a 32-bit integer).

Complexity

Upon the -th evaluation of the generator function, the algorithm compares the generated value with previous values; observe that goes up to at least and at most . Therefore, the time complexity of this algorithm is . Since it stores values, its space complexity is . This is under the usual assumption, present throughout this article, that the size of the function values is constant. Without this assumption, the space complexity is since we need at least distinct values and thus the size of each value is .

Time–space tradeoffs

A number of authors have studied techniques for cycle detection that use more memory than Floyd's and Brent's methods, but detect cycles more quickly. In general these methods store several previously-computed sequence values, and test whether each new value equals one of the previously-computed values. In order to do so quickly, they typically use a hash table or similar data structure for storing the previously-computed values, and therefore are not pointer algorithms: in particular, they usually cannot be applied to Pollard's rho algorithm. Where these methods differ is in how they determine which values to store. Following Nivasch, [12] we survey these techniques briefly.

Any cycle detection algorithm that stores at most M values from the input sequence must perform at least function evaluations. [18] [19]

Applications

Cycle detection has been used in many applications.

Related Research Articles

<span class="mw-page-title-main">Jensen's inequality</span> Theorem of convex functions

In mathematics, Jensen's inequality, named after the Danish mathematician Johan Jensen, relates the value of a convex function of an integral to the integral of the convex function. It was proved by Jensen in 1906, building on an earlier proof of the same inequality for doubly-differentiable functions by Otto Hölder in 1889. Given its generality, the inequality appears in many forms depending on the context, some of which are presented below. In its simplest form the inequality states that the convex transformation of a mean is less than or equal to the mean applied after convex transformation; it is a simple corollary that the opposite is true of concave transformations.

In mathematics, the Radon–Nikodym theorem is a result in measure theory that expresses the relationship between two measures defined on the same measurable space. A measure is a set function that assigns a consistent magnitude to the measurable subsets of a measurable space. Examples of a measure include area and volume, where the subsets are sets of points; or the probability of an event, which is a subset of possible outcomes within a wider probability space.

In mathematics, particularly in functional analysis, a projection-valued measure is a function defined on certain subsets of a fixed set and whose values are self-adjoint projections on a fixed Hilbert space. A projection-valued measure (PVM) is formally similar to a real-valued measure, except that its values are self-adjoint projections rather than real numbers. As in the case of ordinary measures, it is possible to integrate complex-valued functions with respect to a PVM; the result of such an integration is a linear operator on the given Hilbert space.

<span class="mw-page-title-main">Laplace distribution</span> Probability distribution

In probability theory and statistics, the Laplace distribution is a continuous probability distribution named after Pierre-Simon Laplace. It is also sometimes called the double exponential distribution, because it can be thought of as two exponential distributions spliced together along the abscissa, although the term is also sometimes used to refer to the Gumbel distribution. The difference between two independent identically distributed exponential random variables is governed by a Laplace distribution, as is a Brownian motion evaluated at an exponentially distributed random time. Increments of Laplace motion or a variance gamma process evaluated over the time scale also have a Laplace distribution.

Variational Bayesian methods are a family of techniques for approximating intractable integrals arising in Bayesian inference and machine learning. They are typically used in complex statistical models consisting of observed variables as well as unknown parameters and latent variables, with various sorts of relationships among the three types of random variables, as might be described by a graphical model. As typical in Bayesian inference, the parameters and latent variables are grouped together as "unobserved variables". Variational Bayesian methods are primarily used for two purposes:

  1. To provide an analytical approximation to the posterior probability of the unobserved variables, in order to do statistical inference over these variables.
  2. To derive a lower bound for the marginal likelihood of the observed data. This is typically used for performing model selection, the general idea being that a higher marginal likelihood for a given model indicates a better fit of the data by that model and hence a greater probability that the model in question was the one that generated the data.
<span class="mw-page-title-main">Interior-point method</span> Algorithms for solving convex optimization problems

Interior-point methods are algorithms for solving linear and non-linear convex optimization problems. IPMs combine two advantages of previously-known algorithms:

In statistics and probability theory, a point process or point field is a collection of mathematical points randomly located on a mathematical space such as the real line or Euclidean space. Point processes can be used for spatial data analysis, which is of interest in such diverse disciplines as forestry, plant ecology, epidemiology, geography, seismology, materials science, astronomy, telecommunications, computational neuroscience, economics and others.

<span class="mw-page-title-main">Inverse Gaussian distribution</span> Family of continuous probability distributions

In probability theory, the inverse Gaussian distribution is a two-parameter family of continuous probability distributions with support on (0,∞).

In mathematics, the structure tensor, also referred to as the second-moment matrix, is a matrix derived from the gradient of a function. It describes the distribution of the gradient in a specified neighborhood around a point and makes the information invariant to the observing coordinates. The structure tensor is often used in image processing and computer vision.

<span class="mw-page-title-main">Young stellar object</span> Star in its early stage of evolution

Young stellar object (YSO) denotes a star in its early stage of evolution. This class consists of two groups of objects: protostars and pre-main-sequence stars.

Covariance matrix adaptation evolution strategy (CMA-ES) is a particular kind of strategy for numerical optimization. Evolution strategies (ES) are stochastic, derivative-free methods for numerical optimization of non-linear or non-convex continuous optimization problems. They belong to the class of evolutionary algorithms and evolutionary computation. An evolutionary algorithm is broadly based on the principle of biological evolution, namely the repeated interplay of variation and selection: in each generation (iteration) new individuals are generated by variation of the current parental individuals, usually in a stochastic way. Then, some individuals are selected to become the parents in the next generation based on their fitness or objective function value . Like this, individuals with better and better -values are generated over the generation sequence.

<span class="mw-page-title-main">Generalized Pareto distribution</span> Family of probability distributions often used to model tails or extreme values

In statistics, the generalized Pareto distribution (GPD) is a family of continuous probability distributions. It is often used to model the tails of another distribution. It is specified by three parameters: location , scale , and shape . Sometimes it is specified by only scale and shape and sometimes only by its shape parameter. Some references give the shape parameter as .

In mathematics, in particular in measure theory, a content is a real-valued function defined on a collection of subsets such that

In mathematics, the Prékopa–Leindler inequality is an integral inequality closely related to the reverse Young's inequality, the Brunn–Minkowski inequality and a number of other important and classical inequalities in analysis. The result is named after the Hungarian mathematicians András Prékopa and László Leindler.

In probability and statistics, the Tweedie distributions are a family of probability distributions which include the purely continuous normal, gamma and inverse Gaussian distributions, the purely discrete scaled Poisson distribution, and the class of compound Poisson–gamma distributions which have positive mass at zero, but are otherwise continuous. Tweedie distributions are a special case of exponential dispersion models and are often used as distributions for generalized linear models.

In mathematics, the spectral theory of ordinary differential equations is the part of spectral theory concerned with the determination of the spectrum and eigenfunction expansion associated with a linear ordinary differential equation. In his dissertation, Hermann Weyl generalized the classical Sturm–Liouville theory on a finite closed interval to second order differential operators with singularities at the endpoints of the interval, possibly semi-infinite or infinite. Unlike the classical case, the spectrum may no longer consist of just a countable set of eigenvalues, but may also contain a continuous part. In this case the eigenfunction expansion involves an integral over the continuous part with respect to a spectral measure, given by the Titchmarsh–Kodaira formula. The theory was put in its final simplified form for singular differential equations of even degree by Kodaira and others, using von Neumann's spectral theorem. It has had important applications in quantum mechanics, operator theory and harmonic analysis on semisimple Lie groups.

In combinatorics, a Davenport–Schinzel sequence is a sequence of symbols in which the number of times any two symbols may appear in alternation is limited. The maximum possible length of a Davenport–Schinzel sequence is bounded by the number of its distinct symbols multiplied by a small but nonconstant factor that depends on the number of alternations that are allowed. Davenport–Schinzel sequences were first defined in 1965 by Harold Davenport and Andrzej Schinzel to analyze linear differential equations. Following Atallah (1985) these sequences and their length bounds have also become a standard tool in discrete geometry and in the analysis of geometric algorithms.

<span class="mw-page-title-main">Poisson distribution</span> Discrete probability distribution

In probability theory and statistics, the Poisson distribution is a discrete probability distribution that expresses the probability of a given number of events occurring in a fixed interval of time if these events occur with a known constant mean rate and independently of the time since the last event. It can also be used for the number of events in other types of intervals than time, and in dimension greater than 1.

In differential geometry, a fibered manifold is surjective submersion of smooth manifolds YX. Locally trivial fibered manifolds are fiber bundles. Therefore, a notion of connection on fibered manifolds provides a general framework of a connection on fiber bundles.

(Stochastic) variance reduction is an algorithmic approach to minimizing functions that can be decomposed into finite sums. By exploiting the finite sum structure, variance reduction techniques are able to achieve convergence rates that are impossible to achieve with methods that treat the objective as an infinite sum, as in the classical Stochastic approximation setting.

References

  1. Joux, Antoine (2009), Algorithmic Cryptanalysis, CRC Press, p. 223, ISBN   9781420070033 .
  2. 1 2 Joux (2009 , p. 224).
  3. 1 2 Knuth, Donald E. (1969), The Art of Computer Programming, vol. II: Seminumerical Algorithms, Addison-Wesley, p. 7, exercises 6 and 7
  4. Handbook of Applied Cryptography, by Alfred J. Menezes, Paul C. van Oorschot, Scott A. Vanstone, p. 125, describes this algorithm and others
  5. Floyd, R.W. (1967), "Nondeterministic Algorithms", J. ACM, 14 (4): 636–644, doi: 10.1145/321420.321422 , S2CID   1990464
  6. The Hash Function BLAKE, by Jean-Philippe Aumasson, Willi Meier, Raphael C.-W. Phan, Luca Henzen (2015), p. 21, footnote 8
  7. Joux (2009), Section 7.1.1, Floyd's cycle-finding algorithm, pp. 225–226.
  8. 1 2 3 4 Brent, R. P. (1980), "An improved Monte Carlo factorization algorithm" (PDF), BIT Numerical Mathematics , 20 (2): 176–184, doi:10.1007/BF01933190, S2CID   17181286 .
  9. Joux (2009), Section 7.1.2, Brent's cycle-finding algorithm, pp. 226–227.
  10. 1 2 "Archived copy". Archived from the original on 2016-04-14. Retrieved 2017-02-08.{{cite web}}: CS1 maint: archived copy as title (link)
  11. 1 2 "Hakmem -- Flows and Iterated Functions -- Draft, Not Yet Proofed". Archived from the original on 2020-03-18. Retrieved 2024-05-02.
  12. 1 2 3 4 Nivasch, Gabriel (2004), "Cycle detection using a stack", Information Processing Letters , 90 (3): 135–140, doi:10.1016/j.ipl.2004.01.016 .
  13. Schnorr, Claus P.; Lenstra, Hendrik W. (1984), "A Monte Carlo factoring algorithm with linear storage", Mathematics of Computation , 43 (167): 289–311, doi:10.2307/2007414, hdl: 1887/3815 , JSTOR   2007414 .
  14. 1 2 Teske, Edlyn (1998), "A space-efficient algorithm for group structure computation", Mathematics of Computation , 67 (224): 1637–1663, Bibcode:1998MaCom..67.1637T, doi: 10.1090/S0025-5718-98-00968-5 .
  15. Sedgewick, Robert; Szymanski, Thomas G.; Yao, Andrew C.-C. (1982), "The complexity of finding cycles in periodic functions", SIAM Journal on Computing , 11 (2): 376–390, doi:10.1137/0211030 .
  16. van Oorschot, Paul C.; Wiener, Michael J. (1999), "Parallel collision search with cryptanalytic applications", Journal of Cryptology , 12 (1): 1–28, doi: 10.1007/PL00003816 , S2CID   5091635 .
  17. 1 2 Quisquater, J.-J.; Delescaille, J.-P., "How easy is collision search? Application to DES", Advances in Cryptology – EUROCRYPT '89, Workshop on the Theory and Application of Cryptographic Techniques, Lecture Notes in Computer Science, vol. 434, Springer-Verlag, pp. 429–434, doi: 10.1007/3-540-46885-4_43 .
  18. 1 2 Fich, Faith Ellen (1981), "Lower bounds for the cycle detection problem", Proc. 13th ACM Symposium on Theory of Computing , pp. 96–105, doi:10.1145/800076.802462, S2CID   119742106 .
  19. Allender, Eric W.; Klawe, Maria M. (1985), "Improved lower bounds for the cycle detection problem", Theoretical Computer Science , 36 (2–3): 231–237, doi: 10.1016/0304-3975(85)90044-1 .
  20. Pollard, J. M. (1975), "A Monte Carlo method for factorization", BIT, 15 (3): 331–334, doi:10.1007/BF01933667, S2CID   122775546 .
  21. Pollard, J. M. (1978), "Monte Carlo methods for index computation (mod p)", Mathematics of Computation , 32 (143), American Mathematical Society: 918–924, doi:10.2307/2006496, JSTOR   2006496, S2CID   235457090 .
  22. 1 2 Kaliski, Burton S. Jr.; Rivest, Ronald L.; Sherman, Alan T. (1988), "Is the Data Encryption Standard a group? (Results of cycling experiments on DES)", Journal of Cryptology , 1 (1): 3–36, doi: 10.1007/BF00206323 , S2CID   17224075 .
  23. Joux (2009), Section 7.5, Collisions in hash functions, pp. 242–245.
  24. Van Gelder, Allen (1987), "Efficient loop detection in Prolog using the tortoise-and-hare technique", Journal of Logic Programming, 4 (1): 23–31, doi: 10.1016/0743-1066(87)90020-3 .
  25. Auguston, Mikhail; Hon, Miu Har (1997), "Assertions for Dynamic Shape Analysis of List Data Structures", AADEBUG '97, Proceedings of the Third International Workshop on Automatic Debugging, Linköping Electronic Articles in Computer and Information Science, Linköping University, pp. 37–42.