Non-linear Fractal Interpolating Functions

R. Kobes, H. Letkeman
Department of Physics, University of Winnipeg,
Winnipeg, Manitoba, R3B 2E9 Canada
randy@theory.uwinnipeg.ca

Abstract:

We consider two non-linear generalizations of fractal interpolating functions generated from iterated function systems. The first corresponds to fitting data using a Kth-order polynomial, while the second relates to the freedom of adding certain arbitrary functions. An escape-time algorithm that can be used for such systems to generate fractal images like those associated with Julia or Mandelbrot sets is also described.

Introduction

An Iterated Function System (IFS) can be used to construct a fractal interpolating function for a given set of data [1,2]. The simplest such system defines an IFS

 
Wn$\displaystyle \left(\vphantom{\begin{array}{c}t\\ x\end{array}}\right.$$\displaystyle \begin{array}{c}t\\ x\end{array}$ $\displaystyle \left.\vphantom{\begin{array}{c}t\\ x\end{array}}\right)$ = $\displaystyle \left(\vphantom{ \begin{array}{cc} a_n & 0\\ c_n & 0\end{array} }\right.$$\displaystyle \begin{array}{cc} a_n & 0\\ c_n & 0\end{array}$ $\displaystyle \left.\vphantom{ \begin{array}{cc} a_n & 0\\ c_n & 0\end{array} }\right)$$\displaystyle \left(\vphantom{\begin{array}{c}t\\ x\end{array}}\right.$$\displaystyle \begin{array}{c}t\\ x\end{array}$ $\displaystyle \left.\vphantom{\begin{array}{c}t\\ x\end{array}}\right)$ + $\displaystyle \left(\vphantom{\begin{array}{c}e_n\\ f_n\end{array}}\right.$$\displaystyle \begin{array}{c}e_n\\ f_n\end{array}$ $\displaystyle \left.\vphantom{\begin{array}{c}e_n\\ f_n\end{array}}\right)$. (1)

with coefficients an, cn, en, and fn determined from discrete data points (ti, xi), i = 0, 1,..., N. Such an IFS interpolates the data set in the sense that, under certain assumptions on the coefficients [1], the attractor of the IFS is a graph that passes through the data points. In this particular case, the IFS of Eq. (1) uses a linear (in t) function between adjacent points to interpolate the data.

Various generalizations of fractal interpolating functions have been considered, including those for higher dimensional functions, the use of hidden variables, and extensions to certain non-linear distortions [3,4,5,6,7]. In this note we describe a generalization whereby the transformation incorporates a Kth-order polynomial (in t) interpolation. We also discuss certain classes of non-linear functions that can arise in such interpolating functions, and show how such functions can, with the use of a particular escape-time algorithm, be used to generate fractal images.

   
Linear Interpolating Functions

We first review the construction of a linear fractal interpolating function. Suppose we have data points (ti, xi), i = 0...N, describing a function x(t). Consider the IFS of Eq. (1) and impose, for n = 1, 2,..., N,

 
Wn$\displaystyle \left(\vphantom{\begin{array}{c}t_0\\ x_0\end{array}}\right.$$\displaystyle \begin{array}{c}t_0\\ x_0\end{array}$ $\displaystyle \left.\vphantom{\begin{array}{c}t_0\\ x_0\end{array}}\right)$ = $\displaystyle \left(\vphantom{\begin{array}{c}t_{n-1}\\ x_{n-1}\end{array}}\right.$$\displaystyle \begin{array}{c}t_{n-1}\\ x_{n-1}\end{array}$ $\displaystyle \left.\vphantom{\begin{array}{c}t_{n-1}\\ x_{n-1}\end{array}}\right)$,        Wn$\displaystyle \left(\vphantom{\begin{array}{c}t_N\\ x_N\end{array}}\right.$$\displaystyle \begin{array}{c}t_N\\ x_N\end{array}$ $\displaystyle \left.\vphantom{\begin{array}{c}t_N\\ x_N\end{array}}\right)$ = $\displaystyle \left(\vphantom{\begin{array}{c}t_{n}\\ x_{n}\end{array}}\right.$$\displaystyle \begin{array}{c}t_{n}\\ x_{n}\end{array}$ $\displaystyle \left.\vphantom{\begin{array}{c}t_{n}\\ x_{n}\end{array}}\right)$. (2)

This leads to determination of the coefficients as

 
an = $\displaystyle {\frac{t_n - t_{n-1}}{t_N-t_0}}$,     en = $\displaystyle {\frac{t_{n-1}t_N - t_{n}t_0}{t_N-t_0}}$
cn = $\displaystyle {\frac{x_n - x_{n-1}}{t_N-t_0}}$,     fn = $\displaystyle {\frac{x_{n-1}t_N - x_nt_{0}}{t_N-t_0}}$,
(3)

whereby the transformation of Eq. (1) can be written as
 
Wn(t) $\displaystyle \equiv$ t$\scriptstyle \prime$ = $\displaystyle {\frac{(t-t_0)}{(t_N-t_0)}}$ tn + $\displaystyle {\frac{(t-t_N)}{(t_0-t_N)}}$ tn - 1  
Wn(x) $\displaystyle \equiv$ x$\scriptstyle \prime$ = $\displaystyle {\frac{(t^\prime-t_{n-1})}{(t_n-t_{n-1})}}$ xn + $\displaystyle {\frac{(t^\prime-t_n)}{(t_{n-1}-t_n)}}$ xn - 1 (4)

Thus, Wn(x) $ \equiv$ x$\scriptstyle \prime$ is determined by a linear (in t) interpolating function constructed between the points ( tn - 1, xn - 1) and (tn, xn). Assuming that the corresponding IFS transformation is contractive, so that the distance d (t, x) between any two points in the range of interest satisfies d (Wn(t), Wn(x))$ \le$sn d (t, x),where 0 < sn$ \le$1 is the contractivity factor, graphs of fractal interpolating functions can be made by applying the standard random iteration algorithm to the IFS:
  • initialize (t, x) to a point in the interval of interest
  • for a set number of iterations
    • randomly select a transformation Wn(t, x)
    • plot ( t$\scriptstyle \prime$, x$\scriptstyle \prime$) = Wn(t, x)
    • set (t, x) = ( t$\scriptstyle \prime$, x$\scriptstyle \prime$)
  • end for

A generalization of this type of fractal interpolating function can be found by considering an IFS of the form

 
Wn(t) = ant + en  
Wn(x) = cnt + fn + gn(x) (5)

where gn(x) is, at this stage, an arbitrary function (Ref. [1] has considered the cases gn(x) = dnx and gn(x) = dnx2). Imposing the conditions (2) leads to the same expressions for an and en as in Eq. (3) but cn and fn given by

cn = $\displaystyle {\frac{x_n - x_{n-1}}{t_N-t_0}}$ - $\displaystyle {\frac{g_n(x_N)-g_n(x_0)}{t_N-t_0}}$,        fn = $\displaystyle {\frac{x_{n-1}t_N - x_nt_{0}}{t_N-t_0}}$ - $\displaystyle {\frac{g_n(x_0)t_N-g_n(x_N)t_0}{t_N-t_0}}$. (6)

The resulting transformation then has the same form as Eq. (4) with the following replacements:
xn $\displaystyle \to$ xn + gn(x) - gn(xN)  
xn - 1 $\displaystyle \to$ xn - 1 + gn(x) - gn(x0). (7)

   
Quadratic and Higher Order Interpolating Functions

The interpolating function of the last section used a linear (in t) approximation between adjacent points. In this section we indicate how a quadratic, and then arbitrary Kth-order polynomial, approximation may be constructed. Let us consider a transformation of the form
Wn(t) = ant + en  
Wn(x) = cnt + dnt2 + fn (8)

and impose the conditions, for n = 2, 3,..., N,

 
Wn$\displaystyle \left(\vphantom{\begin{array}{c}t_0\\ x_0\end{array}}\right.$$\displaystyle \begin{array}{c}t_0\\ x_0\end{array}$ $\displaystyle \left.\vphantom{\begin{array}{c}t_0\\ x_0\end{array}}\right)$ = $\displaystyle \left(\vphantom{\begin{array}{c}t_{n-2}\\ x_{n-2}\end{array}}\right.$$\displaystyle \begin{array}{c}t_{n-2}\\ x_{n-2}\end{array}$ $\displaystyle \left.\vphantom{\begin{array}{c}t_{n-2}\\ x_{n-2}\end{array}}\right)$,        Wn$\displaystyle \left(\vphantom{\begin{array}{c}t_m\\ x_m\end{array}}\right.$$\displaystyle \begin{array}{c}t_m\\ x_m\end{array}$ $\displaystyle \left.\vphantom{\begin{array}{c}t_m\\ x_m\end{array}}\right)$ = $\displaystyle \left(\vphantom{\begin{array}{c}t_{n-1}\\ x_{n-1}\end{array}}\right.$$\displaystyle \begin{array}{c}t_{n-1}\\ x_{n-1}\end{array}$ $\displaystyle \left.\vphantom{\begin{array}{c}t_{n-1}\\ x_{n-1}\end{array}}\right)$,        Wn$\displaystyle \left(\vphantom{\begin{array}{c}t_N\\ x_N\end{array}}\right.$$\displaystyle \begin{array}{c}t_N\\ x_N\end{array}$ $\displaystyle \left.\vphantom{\begin{array}{c}t_N\\ x_N\end{array}}\right)$ = $\displaystyle \left(\vphantom{\begin{array}{c}t_{n}\\ x_{n}\end{array}}\right.$$\displaystyle \begin{array}{c}t_{n}\\ x_{n}\end{array}$ $\displaystyle \left.\vphantom{\begin{array}{c}t_{n}\\ x_{n}\end{array}}\right)$. (9)

The point tm is determined as

 
tm = $\displaystyle {\frac{(t_{n-1}-t_{n-2})}{(t_n-t_{n-2})}}$ tN + $\displaystyle {\frac{(t_{n-1}- t_{n})}{(t_{n-2}-t_{n})}}$ t0 (10)

with corresponding point xm. The coefficients of the IFS are determined as
an = $\displaystyle {\frac{t_n-t_{n-2}}{t_N-t_0}}$  
en = $\displaystyle {\frac{t_Nt_{n-2}-t_0t_n}{t_N-t_0}}$  
cn = $\displaystyle {\frac{x_n(t_0^2-t_m^2) +x_{n-1}(t_N^2-t_0^2)+x_{n-2}(t_m^2-t_N^2)}{ (t_N-t_0)(t_N-t_m)(t_m-t_0)}}$  
dn = $\displaystyle {\frac{x_n(t_m-t_0) +x_{n-1}(t_0-t_N)+x_{n-2}(t_N-t_m)}{ (t_N-t_0)(t_N-t_m)(t_m-t_0)}}$  
fn = $\displaystyle {\frac{x_nt_mt_0(t_m-t_0) +x_{n-1}t_Nt_0(t_0-t_N)
+x_{n-2}t_Nt_m(t_N-t_m)}{ (t_N-t_0)(t_N-t_m)(t_m-t_0)}}$ (11)

With this, the transformation can be written as
 
Wn(t) $\displaystyle \equiv$ t$\scriptstyle \prime$ = $\displaystyle {\frac{(t-t_0)}{(t_N-t_0)}}$ tn + $\displaystyle {\frac{(t-t_N)}{(t_0-t_N)}}$ tn - 2  
Wn(x) $\displaystyle \equiv$ x$\scriptstyle \prime$ = $\displaystyle {\frac{(t^\prime-t_{n-1})(t^\prime-t_{n-2})}{(t_n-t_{n-1})(t_n-t_{n-2})}}$ xn + $\displaystyle {\frac{(t^\prime-t_{n})(t^\prime-t_{n-2})}{(t_{n-1}-t_{n})(t_{n-1}-t_{n-2})}}$ xn - 1  
  + $\displaystyle {\frac{(t^\prime-t_{n-1})(t^\prime-t_{n})}{(t_{n-2}-t_{n-1})(t_{n-2}-t_{n})}}$ xn - 2 (12)

which thus uses a quadratic (in t$\scriptstyle \prime$) interpolating function between the points (tn, xn), ( tn - 1, xn - 1), and ( tn - 2, xn - 2).

As in the previous section, including an arbitrary function gn(x) in the IFS transformation via

 
Wn(t) = ant + en  
Wn(x) = cnt + dnt2 + fn + gn(x) (13)

is straightforward. The conditions (9) leads to determination of the point tm of Eq. (10) as before, together with the accompanying point xm. The transformation itself has the same form as Eq. (12) with the following replacements:
xn $\displaystyle \to$ xn + gn(x) - gn(xN)  
xn - 1 $\displaystyle \to$ xn - 1 + gn(x) - gn(xm)  
xn - 2 $\displaystyle \to$ xn - 2 + gn(x) - gn(x0). (14)

From these considerations, the pattern to constructing a Kth-order fractal interpolating function is apparent. Start with a transformation of the form

 
Wn(t) = ant + en  
Wn(x) = Bn(0) + Bn(1)t + Bn(2)t2 +...+ Bn(K)tK (15)

and impose the conditions, for n = K, K + 1,..., N,

 
Wn$\displaystyle \left(\vphantom{\begin{array}{c}t_0\\ x_0\end{array}}\right.$$\displaystyle \begin{array}{c}t_0\\ x_0\end{array}$ $\displaystyle \left.\vphantom{\begin{array}{c}t_0\\ x_0\end{array}}\right)$ = $\displaystyle \left(\vphantom{\begin{array}{c}t_{n-K}\\ x_{n-K}\end{array}}\right.$$\displaystyle \begin{array}{c}t_{n-K}\\ x_{n-K}\end{array}$ $\displaystyle \left.\vphantom{\begin{array}{c}t_{n-K}\\ x_{n-K}\end{array}}\right)$,        Wn$\displaystyle \left(\vphantom{\begin{array}{c}t_{m1}\\ x_{m1}\end{array}}\right.$$\displaystyle \begin{array}{c}t_{m1}\\ x_{m1}\end{array}$ $\displaystyle \left.\vphantom{\begin{array}{c}t_{m1}\\ x_{m1}\end{array}}\right)$ = $\displaystyle \left(\vphantom{\begin{array}{c}t_{n-K+1}\\ x_{n-K+1}\end{array}}\right.$$\displaystyle \begin{array}{c}t_{n-K+1}\\ x_{n-K+1}\end{array}$ $\displaystyle \left.\vphantom{\begin{array}{c}t_{n-K+1}\\ x_{n-K+1}\end{array}}\right)$,    ...,    Wn$\displaystyle \left(\vphantom{\begin{array}{c}t_N\\ x_N\end{array}}\right.$$\displaystyle \begin{array}{c}t_N\\ x_N\end{array}$ $\displaystyle \left.\vphantom{\begin{array}{c}t_N\\ x_N\end{array}}\right)$ = $\displaystyle \left(\vphantom{\begin{array}{c}t_{n}\\ x_{n}\end{array}}\right.$$\displaystyle \begin{array}{c}t_{n}\\ x_{n}\end{array}$ $\displaystyle \left.\vphantom{\begin{array}{c}t_{n}\\ x_{n}\end{array}}\right)$ (16)

which determine the coefficients Bn(i) of Eq. (15), as well as the K - 1 intermediate points tmj, j = 1, 2,..., K - 1 and corresponding xmj points. The resulting transformation will be of the form given by Lagrange's formula for a Kth-order polynomial interpolating function constructed from K + 1 points:
Wn(t) $\displaystyle \equiv$ t$\scriptstyle \prime$ = $\displaystyle {\frac{(t-t_0)}{(t_N-t_0)}}$ tn + $\displaystyle {\frac{(t-t_N)}{(t_0-t_N)}}$ tn - K  
Wn(x) $\displaystyle \equiv$ x$\scriptstyle \prime$ = $\displaystyle {\frac{(t^\prime-t_{n-1})(t^\prime-t_{n-2})\cdots(t^\prime-t_{n-K})}{(t_n-t_{n-1})(t_n-t_{n-2})\cdots(t_n-t_{n-K})}}$ xn +  
  + $\displaystyle {\frac{(t^\prime-t_{n})(t^\prime-t_{n-2})\cdots(t^\prime-t_{n-K})}{(t_{n-1}-t_{n})(t_{n-1}-t_{n-2})\cdots(t_{n-1}-t_{n-K})}}$ xn - 1 +...+  
  + $\displaystyle {\frac{(t^\prime-t_{n})(t^\prime-t_{n-1})\cdots(t^\prime-t_{n-K+1})}{(t_{n-K}-t_{n})(t_{n-K}-t_{n-1})\cdots(t_{n-K}-t_{n-K+1})}}$ xn - K (17)

The inclusion of an arbitrary function gn(x) in the transformation Wn(x) of Eq. (15), as was done for the linear and quadratic transformations of Eqs. (5) and (13) respectively, is straightforward. As might be expected, the use of these higher-order interpolating functions can increase the accuracy of the interpolation significantly, at least for smooth functions. As an informal example, consider the interpolation of the function x4exp(- x)cos(x) - a comparison of the relative error using a linear vs. a quadratic interpolation function appears in Fig. 1, where it can be seen that the a quadratic interpolation offers an order-of-magnitude improvement in the fit over most of the range.


  
Figure: Comparison of interpolations of x4exp(- x)cos(x) using a linear (boxes) and a quadratic (stars) function.
\includegraphics[width=3.375in]{error}

   
Escape-Time Algorithm

In this section we describe an algorithm for generating fractal images like those for Julia or Mandelbrot sets from IFS interpolating functions. As well as allowing one to associate fractal images with discrete data sets, this will also offer another comparison of the detail available from linear vs. higher-order fractal interpolation functions.

Suppose we have an IFS transformation Wn(t, x), generated by some data points (xi, yi), i = 0, 1,..., N, which includes a non-linear function gn(x), as was done for the linear and quadratic transformations of Eqs. (5) and (13) respectively. We now continue the real variable x of this transformation to complex values: x$ \to$z = (zR, zI), in effect defining a complex map:

tn + 1 = antn + en  
zn + 1 = Bn(0) + Bn(1)tn + Bn(2)tn2 +...+ Bn(K)tnK + gn(zn) (18)

for n = 0, 1,..., N. We can then, in analogy with the algorithm used for Julia sets, define the following escape-time algorithm to generate a fractal pattern:
  • for each pixel in a region of interest
    • initialize t
    • initialize z = (zR, zI) to the pixel coordinates
    • for n = 0, 1,..., N
      • calculate ( t$\scriptstyle \prime$, z$\scriptstyle \prime$) = Wn(t, z)
      • break if $ \sqrt{z_R^{\prime\,2} + z_I^{\prime\,2}}$ exceeds a maximum
      • set (t, z) = ( t$\scriptstyle \prime$, z$\scriptstyle \prime$)
    • end for
    • plot the pixel
  • end for
where the pixel is plotted using a coloring algorithm based upon, amongst perhaps other factors, the number of iterations attained when the break condition was met [8].

The arbitrariness of the function gn(x) and the data set (ti, xi) used to fix the IFS leads to a wide variety of possible fractal images. An interesting class of functions gn(x) to consider are those for which, when continued to the complex plane x$ \to$z = (zR, zI), lead to a map having a fixed point zI* = 0 = Imgn(zR, zI*). In such a case one could augment the usual condition of the escape-time algorithm to cease iteration: $ \sqrt{z_R^2 + z_I^2}$ > $ \Lambda$, where $ \Lambda$ is some suitably large value, to also cease iteration when | zI| < $ \epsilon$, where $ \epsilon$ is a suitably small value. The coloring algorithm used to plot a pixel, which depends on the number of iterations attained when this break-out condition was met (if at all), will then lead to structure in the region where the break-out condition on the magnitude of z is not met.

We give two examples of fractal images generated this way for the choice gn(x) = - 0.4x + 0.5x2 + 0.2x3, with the data generated from the logistic map xn + 1 = 3.4xn(1 - xn), with n = 0, 1,...60. The first one, appearing on the left in Fig. 2, corresponds to the generalization Eq. (5) of a linear (in t) fractal interpolating function, while the second image on the right of Fig. 2 corresponds to the generalization Eq. (13) of a quadratic (in t) interpolating function. A coloring algorithm that simply mapped a color to the number of iterations attained when the break-out condition became satisfied was used in both cases.


  
Figure 2: Fractal images from a t-linear (top) and t-quadratic (bottom) interpolating function
\includegraphics[width=2.75in]{ifs1} \includegraphics[width=2.75in]{ifs2}

These figures illustrate, in the interior of the fractal object, the richer structure arising from the quadratic over the linear interpolation function. In this region the break-out condition | zI| < $ \epsilon$ is satisfied, which numerically for $ \epsilon$ $ \sim$ 10-5 is attained after a relatively small number (10-30) of iterations.

   
Conclusions

We considered two non-linear generalizations of fractal interpolating functions constructed from iterated function systems. One - using a Kth-order interpolating polynomial - can potentially improve the accuracy of fractal interpolating functions. The other - use of certain arbitrary functions in the IFS - can, together with an appropriate escape-time algorithm, generate fractal images of the type associated with Julia or Mandelbrot sets.
This work was supported by the Natural Sciences and Engineering Research Council of Canada.

Bibliography

1
M. F. Barnsley, Fractals Everywhere (Academic Press, San Diego, CA, 1993)

2
H. O. Peitgen, H. Jürgens, and D. Saupe, Chaos and Fractals - New Frontiers of Science (Springer, New York, 1992).

3
M. F. Barnsley, J. Elton, D. Hardin, and P. Massopust, Hidden Variable Fractal Interpolation Functions, SIAM J. Math. Anal. 20 (1989), 1218-1242.

4
P. R. Massopust, Fractal Functions, Fractal Surfaces and Wavelets (Academic Press, San Diego, CA, 1994).

5
L. M. Kocic and A. C. Simoncelli, Fractal Interpolation in Creating Prefractal Images, Visual Mathematics 2, No. 2 (2000).

6
H. O. Peitgen and D. Saupe, The Science of Fractal Images (Springer, New York, 1988).

7
E. Gröller, Modeling and Rendering of Nonlinear Iterated Function Systems, Institute of Computer Graphics, Technical University of Vienna Report TR-186-2-94-12 (1994).

8
J. Barrallo and D. M. Jones, Coloring Algorithms for Dynamical Systems in the Complex Plane, Visual Mathematics 1, No. 4 (1999).



VisMath HOME