Rational Approximation on Exponential Meshes

: Error estimates of pointwise approximation, that are not possible by polynomials, are obtained by simple rational operators based on exponential-type meshes, improving previous results. Rational curves deduced from such operators are analyzed by Discrete Fourier Transform and a CAGD modeling technique for Shepard-type curves by truncated DFT and the PIA algorithm is developed.


Introduction
In [1], simple rational operators r n and q n were constructed, reaching the following pointwise approximation error estimates ∀x ∈ [0, 1] and ∀ f ∈ C([0, 1]) and Here, C denotes a positive constant independent of n and ω( f ) is the usual modulus of continuity of f .
Moreover, in [1] the authors proved that the exponents α and β in (1) and (2), respectively, cannot be equal to 1 (see remark to Corollary 2.5 proposed by Totik [1]). A crucial key for such estimates is the careful choice of node meshes of algebraic-type.
Zeros of orthogonal polynomials as nodes for rational approximation were examined in [3], but corresponding results are weaker than (1) and (2).
It was an open problem to construct simple rational operators, reaching pointwise approximation error estimates in terms of functions vanishing at any fixed point in [0, 1] faster than in (1) and (2), like or On the other hand, it was an unsolved question to get results of convergence and pointwise estimates of the approximation error by the above rational operators, based on exponential-type meshes. In such cases, the main difficulty is to balance the rational nature of above operators with the exponential mesh behavior.
The present paper aims at giving a first successful answer to the above problems, considering easy rational operators based on exponential-type meshes.
In Section 2, convergence statements and pointwise approximation error estimates, involving functions of type (3) or (4), are given in  Then, in Section 3, we use above operators to construct rational curves by Discrete Fourier Transform (DFT) and design a truncated progressive iterative approximation technique, useful in CAGD modeling.
Recently, parametric Shepard-type curves were studied and properties and algorithms, interesting in CAGD and image processing, were deduced in [4,5]. In Section 3, we study how to transform Shepard-type curves into a different domain using Discrete Fourier Transform through the well-known Fast Fourier Transform algorithm. Theorem 6 presents the inner structure of Shepard-type curves in the form of other Shepard-type curves, known as base Shepard-type curves in the transform form. These base Shepard-type curves have as control points twiddle factors of DFT; hence, the geometry of these curves is determined by the geometry of the corresponding star polygons. The modeling power of such basic Shepard-type curves overcomes the analogous Bézier and B-spline ones, as in original Shepard-type curves case (cf. [4]). Then, in Section 3.1, by truncated DFT and the progressive iterative approximation (PIA in short) algorithm introduced in [4], we deduce a method to construct new Shepard-type curves having good shaping behavior. Indeed, the number of base Fourier functions and the iteration level of PIA algorithm handle modeling the outline of Shepard-type curve, getting as a limit case the Shepard-type interpolating curve.
The proofs of the above results are included in Section 4. They are based on smart choices of nodes meshes, careful estimates for our operators and suitable matrix formulations of DFT and the PIA algorithm. Finally, Section 5 shows numerical experiments, illustrating achieved results.

Rational Approximation on Exponential-Type Meshes
For n ∈ N, we introduce the nodes matrix X = (x n,k = x k , k = 0, . . . , n) ⊆ [a, b], a, b ∈ R. Then, for any function f ∈ C([a, b]), we consider the Shepard operator S n as with x ∈ [a, b] and s > 2 even. From the definition, we can deduce that S n is a positive, linear operator, preserving constants; S n ( f ) is a rational function having degree (sn, sn), interpolating f at x k , k = 0, . . . , n and is stable in Fejér sense, i.e., min If the mesh is equispaced on [a, b], and we assume [a, b] = [0, 1] for example, then S n is a symmetric operator, i.e., This operator is widely used in classical approximation theory and problems of scattered data interpolation (see, e.g., [6][7][8][9]).
When the nodes' mesh is equispaced, direct and converse results for the operator S n exist (e.g., [7,10]). For nodes mesh of algebraic type, pointwise estimates were obtained in [1].
The case of nodes exponentially thicker near one fixed point was an open problem. Now we prove that, by nodes exponentially dense, S n operator achieves results of convergence and estimates of pointwise approximation error, involving functions of type (3) and (4).
To this aim, we assume [a, b] = [0, 1] and denote by || f || the usual supremum norm on [0, 1] of f ∈ C([0, 1]). Moreover, C denotes a positive constant that can assume different values, even when it appears more times in the same formula. Then let be the nodes' matrix. This mesh is exponentially thick near 0. Then, we consider operator S n (X) defined by (5) based on the nodes matrix X in (7). We have 1]) and S n (X) be the operator defined by (5)- (7). Then where Remark 1. Equation (8) allows us to deduce the uniform convergence of S n (X; f ) to f , as n → ∞, ∀ f ∈ Lip 1 ([0, 1]). Moreover, estimate (9) gives an answer to the problems presented in Introduction. The thickness of the exponential-type mesh near 0 affects estimate of the error (due to the presence of functions χ 1 and χ 1 at the r.h.s. in (9)).
As remarked in [1], the left-end point 0 is not a special point and we can get analogous results for the right-endpoint 1 case, both endpoints case and any interior point case. For example, let with x k given in (7). Now consider the matrix Y = (y k , k = 0, . . . , n, n ∈ N). This mesh is exponentially finer near 1. Moreover, for any f ∈ C([0, 1]), consider the operator S n (Y) based on the matrix Y. We can state 1]) and S n (Y) be the operator defined by (5) and (10). Then with
Combining the above results, we can consider the matrix Z = (z k , k = 0, . . . , n, n ∈ N) , with n even, x n/2,k as in (7). Then, denote by S n (Z) the operator S n based on the matrix Z. When n is odd, we can replace S n (Z) by S n+1 (Z). We have 1]) and S n (Z) be the operator defined by (5) and (13). Then In

Remark 3.
From (14), we have the uniform convergence of S n (Z; Moreover, the mesh thickness near 0 and 1 affects estimate of the error (due to functions χ 3 (x) and χ 3 at the r.h.s. of Equation (15)). Hence, Theorem 3 successfully answers the problems posed in the introduction. For example, Additionally for n even, let This mesh is thicker near the inner point 1/2. Denote by S n (T) the operator S n based on the matrix T = (t k , k = 0, . . . , n, n ∈ N) ⊆ [0, 1]. If n is odd, we can replace S n (T) by S n+1 (T).
We prove 1]) and S n (T) be the operator defined by (5) and (16). Then

Shepard-Type Curves and Discrete Fourier Transform
Let us first recall some properties of parametric curves of Shepard-type (see [4]). Let A n (t) = [A n,0 (t), A n,1 (t), . . . , A n,n (t)], where where C is any fixed positive constant, and s > 2 is even.
Let P = [P 0 , P 1 , . . . , P n ] T , P i ∈ R d , i = 0, 1, . . . , n, d ≥ 2, be the control vector. Then the parametric curve of Shepard-type S n [P, t] is defined as (20) In [4], some properties of such curves were studied, that are interesting in CAGD. In particular, S n [P] is a rational curve that reproduces points and lies within the convex hull of the control polygon P. Such a curve satisfies the pseudo-local control property; that is, each function A n,j (t), 0 ≤ j ≤ n, reaches its maximum value close to 1 at t = t j ; therefore the point P j strongly affects the shape of the curve close to t = t j . From (6), Shepard-type curves are symmetric, i.e., S n [P, t] = S n [P, 1 − t], withP = [P n , P n−1 , . . . , P 0 ] T . The choice 0 < n s λ ≤ C makes S n [P] a near-interpolating curve of the original control polygon. This fixes the flat spots artifact that affects the original Shepard operator (cf. [4]).
In [4], we introduced and investigated a PIA technique for curves of Shepard-type. The PIA process starts with an initial Shepard-type curve; then through iterations, it adjusts the control points at each iteration, resulting in a sequence of fitting curves. The limiting of the curves at different iterations is the Shepard-type curve interpolating the data. In details, given the control vector P and the basis A n,i (t), i = 0, . . . , n, defined by (19), the initial curve is generated as with P 0 i = P i , i = 0, . . . , n. Then the remaining curves of the sequence γ k+1 (t), for k ≥ 0, are computed as and ∆ k i the adjusting vectors given by The iterative process can be also expressed in matrix form as We remark that B is a symmetric, stochastic, diagonally dominant matrix (see [4]). We say that the γ 0 curve satisfies the PIA property iff lim k γ k (t i ) = P i , i = 0, . . . , n. We can state (cf. Thanks to the PIA property, we can construct a sequence of control polygons that converge to the control polygon generating a Shepard-type interpolating curve. In addition, we can use k as a shape parameter so to model several outlines; the extreme cases are the curve of Shepard-type and the global interpolating curve of Shepard-type. As remarked in [4], the rate of convergence of such a procedure is faster than the Bézier case. Then we recall some preliminaries on Fourier analysis, namely the resulting Fourier Transform. The Fourier Transform is the decomposition of a function into components at higher and higher frequencies. Discrete Fourier Transform (DFT) is the discrete counterpart of the Fourier Transform, yielding an estimate starting from a finite sample. DFT maps an ordered set of n + 1 complex number to a different one. Let x(k) ≡ x 0 , x 1 , x 2 , . . . , x n be a series of n + 1 complex numbers. We assume a periodicity condition that outside the range 0, n the series is extended n + 1− periodic, i.e., x L = x L+n+1 , ∀L. We denote DFT of this series x(L). The forward transform is defined as x(i)w Li , L = 0, 1, . . . , n, with w = exp(−j 2π n+1 ) known as twiddle factor and j the solution of j 2 = −1 (complex unity). In practical applications, DFT is computed by the well-known Fast Fourier Transform algorithm. Now we are ready to analyze Shepard-type curves by DFT.
Remark 6. From Theorem 6, S n [P, t] curve can be expressed as weighted average of base Shepard-type curves, by DFT of the original control points. It can be easily proved thatÃ behaves as a basis for Shepard-type curve by polygon points in transformed form. Thus, the original Shepard-type curves are generated by these base Shepard-type curves.
The control points that are obtained for them are w −L , w −2L , w −3L , . . . , w −(n+1)L , so the corresponding polygon is a regular or star one. The relative base curves of Shepard-type lie in the convex-hull of these regular polygons or star polygons, (cf. [11]).
Analogous Fourier analysis for Bézier, B-splines and rational B-splines curves was presented in [12,13], but at lower modeling power than our curves (see Section 5).

Shepard-Type Curves via DFT and PIA
Consider the global interpolating Shepard-type curve (cf. [4]) with Q k ∈ R d , d ≥ 2, such that Now let us consider only the first k, k ≤ n, Fourier basis functions in W, or in matrix notation W (k) ∈ R k×n+1 , with index (k) denoting the truncation procedure. This choice is made in analogy with truncation occurring in some statistical contexts involving FT. Then the truncated interpolating curve of Shepard-type is introduced: Obviously, if k = n, W (n) H W (n) = I and we get back (25). Hence, varying 0 ≤ k ≤ n in (27), we get different curves approaching the interpolating Shepard-type curve.
Representation (27) suggests to use the above PIA technique (see Theorem 5,(22) and consequent remark) to construct a method giving a pencil of curves of Shepard-type, modeling the original data points.
We sketch the procedure. In (27), we replace V by with r being the iteration level. From [4], lim r V r = V and convergence rate is fast, so a few iterations are enough to go close to V. Then playing on two shape parameters, the number k of basis Fourier functions in (27) and the iteration level of the PIA algorithm-the index r in (28), the designer can get intermediate contours not reachable by original PIA format (see Section 5).
Proof of Theorem 2. Analogously to Theorem 1, we have Then, as in the proof of Theorem 1, Hence, Then we assume 1 − x = o(1). By |x − y j | ≤ |x − y k |, k = j + 2, . . . , n and f ∈ Lip 1 ([0, 1]), Other cases can be handled as in Theorem 1 and we deduce the statement.
Proof of Theorem 3. Proceeding as in Theorems 1 and 2, we obtain Working as before (cf. [1]), the assertion is proved.
Proof of Theorem 4. Working as in the proof of Theorems 1 and 2, we get Following the above demonstrations (cf. [1]), we derive the statement.
Proof of Theorem 6. Let S n [P, t] be the Shepard-type curve defined by (20). The control points P i , i = 0, 1, . . . , n, can be considered as points on complex plane. Therefore, each P i , i = 0, 1, . . . , n, can be equivalently seen as the inverse transform of points P i , with P i being the Fourier Transform of points P i , i.e., with w = exp −j 2π n+1 the twiddle factor. Thus, equation (20) can be written as As P L 's are independent of i, we get and by (23) the assertion follows.

Numerical Experiments
In this section, we show some examples of Shepard-type curves in the transformed domain.
From such figures, the modeling outperformance of basic Shepard-type curves with respect to Bézier case is clear (cf. [12]).

Real part
Imaginary part

Example 3
Let n = 6, s = 4 and λ = 5 × 10 −5 . Figure 4 shows the corresponding basic Shepard-type curves satisfying properties mentioned in Section 3. Comparing Figure 4 with Figure 1 in [13] for the corresponding base B-spline curves, we can note that basic Shepard-type curves show superior shaping capabilities.
We extract from the spiral a sample of 11 (n = 10) control points as (x(s i ), y(s i )), s i = i n , i = 0, . . . , n.
They form the control polygon in Figure 5. We start from these control points and fit the spiral by a sequence of three curves that are defined by the truncated iterative procedure introduced in Section 3.1 with k = 9, 10, 11, respectively, in (27), r = 1, s = 4 and λ = 10 −5 (see Figure 5, left). From Figure 5 (left), we can see the modeling behavior of our procedure; indeed for k = 11 and r = 1, we find back the curve obtained by one iteration of the PIA algorithm, while the choices k =9, 10 and r = 1 give intermediate curves approximating the given data points. Analogous curves for the choice r = 2, s = 4 and λ = 10 −5 are shown in Figure 5 (right). In this way, the designer has different choices for outlines modeling the spiral.
Then we start from these control points and fit the spiral by a sequence of three curves defined by our truncated iterative procedure with s = 4, λ = 3 × 10 −6 , k = 17, 18, 19 in (27) and r = 1 (Figure 6, left) and r = 3 ( Figure 6, right), respectively. Comparison of Figure 6 (left) with Figure 6 (right) shows the influence of shape parameters in slight changes of our curves, giving new contours fitting the spiral.

Example 6
Consider the Golden spiral curve in a sequence of 15 (n = 14) points defined by (x(t i ), y(t i )) = aϕ 2τ i /π cos τ i , aϕ 2τ i /π sin τ i , ϕ = 1 2 , τ i = −5 + 10 i n , i = 0, . . . , n, forming the control polygon in Figure 7. Actually, since the curve lacks periodicity, we virtually extend the control points at the right and/or left with mirror points. As is well known, the enlarged function is now even, then the Fourier Transform reverts to the Cosine Fourier Transform [14], with the favorable consequence that the coefficients and the inverse transform keep real. From the computational point of view this is obtained by simply computing the Discrete Cosine Fourier Transform of the original 15 points. We fit these points by 6 sequences of curves defined by the above procedure with s = 4, λ = 10 −5 , k = 10, 11, . . . , 15 and r = 1 (see Figure 7), from which one can see several contours approximating Golden Spiral. Funding: This research received no external funding.