ЦЕЛЬ

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

ПОСТАНОВКА ЗАДАЧИ

Разработать алгоритм построения полинома Ньютона и полинома Лагранжа для заданной функции: \[ \ln{x} + \sqrt{1+x}\] В качестве заданных узлов интерполяции возьмем точки, принадлежащие интервалу [1,16] с шагом h = 2.5. Таким образом, у нас заданы 7 точек:

x y
1.0 1.41
3.5 3.37
6.0 4.44
8.5 5.22
11.0 5.86
13.5 6.41
16.0 6.90

Примечание: вектор значений \(х\) и вектор значений \(у\) хранятся в переменных \(х\) и \(у\) соответственно.

ОПИСАНИЕ МЕТОДА

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

Первая интерполяционная формула Ньютона

Метод Ньютона заключается в том, чтобы построить полином вида: \[ P_n(x) = y_0 + q \Delta y_0 + \frac{q(q-1)}{2!} \Delta^2 y_0 + \ldots + \frac{q(q-1)\ldots(q-n+1)}{n!} \Delta^n y_0\] \(\Delta^ky_0\) — конечные разности для первого узла
\(q=\frac{x-x_0}h\), где \(x_0\) - значение функции в первом узле, \(h\) - шаг интерполяции

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

И проверим результаты ее работы. Найдем вектор первых конечных разностей:

x y d1
1.0 1.4142 1.9599
3.5 3.3741 1.0634
6.0 4.4375 0.7848
8.5 5.2223 0.6397
11.0 5.8620 0.5486
13.5 6.4106 0.4851
16.0 6.8957 NA

Аналогично найдем другие конечные разности и составим таблицу:

Node # x y d1 d2 d3 d4 d5 d6
1 1 1.4142 1.9599 -0.8964 0.6178 -0.4842 0.4044 -0.3509
2 3.5 3.3741 1.0634 -0.2787 0.1336 -0.0797 0.0535
3 6 4.4375 0.7848 -0.145 0.0539 -0.0262
4 8.5 5.2223 0.6397 -0.0911 0.0277
5 11 5.862 0.5486 -0.0635
6 13.5 6.4106 0.4851
7 16 6.8957

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

y d1 d2 d3 d4 d5 d6
1.414214 1.95987 -0.8964423 0.6177772 -0.4841508 0.4044186 -0.3508977

Когда мы подставим в формулу Ньютона полученных конечные разности, то получим формулу вида:

\[\begin{eqnarray*} P_n(x) & = & 1.4142136 + (1.9598697)q + (-0.8964423) \frac{q(q-1)}{2!} \\ & + & (0.6177772) \frac{q(q-1)(q-2)}{3!} + (-0.4841508) \frac{q(q-1)(q-2)(q-3)}{4!} \\ & + & (0.4044186) \frac{q(q-1)(q-2)(q-3)(q-4)}{5!} \\ & + & (-0.3508977) \frac{q(q-1)(q-2)(q-3)(q-4)(q-5)}{6!} \end{eqnarray*}\] Упростим это выражение с помощью Wolphram Simplify и получим:

\[P_6(q) = 1.4142+2.87445 q-1.28106 q^2+0.451623 q^3-0.0953007 q^4+0.0106804 q^5-0.000487361 q^6\] Так как мы интерполировали относительно первого узла, где \(х = 1\), то в формуле для q у нас \(x_0=1\). Из условия также имеем \(h = 2.5\) То есть нам остается только подставить: \[q =\frac{x-1}{2.5}\] Подставим: \[\begin{eqnarray*} P_6(x) & = & 1.4142+2.87445 \frac{x-1}{2.5}-1.28106 (\frac{x-1}{2.5})^2\\ &+&0.451623 (\frac{x-1}{2.5})^3-0.0953007 (\frac{x-1}{2.5})^4\\ &+&0.0106804 (\frac{x-1}{2.5})^5-0.000487361 (\frac{x-1}{2.5})^6 \end{eqnarray*}\] Аналогично, упростим данное выражение с помощью Wolphram Simplify и получим:

\[P_6(x)=0.0279955+1.65675 x-0.307443 x^2+0.0397963 x^3-0.00301648 x^4+0.000121345 x^5-1.99623*10^{-6} x^6\] Теперь проверим, что мы действительно нашли то, что нужно и полученный полином в узлах совпадает с исходной функцией, то есть, что: \[P_6(x_i) = y_i\]

x y P(x)
1 1.414214 1.414202
3.5 3.374083 3.374113
6 4.437511 4.437633
8.5 5.222273 5.222567
11 5.861997 5.862524
13.5 6.410576 6.411319
16 6.895694 6.896472

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

Теперь построим графики обеих функций.