An extension of the Floater–Hormann family of barycentric rational ...

Report 5 Downloads 36 Views
An extension of the Floater–Hormann family of barycentric rational interpolants Georges Klein University of Oxford

Fox Prize Meeting, Edinburgh, June 24, 2013

Introduction Extended FH Applications

Outline

1

Introduction: linear barycentric rational interpolation

2

Extended Floater–Hormann interpolation

3

Applications of extended Floater–Hormann interpolation

G. Klein

Extended Floater–Hormann interpolation

2/42

Introduction Extended FH Applications

Floater–Hormann interpolation Convergence and stability Differentiation

One-dimensional interpolation

Given: a ≤ x0 < x1 < . . . < xn ≤ b, f (x0 ), f (x1 ), . . . , f (xn ),

n + 1 distinct nodes and corresponding values.

We study functions g from a finite-dimensional linear subspace of (C [a, b], k · k∞ ) which interpolate f between the nodes, g (xi ) = f (xi ) = fi ,

i = 0, . . . , n.

Main interest: well-conditioned equispaced interpolation.

G. Klein

Extended Floater–Hormann interpolation

3/42

Introduction Extended FH Applications

Floater–Hormann interpolation Convergence and stability Differentiation

The polynomial in barycentric form There exists a unique polynomial pn of degree ≤ n, which interpolates the data and can be written in Lagrange form pn (x) =

n X

`i (x)fi

with

`i (x) =

i=0

Y x − xk , xi − xk

k6=i

in the first barycentric form (backward stable (Higham ’04)) pn (x) = L(x)

n X i=0

λi fi , x − xi

.Q Qn with λi = 1 k6=i (xi − xk ), and L(x) = k=0 (x − xk ), and in the second barycentric form (forward stable (Higham ’04)) Pn λi i=0 x−xi fi . pn (x) = Pn λ i

i=0 x−xi

G. Klein

Extended Floater–Hormann interpolation

4/42

Introduction Extended FH Applications

Floater–Hormann interpolation Convergence and stability Differentiation

From polynomial to linear rational interpolation Lemma Let {xi }ni=0 be a set of n + 1 distinct nodes, {fi }ni=0 corresponding real numbers and let {vi }ni=0 be any non-zero real numbers. Then (a) the barycentric rational function n X

R(x) =

vi fi x − xi

i=0 n X i=0

, vi x − xi

interpolates fk at xk : limx→xk R(x) = fk ; (b) conversely, every rational interpolant of the fi may be written in barycentric form for some weights vi . G. Klein

Extended Floater–Hormann interpolation

5/42

Introduction Extended FH Applications

Floater–Hormann interpolation Convergence and stability Differentiation

Construction of Floater–Hormann interpolation for arbitrary nodes - Given n + 1 nodes, x0 , . . . , xn , and function values, f0 , . . . , fn , choose an integer d ∈ {0, 1, . . . , n}, “blending parameter”, - for i = 0, . . . , n − d, define pi (x), the polynomial of low degree ≤ d interpolating fi , fi+1 , . . . , fi+d . The d-th interpolant of the family is a “blend” of the pi (x), n−d X

rn (x) =

λi (x)pi (x)

i=0 n−d X

,

with

λi (x) =

(−1)i . (x − xi ) . . . (x − xi+d )

λi (x)

i=0

Notice that for d = n, rn simplifies to pn . G. Klein

Extended Floater–Hormann interpolation

6/42

Introduction Extended FH Applications

Floater–Hormann interpolation Convergence and stability Differentiation

Linear barycentric rational form For its evaluation, we write rn in linear barycentric form n−d X

rn (x) =

n X

λi (x)pi (x)

i=0 n−d X

= λi (x)

i=0 n X i=0

i=0

wi fi x − xi . wi x − xi

For equispaced nodes, the weights wi do not depend on f , and oscillate in sign with absolute values 1, 1, . . . , 1, 1, 1 4, 1 4 8, 8,

1 2 , 1, 1, . . . , 1, 1, 3 4 , 1, 1, . . . , 1, 1, 7 8 , 1, 1, . . . , 1, 1,

··· G. Klein

d = 0, 1 2, 3 4, 7 8,

d = 1, 1 4, 4 1 8, 8,

d = 2, d = 3, d ≥ 4.

Extended Floater–Hormann interpolation

7/42

Introduction Extended FH Applications

Floater–Hormann interpolation Convergence and stability Differentiation

Algebraic convergence and polynomial reproduction

Theorem (Floater–Hormann ’07) Let 1 ≤ d ≤ n and f ∈ C d+2 [a, b], h =

max (xi+1 − xi ), then

0≤i≤n−1

kf − rn k∞ ≤ Khd+1 , where K depends only on d, b − a and derivatives of f ; the analytic rational function rn has no real poles; rn reproduces polynomials of degree at least d.

G. Klein

Extended Floater–Hormann interpolation

8/42

Introduction Extended FH Applications

Floater–Hormann interpolation Convergence and stability Differentiation

Lebesgue function and constant The Lebesgue constant associated with linear barycentric interpolation,  X n X n wi |wi | , Λn,d = max Λn,d (x) = max a≤x≤b a≤x≤b |x − xi | x − xi i=0

i=0

is the condition number of the interpolation scheme. 4

15

3

10

2

5

1 −1

0

1

−1

0

1

Figure: Lebesgue function for Floater–Hormann interpolation with equispaced nodes in [−1, 1] with d = 2 and d = 5 and n = 40. G. Klein

Extended Floater–Hormann interpolation

9/42

Introduction Extended FH Applications

Floater–Hormann interpolation Convergence and stability Differentiation

Lebesgue function with n = 40 and d = 0, . . . , 40

(loading movie...)

G. Klein

Extended Floater–Hormann interpolation

10/42

Introduction Extended FH Applications

Floater–Hormann interpolation Convergence and stability Differentiation

Lebesgue constant Theorem (Bos–De Marchi–Hormann–K. ’12) Let 0 ≤ d ≤ n and the nodes xi , i = 0, . . . , n, be equispaced. Then n  2d−2 log − 1 ≤ Λn,d ≤ 2d−1 (2 + log n). d +1 d 108

8 d=3

7

107

6

106

5

105

d=2

n = 100

104

4 d=1

3

103 102

2 1

101

0 0

100 0

20

40

60

80

100 120 140 160 180 200

5

10

15

20

25

Figure: Logarithmic growth with n (left) and exponential with d (right). G. Klein

Extended Floater–Hormann interpolation

11/42

Introduction Extended FH Applications

Floater–Hormann interpolation Convergence and stability Differentiation

Central question: How to choose d? Answer: It depends on the nodes and on the regularity of f . FH ’07: large d for very smooth functions, small d for less smooth functions; Numerical Recipes ’07: “It is wise to start with small values of d before trying larger values”; BMHK ’12: not too large, the Lebesgue constant grows exponentially with d; Cross validation (no extra knowledge on the data): remove a few well chosen data values, compute the interpolant for various values of d, compare the errors at the removed data points and keep the d that gives the smallest error (funqui); EFH: wider range of good choices for d. G. Klein

Extended Floater–Hormann interpolation

12/42

Introduction Extended FH Applications

Floater–Hormann interpolation Convergence and stability Differentiation

Interpolation with FH and d = 15

(loading movie...)

G. Klein

Extended Floater–Hormann interpolation

13/42

Introduction Extended FH Applications

Floater–Hormann interpolation Convergence and stability Differentiation

Interpolation of perturbed data with FH and d = 15

(loading movie...)

G. Klein

Extended Floater–Hormann interpolation

14/42

Introduction Extended FH Applications

Floater–Hormann interpolation Convergence and stability Differentiation

Interpolation of f (x) = sin(20x) perturbed data 2

10

exact data FH funqui spline cheb

1

10

FH funqui spline cheb

0

10

−5

10

0

10

−10

10

−1

10

50

100

150

200

50

100

150

200

Figure: Errors with 8 ≤ n ≤ 200 using FH with d = 15, funqui, B-spline of degree 15, and polynomial interpolation in Chebyshev points.

G. Klein

Extended Floater–Hormann interpolation

15/42

Introduction Extended FH Applications

Floater–Hormann interpolation Convergence and stability Differentiation

Differentiation of barycentric functions Proposition (Schneider–Werner ’86) Let rn be a rational function given in its barycentric form with non vanishing weights. Assume that x is not a pole of rn . Then for k ≥1 n X (k)

rn (x) = k!

j=0

wj rn [(x)k , xj ] x − xj n X j=0

(k) rn (xi )

= −k!

n X

wj x − xj

,

x not a node,

! . wi , wj rn [(xi )k , xj ]

i = 0, . . . , n.

j=0 j6=i

G. Klein

Extended Floater–Hormann interpolation

16/42

Introduction Extended FH Applications

Floater–Hormann interpolation Convergence and stability Differentiation

Differentiation matrices Define the matrix D(1) (Baltensperger–Berrut–No¨el ’99):  wj 1   , i 6= j,    wi xi − xj n (1) X Dij = (1)  − Di` , i = j,     `=0 `6=i

and D(k) , k > 1, with the “hybrid formula” (Tee ’06),   k  wj (k−1) (k−1)   D − D , i= 6 j,  x − x w ii ij   i j i n (k) X Dij = (k)  − Di` , i = j.     `=0 `6=i

G. Klein

Extended Floater–Hormann interpolation

17/42

Introduction Extended FH Applications

Floater–Hormann interpolation Convergence and stability Differentiation

Evaluation of derivatives of rn at the nodes

With f = (f0 , . . . , fn )| , the product D(k) · f returns

| (k) (k) rn (x0 ), . . . , rn (xn ) ,

the vector of the k-th derivative of rn at the nodes.

G. Klein

Extended Floater–Hormann interpolation

18/42

Introduction Extended FH Applications

Construction Convergence and stability Examples

Extended Floater–Hormann interpolation

G. Klein

Extended Floater–Hormann interpolation

19/42

Introduction Extended FH Applications

Construction Convergence and stability Examples

Extended Floater–Hormann interpolation For given d, Λn (x) typically has at most d large local maxima at the ends the interval if the nodes are equispaced. 15 10 5 −1

0

1

Idea: construct 2d new data values e f−d , . . . , e f−1 and e fn+1 , . . . , e fn+d out of the given data f0 , . . . , fn . We choose positive integers e n  n and de ≤ e n, and compute (k) (k) e ren (x0 ) and ren (xn ), k = 1, . . . , d. G. Klein

Extended Floater–Hormann interpolation

20/42

Introduction Extended FH Applications

Construction Convergence and stability Examples

Extended Floater–Hormann interpolation We set  e d  X  (xi − x0 )k (k)   f + r (x ) , −d ≤ i ≤ −1,  0 0 e n  k!   k=1  e 0 ≤ i ≤ n, fi := fi ,    e d  X  (xi − xn )k  (k)  , n + 1 ≤ i ≤ n + d. f + r (x )  n n e  n k! k=1

The extension of the Floater–Hormann family then is n+d X

e rn (x) :=

i=−d n+d X

wi e fi x − xi

i=−d G. Klein

. wi x − xi

Extended Floater–Hormann interpolation

21/42

Introduction Extended FH Applications

Construction Convergence and stability Examples

Convergence of extended Floater–Hormann interpolation Notation: e D := min{d, d}. Theorem Suppose n, d, e n < n and de ≤ e n are positive integers and assume e d+2 that f ∈ C [a − dh, b + dh] ∩ C 2d+1 ([a, a + e nh] ∪ [b − e nh, b]) is sampled at n + 1 equispaced nodes in [a, b]. Then kf − e rn k∞ ≤ KhD+1 . e rn is analytic, has no real poles and reproduces polynomials of degree ≤ D.

G. Klein

Extended Floater–Hormann interpolation

22/42

Introduction Extended FH Applications

Construction Convergence and stability Examples

A generalised Hermite formula (Platte’s method) Let L be a general linear approximation, n X L(x) = bi (x)fi . i=0

Suppose f is analytic inside a simple, closed, rectifiable curve C contained in a closed simply connected region around the nodes, then the error satisfies Z f (s) 1 f (x) − L(x) = R(x, s) ds, 2πi C s − x where R(x, s) = 1 − (s − x)

n X bi (x) s − xi i=0

is the relative error in the approximation of 1/(s − x) by L. G. Klein

Extended Floater–Hormann interpolation

23/42

Introduction Extended FH Applications

Construction Convergence and stability Examples

Interpolation of analytic functions −1 0.6

original FH

extended FH

−2 −3

0.4

−4 0.2

−5 −6

0 −7 −8

−0.2

−9 −0.4 −0.6 −1.5

−10 spline 15 −1

−11

funqui −0.5

0

0.5

1

1.5

−12

Figure: Comparison of the relative error of the prototype function 1/(s − x) with FH with d = 15, EFH with d = 15, de = e n = 15, B-spline of degree 15 and funqui for n = 100 on a log10 scale. G. Klein

Extended Floater–Hormann interpolation

24/42

Introduction Extended FH Applications

Construction Convergence and stability Examples

Lebesgue constant for extended FH interpolation

Theorem en associated For positive integers n and d, the Lebesgue constant Λ with extended Floater–Hormann interpolation is bounded as en ≤ 2 + log(n + 2d). Λ

Logarithmic growth with n and d.

G. Klein

Extended Floater–Hormann interpolation

25/42

Introduction Extended FH Applications

Construction Convergence and stability Examples

Lebesgue function for extended FH interpolation

(loading movie...)

G. Klein

Extended Floater–Hormann interpolation

26/42

Introduction Extended FH Applications

Construction Convergence and stability Examples

Lebesgue functions EFH

2nd kind

1st kind

3.5

3.5

3.5

3

3

3

2.5

2.5

2.5

2

2

2

1.5

1.5

1.5

1 −1

0

1

1 −1

0

1

1 −1

0

1

Figure: Lebesgue function with n = 50 for extended FH interpolation in equispaced nodes with d = 3 (left), for polynomial interpolation in Chebyshev points of the second kind (center) and first kind (right).

G. Klein

Extended Floater–Hormann interpolation

27/42

Introduction Extended FH Applications

Construction Convergence and stability Examples

Perturbed data with EFH and d = 15, de = e n = 15

(loading movie...)

G. Klein

Extended Floater–Hormann interpolation

28/42

Introduction Extended FH Applications

Construction Convergence and stability Examples

Interpolation of f (x) = sin(20x) perturbed data 2

10

exact data EFH FH funqui spline cheb

1

10

0

10

EFH FH funqui spline cheb

0

10

−5

10

−10

10

−1

10

50

100

150

200

50

100

150

200

Figure: Errors with 8 ≤ n ≤ 200 using EFH with d = 15, de = e n = 15, FH with d = 15, funqui, B-spline of degree 15, and polynomial interpolation in Chebyshev points.

G. Klein

Extended Floater–Hormann interpolation

29/42

Introduction Extended FH Applications

Construction Convergence and stability Examples

Interpolation of f (x) = sin(20x)

EFH FH spline funqui

−5

10

−10

10

5

10

15

20

25

30

Figure: Errors with 1 ≤ d ≤ 30 and n = 200 using EFH with de = e n = 15, FH and B-spline of degree d.

G. Klein

Extended Floater–Hormann interpolation

30/42

Introduction Extended FH Applications

Differentiation Integration

Applications

G. Klein

Extended Floater–Hormann interpolation

31/42

Introduction Extended FH Applications

Differentiation Integration

Approximation of derivatives via finite differences

Extended rational finite difference methods (ERFD) for the approximation at the nodes xi , 0 ≤ i ≤ n, of the k-th derivative of a sufficiently smooth function, (k)

f (k) (xi ) ≈ e rn (xi ) =

n+d X

(k) (k) Dij e fj =: e fi .

j=−d (k)

The weights Dij are the elements from the (d + 1)st to the (n + d + 1)st row of the (n + 2d + 1) × (n + 2d + 1) differentiation matrix D (k) .

G. Klein

Extended Floater–Hormann interpolation

32/42

Introduction Extended FH Applications

Differentiation Integration

Convergence of ERFD

Theorem e are positive Suppose n, d, e n < n, de ≤ e n and k ≤ D = min{d, d} integers and assume that e f ∈ C d+1+k [a − dh, b + dh] ∩ C 2d+1 ([a, a + e nh] ∪ [b − e nh, b]) is sampled at n + 1 equispaced nodes a = x0 < x1 < . . . < xn = b. Then (k)

|e rn (xj ) − f (k) (xj )| ≤ KhD+1−k ,

G. Klein

−d ≤ j ≤ n + d.

Extended Floater–Hormann interpolation

33/42

Introduction Extended FH Applications

Differentiation Integration

Approximation of derivatives at intermediate points For the approximation of the k-th derivative of a function f at intermediate points x ∈ [a, b], there are several possibilities, among them: interpolation of the derivatives at the nodes: en(k) (x) := R

n+d X i=−d

wi (k) e rn (xi ) x − xi

 n+d X i=−d

wi . x − xi

approximate e rn in Chebfun and differentiate the resulting chebfun. Usually not much loss in accuracy with such cheaper methods.

G. Klein

Extended Floater–Hormann interpolation

34/42

Introduction Extended FH Applications

Differentiation Integration

2nd derivative with FH and d = 15

(loading movie...)

G. Klein

Extended Floater–Hormann interpolation

35/42

Introduction Extended FH Applications

Differentiation Integration

2nd derivative with EFH and d = 15, de = e n = 15

(loading movie...)

G. Klein

Extended Floater–Hormann interpolation

36/42

Introduction Extended FH Applications

Differentiation Integration

2nd derivative of f (x) = 1/(1 + 100x 2 ) 5

10

EFH FH funqui spline cheb

0

10

−5

10

50

100

150

200

250

300

Figure: Errors with 20 ≤ n ≤ 300 using EFH with d = 15, de = e n = 15, FH with d = 15, funqui, B-spline of degree 15, and polynomial interpolation in Chebyshev points.

G. Klein

Extended Floater–Hormann interpolation

37/42

Introduction Extended FH Applications

Differentiation Integration

Integration of rational interpolants Linear interpolation formulas leads to a linear quadrature rules. For linear barycentric rational interpolant, e rn in particular, we have: Z

b

Z

b

f (x) dx ≈ a

a

Z e rn (x) dx = =

b

Pn+d

a n+d X

wk e k=−d x−xk fk wj j=−d x−xj

Pn+d

dx

ωk e fk ,

k=−d

where Z ωk := a

b

wk x−xk Pn+d wj j=−d x−xj

G. Klein

dx.

Extended Floater–Hormann interpolation

38/42

Introduction Extended FH Applications

Differentiation Integration

Convergence of EFH integration

Theorem Suppose n, d, e n < n and de ≤ e n are positive integers and assume e d+2 that f ∈ C [a − dh, b + dh] ∩ C 2d+1 ([a, a + e nh] ∪ [b − e nh, b]) is sampled at n + 1 equispaced nodes in [a, b]. Let the integral of e rn be approximated by a method of order at least O(hd+2 ). Then for any x ∈ [a, b], Z x Z x e f (y ) dy − rn (y ) dy ≤ KhD+2 log(n). a

a

G. Klein

Extended Floater–Hormann interpolation

39/42

Introduction Extended FH Applications

Differentiation Integration

Antiderivative of f (x) = 1/(1 + 16 cos(3x)2 ) 0

10

EFH FH funqui spline cheb

−10

10

50

100

150

200

250

300

Figure: Errors with 20 ≤ n ≤ 300 using EFH with d = 15, de = e n = 15, FH with d = 15, funqui, B-spline of degree 15, and polynomial interpolation in Chebyshev points.

G. Klein

Extended Floater–Hormann interpolation

40/42

Introduction Extended FH Applications

Summary

Floater–Hormann interpolation Stability and choice of d Extended Floater–Hormann interpolation Enhanced stability, convergence, parameter choice Applications: differentiation, integration

G. Klein

Extended Floater–Hormann interpolation

41/42

Introduction Extended FH Applications

Thank you for your attention!

G. Klein

Extended Floater–Hormann interpolation

42/42