Метод дробления шага. Правило Рунге
Тема 4. Численное ИНТЕГРИРОВАНИЕ В прикладных исследованиях часто возникает необходимость вычисления значения определенного интеграла , К сожалению, в подавляющем большинстве случаев получить значение интеграла с помощью формулы Ньютона-Лейбница или других аналитических методов не удается. Так, интеграл широко используется при исследовании процессов массо- и теплообмена, в статистической физике и теории вероятностей. Однако его значение не может быть выражено в виде конечной комбинации элементарных функций. При разработке численных методов нахождения значения определенных интегралов наиболее распространенным подходом является аппроксимация. Подынтегральную функцию f(x) на отрезке [a, b] приближают некоторой функцией F(x), которая легко интегрируется аналитически. Затем полагают . Если функция f(x) задана аналитически, то ставится вопрос об оценке погрешности такой замены.
Формула трапеций
1.1. Вывод формулы. Пусть функция f(x) интерполируется на [a, b] многочленом 1-й степени: F(x) = P1 (x) = y0 + (х – x0) ×D1[x0, x1], где x0 = a, x1= b. Обозначим h = x1–x0. Или Тогда
I*= .
Сделаем замену переменных t = x – x0, получим
I*= = =
После интегрирования и преобразования получим следующую формулу:
» [f(a) + f(b)] = [f(x0) + f(x1)] = I*, (4.1)
называемой формулой трапеций. Геометрически она означает замену площади криволинейной трапеции обычной трапецией (рис. 4.1), отсюда название формулы.
|
|
|
1.2. Погрешность формулы. Воспользуемся известной оценкой (3.9) погрешности интерполирования многочленом P1(x):
R1(z) = .
Имеем: R1 = f(x) – F(x) Þ D = I – I* = = =
= = .
Здесь мы опять сделали замену переменных t = x – x0. Тогда x = t + x0, а x – x1 = t + (x0 – x1) = t – (x1 – x0) = t – h. Отсюда
= = – f’’(x)× .
Итак,
DТРАП = – f’’(x)× , xÎ(a, b). (4.2)
Знак “–” означает, что при f’’ > 0 формула (4.1) дает значение интеграла с избытком, а при f’’ < 0 – с недостатком. На рис. 4.1 функция f(x) вогнута, т.е. f’’ < 0, следовательно, I* < I. При f’’ = 0, т.е. если f(x) = P1(x), формула (4.1) не содержит погрешности и абсолютно точна.
Формула Симпсона
2.1. Вывод формулы. Пусть функция f(x) интерполируется на [a, b] многочленом 2-й степени:
F(x) = P1 (x) = y0 + (х – x0) ×( D1[x0, x1] + (х – x1) ×D2[x0, x1, x2] ),
где x0 = a, x1 = (a + b)/2, x2 = b, h = xi+1 – xi. Преобразуем формулу, используя явное выражение для D1[x0, x1] = и D2[x0, x1, x2] = :
P1 (x) = y0 + ×(х – x0) + ×(х – x0)×(х – x1).
Сделаем замену переменных t = x – x0. Тогда x = t + x0, а x – x1 = t – h. После интегрирования и преобразования получим следующую формулу:
» [ f (x0) + 4f (x1) + f (x2) ] = I*, (4.3)
|
|
называемой формулой Симпсона. Ее геометрический смысл представлен на рис. 4.2.
|
Оценка погрешности формулы Симпсона.
Метод вывода оценки погрешности, применённый для формулы трапеции, не годится для формулы Симпсона. Поэтому без доказательства примем следующую оценку погрешности формулы (4.3):
DСИМП = I – I* = , xÎ[ x0 , x2]. (4.4)
Знак “–” означает, что при f (4) > 0 формула (4.3) дает значение интеграла с избытком, а при f (4) < 0 – с недостатком. На рис. 4.3 функция f(x) представляет собой многочлен 4-й степени с постоянным значением f (4) < 0. Из рисунка видно, что площадь фигуры под графиком f(x) (окрашенная область) больше, чем под F(x). Следовательно, I* < I. При f (4) = 0, т.е. если f(x) = P3(x), формула (4.3) не содержит погрешности и абсолютно точна. Это видно на рис. 4.2 – f(x) представляет собой многочлен 3-й степени и площади под f(x) и под F(x) абсолютно одинаков, хотя сами функции не совпадают.
|
3. Составные формулы численного интегрирования
Описанные формулы численного интегрирования легко обобщить на случай применения многочленов более высоких степеней и тем самым добиться улучшения точности решения задачи (подобный метод широко применяют, при этом, в частности, получаются так называемые квадратурные формулы Ньютона-Котеса). Однако, более перспективен другой способ, основанный на кусочной аппроксимации подынтегральной функции многочленами невысокой степени – 1-й или 2-й.
|
|
Идея метода: разбить отрезок [a, b] на n частей и на каждой части применить формулу трапеций или Симпсона. Из формул (4.2) и (4.4) видно, что при h < 1 погрешность этих формул уменьшается.
3.1. Составная формула трапеций. Разобьём отрезок [a, b] на n частей равноотстоящими точками x0 = a, x1 = a + h, …, xn = b, и на каждой части [xi, xi+1] применим формулу трапеций:
» ×( yi + yi+1).
Получим
» ×[( y0 + y1) + ( y1 + y2) + … + ( yn–1 + yn)] Þ
» h× – (4.5)
составная формула трапеций. Её геометрический смысл – подынтегральная функция заменяется ломанной.
3.2. Оценка погрешности. Очевидно, погрешность формулы (4.5) равна сумме погрешностей на каждом i-м интервале. Воспользуемся формулой (4.2), получим = – , где xiÎ[xi, xi+1].
Обозначим m = – среднее арифметическое значений второй производной, а также m = , M = . По известному свойству среднего арифметического имеем: m £ m£ M.
|
|
Из математического анализа известна 2-я теорема Больцано-Коши: если функция y(x) непрерывна на отрезке, то она принимает все свои значения, от минимального до максимального. Иными словами, для любого р, такого что min y(x) £ p £ max y(x), уравнение p = y(x) всегда имеет решение при непрерывной функции у(х). Для разрывной функции это, вообще говоря, неверно (см. рис. 4.4).
|
Если f’’(x) непрерывна на [a, b], то к ней применима указанная теорема при р = m. Значит, существует такая точка xÎ[a, b], что m = f’’(x). Следовательно,
= n×m = n×f ’’(x). Тогда
= – = – .
Следовательно,
= – . (4.6)
3.3. Составная формула Симпсона. Разобьём отрезок [a, b] на n частей, где n – чётно. Применим формулу Симпсона к каждому двойному отрезку [x2i–2, x2i]. Тогда
» ×[(y0 + 4y1 + y2) + (y2 + 4y3 + y4) + … + (yn–2 + 4yn–1 + yn)].
Заметим, что y2, y4,… yn–2 повторяются дважды. Отсюда
» ×[( y0 + yn) + 4Sнечётн + 2Sчётн ] (4.7)
составная формула Симпсона.
Здесь Sнечётн = y1 + y3 + … + yn–1, Sчётн = y2 + y4 + … + yn–2.
3.4. Оценка погрешности. Для оценки погрешности формулы (4.7) воспользуемся оценкой (4.4) погрешности простой формулы. Применим тот же метод, что и для составной формулы трапеции. Получим
= – . (4.8)
Задание. Вывести формулу (4.8) самостоятельно.
Метод дробления шага. Правило Рунге
Формулы (4.6) и (4.8) для практики неудобны, т.к. значения f(k) (x) обычно неизвестны. Правило Рунге позволяет найти достаточно точные оценки погрешности, используя значения I*, вычисленные при различных h.
Согласно полученным оценкам, погрешность составных формул численного интегрирования имеет вид
I – Ih » C×hm, (4.9)
где I и Ih – точное значение интеграла и значение, найденное по составной формуле при шаге h, соответственно; С – некоторая константа; m – порядок точности численного метода (согласно формулам (4.6) и (4.8), m = 2 для формулы трапеции и m = 4 для формулы Симпсона).
Уменьшим шаг вдвое. Будем иметь:
I – Ih/2 » 2–m×C×hm. (4.10)
Вычтем (4.10) из (4.9). Получим: Ih/2 – Ih » 2–m×C×hm(2m – 1) = (I – Ih/2)(2m – 1).
Следовательно,
. (4.11)
Формула (4.11) – оценка погрешности вычисления Ih/2 по правилу Рунге. Если значение правой части приближенного равенства (4.11) окажется меньше предельно допустимой погрешности e, то можно считать Ih/2 ответом задачи. В противном случае – требуемая точность не достигнута и необходимо использовать более мелкий шаг, т.е. продолжить процесс дробления шага.
Замечание. При приближённом вычислении интегралов, абсолютная величина | I | которых мала, использование значения абсолютной погрешности D для контроля точности необоснованно, более приемлема относительная погрешность d. Поэтому при использовании правила Рунге рекомендуется следующий критерий останова процесса дробления шага:
если | Ih/2 | ³ 1, то заканчивать при ,
если | Ih/2 | < 1, то заканчивать при .
Пример использования правила Рунге см. «Вычислительная практика» Пример 4.4. Там же в разделе 5 приведены тестовые примеры для программ численного интегрирования.
Задание. Обобщите правило Рунге на случай, когда шаг h на каждой итерации уменьшается не в 2, а в s раз.
Дата добавления: 2018-04-15; просмотров: 1919; Мы поможем в написании вашей работы! |
Мы поможем в написании ваших работ!