Preserving convexity through rational cubic spline fractal interpolation function

Journal of Computational and Applied Mathematics 263 (2014) 262–276

Contents lists available at ScienceDirect

Journal of Computational and Applied Mathematics journal homepage: www.elsevier.com/locate/cam

Preserving convexity through rational cubic spline fractal interpolation function P. Viswanathan a , A.K.B. Chand a,∗ , R.P. Agarwal b,c a

Department of Mathematics, Indian Institute of Technology Madras, Chennai - 600036, India

b

Department of Mathematics, Texas A&M University - Kingsville, TX 78363, USA

c

Department of Mathematics, Faculty of Science, King Abdulaziz University, 21589 Jeddah, Saudi Arabia

article

info

Article history: Received 15 July 2013 Received in revised form 7 November 2013 MSC: 65D05 65D07 65D10 41A20 28A80 Keywords: Spline Rational spline Rational fractal interpolation Convergence Convexity

abstract We propose a new type of C 1 -rational cubic spline Fractal Interpolation Function (FIF) for convexity preserving univariate interpolation. The associated Iterated Function System P (x) (IFS) involves rational functions of the form Qn (x) , where Pn (x) are cubic polynomials n determined through the Hermite interpolation conditions of the FIF and Qn (x) are preassigned quadratic polynomials with two shape parameters. The rational cubic spline FIF converges to the original function Φ as rapidly as the rth power of the mesh norm approaches to zero, provided Φ (r ) is continuous for r = 1 or 2 and certain mild conditions on the scaling factors are imposed. Furthermore, suitable values for the rational IFS parameters are identified so that the property of convexity carries from the data set to the rational cubic FIFs. In contrast to the classical non-recursive convexity preserving interpolation schemes, the present fractal scheme is well suited for the approximation of a convex function Φ whose derivative is continuous but has varying irregularity. © 2013 Elsevier B.V. All rights reserved.

1. Background and preliminaries Suppose a set of data points D = {(xn , yn ) ∈ I × R : n = 1, 2, . . . , N } is given, where x1 < x2 < · · · < xN and I = [x1 , xN ]. The problem of interpolation in numerical analysis and approximation theory deals with the construction of a continuous function S : I → R satisfying S (xn ) = yn for n = 1, 2, . . . , N. The interpolants produced by these traditional methods are smooth, sometimes infinitely (piecewise) differentiable. As a consequence, these methods become unsuitable for interpolating a highly irregular data or approximating a function whose derivative of a certain order is irregular in a dense subset of the interpolation interval. This served as a motivation for the development of a new interpolation technique using fractal methodology. The theory of fractals and fractal interpolation functions has evolved beyond its mathematical framework and has become a powerful tool in the applied sciences as well as engineering [1–3]. The realm of applications of fractals and FIFs includes but not limited to geometric design [4], structural mechanics [5], physics and chemistry [6,7], signal processing and decoding [8,9], and applied wavelet theory [10]. The reason for the vast applications of FIFs is attributable to their ability to produce complicated mathematical structures with a simple recursive procedure. A FIF is constructed as a fixed point of the Read–Bajraktarević operator defined on a suitable function space. In what follows, we shall recall the precise definition of a FIF and its construction proposed by Barnsley [2].

∗

Corresponding author. Tel.: +91 44 22574629; fax: +91 44 22574602. E-mail addresses: [email protected] (P. Viswanathan), [email protected] (A.K.B. Chand), [email protected] (R.P. Agarwal).

0377-0427/$ – see front matter © 2013 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.cam.2013.11.024

P. Viswanathan et al. / Journal of Computational and Applied Mathematics 263 (2014) 262–276

263

Set In = [xn , xn+1 ] for n ∈ J = {1, 2, . . . , N − 1}. Let Ln : I → In , n ∈ J, be contraction homeomorphisms such that: Ln (x1 ) = xn ,

Ln (xN ) = xn+1 .

(1.1)

Let −1 < αn < 1, n ∈ J. Further, let K = I × [a, b] for some −∞ < a < b < +∞, and N − 1 continuous mappings Fn : K → [a, b] be given satisfying: Fn (x1 , y1 ) = yn ,

Fn (xN , yN ) = yn+1 , n ∈ J , |Fn (x, yl ) − Fn (x, ym )| ≤ |αn ||yl − ym |, x ∈ I , yl , ym ∈ [a, b].

(1.2)

Define functions wn : K → K , wn (x, y) = (Ln (x), Fn (x, y)) ∀ n ∈ J. Consider the collection I ≡ {K ; wn : n ∈ J }, which is termed as an IFS. Associated with the collection of functions in I, there is a set valued mapping w from the hyperspace H (K ) of nonempty compact subsets of K into itself. More precisely, w(A) = ∪n∈J wn (A) for A ∈ H (K ), where wn (A) = {wn (a) : a ∈ A}. There exists a metric h, called the Hausdorff metric, which completes H (K ). This metric is defined as h(A, B) = max{maxx∈A miny∈B d(x, y), maxy∈B minx∈A d(y, x)} ∀A, B ∈ H (K ). Here d is a metric that is equivalent to the Euclidean metric on R2 with respect to which each wn is a contraction. It is well known [2] that w is a contraction on the complete metric space (H (K ), h). Consequently, by the Banach Fixed Point Theorem there exists a unique set G such that G = limn→∞ w n (A0 ) and w(G) = G, where A0 ∈ H (K ) is arbitrary. Here w n denotes the n-fold composition of w , and the limit is taken in the Hausdorff metric: w n (A0 ) → G ⇔ h(G, w n (A0 )) → 0. Such a set G is called an attractor or a deterministic fractal. The next proposition relates G with a function interpolating the data set. Proposition 1.1 (Barnsley [1]). The IFS {K ; wn : n ∈ J } defined above admits a unique attractor G. Further, G is the graph of a continuous function S : I → R which obeys S (xn ) = yn for n = 1, 2, . . . , N. The function S whose graph is the attractor of an IFS as described in Proposition 1.1 is called a FIF. Now we provide some excerpts from the proof of the above proposition that yield a functional equation corresponding to the interpolant S. Let G := {H : I → R | H is continuous, H (x1 ) = y1 and H (xN ) = yN }. Then G endowed with the uniform metric τ (H , H ∗ ) := max{|H (x) − H ∗ (x)| : x ∈ I } is a complete metric space. Define the Read–Bajraktarević operator T on (G, τ ) as follows:

 1  −1 (TH )(x) = Fn L− n (x), H ◦ Ln (x) ,

x ∈ In , n ∈ J .

(1.3)

Due to the conditions on the maps Ln and Fn , n ∈ J, it follows that TH is continuous on I. Furthermore, the map T is a contraction on the metric space (G, τ ), i.e., τ (TH , TH ∗ ) ≤ |α|∞ τ (H , H ∗ ), where |α|∞ := max{|αn | : n ∈ J } < 1. Since (G, τ ) is complete, T possesses a unique fixed point S, i.e., there is S ∈ G such that (TS )(x) = S (x) ∀ x ∈ I. This function S is the FIF corresponding to the IFS I. Therefore, from (1.3) it follows that S satisfies the functional equation:

 1  −1 S (x) = Fn L− n (x), S ◦ Ln (x) ,

x ∈ In , n ∈ J .

(1.4)

The most extensively studied FIFs so far in the literature stem from the IFS K ; wn (x, y) ≡ Ln (x) = an x + bn , Fn (x, y) = αn y + Rn (x) : n ∈ J ,









(1.5)

where Rn : I → R are suitable continuous functions, generally polynomials, such that the conditions prescribed in (1.2) are satisfied. The parameter αn is called a scaling factor of the transformation wn , and α = (α1 , α2 , . . . , αN −1 ) is the scale vector of the IFS. The main differences of a FIF with the traditional interpolation techniques consist: (i) in the definition in terms of a functional equation that implies a self similarity in small scales, (ii) in the constructive way through iterations that is used to compute the interpolant instead of analytic formulae, (iii) in the presence of free parameters αn that replace the unicity of the traditional interpolant for a fixed set of interpolation data with unicity of the interpolant for a fixed data set and a fixed choice of scale vector, and (iv) in the fact that the fractal dimension of a FIF is, in general, non-integer. To obtain the actual fractal interpolant, one needs to continue the iterations indefinitely. However, a small number of iterations usually yields a close approximation. Next we recall the following result of Barnsley and Harrington [11] for the construction of an IFS that generates smooth interpolants corresponding to a finite set of data points. Proposition 1.2 (Barnsley and Harrington [11]). Let {(xn , yn ) : n = 1, 2, . . . , N } be a given set of data points, where x1 < x2 < · · · < xN . Let Ln (x) = an x + bn , n ∈ J, be the affine functions satisfying (1.1) and Fn (x, y) = αn y + Rn (x), n ∈ J, satisfy (1.2). p Suppose that for some integer p ≥ 0, |αn | < an and Rn ∈ C p [x1 , xN ], n ∈ J. Let Fn,k (x, y) =

αn y + R(nk) (x) akn

,

(k)

y1,k =

R 1 ( x1 ) ak1 − α1

(k)

,

y N ,k =

RN −1 (xN ) akN −1 − αN −1

,

k = 1, 2, . . . , p.

If Fn−1,k (xN , yN ,k ) = Fn,k (x1 , y1,k ) for n = 2, 3, . . . , N − 1 and k = 1, 2, . . . , p, then the IFS {I × R; (Ln (x), Fn (x, y)), n ∈ J } determines a FIF S ∈ C p [x1 , xN ], and S (k) is the FIF determined by {I × R; (Ln (x), Fn,k (x, y)), n ∈ J } for k = 1, 2, . . . , p.

264

P. Viswanathan et al. / Journal of Computational and Applied Mathematics 263 (2014) 262–276

These differentiable FIFs (or fractal splines) cannot actually be fractals. However, the name fractal interpolation function is retained because of the flavor of the scaling in its definition and because of the fact that certain derivative of this function is typically a fractal [11]. Further, the graph of a fractal spline is a union of its transformed copies, and hence possesses self-similarity, a characteristic feature of the fractal sets. Since most of the classical splines can be generalized through the fractal splines, they constitute an advance in the technique of approximation. Note that classical interpolation splines are generally infinitely differentiable except perhaps at finite number of knots, whereas fractal splines of C r -continuity include interpolants whose rth derivative can be non-differentiable in a dense set of points of the interpolation interval. In addition, if the experimental data are approximated by a C r -FIF S, then one can aptly use the fractal dimension of S (r ) as a quantitative parameter for the analysis of underlying physical phenomenon. Due to their simplicity, versatility, and flexibility, the smooth FIFs (where Rn are polynomial functions) and their extensions (such as coalescence hidden variable FIFs, fractal interpolation surfaces) have received a great deal of attention in the literature [12–18]. Both the classical polynomial spline interpolation and the fractal polynomial spline interpolation may, in general, ignore the intrinsic form implied by the given data points. Consequently, interpolants produced by these methods may have undesirable inflections or oscillations. To obtain a valid physical interpretation of the underlying process, it is important to develop interpolation schemes that honor the properties inherent in the data, particularly when the data is produced by some scientific phenomena. The problem of searching a sufficiently smooth function that preserves qualitative shape properties inherent in the input data is generally called a shape preserving interpolation problem. Various shape properties are mathematically expressed in terms of conditions such as positivity, monotonicity, and convexity. Owing to the wide applicability, overtime, much effort has been expended by researchers to contribute to the topic of shape preserving interpolation using the traditional polynomial and rational functions (see, for instance, [19–23] and references quoted therein). In the context of shape preserving interpolation, the rational functions are more effective than the polynomial functions. In fact, one of the strong reasons for the importance of the rational interpolation is the perennial and concrete problem of interpolating a set of data efficiently with a function that preserves the inherent shapes. In spite of the promising versatility and flexibility of smooth FIFs, their shape preserving aspects are not well explored hitherto. The present article seeks to fill this gap by picking convexity as an example of shape. On one hand, the smooth fractal interpolants with rational Rn are not well studied so far in the literature. On the other hand, rational functions turn out to be an appropriate medium for an initial exposition of the shape preservation to the fractal interpolation theory. Motivated by these facts, we abandon the common assumption of polynomiality in the FIF and construct a C 1 -rational FIF P (x) with Rn (x) = Qn (x) , where Qn (x) ̸= 0 is a preassigned polynomial involving two shape control parameters. To explicate n the effectiveness of the proposed interpolation, we carry out its convergence analysis. Suitable values of the parameters involved in the IFS are identified so that the corresponding FIF remains convex for a given set of convex data. 2. C 1 -rational cubic spline FIF In this section, we apply Barnsley’s result on smooth FIFs (see Proposition 1.2) to construct a new kind of rational cubic spline FIF with C 1 -continuity. Further, we carry out the convergence analysis of the constructed rational FIF. Let {(xn , yn ) : n = 1, 2, . . . , N } be a given set of data points, where x1 < x2 < · · · < xN . Let yn and dn , respectively, be the function value and thederivative value (given or estimated by some standard methods) at the knot point xn . Consider the IFS {K ; Ln (x), Fn (x, y) : n ∈ J }, Ln (x) = an x + bn , Fn (x, y) = αn y + Rn (x), where Rn ∈ C 1 [x1 , xN ] is chosen such that An (1 − θ )3 + Bn θ (1 − θ )2 + Cn θ 2 (1 − θ ) + Dn θ 3 Pn (θ ) = , Rn (x) ≡ R∗n (θ ) = Qn (θ ) (1 − θ )2 vn + 2un vn θ (1 − θ ) + θ 2 un

θ =

and un , vn are free shape parameters. Adapting the procedure for construction of a FIF given in the previous section, we obtain the following functional equation for the desired interpolant: x−x1 x N −x 1

S Ln (x) = Fn x, S (x) = αn S (x) + Rn (x).









(2.1)

In accordance with the principle of construction of a C 1 -spline FIF (see Proposition 1.2), we impose condition on the scaling factor to be |αn | < an , n ∈ J. Adhering to the notation in Proposition 1.2, let y1,1 = d1 , yN ,1 = dN , Fn,1 (x1 , d1 ) = dn , and Fn,1 (xN , dN ) = dn+1 , n ∈ J. Then, S turns out to be differentiable. The derivative S ′ is the fractal function corresponding to   the IFS {I × R; Ln (x), Fn,1 (x, y) : n ∈ J }, and interpolates {(xn , dn ) : n = 1, 2, . . . , N }. Further, S ′ fulfills the functional equation: S ′ Ln (x) = Fn,1 x, S ′ (x) =









αn S ′ (x) + R′n (x) an

.

(2.2)

The constants An , Bn , Cn , and Dn appearing in the expression are evaluated based on the Hermite interpolatory conditions, viz., S (xn ) = yn , S (xn+1 ) = yn+1 , S ′ (xn ) = dn , and, S ′ (xn+1 ) = dn+1 (these are equivalent to Fn (x1 , y1 ) = yn , Fn (xN , yN ) = yn+1 , Fn,1 (x1 , d1 ) = dn , and Fn,1 (xN , dN ) = dn+1 , and hence satisfy the hypotheses of the Barnsley–Harrington theorem).

P. Viswanathan et al. / Journal of Computational and Applied Mathematics 263 (2014) 262–276

265

Substituting the Hermite interpolation conditions in (2.1)–(2.2), letting hn = xn+1 − xn , and performing some routine algebra, we obtain the desired rational cubic spline FIF as follows: S Ln (x) = αn S (x) +





Pn (θ ) Qn (θ )

,

(2.3)

 Pn (θ ) = {yn − αn y1 }vn (1 − θ )3 + {yn+1 − αn yN }un θ 3 + (2un vn + vn )yn + vn hn dn   − αn [(2un vn + vn )y1 + vn (xN − x1 )d1 ] θ (1 − θ )2 + (2un vn + un )yn+1  2 − un hn dn+1 − αn [(2un vn + un )yN − un (xN − x1 )dN ] θ (1 − θ ), Qn (θ ) = (1 − θ )2 vn + 2un vn θ (1 − θ ) + θ 2 un . The application of our interpolation scheme requires knowledge of the derivative parameters dn (n = 1, 2, . . . , N ). In the absence of other information, the derivative parameters must be estimated from the data. Methods that associate derivatives with the data points involve estimates based on nearby slopes or data differences. Depending on the applications, various schemes based on linear combination (e.g., arithmetic mean method) or multiplicative combination (e.g., geometric mean method) of chord-slopes are developed in the literature (see, for instance, [24]). In this article, we use the following arithmetic mean method: dn =

hn ∆n−1 + hn−1 ∆n hn−1 + hn

 d1 =

1+

 dN =

1+

h1

 ∆1 −

h2

hN −1 hN −2

,

yn+1 − yn

∆n =

h1 y3 − y1 h2 x 3 − x 1

 ∆N −1 −

hn

, n = 2, 3, . . . , N − 1;

;

hN −1 yN − yN −2 hN −2 x N − x N −2

.

This is a good time to draw the reader’s attention in regard to our choice of Rn in the IFS. The choice of cubic polynomial in the numerator gives four degrees of freedom to impose the Hermite interpolation conditions, and consequently provides C 1 -continuity. Prescribing the polynomial in the denominator on a priori grounds avoids the possibility of nonlinearity in the system governing the coefficients of the rational expression. Presence of the shape parameters un and vn gives tension property to the interpolant. That is to say that for increasing values of the shape parameters un and vn , indeed with scaling parameters αn approaching to zero, the rational fractal interpolant approaches to a piecewise linear interpolant. This property ensures that in the limiting configuration, the interpolant preserves the shape inherent in the data set. Consequently, the study of convexity property of the fractal interpolant reduces to the development of suitable conditions on the scaling factors and the shape parameters that answers the question: to what extent the scaling parameters can be increased from zero and the shape parameters can be lowered from the theoretical limiting value infinity so that the resulting fractal interpolant preserves convexity without compromising for the C 1 -continuity? Of course, various other choices for the degrees of numerator/denominator and the number of shape parameters can be considered depending on the application. Remark 2.1. If αn = 0, for all n ∈ J, then the rational spline FIF defined above reduces to the classical rational cubic spline C defined as follows. For x ∈ [xn , xn+1 ], C (x) =

Un (Θ ) Vn (Θ )

,

Θ=

x − xn xn+1 − xn

,

(2.4)

Un (Θ ) = vn yn (1 − Θ )3 + [(2un vn + vn )yn + vn hn dn ]Θ (1 − Θ )2

+ [(2un vn + un )yn+1 − un hn dn+1 ]Θ 2 (1 − Θ ) + un yn+1 Θ 3 , Vn (Θ ) = (1 − Θ )2 vn + 2un vn Θ (1 − Θ ) + Θ 2 un . This classical rational spline is studied in detail in [23]. If for all n ∈ J , αn = 0 and un = vn = 1, then the rational cubic spline FIF reduces to the standard cubic Hermite interpolant. Remark 2.2. Taking αn = 0 and un = vn = C ∗ (x) =

Wn (Θ ) Zn (Θ )

,

Θ=

x − xn xn+1 − xn

rn −1 , 2

our rational cubic spline FIF reduces to

,

(2.5)

Wn (Θ ) = yn (1 − Θ )3 + (rn yn + hn dn )Θ (1 − Θ )2 + (rn yn+1 − hn dn+1 )Θ 2 (1 − Θ ) + yn+1 Θ 3 , Zn (Θ ) = 1 + (rn − 3)Θ (1 − Θ ). In this special case, rn > −1 ensures the strict positivity of the denominator. The aforementioned one parameter family of C 1 -rational cubic spline is introduced and investigated by Delbourgo and Gregory [25].

266

P. Viswanathan et al. / Journal of Computational and Applied Mathematics 263 (2014) 262–276

Having introduced a new interpolation scheme, it is natural to query on its effectiveness in the approximation of a function, and it is to this that we now turn. The following theorem establishes that for suitable values of the scaling factors and for Φ ∈ C r , ∥Φ − S ∥∞ = O(hr ), r = 1, 2. Theorem 2.1. Let {(xn , yn ) : n = 1, 2, . . . , N } be a given set of interpolation data generated by the original function Φ . Suppose that dn is the derivative value of the original function Φ at the knot xn or the derivative value estimated from the data points yn = Φ (xn ) using some linear approximation methods. Let S and C , respectively, be the corresponding rational cubic fractal P (α ,x) spline and its classical counterpart. Further, let the rational functions Rn (αn , x) = nQ (nx) (where the dependence with un and vn

n    ∂ Rn  is suppressed for convenience) involved in the IFS that generate the FIF S satisfy  ∂α (ξ , x) ≤ D0 ∀ (ξ , x) ∈ (−an , an )× In , n ∈ J, n and h := max{hn : n ∈ J } be the mesh norm. Then for suitable constants c and c ∗ :

∞ (∥C ∥∞ + D0 ). (i) For Φ ∈ C 1 [x1 , xN ], we have ∥Φ − S ∥∞ ≤ c h ∥Φ ′ ∥∞ + 1−|α| ∞ |α|∞ 2 ∗ 2 ′′ (ii) For Φ ∈ C [x1 , xN ], we have ∥Φ − S ∥∞ ≤ c h ∥Φ ∥∞ + 1−|α| (∥C ∥∞ + D0 ).

|α|

∞

Proof. By triangle inequality, we have

∥Φ − S ∥∞ ≤ ∥Φ − C ∥∞ + ∥C − S ∥∞ .

(2.6)

In Refs. [13,18], the following upper bound for the uniform distance between the classical interpolant C and its fractal perturbation S is established:

∥S − C ∥∞ ≤

|α|∞ 1 − |α|∞

(∥C ∥∞ + D0 ).

(2.7)

Now we turn our attention to the first summand occurring in the right hand side of (2.6), viz., the uniform error bound for the classical rational cubic spline C . Since the interpolant C is local, without loss of generality, it is enough to consider the uniform error bound in the subinterval [xn , xn+1 ]. For Φ ∈ C 1 and for a fixed x ∈ [xn , xn+1 ], let us consider the error E (Φ ; x) = Φ (x) − C (x) as a linear functional which operates on Φ and which annihilates all elements of P1 , the space of polynomials of degree at most one. By the Peano Kernel Theorem [26]: L[Φ ] = E (Φ ; x) = Φ (x) − C (x) =



xn+1

Φ ′ (τ )Lx [(x − τ )0+ ]dτ .

(2.8)

xn

The kernel function is given by Lx [(x − τ )0+ ] =



r (τ , x) if xn < τ < x, s(τ , x) if x < τ < xn+1 .

Here the notation Lx is used to emphasize that the functional L is applied to the expression (x − τ )n+ considered as a function of x. Also, (x − τ )n+ is the truncated power function:

(x − τ )+ := n



(x − τ )n if τ < x, 0 if τ > x,

vn (1 − Θ )2 (1 + 2un Θ ) , (1 − Θ )2 vn + 2un vn Θ (1 − Θ ) + Θ 2 un −Θ 2 un (1 + 2vn (1 − Θ )) s(τ , x) = . (1 − Θ )2 vn + 2un vn Θ (1 − Θ ) + Θ 2 un

r (τ , x) =

Using the expression of kernel function, (2.8) is rewritten in successive self-explanatory stages as:

   xn+1    ′ 0 |Φ (x) − C (x)| =  Φ (τ )Lx [(x − τ )+ ]dτ ,  xn   xn+1 ≤ ∥Φ ′ ∥ |Lx [(x − τ )0+ ]dτ |, xn

≤ ∥Φ ∥ ′



x

r (τ , x)dτ − xn

xn+1

 x

 s(τ , x)dτ ,

P. Viswanathan et al. / Journal of Computational and Applied Mathematics 263 (2014) 262–276

267

where ∥ · ∥ denotes the uniform norm on [xn , xn+1 ]. Performing the required integration, we obtain the following local error bound:

  vn Θ (1 − Θ )2 (1 + 2un Θ ) + un Θ 2 (1 − Θ ) 1 + 2vn (1 − Θ ) . |Φ (x) − C (x)| ≤ ∥Φ ∥hn (1 − Θ )2 vn + 2un vn Θ (1 − Θ ) + Θ 2 un ′

Let c := maxn∈J max0≤Θ ≤1

(2.9)

vn Θ (1−Θ )2 (1+2un Θ )+un Θ 2 (1−Θ )(1+2vn (1−Θ )) . (1−Θ )2 vn +2un vn Θ (1−Θ )+Θ 2 un

The local error bound (2.9) for the classical rational interpolant C along with the aforementioned notation suggests that

∥Φ − C ∥∞ ≤ c h ∥Φ ′ ∥∞ .

(2.10)

The triangle inequality (2.6) coupled with (2.7) and (2.10) proves the first part of the theorem. Next we assume that Φ ∈ C 2 [x1 , xN ]. For x ∈ [xn , xn+1 ], invoking the Peano Kernel Theorem once again, L[Φ ] = Φ (x) − C (x) =

xn+1



Φ ′′ (τ )Lx [(x − τ )+ ]dτ ,

(2.11)

xn

where Lx [(x − τ )+ ] =



r (τ , x) if xn < τ < x, s(τ , x) if x < τ < xn+1 ,

r (τ , x) = (x − τ ) − s(τ , x) =

Θ 2 (1 − Θ )[(2un vn + un )(xn+1 − τ ) − un hn ] + Θ 3 un (xn+1 − τ ) , (1 − Θ )2 vn + 2un vn Θ (1 − Θ ) + Θ 2 un

−Θ 2 (1 − Θ )[(2un vn + un )(xn+1 − τ ) − un hn ] + Θ 3 un (xn+1 − τ ) . (1 − Θ )2 vn + 2un vn Θ (1 − Θ ) + Θ 2 un

(2.12)

(2.13)

Therefore we have

|Φ (x) − C (x)| ≤ ∥Φ ′′ ∥



xn+1

|Lx [(x − τ )+ ]| dτ .

(2.14)

xn

The integral involved in (2.14) is expressed as xn+1



|Lx [(x − τ )+ ]| dτ =

xn



x

|r (τ , x)| dτ +

xn+1



xn

|s(τ , x)| dτ .

x

To study the sign changes of the kernel function, we observe the following properties of r (τ , x) and s(τ , x): Roots of r (x, x), s(x, x) in [0, 1] are Θ = 0 and Θ = 1, 2

n un Θ Root of r (τ , x) is τ ∗ = x − 2h , and τ ∗ ∈ [xn , x], 1+2un Θ h (1−Θ ) Root of s(τ , x) is τ∗ = xn+1 − 1+2nv (1−Θ ) , and τ∗ ∈ [x, xn+1 ]. n

The expressions for r (τ , x) and s(τ , x) given in (2.12) and (2.13) are simplified as

(1 − Θ )2 vn (1 + 2un Θ )(τ ∗ − τ ) , Vn (Θ ) (2un vn Θ 2 (1 − Θ ) + un Θ 2 )(τ − τ∗ ) s(τ , x) = , Vn (Θ )

r (τ , x) =

respectively, where Vn (Θ ) is given in (2.4). With these observations, we get x



|r (τ , x)| dτ = xn

τ∗



r (τ , x) dτ −

xn



x

τ∗

r (τ , x) dτ ,

Θ 2 (1 − Θ )2 h2n vn 2Θ 4 (1 − Θ )2 h2n u2n vn + , 2Vn (Θ )[1 + 2un Θ ] Vn (Θ )[1 + 2un Θ ]  xn+1  τ∗  xn+1 |s(τ , x)| dτ = − s(τ , x) dτ + s(τ , x) dτ ,

=

x

x

=

τ∗

2Θ 2 (1 − Θ )4 h2n vn2 un

Vn (Θ )[1 + 2vn (1 − Θ )]

+

Θ 2 (1 − Θ )2 h2n un . 2Vn (Θ )[1 + 2vn (1 − Θ )]

(2.15)

268

P. Viswanathan et al. / Journal of Computational and Applied Mathematics 263 (2014) 262–276

Using the above estimates, (2.14)–(2.15) give

|Φ (x) − C (x)| ≤ cn∗ h2n ∥Φ ′′ ∥ where cn∗ := max{w(un , vn , Θ ) : 0 ≤ Θ ≤ 1}, w(un , vn , Θ ) =

Θ 2 (1 − Θ )2 un + 4Θ 2 (1 − Θ )4 un vn2 Θ 2 (1 − Θ )2 vn + 4Θ 4 (1 − Θ )2 u2n vn + . 2Vn (Θ )[1 + 2un Θ ] 2Vn (Θ )[1 + 2vn (1 − Θ )]

Since the above inequality is true for x ∈ [xn , xn+1 ], n ∈ J, the desired error bound is given by

∥Φ − C ∥∞ ≤ c ∗ h2 ∥Φ ′′ ∥∞ ,

where c ∗ := max{cn∗ : n ∈ J }.

(2.16)

The inequality (2.6) along with (2.7) and (2.16) establishes the second part of the theorem. Consequences: By the hypotheses of the theorem on differentiability of FIFs: |αn | ≤ κ an = |α|∞ 1−|α|∞

≤

κh . (xN −x1 )−κ h

 κ hn xN −x1

for all n ∈ J, and therefore

Hence the above theorem states, in particular, that

(i) If Φ ∈ C [x1 , xN ] and |αn | ≤ κ an for all n ∈ J, then ∥Φ − S ∥∞ = O(h). (ii) If Φ ∈ C 2 [x1 , xN ] and |αn | ≤ κ a2n for all n ∈ J, then ∥Φ − S ∥∞ = O(h2 ). 1

Remark 2.3. Since the rational cubic spline FIF does not possess a closed form expression, the standard methods such as a Taylor series analysis, the Cauchy remainder formula, and the Peano Kernel theorem cannot be easily adapted for analyzing its convergence. Instead, we have used the error bound for the classical rational cubic spline via the triangle inequality ∥Φ − S ∥∞ ≤ ∥Φ − C ∥∞ + ∥C − S ∥∞ to establish that S has the same order of convergence as that of C . As a consequence, it is possible to approximate any regular function by using a rational cubic spline FIF with arbitrary accuracy. Application of triangle inequality to obtain the convergence does not imply that the error committed by the rational cubic spline FIF in approximating a smooth function will be always greater than the error committed by its classical counterpart (see Section 4 for a numerical illustration). Furthermore, we feel that any possible loss of precision is to be counteracted with the generality offered by the method. 3. Identification of parameters for convex rational FIFs The goal of this section is to identify suitable values for the parameters involved in the rational IFS (cf. Section 2) so that the resulting FIF is convex for given convex data. For the sake of definiteness, we assume a strictly convex set of y −y data {(xn , yn ) : n = 1, 2, . . . , N } so that ∆1 < ∆2 < · · · < ∆N −1 , where ∆n := n+h1 n , n ∈ J. To have a convex n interpolant S that avoids possibility of having straight line segments, it is necessary that the derivative parameters satisfy d1 < ∆1 < · · · < dn < ∆n < dn+1 < · · · < ∆N −1 < dN . Recall that the convexity of a univariate function S is derived generally using the following standard results: (i) a twice differentiable function S is convex on I if and only if S ′′ (x) ≥ 0, or more generally (ii) if there exists a partition a = x1 < x2 < · · · < xN = b of I = [a, b] such that S is continuous on [a, b], S is of class C 2 on each subinterval (xn , xn+1 ), ′ + and S has one sided derivative at x2 , x3 , . . . , xN −1 satisfying S ′ (x− n ) ≤ S (xn ) for n = 2, 3, . . . , N − 1, then S is convex on I. However, the regularity of the developed rational cubic spline FIF S is, in general, restricted to C 1 , and S ′ is a fractal function that may be nowhere differentiable. Therefore S ′′ is not even piecewise C 2 , and the convexity of S cannot be deduced from these results. Instead, we employ the following results: (i) a differentiable univariate function is convex on an interval if and only if its derivative is monotonic increasing on that interval, (ii) let S ∈ C [a, b]. If for each x ∈ (a, b) one of the one sided derivative S ′ (x+ ) or S ′ (x− ) exists and is nonnegative (possibly +∞), then S is monotone increasing on [a, b] (see [27]). Owing to these results, for the convexity of S ∈ C 1 [x1 , xN ], it suffices to show that S ′′ (x+ ) or S ′′ (x− ) exists and is nonnegative for each point in (x1 , xN ). By the principle of construction of twice differentiable FIFs (see [13]), we take |αn | < a2n for all n ∈ J. Informally, S ′′ (Ln (x)) =

αn a2n

S ′′ (x) + Tn (x),

(3.1)

where Tn (x) = Tn x1 + θ (xN − x1 )



=



 1 5 4  3 K1n (1 − θ ) + K2n θ (1 − θ ) hn (1 − θ )2 vn + 2un vn θ (1 − θ ) + θ 2 un  + K3n θ 2 (1 − θ )3 + K4n θ 3 (1 − θ )2 + K5n θ 4 (1 − θ ) + K6n θ 5 , 

K1n =

v 

4un n3

∆n −

− dn+1 −

αn hn

αn hn



αn

(yN − y1 ) − dn − d1 (xN − x1 ) hn 

dN (xN − x1 )

,



 +

v

2un n2

∆n −

αn hn

(yN − y1 )

P. Viswanathan et al. / Journal of Computational and Applied Mathematics 263 (2014) 262–276





269



αn αn (yN − y1 ) − dn − d1 (xN − x1 ) + 4un vn2 ∆n − (yN − y1 ) hn hn hn      αn α α n n − dn+1 − dN (xN − x1 ) + 6un vn2 ∆n − (yN − y1 ) − dn − d1 (xN − x1 ) ,

K2n = 8un vn3 ∆n −

αn



hn

hn

 K3n =

v 

6u2n n

d n +1 −

αn hn

 dN (xN − x1 ) − ∆n −



αn hn

hn

 (yN − y1 )

 +



v

∆n −

4un n3

αn hn



(yN − y1 ) 

αn αn αn − d n − d 1 ( xN − x1 ) + 12un vn2 ∆n − (yN − y1 ) − dn − d1 (xN − x1 ) hn hn hn    αn αn + 2un vn2 ∆n − (yN − y1 ) − dn+1 − dN (xN − x1 ) , hn



K4n =

hn



αn

∆n −

 v 

8u3n n

d n +1 −

αn

− ∆n −

hn

v 

4u3n n

αn hn

hn

 dN (xN − x1 ) − ∆n −

(yN − y1 )

d n +1 −

− ∆n −

αn

hn



 K6n =



4u3n n

hn

hn

K5n =



αn

αn (yN − y1 ) − dn − d1 (xN − x1 ) + v dn+1 − dN (xN − x1 ) hn hn     αn αn αn 2 + 12un vn dn+1 − dN (xN − x1 ) − ∆n − (yN − y1 ) − ∆n − (yN − y1 ) hn hn hn    αn αn + 2u2n vn dn − d1 (xN − x1 ) − ∆n − (yN − y1 ) , v 

6un n2

αn hn

αn hn

 (yN − y1 )



αn

+ 6u2n vn dn+1 − 

dN (xN − x1 ) − ∆n −

αn hn

hn

 +

v 

4u2n n

dn −

dN (xN − x1 ) − ∆n −

 (yN − y1 )

αn hn

 +

v

2u2n n

dn −

αn hn

d1 (xN − x1 )

 (yN − y1 ) αn hn

,

d1 (xN − x1 )

 (yN − y1 )

.

We shall formalize the things in the following theorem that provides a set of sufficient conditions on the rational IFS parameters so as to ensure the convexity of the corresponding rational cubic spline FIF. Theorem 3.1. Suppose {(xn , yn ) : n = 1, 2, . . . , N } is a set of strictly convex data, and S is the corresponding rational cubic spline FIF described in (2.3). Assume that the derivative parameters at the knots satisfy d1 < ∆1 < · · · < dn < ∆n < dn+1 < · · · < ∆N −1 < dN . Then the following conditions on the scaling factors and the shape parameters for each n ∈ J are sufficient for the convexity of S on the interpolation interval I = [x1 , xN ]: 0 ≤ αn < min



a2n

,

hn (∆n − dn )

,

hn (dn+1 − ∆n )

yN − y1 − d1 (xN − x1 ) dN (xN − x1 ) − (yN − y1 )

 ∆n − un ≥

αn hn



(yN − y1 ) − dn −



αn hn

d1 (xN − x1 )

 ,



2 dn+1 −

αn hn

dN (xN − x1 ) − ∆n −

αn hn

(yN − y1 )

 d n +1 −

αn hn

vn ≥  2 ∆n −



dN (xN − x1 ) − ∆n −

hn

(yN − y1 ) − dn −

hn

(yN − y1 )  .

 αn

αn

αn hn

d1 (xN − x1 )



,

270

P. Viswanathan et al. / Journal of Computational and Applied Mathematics 263 (2014) 262–276

Proof. Recall that for n ∈ J, the maps Ln : [x1 , xN ] → [xn , xn+1 ] satisfy Ln (x1 ) = xn and Ln (xN ) = xn+1 . Therefore we obtain: S (x1 ) = +

′′

S ′′ (x+ m) =



K11

1−

v13 h1 αm a2m

α1

−1

a21

S ′′ (x+ 1 )+

,

K1m

vm3 hm

S ( xN ) = ′′

,



K 6 ,N − 1

−

u3N −1 hN −1

αN −1

1−

−1 ,

a2N −1

m = 2, 3, . . . , N − 1.

(3.2)

For 0 ≤ αn < a2n , it follows from (3.2) that if K1n ≥ 0 (n ∈ J) and K6,N −1 ≥ 0, then the second derivatives (right-handed) at the knots xn , n ∈ J, and the second derivative (left-handed) at xN are nonnegative. For a knot point xm , m ∈ J we have: S ′′ Ln (xm )+ =





αm a2m

S ′′ (x+ m ) + Tn (xm ).

(3.3)

Whence, with the assumption K1n ≥ 0 for all n ∈ J, we have S ′′ Ln (xm )+ ≥ 0, provided Tn (xm ) ≥ 0. Note that Tn (xm ) ≥ 0 is satisfied if the coefficients Kjn ≥ 0 for j = 1, 2, . . . , 6. From the Three Chords Lemma for convex functions [28], it follows y −y that for a strictly convex interpolant, the end point derivatives should necessarily satisfy d1 < xN −x1 < dN . In view of this,





N

  the following condition on αn implies ∆n − αh n (yN − y1 ) − dn − αh n d1 (xN − x1 ) > 0: n

αn <

1

n

hn (∆n − dn ) . (yN − y1 ) − d1 (xN − x1 )

(3.4)

Note that



αn

K1n ≥ 0 ⇐⇒ 2vn ∆n −

hn

 − dn+1 − 

− ∆n −

hn

αn

⇐⇒ 2vn ∆n − 

αn

hn

αn hn

 (yN − y1 ) − dn −



αn

d1 (xN − x1 )

hn

+ ∆n −

αn hn

(yN − y1 )

 d N ( xN − x1 )

≥0 

(yN − y1 ) − dn −



αn

d1 (xN − x1 )

hn

≥ dn+1 −

αn hn

d N ( xN − x1 )

 (yN − y1 ) .

(3.5)

From (3.4) and (3.5), it follows that K1n ≥ 0 whenever

 dn+1 −

αn hn

 αn

dN (xN − x1 ) − ∆n −

vn ≥ 

hn

(yN − y1 )  .



2 ∆n −

αn hn

(yN − y1 ) − dn − 

αn hn

(3.6)

d1 (xN − x1 )

Since K2n = 2K1n + 6un vn2 ∆n − αh n (yN − y1 )− dn − αh n d1 (xN − x1 )



n

n



, (3.4)–(3.6) give K2n ≥ 0. Now correlating the expression

for K3n with that of K1n , the additional terms involved suggest that K3n ≥ 0 if dn+1 − αh n dN (xN − x1 ) − ∆n − αh n (fN − f1 ) > 0, n n which can be verified from the following condition on the scaling factor:



αn <

hn (dn+1 − ∆n ) dN (xN − x1 ) − (yN − y1 )

.

(3.7)

In view of (3.4) and (3.7),

 K4n ≥ 0

if 2un dn+1 −

 − ∆n −

αn hn

αn hn

 dN (xN − x1 ) − ∆n −

 (yN − y1 ) ≥ 0,



αn hn

 (yN − y1 )

+ dn −

αn hn

d1 (xN − x1 )

P. Viswanathan et al. / Journal of Computational and Applied Mathematics 263 (2014) 262–276

271

Table 1 A convex data generated by a function Φ with Φ ′′ having irregularity in a dense subset of unit interval. x

0

0.5

0.75

1

y

0

0.0783

0.1918

0.55

 i.e., if 2un dn+1 −

 − dn −

αn hn

αn hn

 dN (xN − x1 ) − ∆n −

αn hn

 (yN − y1 )

≥ ∆n −

αn hn

(yN − y1 )

 d 1 ( xN − x1 ) .

(3.8)

Hence, K4n ≥ 0 if

 ∆n − un ≥

αn hn

(yN − y1 ) − dn −

 2 dn+1 −

 αn hn

d1 (xN − x1 )

 .

 αn hn

dN (xN − x1 ) − ∆n −

αn hn

(3.9)

(yN − y1 ) 



  From (3.9) via (3.8), we obtain K6n ≥ 0. Finally, since K5n = 2K6n + 6u2n vn dn+1 − αh n dN (xN − x1 ) − ∆n − αh n (yN − y1 ) , in n

n

aid of earlier assumptions we have K5n ≥ 0. Thus the conditions on the scaling factors and the shapeparameters  prescribed in the theorem ensure Kjn ≥ 0 for all j = 1, 2, . . . , 6. This in turn ensures the nonnegativity of S ′′ Ln (xm )+ for n, m ∈ J, ′′ + ′′ − and S ′′ (x− N ). Consequently, the nonnegativity of S (x ) and S (x ) follows from the iterative nature of a fractal function. The convexity of S can now be deduced by using the result quoted at the beginning of this section.  Remark 3.1. If the given set of data is not strictly convex but ∆n − dn = 0 or dn+1 − ∆n = 0, then we take αn = 0. Now for S ′′ (x) ≥ 0 (see the expressions for the coefficients Kjn ), we take dn = dn+1 = ∆n . Thus we get S Ln (x) = i.e., the interpolant reduces to a straight line segment on the interval [xn , xn+1 ].





(xN −x)yn +(x−x1 )yn+1 , x N −x 1

Remark 3.2. When αn = 0 ∀ n ∈ J, Theorem 3.1 recaptures the sufficient conditions for the convexity of the classical rational cubic spline C (see Remark 2.1) described in [29]. Remark 3.3. Taking all αn = 0 and un = vn = rn 2−1 in Theorem 3.1, we deduce the convexity conditions for the cubic rational spline introduced by Delbourgo et al. [25] as a special case of our result. −∆

n +1 , ∀ n ∈ J, obviously satisfy the convexity conditions Remark 3.4. The choices αn = 0 and un = vn = 2(d∆n −−dn∆ ) + 2(n∆ n n −dn ) n+1 prescribed in the preceding theorem. For these choices of spline parameters, our rational cubic spline FIF reduces to the rational functions with quadratic numerator and linear denominator studied in detail in [30]. Therefore, our scheme includes the convexity preserving rational interpolation by Delbourgo as a very special case.

d

4. Demonstration In this section, we implement our rational convexity preserving interpolation scheme on the data set reported in Table 1. This convex data is taken from a smooth function Φ displayed in Fig. 1(a) whose second derivative Φ ′′ is nondifferentiable in a dense subset of [0, 1] (see Fig. 1(b)). The derivatives at the knots are estimated using the arithmetic mean method (see Section 2). For the present data, we obtain: d1 = −0.0417, d2 = 0.3549, d3 = 0.9434, and d4 = 1.9222, which satisfy the y −y necessary convexity conditions: d1 < ∆1 < d2 < ∆2 < d3 < ∆3 < d4 ; d1 < x4 −x1 < d4 . 4 1 For α1 = 0.24, α2 = −0.06, α3 = 0.03 and the shape parameters as shown in Table 2, the iterations of the rational IFS produce the rational cubic spline FIF (see Fig. 2(a)). Since the scaling factors are selected merely to satisfy |αn | < a2n , n = 1, 2, and 3, the fractal curve in Fig. 2(a) has undesired inflections, and hence fails to be convex. This illustrates the need of an automated algorithm for the selection of the rational IFS parameters that ensure convex fractal curves for prescribed convex data. Next we employ Theorem 3.1 to get suitable values of the scaling factors and the shape parameters that generate convex fractal curves. The computed values of the scaling factors prescribed by Theorem 3.1 are: 0 ≤ α1 < 0.0723, 0 ≤ α2 < 0.0419, 0 ≤ α3 < 0.0625. The rational cubic spline FIF in Fig. 2(b) is generated with the scaling factors and shape parameters (see Table 2) according to Theorem 3.1. To recapture the classical rational cubic spline, we assume zero scaling on each subinterval. The specific values of the shape parameters that are taken to generate the classical rational cubic spline (Fig. 2(c)) is given in Table 2. Thus the presence of the scaling factors in the proposed fractal scheme provides more flexibility and diversity for the choice of a suitable interpolant.

272

P. Viswanathan et al. / Journal of Computational and Applied Mathematics 263 (2014) 262–276

(a) Original function Φ .

(b) Second derivative of Φ .

Fig. 1. The nature of original function for the convex data in Table 1.

Table 2 Parameters corresponding to Fig. 2(a)–(f). Figure no

Choice of rational IFS parameters

αn Fig. 2(a)

un

vn

αn Fig. 2(b)

un

vn

αn Fig. 2(c)

un

vn

αn Fig. 2(d)

un

vn

αn Fig. 2(e)

un

vn

αn Fig. 2(f)

un

vn

0.024 0.1 10

−0.06 0.5 0.5

0.03 0.1 1

0.0722 365.6640 0.0007

0.0418 0.0004 774.8217

0.0624 1.1631 0.2150

0 0.51 0.51

0 0.1013 2.4693

0 0.51 0.5

0.02 1.1390 0.4107

0.02 0.0682 3.6669

0.02 0.5823 0.4294

0.05 1.1390 0.4107

0.02 0.0682 3.6669

0.02 0.5823 0.4294

−0.001 5 5

−0.001 5 5

−0.001 5 5

At first glance, it might appear that the convexity preserving algorithm presented here has a drawback due to its global nature, which might limit its applications. However, numerical experiments conducted on different data sets suggest that the portions of a fractal curve pertaining to other subintervals are not extremely sensitive towards the changes of parameters in a particular subinterval. For instance, Fig. 2(d)–(e) are generated with the same sets of parameters (see Table 2) except for the scaling factor in the first subinterval and the corresponding curves differ only in the first subinterval. Thus our scheme is local or global depending on the magnitude of the scaling factor in each subinterval. In Fig. 2(f), we generate a convex fractal curve with negative scalings in all the subintervals on a trial and error basis, whence the conditions prescribed by Theorem 3.1 are sufficient but not necessary. This provides a scope for further research in developing the most appropriate sufficient conditions on the rational IFS parameters where negative values of the scaling factors are allowed or more generally in developing necessary and sufficient conditions. Recall that the classical rational cubic spline C (cf. Remark 2.1) can also be used to generate a large number of convex curves by varying the shape parameters involved. However, the second derivatives of convex interpolants produced by the classical method are differentiable except possibly at the knots. Consequently, the convexity preserving classical cubic spline interpolants fail to perform well for a data set generated from a function whose second derivative is irregular in a dense set of points. In this capacity, the present rational fractal scheme seems to be more effective, for instance, in the approximation of motion of a body having very irregular acceleration like the motion of a pendulum on a cart [31,32]. In fact, the proposed convexity preserving fractal interpolation scheme works equally good for a data set generated from a function Φ with varying irregularity in Φ ′′ . To illustrate this point, we proceed as follows. Let Si , i = 1, 2, 3, . . . , 6, denote the rational FIFs corresponding to the function Φ that are generated in Fig. 2(a)–(f). The corresponding (one-sided) second derivatives Si′′ , i = 1, 2, 3, . . . , 6 are generated in Fig. 3(a)–(f). The uniform errors

P. Viswanathan et al. / Journal of Computational and Applied Mathematics 263 (2014) 262–276

(a) Nonconvex rational cubic FIF S1 .

(b) Convex rational cubic FIF S2 .

(c) Classical rational cubic spline C = S3 .

(d) Convex rational cubic FIF S4 .

(e) Convex rational cubic FIF S5 .

(f) Convex rational cubic FIF S6 .

273

Fig. 2. Rational cubic spline FIFs for the data set in Table 1 with IFS parameters as in Table 2.

(rounded off to five decimal places) ∥Φ − Si ∥∞ and ∥Φ ′′ − Si′′ ∥∞ calculated with 26 244 sample points (i.e., with 8 iterations of the IFS code) are displayed in Table 3. The fractal interpolant S2 provides a better convex approximant to Φ with a minimum uniform error in the second derivative. As the second derivative of the classical convex rational spline S3 is smooth except at a finite number of points, it remains unsuitable for approximating the nondifferentiable function Φ ′′ (see Figs. 1(b) and 3(c)). This is reflected in the calculated uniform error of the second derivative ∥Φ ′′ − S3′′ ∥∞ (see Table 3). In comparison with the classical convex rational spline S3 , the rational cubic spline FIF S1 gives a better approximant for the second derivative of Φ , however, S1 is not convex. Thus, in general, the diversity offered by the present scheme can be exploited in data visualization, engineering design environments by considering various desirable properties such as convexity, fairness (visual pleasantness), approximation order, localness, and fractality in the second derivative and trade-

274

P. Viswanathan et al. / Journal of Computational and Applied Mathematics 263 (2014) 262–276

(a) The graph of fractal function S1′′ .

(c) The graph of C ′′ = S3′′ .

(b) The graph of fractal function S2′′ .

(d) The graph of fractal function S4′′ .

(e) The graph of fractal function S5′′ .

(f) The graph of fractal function S6′′ . Fig. 3. Second derivatives of the rational cubic spline FIFs in Fig. 2.

Table 3 Uniform errors ∥Φ − Si ∥∞ and ∥Φ ′′ − Si′′ ∥∞ . i

∥Φ − Si ∥∞

∥Φ ′′ − Si′′ ∥∞

1 2 3 4 5 6

0.02689 0.00646 0.01335 0.00844 0.00686 0.02102

712.35455 460.89916 769.16695 769.16986 769.16986 730.95933

P. Viswanathan et al. / Journal of Computational and Applied Mathematics 263 (2014) 262–276

275

offs between them. Further, it is felt that the box counting dimension of the second derivative, which can be estimated using Matlab package ‘boxcount’, can be used as a quantifying parameter to make comparison between various convex rational spline FIFs. It is worthwhile to mention here that shape preserving interpolation schemes that introduce fractality in the derivatives, namely subdivision schemes, are available in classical (non-fractal) numerical analysis as well. However, up to our knowledge, the fractality cannot be controlled in terms of the parameters involved in those schemes. The fractal dimension, which is a quantification of the irregularity, can be measured in terms of the scaling factors involved in the IFS. As |αn | increases the fractal dimension increases [1]. By taking the scaling factors in other subintervals as zero, fractality can be restricted to a small portion of the domain. The reader is invited to refer [18] for a detailed comparative study of the two traditions-recursive subdivision and fractal interpolation. 5. Conclusions In the present article, we have proposed a new kind of C 1 -rational cubic spline FIF that has cubic numerator and quadratic denominator with two families of shape parameters. The shape parameters and the scaling factors are constrained so that the rational cubic FIF preserves the convexity property inherent in a prescribed data set. The parameter identification problem discussed here combines the three ingredients: (i) self-similarity of the interpolating functions and derivatives (ii) convexity (iii) irregularity of the first/second derivative that can be measured via Hausdorff/box counting dimension. Though it may appear at first glance that the globality limits its applications, it is observed through several examples that the changes in the parameters of a particular subinterval do not produce a major impact on the shape of the curve in other subintervals. More precisely, since the classical rational interpolation that emerges as a special case of the proposed FIF is completely local, the proposed scheme is local or global depending on the values of the scaling factors. Besides adding a layer of flexibility to the corresponding classical counterpart, the novelty of the present convexity scheme lies in the fact that it is well suited for the data set generated by a function with irregular derivative arising in connection with nonlinear and nonequilibrium phenomena, for example in the approximation of the motion of a pendulum on a cart from its sampled data. The computational complexity to achieve convexity for the rational FIF seems to be equivalent to that for its classical counterpart. With some mild conditions on the scaling factors, it is shown that the proposed FIF converges to the data defining function Φ as rapidly as the rth power of the mesh norm approaches to zero, provided that Φ (r ) is continuous for r = 1 or 2. Owing to the implicit and recursive nature of a FIF, we imposed nonnegativity constraints on the scaling factors to develop easily calculable conditions for convexity. However, test examples demonstrate the possibility of preserving convexity with negative scalings as well; providing a scope for the extension of the present work. Another question that we leave open for further study is on the ‘‘optimal’’ choice of parameters involved in the IFS that generates convex fractal curves. However, we note that the common practice in traditional non-recursive shape preserving interpolation is to select a preferable solution by minimizing a choice functional subject to the constraints that arise from the shape preservation requirements. Widely used functional in classical shape preservation is the Holladay functional or some approximations of this. We surmise that certain differential evolution optimization algorithm/genetic algorithm can be used to solve this type of constrained optimization problem. Acknowledgments The first author is partially supported by the Council of Scientific and Industrial Research India (Grant No. 09/084 (0531)/2010-EMR-I) for this work. The second author is thankful to the Science and Engineering Research Council, Department of Science and Technology India (Project No. SR/S4/MS: 694/10). References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16]

M.F. Barnsley, Fractal functions and interpolation, Constr. Approx. 2 (1986) 303–329. M.F. Barnsley, Fractals Everywhere, Academic Press, Orlando, Florida, 1988. K.J. Falconer, Fractal Geometry: Mathematical Foundations and Applications, John Wiley, Chichester, 1990. C.E. Zair, E. Tosan, Computer aided geometric design with IFS Techniques, in: M.M. Novak, T.G. Dewey (Eds.), Fractal Frontiers, World Scientific Publishing, 1997, pp. 443–452. H. Xie, Fractals in Rock Mechanics, A. A. Balkema Publishers, Netherlands, 1973. L. Pietronero, Fractals in physics: applications and theoretical developments, Physica A 191 (1992) 85–94. K.S. Birdi, Fractals in Chemistry, Geochemistry, and Biophysics, Plenum, New York, 1995. M. Ali, T.G. Clarkson, Using linear fractal interpolation function to compress video images, Fractals 2 (3) (1994) 417–421. J.L. Véhel, K. Dauodi, E. Lutton, Fractal modelling of speech signals, Fractals 2 (3) (1994) 379–382. P.R. Massopust, Fractal Functions, Fractal Surfaces, and Wavelets, Academic Press, Orlando, 1995. M.F. Barnsley, A.N. Harrington, The calculus of fractal interpolation functions, J. Approx. Theory 57 (1) (1989) 14–34. P. Bouboulis, L. Dalla, Hidden variable vector valued fractal interpolation functions, Fractals 13 (3) (2005) 227–232. A.K.B. Chand, G.P. Kapoor, Generalized cubic spline fractal interpolation functions, SIAM J. Numer. Anal. 44 (2) (2006) 655–676. A.K.B. Chand, G.P. Kapoor, Cubic spline coalescence fractal interpolation function through moments, Fractals 15 (2007) 41–53. A.K.B. Chand, M.A. Navascués, Generalized Hermite fractal interpolation, Rev. R. Acad. Cienc. Zaragoza 64 (2) (2009) 107–120. A.K.B. Chand, Natural cubic spline coalescence hidden variable fractal interpolation surfaces, Fractals 20 (2012) 117–131.

276

P. Viswanathan et al. / Journal of Computational and Applied Mathematics 263 (2014) 262–276

[17] A.K.B. Chand, P. Viswanathan, Cubic Hermite and cubic spline fractal interpolation functions, AIP Conf. Proc. 1479 (2012) 1467–1470. [18] A.K.B. Chand, P. Viswanathan, A constructive approach to cubic Hermite fractal interpolation function and its constrained aspects, BIT 53 (4) (2013) 841–865. [19] F.N. Fritsch, R.E. Carlson, Monotone piecewise cubic interpolation, SIAM J. Numer. Anal. 17 (2) (1980) 238–246. [20] P. Costantini, On monotone and convex spline interpolation, Math. Comp. 46 (1986) 203–214. [21] J.A. Gregory, R. Delbourgo, Shape preserving piecewise rational interpolation, SIAM J. Stat. Comput. 6 (4) (1985) 967–976. [22] M.Z. Hussain, M. Sarfraz, Positivity preserving interpolation of positive data by rational cubics, J. Comput. Appl. Math. 218 (2) (2008) 1–13. [23] M. Tian, Monotonicity-preserving piecewise rational cubic interpolation, Int. J. Math. Anal. 5 (3) (2011) 99–104. [24] R. Delbourgo, J.A. Gregory, Determination of derivative parameters for a monotonic rational quadratic interpolant, IMA J. Numer. Anal. 5 (1) (1985) 397–406. [25] R. Delbourgo, J.A. Gregory, Shape preserving piecewise rational interpolation, SIAM J. Stat. Comput. 6 (5) (1985) 967–976. [26] M.J.D. Powell, Approximation Theory and Methods, Cambridge University Press, 1981. [27] A.D. Millar, R. Vyborny, Some remarks on functions with one sided derivative, Amer. Math. Monthly 93 (6) (1986) 471–475. [28] L.D. Berkovitz, Convexity and Optimization in Rn , John Wiley & Sons, 2002. [29] M. Tian, L.S. Long, Convexity-preserving piecewise rational cubic interpolation, J. Shandong Univ. Nat. Sci. 42 (10) (2007) 1–5. [30] R. Delbourgo, Shape preserving interpolation to convex data by rational functions with quadratic numerator and linear denominator, IMA J. Numer. Anal. 9 (1) (1989) 123–126. [31] L. Fantoni, R. Lazano, Nonlinear Control for Underactuated Mechanical Systems, Springer Verlag, 2001. [32] W.M. Siebert, Circuits, Signals, and Systems, MIT Press, Cambridge, 1986.