Higher order lacunary interpolation

              Rev.Acad.Canar.Cienc., IX (Núm. 1), 15-22 (1997)
HIGHER ORDER LACUNARY INTERPOLATION*
GIOVANNA PITTALUGA - LAURA SACRIPANTE -EZIO VENTURINO
ABSTRACT - A generalization of a previously analyzed lacunary interpolation problem is
considered. We prove the unconditional solvability of the problem, using elements in an
appropriate class of "defi cient" splines. Also, the error analysis presented here sharpens the
estimates obtained in the earlier study. Moreover the convergence rates obtained here are
shown to be optima!.
AMS{MOS) Subject Classifications: 41Al5, 65005
Key words: lnterpolation , spline.
l. Introduction. In this paper we considera lacunary interpolation problem, which is
a generalization of a previous investigation [7], the so-called (0,4) lacunary interpolation
problem and improve the results therein. We assume that information is provided on
a function and its q-th order derivative at a set of equispaced nodes. Our task is the
reconstruction of the function by means of a suitably defined spline. It turns out that
in such instance, the spline is " deficient", in the sense that it is impossible to ensure
continuity of ali derivatives up to the maximum possible order, since there is a constraint
relating the number of the conditions, of interpolatory or continuity type, and the number
of coefficients in each polynomial are of the spline. Solution of lacunary interpolation
problems is usually obtained using special classes of lacunary splines [1,2,3,4,6].
This study can havc application in the solution of boundary value problems governed by
ordinary differential equations of higher order. The fourth order case already dealt with
has a direct practica! relevance, si nce it corresponds to the so-called cantilever problem
[5] and is related to the calculation of the deformations of beams. The method presented
here assumes that the underlying differential equations governing the two point boundary
value problem is solved by means of a finite difference approximation scheme. At the
expense of just one more function evaluation, it is possible to obtain data also on the
highest order derivative, in adclition to the value of the function. These are the data we
assume given for the problem here at hand. The higher order case considered here in our
opinion possesses a mathematical interest, in itself worthy of investigation.
The paper is organized as follows: after giving the necessary definitions and stating
the problem, we reduce its solution to a linear algebraic system. The investigation of
the structure of the matrix is performecl in section 3. In section 4, we analyze a solution
scheme and show that the system is always solvable. Section 5 contains the error analysis.
*Work supported by the Consiglio Nazion ale delle Ri cerche of ltaly a nd by the Ministero dell'Unive rsita
e della Ricerca Scientifica e Tecnologica of ltaly.
IS
© Del documento, de los autores. Digitalización realizada por ULPGC. Biblioteca Universitaria, 2017
We would like to emphasize that here the resolution scheme is directly used as a model
to follow for the derivation of the stability result. In this way the use of the 2- norm as
in [7] is avoided. The direct use of the supremum norm indeed allows the recovering of
a factor of 1/2, lost in the previous analysis. As a consequence, it turns out that we can
recover the optima! order of convergence for the spline, thereby sharpening the results
of the former investigation.
2. The problem. We assume to work on the normalized interval [a, b] = [O, 1],
partitioned by means of the nodes Xk = kh, k = O, ... , n, where h = l/n. We want
to determine the spline function s( x) of degree q + 2 in each subinterval, satisfying the
following conditions:
i) s(x) E 5~~~+2 i.e. its "deficiency" is 2, as emphasized by the upper index;
ii) s(x) E Cq[O, 1] i.e. it is continuous up to order q included.
Let Sk(x) be the restriction of s(x) to the interval 6..k = [xk-1,Xk], k = l, ... ,n.
We assume that the number of subintervals in the partition is larger than the order of the
highest known derivative, thus q < n. In this situation there are n( q + 3) free parameters
to be determined, i.e. the coefficients of the polynomials making up the spline. On the
other hand the number of interpolatory conditions is 2( n + l) .
Explicitly, the latter are
(2. 1)
(i) (i)
sk (xk_¡) = Ík-1
(i) (i)
Sn (xn) = Ín (xn)
i = o, q
i = o, q
k = l , ... , n
where for easeness we use the shorthand notation Jli) = ¡<i>(xk),
We need also to satisfy the following continuity conditions
(2.2) i =O, .. . , q k=l , ... ,n - l.
Their number is then ( q + l) ( n - l ). The free parameters that are still undetermined are
then n(q +3)-2(n + 1)-( q + l)(n- l) = q- l and equal the number of extra (boundary)
conditions that need to be specified to guarantee a unique solution for the problem. For
easeness but without loss of generality, in the sense that specification of conditions at the
other endpoint would result only in marginal modifications of the matrix of the system
and of the relative solution scheme, we choose them to be "initial" conditions, in the
form
(2 .3) _(i)( . ) - ,(i)
S¡ Xo - JO i=l, ... ,q-l.
We introduce now the "unknowns" of the problem, namely the moments of s( x) and of
its derivatives at the breakpoints:
16
© Del documento, de los autores. Digitalización realizada por ULPGC. Biblioteca Universitaria, 2017
M(i) (i)( )
k-1 = S Xk-1
( i) (i) ( ) Nlk-l = sk Xk-1
M~i) = s~\xn)
i =o, ... ,q'
i = q + 1, q + 2 ,
i = O,q.
Then the restriction of the spline to the k-th subinterval can be expressed as follows
( )_ 'q° +2' . M(j) (x - Xk-1) 1
Sk X - L k-1- -,
j=O J.
k = 1, ... ,n
Let us recall that
M(i) -f(i)
k-1 - k-1
M (i) _¡·(i)
o - o
i = O,q k = 2, ... , n + 1 ,
i =O, ... , q
Differentiating,
q+2-i ·
(i)(·)- "M(i+j)(X-Xk-1)1
sk x - L k-1 ·1
j=O J.
i =O, ... ,q.
If we now impose the interpolation conditions that are not yet implicitly satisfied, i.e.
(2.1 ), and the continuity conditions (2.2) we obtain a square linear algebraic system of
( n - 1 )( q + 1) + 2 equations in as many unknowns. Explicitly, the latter are
mT = [M<,+2 ) M(,+l) M<,+2 ) M(v+i) M(v-i) M(i)
O ,o ,1 ,1 ,1 ,···,1,···
Notice that they are written in reverse order, with respect to the derivative order. This
will enhance the resulting structure of the matrix for the investigation of the next section.
Let us write explicitly the i-th equation of the system. We need to distinguish three
different cases:
for i = q it is
2 . '°' M<v+i) !:!._ - M(v) - M(v)
¿ •-1 -1 - k k-1
j=l J.
k = 1, ... ,n,
while for i = q-1(-1)1 it is
~ M(q+i) hq+j-i - M(i) = - ~ M(i+i) hi
¿ o ( · -)1 1 ¿ o -1 ,
j=I q + J - i . j=O J.
17
© Del documento, de los autores. Digitalización realizada por ULPGC. Biblioteca Universitaria, 2017
k = 2, ... , n - l ,
j=O
j#q-i
and finally for i = O it is given by
..¿--. M(q+j) ~ = M(o) - ~ M(j) hi
0 O ( + .), 1 0 O .! ,
j=I q J . j=O J
q+2 hi hq
~ M(j) - - M(O) - M(O) - M(q) - 0 k-I ·¡ - k k-I k-I t
j=I J. q.
k = 2, ... ,n .
j#q
3. Matrix structure. Let the system be written in compact form as Am = b. Here
the matrix A = [ A; ,j J, i, j = 1, ... , q + l denotes a block Hessemberg matrix, w here
however each block has dimensions that may differ from those of other blocks.
More precisely,
A(n,n)
',) i=l,q+l, J = 1,2,
A(n-I ,n)
•,J i = 2, ... , q, J = 1,2,
A(n-I,n-I)
' ,) i = 2, . . . ,q, J = 3, ... , q + l ,
A(n,n-I)
1,J i=q+l J = 3, ... , q + l
Moreover, !et us observe that A;,i = O for j > i + 1, i = 1, ... , q - l , while ali the other
blocks are either diagonal, subdiagonal or bidiagonal as explicitly stated below
hi-j+2
A;i =----!
(i-j+2)! n
hi-j+2
A;i = -----E
(i-j+2)!
and for the cases j = 3, ... , q + l
i=l,q+l
i = 2, ... ,q
i= j-1
i=J, ... ,q
i = q + l
where, by denoting by b;,j the Kronecker symbol, we write
E =( Ói,j )n- 1,n ,
F = (bi - 1,j)n- I,n - l ,
D =( Ói-1,j )n,n-I ·
IX
j = 1,2 ,
j = 1,2 ,
© Del documento, de los autores. Digitalización realizada por ULPGC. Biblioteca Universitaria, 2017
4. The algoritm. It is possible to easily solve the system by "almost" direct forward
substitution. Let us assume to be at the k-th stage of the procedure, k = l , ... , n - l.
Then we proceed as follows:
1) We solve the 2 by 2 system obtained from the k-th equations of the füst and last
horizontal blocks, getting the unknowns Mk"...~2 ), Mk"...~1). Then the k-th columns in the
first two vertical blocks can be eliminated;
2) from the k-th equations of the j-th horizontal blocks, j = 2, ... , q, we can
immediately determine the unknowns M?J,i = q-1(-1)1, since Mii)' s = 1, . .. ,k-
1, i = q -1(-1)1 are already known from the previous steps;
3) at the end of the ( n - l )-th step of the procedure, it is sufficient to sol ve the
last equations of the first and last horizontal blocks, thus determining the remaining
unknowns M(q+2 ) M(q+l).
n-1 , n-1
It is fundamental to not'e that in the implementation of the above procedure, the structure
of the matrix is left unchanged. In particular after the deletion at each stage of the ( q + l)
columns relative to the unknowns determined at that step we find an analogous matrix
to the one present before the elimination took place, but of smaller dimensions.
The 2 by 2 matrix B of the linear system that needs to be sol ved in step 1) of each stage
of the procedure to obtain the unknowns M?+2 ) M?+ 1), j = O, ... , n - l, and finally in
step 3) has the form
Easily, detB =/- O, so that in view of this and the fact that the procedure never breaks
clown, we have the following important result.
THEOREM l. The linear algebra.ic system for the calculation of the moments of the
spline function is always nonsingular.
5. Convergence results and stability estimates. Let us assume that f(x) E C(q+J)
over [O, 1] and Jet us denote by T( x) its truncated Taylor expansion, i.e. the Taylor
polynomial of degree q + 2
q+2 .
~ ( ·) (x - Xk- 1)1
T(x) = L., f 1 (xk-il 1
J=Ü J.
Then
f(x) - T(x) = Khq+J .
Here and in what follows, J( and J( k ,i denote suitable constants. Let us define the errors
we need to estimate:
p =o, ... ,q + 2,
I ')
© Del documento, de los autores. Digitalización realizada por ULPGC. Biblioteca Universitaria, 2017
(5.1) e(p) - ¡<P) - M(p) k - k k p= 1, ... ,q-l,q+l,q+2 .
If x E É>k, denoting by ek(x) the restriction of e(x) to É>k, we have
(5.2)
Since
ek(x) = f(x) - T(x) + T(x)- sk(x).
q+2 )'
" (j) (x - Xk-1 J T(x) - Sk(x) = L ek-l . , ,
j=l J.
#q
by imposing the continuity at the point Xk , we obtain
k = l ,
(5.3)
k = 2, ... , n .
Pr_oceeding in a similar fashion , by imposing the continuity of the i-th derivative,
i = 1, ... , q - l at the breakpoints, we get
'q°+'2 . hi-i . e(J) ___ - e( •) =K hq+3-i
L o (j _ i)! 1 1 ,,
j=q+l
k = l '
(5.4) q+2 . . '°' e(j) ~ - e(i) =K hq+3-i L k-1 ( · _ ·)' k k ,,
j=i J i .
k = 2, .. . , n - l .
#q
But, in view of the interpolation conditions,
k = l, .. . ,n,
so that
(5.5)
q+2 .
" (j) hJ-q , 3
L ek-(. 1- -)' = !'l.k ,qh
j=q+l J q .
The "error equation" Ae =€is given by the system whose equations are (5.5), (5.4), (5.3) .
Here, e represents the error vector, A is the same matrix discussed in section 3, € denotes
the consistency error, in this case the interpolation error. The latter has components of
20
© Del documento, de los autores. Digitalización realizada por ULPGC. Biblioteca Universitaria, 2017
order €~i) = O(hq+J-i) in the equation corresponding to the i-th derivative, where the
notation for the consistency error mimics the one used for the discretization error ( 5.1).
Once again, we denote the unknown vector using a reverse order:
eT = [ (q+2) (q+l) (q+2) (q+l) (q-1) (1)
e0 , e 0 , e 1 , e1 , e1 , ... , e 1 , ...
(q+2) (q+l) (q-1) (1) l
· · · , en-1 , en-1 , en-1 , · ·' , en-1
so that the matrix of the system is the same as the one for the moments, analyzed in
the previous section. To determine the errors e~), the same earlier considerations hold.
Step 1) of the algorithm applied to the error equation leads then to the solution of the
following system, initially in the unknowns e~q+2> , e~q+l)
(5.6)
Easily, then we have the following important stability result
(5.7)
Direct solution however allow.s us to obtain a sharper estímate for e~q+l) than the one
that can be obtained by a crude application of (5.7) to (5.6):
(5.8) e~g+2l = O( h)
Iteration of the calculations as prescribed by the algorithm leads to the convergence
result. Let us denote by R the matrix corresponding to all the elementary operations
described by the algorithm. The solution of the error equation is equivalent then to the
solution of the equation e = R€, with I = RA, the identity matrix.
In view of the above remarks, it is indeed suffi.cient to note that A-1 = R, and that
the algorithm when considering the equations related to the i-th derivative, consists in
forward substitution, i.e. in moving to the right hand side at most q terms of the same
order as the one of the right hand side, i.e. O(hq+3-i). The previous statement is easily
established by induction. Since q < n, no loss of accuracy ensues for the unknown being
determined at this stage, the right hand side remaining of order O(hq+J-i). Finally, for
step 3) of the algorithm, the analysis done earlier for the system (5.6) applies, the only
change being in the name of the unknowns. It thus leads to the same conclusion, i.e.
estímate (5.8) for the unknowns e\2l, e\1>.
The convergence of the algorithm is thus ensured, and the rate equals the one of the
consistency error. In summary
THEOREM 2. For the error in the calculation of the moments, on top of the above
estimations for e~q+2> and e~q+l) we have
i=q-1(-1)1,
21
© Del documento, de los autores. Digitalización realizada por ULPGC. Biblioteca Universitaria, 2017
and, iterating the procedure,
k=l, .. . ,n- 1 i = 1, ... , q - 1, q + 1, q + 2 .
As a consequence, since the problem is linear and we use a consistent method, in view of
estimate (5.7) we have also the stability result
THEOREM 3. The norm of the inverse matrix of the system satisfies the following
estimate
On using the triangular inequality on the error representation (5.2), we have also the
error for the spline at arbitrary points in [0,1].
THEOREM _4. For the error of the spline function and its derivat.ives the following
estimate holds
x E [O, 1] p =o, ... ,q + 2.
References
[1] T.Fawzy and L.L.Schumaker, A Piecewise Polynomial Lacunary Interpolation
Methocl, J. Approx. Thcory 48 ( 1986) 407-426.
[2] A.Meier and A.Sharma., La.cuna.ry interpolation by splines, SIAM J.Numer. Anal. 10
(1973) 433-442.
[3] A.Sa.xena., Interpolation by a.lmost c¡ua.rtic splines, A cta Math.Hungar. 51 (1988)
283-292.
[4] A.Sa.xena, Interpolation by c¡ua.rtic splines, Ganita 38 (1989) 76-90.
[5] G.Stra.ng, " Introcluction to appliecl mathema.tics", vVellesley-Cambriclge Press, 1986.
[6] A.K.Va.rma, La.cunary interpola.tion by splines II, Acta Math.Hungar. 31 (1978) 193-
203.
[7] E.Venturino, On the (0,4) la.cunary interpolation problem a.nd sorne related c¡uestions,
(submittecl for publica.tion).
Giovanna. Pitta.luga - Laura Sacripante, Dipa.rtimento di Matematica., Universita di
Torino, Via Ca.rlo Alberto 10, 10123 Torino, Ita.ly.
Ezio Venturino, Dipartimento di Ma.tema.tica, Cittá Universitaria., Viale A. Doria 6, 95125
Ca.ta.nia, Ita.ly.
22
© Del documento, de los autores. Digitalización realizada por ULPGC. Biblioteca Universitaria, 2017