Аппроксимация табличной функции методом Ньютона

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

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

				...
				...

Вводим некоторое , .

Необходимо найти приближённое значение этой функции в некоторой точке x.

Особенность метода – может применяться только для равноотстоящих узлов:

x_i+1 – x_i=h-const (шаг интерполяции).

Если z ближе к x₀ – I многочлен Ньютона, а если z ближе к x_n – II многочлен Ньютона.

I многочлен Ньютона:

II многочлен Ньютона:

Cтроим таблицу:

			...
			...
			...
			...
...	...	...	...	...
			...
			...

Исходные данные:

Функция для проверки: y=tg(cos(x)).

Алгоритм решения:

1) Создать функцию для вычисления факториала:

fakt(int n)

2) Входные данные:

float X[N] = {..,..,..,..,..,..,..,.., };

float Y[N] = {..,..,..,..,..,..,..,.., };

3) Рассчитать h:

h=x₀ – x₁

4) Формирование матрицы конечных разностей:

4.1. заполняем 1 столбец

н. ц. i = 0, n-1

a_i0=y_i

к. ц. i

4.2. заполняем остальные столбцы

н. ц. j=1, n-1

н. ц. i=0, n-j

a_ij = a_i+1j-1-a_ij-1

к. ц. i

к. ц. j

5) Ввод z.

6) Начальные установки:

N=0

k=1

7) Проверка условия:

если z< x[N/2]

7.1. t=(z – x[0])/h

7.2. н. ц. j=0, N-1

N=(N+k*a_0j)/fakt(j)

k=k*(t-j)

к. ц. j

иначе

7.3. t=(z – x[N])/h

н. ц. j=0, N-1

N=(N+k*a_{n-j-1 j})/fakt(j)

k=k*(t+j)

к. ц. j

8) Вывод результата: N.

9) Проверка:

y=f(z)

y ≈ N

Текст программы:

#include< conio.h>

#include< math.h>

#include< stdio.h>

#include< iostream.h>

int n=8, pr=0;

float *X, *Y, h, **A;

void init(int q)

{

if(q)

{X=new float [n];

Y=new float [n];

A=new float *[n];

for(int i=0; i< n; i++) A[i]=new float [n-i];

}

else

{delete[] X;

delete[] Y;

for(int i=0; i< n; i++) delete[] A[i];

delete[] A; }

}

void example(void)

{

n=8;

init(1);

float x[8]={0.1, 0.3, 0.5, 0.7, 0.9, 1.1, 1.3, 1.5},

y[8]={1.540, 1.414, 1.204, 0.960, 0.716, 0.487, 0.274, 0.017};

for(int i=0; i< 8; i++) {X[i]=x[i]; Y[i]=y[i]; }

cout< < " \nVvod x[8]: ";

for(i=0; i< n; i++) printf(" %.1f ", X[i]);

cout< < " \n\nVvod y[8]: ";

for(i=0; i< n; i++) printf(" %.3f ", Y[i]);

cout< < endl;

pr=1;

}

float f(float x)

{

float t=1;

for(int i=1; i< =x; i++)

t*=i;

return t;

}

void main(void)

{

clrscr(); textmode(2);

cout< < " Vvod 0 "; gotoxy(9, 1); cin> > n;

int i, j;

if(n< =0) example();

else

{

init(1);

cout< < " Vvod x[" < < n< < " ] ";

for(i=0; i< n; i++)

cin> > X[i];

cout< < " Vvod y[" < < n< < " ] ";

for(i=0; i< n; i++)

cin> > Y[i];

}

h=X[1]-X[0];

for(i=0; i< n; i++)

A[i][0]=Y[i];

for(j=1; j< n; j++)

for(i=0; i< n-j; i++)

A[i][j]=A[i+1][j-1]-A[i][j-1];

float z;

cout< < " \nVvod x: "; cin> > z;

float N=0, K=1, t;

if(z< X[n/2])

{t=(z-X[0])/h;

for(j=0; j< n; j++)

{N+=K*A[0][j]/f(j);

K*=(t-j); }

}

else

{t=(z-X[n-1])/h;

for(j=0; j< n; j++)

{N+=K*A[n-j-1][j]/f(j);

K*=(t+j); }

}

cout< < " \nf(" < < z< < ")=" < < N;

if(pr) cout< < " \n\nProverka: tg(cos(" < < z< < "))=" < < sin(cos(z))/cos(cos(z));

init(0);

getch();

}

Скриншот результата программы (при x=0, 4):

Результаты работы программы и проверка:

Результаты работы программы и проверка (проверка говорит о том, что решение правильное):

f(0, 4) = 1, 316354

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

Скриншот результата программы (при x=0, 8):

Результаты работы программы и проверка:

Результаты работы программы и проверка (проверка говорит о том, что решение правильное):

f(0, 8) = 0, 836674