Z uwagi na liczne zastosowania zadań interpolacyjnych Lagrange’a i Hermite’a w aproksymacji funkcji i w konstrukcji różnych metod numerycznych (np.
rozwiązywania równań nieliniowych i obliczania całek), duże znaczenie ma wzór opisujący resztę interpolacyjną.
Twierdzenie. Jeśli funkcja f jest klasy Cn+1 w przedziale A ⊂ R i h(x) oznacza wielomian interpolacyjny Hermite’a funkcji f dla węzłów x0, . . . , xn∈ A, to dla każdego x ∈ A istnieje liczba ξ ∈ A, taka że
f(x) − h(x) =f(n+1)(ξ) (n + 1)!pn+1(x).
Dowód. Jeśli x = xidla pewnego i ∈ {0, . . . , n}, to pn+1(x) = 0i dowodzona równość jest oczywista (z dowolnym ξ ∈ A). Dla ustalonego x ∈ A \ {x0, . . . , xn} uporządkujmy ciąg x0, . . . , xn, xtak, aby otrzymać ciąg niemalejący
x(0)0 6· · · 6 x(0)n+1. Określamy funkcję gx(s)def= f(s) − h(s) − zpn+1(s),
z parametrem z = z(x), który dobierzemy za chwilę. Korzystając z tego, że pn+1(x)6= 0, bierzemy
z = f(x) − h(x) pn+1(x) ,
i w ten sposób dostajemy gx(x) = 0. Tak określona funkcja spełnia warunek gx(x(0)i ) = 0dla i = 0, . . . , n + 1, tzn. ma co najmniej n + 2 miejsca zerowe8. Funkcja gxjest klasy Cn+1(A). Jej pochodna rzędu k 6 n + 1 ma co najmniej n + 2 − kmiejsca zerowe, które tworzą ciąg niemalejący x(k)0 6. . . 6 x(k)n+1−k. Istotnie, jeśli x(0)i < x(0)i+1, to (z twierdzenia Rolle’a) funkcja gx, która na końcach przedziału [x(0)i , x(0)i+1]przyjmuje tę samą wartość 0, osiąga w tym przedziale maksimum lub minimum, w punkcie x(1)i będącym miejscem zerowym funkcji gx′. Jeśli zaś funkcja gxma miejsce zerowe o krotności r > 1
(w węźle x(0)i =· · · = x(0)i+r−1), to jej pochodna ma w tym punkcie miejsce zerowe o krotności r − 1 (zatem mamy podciąg x(1)i =· · · = x(1)i+r−2). Korzystając
z indukcji, stosujemy to rozumowanie do kolejnych pochodnych. Wynika z niego, że pochodna rzędu n + 1 funkcji gxma w przedziale A co najmniej jedno miejsce zerowe, ξ = x(n+1)0 .
Podstawiając s = ξ, dostajemy
0 = g(n+1)x (ξ) = f(n+1)(ξ) − h(n+1)(ξ) − zp(n+1)n+1 (ξ).
Pochodna rzędu n + 1 wielomianu h(s) (stopnia n) jest równa 0, zaś pochodna wielomianu pn+1(s)(stopnia n + 1), którego współczynnik (w bazie potęgowej) przy sn+1jest równy 1, jest dla każdego s równa (n + 1)!. Zatem
z = f(n+1)(ξ) (n + 1)!.
Dowód zakończymy, wstawiając to do definicji funkcji gxi biorąc s = x. ✷ Przypuśćmy, że węzły x0, . . . , xnsą parami różne i weźmy dowolną liczbę
x /∈ {x0, . . . , xn}. Rozważmy wielomian interpolacyjny Lagrange’a h(s) funkcji f dla
8Jeśli pewien węzeł ma krotność r > 1, to liczymy go r razy; funkcja gxprzyjmuje w tym punkcie wartość 0 razem z pochodnymi rzędu 1, . . . , r − 1; oczywiście ciąg x(0), . . . , x(0)n+1nie jest wtedy ściśle rosnący.
węzłów x0, . . . , xni wielomian hn+1(s)stopnia co najwyżej n + 1, taki że hn+1(s) = f(s)dla każdego s ∈ {x0, . . . , xn, x}. Wtedy
hn+1(s) = h(s) + f[x0, . . . , xn, x]pn+1(s),
gdzie pn+1(s) = (s − x0)· . . . · (s − xn). Powyższa równość ma miejsce dla dowolnej funkcji f określonej w punktach x0, . . . , xn, x. W szczególności, dla s = x mamy
hn+1(x) = f(x) = h(x) + f[x0, . . . , xn, x]pn+1(x).
W przypadku, gdy funkcja f jest klasy Cn+1 w przedzialezawierającym wszystkie te punkty, na podstawie udowodnionego twierdzenia
f(x) = h(x) +f(n+1)(ξ) (n + 1)!pn+1(x), a ponieważ pn+1(x)6= 0, zachodzi równość
f[x0, . . . , xn, x] =f(n+1)(ξ) (n + 1)!
dla pewnego ξ należącego do tego przedziału. Stąd wynika fakt, wcześniej podany bez dowodu:
Wniosek. Jeśli funkcja f jest klasy Ck w otoczeniu punktu xi, to
xi+1,...,xlimi+k→xi
f[xi, . . . , xi+k] =f(k)(xi) k! .
Dowód. Wystarczy dla każdego układu liczb xi, . . . , xi+knależących do
dostatecznie małego otoczenia punktu xi przyjąć najkrótszy przedział zawierający te liczby i zauważyć, że potrzebną liczbę ξ możemy znaleźć w tym przedziale. ✷ Fakt ten jest podstawą definicji różnic dzielonych także dla węzłów o krotnościach większych niż 1. Wzór na resztę interpolacyjną w szczególnym przypadku, gdy x0=· · · = xn, jest wzorem Taylora (z resztą w postaci Lagrange’a). Inny przypadek szczególny, dla dwóch węzłów jednokrotnych, wykorzystaliśmy już w analizie metody siecznych. Dalsze zastosowania nastąpią.
Zadania i problemy
1. Rozwiąż zadanie interpolacyjne Lagrange’a (tj. skonstruuj bazę Newtona i oblicz współczynniki wielomianu interpolacyjnego w tej bazie) dla węzłów i wartości funkcji podanych w tabelce:
xi −2 −1 0 2 3 f(xi) 17 1 1 1 37
xi −2 −1 0 1 2 3
f(xi) −63 −6 −1 0 21 188 2. Rozwiąż zadanie interpolacyjne Hermite’a dla węzłów i wartości funkcji
i pochodnych podanych w tabelce:
xi 0 1 3
f(xi) −1 −2 80 f′(xi) −2 f′′(xi) 12
xi −2 −1 0 2
f(xi) −64 −7 −2 20
f′(xi) 1
f′′(xi) −2
3. Znajdź wyrażenia opisujące współczynniki kombinacji liniowej wartości funkcji f[a, b, c, d] = Af(a) + Bf(b) + Cf(c) + Df(d)
oraz wartości funkcji i pochodnej
f[a, b, b, c] = Af(z) + B0f(b) + B1f′(b) + Cf(c) w punktach a, b, c, d parami różnych.
4. Przejście między bazą Newtona i potęgową może być dokonane za pomocą schematu Hornera; zobaczmy, w jaki sposób. Rozważmy dzielenie z resztą wielomianu w(x) = anxn+· · · + a1x + a0 przez dwumian (x − x0):
w(x) = Xn
i=0
aixi=
Xn−1 i=0
ci+1xi
(x − x0) + c0= cnxn+ Xn−1
i=0
(ci− ci+1x0)xi. Stąd cn= anoraz ci= ci+1x0+ ai dla i < n. Zauważamy, że liczby cisą kolejno nadawanymi wartościami zmiennej w podczas obliczania za pomocą (zwykłego) schematu Hornera wartości w(x0). Otrzymujemy resztę z dzielenia,
b0= c0= w(x0), i współczynniki w bazie potęgowej ilorazu c(x). Jeśli określimy bazę Newtona, której elementami są wielomiany
q0(x) = 1, qi(x) =
Yi−1 k=1
(x − xk) dla i = 1, . . . , n − 1,
to współczynniki wielomianu c(x) w tej bazie są identyczne ze współczynnikami b1, . . . , bnwielomianu w w bazie Newtona {p0, . . . , pn}. Aby otrzymać b1, można następnie podzielić wielomian c(x) przez (x − x1), i tak dalej, rekurencyjnie.
Mamy zatem algorytm przejścia od bazy potęgowej do bazy Newtona:
for ( j = 0; j < n; j++ ) for ( i = n − 1; i > j; i-- )
a[i] += a[i + 1]∗ xj;
(algorytm ten zastępuje w tablicy a współczynniki wielomianu w bazie potęgowej współczynnikami w bazie Newtona). Algorytm przejścia w drugą stronę wykonuje operacje przeciwne w odwrotnej kolejności:
for ( j = n − 1; j > 0; j-- ) for ( i = j; i < n; i++ )
a[i] -= a[i + 1]∗ xj;
5. Algorytm Aitkena. Wartość w dowolnym punkcie x wielomianu interpolacyjnego Lagrange’a określonego przez węzły x0, . . . , xni wartości w tych węzłach y0, . . . , yn, może być obliczona (kosztem O(n2)operacji) bezpośrednio na podstawie tych danych. Przypuśćmy, że dwa wielomiany stopnia co najwyżej k − 1, p(k−1)i (x) i p(k−1)i+1 (x), są rozwiązaniami zadań interpolacyjnych Lagrange’a odpowiednio dla węzłów xi, . . . , xi+k−1oraz xi+1, . . . , xi+k. Wtedy wielomian
p(k)i (x) = xi+k− x
xi+k− xip(k−1)i (x) + x − xi
xi+k− xip(k−1)i+1 (x)
ma stopień co najwyżej k i jest rozwiązaniem zadania dla węzłów xi, . . . , xi+k. Istotnie, dla x = xipierwszy z ułamków w powyższym wzorze ma wartość 1, a drugi 0 (zatem p(k)i (xi) = p(k−1)i (xi) = yi, dla x = xi+kułamki mają wartości 0 i 1 (skąd wynika p(k)i (xi+k) = p(k−1)i+1 (xi+k) = yi+k), zaś dla każdego x ∈ R (w tym dla x∈ {xi+1, . . . , xi+k−1}) suma ułamków jest równa 1, a więc p(k)i (xj) = yjdla j = i + 1, . . . , i + k − 1.
Wartości wielomianów stopnia 0 w punkcie x, p(0)i (x) = yi, mamy za darmo.
Możemy zatem wpisać liczby y0, . . . , yndo tablicy p, i wykonać algorytm for ( k = 1; k 6 n; k++ )
for ( i = 0; i 6 n − k; i++ )
p[i] = (xi+k− x)∗ p[i] + (x − xi)∗ p[i + 1]/(xi+k− xi);
otrzymując w ten sposób wartość wielomianu interpolacyjnego w zmiennej p[0].
6. Zachodzi wzór Leibniza: jeśli funkcje g i h są odpowiednio gładkie i f(x) = g(x)h(x), to
f[x0, . . . , xk] = Xk
l=0
g[x0, . . . , xl]h[xl, . . . , xk].
Aby go udowodnić, określamy bazy Newtona:
pl(x) = Yl−1
j=0
(x − xj), qm(x) = Yk j=k−m+1
(x − xj).
Niech F(x) oznacza iloczyn wielomianów interpolacyjnych Hermite’a funkcji g i h opartych na węzłach x0, . . . , xk:
F(x) = Xk
l=0
g[x0, . . . , xl]pl(x) Xk m=0
h[xm, . . . , xk]qk−m(x).
Dla dowolnego węzła xi o krotności r funkcje F i f mają w tym węźle takie same wartości i pochodne do rzędu r − 1 włącznie. Możemy napisać
F(x) = X
06l6m6k
g[x0, . . . , xl]h[xm, . . . , xk]pl(x)qk−m(x) + X
06m<l6k
. . .
(składniki obu sum są opisane tym samym wzorem). Dla m < l iloczyn
wielomianów pl(x)qk−m(x)ma w węźle xio krotności r miejsce zerowe o krotności co najmniej r, a zatem pierwsza suma (której wszystkie składniki mają stopień nie większy niż k) reprezentuje wielomian interpolacyjny Hermite’a funkcji f, oparty na węzłach x0, . . . , xk. Współczynnik tego wielomianu odpowiadający
wielomianowi pkw bazie Newtona i jednocześnie xkw bazie potęgowej jest równy f[x0, . . . , xk]. Składniki stopnia k w pierwszej sumie mają l = m, a pozostałe mają stopień niższy. Pozostaje zauważyć, że współczynnik przy xk(w bazie potęgowej) pierwszej sumy jest równy wyrażeniu po prawej stronie wzoru Leibniza. ✷