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

Страницы работы

Фрагмент текста работы

РАСЧЕТНО-ГРАФИЧЕСКОЕ ЗАДАНИЕ 1

Задание 1. Построить алгоритм для интерполяции функции  на отрезке  интерполяционным линейным сплайном  с равномерным шагом. Разработать программу, которая реализует этот алгоритм и выводит на печать значение сплайна  и значение функции , значение производной сплайна и производной функции при . Размерность интерполяционной таблицы равна 21.

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

Варианты заданий

Номер

варианта

1

-1

1

0.15

2

-2

2

0.1

3

-1.5

1.5

0.2

4

-1.2

1.2

0.3

5

-1.6

1.6

0.4

6

-1

1

0.25

7

-2

2

-0.1

8

-1.5

1.5

-0.2

9

-1.2

1.2

-0.3

10

-1.6

1.6

-0.4

11

-1

1

-0.15

12

-2

2

-0.5

13

-1.5

1.5

0.1

14

-1.2

1.2

0.2

15

-1.6

1.6

0.6

16

-1

1

-0.25

17

-2

2

0.3

18

-1.5

1.5

-0.1

19

-1.2

1.2

-0.2

20

-1.6

1.6

-0.6

21

-1

1

0.35

22

-2

2

-0.3

23

-1.5

1.5

0.4

24

-1.2

1.2

0.1

25

-1.6

1.6

0.2

26

-1

1

-0.35

27

-2

2

0.5

28

-1.5

1.5

-0.4

29

-1.2

1.2

-0.1

30

-1.6

1.6

-0.2

ПРИМЕР ВЫПОЛНЕНИЯ ЗАДАНИЯ 1

Построим алгоритм для интерполяции функции  на отрезке  интерполяционным линейным сплайном с равномерным шагом. Будем опираться на материал, изложенный в п.п. 2.1.

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

Узлы сетки  находятся по формуле . Значения  находятся по формуле . Таким образом, по формуле функции  на отрезке  построена интерполяционная таблица  размерности 21 с равномерным шагом. Так как все узлы сетки различны, то существует единственный интерполяционный линейный сплайн , удовлетворяющий этой таблице. Для построения сплайна  сначала найдем его коэффициенты по формулам

,   ,   ,   .

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

.


Программа

#include <stdio.h>

#include <conio.h>

#include <math.h>

float A, B;          // границы отрезка [а, b]

float Xi[21], Yi[21]; // интерполяционная таблица

float Ai[21], Bi[21]; // коэффициенты сплайна ai, bi

float h;             // шаг сетки

int   n;             // число отрезков

// Функция sin(x)

float fn (float x)  { return (sin (x));}

// Вычисление значений сплайна L(x)

float L (float x)

{

int i=int ((x-A)/h);   // номер интервала

return Ai[i]+Bi[i]*(x-Xi[i]);

}

void main ()

{

// Инициализация начальных значений

A=-2, B=2;          // отрезок [-2,2]

n=20;               // число отрезков n

h=float (B-A)/n;    // шаг h=(b-a)/n

// Создание интерполяционной таблицы

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

{

Xi[i]=A+i*h;

Yi[i]=fn (A+i*h);

};

// Расчет коэффициентов сплайна L(x)

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

{

Ai[i]=Yi[i];

Bi[i]=(Yi[i+1]-Yi[i])/h;

};

Ai[n]=Yi[n],

Bi[n]=Bi[n-1];

printf ("Значение функции f в точке x=%.2g равно %.2g\n", 0.1, fn (0.1));

printf ("Значение сплайна L в точке x=%.2g равно %.2g\n", 0.1, L (0.1));

getch ();

}

Результаты работы программы

Значение функции f в точке x=0.1 равно 0.1

Значение сплайна L в точке x=0.1 равно 0.099

ПРИМЕР ВЫПОЛНЕНИЯ ЗАДАНИЯ 2

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

Как и в задании 1, число отрезков , шаг, где , . Введем сетку на отрезке : .

Узлы сетки  находятся по формуле . Значения  находятся по формуле . Таким образом, по формуле функции  на отрезке  построена интерполяционная таблица  размерности 21 с равномерным шагом. Так как все узлы сетки различны, то существует единственный интерполяционный кубический сплайн , удовлетворяющий этой таблице.

Для построения сплайна  сначала найдем его коэффициенты. Коэффициенты  находятся по явным формулам . Для нахождения коэффициентов , решается система линейных уравнений с трехдиагональной матрицей:

  ,

Эта система линейных уравнений имеет единственное решение, которое находится методом прогонки. Запишем формулы метода прогонки для данного примера. Сначала в порядке возрастания находятся прогоночные коэффициенты  и :

,   ,   ,   ,   , где   .

Затем в порядке убывания находятся неизвестные :

,   ,  

После нахождения коэффициентов  и , находим по явным формулам коэффициенты  и :

,    ,    ,

,   

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

Значение производной сплайна при  находим по формуле

.

Программа

#include <stdio.h>

#include <conio.h>

#include <math.h>

float A,B;                     // границы отрезка

float Xi[21],Yi[21];           // интерполяционная таблица

float Ai[21],Bi[21],Ci[21],Di[21];   // коэффициенты ai, bi, ci, di

float Alpha[21],Betta[21];     // прогоночные коэффициенты

float h;                       // шаг сетки

float n;                       // число отрезков

// Функция sin(x)

float fn (float x) { return (sin(x)); }

// Производная функции sin(x)

float fsh (float x) { return (cos(x)); }

// Вычисление значений сплайна S(x)

float Sx (float x)

{

int i=int((x-A)/h);

return Ai[i]+Bi[i]*(x-Xi[i])+Ci[i]*(x-Xi[i])*(x-Xi[i])+Di[i]*(x-Xi[i])*

(x-Xi[i])*(x-Xi[i]);

}

// Вычисление значений производной сплайна

float Sshx (float x)

{

int i=int((x-A)/h);

return Bi[i]+2*Ci[i]*(x-Xi[i])+3*Di[i]*(x-Xi[i])*(x-Xi[i]);

}

void main()

{

clrscr();

float Ti=0;       // вспомогательная  переменная ti

// Инициализация начальных значений

A=-2,B=2;         // отрезок [-2,2]

n=20;             // число отрезков n

h=float(B-A)/n;   // шаг

// Создание интерполяционной таблицы

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

{

Xi[i] = A+i*h;

Yi[i] = fn(Xi[i]);

}

// Находим коэффициенты сплайна ai, bi, ci, di

// Расчет коэффициентов ai

for(i=0;i<=n;i++) Ai[i]=Yi[i];

// Решение СЛУ с трехдиагональной матрицей методом прогонки

//                           для нахождения коэффициентов ci

// Вычисление прогоночных коэффициентов (прямая прогонка)

Alpha[1] = 0;

Betta[1]  = 0;

for(i=1;i<=n-1;i++)

{

Alpha[i+1] = (-1)/(4+Alpha[i]);

Ti         = 3*(Yi[i+1]-2*Yi[i]+Yi[i-1])/(h*h) ;

Betta[i+1] = (Ti-Betta[i])/(4+Alpha[i]);

}

// Вычисление коэффициентов ci

// Обратная прогонка

Ci[n]=0;

for(i=n-1;i>=0;i--) { Ci[i] = Alpha[i+1]*Ci[i+1]+Betta[i+1]; }

// Находим коэффициенты bi,di

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

{

Bi[i] = ((Yi[i+1]-Yi[i])/h)-(h*(Ci[i+1]+2*Ci[i])/3);

Di[i] = (Ci[i+1]-Ci[i])/(3*h);

}

Bi[n] = Bi[n-1]+2*Ci[n-1]*h+3*Di[n-1]*h*h;

Di[n] = Di[n-1];

printf("Значение функции f в точке x=%.2g равно %.4g\n",0.1,fn(0.1));

printf("Значение сплайна S в точке x=%.2g равно %.4g\n",0.1,Sx(0.1));

printf("Значение производной функции f' в точке x=%.2g равно

%.4g\n",0.1,fsh(0.1));

printf("Значение производной сплайна S' в точке x=%.2g равно

%.4g\n",0.1,Sshx(0.1));

getch();

}

Результаты работы программы

Значение функции f в точке x=0.1 равно 0.09983

Значение сплайна S в точке x=0.1 равно 0.09983

Значение производной функции f' в точке x=0.1 равно 0.995

Значение производной сплайна S' в точке x=0.1 равно 0.995

РАСЧЕТНО-ГРАФИЧЕСКОЕ ЗАДАНИЕ 2

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

Замечание. Для каждого варианта функция  выбирается та же, что и в расчетно-графическом задании 1.

Варианты заданий

Номер

варианта

g(x)

a

b

1

-0.2

0.2

2

0

0.5

3

-0.2

0.2

4

-0.18

0.22

5

-0.2

0.2

6

-0.2

0.2

7

0

0.6

8

-0.4

0

9

0

0.5

10

-0.5

0.5

11

-0.2

0.2

12

0.5

1.5

13

-0.5

0.5

14

-0.5

0.5

15

-0.2

0.2

16

-0.5

0.5

17

0.5

1.5

18

0.5

0.7

19

-1

0.2

20

0.73

0.89

21

-0.2

0.2

22

-0.2

0.2

23

-0.2

0.2

24

0

0.2

25

-0.5

0.5

26

0.4

0.9

27

0.6

0.9

28

-0.5

0.5

29

0.5

1

30

0.4

0.8

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

Варианты заданий

Номер варианта

1

20

-1.35; -1.66; -1.96; -2.26; -2.56; -2.86; -3.16; -3.45; -3.75; -4.04; -4.33;

-4.61;-4.90; -5.18; -5.45; -5.73; -6.00; -6.26; -6.52; -6.78

2

21

0.75; 0.44; 0.14; -0.16; -0.46; -0.76; -1.06; -1.35; -1.65; -1.94; -2.23; -2.51; -2.80; -3.08; -3.35; -3.63; -3.90; -4.16; -4.42; -4.68; -4.93

3

17

1.81; 2.56; 3.31; 4.06; 4.80; 5.53; 6.25; 6.95; 7.64; 8.32; 8.98; 9.62; 10.24; 10.84; 11.41; 11.96; 12.48

4

19

1.07; 1.14; 1.21; 1.28; 1.35; 1.42; 1.48; 1.55; 1.61; 1.67; 1.72; 1.77; 1.82; 1.87; 1.91; 1.95; 1.98; 2.01; 2.04

5

18

0.69; 0.32; -0.04; -0.40; -0.76; -1.11; -1.46; -1.80; -2.13; -2.46; -2.77;

-3.08;-3.38; -3.67; -3.95; -4.21; -4.46; -4.70

6

20

1.36; 1.42; 1.48; 1.53; 1.59; 1.65; 1.71; 1.77; 1.82; 1.88; 1.94; 2.00; 2.05; 2.11; 2.17; 2.22; 2.28; 2.33; 2.38; 2.44

7

21

0.93; 0.65; 0.38; 0.11; -0.16; -0.43; -0.69; -0.94; -1.19; -1.44; -1.67; -1.91;

-2.13; -2.34; -2.55; -2.75; -2.93; -3.11; -3.27; -3.43; -3.57

8

19

1.56; 2.07; 2.58; 3.09; 3.58; 4.08; 4.56; 5.04; 5.51; 5.97; 6.41; 6.85; 7.27; 7.67; 8.06; 8.43; 8.78; 9.12; 9.43

9

18

1.78; 2.05; 2.33; 2.59; 2.85; 3.09; 3.33; 3.54; 3.75; 3.93; 4.09; 4.23; 4.35; 4.45; 4.52; 4.57; 4.60; 4.60

10

17

-2.07; -2.63; -3.19; -3.74; -4.28; -4.80; -5.31; -5.80; -6.27; -6.72; -7.14;

-7.53; -7.89; -8.23; -8.53; -8.79; -9.02

11

20

2.21; 2.91; 3.61; 4.29; 4.96; 5.62; 6.25; 6.86; 7.45; 8.01; 8.53; 9.02; 9.47; 9.89; 10.26; 10.59; 10.88; 11.12; 11.31; 11.45

12

21

1.42; 1.44; 1.46; 1.48; 1.50; 1.51; 1.53; 1.55; 1.57; 1.59; 1.61; 1.63; 1.64; 1.66; 1.68; 1.70; 1.72; 1.73; 1.75; 1.77; 1.78

13

19

-1.63; -1.66; -1.69; -1.72; -1.75; -1.78; -1.81; -1.83; -1.86; -1.89; -1.92;

-1.95; -1.98; -2.01; -2.03; -2.06; -2.09; -2.12; -2.14

14

18

0.11; -0.52; -1.16; -1.78; -2.40; -3.02; -3.62; -4.22; -4.80; -5.36; -5.91;

-6.45; -6.97; -7.46; -7.94; -8.40; -8.83; -9.24

15

17

2.39; 3.02; 3.66; 4.28; 4.90; 5.52; 6.12; 6.72; 7.30; 7.86; 8.41; 8.95; 9.47; 9.96; 10.44; 10.90; 11.33

16

20

1.27; 1.22; 1.17; 1.13; 1.08; 1.03; 0.98; 0.93; 0.88; 0.84; 0.79; 0.74; 0.69; 0.65; 0.60; 0.55; 0.51; 0.46; 0.42; 0.37

17

21

0.68; 0.11; -0.46; -1.03; -1.60; -2.16; -2.72;

-3.27; -3.82; -4.36; -4.89; -5.42; -5.94; -6.44;

-6.94; -7.43; -7.91; -8.38; -8.84; -9.28; -9.71

18

19

-1.82; -2.39; -2.96; -3.53; -4.10; -4.66; -5.22; -5.77; -6.32; -6.86; -7.39;

-7.92; -8.44; -8.94; -9.44; -9.93; -10.41; -10.88; -11.34

19

18

2.76; 4.14; 5.51; 6.87; 8.21; 9.55; 10.86; 12.15; 13.42; 14.66; 15.86; 17.04; 18.18; 19.28; 20.34; 21.35; 22.32; 23.24

20

17

-1.57; -1.90; -2.23; -2.55; -2.86; -3.16; -3.44; -3.70; -3.95; -4.17; -4.37;

-4.54; -4.68; -4.80; -4.89; -4.95; -4.98

21

20

-2.11; -2.88; -3.65; -4.41; -5.15; -5.88; -6.59; -7.28; -7.95; -8.59; -9.19;

-9.77; -10.31;-10.82; -11.28; -11.71; -12.09; -12.43; -12.72; -12.97

22

21

-1.82; -2.61; -3.40; -4.17; -4.94; -5.70; -6.44; -7.16; -7.86; -8.55; -9.21;

-9.84; -10.45; -11.02; -11.57; -12.08; -12.56; -13.00; -13.40; -13.77;

-14.09

23

19

0.74; 0.66; 0.58; 0.50; 0.42; 0.34; 0.26; 0.18; 0.10; 0.02; -0.06; -0.13;

-0.21;-0.29; -0.36; -0.44; -0.51; -0.58; -0.65

24

18

0.43; 0.11; -0.21; -0.52; -0.83; -1.12; -1.40; -1.67; -1.93; -2.17; -2.39;

-2.59; -2.77; -2.93; -3.07; -3.19; -3.28; -3.34

25

17

0.84; 0.78; 0.71; 0.65; 0.58; 0.52; 0.45; 0.39; 0.33; 0.26; 0.20; 0.14; 0.07; 0.01; -0.05; -0.11; -0.17

26

20

-0.59; -0.53; -0.48; -0.42; -0.36; -0.30; -0.25; -0.19; -0.13; -0.08; -0.02; 0.03; 0.09; 0.14; 0.20; 0.25; 0.31; 0.36; 0.41; 0.46

27

21

-0.72; -0.78; -0.85; -0.91; -0.98; -1.05; -1.11; -1.18; -1.24; -1.31; -1.37;

-1.43; -1.50; -1.56; -1.62; -1.69; -1.75; -1.81; -1.87; -1.93; -1.99

28

19

-0.58; -0.51; -0.44; -0.37; -0.31; -0.24; -0.17; -0.10; -0.03; 0.03; 0.10; 0.17; 0.23; 0.30; 0.36; 0.42; 0.49; 0.55; 0.61

29

18

-0.40; 0.05; 0.49; 0.94; 1.38; 1.82; 2.26; 2.70; 3.14; 3.57; 4.00; 4.42; 4.84; 5.25; 5.66; 6.06; 6.46; 6.85

30

17

-1.30; -1.75; -2.19; -2.64; -3.08; -3.52; -3.96; -4.40; -4.84; -5.27; -5.70;

-6.12; -6.54; -6.95; -7.36; -7.76; -8.16

ПРИМЕР ВЫПОЛНЕНИЯ ЗАДАНИЯ 1

Задание. Построить алгоритм для аппроксимации функции  на отрезке  регрессионной функцией  с использованием метода наименьших квадратов. Разработать программу, реализующую этот алгоритм. Размерность равномерной сетки равна 21.

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

Сначала построим сетку и таблицу данных. Так как сетка равномерная и число узлов сетки равно 21, то число отрезков  равно 20 (). Обозначим через  шаг, который вычисляется по формуле . В нашем случае , .

Введем сетку на отрезке : . Узлы сетки находятся по формуле ,  . Значение находится по формуле . Таким образом, по формуле функции  на отрезке  построена таблица данных

Похожие материалы

Информация о работе