Федоренко Р.П. Введение в вычислительную физику (Федоренко Р.П. Введение в вычислительную физику.djvu), страница 11

Вычисление значений 5л для нескольких Ф открывает возможность в качестве «бесплатного приложения» получить значительно более точное значение интеграла. Это достигается простой процедурой экстраполяции полученных значений. Идея очень проста. Пусть для величины Я имеются приближенный метод вычисления с малым параметром т и теоретическая оценка 5=5, + Стг+ о(тг). Посюянная С, конечно, неизвестна, порядок же погрешности р иавестен. 3 41 вычислеиие опгеаелеииых иитнтллов Если вычислены две величины Ю, и Л„мы имеем два следующих соотношения: Я= Лм + С(гт)г+ о(тл), Я= Б, + Стг + о(тг), Из них можно найти о с точностью до о(тл), исключив главный член погрешности. Умножая второе уравнение на гг и вычитая из него первое, получаем 5 = (гав, — Бы)/(гг — 1) + о(тг) В табл, 6 представлены результаты: значения Я„вычисленные по формуле Симпсона, и значения Б, полученные экстраполяцией Ричардсона по значениям Яп Явл и Зол, Бом. Для всех величин приведена погрешность (абсолютная).

Выше было указано, что погрешность формулы Симпсона не превосходит о(т'), поэтому экстраполяция проводилась сначала при р= 3. Результат оказался обескураживающим: экстраполяция не только ие повысила точности, но и дала заметно худший результат. В чем же дело? Причина кроется в очень важном обстоятельстве, которое никогда не следует забывать, применяя подобную экстраполяцию. Обычно порядок погрешности р устанавливается на основе теоретиче- Таблица б скнх оценок, которые, как правило, дают завышенное значение.

Завышение оценки может влиять как на коэффициент при тг, так и на саму степень. Фактически точность может быть более высокой. Те же соображения могут быть использованы в условиях, когда не только коэффициент С в оценке погрешности, но и степень р не известны. Разумеется, теперь надо иметь три приближенных значения, чтобы исключить член Стг с двумя неизвестными. Формулы выводятся просто: б,=я+ Се" + о( ), бм = и+ г С" + о("), б„= я+2 г С. + о(. ), основы вычяслитзльноя млтвмлтихн Пренебрегая членами о(тг), имеем (Я, Ям)/(Я, Я,) 2 . Используя числа табл. 6, получаем 2г 14.64 = 21~~. Поскольку порядок точности формул обычно является целым числом, положим р= 4, относя разницу на счет величин О(тг).

Если теперь проделать эксграполяцию при р= 4, получим числа, представленные в той же табл. 6. Эффект уточнения явно виден и не нуждается в комментариях. Однако читатель должен усвоить и другой урок: применение этой экстраполяции требует определенной осторожности, особенно в снежных задачах, в которых нет полной ясности с главным членом асимптотикн погрешности.

Что касается формулы Симпсона, то более аккуратный анализ показывает, что она на самом деле имеет четвертый порядок точности, Это дополнительное повышение ее точности — следствие симметричности формулы. Симметричность часто и в других случаях приводит к повышению точности. Например, формула численного дифференцирования (односторонняя, несимметричная) и формула центральной разности (снмметричная) имеют первый и второй порядки точности соответственно, хотя они одинаковы по трудоемкости (каждая гстонть двух вычислений функции).

Несколько слов о чувствительности формулы Симпсона по отношению к возможной неточности априорного предположения о функции /. Пусть на самом деле функция /(1) лишь кусочно-гладкая, для простоты, имеет один разрыв на [а, Ь[, а в остальном гладкая. Нетрудно понять, что если разрыв попал на четный узел сетки, точность формулы сохраняется. Прн произвольном положении разрыва испортится только значение интеграла по той паре шагов сетки, в которую попал разрыв.

Этот интеграл и вычисленный по Снмпсону не имеют между собой ничего общего кроме того, что оба суть О(т). Следовательно„точность вычисления интеграла станет О(т) вмеето ожидаемой О(тз). Вычисление интегралов с особенностями. Рассмотрим простой пример — приближенное вычисление интеграла пв г о Подынтегральная функция обращается в бесконечность при 1'= О, но интеграл существует. Попытка его прямого вычисления но формуле Симпсона сразу же приведет к неудаче: первое слагаемое Уо 53 вычисления опгвдвязнных интепъяов- обращается в бесконечность.

Грамотный студент легко находит выход! надо вычислить интеграл ! 1 — "'! дг, е>0; т7 ь при достаточно малом е он приближает нужное значение. И даже оценка связанной с этим погрешности легко вычисляется: зто величина с — в!! = 2!/е. а Некоторые даже догадаются взять в качестве приближенного решения сумму ! 2г!е + ~ — (Й, с где интеграл по [з, 11 вычисляется, например, по формуле Симпсона.

В вычислительной математике многие проблемы допускают решения на уровне студенческой грамотности. Нередко на основе подобных соображений работы пишутся людьми далеко не студенческого возраста„знающими о проблемах этой нзуки с чужих слов и не «чувствующих» такого важнейшего фактора, как число операций (речь идет, конечно, о более сложных задачах). И здесь приведенный выше, вполне правильный ответ не устраивает профессиональных вычислителей именно в силу нерационально большого обьема вычислений. Аккуратное вычисление интеграла с особенностью может быть выполнено гораздо более экономными средствами.

Это достигается с помощью приема рееуляриэации, или выделения особенности. Поясним его в более общей ситуации. Пусть требуется вычислить где у'(!) — гладкая функция. Регуляризация состоит в том, что проделывается тождественное преобразование ! ! ! ~ у[,),-!!г в!, ~ [у[,) р[,)],-!п,~с+ ~ р[,),-г!г в!, о о о Функция р(!) выбирается такой, чтобы первый интеграл правой части не содержал особенности и при небольшом объеме вычислений достаточно точно определялся хотя бы по формуле Симпсона. Второй интеграл особенность содержит, но вычисляется аналитически.

основы вычислительной млтвылтнкн 54 В данном случае цель будет достигнута, если в качестве р(!) взять отрезок ряда Тейлора У(г) в точке ! = О. Это приводит к вычислению ! ! ! 1 У!!! У!о1 1~ !о! ж+ЛО)1 ! 41 +!'(0)1 ~ с! /! !!! о о о В примере, с которого мы начали (у(1) = соз !), приходим к вычислению Второе слагаемое есть 2, первое вычислим по формуле Симпсона: сначала с шагом 0.5, что даст значение 1.807967, затем с шагом 0.25, что даст значение 1.808850. Эти вычисления *стоили» всего четырех вычислений подынтегральной функции. Поучительно сравнить их с тем, сколько вычислений этой функции потребуется при «студенческом» рецепте (для достижения такой же точности).

Вычисление интегралов от быстроосциллнрующих функций. Начнем с простого примера. Пусть требуется вычислить ') е !з!пЫ Й о при большом значении я, например к 100. Интегралы типа ~ У(!) зш х! а1, где ~(!) — гладкая функция, часто приходится вычислять в некоторых разделах физики.

Сложность задачи состоит в том, что подынтегральная функция совершает большое число колебаний. Вычисление интеграла по стандартной формуле Симпсона, конечно, возможно, но требует сетки с очень малым шагом: каждая волна должна быть описана некоторым числом узлов сетки, а волн мною. Дело осложняется еще и тем, что вычисление должно проводиться с высокой точностью, так как результат есть сумма большого числа близких величин с противоположными знаками (интегралов от отдельных волн подынтегральной функции), происходит сильное сокращение знаков и для обеспечения точности остатка (результата) отдельные слагаемые должны вычисляться с существенно более высокой точностью. Дая вычисления подобных интегралов используется следующий прием: гладкая функция Яг) аппроксимируется некоторой другой гладкой функцией у(!), такой, чтобы интеграл от ,г(!) з1п к!' вычисаялся аналитически.

55 Я 4! ВБ1ЧИСЛЕНИЕ ОПРЕДЕЛЕННЫХ ННТЕГГАЛОЕ Итак, дело сводится к тождественному преобразованию л л л ~ 1(!) Е1п и ы! = ~ 7(!) зш ! ! 5! + ][1(г) — 1(!)] 51п к! (!. о а о Второе слагаемое является малым и отбрасывается. Правда, если оценить.ртбрасываемую величину, опираясь только на оценку типа ]У(!) — у'(!) ] ц е, т,е. в дан- ном случае величиной яе, ничего хорошего (даже если е — точная оценка погрешности аппроксимации) не получится, так как величина яе может оказаться значительно большей интере- ! 4 1 4 2 4 2 4 сующего нас интеграла. На о 6 самом деле погрешность существенно меньше. Это ведь Рис 7 интеграл от гладкой функции, не превосходящей с, умноженной на быстроосциллирующую функцию.

Естественно ожидать, что погрешнос!ь будет во столько раз меньше результата, во сколько раз ]У вЂ” 7] меньше У. При у(!) = е ' интеграл вычисляется аналитически, Поучительно вычислить его приближенно, заменив функцию е ' ее интерполяционным полиномом всего лишь второй степени. Интересно, совпадает ли точность результата с ожиданиями? На этой идее построена формула механической квадратуры, аналогичная формуле'Симпсона.

Интервал интегрирования разбивается на четное число шагов длиной т, на каждой паре шагов функция заменяется ингерполяционным полиномом второй степени 2;2(!), интегралы от Ьз(!) 51п й! вычисляются аналитически, полученные выражбния суммируются. Вычисление многомерных интегралов. Метод Моите-Карло.

Рассмотрим задачу вычисления интегрзла по многомерному кубу 11 1 ~ ''' ~ ! (Х!' ХМ "'~ Хл) иХ! икз" о о о Нетрудно и здесь построить формулы механических квадратур, аналогичные, например, формуле Симпсона. Проще всего такие формулы получить, используя процедуру повторного интегрирования, т.е. заменяя многомерный интеграл иа равный ему ! 1 1 ]с с(х, ~ а'хз ... ]г У(х„х, ..., х„) !5х„. о о о 1ч. г ОСНОВЫ ВЫЧИСЛИТЕЛЬНОЙ МАТЕМАТИЕЯ 56 В двумерном случае для «элементарной ячейки» имеем коэффициенты (аналогичные коэффициентам 1, 4, 1 в одномерной формуле Симпсона), изображенные на рис. 7О.