3.3. Интерполяционный полином Ньютона

По данным табл. 3.1 построим интерполяционный полином степени в виде, предложенном Ньютоном

Равносильный вариант полинома можно записать при симметричной перенумерации узлов исходной таблицы узлов и значений функции

Коэффициенты полиномов (3.6) и (3.7) определяются из условий Лагранжа

Полагаем тогда в формуле (3.6) все слагаемые, кроме обращаются в нуль, следовательно,

Затем полагаем тогда по условию (3.8)

откуда находим коэффициент

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

При полином (3.6) принимает значение

из условия Лагранжа (3.8) определяем искомый коэффициент

где

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

Аналогичным образом при находим коэффициент полинома Ньютона

где

Для коэффициента методом математической индукции запишем следующее выражение:

Таблица 3.3

Полученные результаты сведем в табл. 3.3.

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

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

Заметим, что добавление новых узлов в табл. 3.3 не изменит уже вычисленных коэффициентов, таблица будет дополнена новыми строками и столбцами разделенных разностей.

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

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

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

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

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

Значения полинома Ньютона при циклическом изменении значений аргумента (переменная определяются с помощью отдельной подпрограммы, реализующей схему Горнера (строки 300-390 программы подпрограммы программ Наряду со значениями полинома вычисляются первая и вторая производные при соответствующих аргументах (переменные Операторы для производных записаны путем дифференцирования по переменной правой части оператора для переменной Таким же способом записаны операторы инициализации вне цикла для выполнения схемы Горнера

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

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

Погрешность полиномиальной аппроксимации функции

определяется соотношением [1]

где

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

Если узлы расположены равномерно с шагом то наименьшая погрешность аппроксимации будет в интервалах» примыкающих к центральному узлу, за счет минимальной величины произведения в правой части оценки (3.15). Особенно резко увеличивается погрешность при экстраполяции. В центральном интервале (при четном количестве узлов) получена следующая оценка погрешности [1]:

(см. скан)

(см. скан)

(см. скан)