Міністерство освіти і науки України
Вінницький національний технічний університет
Інститут АЕКСУ
Факультет ФЕЛТ
Кафедра Електроніки
КУРСОВА РОБОТА
З дисципліни “Обчислювальна математика”
Алгоритм та програма знаходження значення функції за допомогою інтерполяційної формули Бесселя
2004

ЗАВДАННЯ
на курсову роботу з дисципліни
“Обчислювальна математика”
студентові гр. МП-02 ТвердохлібуА.М. пропонується розробити
алгоритм та програму мовою програмування Турбо Паскаль на знаходження значення функції  за допомогою інтерполяційної формули Бесселя
Основні вхідні дані:
1. Кількість вузлів таблично заданої функції – не більше 100.
2. Похибка обчислень – не більше 0,001.
Основні вихідні дані:
1. Пояснювальна записка до курсової роботи.
2. Виконуваний файл програми.

АНОТАЦІЯ
В даній курсовій роботі розроблений ефективний алгоритм та програма мовою Турбо Паскаль знаходження значення функції за допомогою інтерполяційної формули Бесселя. Розроблений алгоритм є досить непоганим за розміром пам’яті, необхідної для збереження даних, котрі обчислюються в ході виконання алгоритму, та за кількістю арифметичних операцій для обчислення за основною формулою.

ВСТУП
Задача знаходження значення функції у міжвузловій точці за допомогою інтерполяційної формули Бесселя має важливе значення при вирішенні як наукових, так і практичних задач, оскільки дає можливість знаходження значення функції у будь-якій точці, в якій це потрібно. В багатьох випадках функція не має аналітичного вигляду, тобто він невідомий, а задана лише декількома точками та значеннями функції в цих точках. Тому для отримання значення функції в проміжних точках застосовуються інтерполяційна формули Гауса (1-а та 2-а), інтерполяційна формула Стірлінга та Бесселя. Останні дві формули є похідними від першої та другої інтерполяційних формул Гауса. Кожна з цих формул має свої переваги та недоліки, що полягають у кількості обчислювальних операцій та в похибці обчислень.
Серед сучасного програмного забезпечення є багато програм чисельного аналізу, до яких можна віднести всесвітньо відомі пакети програм MathCad та MatLab. Вони, як правило, мають зручний інтерфейс та є багатофункціональними. Але їх недоліком є те, що задачі чисельного аналізу певного класу (наприклад, знаходження першої похідної) вирішуються за допомогою лише деякого одного методу. Крім того, вони займають багато дискової пам’яті та вимагають певного часу для того, щоб навчитися ними користуватися. Тому в даній курсовій роботі була поставлена задача розробити програму знаходження значення функції у міжвузловій точці за допомогою інтерполяційної формули Бесселя, яка займала б небагато пам’яті та була б простою у користуванні.
Курсова робота складається з трьох основних розділів. В першому розділі наведені основні теоретичні відомості про метод знаходження значення функції у міжвузловій точці за допомогою інтерполяційної формули Бесселя та приклад його застосування. У другому розділі розроблено алгоритм за даним методом. Третій розділ містить загальний опис програми, лістинг програми та результати тестування.

ТЕХНІЧНЕ ЗАВДАННЯ
1.Основою для проведення роботи є навчальний план кафедри Електроніки ВНТУ.
Замовник - кафедра Електроніки ВНТУ.
Виконавець – студент гр.. МП – 02 Твердохліб А.М.
2. Мета роботи.
Метою роботи є розробка ефективного алгоритму та програми мовою Турбо Паскаль знаходження значення функції за допомогою інтерполяційної формули Бесселя.
3. Етапи виконання роботи.
4.1 Кількість вузлів таблично заданої функції – не більше 100.
4.2 Похибка обчислень – не більше 0,001.
4.3 Алгоритм повинен бути оптимізований за часом виконання та розміром.
4.4 Програма повинна бути розроблена за принципами структурного та модульного програмування.
5. Спосіб реалізації результатів.
Робота повинна закінчуватись передачею замовнику пояснювальної записки до курсової роботи та виконуваного файла програми.
6. Техніко економічне обґрунтування.
В результаті виконання роботи планується розробити ефективний алгоритм та програму знаходження значення функції за допомогою інтерполяційної формули Бесселя з мінімальним часом виконання та розміром.
7. Порядок розгляду і приймання роботи:
7.1 Курсова робота приймається комісією в складі двох викладачів за участю керівника роботи.
7.2 Програма перевіряється шляхом тестування на комп’ютері тестових завдань, розроблених замовником і виконавцем роботи.
8. Додаткові відомості.
Дане ТЗ може змінюватись і корегуватись за спільною домовленістю замовника та виконавця.

ТЕОРЕТИЧНІ ВІДОМОСТІ. ІНТЕРПОЛяція ФУНКЦІЙ

1. ПОСТАНОВКА ЗАДАЧІ ІНТЕРПОЛЯЦІЇ

Нехай деяка функція у=f(х) задана таблицею (табл.1), тобто при значеннях аргументу х=х₀, х₁, ... , х_n функція f(х) приймає відповідні значення
у₀, у₁,... , у_n.
Таблиця 1

Таблиця експериментальних значень

Також нехай необхідно визначити значення у=f(х), (х_i-1< <х_i). Величина х=  потрапляє між двома табличними значеннями, тому для обчислення значення функції необхідно запропонувати деякий характер її зміни між відомими експериментальними даними.
Інтерполяцію можна розглядати як процес визначення для даного аргументу х значення функції у=f(х) по її декількох відомих значеннях. При цьому розрізняють інтерполяцію у вузькому смислі, коли х знаходиться між x₀ і x_n, і екстраполювання, коли х знаходиться поза відрізком інтерполяції [x₀ , x_n].
Задача інтерполяції полягає в наступному. На відрізку [а, b] задані n+1 точки х₀, х₁, ... , х_n, що називаються вузлами інтерполяції, і значення деякої функції f(x) у цих точках.
f(x₀) = y₀;
f(x₁) = y₁;(1)
f(x_n) = y_n
Потрібно побудувати функцію Р_n(х) (інтерполюючу функцію), яка б задовольняла таким умовам:
P_n(x₀) = y₀;
P_n(x₁) = y₁;(2)
P_n(x_n) = y_n
тобто інтерполююча функція Р_n(х) повинна приймати ті ж значення, що і функція f(х), яку ми визначаємо (що інтерпелюється), для вузлових значень аргументу х₀, х₁, ... , х_n.
Геометрично це означає, що потрібно знайти криву y=P_n(х) деякого визначеного типу, що проходить через задану систему точок М_i (х_i,у_i) (i=0,1,2,..,n). Очевидно, можна побудувати множину неперервних функцій, що будуть проходити через задані вузлові точки.[1]
Заміна функції f(х) її інтерполяційним багаточленом Р_n(x) може знадобитися не тільки тоді, коли відома лише таблиця її значень, але і коли аналітичний вираз для f(х) відомо, проте є занадто складним і незручним для подальших математичних перетворень (наприклад, для інтегрування, диференціювання та ін.). Іноді розглядаються задачі тригонометричної інтерполяції (інтерполююча функція – тригонометричний поліном). Інтерполюючою може бути також раціональна функція.
У загалі залежність, якої підпорядковується функція, може бути апроксимована багаточленом ступеня n:
Р_n(x) = y = a₀ + a₁ ∙ x + a₂ ∙ x² + ... + a_n ∙ xⁿ.                               (3)
Таку задачу називають задачею параболічної інтерполяції (або інтерполюванням).
Загалом є багато інтерполяційних формул та методів. До них відносяться такі: інтерполяційні формули Гауса (дві), Стерлінга та Бесселя (які є похідними від формул Гауса), Ньютона (дві) та багато інших.

2. ПАРАБОЛІЧНА ІНТЕРПОЛЯЦІЯ

Для визначення коефіцієнтів багаточлена (3) необхідно мати n+1 вузлову точку. Аналітичне визначення коефіцієнтів інтерполяційного багаточлена для n+1 точки зводиться до рішення системи лінійних рівнянь n+1 порядку, кожне з яких являє собою вираз (3), записаний для визначеної вузлової точки
y_i = a₀ + a₁ ∙ x_i + a₂ ∙ x_i² + ... + a_n ∙ x_iⁿ,(4)
де i = 1, 2,. . . n+1.
Даним методом побудови інтерполяційного поліному зручно користуватися, маючи персональний комп’ютер і відповідні програми. Даний метод не є єдиним способом побудови інтерполяційного поліному. Інший підхід, яким часто користуються на практиці, називається методом Лагранжа.[2]
3. МЕТОД ЛАГРАНЖА
Нехай при х=х₀, х₁, ... , х_n функція f(х) приймає відповідно значення у₀, у₁,... , у_n. Багаточлен ступеня не вище n, що приймає у вузлових точках задані значення, має вид:
Р_n(х)=у= . (5)

Цей багаточлен (5) називається інтерполяційною формулою Лагранжа і має такі властивості:
1.                При заданій сукупності вузлових точок будова багаточлена можлива тільки єдиним способом.
2.                Багаточлен Лагранжа може бути побудовано при будь-якому розташуванні вузлів інтерполяції (включаючи і нерівномірне).
У розгорнутому виді форма Лагранжа має вид:
Р_n(х)=  +
+ +
+ … +  +
+ … + .(6)
При n=1 формула Лагранжа має вид:
Р(х) = (7)
і називається формулою лінійної інтерполяції.
При n=2 одержимо формулу квадратичної інтерполяції:
Р(х)= . (8)

4. ЗВОРОТНА ІНТЕРПОЛЯЦІЯ
Нехай функція у= f(х) задана таблицею. Задача зворотної інтерполяції полягає в тому, щоб по заданому значенню функції у визначити відповідне значення аргументу х.
Якщо вузли інтерполяції x₀, x₁, x₂, … x_n нерівновіддалені, задача легко вирішується за допомогою інтерполяційної формули Лагранжа (5). Для цього достатньо прийняти у за незалежну змінну, а х вважати функцією. Тоді отримаємо
x = (9)
Розглянемо тепер задачу зворотної інтерполяції для випадку рівновіддалених вузлів інтерполяції. Припустимо, що функція f(х) монотонна і дане значення у знаходиться між y₀=f(x₀) і y₁ = f(x₁).
Замінюючи функцію у=f(x) першим інтерполяційним багаточленом Ньютона, одержимо:
y = y₀ + q Dy₀ + D²y₀ + D³y₀ +…+ Dⁿy₀ .
Звідси
q = D²y₀  – …–Dⁿy₀ ,
тобто q=j(q).
Розмір q визначаємо методом послідовних наближень як границю послідовності:
q = ,
де q_i = j (q_i–1) (i=1, 2,…).
За початкове наближення приймаємо
q₀ =  (10)
Для i-го наближення маємо:
 q_i = q₀ – D²y₀  – …–Dⁿy₀ . (11)
На практиці ітераційний процес продовжують доти, поки не установляться значення, що відповідають необхідній точності, причому q » q_m, де m – останнє зі знайдених наближень. Знайдемо q, визначаємо х по формулі
 = q,
звідки
х = x₀ + q h.(12)

Ми застосували метод ітерації для рішення задачі зворотної інтерполяції, користуючись першою інтерполяційною формулою Ньютона. Аналогічно можна застосувати цей спосіб і до другої формули Ньютона:
y = y_n + qDy_n–1 + D²y_n–2 + D³ y_n–3 + …
+ Dⁿy₀ .
Звідси
q = D^2yn-2  – …–Dⁿy₀ ...
Позначимо q₀ =  – початкове наближення.
Для i-го наближення маємо:
q_i = q₀–D²y_n–2  – …–Dⁿy₀ ...  (13)
Знайдемо
q = ,
визначимо х по формулі
х = x_n + q h .[3], [2]

Далі розглянемо запропоновану мені інтерполяційну формулу Бесселя, яка часто використовується для знаходження значення функції у між вузловій точці. Вона подібна до інтерполяційної формули Стерлінга і обидві вони є похідними від першої та другої інтерполяційних формул Гауса.
ІНТЕРПОЛЯЦІЙНА ФОРМУЛА БЕССЕЛЯ
Часто використовується інтерполяційна формула Бесселя. Для виведення цієї формули скористаємось другою інтерполяційною формулою Гауса:

у скороченому вигляді:

де х=х₀+qh
Візьмемо 2n+2 рівновіддалених вузлів інтерполювання
x_-n, x_-(n-1),..., x₀,..., x_n-1, x_n, x_n+1
з кроком h, і нехай
y_i= f(x_i) (i =-n,…,n+1)
- задані значення функції y= f(x).
Якщо вибрати за початкові значення x= x₀ та y= y₀, то, використовуючи вузли x_k (k= 0, ±1, …,  n), будемо мати:
 (1)
Приймемо тепер за початкові значення х=х₁ і у=у₁ і використаємо вузли х_1+к (к=0, 1,..., n). Тоді

причому відповідно індекси всіх різниць в правій частині формули (1) зростуть на одиницю. Замінивши в правій частині формули (1) q на q-1 і збільшивши індекси всіх різниць на 1 , отримаємо допоміжну інтерполяційну формулу
 (2)
Взявши середнє арифметичне формул (1) і (2), після простих перетворень отримаємо інтерполяційну формулу Бесселя
 (3)
де
Інтерполяційна формула Бесселя (3), як слідує з способа отримання її, представляє собою поліном, що співпадає з даною функцією y= f(x) в 2n+2 точках
x_-n , x_-(n-1),…, x_n , x_n+1.
В частинному випадку, при n=1, нехтуючи різницею ∆³y_-1, отримаємо формулу квадратичної інтерполяції по Бесселю
P(x)=  
Або

де

В формулі Бесселя всі члени, які містять різниці непарного порядку, мають множник q- ; тому при  формула (3) значно спрощується :

Цей спеціальний випадок формули Бесселя називається формулою інтерполювання на середину.    Якщо в формулі Бесселя (3) зробити заміну по формулі
то вона приймає більш симетричний вид

де
ЗАГАЛЬНИЙ ОПИС ПРОГРАМИ
В програмі використано кілька процедур та функцій, в яких використовуються різноманітні позначення, змінні та тому подібне. Тому далі буде дано пояснення (розшифровування) що яким символом, чи їх сукупністю позначено.
PROCEDURE vvod – процедура введення даних, тобто задання початкових умов (кількість вузлів, задаються X та відповідні їм Y, а також X, в яких потрібно знайти значення ф-ї).
PROCEDURE ddd –дана процедура формує трикутну таблицю різниць значень
ф-ї.
PROCEDURE rech- процедура рішення задачі за інтерполяційною формулою Бесселя.
PROCEDURE vivod- це процедура виведення результатів роботи програми, які виводяться файлі f2, який має назву «ANA.NAS».
FUNCTION pos1 (j, q: real): real; - функція, яка обчислює для формули Бесселя чисельник її доданків.
FUNCTION fak (j: integer): real; - функція для обрахунку факторіала в знаменнику доданків формули.
Також в програмі використовуються такі позначення:
kol- позначено кількість точок (вузлів), в яких потрібно знайти значення функції.
zad- назва масива, в який заносимо ікси, значення функції в яких бажаємо знайти.
otv- назва масива для відповідних ігреків.
f1- це файл, який має назву «st1.tab», з нього читаються задані ікси та ігреки.
f2- файл, який має назву «ana.nas», - це файл виведення, з нього читаються результати роботи програми.
n- кількість вузлів інтерполяції.
h- крок задання точок.

Література
1.         Вычислительная техника в инженерных и экономических расчетах
/ А.В.Крушевский, А.В.Беликов, В.Д.Тищенко – Киев: «Высшая школа». Главное изд-во, 1985. – 290 с.
2. Дифференциальное и интегральное исчисления /Н.С.Пискунов- Москва: «Наука». Главная редакция физико-математической литературы, 1978.-576 с.
3.Численые методы в инженерных исследованиях /В.Е.Краскевич, К.Х.Зеленский, В.И.Гречко - Киев: Главное изд-во «Высшая школа», 1986.-263 с.