- Email: [email protected]
On approximation accuracy of the Chaos Game's finite-time activity
Computers & Graphics 26 (2002) 753–764
Technical section
On approximation accuracy of the Chaos Game’s finite-time activity Tomek Martyn* Computer Graphics Laboratory, Institute of Computer Science, Warsaw University of Technology, ul. Nowowiejska 15/19, 00-665, Warsaw, Poland
Abstract This paper describes a method for estimating probability that the Chaos Game generates, within L steps, a sequence of points which approximates the IFS attractor at a given accuracy. The presented algorithm is a hybrid based on methods from the theory of Markov chains used in the pattern-matching context. Unlike the previous approaches, ours is general both in the sense of component IFS mappings and the associated probability distribution. r 2002 Elsevier Science Ltd. All rights reserved. Keywords: Iterated function system; Fractals; IFS attractors; The Chaos Game; Probability
1. Introduction The iterated function system (IFS) is the most popular way to describe fractal shapes. Each IFS specifies implicitly the set of points called the IFS attractor. In order to obtain an explicit image of the attractor you can use one of the known methods of IFS rendering [1]. Because of relatively easy implementation, the method of choice is usually the random iteration algorithm [2] that is known wider as the Chaos Game. As far as the infinite number of iteration steps is considered, the algorithm almost surely (i.e., with probability one) generates an infinite sequence of points which, ignoring a finite number of initial points, approximates the attractor with an arbitrarily small error. The formal proofs of this fact usually do not take into account the results of the algorithm’s finite-time activity. Since in practice only a finite number of points can be generated, knowledge of the asymptotic behaviour of the process is very interesting from the mathematical point of view, but as it concerns the Plato world exclusively, it says nothing about the world’s merely shadows that we are
*Tel.: +48-22-285-31-79; fax: +48-22-285-16-35. E-mail address: [email protected] (T. Martyn).
given. Therefore, knowing all those things about infinity but being a material person living in the ordinary nonperfect world, one could ask: How many points have to be generated to approximate an IFS attractor at a given accuracy? As far as you use the Chaos Game just to have a glance at fractals for pleasure, trying to answer this question may seem another whim of a mathematician, because you are always able to generate more points to refine an image. However, when the Chaos Game is regarded as a raster conversion function for purposes of visualization at pixel size accuracy, or if the stream of points generated by the algorithm is directed to the input of a computational geometry procedure (e.g. for determination of the minimum ball enclosing an IFS attractor), then the stated problem is of key importance. As a previous work concerning, in many aspects, the above problem one should mention the paper by Goodman [3] devoted to the analysis of the Chaos Game in its original version (i.e., for generating Sierpinski’s Triangle with the uniform probability distribution). The approach presented in this paper has its roots in Goodman’s analysis but we extend the scope to the general case in which an IFS can consist of any type of contraction mappings equipped with an arbitrary assigned probability distribution.
0097-8493/02/$ - see front matter r 2002 Elsevier Science Ltd. All rights reserved. PII: S 0 0 9 7 - 8 4 9 3 ( 0 2 ) 0 0 1 3 0 - 9
T. Martyn / Computers & Graphics 26 (2002) 753–764
754
The diameters of these subsets satisfy
2. IFS and the Chaos Game An iterated function system ðX ; w1 ; y; wN Þ; shortly IFS, is defined as a finite set of contraction mappings wi : X -X ; i ¼ 1; y; N; on a complete metric space ðX ; dÞ: Using the Banach fixed point theorem one can prove that each IFS implicitly specifies a compact (i.e., totally bounded and closed) and nonempty set A; which is the unique solution of the equation A¼
N [
wi ðAÞ:
ð1Þ
i¼1
Such a set is called the attractor of the IFS. The assignation of probabilities to IFS mappings leads to the concept of an IFS with probabilities ðX ; wP where pi > 0; i ¼ 1; y; N; 1 ; y; wN ; p1 ; y; pN Þ; and N i¼1 pi ¼ 1: Having such an IFS, we can approximate its attractor by means of the random iteration algorithm [2], known wider as the Chaos Game. The algorithm works in an iterative way and produces a sequence of points fxk gN k¼0 such that x0 is any point of X and xkþ1 ¼ wik ðxk Þ; where wik is chosen independently in each step with the probability pik from the set of the IFS mappings. The sequence fxk gN k¼0 approximates the attractor with an arbitrarily small approximation error e > 0: If the starting point x0 is in the attractor, every point of the sequence belongs to the attractor as well. The initial point can be easily determined as the fixed point x ¼ wi ðxÞ of any IFS mapping. In this case it is possible to achieve the best approximation of an attractor with the least number of iteration steps, which is desirable for computer graphics purposes. Therefore, our attention will focus on the stated case. However, if necessary, the considerations presented in this paper can be easily extended to the general case in which x0 is any point of X :
3. IFS attractor structure Let A be the ðX ; w1 ; y; wN ; p1 ; y; pN Þ: attractor A consists of N Since the IFS mappings wi of each subset satisfies
attractor of an IFS On the basis of (1), the subsets wi ðAÞ; iAf1; y; Ng: are contractive, the diameter
diamðwi ðAÞÞpsðwi ÞdiamðAÞ;
8iAf1; y; Ng;
ð2Þ
where sðÞ denotes the contractivity factor of a given mapping. In turn, each subset wi ðAÞ is made up of N subsets of the form wi 3wj ðAÞ; jAf1; y; Ng; for ! N N [ [ wj ðAÞ ¼ wi 3wj ðAÞ: ð3Þ wi ðAÞ ¼ wi j¼1
j¼1
diamðwi 3wj ðAÞÞpsðwi 3wj ÞdiamðAÞ psðwi Þsðwj ÞdiamðAÞ:
ð4Þ
Generally, by induction, any subset wi1 3wi2 3?3win ðAÞ; ik Af1; y; Ng; of the attractor consists of N subsets y; Ng; whose diawi1 3wi2 3?3win 3winþ1 ðAÞ; inþ1 Af1; Qnþ1 meters are not greater than k¼1 sðwik Þ diamðAÞ: The contractivity factors sðwik ÞA½0; 1Þ; ik Af1; y; Ng; so given e > 0; we are able to decompose the attractor to a finite number of subsets wi1 3?3win ðAÞ; whose diameters do not exceed e: Such a family of the attractor subsets we will call an e-decomposition of the attractor. Choosing at least one arbitrary point from each subset of the e-decomposition we can construct a finite set Ae such that the Hausdorff distance hðA; Ae Þpe; i.e., the set Ae approximates A with an error not greater than e:
4. Estimating probability Let e > 0 be given and let aAA be an arbitrary point of the attractor. Then a belongs to at least one of the e-decomposition subsets; say aAwa1 3?3wan 1 3wan ðAÞ; ai Af1; y; Ng: Assume the initial point for the Chaos Game is x0 AA: Then every point generated by the algorithm also belongs to the attractor, because from (1) the IFS mappings take A into itself. The n-fold composition wa1 3?3wan maps the attractor A to its subset whose diameter is not greater than e and which includes the point a: Therefore, appearing of the sequence wan ; y; wa1 in arbitrary n successive steps of the algorithm is the sufficient condition for at least one of the points produced by the algorithm to be within the distance e to a: Each mapping sequence wan ; y; wa1 is specified by the string of indices an an 1 ya1 (called in the IFS nomenclature an address of the subset wan ; y; wa1 ðAÞ), so given e > 0; we can restrict our considerations to the strings from the set that is induced by the P e-decomposition of the attractor. Let ¼ f1; y; Ng be P the set of the IFS map indices (or symbols), P and let denote the set of all (finite) strings over : Let fp1 ; y;P pN g be a set of IFS probabilities. Denote by DC the set of strings, which is induced by the attractor e-decomposition, i.e., if an an 1 ya1P AD; then diamðwa1 3?3wan ðAÞÞpe: Let foi gLi¼1 ; oi A ; and respectively fxi gLi¼1 ; xi AA; be a string of indices, and a sequence of points, respectively, produced in L initial steps of the algorithm. Then fxi gLi¼1 approximates the attractor with an error not greater than e, if foi gLi¼1 includes all strings from D as its substrings. Just like in the case of Goodman’s analysis of the Chaos Game for Sierpinski’s Triangle [3], the process of generating the strings from D within foi gLi¼1 possesses the Markov property, because knowledge of the present
T. Martyn / Computers & Graphics 26 (2002) 753–764
substring completely determines probability for the next one to arise, without regard to the prior history of the process. The main differences are that, in our general case, probabilities do not have P to be distributed uniformly on the symbols from ; and the strings from D do not have to be of equal lengths. Let us denote by ‘!’ and respectively by ‘g’ P the P prefix and suffix relations on ; i.e., for x; yA P ; y!x and respectively ygx; iff there is some wA such that x ¼ y3w and x ¼ w3y; respectively, where ‘3’ denotes thePconcatenation of strings. Now, defining Dp :¼ fxA : x!s for some sADg and identifying the set S ¼ D,Dp with the space of states, we can consider the action of the Chaos Game within the category of a Markov chain, i.e., a sequence fXi gN i¼0 of random variables with values in S; such that for every sequence fsi gli¼0 ; si AS; we have PðXi ¼ si jXi 1 ¼ si 1 ; y; X1 ¼ s1 ; X0 ¼ s0 Þ ¼ PðXi ¼ si jXi 1 ¼ si 1 Þ:
ð5Þ
(Note that the empty string l; i.e., consisting ofPno symbols, is a prefix and suffix of every string from ). The Markov chain proceeds from state to state as follows: P P Let Xi 1 ¼ a3s; where aA and sAS; and let bA be the symbol chosen in the ith step (i > 0); then X0 ¼ l ( Xi ¼
s3b; a3s3bAS;
if s3bAS; otherwise;
ð6Þ
where probability of a transition from the state to one of the states is equal to probability of choosing the symbol b: Thus, in order to estimate probability that the Chaos Game generates, in L initial steps, the sequence of points that approximates the attractor at a given accuracy, all we have to do is determine probability that the Markov chain visits, in L steps, all states from D: Pð8sAD; (kpL; Xk ¼ sÞ ¼ P
L \[
! fXk ¼ sg
sAD k¼1
¼1 P
!
L [\
fXk asg :
sAD k¼1
ð7Þ Using Bool’s inequality, we obtain P
L [\
! fXk asg p
sAD k¼1
¼ jDj
X sAD
X sAD
P
L [ k¼1
P
L \
! fXk asg
k¼1
!
fXk ¼ sg ;
ð8Þ
755
and hence Pð8sAD; (kpL; Xk ¼ sÞX ! L X [ 1þ P fXk ¼ sg jDj; sAD
ð9Þ
k¼1
where jDj denotes the number of states in D: Thus, in order to achieve our goal, we have to
S determine for each sAD the value P Lk¼1 fXk ¼ sg ; that is probability that the Markov chain fXi gN i¼0 will visit, in L initial steps, the state s at least once.
5. Probability of generating strings It is well-known fact that probability of the occurrence of a substring sAS within a random sequence P depends not only on the probability of foi gLi¼1 over the component symbols to appear but also on the structure of s itself. For example, in the case P when probabilities are distributed uniformly on ; the occurrence of a substring that consists of a single repeated symbol is a less probable event than the occurrence of any string of the P same length but composed of varied symbols of : The problem considered here resembles in many aspects the one of pattern matching. In fact, the solution presented further can be regarded as an adaptation of the method proposed by Knuth et al. [4] for searching for a string within a given text by means of finite automata (for a good discussion on the subjectPsee also [5]). Let a ¼ a1 a2 ; y; an AD; ai A ; and let aðkÞ denote the a string’s prefix of the length kA½0; n : Assuming iXn; it is easy to see that the following identities are satisfied: ð1Þ Xi ¼ a3aðn 1ÞgXi 1 and an is chosen in the ith step; ð2Þ aðn 1ÞgXi 1 3aðn 2ÞgXi 2 and an 1 is chosen in the ði 1Þth step; ^ ð3Þ að1ÞgXi nþ1 3lgXi n and a1 is chosen in the ði n þ 1Þth step:
ð10Þ
From the above sequence of implications we conclude that probability for a given string x to give rise to the string a depends only on the maximal suffix of x; which matches a prefix of a; no matter what the x string’s remaining (i.e., beyond the suffix) initial symbols are. Thus, if x; yAD are strings with suffixes fitting the a prefix at the same positions, then they are indistinctive states from the point of view of conditional probability for the Markov chain to visit the state a: Therefore, we can construct a Markov chain with the space of states consisting of subsets of the space S; which are induced
T. Martyn / Computers & Graphics 26 (2002) 753–764
756
with transition probability specified by
by a suffix equivalence relation. The specification of the process affects the statistical behaviour of this chain to be identical to that one of fXi gN i¼0 in the sense of probability of the string a to arise in the random sequence foi gLi¼1 : Moreover, such an approach will allow us to solve the considered problem by means of standard techniques from the theory of Markov chains used in the context of pattern matching.
PðYi ¼ Qm jYi 1 ¼ Qk Þ 8 > < 1; ¼ 0; > : PðsðXi Þ ¼ mjsðXi 1 Þ ¼ kÞ;
Let s : S-N,f0g be the suffix function of the string a; i.e., ð11Þ
Now, we can introduce a relation, denoted by ‘B’, on S such that xBy :¼ ðsðxÞ ¼ sðyÞÞ:
ð12Þ
The relation satisfies properties of equivalence relation. Hence, we immediately get the partition Q ¼ fQ0 ; Q1 ; y; Qn g of S with respect to B; where Qi ¼ ½x is the equivalence class of x; such that i ¼ sðxÞ: Let fYi gN i¼0 be a sequence of random variables taking values from the set Q in the following way: P Let Yi 1 ¼ Qk and let bA be a symbol chosen in the ith ði > 0Þ step; then Y 0 ¼ Q0 ; ( Y i ¼ Qm ; m ¼
PðYl ¼ Qn Þ ¼ Pð(kpL; Xk ¼ aÞ; if k ¼ n otherwise
n; sðaðkÞ3bÞ;
(a) Q0
β Q1
β
α Q2
Q3
α
(b) State Q0 Q1 Q2 Q3 Q4 Q5 Q6 Q7 Q8
α Q4
χ
β Q6
Q5
α
α, β , χ
χ Q7
Q8
α
α
α
ð15Þ
i.e., probability that fYi gN i¼0 ; after L steps, is in the state Qn is equal to probability that fXi gN i¼0 visits the state a
ð13Þ
α
ð14Þ
From the above definition the sequence fYi gN i¼0 satisfies the properties of a Markov chain. Moreover, the chain is time homogenous, since the transition probability matrix induced by (14) is the same for every Yi ; i.e., it does not depend on the number of steps. The conditional probability PðsðXi Þ ¼ mjsðXi 1 Þ ¼ kÞ is specified by P the sum of probabilities associated with symbols bA that map a state Xi 1 ¼ x3aðkÞ to a state Xi ¼ y3aðkÞ3b ¼ z3aðmÞ; x; y; zAS; according to the transition rules (6). Note that the properties of suffix matching imply that there is usually at most one symbol taking a given state for which sðXi 1 Þ ¼ k to another one such that sðXi Þ ¼ m; see Fig. 1. The exceptions are, from definition, the transitions from the final state Qn ; as well as the transitions to the states with sðX Pi Þ ¼ 0 and the ones for which there is no symbol in to take a state to another. Obviously, in the last case probability is equal to zero. The definition of fYi gN i¼0 distinguishes the state Qn as the unique fixed point of the chain in the sense that for every kXi; Yi ¼ Qn implies Yk ¼ Qn : Since Qn ¼ fag; therefore YL ¼ Qn if and only if Xi ¼ a for some ipL: Putting the above observations in the context of the sequence of identities (10), we conclude that
5.1. The equivalent Markov chain
sðxÞ :¼ maxfk : aðkÞgxg; xAS:
if k ¼ n ¼ m; if k ¼ nam; if kan:
(c) α Q1 Q1 Q3 Q1 Q5 Q1 Q5 Q1 Q8
β Q0 Q2 Q0 Q4 Q0 Q6 Q0 Q0 Q8
χ Q0 Q0 Q0 Q0 Q0 Q0 Q7 Q8 Q8
Q0 Q1 Q2 Q3 Q4 Q5 Q6 Q7 Q8
Q0 β+χ χ β+χ χ β+χ χ β β 0
Q1 α α 0 α 0 α 0 α 0
Q2 0 β 0 0 0 0 0 0 0
Q3 0 0 α 0 0 0 0 0 0
Q4 0 0 0 β 0 0 0 0 0
Q5 0 0 0 0 α 0 α 0 0
Q6 0 0 0 0 0 β 0 0 0
Q7 0 0 0 0 0 0 χ 0 0
Q8 0 0 0 0 0 0 0 χ 1
Fig. 1. The Markov chain induced by the string abababww: (a) the graph of transitions (the transitions to the state Q0 have been omitted); (b) the transition table; (c) the transition probability matrix—the component symbols are identified with the associated values of probabilities.
T. Martyn / Computers & Graphics 26 (2002) 753–764
within L steps at least once. As we stated above, fYi gN i¼0 is a time homogenous Markov chain, so in order to compute this probability we can use a standard technique from the relevant theory of such processes. Having a ðn þ 1Þ ðn þ 1Þ transition probability matrix M; the probability distribution for the variable YL is specified by the vector y ¼ x ML ; x; yARnþ1 ; where the ith coordinate of the vector y is probability of the event YL ¼ Qiþ1 ; and x specifies an initial distribution for Y0 ; in our case x ¼ ½1; 0; y; 0 : Because the exponent L as the number of the Chaos Game steps is supposed to be a relatively large number, it is convenient to use the following algorithm of fast matrix involution that raises a n n matrix M to a given power L in time Oðn3 log LÞ: out’identity matrix repeat { if L mod 2=1 then out’out M M’M M L’Lb1 } until L ¼ 0
757
iSþ 1:
Because, from Proposition 1, the set 1q aðp ðkÞÞ is composed of all the proper qAN suffixes aðiÞ of aðkÞ; we can determine the value pðk þ 1Þ on the basis of the values p1q ðkÞAfpðkÞ; pðk
1Þ; y; pð1Þg computedSin the previous steps. We do so by searching the set qAN aðp1q ðkÞÞ for the prefix aðiÞ such that aiþ1 ¼ akþ1 ; where iAf0; y; k 1g is maximal. If there is such a string aðiÞ in the set, then pðk þ 1Þ ¼ i þ 1; otherwise pðk þ 1Þ ¼ 0: Furthermore, note that for any i; mAN such that iom; the length of the prefix aðp1i ðkÞÞ is not greater than the length of aðp1m ðkÞÞ: S Thus, in practice, we do not have to deal with the set qAN aðp1q ðkÞÞ of prefixes directly to compute pðk þ 1Þ; but we just search a finite set f1; y; mg; mpk; for the minimal q such that the ðk þ 1Þth and ðp1q ðkÞ þ 1Þth symbols of the string a are equal. As a result we have the following algorithm for the determination of the value pðkÞ on the basis of values pðiÞ; iAf1; y; k 1g:
5.2. Determination of the transition probability matrix
let a1 yak ADp be a prefix of aAD; then if ko2 then pðkÞ’0 else { i’pðk 1Þ while i > 0 and aiþ1 aak do i’pðiÞ if aiþ1 ¼ ak then i’i þ 1 pðkÞ’i }
The last remaining problem to solve is the determination of the transition probability matrix. In order to make it in an efficient way our approach takes advantage of the so-called prefix function introduced by Knuth et al. for their method of pattern matching [4]. For a given string a ¼ a1 a2 yan AD; the prefix function p : f0; 1; y; ng-f0; 1; y; n 1g is defined as
Proposition 2. Let a ¼ a1 a2 yan AD, and let s and p be the suffix and prefix function ofP a, respectively. Then, for every kAf0; y; n 1g and bA 8 if k ¼ 0 and baa1 ; > < 0; sðaðkÞ3bÞ ¼ k þ 1; if b ¼ akþ1 ; ð17Þ > : sðaðpðkÞÞ3bÞ; if k > 0 and baakþ1 :
where ‘b’ denotes the right-bit-shift operation in the C language convention.
pð0Þ ¼ pð1Þ ¼ 0; pðqÞ :¼ maxfk : koq and aðkÞgaðqÞg;
ð16Þ
i.e., for q > 0; pðqÞ is the maximum length of a prefix of the string a, which is simultaneously a proper suffix of aðqÞAS: Proposition 1. Let a ¼ a1 a2 yan AD be a string with the prefix function p: Denote by p1q the q-fold S composition of p: Then, for a given kAf1; y; ng; the set qAN aðp1q ðkÞÞ consists of all the strings, which being of the form aðiÞ, iok, are simultaneously the proper suffixes of the prefix aðkÞ. Proof. See e.g. [5].
&
It is easy to see that aði þ 1Þgaðk þ 1Þ if and only if aðiÞgaðkÞ and aiþ1 ¼ akþ1 : Hence, if aðiÞ is the proper suffix of aðkÞ such that the length of aðiÞ is maximal in the sense of the condition aiþ1 ¼ akþ1 ; then pðk þ 1Þ ¼
Proof. The first and second equalities are obvious. To prove the last one, we will show that if kAf0; y; n 1g and baakþ1 then aðiÞgaðkÞ3b3aðiÞgaðpðkÞÞ3b: Hence, from the suffix function definition, we will get the thesis. First of all, let us note that the suffix relation is transitive. Furthermore, it is not difficult to see that if xgz; ygz; and x is longer than y, then ygx: Let aðiÞgaðkÞ3b: Since baakþ1 ; we have aðiÞaaðkÞ3b; hence aðiÞ is a proper suffix of aðkÞ3b: From the definition of the function p; aðpðkÞÞ is the proper suffix of aðkÞ so that the length of aðpðkÞÞ is maximal. Therefore aðpðkÞÞ3b is the maximal proper suffix of the string aðkÞ3b: Since aðiÞ is a proper suffix of aðkÞ3b; using our observations on the suffix relation, we conclude that aðiÞgaðpðkÞÞ3b: Conversely, suppose that aðiÞgaðpðkÞÞ3b: As we note above aðpðkÞÞ3bgaðkÞ3b; and the suffix relation is transitive. Hence aðiÞgaðkÞ3b: &
758
T. Martyn / Computers & Graphics 26 (2002) 753–764
Using the above observations, we can construct an efficient algorithm for the successive determination of the transitions (13) and their probabilities (14). Given a state Qk ; in order to compute transitions, the algorithm uses information on transitions determined for one of the previous states Qi ; iok: Moreover, the algorithm works in on-line manner, i.e., to make the computation concerning the state Qk ; it does not require knowledge of the further string suffix symbols am nor information about the states Qm ; kom: Consequently, the transition probability matrices can be determined simultaneously together with constructing the individual strings from the set of the states D during the e-decomposition process. As a result, we have the following algorithm for the determination of the vector M½k; : that specifies probabilities of the transitions from the state Qk to others. Given a string symbol akþ1 ; the algorithm calculates M½k; : using simultaneously determined vaP lues of the transition function s½k; b ; bA ; until the input symbols are exhausted: let Qk be the current state; then for each state Qm do M½k; m ’0 P if there is a symbol akþ1 A on the input then { calculate theP value p½k for each bA do { if b ¼ akþ1 then s½k; b ’k þ 1 else if k ¼ 0 then s½k; b ’0 else s½k; b ’s½p½k ; b M½k; s½k; b
’M½k; s½k; b
þ p½b } } else M½k; k ’1 where p½b is probability associated with the symbol b (or equivalently with the IFS mapping wb ).
6. The algorithm—summary Let ðX ; w1 ; y; wN ; p1 ; y; pN Þ be an IFS with the attractor A: The heart of the algorithm is the following procedure that determines the sum of
S P recurrently L probabilities during the esAD P k¼1 fXk ¼ sg decomposition process: procedure Descent(MapComp, k) { P for each aA do { Calculate M[k; . ] with the symbol a on the input. Determine a new composition of IFS mappings and save it in the local variable NewMapComp. NewMapComp’MapComp3wa if diamðNewMapCompðAÞÞpe then { Here is the end of the current string,
so calculate M½k þ 1; : for the empty input. Increment the global variable NumStrings that represents the number of strings in the set of states D: NumStrings’NumStrings+1 Finally, calculate probability the Markov chain visits the string (state), and then use that value to increment the global variable Prob Prob’Prob þ ML ½0; k þ 1 } else Descent(NewMapComp, k þ 1) } } P where ¼ f1; y; Ng; and L is the number of the Chaos Game steps. Now, armed with the above procedure, we calculate the estimated probability (9) as follows: let map be equal to the identity mapping Descent(map, 0) Probability’1+Probability NumStrings.
7. Results and open problems Two classic IFS attractors were tested for the sufficient number of the Chaos Game steps to obtain the approximation of an attractor at given accuracy and probability. The first remark is that the results obtained by means of our approach appear to substantiate the recommendation by Barnsley to continue iterating millions of times [6, p. 28], especially when one looks at the numbers of steps calculated for Barnsley’s Fern (see Fig. 4). Also in the case of Sierpinski’s Triangle W with equal probabilities the results agree in a sense with ones determined by Goodman [3]. In the latter case, the main reason of the existing differences in number of steps follows from the fact that, for each string, our method calculates probability of the string to occur in the accurate way, whereas Goodman estimates it on the basis of the worst case in which a string consists of a single repeated symbol. As a consequence, our method is more accurate which manifests itself in a fewer number of steps than presented in [3]. For example, Goodman estimates that for the approximation with an error not greater than diamðDÞ2 9 (the length of strings in D is 9), ‘‘the chance that more than 656,829 are required is less 8.44%’’, whereas given such the approximation error and probability 100%–8.44%=91.56%, our algorithm results in 221,825 steps. Fig. 2 presents the results that our method supplied for Sierpinski’s Triangle. The selected values of e suits visualization at 100 100, 200 200, and 400 400
T. Martyn / Computers & Graphics 26 (2002) 753–764
759
Numer of steps
(a) ε = diam(A ) /100 275 000 250 000 225 000 200 000 175 000 150 000 125 000 100 000 75 000 50 000 25 000 0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
14709
14965
15256
15591
15989
16476
17105
17993
19515
[10/24, 7/24, 7/24]
18230
18621
19067
19585
20204
20967
21963
23385
25863
[1/2, 1/3, 1/6]
138593 144109 150529 158163 167499 179377 195422 219419 264102
[1/3, 1/3, 1/3]
Probability
(b) ε = diam(A ) /200
Numer of steps
1 650 000 1 500 000 1 350 000 1 200 000 1 050 000 900 000 750 000 600 000 450 000 300 000 150 000 0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
51258
52028
52901
53910
55104
56566
58452
61113
65667
[10/24, 7/24, 7/24]
71026
72398
73963
75780
77945
80615
84091
89048
97659
[1/2, 1/3, 1/6]
857834 890499 928507 973678 1028898 1099126 1193951 1335737 1599805
[1/3, 1/3, 1/3]
Probability
(c) ε = diam(A ) /400 10 000 000 9 000 000 Numer of steps
8 000 000 7 000 000 6 000 000 5 000 000 4 000 000 3 000 000 2 000 000 1 000 000 0 [1/3, 1/3, 1/3]
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
175262 177575 180198 183227 186811 191198 196856 204834 218483
[10/24, 7/24, 7/24] 272901 277702 283172 289522 297077 306389 318491 335722 365587 [1/2, 1/3, 1/6]
5288610 5482463 5707940 5975800 6303110 6719206 7280805 8120229 9683557 Probability
Fig. 2. The sufficient number of the Chaos Game steps to produce an approximation of Sierpinski’s Triangle at a given probability and error e: The results are listed for three different IFS probability vectors. The error was measured using the maximum metric.
screen resolution, respectively. In order to get some view how IFS probabilities influence the number of the Chaos Game’s steps, three different distributions of probabilities have been tested. In this case, as one might have expected, the uniform distribution results in the
least number of steps. This is no surprise, however, the noteworthy thing is that the graphs of the presented results look very similar to each other, with no regard of approximation accuracy. (The same can be noticed for Barnsley’s fern; Fig. 4).
3.485 3.477 3.468 3.458 3.433 3.417 3.396 3.365 3.361 3.356 3.351 3.345 3.33
[1/2, 1/3, 1/6]
[1/2, 1/3, 1/6]
Fig. 3. The values of the function (18) (see text) evaluated for Sierpinski’s Triangle.
Probability
6.165 6.157 6.147 6.137 6.113 6.098 6.079 6.053 6.05 6.046 6.042 6.038 6.028 6.023 6.016 6.009 6.008 6.007 6.006 6.005 6.004 6.003 6.002 6.001
3.642 3.63
3.419 3.413 3.406 3.399 3.38 3.368 3.352 3.327 3.324 3.32 3.316 3.311 3.299 3.291 3.281 3.264 3.262 3.259 3.256 3.253 3.245 3.239 3.231 3.219
0.100 0.200 0.300 0.400 0.600 0.700 0.800 0.900 0.910 0.920 0.930 0.940 0.960 0.970 0.980 0.990 0.991 0.992 0.993 0.994 0.996 0.997 0.998 0.999
Probability
6.045 6.041 6.03 6.024 6.017 6.009 6.008 6.008 6.007 6.006 6.004 6.003 6.002 6.001
5 [10/24, 7/24, 7/24] 3.842 3.836 3.829 3.821 3.801 3.787 3.77 3.744 3.74 3.736 3.731 3.726 3.713 3.705 3.693 3.676 3.673 3.671 3.668 3.664 3.656 3.6
[1/3, 1/3, 1/3]
6.25 6 5.75 5.5 5.25 5 4.75 4.5 4.25 4 3.75 3.5 3.25 3
6.19 6.179 6.168 6.156 6.127 6.11 6.088 6.058 6.054 6.05
3.276 3.272 3.262 3.255 3.245 3.229
3.73 3.717 3.697 3.694 3.691 3.688 3.684 3.674 3.667 3.659 3.645
3.32 3.307 3.286 3.283 3.28
0.100 0.200 0.300 0.400 0.600 0.700 0.800 0.900 0.910 0.920 0.930 0.940 0.960 0.970 0.980 0.990 0.991 0.992 0.993 0.994 0.996 0.997 0.998 0.999
[10/24, 7/24, 7/24] 3.896 3.888 3.879 3.869 3.845 3.829 3.808 3.776 3.772 3.767 3.761 3.755 3.74
[1/3, 1/3, 1/3]
6.25 6 5.75 5.5 5.25 5 4.75 4.5 4.25 4 3.75 3.5 3.25 3
ε = diam(A) /200
(b) ε = diam(A) /400
(a)
Numer of steps
Numer of steps
760 T. Martyn / Computers & Graphics 26 (2002) 753–764
T. Martyn / Computers & Graphics 26 (2002) 753–764
761
(a) = diam(A ) /100 50 000 000 47 500 000 45 000 000 42 500 000 40 000 000 37 500 000 35 000 000 Numer of steps
32 500 000 30 000 000 27 500 000 25 000 000 22 500 000 20 000 000 17 500 000 15 000 000 12 500 000 10 000 000 7 500 000 5 000 000 2 500 000 0 0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
[0.85, 0.07, 0.07, 0.01]
32793469 33605931 34530466 35602082 36875166 38441062 40471503 43353455 43828176
[0.73, 0.11, 0.13, 0.03]
24231552 2463955
[0.714, 0.1099, 0.119, 0.068] 1770396
1799714
2510400
2564259
2628294
2707148
2809577
2955377
3208348
1833099
1871829
1917896
1974652
2048414
2153465
2335844
0.7
0.8
0.9
Probability
(b) = diam(A ) /200 400 000 000 375 000 000 350 000 000 325 000 000 300 000 000 Numer of steps
275 000 000 250 000 000 225 000 000 200 000 000 175 000 000 150 000 000 125 000 000 100 000 000 75 000 000 50 000 000 25 000 000 0
0.1
0.2
0.3
0.4
0.5
0.6
[0.85, 0.07, 0.07, 0.01]
231062775 239330767 248802636 259859264 273093458 289502282 310962523 341705842 395301940
[0.73, 0.11, 0.13, 0.03]
11409316 11580828 11775951 12002085 12270738 13030087 13030087 13639459 14693869
[0.714, 0.1099, 0.119, 0.068] 8867802
8999016
9148292
9321285
9526792
9779595 10107517 10573329 11378714
Probability
Fig. 4. The sufficient number of the Chaos Game steps to produce an approximation of Barnsley’s Fern at a given probability and error e: The first and respectively second rows of the tables present results for IFS probability vectors proposed by Barnsley and Hepting et al., respectively. The third row displays number of steps computed for the probability vector estimated by means of the method presented in this paper. The error was measured using the weighted maximum metric dðx; yÞ ¼ jx yj; jxj :¼ maxð13=28xð1Þ ; xð2Þ Þ:
Taking advantage of this fact allows computing to be done efficiently by estimating results for high approximation accuracy on the basis of the outcomes obtained for a lower one. It is of the key importance especially if we search for the minimum number of steps required to attain given approximation accuracy with given probability fixed. There is also another thing to notice, namely, a conjecture between the number of steps for corresponding probabilities in the neighbouring tables in Fig. 2. Assuming fixed IFS probabilities, let nðp; eÞ denote the
sufficient number of steps determined by means of our method to obtain the attractor approximation with an error pe and probability p: The subject of interest is the function f : ½0; 1Þ ð0; 1Þ-R defined by f ðp; eÞ :¼
nðp; eÞ : nðp; 2eÞ
ð18Þ
First of all, numerical experiments allow us to believe that for a fixed e > 0; f is a decreasing function of probability p; see Fig. 3. Furthermore, the values of the function seem to converge to a certain limiting value as p
762
T. Martyn / Computers & Graphics 26 (2002) 753–764
tends to one. Moreover the limit appears not to depend on the e chosen. Looking at the ½13; 13; 13 rows of the tables (a) and (b) listed in Fig. 3, it is hard to resist from the conclusion that for the equal IFS probabilities, the function f ‘‘tries to mimic’’ the ratio of the Ddimensional Hausdorff measures H D ðÞ of Sierpinski’s Triangle W and its double reduced copy 12W, for H D ðnÞ ¼ 3; f ðp; eÞE D H ð1=2nÞ
ð19Þ
where D ¼ log 3=log 2 is the Hausdorff dimension of the fractal. (For the scaling property of the Hausdorff measures see e.g. [7, p. 27].) Although the results shown in the tables cannot be considered as a proof that for every e > 0; limp-1 f ðp; eÞ ¼ 3; they suggest that such a hypothesis is probable. This is so especially if we consider the function nðp; eÞ as a measure of a set generated by the Chaos Game to obtain the attractor approximation of parameters p and e; then nðp; 2eÞ is the corresponding measure of the set scaled by the factor 12: However, even if our supposition turns out true, one should ask what the assumed limits are in the remaining cases in which the IFS probabilities are not equal. And
again one could speculate in a similar fashion that they are ratios of the corresponding measures induced by the distribution of IFS probabilities (see e.g. [8]). However, the problem appears to deserve a deeper and more formal analysis, so we leave it for the further research. The second example is Barnsley’s Fern (see appendix) that was tested for the original IFS probabilities given in [2] and the ones proposed in [1]. The selected values of e suits visualization at 47 100, 93 200. The results presented in Fig. 4 can be regarded as a confirmation of the thesis by Hepting et al. [1] that their IFS probabilities are really improved in the sense of minimization of the sufficient number of the algorithm’s steps as computed by means of our method. However, it should be noted that the problem whether the number of the Chaos Game steps specified in this way could be considered as an objective measure of the quality of IFS probabilities (with respect to given approximation accuracy) remains open. Nevertheless, we tried to determine an IFS probability distribution that would result in a less number of steps than computed with our algorithm for the distribution proposed in [1]. We found that, in the context of our method, the vector [0.714,
Fig. 5. Rendering of p-balanced measures supported by Barnsley’s Fern. The measures are generated by the following IFS probability vectors: (a) Barnsley’s [0.85, 0.07, 0.07, 0.01], (b) Hepting et al. [0.73, 0.11, 0.13, 0.03], (c) ours improved [0.714, 0.099, 0.119, 0.068]. The same number of the Chaos Game steps was performed to render each fern.
T. Martyn / Computers & Graphics 26 (2002) 753–764
0.099, 0.119, 0.068] seems to be a better candidate for the Fern’s probability distribution than the ones proposed previously (see Fig. 4). The results of visualization of the p-balanced measures induced by the three discussed probability distributions are displayed in Fig. 5. On the basis of the presented images it seems that the new probability vector offers more uniform distribution of points than the previous ones. There are also many other problems open for further research. Some of them follow: *
*
*
Our method estimates probability for the Chaos Game accuracy as the sum of probabilities of the occurrence of each string sAD separately, paying no attention on probabilities of the string occurrence in common. How accurate is such estimation in general case? Can our algorithm be improved by taking into consideration probability of the co-occurrence of strings without a significant increase of the algorithm’s computational complexity? Our algorithm works on the basis of the IFS specification only, regardless of the resulting attractor’s geometry. However, in the case of attractor’s subsets overlapping, there can be more than one string (address) assigned to a given subset [2]. In such cases it is sufficient to consider only one representative for each of the overlapping subsets. Taking advantage of this fact would improve our algorithm both in the accuracy of estimation and time complexity. If the number of the Chaos Game steps that is determined by our method turns out to be a correct measure of the quality of IFS probabilities, a research concerning an efficient optimization algorithm should be taken. At this stage it seems that further work on application of heuristic approaches such as genetic algorithms could prove to be fruitful both in the practical and theoretical sense.
Game will produce, within L of its steps, a sequence that approximates the IFS attractor at a given value of accuracy. As a result, we can determine the sufficient number of the Chaos Game steps to approximate the attractor for given accuracy and probability. However, whether this number can be identified with a measure of the quality of the probability vector assigned to the IFS, remains an open problem. We also outlined a number of other noteworthy problems and questions that have arisen from the analysis of our method and results, which as we believe can be a stimulus for further research.
Acknowledgements I would like to thank Dr. Andrzej Paja) k for fruitful discussions on Markov chains during the development of this method, and to Krzysiek Gracki for his help in production of the figures and patient listening to my tales on the fractals’ and other strange worlds. The shape of this paper could not be reached without suggestions by Prof. Jan Zabrodzki. My special thanks go to the C&G reviewers whose encouraging and valuable comments help me to improve the original version of this paper.
Appendix Here is the IFS code for the Barnsley Fern on the basis of which computation was done. " # " # 0:85 0:04 0:0 w1 ðxÞ ¼ xþ ;
0:04 0:85 1:6 " # " # 0:2 0:26 0:0 xþ ; w2 ðxÞ ¼ 0:23 0:22 1:6 "
8. Summary The IFS representation is a powerful and flexible tool for description of the fractal shapes. The most popular algorithm for rendering IFS attractors is the Chaos Game. As far as the infinite number of iteration steps is considered, the algorithm produces an infinite sequence of points that, with probability one, approximates the IFS attractor at any accuracy. However, in real applications only a finite number of points can be generated, so it is often necessary to know how many iteration steps should be performed to produce an approximation at a given accuracy. The method presented in this paper makes it possible to answer this question. More precisely, given an IFS with probabilities, our algorithm estimates probability that the Chaos
763
# " # 0:28 0:0 xþ ; w3 ðxÞ ¼ 0:24 0:44 " # " # 0:0 0:0 0:0 xþ : w4 ðxÞ ¼ 0:0 0:16 0:0
0:15 0:26
References [1] Hepting D, Prusinkiewicz P, Saupe D. Rendering methods for iterated function systems. Proceedings of Fractals ’90, 1990, IFIP. [2] Barnsley MF. Fractals everywhere, 2nd ed. Boston: Academic Press, 1993. [3] Goodman GS. A probabilist looks at the Chaos Game. Proceedings of Fractals ’90, 1990, IFIP.
764
T. Martyn / Computers & Graphics 26 (2002) 753–764
[4] Knuth DE, Morris JH, Pratt VR. Fast pattern matching in strings. SIAM Journal on Computing 1977;6(2):323–50. [5] Cormen ThH, Leiserson ChE, Rivest RL. Introduction to algorithms. MA: MIT, 1990. [6] Barnsley MF. The desktop fractal design handbook. New York: Academic Press, 1989.
[7] Falconer K. Fractal geometry. mathematical foundations and applications. Chichester: Wiley, 1990. [8] Demko S, Hodges L, Naylor B. Construction of fractal objects with iterated function systems. Computer Graphics (SIGGRAPH ’85) 1985;19(3):271–8.
Technical section
On approximation accuracy of the Chaos Game’s finite-time activity Tomek Martyn* Computer Graphics Laboratory, Institute of Computer Science, Warsaw University of Technology, ul. Nowowiejska 15/19, 00-665, Warsaw, Poland
Abstract This paper describes a method for estimating probability that the Chaos Game generates, within L steps, a sequence of points which approximates the IFS attractor at a given accuracy. The presented algorithm is a hybrid based on methods from the theory of Markov chains used in the pattern-matching context. Unlike the previous approaches, ours is general both in the sense of component IFS mappings and the associated probability distribution. r 2002 Elsevier Science Ltd. All rights reserved. Keywords: Iterated function system; Fractals; IFS attractors; The Chaos Game; Probability
1. Introduction The iterated function system (IFS) is the most popular way to describe fractal shapes. Each IFS specifies implicitly the set of points called the IFS attractor. In order to obtain an explicit image of the attractor you can use one of the known methods of IFS rendering [1]. Because of relatively easy implementation, the method of choice is usually the random iteration algorithm [2] that is known wider as the Chaos Game. As far as the infinite number of iteration steps is considered, the algorithm almost surely (i.e., with probability one) generates an infinite sequence of points which, ignoring a finite number of initial points, approximates the attractor with an arbitrarily small error. The formal proofs of this fact usually do not take into account the results of the algorithm’s finite-time activity. Since in practice only a finite number of points can be generated, knowledge of the asymptotic behaviour of the process is very interesting from the mathematical point of view, but as it concerns the Plato world exclusively, it says nothing about the world’s merely shadows that we are
*Tel.: +48-22-285-31-79; fax: +48-22-285-16-35. E-mail address: [email protected] (T. Martyn).
given. Therefore, knowing all those things about infinity but being a material person living in the ordinary nonperfect world, one could ask: How many points have to be generated to approximate an IFS attractor at a given accuracy? As far as you use the Chaos Game just to have a glance at fractals for pleasure, trying to answer this question may seem another whim of a mathematician, because you are always able to generate more points to refine an image. However, when the Chaos Game is regarded as a raster conversion function for purposes of visualization at pixel size accuracy, or if the stream of points generated by the algorithm is directed to the input of a computational geometry procedure (e.g. for determination of the minimum ball enclosing an IFS attractor), then the stated problem is of key importance. As a previous work concerning, in many aspects, the above problem one should mention the paper by Goodman [3] devoted to the analysis of the Chaos Game in its original version (i.e., for generating Sierpinski’s Triangle with the uniform probability distribution). The approach presented in this paper has its roots in Goodman’s analysis but we extend the scope to the general case in which an IFS can consist of any type of contraction mappings equipped with an arbitrary assigned probability distribution.
0097-8493/02/$ - see front matter r 2002 Elsevier Science Ltd. All rights reserved. PII: S 0 0 9 7 - 8 4 9 3 ( 0 2 ) 0 0 1 3 0 - 9
T. Martyn / Computers & Graphics 26 (2002) 753–764
754
The diameters of these subsets satisfy
2. IFS and the Chaos Game An iterated function system ðX ; w1 ; y; wN Þ; shortly IFS, is defined as a finite set of contraction mappings wi : X -X ; i ¼ 1; y; N; on a complete metric space ðX ; dÞ: Using the Banach fixed point theorem one can prove that each IFS implicitly specifies a compact (i.e., totally bounded and closed) and nonempty set A; which is the unique solution of the equation A¼
N [
wi ðAÞ:
ð1Þ
i¼1
Such a set is called the attractor of the IFS. The assignation of probabilities to IFS mappings leads to the concept of an IFS with probabilities ðX ; wP where pi > 0; i ¼ 1; y; N; 1 ; y; wN ; p1 ; y; pN Þ; and N i¼1 pi ¼ 1: Having such an IFS, we can approximate its attractor by means of the random iteration algorithm [2], known wider as the Chaos Game. The algorithm works in an iterative way and produces a sequence of points fxk gN k¼0 such that x0 is any point of X and xkþ1 ¼ wik ðxk Þ; where wik is chosen independently in each step with the probability pik from the set of the IFS mappings. The sequence fxk gN k¼0 approximates the attractor with an arbitrarily small approximation error e > 0: If the starting point x0 is in the attractor, every point of the sequence belongs to the attractor as well. The initial point can be easily determined as the fixed point x ¼ wi ðxÞ of any IFS mapping. In this case it is possible to achieve the best approximation of an attractor with the least number of iteration steps, which is desirable for computer graphics purposes. Therefore, our attention will focus on the stated case. However, if necessary, the considerations presented in this paper can be easily extended to the general case in which x0 is any point of X :
3. IFS attractor structure Let A be the ðX ; w1 ; y; wN ; p1 ; y; pN Þ: attractor A consists of N Since the IFS mappings wi of each subset satisfies
attractor of an IFS On the basis of (1), the subsets wi ðAÞ; iAf1; y; Ng: are contractive, the diameter
diamðwi ðAÞÞpsðwi ÞdiamðAÞ;
8iAf1; y; Ng;
ð2Þ
where sðÞ denotes the contractivity factor of a given mapping. In turn, each subset wi ðAÞ is made up of N subsets of the form wi 3wj ðAÞ; jAf1; y; Ng; for ! N N [ [ wj ðAÞ ¼ wi 3wj ðAÞ: ð3Þ wi ðAÞ ¼ wi j¼1
j¼1
diamðwi 3wj ðAÞÞpsðwi 3wj ÞdiamðAÞ psðwi Þsðwj ÞdiamðAÞ:
ð4Þ
Generally, by induction, any subset wi1 3wi2 3?3win ðAÞ; ik Af1; y; Ng; of the attractor consists of N subsets y; Ng; whose diawi1 3wi2 3?3win 3winþ1 ðAÞ; inþ1 Af1; Qnþ1 meters are not greater than k¼1 sðwik Þ diamðAÞ: The contractivity factors sðwik ÞA½0; 1Þ; ik Af1; y; Ng; so given e > 0; we are able to decompose the attractor to a finite number of subsets wi1 3?3win ðAÞ; whose diameters do not exceed e: Such a family of the attractor subsets we will call an e-decomposition of the attractor. Choosing at least one arbitrary point from each subset of the e-decomposition we can construct a finite set Ae such that the Hausdorff distance hðA; Ae Þpe; i.e., the set Ae approximates A with an error not greater than e:
4. Estimating probability Let e > 0 be given and let aAA be an arbitrary point of the attractor. Then a belongs to at least one of the e-decomposition subsets; say aAwa1 3?3wan 1 3wan ðAÞ; ai Af1; y; Ng: Assume the initial point for the Chaos Game is x0 AA: Then every point generated by the algorithm also belongs to the attractor, because from (1) the IFS mappings take A into itself. The n-fold composition wa1 3?3wan maps the attractor A to its subset whose diameter is not greater than e and which includes the point a: Therefore, appearing of the sequence wan ; y; wa1 in arbitrary n successive steps of the algorithm is the sufficient condition for at least one of the points produced by the algorithm to be within the distance e to a: Each mapping sequence wan ; y; wa1 is specified by the string of indices an an 1 ya1 (called in the IFS nomenclature an address of the subset wan ; y; wa1 ðAÞ), so given e > 0; we can restrict our considerations to the strings from the set that is induced by the P e-decomposition of the attractor. Let ¼ f1; y; Ng be P the set of the IFS map indices (or symbols), P and let denote the set of all (finite) strings over : Let fp1 ; y;P pN g be a set of IFS probabilities. Denote by DC the set of strings, which is induced by the attractor e-decomposition, i.e., if an an 1 ya1P AD; then diamðwa1 3?3wan ðAÞÞpe: Let foi gLi¼1 ; oi A ; and respectively fxi gLi¼1 ; xi AA; be a string of indices, and a sequence of points, respectively, produced in L initial steps of the algorithm. Then fxi gLi¼1 approximates the attractor with an error not greater than e, if foi gLi¼1 includes all strings from D as its substrings. Just like in the case of Goodman’s analysis of the Chaos Game for Sierpinski’s Triangle [3], the process of generating the strings from D within foi gLi¼1 possesses the Markov property, because knowledge of the present
T. Martyn / Computers & Graphics 26 (2002) 753–764
substring completely determines probability for the next one to arise, without regard to the prior history of the process. The main differences are that, in our general case, probabilities do not have P to be distributed uniformly on the symbols from ; and the strings from D do not have to be of equal lengths. Let us denote by ‘!’ and respectively by ‘g’ P the P prefix and suffix relations on ; i.e., for x; yA P ; y!x and respectively ygx; iff there is some wA such that x ¼ y3w and x ¼ w3y; respectively, where ‘3’ denotes thePconcatenation of strings. Now, defining Dp :¼ fxA : x!s for some sADg and identifying the set S ¼ D,Dp with the space of states, we can consider the action of the Chaos Game within the category of a Markov chain, i.e., a sequence fXi gN i¼0 of random variables with values in S; such that for every sequence fsi gli¼0 ; si AS; we have PðXi ¼ si jXi 1 ¼ si 1 ; y; X1 ¼ s1 ; X0 ¼ s0 Þ ¼ PðXi ¼ si jXi 1 ¼ si 1 Þ:
ð5Þ
(Note that the empty string l; i.e., consisting ofPno symbols, is a prefix and suffix of every string from ). The Markov chain proceeds from state to state as follows: P P Let Xi 1 ¼ a3s; where aA and sAS; and let bA be the symbol chosen in the ith step (i > 0); then X0 ¼ l ( Xi ¼
s3b; a3s3bAS;
if s3bAS; otherwise;
ð6Þ
where probability of a transition from the state to one of the states is equal to probability of choosing the symbol b: Thus, in order to estimate probability that the Chaos Game generates, in L initial steps, the sequence of points that approximates the attractor at a given accuracy, all we have to do is determine probability that the Markov chain visits, in L steps, all states from D: Pð8sAD; (kpL; Xk ¼ sÞ ¼ P
L \[
! fXk ¼ sg
sAD k¼1
¼1 P
!
L [\
fXk asg :
sAD k¼1
ð7Þ Using Bool’s inequality, we obtain P
L [\
! fXk asg p
sAD k¼1
¼ jDj
X sAD
X sAD
P
L [ k¼1
P
L \
! fXk asg
k¼1
!
fXk ¼ sg ;
ð8Þ
755
and hence Pð8sAD; (kpL; Xk ¼ sÞX ! L X [ 1þ P fXk ¼ sg jDj; sAD
ð9Þ
k¼1
where jDj denotes the number of states in D: Thus, in order to achieve our goal, we have to
S determine for each sAD the value P Lk¼1 fXk ¼ sg ; that is probability that the Markov chain fXi gN i¼0 will visit, in L initial steps, the state s at least once.
5. Probability of generating strings It is well-known fact that probability of the occurrence of a substring sAS within a random sequence P depends not only on the probability of foi gLi¼1 over the component symbols to appear but also on the structure of s itself. For example, in the case P when probabilities are distributed uniformly on ; the occurrence of a substring that consists of a single repeated symbol is a less probable event than the occurrence of any string of the P same length but composed of varied symbols of : The problem considered here resembles in many aspects the one of pattern matching. In fact, the solution presented further can be regarded as an adaptation of the method proposed by Knuth et al. [4] for searching for a string within a given text by means of finite automata (for a good discussion on the subjectPsee also [5]). Let a ¼ a1 a2 ; y; an AD; ai A ; and let aðkÞ denote the a string’s prefix of the length kA½0; n : Assuming iXn; it is easy to see that the following identities are satisfied: ð1Þ Xi ¼ a3aðn 1ÞgXi 1 and an is chosen in the ith step; ð2Þ aðn 1ÞgXi 1 3aðn 2ÞgXi 2 and an 1 is chosen in the ði 1Þth step; ^ ð3Þ að1ÞgXi nþ1 3lgXi n and a1 is chosen in the ði n þ 1Þth step:
ð10Þ
From the above sequence of implications we conclude that probability for a given string x to give rise to the string a depends only on the maximal suffix of x; which matches a prefix of a; no matter what the x string’s remaining (i.e., beyond the suffix) initial symbols are. Thus, if x; yAD are strings with suffixes fitting the a prefix at the same positions, then they are indistinctive states from the point of view of conditional probability for the Markov chain to visit the state a: Therefore, we can construct a Markov chain with the space of states consisting of subsets of the space S; which are induced
T. Martyn / Computers & Graphics 26 (2002) 753–764
756
with transition probability specified by
by a suffix equivalence relation. The specification of the process affects the statistical behaviour of this chain to be identical to that one of fXi gN i¼0 in the sense of probability of the string a to arise in the random sequence foi gLi¼1 : Moreover, such an approach will allow us to solve the considered problem by means of standard techniques from the theory of Markov chains used in the context of pattern matching.
PðYi ¼ Qm jYi 1 ¼ Qk Þ 8 > < 1; ¼ 0; > : PðsðXi Þ ¼ mjsðXi 1 Þ ¼ kÞ;
Let s : S-N,f0g be the suffix function of the string a; i.e., ð11Þ
Now, we can introduce a relation, denoted by ‘B’, on S such that xBy :¼ ðsðxÞ ¼ sðyÞÞ:
ð12Þ
The relation satisfies properties of equivalence relation. Hence, we immediately get the partition Q ¼ fQ0 ; Q1 ; y; Qn g of S with respect to B; where Qi ¼ ½x is the equivalence class of x; such that i ¼ sðxÞ: Let fYi gN i¼0 be a sequence of random variables taking values from the set Q in the following way: P Let Yi 1 ¼ Qk and let bA be a symbol chosen in the ith ði > 0Þ step; then Y 0 ¼ Q0 ; ( Y i ¼ Qm ; m ¼
PðYl ¼ Qn Þ ¼ Pð(kpL; Xk ¼ aÞ; if k ¼ n otherwise
n; sðaðkÞ3bÞ;
(a) Q0
β Q1
β
α Q2
Q3
α
(b) State Q0 Q1 Q2 Q3 Q4 Q5 Q6 Q7 Q8
α Q4
χ
β Q6
Q5
α
α, β , χ
χ Q7
Q8
α
α
α
ð15Þ
i.e., probability that fYi gN i¼0 ; after L steps, is in the state Qn is equal to probability that fXi gN i¼0 visits the state a
ð13Þ
α
ð14Þ
From the above definition the sequence fYi gN i¼0 satisfies the properties of a Markov chain. Moreover, the chain is time homogenous, since the transition probability matrix induced by (14) is the same for every Yi ; i.e., it does not depend on the number of steps. The conditional probability PðsðXi Þ ¼ mjsðXi 1 Þ ¼ kÞ is specified by P the sum of probabilities associated with symbols bA that map a state Xi 1 ¼ x3aðkÞ to a state Xi ¼ y3aðkÞ3b ¼ z3aðmÞ; x; y; zAS; according to the transition rules (6). Note that the properties of suffix matching imply that there is usually at most one symbol taking a given state for which sðXi 1 Þ ¼ k to another one such that sðXi Þ ¼ m; see Fig. 1. The exceptions are, from definition, the transitions from the final state Qn ; as well as the transitions to the states with sðX Pi Þ ¼ 0 and the ones for which there is no symbol in to take a state to another. Obviously, in the last case probability is equal to zero. The definition of fYi gN i¼0 distinguishes the state Qn as the unique fixed point of the chain in the sense that for every kXi; Yi ¼ Qn implies Yk ¼ Qn : Since Qn ¼ fag; therefore YL ¼ Qn if and only if Xi ¼ a for some ipL: Putting the above observations in the context of the sequence of identities (10), we conclude that
5.1. The equivalent Markov chain
sðxÞ :¼ maxfk : aðkÞgxg; xAS:
if k ¼ n ¼ m; if k ¼ nam; if kan:
(c) α Q1 Q1 Q3 Q1 Q5 Q1 Q5 Q1 Q8
β Q0 Q2 Q0 Q4 Q0 Q6 Q0 Q0 Q8
χ Q0 Q0 Q0 Q0 Q0 Q0 Q7 Q8 Q8
Q0 Q1 Q2 Q3 Q4 Q5 Q6 Q7 Q8
Q0 β+χ χ β+χ χ β+χ χ β β 0
Q1 α α 0 α 0 α 0 α 0
Q2 0 β 0 0 0 0 0 0 0
Q3 0 0 α 0 0 0 0 0 0
Q4 0 0 0 β 0 0 0 0 0
Q5 0 0 0 0 α 0 α 0 0
Q6 0 0 0 0 0 β 0 0 0
Q7 0 0 0 0 0 0 χ 0 0
Q8 0 0 0 0 0 0 0 χ 1
Fig. 1. The Markov chain induced by the string abababww: (a) the graph of transitions (the transitions to the state Q0 have been omitted); (b) the transition table; (c) the transition probability matrix—the component symbols are identified with the associated values of probabilities.
T. Martyn / Computers & Graphics 26 (2002) 753–764
within L steps at least once. As we stated above, fYi gN i¼0 is a time homogenous Markov chain, so in order to compute this probability we can use a standard technique from the relevant theory of such processes. Having a ðn þ 1Þ ðn þ 1Þ transition probability matrix M; the probability distribution for the variable YL is specified by the vector y ¼ x ML ; x; yARnþ1 ; where the ith coordinate of the vector y is probability of the event YL ¼ Qiþ1 ; and x specifies an initial distribution for Y0 ; in our case x ¼ ½1; 0; y; 0 : Because the exponent L as the number of the Chaos Game steps is supposed to be a relatively large number, it is convenient to use the following algorithm of fast matrix involution that raises a n n matrix M to a given power L in time Oðn3 log LÞ: out’identity matrix repeat { if L mod 2=1 then out’out M M’M M L’Lb1 } until L ¼ 0
757
iSþ 1:
Because, from Proposition 1, the set 1q aðp ðkÞÞ is composed of all the proper qAN suffixes aðiÞ of aðkÞ; we can determine the value pðk þ 1Þ on the basis of the values p1q ðkÞAfpðkÞ; pðk
1Þ; y; pð1Þg computedSin the previous steps. We do so by searching the set qAN aðp1q ðkÞÞ for the prefix aðiÞ such that aiþ1 ¼ akþ1 ; where iAf0; y; k 1g is maximal. If there is such a string aðiÞ in the set, then pðk þ 1Þ ¼ i þ 1; otherwise pðk þ 1Þ ¼ 0: Furthermore, note that for any i; mAN such that iom; the length of the prefix aðp1i ðkÞÞ is not greater than the length of aðp1m ðkÞÞ: S Thus, in practice, we do not have to deal with the set qAN aðp1q ðkÞÞ of prefixes directly to compute pðk þ 1Þ; but we just search a finite set f1; y; mg; mpk; for the minimal q such that the ðk þ 1Þth and ðp1q ðkÞ þ 1Þth symbols of the string a are equal. As a result we have the following algorithm for the determination of the value pðkÞ on the basis of values pðiÞ; iAf1; y; k 1g:
5.2. Determination of the transition probability matrix
let a1 yak ADp be a prefix of aAD; then if ko2 then pðkÞ’0 else { i’pðk 1Þ while i > 0 and aiþ1 aak do i’pðiÞ if aiþ1 ¼ ak then i’i þ 1 pðkÞ’i }
The last remaining problem to solve is the determination of the transition probability matrix. In order to make it in an efficient way our approach takes advantage of the so-called prefix function introduced by Knuth et al. for their method of pattern matching [4]. For a given string a ¼ a1 a2 yan AD; the prefix function p : f0; 1; y; ng-f0; 1; y; n 1g is defined as
Proposition 2. Let a ¼ a1 a2 yan AD, and let s and p be the suffix and prefix function ofP a, respectively. Then, for every kAf0; y; n 1g and bA 8 if k ¼ 0 and baa1 ; > < 0; sðaðkÞ3bÞ ¼ k þ 1; if b ¼ akþ1 ; ð17Þ > : sðaðpðkÞÞ3bÞ; if k > 0 and baakþ1 :
where ‘b’ denotes the right-bit-shift operation in the C language convention.
pð0Þ ¼ pð1Þ ¼ 0; pðqÞ :¼ maxfk : koq and aðkÞgaðqÞg;
ð16Þ
i.e., for q > 0; pðqÞ is the maximum length of a prefix of the string a, which is simultaneously a proper suffix of aðqÞAS: Proposition 1. Let a ¼ a1 a2 yan AD be a string with the prefix function p: Denote by p1q the q-fold S composition of p: Then, for a given kAf1; y; ng; the set qAN aðp1q ðkÞÞ consists of all the strings, which being of the form aðiÞ, iok, are simultaneously the proper suffixes of the prefix aðkÞ. Proof. See e.g. [5].
&
It is easy to see that aði þ 1Þgaðk þ 1Þ if and only if aðiÞgaðkÞ and aiþ1 ¼ akþ1 : Hence, if aðiÞ is the proper suffix of aðkÞ such that the length of aðiÞ is maximal in the sense of the condition aiþ1 ¼ akþ1 ; then pðk þ 1Þ ¼
Proof. The first and second equalities are obvious. To prove the last one, we will show that if kAf0; y; n 1g and baakþ1 then aðiÞgaðkÞ3b3aðiÞgaðpðkÞÞ3b: Hence, from the suffix function definition, we will get the thesis. First of all, let us note that the suffix relation is transitive. Furthermore, it is not difficult to see that if xgz; ygz; and x is longer than y, then ygx: Let aðiÞgaðkÞ3b: Since baakþ1 ; we have aðiÞaaðkÞ3b; hence aðiÞ is a proper suffix of aðkÞ3b: From the definition of the function p; aðpðkÞÞ is the proper suffix of aðkÞ so that the length of aðpðkÞÞ is maximal. Therefore aðpðkÞÞ3b is the maximal proper suffix of the string aðkÞ3b: Since aðiÞ is a proper suffix of aðkÞ3b; using our observations on the suffix relation, we conclude that aðiÞgaðpðkÞÞ3b: Conversely, suppose that aðiÞgaðpðkÞÞ3b: As we note above aðpðkÞÞ3bgaðkÞ3b; and the suffix relation is transitive. Hence aðiÞgaðkÞ3b: &
758
T. Martyn / Computers & Graphics 26 (2002) 753–764
Using the above observations, we can construct an efficient algorithm for the successive determination of the transitions (13) and their probabilities (14). Given a state Qk ; in order to compute transitions, the algorithm uses information on transitions determined for one of the previous states Qi ; iok: Moreover, the algorithm works in on-line manner, i.e., to make the computation concerning the state Qk ; it does not require knowledge of the further string suffix symbols am nor information about the states Qm ; kom: Consequently, the transition probability matrices can be determined simultaneously together with constructing the individual strings from the set of the states D during the e-decomposition process. As a result, we have the following algorithm for the determination of the vector M½k; : that specifies probabilities of the transitions from the state Qk to others. Given a string symbol akþ1 ; the algorithm calculates M½k; : using simultaneously determined vaP lues of the transition function s½k; b ; bA ; until the input symbols are exhausted: let Qk be the current state; then for each state Qm do M½k; m ’0 P if there is a symbol akþ1 A on the input then { calculate theP value p½k for each bA do { if b ¼ akþ1 then s½k; b ’k þ 1 else if k ¼ 0 then s½k; b ’0 else s½k; b ’s½p½k ; b M½k; s½k; b
’M½k; s½k; b
þ p½b } } else M½k; k ’1 where p½b is probability associated with the symbol b (or equivalently with the IFS mapping wb ).
6. The algorithm—summary Let ðX ; w1 ; y; wN ; p1 ; y; pN Þ be an IFS with the attractor A: The heart of the algorithm is the following procedure that determines the sum of
S P recurrently L probabilities during the esAD P k¼1 fXk ¼ sg decomposition process: procedure Descent(MapComp, k) { P for each aA do { Calculate M[k; . ] with the symbol a on the input. Determine a new composition of IFS mappings and save it in the local variable NewMapComp. NewMapComp’MapComp3wa if diamðNewMapCompðAÞÞpe then { Here is the end of the current string,
so calculate M½k þ 1; : for the empty input. Increment the global variable NumStrings that represents the number of strings in the set of states D: NumStrings’NumStrings+1 Finally, calculate probability the Markov chain visits the string (state), and then use that value to increment the global variable Prob Prob’Prob þ ML ½0; k þ 1 } else Descent(NewMapComp, k þ 1) } } P where ¼ f1; y; Ng; and L is the number of the Chaos Game steps. Now, armed with the above procedure, we calculate the estimated probability (9) as follows: let map be equal to the identity mapping Descent(map, 0) Probability’1+Probability NumStrings.
7. Results and open problems Two classic IFS attractors were tested for the sufficient number of the Chaos Game steps to obtain the approximation of an attractor at given accuracy and probability. The first remark is that the results obtained by means of our approach appear to substantiate the recommendation by Barnsley to continue iterating millions of times [6, p. 28], especially when one looks at the numbers of steps calculated for Barnsley’s Fern (see Fig. 4). Also in the case of Sierpinski’s Triangle W with equal probabilities the results agree in a sense with ones determined by Goodman [3]. In the latter case, the main reason of the existing differences in number of steps follows from the fact that, for each string, our method calculates probability of the string to occur in the accurate way, whereas Goodman estimates it on the basis of the worst case in which a string consists of a single repeated symbol. As a consequence, our method is more accurate which manifests itself in a fewer number of steps than presented in [3]. For example, Goodman estimates that for the approximation with an error not greater than diamðDÞ2 9 (the length of strings in D is 9), ‘‘the chance that more than 656,829 are required is less 8.44%’’, whereas given such the approximation error and probability 100%–8.44%=91.56%, our algorithm results in 221,825 steps. Fig. 2 presents the results that our method supplied for Sierpinski’s Triangle. The selected values of e suits visualization at 100 100, 200 200, and 400 400
T. Martyn / Computers & Graphics 26 (2002) 753–764
759
Numer of steps
(a) ε = diam(A ) /100 275 000 250 000 225 000 200 000 175 000 150 000 125 000 100 000 75 000 50 000 25 000 0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
14709
14965
15256
15591
15989
16476
17105
17993
19515
[10/24, 7/24, 7/24]
18230
18621
19067
19585
20204
20967
21963
23385
25863
[1/2, 1/3, 1/6]
138593 144109 150529 158163 167499 179377 195422 219419 264102
[1/3, 1/3, 1/3]
Probability
(b) ε = diam(A ) /200
Numer of steps
1 650 000 1 500 000 1 350 000 1 200 000 1 050 000 900 000 750 000 600 000 450 000 300 000 150 000 0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
51258
52028
52901
53910
55104
56566
58452
61113
65667
[10/24, 7/24, 7/24]
71026
72398
73963
75780
77945
80615
84091
89048
97659
[1/2, 1/3, 1/6]
857834 890499 928507 973678 1028898 1099126 1193951 1335737 1599805
[1/3, 1/3, 1/3]
Probability
(c) ε = diam(A ) /400 10 000 000 9 000 000 Numer of steps
8 000 000 7 000 000 6 000 000 5 000 000 4 000 000 3 000 000 2 000 000 1 000 000 0 [1/3, 1/3, 1/3]
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
175262 177575 180198 183227 186811 191198 196856 204834 218483
[10/24, 7/24, 7/24] 272901 277702 283172 289522 297077 306389 318491 335722 365587 [1/2, 1/3, 1/6]
5288610 5482463 5707940 5975800 6303110 6719206 7280805 8120229 9683557 Probability
Fig. 2. The sufficient number of the Chaos Game steps to produce an approximation of Sierpinski’s Triangle at a given probability and error e: The results are listed for three different IFS probability vectors. The error was measured using the maximum metric.
screen resolution, respectively. In order to get some view how IFS probabilities influence the number of the Chaos Game’s steps, three different distributions of probabilities have been tested. In this case, as one might have expected, the uniform distribution results in the
least number of steps. This is no surprise, however, the noteworthy thing is that the graphs of the presented results look very similar to each other, with no regard of approximation accuracy. (The same can be noticed for Barnsley’s fern; Fig. 4).
3.485 3.477 3.468 3.458 3.433 3.417 3.396 3.365 3.361 3.356 3.351 3.345 3.33
[1/2, 1/3, 1/6]
[1/2, 1/3, 1/6]
Fig. 3. The values of the function (18) (see text) evaluated for Sierpinski’s Triangle.
Probability
6.165 6.157 6.147 6.137 6.113 6.098 6.079 6.053 6.05 6.046 6.042 6.038 6.028 6.023 6.016 6.009 6.008 6.007 6.006 6.005 6.004 6.003 6.002 6.001
3.642 3.63
3.419 3.413 3.406 3.399 3.38 3.368 3.352 3.327 3.324 3.32 3.316 3.311 3.299 3.291 3.281 3.264 3.262 3.259 3.256 3.253 3.245 3.239 3.231 3.219
0.100 0.200 0.300 0.400 0.600 0.700 0.800 0.900 0.910 0.920 0.930 0.940 0.960 0.970 0.980 0.990 0.991 0.992 0.993 0.994 0.996 0.997 0.998 0.999
Probability
6.045 6.041 6.03 6.024 6.017 6.009 6.008 6.008 6.007 6.006 6.004 6.003 6.002 6.001
5 [10/24, 7/24, 7/24] 3.842 3.836 3.829 3.821 3.801 3.787 3.77 3.744 3.74 3.736 3.731 3.726 3.713 3.705 3.693 3.676 3.673 3.671 3.668 3.664 3.656 3.6
[1/3, 1/3, 1/3]
6.25 6 5.75 5.5 5.25 5 4.75 4.5 4.25 4 3.75 3.5 3.25 3
6.19 6.179 6.168 6.156 6.127 6.11 6.088 6.058 6.054 6.05
3.276 3.272 3.262 3.255 3.245 3.229
3.73 3.717 3.697 3.694 3.691 3.688 3.684 3.674 3.667 3.659 3.645
3.32 3.307 3.286 3.283 3.28
0.100 0.200 0.300 0.400 0.600 0.700 0.800 0.900 0.910 0.920 0.930 0.940 0.960 0.970 0.980 0.990 0.991 0.992 0.993 0.994 0.996 0.997 0.998 0.999
[10/24, 7/24, 7/24] 3.896 3.888 3.879 3.869 3.845 3.829 3.808 3.776 3.772 3.767 3.761 3.755 3.74
[1/3, 1/3, 1/3]
6.25 6 5.75 5.5 5.25 5 4.75 4.5 4.25 4 3.75 3.5 3.25 3
ε = diam(A) /200
(b) ε = diam(A) /400
(a)
Numer of steps
Numer of steps
760 T. Martyn / Computers & Graphics 26 (2002) 753–764
T. Martyn / Computers & Graphics 26 (2002) 753–764
761
(a) = diam(A ) /100 50 000 000 47 500 000 45 000 000 42 500 000 40 000 000 37 500 000 35 000 000 Numer of steps
32 500 000 30 000 000 27 500 000 25 000 000 22 500 000 20 000 000 17 500 000 15 000 000 12 500 000 10 000 000 7 500 000 5 000 000 2 500 000 0 0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
[0.85, 0.07, 0.07, 0.01]
32793469 33605931 34530466 35602082 36875166 38441062 40471503 43353455 43828176
[0.73, 0.11, 0.13, 0.03]
24231552 2463955
[0.714, 0.1099, 0.119, 0.068] 1770396
1799714
2510400
2564259
2628294
2707148
2809577
2955377
3208348
1833099
1871829
1917896
1974652
2048414
2153465
2335844
0.7
0.8
0.9
Probability
(b) = diam(A ) /200 400 000 000 375 000 000 350 000 000 325 000 000 300 000 000 Numer of steps
275 000 000 250 000 000 225 000 000 200 000 000 175 000 000 150 000 000 125 000 000 100 000 000 75 000 000 50 000 000 25 000 000 0
0.1
0.2
0.3
0.4
0.5
0.6
[0.85, 0.07, 0.07, 0.01]
231062775 239330767 248802636 259859264 273093458 289502282 310962523 341705842 395301940
[0.73, 0.11, 0.13, 0.03]
11409316 11580828 11775951 12002085 12270738 13030087 13030087 13639459 14693869
[0.714, 0.1099, 0.119, 0.068] 8867802
8999016
9148292
9321285
9526792
9779595 10107517 10573329 11378714
Probability
Fig. 4. The sufficient number of the Chaos Game steps to produce an approximation of Barnsley’s Fern at a given probability and error e: The first and respectively second rows of the tables present results for IFS probability vectors proposed by Barnsley and Hepting et al., respectively. The third row displays number of steps computed for the probability vector estimated by means of the method presented in this paper. The error was measured using the weighted maximum metric dðx; yÞ ¼ jx yj; jxj :¼ maxð13=28xð1Þ ; xð2Þ Þ:
Taking advantage of this fact allows computing to be done efficiently by estimating results for high approximation accuracy on the basis of the outcomes obtained for a lower one. It is of the key importance especially if we search for the minimum number of steps required to attain given approximation accuracy with given probability fixed. There is also another thing to notice, namely, a conjecture between the number of steps for corresponding probabilities in the neighbouring tables in Fig. 2. Assuming fixed IFS probabilities, let nðp; eÞ denote the
sufficient number of steps determined by means of our method to obtain the attractor approximation with an error pe and probability p: The subject of interest is the function f : ½0; 1Þ ð0; 1Þ-R defined by f ðp; eÞ :¼
nðp; eÞ : nðp; 2eÞ
ð18Þ
First of all, numerical experiments allow us to believe that for a fixed e > 0; f is a decreasing function of probability p; see Fig. 3. Furthermore, the values of the function seem to converge to a certain limiting value as p
762
T. Martyn / Computers & Graphics 26 (2002) 753–764
tends to one. Moreover the limit appears not to depend on the e chosen. Looking at the ½13; 13; 13 rows of the tables (a) and (b) listed in Fig. 3, it is hard to resist from the conclusion that for the equal IFS probabilities, the function f ‘‘tries to mimic’’ the ratio of the Ddimensional Hausdorff measures H D ðÞ of Sierpinski’s Triangle W and its double reduced copy 12W, for H D ðnÞ ¼ 3; f ðp; eÞE D H ð1=2nÞ
ð19Þ
where D ¼ log 3=log 2 is the Hausdorff dimension of the fractal. (For the scaling property of the Hausdorff measures see e.g. [7, p. 27].) Although the results shown in the tables cannot be considered as a proof that for every e > 0; limp-1 f ðp; eÞ ¼ 3; they suggest that such a hypothesis is probable. This is so especially if we consider the function nðp; eÞ as a measure of a set generated by the Chaos Game to obtain the attractor approximation of parameters p and e; then nðp; 2eÞ is the corresponding measure of the set scaled by the factor 12: However, even if our supposition turns out true, one should ask what the assumed limits are in the remaining cases in which the IFS probabilities are not equal. And
again one could speculate in a similar fashion that they are ratios of the corresponding measures induced by the distribution of IFS probabilities (see e.g. [8]). However, the problem appears to deserve a deeper and more formal analysis, so we leave it for the further research. The second example is Barnsley’s Fern (see appendix) that was tested for the original IFS probabilities given in [2] and the ones proposed in [1]. The selected values of e suits visualization at 47 100, 93 200. The results presented in Fig. 4 can be regarded as a confirmation of the thesis by Hepting et al. [1] that their IFS probabilities are really improved in the sense of minimization of the sufficient number of the algorithm’s steps as computed by means of our method. However, it should be noted that the problem whether the number of the Chaos Game steps specified in this way could be considered as an objective measure of the quality of IFS probabilities (with respect to given approximation accuracy) remains open. Nevertheless, we tried to determine an IFS probability distribution that would result in a less number of steps than computed with our algorithm for the distribution proposed in [1]. We found that, in the context of our method, the vector [0.714,
Fig. 5. Rendering of p-balanced measures supported by Barnsley’s Fern. The measures are generated by the following IFS probability vectors: (a) Barnsley’s [0.85, 0.07, 0.07, 0.01], (b) Hepting et al. [0.73, 0.11, 0.13, 0.03], (c) ours improved [0.714, 0.099, 0.119, 0.068]. The same number of the Chaos Game steps was performed to render each fern.
T. Martyn / Computers & Graphics 26 (2002) 753–764
0.099, 0.119, 0.068] seems to be a better candidate for the Fern’s probability distribution than the ones proposed previously (see Fig. 4). The results of visualization of the p-balanced measures induced by the three discussed probability distributions are displayed in Fig. 5. On the basis of the presented images it seems that the new probability vector offers more uniform distribution of points than the previous ones. There are also many other problems open for further research. Some of them follow: *
*
*
Our method estimates probability for the Chaos Game accuracy as the sum of probabilities of the occurrence of each string sAD separately, paying no attention on probabilities of the string occurrence in common. How accurate is such estimation in general case? Can our algorithm be improved by taking into consideration probability of the co-occurrence of strings without a significant increase of the algorithm’s computational complexity? Our algorithm works on the basis of the IFS specification only, regardless of the resulting attractor’s geometry. However, in the case of attractor’s subsets overlapping, there can be more than one string (address) assigned to a given subset [2]. In such cases it is sufficient to consider only one representative for each of the overlapping subsets. Taking advantage of this fact would improve our algorithm both in the accuracy of estimation and time complexity. If the number of the Chaos Game steps that is determined by our method turns out to be a correct measure of the quality of IFS probabilities, a research concerning an efficient optimization algorithm should be taken. At this stage it seems that further work on application of heuristic approaches such as genetic algorithms could prove to be fruitful both in the practical and theoretical sense.
Game will produce, within L of its steps, a sequence that approximates the IFS attractor at a given value of accuracy. As a result, we can determine the sufficient number of the Chaos Game steps to approximate the attractor for given accuracy and probability. However, whether this number can be identified with a measure of the quality of the probability vector assigned to the IFS, remains an open problem. We also outlined a number of other noteworthy problems and questions that have arisen from the analysis of our method and results, which as we believe can be a stimulus for further research.
Acknowledgements I would like to thank Dr. Andrzej Paja) k for fruitful discussions on Markov chains during the development of this method, and to Krzysiek Gracki for his help in production of the figures and patient listening to my tales on the fractals’ and other strange worlds. The shape of this paper could not be reached without suggestions by Prof. Jan Zabrodzki. My special thanks go to the C&G reviewers whose encouraging and valuable comments help me to improve the original version of this paper.
Appendix Here is the IFS code for the Barnsley Fern on the basis of which computation was done. " # " # 0:85 0:04 0:0 w1 ðxÞ ¼ xþ ;
0:04 0:85 1:6 " # " # 0:2 0:26 0:0 xþ ; w2 ðxÞ ¼ 0:23 0:22 1:6 "
8. Summary The IFS representation is a powerful and flexible tool for description of the fractal shapes. The most popular algorithm for rendering IFS attractors is the Chaos Game. As far as the infinite number of iteration steps is considered, the algorithm produces an infinite sequence of points that, with probability one, approximates the IFS attractor at any accuracy. However, in real applications only a finite number of points can be generated, so it is often necessary to know how many iteration steps should be performed to produce an approximation at a given accuracy. The method presented in this paper makes it possible to answer this question. More precisely, given an IFS with probabilities, our algorithm estimates probability that the Chaos
763
# " # 0:28 0:0 xþ ; w3 ðxÞ ¼ 0:24 0:44 " # " # 0:0 0:0 0:0 xþ : w4 ðxÞ ¼ 0:0 0:16 0:0
0:15 0:26
References [1] Hepting D, Prusinkiewicz P, Saupe D. Rendering methods for iterated function systems. Proceedings of Fractals ’90, 1990, IFIP. [2] Barnsley MF. Fractals everywhere, 2nd ed. Boston: Academic Press, 1993. [3] Goodman GS. A probabilist looks at the Chaos Game. Proceedings of Fractals ’90, 1990, IFIP.
764
T. Martyn / Computers & Graphics 26 (2002) 753–764
[4] Knuth DE, Morris JH, Pratt VR. Fast pattern matching in strings. SIAM Journal on Computing 1977;6(2):323–50. [5] Cormen ThH, Leiserson ChE, Rivest RL. Introduction to algorithms. MA: MIT, 1990. [6] Barnsley MF. The desktop fractal design handbook. New York: Academic Press, 1989.
[7] Falconer K. Fractal geometry. mathematical foundations and applications. Chichester: Wiley, 1990. [8] Demko S, Hodges L, Naylor B. Construction of fractal objects with iterated function systems. Computer Graphics (SIGGRAPH ’85) 1985;19(3):271–8.