Методи парабол (Сімпсона) і вищих ступенів (Ньютона)
Матеріал із MachineLearning.
Зміст
У цій статті описується спосіб обчислення свого інтеграла гладкої функції за допомогою квадратурних формул. Формули Ньютона-Котеса мають такі особливості:
- Необхідною умовою збіжності даного методу є існування та обмеженість похідної функції (порядок похідної залежить від обраної формули)
- Формули Ньютона-Котеса володіють високим порядком точності
- Формули порядку є чисельно стійкими
Постановка математичного завдання
Завдання чисельного інтегрування полягає у знаходженні наближеного значення інтегралу
де - задана та інтегрована на відрізку функція. На відрізку вводиться сітка і як наближене значення інтеграла розглядається число
де - значення функції у вузлах, де -вагові множники, що залежать тільки від вузлів, але не залежать від вибору. Формула (2) називаєтьсяквадратурною формулою. Завдання чисельного інтегрування за допомогою квадратур полягає у знаходженні таких вузлів і таких ваг, щоб похибка квадратурної формули
була мінімальною за модулем для функції із заданого класу (величина залежить від гладкості). Похибка залежить як від розташування вузлів, і від вибору вагових коефіцієнтів. Введемо нарівномірну сітку з кроком, тобто. безліч точок , і представимо інтеграл (1) у вигляді суми інтегралів за частковими відрізками:
Для побудови формули чисельного інтегрування на всьому відрізку достатньо побудувати квадратурну формулу для інтегралу
на частковому відрізку та скористатися властивістю (3).
Побудова квадратурних формул
З огляду на вищевикладене, обчислення наближеного значення інтеграла проводиться за допомогою квадратурної формули
Цю формулу за допомогою заміни можна привести до стандартного вигляду
У загальному випадку вузли та ваги невідомі та підлягають визначенню.
Розглянемо випадок, коли вузли задані та потрібно знайти ваги квадратурної формули . Користуватимемося вимогою: формула (5) повинна бути точною для будь-якого полінома ступеня . Для того, щоб поліном ступеня задовольняв цю вимогу, достатньо зажадати, щоб квадратурна формула була точною для будь-якого одночлена ступеня . Враховуючи, що , отримуємо рівняння
Ця система має єдине рішення, оскільки її визначником є визначник Вандермонда, відмінний від нуля якщо немає вузлів, що збігаються, .
Так, вважаючи, маємо систему, рішенням якої є ваги формули Сімпсона: . Таким чином, формула Сімпсона є точною для полінома другого ступеня. Однак, з симетрії, вона є точною і для всіх поліномів третього ступеня:
так як вона точна для :
Формули трикутника та трапеції точні для лінійної функції, тобто. для полінома першого ступеня, у чому легко переконатись безпосередньо. У загальному випадку як можна вибрати інтерполяційний поліном Лагранжа
де – інтерполяційний коефіцієнт Лагранжа. З рівності
видно, що формула (5) правильна для полінома ступеня, якщо вагові коефіцієнти визначаються за формулою
Формули такого типу називають квадратурнимиформулами Котеса.
Виклад методу
Застосування квадратурних формул
Повернемося до розгляду інтегралу (1). Як було показано вище, даний інтеграл заміною зводиться до інтегралу на одиничному відрізку, а отже легко узагальнити формулидля наближеного обчислення інтеграла на одиничному відрізку довільний. Застосуємо апарат квадратурних формул. Нехай встановлено рівномірне розбиття відрізка з кроком, де позначимо. Нехай також обрано деяку квадратурну формулу Ньютона-Котеса (тобто обрано ступінь полінома, а отже кожен поліном будується по точці сітки). Ми також вважаємо, що даний набір з точки можна розбити на піднабори по крапках з останніми крайніми, тобто
Тоді, підсумовуючи значення квадратур кожному піднаборі, отримаємо наближене значення шуканого інтеграла. Якщо позначити вагові множники, то наближене значення інтеграла можна записати у вигляді подвійної суми.
Цей алгоритм природним чином узагальнюється на випадок, коли , де , але в кожному відрізку задана рівномірна сітка. Тоді шуканий інтеграл дорівнює
але в кожному зчасткових відрізківнаближене значення інтеграла обчислюється з допомогою квадратурних формул.
Приклади квадратурних формул
Наведемо приклади квадратурних формул Котеса на рівномірній сітці з кроком, де позначимо:
Аналіз методу
Похибка квадратурної формули
Нехай функція має безперервну похідну на відрізку , тобто , - точки, в яких задано значення функції. Нехай використовується квадратурна формула порядку. Введемо функцію
Тоді формула для похибки має такий вигляд
Звідси випливає оцінка для похибки
при , де 0" alt= "M_0" /> - постійна, і при
Якщо не змінює знаку на відрізку, то через теорему про середнє маємо
Чисельна стійкість квадратурних формул
Зазначимо, що формули Ньютона-Котеса рідко використовуються через їх чисельну нестійкість, що призводить до різкого зростання обчислювальної похибки.Причиною такої нестійкості є те, що коефіцієнти формул Ньютона-Котеса при великих мають різні знаки, а саме при існують як позитивні, так і негативні коефіцієнти.
Розглянемо квадратурну суму
Припустимо, значення функції , заданої на відрізку , обчислюються з деякою похибкою, тобто замість точного значення отримуємо наближене значення . Тоді замість отримаємо суму
Оскільки квадратурна формула точна для , маємо
і залежить від .
Припустимо тепер, що це коефіцієнти неотрицательны. Тоді з ( 11 ) та ( 12 ) отримаємо оцінку
яка означає, що за великих похибка у обчисленні квадратурної суми ( 10 ) має той самий порядок, що й похибка у обчисленні функції. І тут кажуть, що сума ( 10 ) обчислюється стійко.
Якщо коефіцієнти мають різні знаки, то може виявитися, що сума не є рівномірно обмеженою і, отже, похибка у обчисленні необмежено зростає зі зростанням . У цьому випадку обчислення за формулою ( 10 ) будуть нестійкі і користуватися такою формулою за великих не можна.
Числовий експеримент
Наведемо приклади обчислення інтегралів із застосуванням формул Ньютона-Котеса. При реалізації використовувався мову C++, наведено нижче код функції, що повертає наближене значення інтеграла.
Вихідний код функції
На вхід функція приймає 4 параметри
double a - лівий кінець досліджуваного відрізка
double b - правий кінець досліджуваного відрізка
int Degree - ступінь використовуваного полінома
int Ndivisions – кількість відрізків, на які розбивається вихідний. Збігається зі значенням у формулі ( 7 )
f - інтегрована функція
Обчислимо за формулами Котесу ступеніввід 1 до 9 значення інтегралу
Обчислення будемо проводити не розбиваючи відрізок на часткові, тобто використовуючи точку, де -ступінь полінома. Результати наведені нижче у таблиці. Округлення проводилося з точністю 6 знаків після коми.
Рекомендації програмісту
Автоматичний вибір кроку інтегрування за допомогою апостеріорної оцінки похибки методом Рунґе
Величина похибки чисельного інтегрування залежить як від кроку сітки, і від гладкості підінтегральної функції. У величину похибки, крім того, входить також величина, яка може сильно змінюватися на відрізку і заздалегідь невідома. Для зменшення величини похибки можна подрібнити сітку на заданому відрізку. Але при цьому необхідно апостеріорно оцінювати похибку. Таку оцінку похибки можна здійснити методом Рунге. Розглянемо застосування якоїсь квадратурної формули на частковому відрізку. Позначимо за значення інтеграла на всьому відрізку, за значення інтеграла на частковому відрізку, за наближене значення інтеграла на всьому відрізку, отримане за допомогою заданої квадратурної формули і рівномірному сітки з кроком, а за наближене значення інтеграла на частковому відрізку. Нехай ця квадратурна формула на даному частковому відрізку має порядок точності, тобто
де з – деяка константа. Тоді
Нехай використовується складова квадратурна формула
причому на всіх часткових відрізках використовуються квадратурні формули з одним і тим же порядком точності (або ж, зокрема, та сама формула). Проведемо на кожному частковому відрізку обчислення двічі - один раз з кроком, другий раз з кроком та апостеріорно оцінимо похибку за правилом Рунґе (14). Якщо для заданого
тобто буде досягнуто заданої точності
Якщо ж у якомусь із часткових відрізків оцінка ( 15 ) нічого очікувати виконуватися, то крок цьому відрізку треба подрібнити ще двічі і знову оцінити похибку. Подрібнення сітки на даному відрізку слід проводити доти, доки не буде досягнуто оцінки виду ( 15 ). Зауважимо, що для деякої функції таке подрібнення може тривати надто довго. Тому у відповідній програмі слід передбачити обмеження зверху на кількість подрібнень.
Таким чином, автоматичний вибір кроку інтегрування призводить до того, що інтегрування ведеться з великим кроком на ділянках плавної зміни функції та з дрібним кроком - на ділянках швидкої зміни. Це дозволяє за заданої точності
Висновок
Формули Сімпсона і Ньютона-Котеса є хорошим апаратом для обчислення певного інтеграла достатню кількість разів безперервно диференційованої функції. На прикладі формул -го порядку бачимо, що з обмеженої похідної -го порядку точність формул є , де - крок сітки, а n" alt= "p>n" />. Тобто при виконанні умов застосування даного методу, формули дають Проте при великих , зокрема вже при , обчислення наближеного значення інтеграла стає чисельно нестійким, що робить їх непридатними.