On Self‐Affine and Self‐Similar Graphs of Fractal Interpolation Functions Generated from Iterated Function Systems

Abstract

This chapter provides a brief and coarse discussion on the theory of fractal interpolation functions and their recent developments including some of the research made by the authors. It focuses on fractal interpolation as well as on recurrent fractal interpolation in one and two dimensions. The resulting self‐affine or self‐similar graphs, which usually have non‐integral dimension, were generated through a family of (discrete) dynamic systems, the iterated function system, by using affine transformations. Specifically, the fractal interpolation surfaces presented here were constructed over triangular as well as over polygonal lattices with triangular subdomains. A further purpose of this chapter is the exploration of the existent breakthroughs and their application to a flexible and integrated software that constructs and visualises the above‐mentioned models. We intent to supply both a panoramic view of interpolating functions and a useful source of links to assist a novice as well as an expert in fractals. The ideas or findings contained in this paper are not claimed to be exhaustive, but are intended to be read before, or in parallel with, technical papers available in the literature on this subject.

1. Introduction

In the mathematical field of numerical analysis, interpolation is a method of constructing new data points within the range of a discrete set of known data points. Interpolation by fractal (graph of) functions, as defined in Refs. [1, 2], is based on the theory of iterated function systems and can be seen as an alternative to traditional interpolation techniques, aiming mainly at data which present detail at different scales, some degree of self‐similarity or self‐affinity. A fractal interpolation function can be seen as a continuous function whose graph is the attractor of an appropriately chosen iterated function system. This attractor is called a fractal interpolation surface, since the graph of a continuous function of two variables defined over a connected open subset of $\mathbb{R}$² is a topological surface, if the so formed graph belongs to the three‐dimensional space and has Hausdorff‐Besicovitch dimension between 2 and 3. It is called a fractal interpolation volume whenever the graph has dimension greater than three. The key difficulty in constructing fractal interpolation surfaces (or volumes) involves ensuring continuity. Another important element necessary in modelling complicated surfaces of this type is the existence of the contractivity, or vertical scaling, factors.

1.1. Historical background

Massopust introduced in Ref. [3] fractal surfaces constructed as attractors of iterated function systems. He considered the case of a triangular domain with coplanar boundary data. Later on, Geronimo and Hardin in Ref. [4] presented a slightly more general construction of such fractal surfaces. They examined the case when the domain is a polygonal region with arbitrary interpolation points but with identical contractivity factors. In Ref. [5], $\mathbb{R}$^m‐valued multivariate fractal functions were investigated. The latter two constructions use the recurrent iterated‐function‐system formalism. The construction of Wittenbrink [6] either produces discontinuous surfaces (and volumes) or reduces to the case, where the contractivity factors must be constant. Zhao [7] allows the contractivity factors to become a continuous ‘contraction function’ and uses consistent triangulation in order to guarantee continuity. All of the previous constructions are based on triangular subdomains. As it is always possible to construct fractal surfaces as tensor products of univariate continuous fractal functions, Massopust in Ref. [8], Section 9.4 suggests a construction by taking the tensor product of two univariate fractal interpolation functions. The derived function is uniquely determined by its evaluation along a pair of adjacent sides of the rectangular domain.

Two piecewise self‐affine models for representing discrete image data on rectangular lattices by using fractal surfaces are proposed in Ref. [9]. In Ref. [10], the piecewise self‐affine IFS model is extended from $\mathbb{R}$³ to $\mathbb{R}$ⁿ (n is an integer greater than 3), which is called the multi‐dimensional piecewise self‐affine fractal interpolation model. The same methodology is repeated in Ref. [11]. According to Ref. [12], the self‐affine iterated function system model is extended from $\mathbb{R}$³ to $\mathbb{R}$ⁿ (n is an integer greater than 3), which is called the multi‐dimensional self‐affine fractal interpolation model. Vijender and Chand in Ref. [13] proposed a class of affine fractal interpolation surfaces that stitch a given set of surface data arranged on a rectangular grid. The created fractal interpolation surfaces are a blend of the affine fractal interpolation functions constructed along the grid lines of a given interpolation domain.

The study of fractals is a field in science that unifies mathematics, theoretical physics, art and computer science. Therefore, it is not difficult to find applications of fractal interpolation functions in almost every scientific field wherein information available in finite number of grid points has to be modelled with a continuous function. Applications of this theory include geometric design, data visualisation, reverse engineering, physics, geology, image encoding and compression (see Ref. [14]), signal processing and wavelet theory. Fractal interpolation also provides a good representation of economic time series such as the stock market fluctuation and weather data. The reason for this variety of applications lies in the underlying complicated mathematical structure of fractal (graph of) functions, produced with simple recursive construction. It has been noted that, for certain problems, they provide better approximants than their classical non‐recursive counterparts.

1.2. Clarifications and interpretations

Although, in most of the cases, the words function, map and mapping can be used as nouns interchangeably in several parts of mathematics, they generally do not have the same meaning. Originally, the word ‘map’ comes from the medieval Latin Mappa mundi, wherein mappa meant napkin or cloth and mundi the world. ‘Mapping’ as a noun is the process of making maps. In mathematics, we think of a map as a way of sending elements from one set to elements of another set. A mapping shows how the elements are paired. A function is a way or rule of associating elements from one set, the domain, to elements of another set, the codomain. So, a domain is where a function maps from, and a range, as a subset of the codomain, is where it maps to. The definition of a function requires that each element in the domain corresponds to one and only one element of the codomain. Moreover, a function is commonly used as a special type of mapping, namely it is used as a mapping from a set into the set of numbers, that is, into $\mathbb{R}$, or $\mathbb{C}$ or into a field $\mathbb{K}$. So, a map is slightly more general, as it allows a many‐to‐one situation, arbitrary sets, etc. For example, we refer to ‘The open mapping theorem’ and not to the ‘The open function theorem’, and we refer to ‘The contraction mapping theorem’ and so on. A mapping of a set into itself is called a transformation.

An affine map (from the Latin, affinis, “connected with”) between two (vector) spaces consists of a linear map followed by a translation. Similarly, one can define an affine transformation. In a geometric setting, these are precisely the functions that map straight lines to straight lines. Be aware that the term linear function refers to two distinct but related notions. It may be a linear map or a polynomial of degree one or less, including the zero polynomial, because its graph, when there is only one independent variable, is a non‐vertical straight line.

Another common error involves the incorrect use of the notions ‘self‐affine function’ and ‘self‐similar function’. There is no ‘self‐affine’ or ‘self‐similar’ function but a function with a self‐affine or self‐similar graph. Therefore, an object or a set, but not a function, can be self‐affine or self‐similar. Each part of a self‐affine object is an image of the whole object (either strictly or in a statistical sense) scaled differently in different directions. Self‐affine sets form an important class of sets, which include self‐similar sets as a particular case. A self‐similar object is exactly or approximately similar to a part of itself (i.e. the whole has the same shape as one or more of the parts).

The aim of our paper is to review the usage of fractal interpolation functions in order to construct self‐affine graphs generated by iterated function systems. Furthermore, we compare and contrast several constructions presented in Refs. [3–8] by pointing out some of their ambiguities, limitations and restrictions. Particularly, in Section 2, we briefly review the theory of iterated function systems. In Section 3, we revisit the 1D fractal interpolation theory and state the prerequisites of the constructions. Necessary conditions for the attractor of an iterated function system to be the graph of a continuous function interpolating a given set of data are also given. In Section 4, we revisit the two‐dimensional fractal interpolation theory. A comparison to existing methods and some examples of fractal interpolation surfaces constructed by them are also presented. The corresponding algorithms for constructing these surfaces are developed and illustrated through several graphic examples. Finally, Section 5 summarises our conclusions and points out areas of future work.

2. Fractal image generation

In mathematics, an iterated map is a map composed with itself, possibly ad infinitum, in a process called iteration. Iteration means the act of repeating a process with the aim of approaching a desired goal, target or result. More formally, let X be a set and f: X→X be a map. Define f^k as the k‐th iterate of f, where k is a non‐negative integer, by f⁰ = id_X and f^k+1 = f ◦ f^k, where id_X is the identity map on X and f ◦ g denotes map composition.

A contraction mapping, or contraction, on a metric space (X, ρ) is a function f from Χ to itself, that is, a transformation, with the property that there is a non‐negative real number s < 1 such that for all x and y in X, ρ(f(x), f(y)) ≤ s⋅ρ(x, y), where ρ is a distance function between elements of X. The smallest such value of s is called the Lipschitz constant or contractivity factor of f. If the above condition is satisfied for s ≤ 1, then the mapping is said to be non‐expansive. A contraction mapping has at most one fixed point, that is, a point x* in X such that f(x*) = x*. Moreover, the Banach fixed point theorem, also known as the contraction mapping theorem or contraction mapping principle, states that every contraction mapping on a non‐empty, complete metric space has a unique fixed point, and that for any x in X the iterated function sequence x, f(x), f(f(x)), f(f(f(x))), … converges to this fixed point. Furthermore, this fixed point can be found as follows: Start with an arbitrary element x₀ in X and define an iterative sequence by x_n = f(x_n−1) for n = 1, 2, 3, …. This sequence converges and its limit is x*.

2.1. Iterated function systems

An iterated function system, or IFS for short, is defined as a collection of a complete metric space (X, ρ), for example, ($\mathbb{R}$ⁿ, ‖⋅‖) or a subset, together with a finite set of continuous transformations {w_i: X→X, i = 1, 2, …, M}. It is often convenient to write an IFS formally as {X; w₁, w₂, …, w_M} or, somewhat more briefly, as {X; w_1−M}. If w_i are contractions with respective contractivity factors s_i, i = 1, 2, …, M, the IFS is termed hyperbolic.

Hutchinson in Ref. [15] showed that, for the metric space $\mathbb{R}$ⁿ, such a (hyperbolic) system of functions has a unique compact (closed and bounded) fixed set S. One way for constructing a fixed set is to start with an initial point or set S₀ and iterate the actions of the w_i, taking S_n+1 to be the union of the image of S_n under the w_i; then taking S to be the closure of the union of the S_n. Symbolically, the unique fixed (non‐empty compact) set S ⊂ $\mathbb{R}$ has the property

$S={\displaystyle \bigcup _{i=1}^M{w}_i(S)}.$E1

The set S is thus the fixed set of the Hutchinson operator

$H(A)={\displaystyle \bigcup _{i=1}^M{w}_i(A)},$E101

where A is any subset of $\mathbb{R}$ⁿ. The operator H itself is a contraction with contractivity factor s = max{s₁, s₂, …, s_M} (Ref. [2], Theorem 7.1, p. 81 or Ref. [16]). The existence and uniqueness of S, which is called the attractor, or invariant set, of the IFS, are a consequence of the contraction mapping principle as is the fact that ${\displaystyle \underset{k\to \infty }{\lim }}{H}^k(A)=S\equiv {\mathcal{A}}_{\infty }$ for all A in $\mathcal{H}$($\mathbb{R}$ⁿ), where $\mathcal{H}$(X) is the metric space of all non‐empty, compact subsets of X with respect to some metric, for example, the Hausdorff metric. The operator H is also called the collage map to alert us to the fact that H(A) is formed as a union or ‘collage’ of sets. If X is a Euclidean space and the w_i are similitudes, that is, similarity transformations, then the attractor is called a self‐similar set. These sets are usually fractals.

A fractal derived by an IFS is made up of the union of several copies of itself, each copy being transformed by a function (hence ‘function system’). A canonical example is the Sierpiński gasket; see Figure 1. The functions are normally contractions which means they bring points closer together and make shapes smaller. Hence, such a shape is made up of several possibly overlapping smaller copies of itself, each of which is also made up of copies of itself, ad infinitum. This is the source of its self‐similar nature. Note that this infinite process is not dependent upon the starting shape being a triangle—it is just clearer that way. The first few steps starting, for example, from a square also tend towards a Sierpiński gasket; see Figure 2.

Figure 1.

The evolution of the Sierpiński gasket.

Figure 2.

Iterating from a square.

Sometimes each function w_i is required to be a linear, or more generally an affine transformation, and hence represented by a matrix. Formally, a transformation w is affine, if it may be represented by a matrix A and translation t as w(x) = Ax + t, or, if X = $\mathbb{R}$²,

$w{\displaystyle \left(\begin{array}{c}x\\ y\end{array}\right)}={\displaystyle \left(\begin{array}{cc}a& b\\ c& s\end{array}\right)}{\displaystyle \left(\begin{array}{c}x\\ y\end{array}\right)}+{\displaystyle \left(\begin{array}{c}d\\ e\end{array}\right)}.$E2

The code of w is the 6‐tuple (a, b, c, s, d, e) and the code of an IFS is a table whose rows are the codes of w₁, w₂, …, w_M. For the three‐dimensional case, this becomes

$w{\displaystyle \left(\begin{array}{c}x\\ y\\ z\end{array}\right)}={\displaystyle \left(\begin{array}{ccc}a& b& c\\ d& e& g\\ h& k& s\end{array}\right)}{\displaystyle \left(\begin{array}{c}x\\ y\\ z\end{array}\right)}+{\displaystyle \left(\begin{array}{c}l\\ m\\ r\end{array}\right)}.$E3

However, IFS’s may also be built from non‐linear functions, including projections and Möbius transformations.

2.2. Recurrent iterated function systems

An IFS with probabilities, written formally as {X; w₁, w₂, …,w_M; p₁, p₂, …, p_M} or, somewhat more briefly, as {X; w_1−M; p_1−M}, gives to each transformation in H a probability or weight. If the weights of transformations differ, so do the measures on different parts of the attractor. A non‐self‐similar attractor, however, is more easily represented with a recurrent iterated function system, or RIFS for short. Each transformation has, instead of a single weight for the next iteration, a vector of weights for each transformation, {X; w_1−M; p_i,j ∈ [0, 1]; i, j = 1, 2, …, M}, so that the matrix of weights is a recurrent Markov operator for the Hutchinson operator’s transformation (see Ref. [17]). Therefore, the attractor of a RIFS need not exhibits the self‐similarity or self‐tiling properties characteristic of the attractor of an IFS, such as Eq. (1).

The most common algorithm to compute fractals derived by IFS is called the chaos game or random iteration algorithm. It consists of picking a random point in the plane, then iteratively applying one of the functions chosen at random from the function system and drawing the point. An alternative algorithm, the deterministic iteration algorithm, or DIA for short, is to generate each possible sequence of functions up to a given maximum length and then to plot the results of applying each of these sequences of functions to an initial point or shape.

3. Fractal interpolation functions in $\mathbb{R}$

FIFs are suitable for data sets with points linearly ordered with respect to their abscissa. This is often sufficient, for example, when interpolating time series data. In practice, however, there are many cases where the data are suitable for fractal interpolation but define a curve rather than a function, for example, when modelling coastlines or plants. There exist methods for constructing fractal interpolation curves based on the theory of FIFs. These methods use various approaches, such as generalisations to higher dimensions, use of index coordinates or application of reversible transformations.

Based on a theorem of Hutchinson ([15], p. 731) and using IFS theory [16], Barnsley introduced a class of functions (see Ref. [1]) which he called fractal interpolation functions. He basically worked with affine fractal interpolation functions, in the sense that they are obtained using affine transformations. Their main difference from elementary functions is their fractal character. Since their graphs usually have non‐integral dimension, they can be used to approximate image components such as the profiles of mountain ranges, the tops of clouds and horizons over forests, to name but a few. For a short survey on fractal interpolation functions see Ref. [19].

Let the continuous function f be defined on a real closed interval I = [x₀, x_M] and with codomain the metric space ($\mathbb{R}$, |⋅|), where x₀ < x₁ < … < x_M. It is not assumed that these points are equidistant. The function f is called an interpolation function corresponding to the generalised set of data {(x_m, y_m) ∈ K : m = 0, 1, …, M}, if f(x_m) = y_m for all m = 0, 1, …, M and K = I × $\mathbb{R}$. The points (x_m, y_m) ∈ $\mathbb{R}$² are called the interpolation points. We say that the function f interpolates the data and that (the graph of) f passes through the interpolation points. The graph of f is the set of points G(f) = {(x, f(x) : x ∈ X}.

3.1. Affine fractal interpolation

Let us represent our, real valued, set of data points as {(u_n, v_n) : n = 0, 1, …, N; u_n < u_n+1} and the interpolation points as {(x_m, y_m) : m = 0, 1, …, M; M ≤ N}, where u_n is the sampled index and v_n the value of the given point in u_n. Let {$\mathbb{R}$²; w_1−M} be an IFS with affine transformations of the special form (see Eq. (2))

$w_i{\displaystyle \left(\begin{array}{c}x\\ y\end{array}\right)}={\displaystyle \left(\begin{array}{cc}{a}_i& 0\\ c_i& s_i\end{array}\right)}{\displaystyle \left(\begin{array}{c}x\\ y\end{array}\right)}+{\displaystyle \left(\begin{array}{c}{d}_i\\ e_i\end{array}\right)}$E102

constrained to satisfy

$w_i{\displaystyle \left(\begin{array}{c}{x}_0\\ y_0\end{array}\right)}={\displaystyle \left(\begin{array}{c}{x}_{i-1}\\ y_{i-1}\end{array}\right)}\text{and}{w}_i{\displaystyle \left(\begin{array}{c}{x}_M\\ y_M\end{array}\right)}={\displaystyle \left(\begin{array}{c}{x}_i\\ y_i\end{array}\right)}$E4

for every i = 1, 2, …, M. Solving the above equations results in

$a_i=\frac{x_i-x_{i-1}}{x_M-x_0},d_i=\frac{x_M{x}_{i-1}-x_0{x}_i}{x_M-x_0}$E103

$c_i=\frac{y_i-y_{i-1}}{x_M-x_0}-s_i\frac{y_M-y_0}{x_M-x_0},e_i=\frac{x_M{y}_{i-1}-x_0{y}_i}{x_M-x_0}-s_i\frac{x_M{y}_0-x_0{y}_M}{x_M-x_0}$E104

i.e., the coefficients a_i, c_i, d_i, e_i are completely determined by the interpolation points, while the s_i are free parameters satisfying |s_i| < 1 in order to guarantee that the IFS is hyperbolic with respect to an appropriate metric for every i = 1, 2, …, M. The transformations w_i are shear transformations: line segments parallel to the y‐axis are mapped to line segments parallel to the y‐axis contracted by the factor |s_i|. For this reason, the s_i are called vertical scaling (or contractivity) factors.

The IFS {$\mathbb{R}$²; w_1−M} has a unique attractor that is the graph of some continuous function which interpolates the data points. This function is called a fractal interpolation function, or FIF for short, because its graph usually has non‐integral dimension. A section is defined as the function values between interpolation points. It is a function with a self‐affine graph since each affine transformation w_i maps the entire (graph of the) function to its section. The above function is known as affine FIF, or AFIF for short. For example, let {(0, 0), (0.4, 0.5), (0.7, 0.2), (1, 0)} be a given set of data points. Figure 3 shows the graph of an AFIF with s₁ = 0.5, s₂ = −0.2 and s₃ = 0.4. The closeness of fit of a FIF is mainly influenced by the determination of its vertical scaling factors. No direct way to find the optimum values of these factors exists, and various approaches have been proposed in the literature.

Figure 3.

The construction of an affine FIF starting from the unit square and using the DIA.

3.2. Piecewise affine fractal interpolation

The piecewise self‐affine fractal model is a generalisation of the affine fractal interpolation model and has its mathematical roots embedded in RIFS theory. A pair of data points $\{({\tilde{x}}_{i,j},{\tilde{y}}_{i,j}):i=1,2,\mathrm{\dots },M;j=1,2\},$ which are called addresses or address points, is now associated with each interpolation interval. Each pair of addresses defines the domain or address interval. The constraints Eq. (4) become

$w_i{\displaystyle \left(\begin{array}{c}{\tilde{x}}_{i,1}\\ {\tilde{y}}_{i,1}\end{array}\right)}={\displaystyle \left(\begin{array}{c}{x}_{i-1}\\ y_{i-1}\end{array}\right)}\text{and}{w}_i{\displaystyle \left(\begin{array}{c}{\tilde{x}}_{i,2}\\ {\tilde{y}}_{i,2}\end{array}\right)}={\displaystyle \left(\begin{array}{c}{x}_i\\ y_i\end{array}\right)}$E105

subjected to ${\tilde{x}}_{i,2}-{\tilde{x}}_{i,1}>x_i-x_{i-1}$ for every i = 1, 2, …, M. Solving the above equations results in

$a_i=\frac{x_i-x_{i-1}}{{\tilde{x}}_{i,2}-{\tilde{x}}_{i,1}},d_i=\frac{{\tilde{x}}_{i,2}{x}_{i-1}-{\tilde{x}}_{i,1}{x}_i}{{\tilde{x}}_{i,2}-{\tilde{x}}_{i,1}},$E106

$c_i=\frac{y_i-y_{i-1}}{{\tilde{x}}_{i,2}-{\tilde{x}}_{i,1}}-s_i\frac{{\tilde{y}}_{i,2}-{\tilde{y}}_{i,1}}{{\tilde{x}}_{i,2}-{\tilde{x}}_{i,1}},e_i=\frac{{\tilde{x}}_{i,2}{y}_{i-1}-{\tilde{x}}_{i,1}{y}_i}{{\tilde{x}}_{i,2}-{\tilde{x}}_{i,1}}-s_i\frac{{\tilde{x}}_{i,2}{\tilde{y}}_{i,1}-{\tilde{x}}_{i,1}{\tilde{y}}_{i,2}}{{\tilde{x}}_{i,2}-{\tilde{x}}_{i,1}}$E107

for every i = 1, 2, …, M. The function constructed as the attractor of the above‐mentioned IFS is called recurrent fractal interpolation function, or RFIF shortly, corresponding to the interpolation points. A RFIF is a piecewise self‐affine function since each affine transformation w_i maps the part of the (graph of the) function defined by the corresponding address interval to each section.

4. Fractal interpolation functions in $\mathbb{R}$²

Let the discrete data {(x_i, y_j, z_i,j = z(x_i, y_j)) ∈  $\mathbb{R}$³ : i = 0, 1, …, N; j = 0, 1, …, M} be known. Each affine mapping that comprises the hyperbolic IFS {$\mathbb{R}$³; w_{1−N, 1−M}} is given by the following special form of Eq. (3)

$w_{n,m}{\displaystyle \left(\begin{array}{c}x\\ y\\ z\end{array}\right)}={\displaystyle \left(\begin{array}{ccc}{a}_{n,m}& b_{n,m}& 0\\ c_{n,m}& d_{n,m}& 0\\ e_{n,m}& g_{n,m}& s_{n,m}\end{array}\right)}{\displaystyle \left(\begin{array}{c}x\\ y\\ z\end{array}\right)}+{\displaystyle \left(\begin{array}{c}{h}_{n,m}\\ k_{n,m}\\ l_{n,m}\end{array}\right)},$E108

with |s_n,m| < 1 for every n = 1, 2, …, N and m = 1, 2,…, M. The condition

$\Vert \begin{array}{cc}{a}_{n,m}& b_{n,m}\\ c_{n,m}& d_{n,m}\end{array}\Vert <1$E109

ensures that

$u_{n,m}{\displaystyle \left(\begin{array}{c}x\\ y\end{array}\right)}={\displaystyle \left(\begin{array}{cc}{a}_{n,m}& b_{n,m}\\ c_{n,m}& d_{n,m}\end{array}\right)}{\displaystyle \left(\begin{array}{c}x\\ y\end{array}\right)}+{\displaystyle \left(\begin{array}{c}{h}_{n,m}\\ k_{n,m}\end{array}\right)}$E110

is a similitude and the transformed surface does not vanish or flip over.

A formal definition for the fractal interpolation surfaces, as presented in Ref. [4], with some generalisations is given below. Let D be a convex polygon in $\mathbb{R}$² with ℓ vertices and let S = {q₀, q_{1, …,} q_m−1} be m, with m > ℓ, distinct points in D such that q₀, q_1,…, q_ℓ−1 are the vertices of D. Given real numbers z₀, z₁, …, z_m−1, we wish to construct a function f such that f(q_j) = z_j, j = 0, 1, …, m−1 and whose graph is self‐affine or self‐similar. Let us denote by C(D) the linear space of all real‐valued continuous functions defined on D, that is,

$C{\displaystyle \left(D\right)}=\{f:D\to \mathbb{R}│f\hspace{0.17em}\text{c}\text{o}\text{n}\text{t}\text{i}\text{n}\text{u}\text{o}\text{u}\text{s}\}.$E9300

A mapping Φ on C(D) which corresponds to piecing the G(Φ(f)) together using copies of parts of G(f) will be defined.

4.1. Affine fractal interpolation

The basic idea is to decompose D into N non‐degenerate subregions Δ₁, Δ₂, …, Δ_N with vertices chosen from S and define affine maps L_i: D → Δ_i and F_i: D × $\mathbb{R}$ → $\mathbb{R}$, i = 1, 2, …, N such that Φ defined by

$\mathrm{\Phi }(f)(x)=F_i{\displaystyle \left(L_i^{-1}(x),f{\displaystyle \left(L_i^{-1}(x)\right)}\right)}$E5

for x ∈ Δ_i maps an appropriate subset of C(D) onto itself. The main work is in showing that Φ is well defined and contractive on some subset of C(D). If L_i is invertible, G(f) is mapped onto $G{\displaystyle \left(\Phi (f){|}_{{\Delta }_i}\right)}$ by (L_i, F_i). Henceforth, we assume that the set ${{\displaystyle \left\{{\Delta }_i\right\}}}_{i=1}^N$ consists of non‐degenerate convex polygons with extreme points whose interior are non‐intersecting, $L_i^{-1}({\Delta }_i)=D$ and the set of vertices of ${{\displaystyle \left\{{\Delta }_i\right\}}}_{i=1}^N$ equals S. Let k: {1, 2, …, N} × {0, 1, …, ℓ−1} → {0, 1, …, m−1} be such that${\left\{q_{k(i,j)}\right\}}_{j=0}^{\ell -1}$ gives the vertices of ${{\displaystyle \left\{{\Delta }_i\right\}}}_{i=1}^N.$

Let i ∈ {1, 2, …, N}. Since D and Δ_i are non‐degenerate, there exists an invertible map satisfying

$L_i(q_j)=q_k{}_{(i,j)},j=0,1,\dots ,\ell -1.$E6

Also, given any necessary free parameters,

$F_i(q_j,z_j)=z_k{}_{(i,j)},j=0,1,\dots ,\ell -1.$E7

With these definitions for L_i and F_i, if f ∈ C(D) and f(q_j) = z_j, j = 0, 1, …, ℓ−1, then Φ$|$_Δi ∈ C(Δ_i) and Φ(f)(q_{k(i, j)}) = z_{k(i, j)}, j = 0, 1, …, ℓ−1. If Δ_i and Δ_i' are adjacent polygons with common edge $\overline{q_j{q}_{j\text{'}}}$, it remains to be determined if Φ is well defined along $\overline{q_j{q}_{j\text{'}}}$, that is, whether Φ(f) satisfies for all x ∈ $\overline{q_j{q}_{j\text{'}}}$ the ‘join‐up’ condition

$F_i(L_i^{-1}(x),f(L_i^{-1}(x)))=F_{i\text{'}}(L_{i\text{'}}^{-1}(x),f(L_{i\text{'}}^{-1}(x))).$E112

When there is no proof that our construction always satisfies it, the surface may be not continuous and a geometrical and visual artefact, known as ‘discontinuity’, appears. When a surface suffers from discontinuities, a correct visualisation should render aligned horizontal holes over the surface. A surface with discontinuities should not be considered as a fractal interpolation surface; the function f is ambiguous on the common edge points and Φ(f) is not well defined.

Let the number of extreme points of the convex region D be 3. A triangular domain is formed, and the set ${{\displaystyle \left\{{\Delta }_i\right\}}}_{i=1}^N$ contains non‐degenerate triangles whose interior is non‐intersecting. In what follows, we shall deal with a special class of maps, namely affine ones. For all i = 1, 2, …, N the invertible map L_i: $\mathbb{R}$²→ $\mathbb{R}$² is defined as follows

$L_i(x,y)=(a_{i}x+b_i,c_{i}y+d_i)$E113

and let the map F_i: $\mathbb{R}$³→ $\mathbb{R}$ be defined by

$F_i(x,y,z)=e_{i}x+f_{i}y+s_{i}z+g_i.$E8

Then, the corresponding IFS is of the form {$\mathbb{R}$³; w_1−N}, where w_i(x, y, z) = (L_i(x, y), F_i(x, y, z)), i = 1, 2, …, N and can be written in a matrix form

$w_i{\displaystyle \left(\begin{array}{l}x\\ y\\ z\end{array}\right)}={\displaystyle \left(\begin{array}{ccc}{a}_i& 0& 0\\ 0& c_i& 0\\ e_i& f_i& s_i\end{array}\right)}{\displaystyle \left(\begin{array}{l}x\\ y\\ z\end{array}\right)}+{\displaystyle \left(\begin{array}{l}{b}_i\\ d_i\\ g_i\end{array}\right)}.$E114

The real numbers a_i, b_i, c_i and d_i for i = 1, 2, …, N are chosen to ensure that Condition (6) holds, that is, L_i(q₀) = q_{k(i, 0)}, L_i(q₁) = q_{k(i, 1)} and L_i(q₂) = q_{k(i, 2)}. Thus, for i = 1, 2, …, N,

$\begin{array}{l}{a}_i=\frac{x_{k(i,j\text{'})}-x_{k(i,j)}}{x_{j^{\prime }}-x_j},b_i=x_{k(i,j)}-x_j\frac{x_{k(i,j\text{'})}-x_{k(i,j)}}{x_{j^{\prime }}-x_j},j\ne j\text{'}\in \{0,1,2\},x_j\ne x_{j^{\prime }},\\ c_i=\frac{y_{k(i,j\text{'})}-y_{k(i,j)}}{y_{j\text{'}}-y_j},d_i=y_{k(i,j)}-y_j\frac{y_{k(i,j\text{'})}-y_{k(i,j)}}{y_{j^{\prime }}-y_j},j\ne j\text{'}\in \{0,1,2\},y_j\ne y_{j^{\prime }}.\end{array}$E115

After selecting the scaling factors s_i with |s_i|< 1, the values e_i, f_i and g_i are chosen to ensure that Condition (7) holds, that is, F_i(q₀, z₀) = z_{k(i, 0)}, F_i(q₁, z₁) = z_{k(i, 1)} and F_i(q₂, z₂) = z_{k(i, 2)}. That is, s_i ∈ (−1, 1) is chosen and then

$e_i=\frac{z_{k(i,0)}(y_1-y_2)+z_{k(i,2)}(y_0-y_1)+z_{k(i,1)}(y_2-y_0)}{x_0(y_1-y_2)+x_2(y_0-y_1)+x_1(y_2-y_0)}+s_i\frac{z_1(y_0-y_2)+z_0(y_2-y_1)+z_2(y_1-y_0)}{x_0(y_1-y_2)+x_2(y_0-y_1)+x_1(y_2-y_0)},$E116

$f_i=\frac{z_{k(i,1)}(x_0-x_2)+z_{k(i,0)}(x_2-x_1)+z_{k(i,2)}(x_1-x_0)}{y_1(x_0-x_2)+y_0(x_2-x_1)+y_2(x_1-x_0)}+s_i\frac{z_0(x_1-x_2)+z_2(x_0-x_1)+z_1(x_2-x_0)}{y_1(x_0-x_2)+y_0(x_2-x_1)+y_2(x_1-x_0)},$E117

$\begin{array}{l}{g}_i=\frac{x_0(y_1{z}_{k(i,2)}-y_2{z}_{k(i,1)})+x_1(y_2{z}_{k(i,0)}-y_0{z}_{k(i,2)})+x_2(y_0{z}_{k(i,1)}-y_1{z}_{k(i,0)})}{x_0(y_1-y_2)+x_1(y_2-y_0)+x_2(y_0-y_1)}+\\ s_i\frac{x_0(y_2{z}_1-y_1{z}_2)+x_1(y_0{z}_2-y_2{z}_0)+x_2(y_1{z}_0-y_0{z}_1)}{x_0(y_1-y_2)+x_1(y_2-y_0)+x_2(y_0-y_1)},\end{array}$E118

for all i = 1, 2, …, N. We construct the IFS of the form {S; w_1−N}, with the intention to produce G(f) that is continuous and passes through the points (q_j, z_j), q_j ∈ S, j = 0, 1, …, m−1 as the unique attractor $\mathcal{A}$ of the IFS. Since we consider f ∈ C(D), we must check if the ‘join‐up’ condition is satisfied for all possible starting points of the IFS. Figure 4 illustrates the surface graph (Level = 0) drawn with the original set of data {(0, 0, 0), (0.5, 0, 0.2), (1, 0, 0), (0, 0.5, 0.3), (0.5, 0.5, 0.5), (0, 1, 0)} where s = 0.6. See also Example 58 of Ref. [8].

Figure 4.

The generation of an affine fractal surface and its view from above.

4.1.1. Coplanar boundary data, arbitrary contractivity factors

If P is a non‐vertical plane in $\mathbb{R}$³, let Ĉ(D) denote the collection of continuous functions f: D → $\mathbb{R}$ such that (x, f(x)) ∈ P for all x ∈ ∂D.

Theorem 4.1.1. Suppose the points {(q_j, z_j) :q_j ∈ ∂D} are contained in a plane P ⊂ $\mathbb{R}$³. Let Φ be defined by (5), where L_i and F_i, i = 1, 2, …, N are defined by Eqs. (6)–(8). Then, Φ: Ĉ(D) → Ĉ(D) is well defined and contractive in the sup‐norm ‖⋅‖_∞. Furthermore, for every j = 0, 1, …, m−1 and f ∈ Ĉ(D), Φ(f)(q_j) = z_j.

We call the unique attractor of the afore‐mentioned IFS a self‐affine fractal interpolation surface, or SAFIS for short, with coplanar boundary data and arbitrary contractivity factors. Figure 5 illustrates Example 1 constructed in Ref. [4].

Figure 5.

Three views of a SAFIS with coplanar boundary data and arbitrary contractivity factors.

4.1.2. Arbitrary boundary data, identical contractivity factors

The chromatic number of a graph is the smallest number of colours needed to colour its vertices so that no two adjacent vertices share the same colour. Let G be the graph with S as its vertices and its edges correspond to the decomposition of D to ${{\displaystyle \left\{{\Delta }_i\right\}}}_{i=1}^N.$ We assign a colour to each vertex of the graph through the labelling function l = l(j) ∈ {0, 1, 2}, j = 0, 1, …, m−1; see Figure 6.

Figure 6.

Domain for fractal interpolating surfaces over triangular lattice and possible subdomains.

For i = 1, 2, …, N and j = 0, 1, 2 let k(i, j) be determined by the condition k(i, l(j′)) = j′ for all vertices q_j’ of Δ_i. Then, for each of the vertices q_j of Δ_i, Eqs. (6) and (7) become

$L_i{\displaystyle \left(q_l{}_{(j)}\right)}=q_j\text{and}{F}_i{\displaystyle \left(q_l{}_{(j)},z_l{}_{(j)}\right)}=z_j.$E9

Let Ĉ(D) denote the collection of continuous functions f such that f(q_j) = z_j, q_j ∈ ∂D.

Theorem 4.1.2. Suppose the graph associated with ${{\displaystyle \left\{{\Delta }_i\right\}}}_{i=1}^N$ has chromatic number 3. Let L_i and F_i, i = 1, 2, …, N be determined by Eqs. (8) and (9) with s_i = s ∈ (−1, 1). Let Φ be defined by Eq. (5). Then, Φ: Ĉ(D) → Ĉ(D) is well defined and contractive in the sup‐norm ‖⋅‖_∞. Furthermore, for every j = 0, 1, …, m−1 and f ∈ Ĉ(D), Φ(f)(q_j) = z_j.

We call the unique attractor of the afore‐mentioned IFS a self‐affine fractal interpolation surface, or SAFIS for short, with arbitrary boundary data and identical contractivity factors. The colouring of Theorem 4.1.2 is also known as ‘consistent triangulation’ within the context of computational geometry. We will be using this term from now on, for simplicity. Figure 7 illustrates Example 2 constructed in Ref. [4].

Figure 7.

Two views of a SAFIS with arbitrary boundary data and identical contractivity factors.

4.1.3. Arbitrary boundary data, arbitrary contractivity factors

Zhao in Ref. [7] deploys a piecewise linear function γ over the graph G of the IFS. Let G be the graph mentioned earlier and S its vertices. We assign a value s ∈ (−1, 1) to each vertex by defining for all i = 1, 2, …, N the piecewise linear function γ = γ_i: $\mathbb{R}$²→ (−1, 1) as follows

${\gamma }_i(q)={\displaystyle \{}\begin{array}{cc}{\beta }_s(q_{k(i,0),}{q}_{k(i,1)},q_{k(i,2)},q),& q\in {\Delta }_i\\ 0,& \text{otherwise}\end{array},$E119

where β_s is the barycentric interpolation function with the attribute s as coefficient

${\beta }_s(q_a,q_b,q_c,q)=\frac{\Vert (q-q_a)\times (q-q_b)\Vert s_c+\Vert (q-q_b)\times (q-q_c)\Vert s_a+\Vert (q-q_c)\times (q-q_a)\Vert s_b}{\Vert (q_a-q_b)\times (q_c-q_b)\Vert }.$E120

Theorem 4.1.3. Consider the IFS {S; w_1−N} and the corresponding graph G described above that also fulfils the conditions mentioned in Theorem 4.1.2, except for the usage of identical scaling factors over the whole surface. For each transformation w_i, we substitute the scaling factors s_i with the function γ_i described above, so the map F_i becomes for all i = 1, 2, …, N

$F_i(x,y,z)=e_{i}x+f_{i}y+{\gamma }_i(L_i(x,y))z+g_i.$E121

Then, the unique attractor of the IFS mentioned previously is the graph of a continuous function f that passes though the points $(q_j,z_j),q_j\in S$, j = 0, 1, …, m−1.

We call the unique attractor of the afore‐mentioned IFS a self‐similar fractal interpolation surface, or SSFIS for short, with arbitrary boundary data and contractivity factors. The assignment of the scaling factor value s_j, j = 0, 1, …, m−1 to each vertex can be done either during the parameter identification process or can be inferred via a given set of vertical scaling factors {s₁, s₂, …, s_N}, where each s_i corresponds to Δ_i, i = 1, 2, …, N. The technique is identical to the calculation of the vertex normal vectors on the Gouraud and Phong shading models, where the facets of the mesh are the sub‐regions Δ_i, i = 1, 2, …, N. Compare and contrast Figures 6 and 9 of Ref. [7] with Figure 8.

Figure 8.

Two views of a SSFIS with arbitrary boundary data and contractivity factors.

Figure 9.

Domain for fractal interpolating surfaces over polygonal lattice and possible subdomains.

4.2. Recurrent affine fractal interpolation

Let D be a polygonal domain, ${{\displaystyle \left\{{\Delta }_i\right\}}}_{i=1}^N$ a triangulation of D consisting of non‐degenerate triangles with non‐intersecting interiors whose union is D and S = {q₀, q₁, …, q_m−1} be the set of vertices of ${{\displaystyle \left\{{\Delta }_i\right\}}}_{i=1}^N$. Let ${{\displaystyle \left\{{\delta }_k\right\}}}_{k=1}^M$ be M triangles each of which is a union of some subset of ${{\displaystyle \left\{{\Delta }_i\right\}}}_{i=1}^N$; see Figure 9. We order S so that the first n points {q₀, q₁, …, q_n−1} ⊂ S are the vertices of ${{\displaystyle \left\{{\delta }_k\right\}}}_{k=1}^M$.

4.2.1. Coplanar boundary data, arbitrary contractivity factors

We call l: {1, 2, …,N} → {1, 2, …, M} a Δ‐labelling which associates the vertices of Δ_i with the vertices of some δ_k. If P_k are non‐vertical planes in $\mathbb{R}$³, let Ĉ(D) denote the collection of continuous functions f_k = f$|$_Δi: δ_k→ $\mathbb{R}$ such that (x, f_k(x)) ∈ P_k for all x ∈ ∂δ_k.

Theorem 4.2.1. Suppose the points {(q_j, z_j) : q_j ∈ ∂δ_k} are contained in the planes P_k ⊂ $\mathbb{R}$³ for every k. Let Φ be defined by Eq. (5), where L_i and F_i, i = 1, 2, …, N are defined by Eqs. (6)–(8). Then, Φ: Ĉ(D) → Ĉ(D) is well defined and contractive in the sup‐norm ‖⋅‖_∞. Furthermore, for every j = 0, 1, …, m−1 and f ∈ Ĉ(D), Φ(f)(q_j) = z_j.

We call the unique attractor of the afore‐mentioned RIFS a recurrent self‐affine fractal interpolation surface, or RSAFIS for short, with coplanar boundary data and arbitrary contractivity factors.

4.2.2. Arbitrary boundary data, identical contractivity factors

We call l: {0, 1, …, m−1} → {0, 1, …, n−1} a Δ‐labelling associated with ${{\displaystyle \left\{{\delta }_k\right\}}}_{k=1}^M$ and ${{\displaystyle \left\{{\Delta }_i\right\}}}_{i=1}^N$ if {q_l(j), q_l(j'), q_l(j”)} are the vertices of some δ_k whenever {q_j, q_j', q_j”}are the vertices of some Δ_i; see Figure 10.

Figure 10.

Domain for fractal interpolating surfaces over polygonal lattices and possible subdomains.

Theorem 4.2.2. Let l be a Δ‐labelling associated with the triangulations ${{\displaystyle \left\{{\Delta }_i\right\}}}_{i=1}^N$ and ${{\displaystyle \left\{{\delta }_k\right\}}}_{k=1}^M$. Let L_i: $\mathbb{R}$² → $\mathbb{R}$² and F_i: $\mathbb{R}$³→$\mathbb{R}$ be the unique affine maps satisfying (9) with s_i = s (|s| < 1). Let Φ be defined by Eq. (5). Then, Φ: Ĉ(D)→Ĉ(D) is well defined and contractive in the sup‐norm ‖⋅‖_∞. Furthermore, for every j = 0, 1, …, m−1 and f ∈ Ĉ(D), Φ(f)(q_j) = z_j.

Geronimo and Hardin in Ref. [4] have showcased a two‐dimensional multi‐resolution analysis. Based on their work, we have created acceptable colourings of graphs with any desirable density, by solving the problem on a small instance that we call pattern graph. If the solution has the identity of self‐similarity, the density of the graph is then enhanced by interpolating the pattern with a RIFS defined on the metric space of the undirected graphs. With the help of a geometric predicate, we unify the resulting points into a single graph that now has the desired density. The results are illustrated in Figure 11.

Figure 11.

Examples of increasing the density of partially self‐similar coloured graphs.

We call the unique attractor of the afore‐mentioned RIFS a recurrent self‐affine fractal interpolation surface, or RSAFIS for short, with arbitrary boundary data and identical contractivity factors. Figure 12 illustrates Example 3 constructed in Ref. [4].

Figure 12.

Two views of a RSAFIS over polygonal domain with arbitrary boundary data and identical contractivity factors.

4.2.3. Arbitrary boundary data, arbitrary contractivity factors

The map F_i becomes for all i = 1, 2, …, N, $F_i(x,y,z)=e_{i}x+f_{i}y+{\gamma }_i(L_i(x,y))z+g_i.$

Theorem 4.2.3. Suppose the graph associated with ${{\displaystyle \left\{{\Delta }_i\right\}}}_{i=1}^N$ has chromatic number 3. Let L_i and F_i, i = 1, 2, …, N are defined by Eqs. (6)–(8) and Φ be defined by Eq. (5). Then, Φ: Ĉ(D)→Ĉ(D) is well defined and contractive in the sup‐norm ‖⋅‖_∞. Furthermore, Φ(f)(q_j) = z_j, j = 0, 1, …, m−1 and f ∈ Ĉ(D).

We call the unique attractor of the afore‐mentioned RIFS a recurrent self‐similar fractal interpolation surface, or RSSFIS for short, with arbitrary boundary data and contractivity factors. We illustrate this construction in Figure 13.

Figure 13.

Two levels of a RSSFIS over different polygonal domains with arbitrary boundary data and contractivity factors.

5. Conclusions, extensions and future work

We have presented an overview of affine interpolation as well as of recurrent affine interpolation using fractal functions in one and two dimensions. The effectiveness of a self‐affine fractal model is limited to the types of data: only those data that are self‐affine, or nearly so, are well represented. Since most data are not self‐affine, the piecewise or recurrent self‐affine fractal model may be a suitable alternative tool.

Nevertheless, our future work should focus on the parameter identification of fractal interpolation surfaces. More methods must be explored, including the ones for self‐affine FIS’s and the vertex‐based ones, while a better sampling technique for the heights is needed. For the recurrent FIS’s, it is utmost important to connect the domains with the subdomains that they resemble the most instead of determining the connections through an arbitrary colouring scheme. Therefore, a premise for future work is to find the optimal connections between domains and subdomains through a statistical analysis and then, starting from a feasible solution, to perform a heuristic local search with the intention to minimise the distance between current connections, implied by the current colouring, and those of the optimal solution. Moreover, methods to construct FIS, other than those based on the iterated function systems, exist, including the ones that use wavelets and tensor products, in which they should be also studied.

Many areas of fractal functions and their applications are yet to be explored, for instance, calculating the Hausdorff dimension of a general FIF is a challenging open problem. By believing that the field of fractal functions has bright future, the reader, in order to be able to contribute to it, should leave the idea of staying in the comfort of well‐known classical approximation theory and start enjoying the benefits of the more versatile fractal approximants, to supplement the former if not to replace it.