Математический анализ в Maple 9



Интерполяция методом Ньютона

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

Параметрами процедуры RS() являются список А узловых точек и значений интерполируемой функции в этих точках, а также порядок N разделительной разности.

В начале процедуры инициализируются с нулевыми значениями переменные i и S. Далее следует условный оператор while. Проверяемым в этом операторе условием является неравенство i=N. Таким образом, команды из оператора цикла будут выполняться до тех пор, пока значение переменной i не превысит N. Изменение значения индексной переменной на единицу осуществляется последней командой в операторе цикла.



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

Теперь определим оператор, с помощью которого будем строить интерполяционный полином Ньютона.

Внимание!
Здесь еще раз хочется напомнить, что при интерполяции функции полиномом последний не зависит от способа его получения. В этом смысле, что полином Ньютона, что Лагран-жа — все едино. Полином один и тот же. Название характеризует только способ, которым этот полином был получен.

Процедура newton() имеет только один параметр, который представляет собой список с узловыми точками и значениями функции.

Локальная переменная L инициализируется со значением А[1][2]. Это значение интерполируемой функции в первой узловой точке. Далее переменной Ямх присваиваем в качестве значения длину списка А. С помощью оператора цикла, согласно формуле для интерполяционного полинома, формируем соответствующее выражение. В этом операторе, кроме прочего, вызывается и разработанная выше процедура определения разделительных разностей RS(). Наконец, собираем слагаемые при соответствующих степенях аргумента (команда collect(L,x)), и для этого выражения вычисляем функциональный оператор (процедура unapply()), который и возвращается в качестве результата. Теперь процедуру можно использовать.

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

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

Для построения интерполяционного полинома методом Ньютона можно было поступить несколько иначе, благо возможности Maple позволяют это сделать. Итак, будем искать интерполяционный полином в виде L(x) = ao+(x-xo)al+... + (x-xB)...(x-x,_,)a,. Коэффициенты а, нужно определить из тех условий, что в узловых точках полином принимает известные значения. Эти коэффициенты, по сути, являются разделительными разностями, о которых речь шла ранее. Однако здесь для нахождения коэффициентов воспользуемся уникальными возможностями Maple. Для этого создадим процедуру newton2(). Параметры процедуры определим несколько отличным от предыдущих случаев способом. Первые два параметра — списки X и Y с узловыми точками и значениями интерполируемой функции в этих точках. Третьим параметром является переменная интерполирования.

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

На заметку
Несовпадение индексов выше объясняется просто. Дело в том, что индексация коэффициентов начинается с нуля, а индексация элементов в списках X и У — с единицы.

Далее определяется функция L(), имеющая структуру интерполяционного полинома Ньютона. В этой функции присутствуют неизвестные пока что коэффициенты а[ j]. Чтобы найти эти коэффициенты, зададим уравнения, из которых они будут определяться. Уравнения будем записывать в качестве элементов множества EqSys, которое инициализируется как пустое. Уравнения определяются следующим образом. В операторе цикла с помощью индексной переменной 1 перебираются значения из списка X и подставляются в функцию L(), задающую полином. По определению в этих точках значения полинома должны совпадать с соответствующими элементами из списка Y. Полученные уравнения заносятся во множество EqSys.

Внимание!
Не следует забывать, что нижняя граница изменения индексной переменной равна О, в то время как индексация элементов списков начинается с 1.

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

На заметку
Для решения такой системы в среде, отличной от Maple, пришлось бы разрабатывать специальные подпрограммы.

Для того чтобы коэффициентам присвоить значения в соответствии с полученным решением, используется процедура assign(). Результат решения системы уравнений возвращается в виде множества равенств, определяющих искомые коэффициенты. Именно на это множество и ссылается переменная среды % в аргументе команды. После этого в выражении L(x) слагаемые группируются по степеням аргумента, и такое выражение возвращается как результат выполнения процедуры.

Внимание!
В отличие от предыдущих случаев, в том числе и процедуры newtonf), где результат процедуры задавался как действие, в процедуре newton2() результатом является выражение. Здесь переменная интерполирования указывается параметром процедуры.

Определяем теперь базовые списки с узловыми точками и значениями функции.

Получаем вполне ожидаемый результат.

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

Само собой разумеется, что интерполирование полиномами не является единственно возможным способом выполнения интерполяции. Рассмотрим иные методы.