Как найти коэффициенты многочлена лагранжа - Исправление недочетов и поиск решений вместе с Examum.ru

Рис. 2.4. Результат линейной интерполяции и значение в точке На рисунке 2.4 представлен результат линейной интерполяции Вторым

полиномом Ньютона и его значение в точке 5.

Перед нами стоит та же задача: пусть дана функция y = f (x), заданы значения этой функции yi = f (xi ), i=0,1,…,n в некоторых интерполяционных узлах X = {x0 , x1, , xn }, принадлежащих некоторому отрезку [a;b]. Требуется найти полином Ln (x) степени не выше п, принимающий в точках xi значения yi = Ln (xi ), i=0,1,…,n. Отметим, что в этой постановке не требуется чтобы h — шаг интерполяции был постоянный.

Для начала построим полином Ki (x), такой чтобы Ki (x j )= 0, если j ≠ i ; и

Ki (xi )=1, если j = i , j=0,1,…,n. Такой полином будет обращаться в нуль в п

точках x0 , x1, xi−1, xi=1, , xn , следовательно, он	будет	иметь вид
Ki (x)= Ci (x − x0 )(x − x1 ) (x − xi−1 )(x − xi+1 ) (x − xn ),	где Ci	константа,

которую следует подобрать так, чтобы Ki (xi )=1. Нетрудно заметить, что если

Ci (xi − x0 )(xi − x1 ) (xi − xi−1 )(xi − xi+1 ) (xi − xn )=1, то

— 36 —

Ci =

− x

)(x

− x ) (x

− x

i−1

)(x

− x

i+1

) (x

− x

)

Ki (x)=

− x0 )(x − x1 ) (x − xi−1 )(x

− xi+1 ) (x − xn )

(xi − x0 )(xi

− x1 ) (xi − xi−1 )(xi − xi+1 ) (xi − xn )

Теперь построим искомый полином Ln (x) с учетом исходных условий

= L

i=0,1,…,n.

Очевидно,

что

(x)

= ∑n K

(x) y

Этот

полином

i=0

удовлетворяет всем исходным требованиям:

•

степень данного полинома L

(x)=

∑n

(x) y

не превышает п, так как

i=0

степень каждого слагаемого Ki (x) yi не больше п deg Ln (x)≤ n ;

•

значения полинома Ln (x) в узлах интерполяции x j , j=0,1,…,n совпадают

со

значением

функции

y = f (x)

этих

точках.

Так

как

Ln (x j

)= ∑n Ki (x j ) yi = K j (x j ) y j = y j .

i=0

Таким образом, полином Ln (x) будет иметь вид

(x)=

(x − x

)(x − x

) (x − x

i−1

)(x − x

) (x − x

)

∑

(2.12)

− x

)(x

) (x

− x

)

i=0

− x

) (x

− x

i−1

− x

i+1

Теперь следует доказать,

что такой полином единственный.

Доказывать

будем методом от противного. Предположим, что существует другой полином M n (x), отличный от полинома Ln (x), причем deg M n (x)≤ n и M n (xi )= yi , i=0,1,…,n. Составим разность Q(x)= Ln (x)− M n (x). Полученный полином Q(x) обращается в нуль в n+1 точке x0 , x1, , xn и имеет степень degQn (x)≤ n . Следовательно, этот полином тождественно равен нулю Q(x)≡ 0 , то есть Ln (x)≡ M n (x). Таким образом, построенный полином Лагранжа единственный.

В случае, если узлы интерполяции являются равноотстоящими, интерполяционный полином Лагранжа совпадает с интерполяционным полиномом Ньютона.

— 37 —

Рассмотрим два частных случая: выведем формулы линейного и квадратичного интерполяционных полиномов Лагранжа.

Пример 1. Построить формулы линейного и квадратичного интерполяционных полиномов Лагранжа, воспользовавшись общей формулой.

В линейном случае нужно построить полином, проходящий через две точки с координатами (x0 , y0 ) и (x1, y1 ). Это означает, что слагаемых в формуле (2.12) останется два:

L (x)=	(x − x1 )		y	+	(x − x0 )	y .	(2.13)

1	(x	0	− x )	0		(x	− x	0	)	1
		1				1

В квадратичном (параболическом) случае нужно построить полином, проходящий через три точки с координатами (x0 , y0 ), (x1, y1 ) и (x2 , y2 ).

Следовательно,	слагаемых		в	формуле			(2.12)		останется	всего	три:
L		(x)=	(x − x1 )(x − x2 )			y		+	(x − x0 )(x − x2 )		y +	(x − x0 )(x − x1 )		y		. (2.14)
	2		(x	0	− x )(x	0	− x	2	)		0		(x	− x	0	)(x	− x	2	)	1	(x	2	− x	0	)(x	2	− x )	2
				1							1		1								1
		Пример 2. Дано множество точек	X = {− 2,0;1,2} и значения некоторой

функции y = f (x) в этих точках: f (− 2)= 0, f (0)= 2, f (1)=1, f (2)= 0. Требуется построить интерполяционный полином Лагранжа, проходящий через эти точки.

Воспользуемся формулой (2.12):

L3 (x)= (− 2(x−−00)()(−x2−−11)()(x−−22−)2) 0 +

(x + 2)(x −1)(x − 2)	2 +	(x + 2)(x − 0)(x − 2)	1	(x + 2)(x − 0)(x −1)	0 =
(0 + 2)(0 −1)(0 − 2)	(1+ 2)(1− 0)(1− 2)	+ (2 + 2)(2 − 0)(2 −1)

16 x3 − 12 x2 − 23 x + 2.

Ответ: L3 (x)= 16 x3 − 12 x2 − 23 x + 2.

— 38 —

Рис. 2.5. Интерполяционный полином Лагранжа На рисунке 2.5 представлен результат интерполяции интерполяционным

полиномом Лагранжа.

Вычисление коэффициентов полинома Лагранжа

Для вычисления коэффициентов в полиноме Лагранжа существует более простой способ, чем использование классической формулы (2.12). Для этого перепишем полином несколько в другой форме:

(x − x

)(x

− x

) (x − x

i−1

)(x − x

i+1

) (x − x

)

(x)= ∑

− x

)(x

− x ) (x

)(x

) (x

− x

)

i=0

− x

i−1

− x

i+1

= (x − x0 )(x − x1 ) (x − xn )×

, (2.15)

× ∑

(x − x

)(x

− x

)(x

− x

) (x

− x

)(x

− x

) (x

− x

)

i=0

i−1

i+1

Составим следующую матрицу:

− x

− x x

− x

− x0

x − x1

− xn

− x

по главной диагонали у нее сто

ят

сомножители

произведения

(x − x0 )(x − x1 ) (x − xn ),

а по строкам – сомножители знаменателя формулы

— 39 —

(2.15). Обозначив Di (x) произведение элементов по

Di (x)= (xi − x0 )(xi − x1 ) (xi − xi−1 )(x − xi )(xi − xi+1 ) (xi − xn ),

формулу (2.15) в виде

Ln (x)= (x − x0 )(x − x1 ) (x − xn )i∑n Dy(ix).

=0 i

i-й строке: перепишем

(2.16)

Это другой вид записи интерполяционной формулы Лагранжа.

Пример. Дано множество точек	X = {−1,0;3,5;} и значения некоторой
функции y = f (x) в этих точках:	f (−1)=1, f (0)=1, f (3)= 0, f (5)= −3.

Требуется построить интерполяционный полином Лагранжа, проходящий через

эти точки.

Составим

матрицу

для

вычисления

коэффициентов:

x +1

−1 − 0

−1 −3

−1 −5

x +1

−1

− 4

− 6

0 +1

x − 0

0 −3

0 −5

−3

−5

3 +1

3 − 0

x −3

3 −5

x −3

− 2

5 +1

5 − 0

5 −3

x −5

Вычислим коэффициенты полинома:

D0 (x)= (−1)(− 4)(−6)(x +1)= −24(x +1), D1(x)= (1)(−3)(−5)x =15x ,

D2 (x)= 4 3(− 2)(x −3)= −24(x −3),

D3 (x)= 6 5 2(x − 5)= 60(x − 5).

Подставим все коэффициенты в формулу (2.16), получим:

−3

L3 (x)

= (x +1)(x −

0)(x

−3)(x −5)

− 24(x +1)

15x

− 24(x −3)

60(x −5)

= −

−

x +1.

120

Ответ: L3 (x)= − 401 x3 − 301 x2 −1201 x +1.

— 40 —

Источник

From Wikipedia, the free encyclopedia

This image shows, for four points ((−9, 5), (−4, 2), (−1, −2), (7, 9)), the (cubic) interpolation polynomial L(x) (dashed, black), which is the sum of the scaled basis polynomials y₀ℓ₀(x), y₁ℓ₁(x), y₂ℓ₂(x) and y₃ℓ₃(x). The interpolation polynomial passes through all four control points, and each scaled basis polynomial passes through its respective control point and is 0 where x corresponds to the other three control points.

In numerical analysis, the Lagrange interpolating polynomial is the unique polynomial of lowest degree that interpolates a given set of data.

Given a data set of coordinate pairs ${displaystyle (x_{j},y_{j})}$ with the $x_{j}$ are called nodes and the $y_{j}$ are called values. The Lagrange polynomial L(x) has degree and assumes each value at the corresponding node, ${displaystyle L(x_{j})=y_{j}.}$

Although named after Joseph-Louis Lagrange, who published it in 1795,^[1] the method was first discovered in 1779 by Edward Waring.^[2] It is also an easy consequence of a formula published in 1783 by Leonhard Euler.^[3]

Uses of Lagrange polynomials include the Newton–Cotes method of numerical integration, Shamir’s secret sharing scheme in cryptography, and Reed–Solomon error correction in coding theory.

For equispaced nodes, Lagrange interpolation is susceptible to Runge’s phenomenon of large oscillation.

Definition[edit]

Given a set of nodes ${displaystyle {x_{0},x_{1},ldots ,x_{k}}}$ , which must all be distinct, ${displaystyle x_{j}neq x_{m}}$ for indices , the Lagrange basis for polynomials of degree for those nodes is the set of polynomials ${textstyle {ell _{0}(x),ell _{1}(x),ldots ,ell _{k}(x)}}$ each of degree which take values ${textstyle ell _{j}(x_{m})=0}$ if and ${textstyle ell _{j}(x_{j})=1}$ . Using the Kronecker delta this can be written ${textstyle ell _{j}(x_{m})=delta _{jm}.}$ Each basis polynomial can be explicitly described by the product:

${displaystyle {begin{aligned}ell _{j}(x)&={frac {(x-x_{0})}{(x_{j}-x_{0})}}cdots {frac {(x-x_{j-1})}{(x_{j}-x_{j-1})}}{frac {(x-x_{j+1})}{(x_{j}-x_{j+1})}}cdots {frac {(x-x_{k})}{(x_{j}-x_{k})}}\[10mu]&=prod _{begin{smallmatrix}0leq mleq k\mneq jend{smallmatrix}}{frac {x-x_{m}}{x_{j}-x_{m}}}.end{aligned}}}$

Notice that the numerator ${textstyle prod _{mneq j}(x-x_{m})}$ has roots at the nodes ${textstyle {x_{m}}_{mneq j}}$ while the denominator ${textstyle prod _{mneq j}(x_{j}-x_{m})}$ scales the resulting polynomial so that ${textstyle ell _{j}(x_{j})=1.}$

The Lagrange interpolating polynomial for those nodes through the corresponding values ${displaystyle {y_{0},y_{1},ldots ,y_{k}}}$ is the linear combination:

${displaystyle L(x)=sum _{j=0}^{k}y_{j}ell _{j}(x).}$

Each basis polynomial has degree , so the sum has degree , and it interpolates the data because ${textstyle L(x_{m})=sum _{j=0}^{k}y_{j}ell _{j}(x_{m})=sum _{j=0}^{k}y_{j}delta _{mj}=y_{m}.}$

The interpolating polynomial is unique. Proof: assume the polynomial of degree interpolates the data. Then the difference is zero at distinct nodes ${textstyle {x_{0},x_{1},ldots ,x_{k}}.}$ But the only polynomial of degree with more than roots is the constant zero function, so or

Barycentric form[edit]

Each Lagrange basis polynomial ${textstyle ell _{j}(x)}$ can be rewritten as the product of three parts, a function ${textstyle ell (x)=prod _{m}(x-x_{m})}$ common to every basis polynomial, a node-specific constant ${textstyle w_{j}=prod _{mneq j}(x_{j}-x_{m})^{-1}}$ (called the barycentric weight), and a part representing the displacement from ${textstyle x_{j}}$ to :^[4]

${displaystyle ell _{j}(x)=ell (x){dfrac {w_{j}}{x-x_{j}}}}$

By factoring out from the sum, we can write the Lagrange polynomial in the so-called first barycentric form:

${displaystyle L(x)=ell (x)sum _{j=0}^{k}{frac {w_{j}}{x-x_{j}}}y_{j}.}$

If the weights $w_{j}$ have been pre-computed, this requires only operations compared to ${displaystyle {mathcal {O}}(k^{2})}$ for evaluating each Lagrange basis polynomial $ell _{j}(x)$ individually.

The barycentric interpolation formula can also easily be updated to incorporate a new node $x_{k+1}$ by dividing each of the $w_{j}$ , by $(x_{j}-x_{k+1})$ and constructing the new $w_{k+1}$ as above.

For any ${textstyle sum _{j=0}^{k}ell _{j}(x)=1}$ because the constant function is the unique polynomial of degree leq k interpolating the data ${textstyle {(x_{0},1),(x_{1},1),ldots ,(x_{k},1)}.}$ We can thus further simplify the barycentric formula by dividing

${displaystyle {begin{aligned}L(x)&=ell (x)sum _{j=0}^{k}{frac {w_{j}}{x-x_{j}}}y_{j}{Bigg /}ell (x)sum _{j=0}^{k}{frac {w_{j}}{x-x_{j}}}\[10mu]&=sum _{j=0}^{k}{frac {w_{j}}{x-x_{j}}}y_{j}{Bigg /}sum _{j=0}^{k}{frac {w_{j}}{x-x_{j}}}.end{aligned}}}$

This is called the second form or true form of the barycentric interpolation formula.

This second form has advantages in computation cost and accuracy: it avoids evaluation of ell (x) ; the work to compute each term in the denominator $w_{j}/(x-x_{j})$ has already been done in computing ${displaystyle {bigl (}w_{j}/(x-x_{j}){bigr )}y_{j}}$ and so computing the sum in the denominator costs only addition operations; for evaluation points which are close to one of the nodes ${textstyle x_{j}}$ , catastrophic cancelation would ordinarily be a problem for the value ${textstyle (x-x_{j})}$ , however this quantity appears in both numerator and denominator and the two cancel leaving good relative accuracy in the final result.

Using this formula to evaluate L(x) at one of the nodes $x_{j}$ will result in the indeterminate ${displaystyle infty y_{j}/infty }$ ; computer implementations must replace such results by ${displaystyle L(x_{j})=y_{j}.}$

Each Lagrange basis polynomial can also be written in barycentric form:

${displaystyle ell _{j}(x)={frac {w_{j}}{x-x_{j}}}{Bigg /}sum _{m=0}^{k}{frac {w_{m}}{x-x_{m}}}.}$

A perspective from linear algebra[edit]

Solving an interpolation problem leads to a problem in linear algebra amounting to inversion of a matrix. Using a standard monomial basis for our interpolation polynomial ${textstyle L(x)=sum _{j=0}^{k}x^{j}m_{j}}$ , we must invert the Vandermonde matrix ${displaystyle (x_{i})^{j}}$ to solve $L(x_{i})=y_{i}$ for the coefficients $m_{j}$ of L(x) . By choosing a better basis, the Lagrange basis, ${textstyle L(x)=sum _{j=0}^{k}l_{j}(x)y_{j}}$ , we merely get the identity matrix, $delta _{ij}$ , which is its own inverse: the Lagrange basis automatically inverts the analog of the Vandermonde matrix.

This construction is analogous to the Chinese remainder theorem. Instead of checking for remainders of integers modulo prime numbers, we are checking for remainders of polynomials when divided by linears.

Furthermore, when the order is large, Fast Fourier transformation can be used to solve for the coefficients of the interpolated polynomial.

Example[edit]

We wish to interpolate $f(x)=x^{2}$ over the domain at the three nodes :

${displaystyle {begin{aligned}x_{0}&=1,&&&y_{0}=f(x_{0})&=1,\[3mu]x_{1}&=2,&&&y_{1}=f(x_{1})&=4,\[3mu]x_{2}&=3,&&&y_{2}=f(x_{2})&=9.end{aligned}}}$

The node polynomial ell is

${displaystyle ell (x)=(x-1)(x-2)(x-3)=x^{3}-6x^{2}+11x-6.}$

The barycentric weights are

${displaystyle {begin{aligned}w_{0}&=(1-2)^{-1}(1-3)^{-1}={tfrac {1}{2}},\[3mu]w_{1}&=(2-1)^{-1}(2-3)^{-1}=-1,\[3mu]w_{2}&=(3-1)^{-1}(3-2)^{-1}={tfrac {1}{2}}.end{aligned}}}$

The Lagrange basis polynomials are

${displaystyle {begin{aligned}ell _{0}(x)&={frac {x-2}{1-2}}cdot {frac {x-3}{1-3}}={tfrac {1}{2}}x^{2}-{tfrac {5}{2}}x+3,\[5mu]ell _{1}(x)&={frac {x-1}{2-1}}cdot {frac {x-3}{2-3}}=-x^{2}+4x-3,\[5mu]ell _{2}(x)&={frac {x-1}{3-1}}cdot {frac {x-2}{3-2}}={tfrac {1}{2}}x^{2}-{tfrac {3}{2}}x+1.end{aligned}}}$

The Lagrange interpolating polynomial is:

${displaystyle {begin{aligned}L(x)&=1cdot {frac {x-2}{1-2}}cdot {frac {x-3}{1-3}}+4cdot {frac {x-1}{2-1}}cdot {frac {x-3}{2-3}}+9cdot {frac {x-1}{3-1}}cdot {frac {x-2}{3-2}}\[6mu]&=x^{2}.end{aligned}}}$

In (second) barycentric form,

${displaystyle L(x)={frac {displaystyle sum _{j=0}^{2}{frac {w_{j}}{x-x_{j}}}y_{j}}{displaystyle sum _{j=0}^{2}{frac {w_{j}}{x-x_{j}}}}}={frac {displaystyle {frac {tfrac {1}{2}}{x-1}}+{frac {-4}{x-2}}+{frac {tfrac {9}{2}}{x-3}}}{displaystyle {frac {tfrac {1}{2}}{x-1}}+{frac {-1}{x-2}}+{frac {tfrac {1}{2}}{x-3}}}}.}$

Notes[edit]

Example of interpolation divergence for a set of Lagrange polynomials.

The Lagrange form of the interpolation polynomial shows the linear character of polynomial interpolation and the uniqueness of the interpolation polynomial. Therefore, it is preferred in proofs and theoretical arguments. Uniqueness can also be seen from the invertibility of the Vandermonde matrix, due to the non-vanishing of the Vandermonde determinant.

But, as can be seen from the construction, each time a node x_k changes, all Lagrange basis polynomials have to be recalculated. A better form of the interpolation polynomial for practical (or computational) purposes is the barycentric form of the Lagrange interpolation (see below) or Newton polynomials.

Lagrange and other interpolation at equally spaced points, as in the example above, yield a polynomial oscillating above and below the true function. This behaviour tends to grow with the number of points, leading to a divergence known as Runge’s phenomenon; the problem may be eliminated by choosing interpolation points at Chebyshev nodes.^[5]

The Lagrange basis polynomials can be used in numerical integration to derive the Newton–Cotes formulas.

Remainder in Lagrange interpolation formula[edit]

When interpolating a given function f by a polynomial of degree k at the nodes ${displaystyle x_{0},...,x_{k}}$ we get the remainder which can be expressed as^[6]

${displaystyle R(x)=f[x_{0},ldots ,x_{k},x]ell (x)=ell (x){frac {f^{(k+1)}(xi )}{(k+1)!}},quad quad x_{0}<xi <x_{k},}$

where ${displaystyle f[x_{0},ldots ,x_{k},x]}$ is the notation for divided differences. Alternatively, the remainder can be expressed as a contour integral in complex domain as

${displaystyle R(x)={frac {ell (x)}{2pi i}}int _{C}{frac {f(t)}{(t-x)(t-x_{0})cdots (t-x_{k})}}dt={frac {ell (x)}{2pi i}}int _{C}{frac {f(t)}{(t-x)ell (t)}}dt.}$

The remainder can be bound as

${displaystyle |R(x)|leq {frac {(x_{k}-x_{0})^{k+1}}{(k+1)!}}max _{x_{0}leq xi leq x_{k}}|f^{(k+1)}(xi )|.}$

Derivation^[7][edit]

Clearly, is zero at nodes. To find R(x) at a point ${displaystyle x_{p}}$ , define a new function and choose ${textstyle {tilde {R}}(x)=Ccdot prod _{i=0}^{k}(x-x_{i})}$ where is the constant we are required to determine for a given $x_{p}$ . We choose so that F(x) has k+2 zeroes (at all nodes and $x_{p}$ ) between $x_{0}$ and $x_{k}$ (including endpoints). Assuming that f(x) is k+1 -times differentiable, since L(x) and are polynomials, and therefore, are infinitely differentiable, F(x) will be k+1 -times differentiable. By Rolle’s theorem, ${displaystyle F^{(1)}(x)}$ has k+1 zeroes, ${displaystyle F^{(2)}(x)}$ has zeroes… ${displaystyle F^{(k+1)}}$ has 1 zero, say ${displaystyle xi ,,x_{0}<xi <x_{k}}$ . Explicitly writing ${displaystyle F^{(k+1)}(xi )}$ :

${displaystyle F^{(k+1)}(xi )=f^{(k+1)}(xi )-L^{(k+1)}(xi )-{tilde {R}}^{(k+1)}(xi )}$

${displaystyle L^{(k+1)}=0,{tilde {R}}^{(k+1)}=Ccdot (k+1)!}$

(Because the highest power of

)

${displaystyle 0=f^{(k+1)}(xi )-Ccdot (k+1)!}$

The equation can be rearranged as

${displaystyle C={frac {f^{(k+1)}(xi )}{(k+1)!}}}$

Since ${displaystyle F(x_{p})=0}$ we have ${displaystyle R(x_{p})={tilde {R}}(x_{p})={frac {f^{k+1}(xi )}{(k+1)!}}prod _{i=0}^{k}(x_{p}-x_{i})}$

Derivatives[edit]

The dth derivative of a Lagrange interpolating polynomial can be written in terms of the derivatives of the basis polynomials,

${displaystyle L^{(d)}(x):=sum _{j=0}^{k}y_{j}ell _{j}^{(d)}(x).}$

Recall (see § Definition above) that each Lagrange basis polynomial is

${displaystyle {begin{aligned}ell _{j}(x)&=prod _{begin{smallmatrix}m=0\mneq jend{smallmatrix}}^{k}{frac {x-x_{m}}{x_{j}-x_{m}}}.end{aligned}}}$

The first derivative can be found using the product rule:

${displaystyle {begin{aligned}ell _{j}'(x)&=sum _{begin{smallmatrix}i=0\inot =jend{smallmatrix}}^{k}{Biggl [}{frac {1}{x_{j}-x_{i}}}prod _{begin{smallmatrix}m=0\mnot =(i,j)end{smallmatrix}}^{k}{frac {x-x_{m}}{x_{j}-x_{m}}}{Biggr ]}\[5mu]&=ell _{j}(x)sum _{begin{smallmatrix}i=0\inot =jend{smallmatrix}}^{k}{frac {1}{x-x_{i}}}.end{aligned}}}$

The second derivative is

${displaystyle {begin{aligned}ell _{j}''(x)&=sum _{begin{smallmatrix}i=0\ineq jend{smallmatrix}}^{k}{frac {1}{x_{j}-x_{i}}}{Biggl [}sum _{begin{smallmatrix}m=0\mneq (i,j)end{smallmatrix}}^{k}{Biggl (}{frac {1}{x_{j}-x_{m}}}prod _{begin{smallmatrix}n=0\nneq (i,j,m)end{smallmatrix}}^{k}{frac {x-x_{n}}{x_{j}-x_{n}}}{Biggr )}{Biggr ]}\[10mu]&=ell _{j}(x)sum _{0leq i<mleq k}{frac {2}{(x-x_{i})(x-x_{m})}}\[10mu]&=ell _{j}(x){Biggl [}{Biggl (}sum _{begin{smallmatrix}i=0\inot =jend{smallmatrix}}^{k}{frac {1}{x-x_{i}}}{Biggr )}^{2}-sum _{begin{smallmatrix}i=0\inot =jend{smallmatrix}}^{k}{frac {1}{(x-x_{i})^{2}}}{Biggr ]}.end{aligned}}}$

The third derivative is

${displaystyle {begin{aligned}ell _{j}'''(x)&=ell _{j}(x)sum _{0leq i<m<nleq k}{frac {3!}{(x-x_{i})(x-x_{m})(x-x_{n})}}end{aligned}}}$

and likewise for higher derivatives.

Finite fields[edit]

The Lagrange polynomial can also be computed in finite fields. This has applications in cryptography, such as in Shamir’s Secret Sharing scheme.

References[edit]

^
Lagrange, Joseph-Louis (1795). «Leçon Cinquième. Sur l’usage des courbes dans la solution des problèmes». Leçons Elémentaires sur les Mathématiques (in French). Paris. Republished in Serret, Joseph-Alfred, ed. (1877). Oeuvres de Lagrange. Vol. 7. Gauthier-Villars. pp. 271–287. Translated as «Lecture V. On the Employment of Curves in the Solution of Problems». Lectures on Elementary Mathematics. Translated by McCormack, Thomas J. (2nd ed.). Open Court. 1901. pp. 127–149.
^ Waring, Edward (1779). «Problems concerning interpolations». Philosophical Transactions of the Royal Society. 69: 59–67. doi:10.1098/rstl.1779.0008.
^ Meijering, Erik (2002). «A chronology of interpolation: from ancient astronomy to modern signal and image processing» (PDF). Proceedings of the IEEE. 90 (3): 319–342. doi:10.1109/5.993400.
^ Berrut, Jean-Paul; Trefethen, Lloyd N. (2004). «Barycentric Lagrange Interpolation» (PDF). SIAM Review. 46 (3): 501–517. Bibcode:2004SIAMR..46..501B. doi:10.1137/S0036144502417715.
^ Quarteroni, Alfio; Saleri, Fausto (2003). Scientific Computing with MATLAB. Texts in computational science and engineering. Vol. 2. Springer. p. 66. ISBN 978-3-540-44363-6..
^ Abramowitz, Milton; Stegun, Irene Ann, eds. (1983) [June 1964]. «Chapter 25, eqn 25.2.3». Handbook of Mathematical Functions with Formulas, Graphs, and Mathematical Tables. Applied Mathematics Series. Vol. 55 (Ninth reprint with additional corrections of tenth original printing with corrections (December 1972); first ed.). Washington D.C.; New York: United States Department of Commerce, National Bureau of Standards; Dover Publications. p. 878. ISBN 978-0-486-61272-0. LCCN 64-60036. MR 0167642. LCCN 65-12253.
^ «Interpolation» (PDF).

External links[edit]

«Lagrange interpolation formula», Encyclopedia of Mathematics, EMS Press, 2001 [1994]
ALGLIB has an implementations in C++ / C# / VBA / Pascal.
GSL has a polynomial interpolation code in C
SO has a MATLAB example that demonstrates the algorithm and recreates the first image in this article
Lagrange Method of Interpolation — Notes, PPT, Mathcad, Mathematica, MATLAB, Maple at Holistic Numerical Methods Institute
Lagrange interpolation polynomial on www.math-linux.com
Weisstein, Eric W. «Lagrange Interpolating Polynomial». MathWorld.
Lagrange polynomial at ProofWiki
Dynamic Lagrange interpolation with JSXGraph
Numerical computing with functions: The Chebfun Project
Excel Worksheet Function for Bicubic Lagrange Interpolation
Lagrange polynomials in Python

Источник

Методы функциональной интерполяции

Постановка задачи

Пусть на множестве [math]Omega=[a,b][/math] задана сетка [math]Omega_n= {x_i,, i=overline{0,n}}[/math], определяемая [math]n+1[/math] точкой [math]x_0,x_1,ldots,x_n[/math], а на сетке задана сеточная функция [math]y_i=f(x_i),~ i=overline{0,n}[/math]

[math]y_0=f(x_0),~ y_1=f(x_1),~ ldots,~ y_n=f(x_n).qquad mathsf{(4.6)}[/math]

В некоторых случаях [math]y_i=f(x_i),~ i=overline{0,n}[/math] является сеточным представлением заданной формульной функции [math]y=f(x)[/math]. Сеточная функция может задаваться совокупностью пар: [math](x_0,y_0),(x_1,y_1),ldots,(x_n,y_n)[/math].

Требуется найти функцию [math]y=F(x)[/math], принимающую в точках [math]x_0,x_1,ldots, x_n[/math] те же значения, что и функция [math]y_i=f(x_i),~ i=overline{0,n}[/math], то есть [math]F(x_i)=y_i,~ i=overline{0,n}[/math].

Точки [math]x_0,x_1,ldots, x_n[/math] называются узлами интерполяции, а искомая функция [math]y=F(x)[/math] — интерполирующей.

Геометрически это означает, что нужно найти кривую, проходящую через заданное множество точек [math](x_i,y_i),~ i=overline{0,n}[/math] (рис. 4.2).

Одной из целей задачи интерполяции является вычисление значения функции в произвольной точке [math]x_{ast}[/math]. При этом различаются собственно интерполирование, когда точка [math]x_{ast}in [x_0,x_n][/math] и экстраполирование, когда [math]x_{ast}notin [x_0,x_n][/math].

Заметим, что можно провести бесчисленное множество «плавных» кривых, проходящих через заданное множество точек. Поэтому задача интерполяции в общей постановке не имеет единственного решения.

Если в качестве интерполирующей функции выбрать алгебраический многочлен, степень которого связана с числом заданных узлов интерполяции (на единицу меньше), решение задачи является единственным. Покажем это. Воспользуемся сначала кусочным способом. Выделим из отрезка [math][x_0,x_n][/math] частичный отрезок [x_i,x_{i+k}] и рассмотрим сеточную функцию [math]y_i=f(x_i)[/math], заданную в (k+1)-м узле [math]x_i,x_{i+1},ldots,x_{i+k}[/math] (узлы не совпадают):

[math]y_i=f(x_i),~ y_{i+1}=f(x_{i+1}),~ ldots,~ y_{i+k}= f(x_{i+k}).[/math]

(4.7)

В качестве интерполирующей функции выберем алгебраический многочлен k-й степени (степень многочлена на единицу меньше количества узлов):

[math]F(x)= widetilde{f}_k(x,overline{a})= sumlimits_{j=0}^{k} a_jx^j= a_0+a_1x+ ldots+a_kx^k.[/math]

(4.8)

Будем искать неизвестные коэффициенты [math]a_0,a_1,ldots,a_k[/math] из условия интерполяции (4.3) , т.е. [math]delta widetilde{f}_k(x_i,overline{a})=0colon[/math]

[math]widetilde{f}_k(x_i,overline{a})=y_i,~ widetilde{f}_k (x_{i+1}, overline{a})= y_{i+1},~ ldots,~ widetilde{f}_k(x_{i+k}, overline{a})= y_{i+k}.[/math]

(4.9)

Теорема 4.1 (о единственности решения задачи интерполяции). Задача о нахождении интерполяционного многочлена [math]textstyle{widetilde{f}_k(x,overline{a})= sumlimits_{j=0}^{k} a_jx^j}[/math], удовлетворяющего условиям (4.9), на частичном отрезке [math][x_i,x_{i+k}][/math] по заданной сеточной функции (4.7) имеет единственное решение.

Доказательство. Запишем условия интерполяции (4.9) с учетом (4.8) и обозначения [math]y_i= f(x_i)= f_icolon[/math]

[math]begin{aligned}& a_0+a_1x_i+ a_2x_i^2+ ldots+ a_kx_i^k=f_i,\ & a_0+a_1x_{i+1}+ a_2x_{i+1}^2+ ldots+ a_kx_{i+1}^k=f_{i+1},\ & quadvdots\ & a_0+a_1x_{i+k}+ a_2x_{i+k}^2+ ldots+ a_kx_{i+k}^k=f_{i+k}. end{aligned}[/math]

(4.10)

Эта система линейных алгебраических уравнений относительно коэффициентов [math]a_0,a_1, ldots,a_k[/math] имеет единственное решение, так как определитель матрицы системы

[math]C_k(x_i,x_{i+1},ldots,x_{i+k})= begin{pmatrix}1& x_i& x_i^2& cdots& x_i^k\ 1& x_{i+1}& x_{i+1}^2& cdots& x_{i+1}^k\ vdots& vdots& vdots& ddots& vdots\ 1& x_{i+k}& x_{i+k}^2& cdots& x_{i+k}^kend{pmatrix}[/math]

(4.11)

не равен нулю (доказательство последнего факта содержится в курсе линейной алгебры, где этот определитель называется определителем Вандермонда). Следовательно, задача интерполяции также имеет единственное решение.

Полагая [math]i=0,~ k=n[/math], приходим к глобальному способу решения поставленной задачи. А именно, если задана сеточная функция в (n+1)-м узле [math]x_0,x_1,ldots,x_n[/math] и требуется найти алгебраический многочлен с использованием условий интерполяции, то единственным решением задачи интерполяции является интерполяционный многочлен n-й степени:

[math]F(x)= widetilde{f}_n(x,overline{a})= sumlimits_{j=0}^{n} a_jx^j= a_0+a_1x+ ldots+ a_nx^n,[/math]

(4.12)

коэффициенты которого находятся из системы (4.10) при [math]i=0,~ k=n[/math].

Замечания

1. Матрица (4.11) имеет дискретный (точечный) характер, так как ее элементы вычисляются по дискретным значениям [math]x_i,x_{i+1},ldots,x_{i+k}[/math].

2. При решении поставленной задачи предполагается, что исходная сеточная функция задана своими точными значениями, хотя класс задач, для которых используются такие функции, ограничен.

3. Имеются и другие формы записи интерполяционных многочленов. По теореме 4.1 все эти многочлены степени [math]n[/math], удовлетворяющие функциональным условиям интерполяции и построенные по одним и тем же точкам, являются одним многочленом, записанным в разных формах.

4. При большом числе узлов решение системы (4.10) затруднительно. Искомый интерполяционный многочлен можно построить, не решая этой системы. Многочлены могут быть построены так, чтобы в самой структуре формулы многочлена условие интерполяции учитывалось.

5. При решении задачи функциональной интерполяции и в ее приложениях требуется:

а) выбрать наиболее удобную форму и степень интерполяционного многочлена. При этом можно использовать многочлены Лагранжа или Ньютона, а также формулу (4.8);

б) оценить погрешность интерполяции;

в) определить значения функции в точках, не совпадающих с узлами;

г) вычислить значения производных или определенных интегралов с использованием полученных интерполяционных многочленов.

Методика решения задачи интерполяции

1. По заданной сеточной функции составить интерполяционный многочлен определенной степени. При выборе степени многочлена следует руководствоваться желаемой точностью интерполяции.

2. Вычислить значения интерполяционного многочлена в заданных точках [math]x_{ast j},~ j=1,ldots,p[/math], путем их подстановки в формулу многочлена.

Многочлен Лагранжа

Пусть исходная сеточная функция задана в (n+1)-й точках сетки [math]Omega_ncolon, y_i=f(x_i),, i=overline{0,n}[/math], где [math]x_iin [a,b]=[x_0,x_n][/math] — в общем случае неравноотстоящие узлы, определяемые шагами [math]h_{i+1}= x_{i+1}-x_i~ (h_{i+1}=text{var})[/math].

Воспользуемся сначала кусочным способом. Здесь и далее будем использовать обозначение [math]f_i=f(x_i)[/math].

Выделим «окно» или частичный отрезок [math][x_i,x_{i+1}][/math], содержащий только две точки (шаблон [math](x_i,x_{i+1})[/math]). Тогда многочлен Лагранжа, интерполирующий исходную функцию на данном шаблоне, имеет вид (где [math]P_{1i}(x),,P_{1i+1}(x)[/math] — коэффициенты)

[math]L_1(x)= frac{(x-x_{i+1})}{(x_{i}-x_{i+1})}f_i+ frac{(x-x_{i})}{x_{i+1}-x_i}f_{i+1}= P_{1i}(x)cdot f_i+ P_{1i+1}(x)cdot f_{i+1}.[/math]

(4.13)

Действительно, легко убедиться в том, что [math]L_1(x)[/math] — алгебраический многочлен первой степени, который удовлетворяет условиям функциональной интерполяции (4.3) , т.е. [math]L_1(x_i)=f_i,~ L_1(x_{i+1})= f_{i+1}[/math].

Выделим «окно» в виде двойного частичного отрезка [math][x_i,x_{i+2}][/math] c шаблоном [math](x_i,x_{i+1},x_{i+2})[/math]. Тогда многочлен Лагранжа записывается в виде (где [math]P_{2i}(x),, P_{2i+1}(x),, P_{2i+2}(x)[/math] коэффициенты)

[math]begin{aligned}L_2(x)&= frac{(x-x_{i+1})(x-x_{i+2})}{(x_{i}-x_{i+1})(x_{i}-x_{i+2})}f_i+ frac{(x-x_{i})(x-x_{i+2})}{(x_{i+1}-x_{i})(x_{i+1}-x_{i+2})}f_{i+1}+ frac{(x-x_{i})(x-x_{i+1})}{(x_{i+2}-x_{i})(x_{i+2}-x_{i+1})}f_{i+2}=\ &= P_{2i}(x)cdot f_{i}+ P_{2i+1}(x)cdot f_{i+1}+ P_{2i+2}(x)cdot f_{i+2}.end{aligned}[/math]

(4.14)

Легко проверить, что (4.14) — многочлен второй степени и также удовлетворяет условиям функциональной интерполяции:

[math]L_2(x_i)= f_i;qquad L_2(x_{i+1})=f_{i+1},qquad L_2(x_{i+2})=f_{i+2}.[/math]

Обобщив запись многочлена на «окно» для к -кратного частичного отрезка [math][x_i,x_{i+k}][/math] с шаблоном [math](x_i,x_{i+1},ldots,x_{i+k})[/math], можно записать многочлен Лагранжа в виде

[math]L_k(x)= sumlimits_{m=i}^{i+k} frac{(x-x_i)(x-x_{i+1})ldots (x-x_{m-1})(x-x_{m+1})ldots (x-x_k)}{(x_m-x_i)(x_m-x_{i+1})ldots (x_m-x_{m-1})(x_m-x_{m+1})ldots (x_m-x_k)}f_m= sumlimits_{m=i}^{i+k} P_{km}(x)cdot f_m.[/math]

(4.15)

где [math]P_{km}(x)[/math] — коэффициенты Лагранжа, которые для внутренних точек шаблона записываются следующим образом:

[math]P_{km}(x)= frac{(x-x_i)(x-x_{i+1})ldots (x-x_{m-1})(x-x_{m+1})ldots (x-x_k)}{(x_m-x_i)(x_m-x_{i+1})ldots (x_m-x_{m-1})(x_m-x_{m+1})ldots (x_m-x_k)},.[/math]

(4.16)

Легко проверить, что [math]P_{km}(x)[/math] удовлетворяют условию

[math]P_{km}(x_j)= begin{cases}1,& j=m,\ 0,& jne m,end{cases} i leqslant j leqslant i+k.[/math]

(4.16)

Если положить [math]i=0,, k=n[/math], то приходим к глобальному способу решения задачи. Тогда интерполяционный многочлен Лагранжа n-й степени имеет вид

[math]L_n(x)= sumlimits_{i=0}^{n} frac{(x-x_0)(x-x_{1})ldots (x-x_{i-1})(x-x_{i+1})ldots (x-x_n)}{(x_i-x_0)(x_i-x_{1})ldots (x_i-x_{i-1})(x_i-x_{i+1})ldots (x_i-x_n)}f_i= sumlimits_{i=0}^{n} P_{ni}(x)cdot f_i,[/math]

(4.17)

где коэффициенты Лагранжа [math]P_{ni}(x)[/math] во внутренних точках отрезка записываются в форме

[math]P_{ni}(x)= frac{(x-x_0)(x-x_{1})ldots (x-x_{i-1})(x-x_{i+1})ldots (x-x_n)}{(x_i-x_0)(x_i-x_{1})ldots (x_i-x_{i-1})(x_i-x_{i+1})ldots (x_i-x_n)},.[/math]

Очевидно, многочлен [math]L_n(x)[/math], заданный равенством (4.17), является многочленом степени [math]n[/math] и удовлетворяет функциональным условиям интерполяции (4.3): [math]L_n(x_i)=f_i,~ i=overline{0,n}[/math]. Для записи интерполяционного многочлена Лагранжа удобно пользоваться табл. 4.1.

[math]begin{array}{|c|c|c|c|c|c|c|} multicolumn{7}{r}{mathit{Table~4.1}}\hline x-x_0& x_0-x_1& x_0-x_2& cdots& x_0-x_n& D_0& f_0 \hline x_1-x_0& x-x_1& x_1-x_2& cdots& x_1-x_n& D_1& f_1 \hline x_2-x_0& x_2-x_1& x-x_2& cdots& x_2-x_n& D_2& f_2 \hline vdots& vdots& vdots& ddots& vdots& vdots& vdots \hline x_n-x_0& x_n-x_1& x_n-x_2& cdots& x-x_n& D_n& f_n \hline multicolumn{5}{|c|}{begin{matrix}{}\[-8pt]{}end{matrix}Pi_{n+1}(x)= (x-x_0)(x-x_1)(x-x_2)ldots (x-x_n)}& D_i&f_i \hline end{array}[/math]

Здесь [math]D_i[/math] — произведение элементов i-й строки, [math]Pi_{n+1}(x)[/math] — произведение элементов главной диагонали, [math]y_i= f(x_i)= f_i,~ i=overline{0,n}[/math]. Тогда многочлен Лагранжа может быть записан в форме

[math]L_n(x)= Pi_{n+1}(x)cdot sumlimits_{i=0}^{n} frac{f_i}{D_i},.[/math]

(4.18)

Замечания

1. Если заданная сеточная функция такая, что [math]f_i=1~(i=overline{0,n})[/math], то из (4.17) следует, что [math]L_n(x)equiv1[/math], и поэтому справедливо равенство [math]textstyle{sumlimits_{i=0}^{n} P_{ni}(x)=1}[/math] (сумма всех коэффициентов Лагранжа в точке [math]x[/math] равна нулю), которое можно использовать для контроля правильности расчетов.

2. Коэффициенты Лагранжа [math]P_{ni}(x)[/math] для некоторой функции [math]y_i= f(x_i)[/math] определяются лишь узлами сетки и точкой [math]x[/math], в которой необходимо вычислить значение многочлена. Если в некоторой точке [math]x_{ast}[/math] требуется определить значения нескольких интерполируемых функций [math]f_1(x_i),f_2(x_i),ldots~ (i=overline{1,n})[/math], то коэффициенты Лагранжа для всех исходных функций подсчитываются только один раз.

3. При введении дополнительных узлов интерполяции все коэффициенты многочлена Лагранжа необходимо пересчитывать заново, что неудобно на практике. От этого недостатка свободны многочлены Ньютона.

Перейдем к рассмотрению примеров решения задачи интерполяции на основе вышеизложенной методики.

▼ Пример 4.1

Построить многочлен Лагранжа третьей степени для сеточной функции, заданной табл. 4.2. Вычислить значение функции в точке [math]x_{ast}=2,!5[/math]. Для записи многочлена использовать формулу (4.18).

[math]begin{array}{|c|c|c|c|c|} multicolumn{5}{r}{mathit{Table~4.2}}\hline i& 0& 1& 2&3 \hline x_i& 2& 3& 4&5 \hline begin{matrix}{}\[-8pt]{}end{matrix} f(x_i)= f_i & 7& 5& 8&7 \hline end{array}[/math]

Решение. 1. Построим многочлен Лагранжа. Для этого составим табл. 4.3, соответствующую табл. 4.1.

[math]begin{array}{|c|c|c|c|c|c|} multicolumn{6}{r}{mathit{Table~4.3}}\hline x-2&-1&-2&-3&-6cdot(x-2)&7 \hline 1&x-3&-1&-2& 2cdot(x-3)&5 \hline 2& 1& x-4&-1&-2cdot(x-4)&8 \hline 3& 2& 1& x-5& 6cdot(x-5)&7 \hline multicolumn{4}{|c|}{begin{matrix}{}\[-8pt]{}end{matrix} Pi_4(x)= (x-2)(x-3)(x-4)(x-5)}& D_i&f_i \hline end{array}[/math]

По формуле (4.18) получаем:

[math]begin{aligned}L_3(x) &= Pi_4(x) sumlimits_{i=0}^{3} frac{f_i}{D_i}=\ &=-frac{7}{6}(x-3)(x-4)(x-5)+ frac{5}{2}(x-2)(x-4)(x-5)-4(x-2)(x-3)(x-5),+\ &quad +,frac{7}{6} (x-2)(x-3)(x-4)=-frac{3}{2},x^3+16x^2-frac{107}{2},x+62. end{aligned}[/math]

2. Вычислим значение функции в заданной точке: [math]L_3(2,!5)= 4,!8125[/math].

Погрешность интерполяции многочленами Лагранжа

При определении значения [math]f(x),~ xne x_{i}[/math] для функции [math]y_i=f(x_i)~ (i=overline{0,n})[/math] с помощью многочлена Лагранжа возникает погрешность или остаточное слагаемое [math]R_n(x)colon[/math]

[math]f(x)= L_n(x)+ R_n(x).[/math]

(4.19)

Здесь предполагается, что используется глобальный способ интерполяции и что [math]f(x)in C_{n+1}[a,b][/math]. Последнее предположение требуется для применения соответствующих теорем математического анализа, однако, приведенное ниже соотношение для [math]R_n(x)[/math] может использоваться и для сеточных функций.

На основе указанных предположений доказано, что при интерполяции функции [math]y_i=f(x_i)[/math], заданной в общем случае на неравномерной сетке [math]Omega_n[/math], интерполяционным многочленом Лагранжа [math]textstyle{L_n(x)= sumlimits_{i=0}^{n} P_{ni}(x)f_i}[/math] для произвольного значения [math]xin (x_0,x_n)[/math] возникает погрешность

[math]R_n(x)= frac{f^{(n+1)}(xi)}{(n+1)!}cdot omega_n(x),[/math]

(4.20)

где [math]omega_n(x)= (x-x_0)(x-x_1)ldots (x-x_n)[/math] — многочлен (n+1)-й степени, а [math]xiin (a,b)[/math]. Поскольку точно найти [math]R_n(x)[/math] нельзя (из-за неопределенности точки [math]xi[/math]), то при проведении вычислений обычно находятся только приближенные оценки погрешностей интерполяции, которые являются априорными.

Оценка погрешности интерполяции в некоторой произвольной фиксированной точке [math]x_{ast}in[a,b][/math] имеет вид, где [math]M_{n+1}= max_{xin[a,b]}|f^{(n+1)}(x)|[/math]

[math]bigl|f(x_{ast}-L_n(x_{ast})bigr| leqslant frac{M_{n+1}}{(n+1)!}cdot bigl|omega_n(x_{ast})bigr|,[/math]

(4.21)

Оценка максимальной погрешности интерполяции в любой точке [math]xin[a,b][/math], т.е. на всем отрезке [math][a,b]colon[/math]

[math]bigl|f(x-L_n(x)bigr| leqslant frac{M_{n+1}}{(n+1)!}cdot max_{xin[a,b]}bigl|omega_n(x_{ast})bigr|.[/math]

(4.22)

Замечание. Для сеточных функций с фиксированными узлами сетки (узлами интерполяции) также можно проводить оценки погрешности по формулам (4.21), (4.22), однако для этого необходимо численно определять [math]M_{n+1}[/math] с помощью аппарата численного дифференцирования. При этом следует учитывать, что при вычислении производных высокого порядка возникают большие погрешности.

▼ Пример 4.2

С какой точностью можно вычислить значение [math]f(x)= Bigl.{sqrt{x} }Bigr|_{x_{ast}=85}[/math], если вычисления производить на основе интерполяционного многочлена Лагранжа первой и второй степени, а в качестве сеточной функции принять [math](x_i,y_i)= bigl{(16;4),(36;6),(100;10)bigr},~ (i=0,1,2)[/math]?

Решение. Так как требуется вычислить погрешность только в одной точке [math]x_{ast}=85[/math], то необходимо использовать формулу (4.21).

Найдем оценку погрешности в точке [math]x_{ast}=85[/math] для [math]L_1(x)[/math]. В качестве «окна» линейной интерполяции [math](n=1)[/math] выбираем отрезок [math][x_1,x_2]=[36;100][/math]. Определим величину [math]M_2colon[/math]

[math]f'(x)= frac{1}{2}x^{-1!!not{phantom{|}},2};qquad f»(x)=-frac{1}{4}x^{-3!!not{phantom{|}},2};qquad M_2= max_{xin[36;100]} bigl|f»(x)bigr|= frac{1}{4 sqrt{36^3}}= frac{1}{864},.[/math]

Тогда [math]omega_1(x_{ast})= (85-36)(85-100)=-735;~~ |R_1(x)| leqslant frac{M_2}{(n+1)!}|omega_1(x_{ast})|= frac{735}{2cdot864}= 0,!425[/math].

Найдем оценку погрешности в точке [math]x_{ast}[/math] для [math]L_2(x)[/math]. В качестве «окна» квадратичной интерполяции [math](n=2)[/math] выбираем отрезок [math][x_0,x_2]= [16;100][/math]. Определим величину [math]M_3colon[/math]

[math]f»'(x)= frac{3}{8}x^{-5!!not{phantom{|}},2};qquad M_3= max_{xin[16;100]} frac{3}{8 sqrt{16^5}}= frac{1}{2730,!7},.[/math]

Тогда [math]omega_2(x_{ast})= (85-16)(85-36)(85-100)=-50715;~~ |R_2(x)| leqslant frac{M_3}{(n+1)!}|omega_2(x_{ast})|=3,!1[/math].

Погрешность для многочлена [math]L_2(x)[/math] получилась больше погрешности для многочлена [math]L_1(x)[/math], что обусловлено величиной [math]omega_2(x_{ast})[/math].

Сходимости функционального интерполяционного процесса для непрерывных функций

Как отмечалось выше, выбор сетки и соответствующей степени интерполяционного многочлена при интерполяции сеточных функций является одной из важных задач, решить которую можно, рассмотрев проблему сходимости интерполяционных процессов [6,40] для непрерывных функций [math]y=f(x)in C_{n+1}[a,b][/math].

Будем считать, что интерполяция проводится на последовательности сеток

[math]Omega_1= bigl{x_0^{(1)},x_1^{(1)}bigr},quad Omega_2= bigl{x_0^{(2)}, x_1^{(2)},x_2^{(2)}bigr},quad ldots,quad Omega_n= bigl{x_0^{(n)},x_1^{(n)},ldots ,x_n^{(n)}bigr}[/math]

с возрастающим разбиением [math]k[/math] отрезка [math][a,b]colon, k_1=1; k_2=2; ldots; k_n=n[/math] и т.д.

Если при данных разбиениях при возрастающем [math]k=1,2,ldots,n,ldots[/math] определяются значения [math]L_k(x_{ast})[/math] в некоторой промежуточной точке [math]x_{ast}[/math], то реализуется интерполяционный процесс, характеризующийся последовательностью значений многочленов: [math]L_1(x_{ast}), L_2(x_{ast}), ldots, L_n(x_{ast}), ldots[/math].

Интерполяционный процесс для функции [math]f(x)[/math] сходится в точке [math]x_{ast}in [a,b][/math], если существует [math]limlimits_{ktoinfty} L_k(x_{ast})= f(x_{ast})[/math] (поточечная сходимость).

Для отрезка [math][a,b][/math] существует понятие равномерной сходимости в некоторой норме, например max, [math]max_{substack{xin[a,b]\ ntoinfty} bigl|f(x)-L_n(x)bigr|to0}[/math].

Характер сходимости или расходимости интерполяционного процесса зависит как от гладкости и поведения функции [math]f(x)[/math], так и от выбора последовательности сеток [math]Omega_n~ (n=1,2,ldots)[/math]. Так, показано, что если f(x) непрерывна на [math][a,b][/math], то найдется такая последовательность [math]Omega_n~ (n=1,2,ldots)[/math], для которой интерполяционный процесс сходится равномерно на [math][a,b][/math] (теорема Марцинкевича). Однако для дискретных функций, рассматриваемых в данном разделе, эта теорема не применима. Отметим также, что построены расходящиеся интерполяционные процессы и для формульных функций, например [math]f(x)=|x|,~ xin[-1;1][/math] и [math]f(x)= (1+25x^2)^{-1},~ xin[-1;1][/math]. Кроме того, применение многочленов высоких степеней приводит к так называемым «провалам» между узлами интерполяции, часто называемым осцилляциями. Указанные свойства интерполяционных процессов обусловливают нецелесообразность применения интерполяционных многочленов высоких степеней. В связи с этим в вычислительной практике для сеточных функций степень [math]n[/math] не берут выше [math]5div 8~ (nleqslant 5div 8)[/math] и задание частичного отрезка согласуют с выбранной степенью многочлена.

Рассмотрим часто использующиеся на практике линейную и параболическую интерполяцию.

Линейная и параболическая интерполяция с помощью многочлена Лагранжа

В прикладных расчетах часто применяется простейшая кусочная интерполяция, основанная на многочленах первой степени [math]L_1(x)[/math] или второй степени [math]L_2(x)[/math]. В этом случае функциональная интерполяция называется линейной или параболической (квадратичной) соответственно. Рассмотрим возможные способы их реализации и найдем оценки их погрешностей.

Пусть для сеточной функции [math]y_i=f(x_i)[/math], заданной на сетке [math]Omega_n= {x_0,x_1, ldots, x_n}[/math], требуется выполнить интерполяционный процесс для определения значения [math]f(x_{ast})[/math], где [math]x_{ast}ne x_i~ (i=overline{0,n})[/math], и оценить погрешности. Для обоснованного выбора степени интерполяционного многочлена необходимо указать, какую погрешность имеют значения исходной функции [math]f(x_i)[/math] в узлах. Если эта погрешность составляет величину [math]O(h_{i+1}^2)[/math] или [math]O(h_{i+1}^3[/math], а в широком классе вычислительных задач обеспечиваются именно такие погрешности, следует использовать линейную или параболическую интерполяцию.

Методика решения задачи линейной интерполяции

1. По расположению заданной точки [math]x_{ast}[/math] на оси [math]Ox[/math] выбрать из всех частичных отрезков [math][x_0,x_1], [x_1,x_2],ldots, [x_i,x_{i+1}], ldots, [x_{n-1},x_n][/math], заданных своими крайними значениями и в совокупности образующих сетку [math]Omega_n[/math], «окно» интерполяции [math]Oequiv [x_i,x_{i+1}][/math], такое, что [math]x_i< x_{ast}< x_{i+1}[/math].

2. Для отрезка [math][x_i,x_{i+1}][/math] вычислить значения коэффициентов Лагранжа [math]P_{1i}(x_{ast})= frac{x_{ast}-x_{i+1}}{x_i-x_{i+1}}[/math] и [math]P_{1i+1}(x_{ast})= frac{x_{ast}-x_i}{x_{i+1}-x_i}[/math] входящих в формулу (4.13) для многочлена [math]L_1(x)[/math]. Правильность полученных значений [math]P_{1i}(x_{ast}),, P_{1i+1} (x_{ast})[/math] проверить по условию [math]P_{1i}(x_{ast})+ P_{1i+1}(x_{ast})=1[/math], которое должно выполняться.

3. Вычислить искомое значение [math]f(x_{ast})[/math] согласно (4.13):

[math]f(x_{ast})approx L_1(x_{ast})= P_{1i}(x_{ast}) f_i+ P_{1i+1} (x_{ast}) f_{i+1}.[/math]

Как показано ниже, порядок этой аппроксимации равен двум, т.е. [math]O(h_{i+1}^2)[/math].

Геометрическая интерпретация линейной интерполяции при известной формульной функции [math]y=f(x)[/math] (штриховая линия) изображена на рис. 4.3. Здесь прямая [math]AB[/math] соответствует графику функции [math]y=L_1(x)[/math] на отрезке [math][x_i, x_{i+1}][/math]. Приближенное значение функции равно [math]L_1(x_{ast})[/math] ( точка [math]C[/math]), и оно отстоит от точного значения [math]f(x_{ast})[/math] на величину [math]CD approx O(h_{i+1}^2)[/math] вдоль оси [math]Oy[/math].

Методика решения задачи параболической интерполяции

1. Из всей совокупности спаренных частичных отрезков [math][x_0,x_2], [x_1,x_3], ldots, [x_{i-1}, x_{i+1}], [x_i, x_{i+2}], ldots, [x_{n-1},x_n][/math], образующих сетку [math]Omega_n[/math], по заданной величине [math]x_{ast}[/math] выбрать два пересекающихся «окна» [math]O_1equiv [x_{i-1},x_{i+1}],~ O_2equiv [x_i, x_{i+2}][/math] (предполагается, что [math]x_{ast}[/math] принадлежит внутреннему отрезку [math][x_i,x_{i+1}][/math]).

2. Для «окон» [math]O_1[/math] и [math]O_2[/math] вычислить значения коэффициентов Лагранжа: [math]P_{2,i-1}^{(1)}(x_{ast}),, P_{2,i}^{(1)}(x_{ast}),, P_{2,i+1}^{(1)}(x_{ast})[/math] — для «окна» [math]O_1[/math] и [math]P_{2,i}^{(2)}(x_{ast}),, P_{2,i+1}^{(2)}(x_{ast}),, P_{2,i+2}^{(2)}(x_{ast})[/math] — для «окна» [math]O_2[/math]. Путем суммирования проверить правильность полученных значений коэффициентов:

[math]sumlimits_{k=i-1}^{i+1} P_{2,k}(x_{ast})=1,qquad sumlimits_{k=i}^{i+2} P_{2,k}(x_{ast})=1.[/math]

3. Используя значения коэффициентов Лагранжа, вычислить значения [math]L_2^{(1)}(x_{ast}),, L_2^{(2)}(x_{ast})[/math] по формуле (4.14).

Если в расчетах не требуется высокая точность интерполирования, то можно ограничиться выбором одного «окна», например [math]O_2[/math], и тогда [math]f(x_{ast})approx L_2^{(2)}(x_{ast})[/math]. Для достижения повышенной точности интерполяцию провести для двух «окон» и результаты [math]L_2^{(1)}(x_{ast}),, L_2^{(2)}(x_{ast})[/math] усреднить:

[math]f(x_{ast})approx frac{1}{2}bigl[L_2^{(1)}(x_{ast})+ L_2^{(2)}(x_{ast})bigr]equiv L_{2text{sr}}(x_{ast}).[/math]

При этом порядок интерполяции повышается на единицу, т.е. [math]L_{2text{sr}}(x_{ast})[/math] аппроксимирует точное значение с четвертым порядком, поскольку погрешность составляет величину [math]O(H^4)[/math], где [math]H=max{h_i,h_{i+1},h_{i+2}}[/math].

Замечание. Если [math]x_{ast}in [x_0,x_2][/math] или [math]x_{ast}in [x_{n-1},x_n][/math], то выбирается одно «окно» [math][x_0,x_2][/math] или [math][x_{n-1},x_n][/math] соответственно.

Геометрическая интерпретация параболической интерполяции изображена на рис. 4.4. Параболе [math]y= L_2^{(1)}(x_{ast})[/math] соответствует кривая [math]A_1A_2B_1[/math], параболе [math]y= L_2^{(2)}(x_{ast})[/math] — кривая [math]A_2B_1A_2[/math]. Точка [math]C[/math] соответствует значению [math]L_2^{(1)}(x_{ast})[/math], точка [math]D[/math] — значению [math]L_2^{(2)}(x_{ast})[/math] точка [math]E[/math] — значению [math]L_{2text{sr}} (x_{ast})[/math].

Приведем оценки погрешностей линейной и параболической интерполяции. Вначале предположим, что сетка [math]Omega_n[/math] равномерная (это имеет значение только для параболической интерполяции). Формулы (4.13), (4.14) и оценки (4.22), записанные для «окон» интерполяции, упрощаются, если ввести в рассмотрение новую переменную — фазу интерполяции и [math]u=frac{x-x_i}{h}colon[/math]

[math]begin{gathered}L_1(u)= (1-u)cdot f_i+ ucdot f_{i+1},\ L_2(u)= frac{(u-1)(u-2)}{2}cdot f_i-u(u-2)cdot f_{i+1}+ frac{u(u-1)}{2}cdot f_{i+2}. end{gathered}[/math]

Здесь учтено, что [math]x-x_{i+1}= x-(x_i+h)= uh-h= h(u-1),~ x-x_{i+2}= h(u-2)[/math]. Очевидно, что для [math]L_1(u)[/math] величина [math]{u}[/math] изменяется в диапазоне [math]0 leqslant u leqslant 1[/math], а для [math]L_2(u)[/math] — в диапазоне [math]0 leqslant u leqslant 2[/math]. Для получения мажорант в оценке (4.22) необходимо найти [math]max|omega_1(x)|[/math] и [math]max|omega_2(x)|[/math]. Преобразуя зависимости [math]omega_1(x)[/math] и [math]omega_2(x)[/math] к новой переменной [math]{u}[/math], получаем [math]omega_1(u)= h^2ucdot (1-u),[/math] [math]omega_2(u)= h^3ucdot (1-u)(2-u)[/math] и находим максимумы:

[math]Omega^1= max_{0 leqslant u leqslant 1} bigl|omega_1(u)bigr|= frac{h^2}{4},,qquad Omega^2= max_{0 leqslant u leqslant 2} bigl|omega_2(u)bigr|= frac{2 sqrt{3}h^2}{9},.[/math]

Таким образом, реализуются следующие оценки погрешностей линейной и параболической интерполяции, справедливых для соответствующих «окон»:

[math]bigl|f(u)-L_1(u)bigr|leqslant frac{h^2}{8},M_{2,i},[/math]

(4.23)

[math]bigl|f(u)-L_2(u)bigr|leqslant frac{sqrt{3}h^3}{27},M_{3,i},[/math]

(4.24)

где [math]M_{2,i}= max_{[x_i, x_{i+1}]}bigl|f»(x)bigr|,~ M_{3,i}= max_{[x_i, x_{i+2}]}bigl|f»'(x)bigr|[/math]. При этом предполагается, что [math]f(x)[/math] принадлежит классам функций [math]f(x)in C_2[a,b],~ f(x)in C_3[a,b][/math] соответственно для линейной и параболической интерполяции. Таким образом, из оценок (4.23), (4.24) следует, что линейная интерполяция обеспечивает на частичном отрезке [math][x_i, x_{i+1}][/math] второй порядок аппроксимации или погрешности по [math]h[/math], а параболическая (без осреднения) на двойном отрезке [math][x_i, x_{i+2}][/math] — третий порядок. Данные пофешности, как отмечалось во введении и предыдущих разделах, сокращенно записываются как [math]O(h^2)[/math] и [math]O(h^3)[/math]. При реализации алгоритма с осреднением порядок параболической интерполяции становится равным [math]O(h^4)[/math].

Замечания

1. Оценка (4.23) для линейной интерполяции инвариантна по отношению к виду сетки [math]Omega_n[/math] (равномерной или неравномерной). Параболическая интерполяция также сохраняет указанную погрешность при выполнении интерполяции на неравномерной сетке [math]Omega_n[/math]. В некоторых источниках для [math]L_2(x)[/math] при выполнении условия [math]f(x)in C_3[a,b][/math] приведена оценка

[math]max_{[x_i,x_{i+2}]} bigl|R_2(x)bigr| leqslant frac{sqrt{3}}{27},H^3M_{3,i},[/math]

(4.25)

где [math]H= max(h_{i+1}, h_{i+2}),~ h_{i+2}= x_{i+2}-x_{i+1},~ M_{3,i}= max_{[x_i,x_{i+2}]} bigl|f»'(x)bigr|[/math].

2. Если гладкость функции [math]f(x)[/math] не достигает вышеуказанной и класс гладкости понижен на единицу [math]bigl(f(x)in C_2[a,b]bigr)[/math], то порядок параболической интерполяции также понижается на единицу:

[math]max_{[x_i,x_{i+2}]} bigl|R_2(x)bigr| leqslant 0,!1546cdot H^2M_{2,i}.[/math]

(4.26)

3. Повышение класса гладкости функции [math]f(x)[/math] выше [math]C_3[a,b][/math] не приводит к увеличению порядка параболической интерполяции относительно Н , т.е. происходит как бы его «замораживание» или «насыщение». Указанные свойства, вытекающие из оценок (4.25), (4.26), носят общий характер и имеют место в других оценках аппроксимации многочленов, сплайнов, производных и интегралов.

4. Для произвольной степени интерполяционного многочлена при [math]h=text{const},~ f(x)in C_{n-1}[a,b][/math] погрешность функциональной интерполяции на отрезке [math][x_0,x_n][/math] выражается следующим образом:

[math]bigl|f(x)-L_n(x)bigr|_{[a,b]}= max_{xin [a,b]} bigl|f(x)-L_n(x)bigr| leqslant h^{n+1} frac{M_{n+1}}{(n+1)!},Omega^n.[/math]

(4.27)

где [math]Omega^n= max_{0 leqslant u leqslant n} bigl|omega_n(u)bigr|,~ M_{n+1}= max_{[x_0,x_n]} bigl|f^{(n+1)}(x)bigr|,~ omega_n(x)= u(u-1)cdot ldotscdot (u-n)[/math]. Для многочленов с [math]n leqslant 5[/math] величины [math]Omega^n[/math] или их оценки являются такими:

[math]Omega^1=frac{1}{4};quad Omega^2=frac{2 sqrt{3}}{9};quad Omega^3=1;quad Omega^4<3,!7;quad Omega^5<17.[/math]

Из (4.27) следует, что на отрезке [math][x_0,x_n][/math] величина [math]|R_n(x)|[/math] есть [math]O(h^{n+1})[/math] и при уменьшении шага [math]h[/math] в два раза мажоранта уменьшается по крайней мере в [math]2^{n+1}[/math] раз. Отсюда вытекает, что для непрерывных функций можно по заданной точности интерполяции выбирать шаг [math]h[/math]. При этом можно путем использования кусочной интерполяции в некоторых пределах изменять степень интерполяционного многочлена. Мажоранту в оценке (4.27) при кусочной интерполяции можно также снизить путем выбора «окна» интерполяции [math][x_i,x_{i+k}][/math] так, чтобы точка [math]x_{ast}[/math] располагалась как можно ближе к его середине. Это обусловлено тем, что колебания функции [math]omega_n(x)[/math] вблизи середины [math][x_i,x_{i+k}][/math] меньше, чем у его концов.

5. В широком классе задач математической физики применяются расчетные схемы в основном второго (и иногда третьего и выше) порядка точности. При их реализации, как правило, используются встроенные интерполяционные алгоритмы, основанные на многочленах и сплайн-функциях. Степени интерполяционных многочленов при этом должны выбираться из условия соответствия порядков их аппроксимации порядкам точности схем. Если эти порядки одинаковы, то порядок точности схем сохраняется, хотя константа в оценке погрешности схемы изменяется. Если же порядок встроенных интерполяционных алгоритмов хотя бы на единицу выше порядка точности схемы, то вместе с порядком точности схемы сохраняется и указанная константа. Отсюда следует, что необходимо выбирать такую степень интерполяционного многочлена, которая либо обеспечивает равенство порядка аппроксимации порядку точности схемы, либо на единицу превышает последний. Таким образом, использование параболической интерполяции в качестве встроенных алгоритмов или для восполнения численных решений, полученных по схемам второго порядка, позволяет сохранить требуемую точность расчета, а также не дает избыточный порядок и, следовательно, не усложняет алгоритм. Это замечание носит общий характер и относится к любым аппроксимационным алгоритмам, выполняющим функцию восполнения или интерполирования.

Перейдем к рассмотрению примеров решения задач линейной и параболической интерполяции.

▼ Пример 4.3

Дана сеточная функция, являющаяся сеточным представлением формульной функции [math]f(x)=x^3[/math] (табл. 4.4). Найти значение [math]f(x_{ast})[/math] при [math]x_{ast}=2[/math] с помощью линейной и параболической интерполяции.

[math]begin{array}{|c|c|c|c|c|c|} multicolumn{6}{r}{mathit{Table~4.4}}\hline i& 0& 1& 2& 3& 4 \hline x_i& -1& 0& 1& 3& 4 \hline f_i& -1& 0& 1& 27& 64 \hline end{array}[/math]

Решение. Применение линейной интерполяции. Воспользуемся методикой.

1. Выбираем «окно» интерполяции [math][x_i,x_{i+1}]= [1;3]~ (i=2,~ i+1=3)[/math].

2. Вычислим коэффициенты Лагранжа: [math]P_{1,2}= frac{x_{ast}-x_3}{x_2-x_3}= 0,!5;~~ P_{1,3}= frac{x_{ast}-x_2}{x_3-x_2}=0,!5[/math]. Вычисления выполнены верно, так как [math]P_{1,2}+P_{1,3}=1[/math].

3. Определим искомое значение [math]f(2)[/math] по формуле (4.13):

[math]f(2)approx L_1(2)= P_{1,2}cdot f_2+ P_{1,3}cdot f_3=14.[/math]

Применение параболической интерполяции с осреднением. Также воспользуемся соответствующей методикой.

1. Выбираем «окна» интерполяции [math]O_1equiv [x_1;x_3]= [0;3];~ O_2equiv [x_2;x_4]= [1;4][/math].

2. Вычислим коэффициенты Лагранжа [math]P_{21}^{(1)},, P_{22}^{(1)},, P_{23}^{(1)}[/math] для «окна» [math]O_1[/math] и [math]P_{22}^{(2)},, P_{23}^{(2)},, P_{24}^{(2)}[/math] для «окна» [math]O_2[/math]. Получим

[math]P_{21}^{(1)}=-frac{1}{3},quad P_{22}^{(1)}=1,quad P_{23}^{(1)}=frac{1}{3};qquad P_{22}^{(2)}=frac{1}{4},quad P_{23}^{(2)}=1,quad P_{24}^{(2)}=-frac{1}{3},.[/math]

Вычисления выполнены правильно, так как [math]textstyle{sumlimits_{k=1}^{3} P_{2k}^{(1)}=1;~ sumlimits_{k=2}^{4} P_{2k}^{(2)}=1}[/math].

3. Определим значения [math]L_2^{(1)}(x_{ast}),~ L_2^{(2)}(x_{ast})[/math] по формуле (4.14):

[math]L_2^{(1)}(x_{ast})=-frac{1}{3}cdot0+1cdo1+ frac{1}{3}cdot27=10;qquad L_2^{(2)}(x_{ast})= frac{1}{3}cdot1+ 1cdot27-frac{1}{3}cdot64=6.[/math]

Искомое значение получим в результате осреднения: [math]f(2)approx frac{L_2^{(1)}(x_{ast})+ L_2^{(2)}(x_{ast})}{2}=8[/math]. Найдено точное значение [math]f(2)=8[/math], что объясняется тем, что оператор осреднения повышает порядок интерполяции. В этом случае в остаточное слагаемое входит [math]f^{(4)}(xi)equiv0[/math] и оно равно нулю.

Многочлены Ньютона

Разделенные и конечные разности. В практике функционального интерполирования иногда удобнее использовать многочлены Ньютона, степень которых можно последовательно повышать путем добавления очередных слагаемых, имеющих более высокую степень. Такие несимметричные многочлены, альтернативные симметричным многочленам Лагранжа, основаны на разделенных и конечных разностях, вычисляемых по интерполируемой сеточной функции.

Разделенные разности вводятся для функции [math]y_i= f(x_i)=f_i,~ i=overline{0,n}[/math], заданной на неравномерной сетке [math](h_{i+1}= text{var})[/math], а конечные разности — для функции [math]y_i= f(x_i)=f_i,~ i=overline{0,n}[/math], определенной на равномерной сетке [math](h_{i+1}= text{const})[/math].

Выбрав внутри неравномерной или равномерной сетки соответствующие шаблоны интерполяции [math](x_i,x_{i+1}), (x_i,x_{i+1},x_{i+2}), ldots, (x_i,x_{i+1},ldots, x_{i+k})[/math] введем следующие определения разделенных и конечных разностей:

– разделенная разность первого порядка: [math]f(x_i,x_{i+1})= frac{f_{i+1}-f_i}{x_{i+1}-x_i}[/math];

– разделенная разность второго порядка: [math]f(x_i,x_{i+1},x_{i+2})= frac{f(x_{i+1},x_{i+2})-f(x_i,x_{i+1})}{x_{i+2}-x_i}[/math];

– разделенная разность k-го порядка: [math]f(x_i,x_{i+1},ldots, x_{i+k})= frac{f(x_{i+1},x_{i+2},ldots, x_{i+k})-f(x_i,x_{i+1},ldots, x_{i+k-1})}{x_{i-k}-x_i}[/math];

– конечная разность первого порядка: [math]Delta f_i= f_{i+1}-f_i[/math];

– конечная разность второго порядка: [math]Delta^2f_i= Delta (Delta f_i)= Delta f_{i+1}-Delta f_i= f_{i+2}-2f_{i+1}+f_i[/math];

– конечная разность k-го порядка: [math]textstyle{Delta^kf_i= Delta (Delta^{k-1}f_i)= sumlimits_{j=0}^{k} (-1)^jC_k^j f_{i+j}}[/math], где [math]C_k^j= frac{k!}{(k-j)!,j!}[/math].

Последовательность получения разделенных и конечных разностей при [math]k=3[/math] для произвольной функции наглядно представляют табл. 4.5 и 4.6.

Для гладких функций числовые значения [math]f(x_j, x_{j+1}, ldots, x_{j+k})[/math] и [math]Delta^kf_i[/math] при возрастании [math]k[/math] уменьшаются и стремятся к нулю, т.е. [math]f(x_j, x_{j+1}, ldots, x_{j+k})to0[/math] и [math]Delta^kf_ito0[/math] при [math]ktoinfty[/math].

Связь между разделенными и конечными разностями k-го порядка при [math]h=text{const}[/math] устанавливается следующим соотношением:

[math]f(x_i,x_{i+1},ldots,x_{i+k})= frac{Delta^kf_i}{k!,h^k},.[/math]

(4.28)

Действительно, при [math]k=1[/math] для разделенной разности [math]f(x_i,x_{i+1})[/math] получаем [math]f(x_i,x_{i+1})= frac{Delta f_i}{h}= frac{Delta f_i}{1!,h^1}[/math], то есть (4.28) при [math]k=1[/math] справедливо.

Пусть [math]k=2[/math]. Тогда получаем

[math]f(x_i,x_{i+1},x_{i+2})= frac{f(x_{i+1},x_{i+2})-f(x_{i},x_{i+1})}{}= frac{1}{2h}! left(frac{f_{i+2}-f_{i+1}}{h}-frac{f_{i+1}-f_i}{h}right)= frac{Delta f_{i+1}-Delta f_i}{2h^2}= frac{Delta^2f_i}{2h^2},.[/math]

Таким образом, связь (4.28) выполняется и при [math]k=2[/math]. Справедливость (4.28) при произвольном к можно доказать методом математической индукции.

Интерполяционный многочлен Ньютона для неравномерной сетки

Пусть исходная (интерполируемая) сеточная функция [math]y_i=f(x_i),~ i=overline{0,n}[/math], задана на неравномерной сетке [math]Omega_nequiv {x_0,x_1,x_2,ldots,x_n}[/math], характеризующейся шагами [math]h_{i+1}= x_{i+1}-x_i=text{var}[/math].

Воспользуемся сначала кусочным способом. Из всей совокупности узлов выбираем шаблон [math](x_i,x_{i+1},ldots,x_{i+k})[/math], соответствующий некоторому «окну» интерполяции [math][x_i,x_{i+k}][/math].

Тогда для функциональной интерполяции может быть использован многочлен Ньютона, основанный на разделенных разностях:

[math]begin{aligned}N_k(x)= f_i&+ f(x_i,x_{i+1})cdot (x-x_i)+ f(x_i,x_{i+1}, x_{i+2})cdot (x-x_i)cdot (x-x_{i+1})+ ldots+\ &+, f(x_i,x_{i+1},ldots,x_{i+k})cdot (x-x_i)cdot (x-x_{i+1})cdot ldotscdot (x-x_{i+k-1}). end{aligned}[/math]

(4.29)

Действительно, [math]N_k(x)[/math] — многочлен k-й степени, что определяется сомножителями последнего слагаемого (разделенные разности, входящие в качестве одного из сомножителей в эти произведения, есть числа). Кроме того, для многочлена [math]N_k(x)[/math] удовлетворяются функциональные условия интерполяции: [math]N_k(x_j)= f_j,~ j=i,ldots,i+k[/math]. Проверим их справедливость при [math]k=1[/math] (шаблон [math](x_i,x_{i+1})[/math]) и [math]k=2[/math] (шаблон [math](x_i,x_{i+1},x_{i+2}})[/math]). Пусть [math]k=1[/math]. Тогда

[math]N_1(x)= f_i+f(x_i,x_{i+1})cdot (x-x_i),[/math]

(4.30)

и поэтому

[math]N_1(x_i)= f_i,;qquad N_1(x_{i+1})= f_i+ frac{f_{i+1}-f_i}{x_{i+1}-x_i}(x_{i+1}-x_i)= f_{i+1}.[/math]

Таким образом, условия интерполяции для [math]N_1(x)[/math] выполнены, следовательно, многочлен (4.30) может быть использован для линейной интерполяции кусочным способом. Пусть [math]k=2[/math]. Тогда

[math]N_2(x)= f_i+ f(x_i,x_{i+1})(x-x_i)+ f(x_i,x_{i+1},x_{i+2}) (x-x_i)(x-x_{i+1}),[/math]

(4.31)

и поэтому

[math]begin{gathered}N_2(x)=f_i;qquad N_2(x_{i+1})= f_i+ frac{f_{i+1}-f_i}{x_{i+1}-x_i}(x_{i+1}-x_i)+0= f_{i+1},,\ N_2(x_{i+2})= f_i+ frac{f_{i+1}-f_i}{x_{i+1}-x_i}(x_{i+2}-x_i)+ frac{1}{x_{i+2}-x_i}! left(frac{f_{i+2}-f_{i+1}}{x_{i+2}-x_{i+1}}-frac{f_{i+1}-f_{i}}{x_{i+1}-x_{i}}right)!(x_{i+2}-x_i)(x_{i+2}-x_{i+1})=f_{i+2}. end{gathered}[/math]

Таким образом, условия интерполяции для многочлена [math]N_2(x)[/math] также выполнены и он может использоваться для параболической интерполяции кусочным способом.

Для произвольного [math]k[/math] справедливость равенств [math]N_k(x_j)=f_i,~ j=overline{i,i+k}[/math], проверяется методом математической индукции.

Полагая [math]i=0,~ k=n[/math], приходим к глобальному способу. Тогда интерполяционный многочлен Ньютона n-й степени имеет вид

[math]begin{aligned}N_n(x)= f_0&+ f(x_0,x_1)(x-x_0)+ f(x_0,x_1,x_2) (x-x_0)(x-x_1)+ ldots+\ &+,f(x_0,x_1,ldots,x_n)(x-x_0)(x-x_1)cdot ldots+(x-x_{n-1}).end{aligned}[/math]

(4.32)

Замечания

1. Согласно теореме 4.1 многочлен Ньютона (4.32) является тождественным многочлену [math]textstyle{widetilde{f}_n(x,overline{a})= sumlimits_{j=0}^{n} a_jx^j}[/math] с коэффициентами, получаемыми из системы (4.10), либо многочлену Лагранжа, т.е. [math]L_n(x)= N_n(x)= widetilde{f}_n(x,overline{a})[/math], если узлы интерполяции и интерполируемая функция одинаковы.

2. Интерполяционный многочлен Ньютона (4.29) или (4.32) (так же, как и многочлен Ньютона, выражаемый ниже через конечные разности) записан не через значения функции, как это имеет место для многочлена Лагранжа, а через разделенные разности. Поэтому при изменении степени [math]k[/math] в процессе интерполирования у многочлена Ньютона [math]N_k(x)[/math] требуется только добавить или отбросить соответствующее число слагаемых. Это иногда упрощает алгоритм интерполирования.

3. При интерполяции на основе (4.29) или (4.32) узлы интерполяции [math]x_i,x_{i+1}, ldots, x_{i+k}[/math] или [math]x_0,x_1,ldots,x_n[/math] определяющие шаблоны интерполяции, целесообразно выбирать так, чтобы точка [math]x_{ast}[/math] была расположена возможно ближе к середине отрезка [math][x_i,x_{i+k}][/math] или [math][x_0,x_n][/math].

4. Остаточное слагаемое многочлена (4.32) совпадает с остаточным слагаемым многочлена Лагранжа, и оценки (4.21), (4.22), справедливые для точки [math]x_{ast}[/math] и всего отрезка [math][x_0,x_n][/math], сохраняются.

Интерполяционные многочлены Ньютона для равномерной сетки

Сначала рассмотрим решение задачи кусочной интерполяции (применение кусочного способа). Если функция [math]y_i= f(x_i),~ i=overline{0,n}[/math] задана на равномерной сетке [math]Omega_n[/math], характеризующейся [math]h_{i+1}=text{const}[/math] для всех [math]i[/math], то многочлен (4.29), соответствующий шаблону [math](x_i,x_{i+1},ldots, x_{i+k})[/math], путем подстановки в него вместо разделенных разностей их выражений через конечные разности, согласно (4.28), преобразуется к виду

[math]N_{k}^{(I)}(q)= f_i+ frac{Delta f_i}{1!},q+ frac{Delta^2 f_i}{2!},q(q-1)+ ldots+ frac{Delta^k f_i}{k!}q(q-1)cdot ldotscdot (q-k+1),[/math]

(4.33)

где [math]q=frac{x-x_i}{h}[/math] — фаза интерполяции, определенная относительно точки [math]x_i[/math]; [math]Delta^jf_i~(j=1,2,ldots,k)[/math] — конечные разности.

В соответствии с формулой для фазы интерполяции [math]q[/math] точка начала ее отсчета расположена в узле [math]x_i[/math] и входящие в (4.33) конечные разности относятся к этой же точке [math]x_i[/math]. В связи с этим (4.33) удобно применять в начале выделенного шаблона [math](x_i,x_{i+1},ldots, x_{i+k})[/math], когда [math]q>0[/math]. Если [math]x_i=x_0[/math], а [math]x_{ast}<x_0[/math], то этот же многочлен используется и для экстраполяции левее точки [math]x_0~(q<0)[/math]. Поэтому многочлен [math]N_k^{(I)}(q)[/math] называется интерполяционным многочленом для интерполяции вперед (в начале таблицы ) или для экстраполяции назад. В (4.33) этот многочлен обозначен цифрой I, указанной в скобках вверху. Для определенности назовем его первым интерполяционным многочленом Ньютона.

Полагая [math]i=0,~k=n[/math], получаем решение задачи глобальной интерполяции на всем отрезке [math][x_0,x_n]colon[/math]

[math]N_{n}^{(I)}(q)= f_0+ frac{Delta f_0}{1!},q+ frac{Delta^2 f_0}{2!},q(q-1)+ ldots+ frac{Delta^n f_i}{n!}q(q-1)cdot ldotscdot (q-n+1),[/math]

(4.34)

где [math]q=frac{x-x_0}{h}[/math]. Остаточное слагаемое этого многочлена имеет вид

[math]R_n(q)= h^{n+1}cdot frac{q(q-1)cdotldotscdot(q-n)}{(n+1)!}cdot f^{(n+1)}(xi),[/math]

где [math]xiin(a,b)[/math] — некоторое промежуточное значение между узлами [math]x_0,x_1,ldots,x_n[/math] и точкой [math]x[/math].

Если фазу интерполяции определить относительно [math]x_{i+k}[/math] некоторой конечной точки шаблона [math](x_i,x_{i+1},ldots, x_{i+k})[/math], то есть [math]widehat{q}= frac{x-x_{i+k}}{h}[/math], то вместо (4.33) получается второй интерполяционный многочлен Ньютона:

[math]N_{k}^{(II)}(widehat{q})= f_{i+k}+ frac{Delta f_{i+k-1}}{1!},widehat{q}+ frac{Delta^2 f_{i+k-2}}{2!},widehat{q}(widehat{q}+1)+ ldots+ frac{Delta^k f_i}{k!}widehat{q}(widehat{q}+1)cdot ldotscdot (widehat{q}+k-1),[/math]

(4.35)

Данный многочлен удобно применять в конце выделенного шаблона или всей таблицы [math]y_i= f(x_i),~ i=overline{0,n}[/math].

Если [math]i+k=n[/math], а [math]x_{ast}>x_n[/math], то (4.35) используется для экстраполяции правее точки [math]x_n~(widehat{h}>0)[/math]. Поэтому многочлен [math]N_{k}^{(II)}(widehat{q})[/math] называется многочленом для интерполяции назад (в конце таблицы) или экстраполяции вперед.

Полагая [math]i=0,~ i+k=n[/math], получаем решение задачи глобальной интерполяции — второй интерполяционный многочлен Ньютона n-й степени (где [math]widehat{q}= frac{x-x_n}{h}[/math]):

[math]N_{n}^{(II)}(widehat{q})= f_{n}+ frac{Delta f_{n-1}}{1!},widehat{q}+ frac{Delta^2 f_{n-2}}{2!},widehat{q}(widehat{q}+1)+ ldots+ frac{Delta^n f_0}{n!}widehat{q}(widehat{q}+1)cdot ldotscdot (widehat{q}+n-1),[/math]

(4.36)

Остаточное слагаемое многочлена (4.36) имеет вид

[math]R_n(widehat{q})= h^{n+1}cdot frac{widehat{q}(widehat{q}+1)cdot ldotscdot (widehat{q}+n)}{(n+1)!}cdot f^{(n+1)}(xi)[/math], где [math]xiin(a,b)[/math].

Схема выбора узлов интерполяции при изменении степеней интерполяционных многочленов [math]N_{k}^{(I)}(q),~ N_{k}^{(II)}(widehat{q})~(k=1,2,3)[/math] в одном алгоритме показана на рис. 4.5. Если точка [math]x_{ast}[/math] находится в начале отрезка [math][x_0,x_n][/math], например на частичном отрезке [math][x_0,x_1][/math], то применяется первый интерполяционный многочлен Ньютона, а если в конце, например на частичном отрезке [math][x_{n-1},x_n][/math], то — второй (рис. 4.5,а). Если точка [math]x_{ast}[/math] находится вдали от концов отрезка [math][x_0,x_n][/math], то может применяться как первый интерполяционный многочлен, так и второй (рис. 4.5,б).

Замечания

1. Из формул интерполяционных многочленов [math]N_{k}^{(I)}(q)[/math] и [math]N_{k}^{(II)}(widehat{q})[/math] видно, что повышение их степеней в процессе реализации алгоритма не требует пересчета предыдущих слагаемых, входящих в многочлены с меньшими степенями.

2. Для гладких функций при повышении порядка конечных разностей справедливо свойство: [math]Delta^kf_ito0[/math] при [math]ktoinfty[/math], и поэтому, как только очередное слагаемое в рассматриваемых многочленах становится меньше требуемой точности интерполяции, увеличение степени [math]N_{k}^{(I)}(q),~ N_{k}^{(II)}(widehat{q})[/math] следует прекратить. Это замечание справедливо также и для [math]N_k(x)[/math].

▼ Пример 4.4-4.6

Пример 4.4. Построить многочлен Ньютона третьей степени для сеточной функции, заданной табл. 4.7. Вычислить значение функции в точке [math]x_{ast}=2,!5[/math]. Решить задачу интерполяции при включении одного дополнительного значения сеточной функции: [math]f(1)=5[/math].

[math]begin{array}{|c|c|c|c|c|} multicolumn{5}{r}{mathit{Table~4.7}}\hline i& 0& 1& 2& 3\hline x_i& 2&3& 4& 5\hline f(x_i)=f_i& 7&5& 8& 7\hline end{array}[/math]

Решение

Построим многочлен Ньютона (4.32), справедливый для произвольного расположения узлов. Для этого составим табл. 4.8, аналогичную табл. 4.5:

[math]begin{gathered}f(x_0,x_1)= frac{5-7}{3-2}=-2;qquad f(x_1,x_2)= frac{8-5}{4-3}=3;qquad f(x_2,x_3)= frac{7-8}{5-4}=-1;\[2pt] f(x_0,x_1,x_2)= frac{3-(-2)}{4-2}=frac{5}{2};qquad f(x_1,x_2,x_3)= frac{-1-3}{5-3}=-2;qquad f(x_0,x_1,x_2,x_3)= frac{-2-frac{5}{2}}{5-2}=-frac{3}{2},.end{gathered}[/math]

По формуле (4.32) для [math]n=3[/math] имеем

[math]begin{gathered}N_3(x)= f(x_0)+ (x-x_0)f(x_0,x_1)+ (x-x_0)(x-x_1)f(x_0,x_1,x_2)+\ +(x-x_0)(x-x_1)(x-x_2)f(x_0,x_1,x_2,x_3)= 7+(x-2)cdot(-2)+ (x-2)(x-3)cdot frac{5}{2}+\ (x-2)(x-3)(x-4)! left(-frac{3}{2}right)=-frac{3}{2},x^3+ 16x^2-frac{107}{2},x+62. end{gathered}[/math]

Поскольку в данной задаче заданы равностоящие узлы, воспользуемся также формулой (4.34) для первого интерполяционного многочлена Ньютона:

[math]N_{3}^{(I)}(q)= f(x_0)+ frac{Delta f_0}{1!},q+ frac{Delta^2f_0}{2!},q(q-1)+ frac{Delta^3f_0}{3!},q(q-1)(q-2),[/math]

где [math]q=frac{x-x_0}{h}= frac{x-2}{1}=x-2[/math]. Составим табл. 4.9, аналогичную табл. 4.6. Имеем:.

[math]begin{gathered}Delta f_0= 5-7=-2;quad Delta f_1=8-5=3;quad Delta f_2=7-8=-1;quad Delta^2f_0=3-(-2)=5;\[2pt] Delta^2f_1=-1-3=-4;qquad Delta^3f_0=-4-5=-9. end{gathered}[/math]

Поэтому

[math]begin{aligned}N_{3}^{(I)}(q)&= f(x_0)+ frac{(-2)}{1!},q+ frac{5}{2!},q(q-1)+ frac{(-9)}{3!},q(q-1)(q-2)= left.{left(-frac{3}{2},q^3+7q^2-frac{15}{2},q+7right)}right|_{q=frac{x-2}{1}}\ &=-frac{3}{2}cdot (x^3-6x^2+12x-8)+ 7cdot (x^2-4x+4)-frac{15}{2}cdot (x-2)+7=-frac{3}{2},x^3+ 16x^2-frac{107}{2},x+62. end{aligned}[/math]

Запишем также второй интерполяционный многочлен Ньютона (4.36):

[math]begin{aligned}N_{3}^{(II)}(q)&= left.{left(f_3+ frac{Delta f_2}{1!},widehat{q}+ frac{Delta^2f_1}{2!}, widehat{q} (widehat{q}+1)+ frac{Delta^3f_0}{3!}, widehat{q}(widehat{q}+1)(widehat{q}+2)right)!}right|_{widehat{q}=frac{x-x_0}{1}}=\ &=7+ frac{-1}{1!},widehat{q}+ frac{-4}{2!}, widehat{q} (widehat{q}+1)+ frac{-9}{3!}, widehat{q}(widehat{q}+1)(widehat{q}+2)=\ &= left.{left(-frac{3}{2},widehat{q}^3-frac{13}{2},widehat{q}^2-6widehat{q}+7right)!}right|_{widehat{q}=x-5}=-frac{3}{2},x^3+ 16x^2-frac{107}{2},x+62. end{aligned}[/math]

Сравнивая с результатом примера 4.1, можно заключить, что [math]L_3(x)= N_3^{(I)}(x)= N_3^{(II)}(x)[/math]. Это еще раз подтверждает единственность решения задачи интерполяции в классе многочленов, удовлетворяющих условиям теоремы 4.1.

Предположим, что в ходе некоторого эксперимента получен новый результат [math]f(1)=5[/math], дополняющий заданную сеточную функцию. Тогда для решения задачи интерполяции с помощью многочлена Ньютона можно использовать уже полученный результат. Для этого дополним табл. 4.8. Новый узел и соответствующее значение функции поместим в конце табл. 4.10.

По формуле (4.34) имеем

[math]begin{aligned}N_4(x)&= N_3(x)+(x-x_0)(x-x_1)(x-x_3)(x-x_4)f(x_0,x_1,x_2, x_3,x_4)=\[2pt] &=-frac{3}{2},x^3+ 16x^2-frac{107}{2},x+62+ (x-2)(x-3)(x-4)(x-5)! left(-frac{3}{4}right)=\[2pt] &=-frac{3}{2},x^3+ 16x^2-frac{107}{2},x+62-frac{3}{4}cdot (x^4-14x^3+71x^2-154x+120)=\[2pt] &=-frac{3}{4},x^4+ 9x^3-frac{149}{4},x^2+62x-28. end{aligned}[/math]

Легко проверить, что новый узел можно добавить и в начале таблицы, т.е. будет найден тот же многочлен [math]N_4(x)[/math]. На рис. 4.6 изображены полученные интерполяционные многочлены.

Вычислим значение функции в точке [math]x_{ast}=2,!5colon[/math]

[math]N_3^{(I)}(x_{ast})= N_3^{(II)}(x_{ast})= 4,!8125;qquad N_4^{(I)}(x_{ast})= 5,!5156.[/math]

Пример 4.5. Для сеточной функции из примера 4.3 найти линейный и параболический многочлены Ньютона и на их основе подсчитать значение функции в точке [math]x_{ast}=2,!5[/math]. Оценить погрешность интерполяции.

Решение

Так как точка [math]x_{ast}=2,!5[/math] находится вблизи узла [math]x_i=x_2[/math] и шаблон [math](x_2,x_3,x_4)[/math] имеет неравномерность по шагу, то для выполнения интерполяции выбираем многочлен Ньютона (4.29).

1. Определим разделенные разности [math]f(x_i,x_{i+1}),~ f(x_i,x_{i+1},x_{i+2})[/math] (табл. 4.11)

Здесь [math]f(x_1,x_2)= frac{27-1}{3-1}=13,~ f(x_3,x_4)= frac{64-27}{4-3}=37,~ f(x_2,x_3,x_4)= frac{37-13}{4-1}=8.[/math]. Определим [math]N_1(x),~N_2(x)colon[/math]

[math]begin{aligned}N_1(x)&= f_2+f(x_2,x_3)cdot (x-x_2)= 1+13cdot(x-1)=13x-12.\[2pt] N_2(x)&= f_2+ f(x_2,x_3)(x-x_2)+ f(x_2,x_3,x_4)(x-x_2)(x-x_3)=\ &=13x-12+8(x-1)(x-3)= 8x^2-19x+12. end{aligned}[/math]

Сравнивая [math]N_2(x)[/math] с соответствующим многочленом Лагранжа L_2(x) (см. пример 4.3), видим, что они совпадают. Это подтверждает единственность решения задачи о построении интерполяционного многочлена.

2. Вычислим значения [math]Bigl.{N_1(x)}Bigr|_{x=x_{ast}},~ Bigl.{N_2(x)}Bigr|_{x= x_{ast}}colon[/math]

[math]N_1(2,!5)= 13cdot2,!5-12=20,qquad N_2(2,!5)= 8cdot2,!5^2-19cdot 2,!5+12=14,!5.[/math]

Итак, линейная интерполяция приводит к отличию от точного значения на [math]left|frac{20-15,!625}{15,!625}right|cdot 100%=28%[/math], а квадратичная — на [math]left|frac{14,!5-15,!625}{15,!625}right|cdot 100%=7%[/math].

Таким образом, найдена фактическая погрешность, которая для квадратичной интерполяции получилась в 4 раза меньше.

3. Найдем априорную оценку погрешности линейной интерполяции на отрезке [math][x_2,x_3][/math] по формуле (4.23). Для этого необходимо сначала оценить производную [math]M_{2,i}[/math]. При этом, как и в реальных задачах, будем определять [math]f»(x)[/math] численным дифференцированием. С этой целью используем многочлен [math]N_2(x)colon, N’_2(x)=16x-9,~ N»_2(x)=16[/math]. Приближенно можно принять, что [math]M_{2,i}=16[/math]. Тогда получаем: [math]|R_1(x)|leqslant frac{2^2}{8}cdot16=8[/math] или [math]frac{8}{15,!625}cdot 100%= 51,!2%[/math]. Заметим, что фактическая погрешность получилась примерно в два раза меньше априорной.

Пример 4.6. Для сеточной функции, заданной в примере 4.4, построить интерполяционные многочлены первой и второй степени для нахождения значений в точках

[math]x_{ast 1}=2,!5;qquad x_{ast 2}=4,!5;qquad x_{ast 3}=3,!5;qquad x_{ast 4}=5,!5.[/math]

Решение

Для подсчета значения функции в точке [math]x_{ast 1}=2,!5[/math] (согласно рис. 4.5,а) выбираем шаблон [math](x_0,x_1,x_2)= (2;3;4)[/math] для параболической интерполяции [math](i=0,~ k=2)[/math] и шаблон [math](x_0,x_1)= (2;3)[/math] Для линейной интерполяции [math](i=0,~ k=1)[/math]. Тогда по формуле (4.34) имеем

[math]begin{gathered}N_1^{(I)}(q)= f_0+ frac{Delta f_0}{1!},q= 7+frac{-2}{1},q=7-2q,\[2pt] N_2^{(II)}(q)= N_1^{(I)}(q)+ frac{Delta^2f_0}{2!},q(q-1)= 7-2q+ frac{5}{2},q(q-1)= frac{5}{2},q^2-frac{9}{2},q+7. end{gathered}[/math]

где [math]q=frac{x-2}{1}[/math] (см. табл. 4.9). В результате подстановки получаем

[math]N_1^{(I)}(x)=-2x+11;qquad N_2^{(I)}(x)= frac{5}{2},x^2-frac{29}{2},x+26.[/math]

В заданной точке [math]N_1^{(I)}(2,!5)=6;~ N_2^{(I)}(2,!5)=5,!375[/math].

Для подсчета значения в точке [math]x_{ast 2}=4,!5[/math] выбираем шаблон [math](3;4;5)= (x_1,x_2,x_3)[/math] для параболической интерполяции [math](i=1,~ k=2)[/math] и шаблон [math](4;5)=(x_2,x_3)[/math] для линейной интерполяции [math](i=2,~ k=1)[/math]. Тогда по формуле (4.35) имеем

[math]begin{gathered}N_1^{(II)}(widehat{q})= f_3+ frac{Delta f_2}{1!},widehat{q}= 7+frac{-1}{1},widehat{q}= 7-widehat{q},\[2pt] N_2^{(II)}(widehat{q})= N_1^{(II)}(widehat{q})+ frac{Delta^2f_1}{2!}, widehat{q}(widehat{q}+1)= 7-widehat{q}+ frac{-4}{2},widehat{q} (widehat{q}+1)=-2widehat{q}^2-3widehat{q}+7, end{gathered}[/math]

где [math]widehat{q}=frac{x-x_3}{h}=x-5[/math]. В результате подстановки получаем

[math]N_1^{(II)}(x)=-x+12;qquad N_2^{(II)}(x)=-2x^2+17x-28.[/math]

В заданной точке [math]N_1^{(II)}(4,!5)=7,!5;~ N_2^{(II)}(4,!5)=8[/math].

Многочлен [math]N_1^{(II)}(x)[/math] можно использовать для экстраполяции, т.е. для подсчета функции в точке [math]x_{ast 4}=5,!5colon, N_2^{(II)}(5,!5)=5[/math].

Чтобы подсчитать значения в точке [math]x_{ast3}=3,!5[/math], сначала выберем шаблон [math](3;4;5)=(x_1,x_2,x_3)[/math] для параболической интерполяции [math](i=1,~k=2)[/math] и шаблон [math](3;4)=(x_1,x_2)[/math] для линейной интерполяции [math](i=1,~ k=1)[/math]. Согласно рис. 4.5,б, применим формулу (4.33):

[math]begin{gathered}N_1^{(I)}(q)= f_1+frac{Delta f_1}{1!},q= 5+frac{3}{1},q=5+3q,\[2pt] N_2^{(I)}(q)= N_1^{(I)}(q)+ frac{Delta^2f_1}{2!},q(q-1)= 5+3q+frac{-4}{2},q(q-1)=-2q^2+5q+5, end{gathered}[/math]

где [math]q=frac{x-x_1}{h}=x-3[/math]. После подстановки получаем

[math]N_1^{(I)}(x)= 3x-4;~~ N_2^{(I)}(x)=-2x^2+17x-28[/math] и [math]N_1^{(I)}(3,!5)= 6,!5;~~ N_2^{(I)}(3,!5)=7[/math].

Подсчитать значения в точке [math]x_{ast3}=3,!5[/math] (согласно рис. 4.5,б) можно, использовав шаблон [math](2;3;4)= (x_0,x_1,x_2)[/math] для параболической интерполяции [math](i=0,~ k=2)[/math] и шаблон [math](3;4)=(x_1,x_2)[/math] для линейной интерполяции [math](i=1,~ k=1)[/math]. По формуле (4.35) получаем

[math]begin{gathered}N_1^{(II)}(widehat{q})= f_2+ frac{Delta f_1}{1!},widehat{q}= 8+frac{3}{1},widehat{q},\[2pt] N_2^{(II)}(widehat{q})= N_1^{(II)}(widehat{q})+ frac{Delta^2 f_0}{2!},widehat{q}= 8+frac{3}{1},widehat{q}(widehat{q}+1)= frac{5}{2},widehat{q}^2+ frac{11}{2},widehat{q}+8. end{gathered}[/math]

где [math]widehat{q}=frac{x-x_2}{h}=x-4[/math]. После подстановки

[math]N_1^{(II)}(x)=3x-4;~~ N_2^{(II)}(x)= frac{5}{2},x^2-frac{29}{2},x+26[/math] и [math]N_1^{(II)}(3,!5)=6,!5;~~ N_2^{(II)}(3,!5)=5,!875[/math].

Математический форум (помощь с решением задач, обсуждение вопросов по математике).

Кнопка "Поделиться"

Если заметили ошибку, опечатку или есть предложения, напишите в комментариях.

Источник

Этот калькулятор может пригодиться при решении задач на интерполяцию полиномом Лагранжа. В таких задачах обычно требуется интерполировать значение неизвестной функции, соответствующее некоторому значению x, использую формулу интерполяционного многочлена Лагранжа, полученную из известного набора точек со значениями неизвестной функции (x, f(x)).

Калькулятор ниже обладает следующими функциями:

Он находит формулу полинома Лагранжа для заданного набора точек.
Он отображает пошаговый вывод формулы.
Он вычисляет значения интерполяционного многочлена Лагранжа для заданных точек (интерполирует функцию полиномом Лагранжа в заданных точках интерполяции)
Он отображает набор точек, значения в точках интерполяции, полином Лагранжа и все базисные полиномы на графике.

Как пользоваться

Сначала вводите набор точек — одна точка на строку в форме x f(x), значения разделены пробелом. Если вы хотите получить интерполяцию, вводите значения точек интерполяции в следующее поле в виде значений x, разделенных пробелом.

По умолчанию, калькулятор отображает формулу многочлена и его значения в точках интерполяции. Если нужно пошаговое решение, включите опцию «Показать пошаговое решение». Также можно отключить отображение базисных полиномов.

Теория и формулы, как обычно, описаны под калькулятором.

Интерполяционный многочлен Лагранжа (полином Лагранжа)

Набор точек, одна точка на строку, значения разделяются пробелом

Точность вычисления

Знаков после запятой: 2

Файл очень большой, при загрузке и создании может наблюдаться торможение браузера.

Показать решение по шагам

Файл очень большой, при загрузке и создании может наблюдаться торможение браузера.

Показать базисные полиномы

Полином Лагранжа

Файл очень большой, при загрузке и создании может наблюдаться торможение браузера.

Интерполяционный многочлен Лагранжа

Предположим, что у нас есть набор значений, соответствующих неизвестной функции, при этом все x различны:

$(x_{0},y_{0}),ldots ,(x_{j},y_{j}),ldots ,(x_{k},y_{k})$

Сконструируем следующий многочлен (называемые многочленом Лагранжа):

$L(x):=sum _{j=0}^{k}y_{j}ell _{j}(x)$

где $ell _{j}(x)$ — базисный полином Лагранда

$ell _{j}(x):=prod _{begin{smallmatrix}0leq mleq k\mneq jend{smallmatrix}}{frac {x-x_{m}}{x_{j}-x_{m}}}={frac {(x-x_{0})}{(x_{j}-x_{0})}}cdots {frac {(x-x_{j-1})}{(x_{j}-x_{j-1})}}{frac {(x-x_{j+1})}{(x_{j}-x_{j+1})}}cdots {frac {(x-x_{k})}{(x_{j}-x_{k})}}$

Если посмотреть на формулу базисного полинома для любого j, то видно что для всех точек i не равных j, значение этого полинома обращается в ноль, а в самой точке j значение этого полинома j равно единице. Таким образом,

$y_{j}ell _{j}(x_{j})=y_{j} cdot 1=y_{j}$

$L(x_{j})=y_{j}+0+0+dots +0=y_{j}$

что означает, что полином Лагранжа точно интерполирует значение функции в заданных точках.

Стоит заметить, что формула интерполяционного многочлена Лагранжа подвержена воздействию так называемого феномена Рунге. Феномен Рунге связан к увеличением колебаний полинома на краях интервала при использовании полиномов высоких степеней на равноудаленных друг от друга точках. Таким образом, наличие большого количества точек далеко не всегда приводит к улучшению точности интерполяции.

Однако также стоит заметить, что в отличие от некоторых других формул интерполяции, формула Лагранжа не требует того, чтобы точки в наборе были равноудалены друг от друга. Это используется в некоторых способах борьбы с феноменом Рунге, например, при использовании в качестве точек интерполяции узлов Чебышева.

Источник

Интерполяцио́нный многочле́н Лагра́нжа — многочлен минимальной степени, принимающий данные значения в данном наборе точек. Для пар чисел ${displaystyle (x_{0},y_{0}),(x_{1},y_{1})dots (x_{n},y_{n})}$ , где все ${displaystyle x_{i}}$ различны, существует единственный многочлен степени не более , для которого ${displaystyle L(x_{i})=y_{i}}$ .

В простейшем случае это линейный многочлен, график которого — прямая, проходящая через две заданные точки.

Определение

Файл:Lagrangepolys.png

Этот пример показывает интерполяционный многочлен Лагранжа для четырёх точек (-9,5), (-4,2), (-1,-2) и (7,9), а также полиномы y_j l_j(x), каждый из которых проходит через одну из выделенных точек, и принимает нулевое значение в остальных x_i

Лагранж предложил способ вычисления таких многочленов:

${displaystyle L(x)=sum _{j=0}^{n}y_{j}l_{j}(x)}$

где базисные полиномы определяются по формуле:

${displaystyle l_{j}(x)=prod _{i=0,jneq i}^{n}{frac {x-x_{i}}{x_{j}-x_{i}}}={frac {x-x_{0}}{x_{j}-x_{0}}}cdots {frac {x-x_{j-1}}{x_{j}-x_{j-1}}}{frac {x-x_{j+1}}{x_{j}-x_{j+1}}}cdots {frac {x-x_{n}}{x_{j}-x_{n}}},!}$

Легко видеть что ${displaystyle l_{j}(x)}$ обладают такими свойствами:

Отсюда следует, что , как линейная комбинация ${displaystyle l_{j}(x)}$ , может иметь степень не больше , и ${displaystyle L(x_{j})=y_{j}}$ , Q.E.D.

Применения

Полиномы Лагранжа используются для интерполяции, а также для
численного интегрирования.

Пусть для функции известны значения ${displaystyle y_{j}=f(x_{j})}$ в некоторых точках. Тогда мы можем интерполировать эту функцию как

${displaystyle f(x)approx sum _{j=0}^{n}f(x_{j})l_{j}(x)}$

В частности,

${displaystyle int limits _{a}^{b}f(x)dxapprox sum _{j=0}^{n}f(x_{j})int limits _{a}^{b}l_{j}(x)dx}$

Значения интегралов от ${displaystyle l_{j}}$ не зависят от , и их можно вычислить заранее, зная последовательность ${displaystyle x_{i}}$ .

Для случая равномерного распределения по отрезку узлов интерполяции

В указанном случае можно выразить ${displaystyle x_{i}}$ через расстояние между узлами интерполяции h и начальную точку ${displaystyle x_{0}}$ :

${displaystyle x_{j}equiv {x_{0}+jh}}$

и, следовательно,

${displaystyle {x_{i}-x_{j}}equiv (i-j)h}$

Подставив эти выражения в формулу полинома и вынеся h за знаки перемножения в числителе и знаменателе, получим

${displaystyle l_{i}(x)={prod _{j=0,,ineq j}^{n}{(x-x_{j}) over (x_{i}-x_{j})}}={prod limits _{j=0,,ineq j}^{n}(x-x_{0}-jh) over h^{n-1}prod limits _{j=0,,ineq j}^{n}(i-j)}}$

Теперь можно ввести замену переменной

${displaystyle y={{x-x_{0}} over h},!}$

и получить полином от XY, который строится с использованием только целочисленной арифметики. Недостатком данного подхода является факториальная сложность числителя и знаменателя, что требует использования алгоритмов с многобайтным представлением чисел.

cs:Lagrangeova interpolace
nl:Lagrange-polynoom
sr:Лагранжов полином
uk:Многочлен Лагранжа

Источник

Вычисление коэффициентов полинома Лагранжа

Definition[edit]

Barycentric form[edit]

A perspective from linear algebra[edit]

Example[edit]

Notes[edit]

Remainder in Lagrange interpolation formula[edit]

Derivation[7][edit]

Derivatives[edit]

Finite fields[edit]

See also[edit]

References[edit]

External links[edit]

Методы функциональной интерполяции

Постановка задачи

Многочлен Лагранжа

Погрешность интерполяции многочленами Лагранжа

Сходимости функционального интерполяционного процесса для непрерывных функций

Линейная и параболическая интерполяция с помощью многочлена Лагранжа

Методика решения задачи линейной интерполяции

Методика решения задачи параболической интерполяции

Многочлены Ньютона

Интерполяционный многочлен Ньютона для неравномерной сетки

Интерполяционные многочлены Ньютона для равномерной сетки

Как пользоваться

Интерполяционный многочлен Лагранжа (полином Лагранжа)

Интерполяционный многочлен Лагранжа

Определение

Применения

Для случая равномерного распределения по отрезку узлов интерполяции

Derivation^[7][edit]