Учебник - Практические занятия по вычислительной математике - Аристова (1238839), страница 29
Текст из файла (страница 29)
Эти вычисления «стоили» всего четырех вычислений подынтегральной функции. Поучительно сравнить их с тем,сколько вычислений этой функции потребуется при «студенческом» рецепте (для достижения такой же точности).Рассмотрим пример вычисления интеграла второго типа:2e x dx.0Можно, как уже отмечалось, свести его к интегралу первого типа.Но мы воспользуемся универсальным приемом выделения особенности.Особенность состоит в том, что верхний предел интегрирования – бесконечность.
Представим интеграл в виде суммы двух интегралов:I = I1 + I2, где I1 – интеграл по конечному отрезку [a, A]; I2 – интеграл по[A, ∞]. Вычисление I1 при заданном A затруднений не вызывает.Выберем теперь A так, чтобы в пределах допустимой погрешностивторым интегралом можно было пренебречь, т.
е. так, чтобы |I2| ≤ ε / 2.Например, учитывая, что при A ≥ 12e x dx A22x e x dx 1 e A ,A2и требуя, чтобы выполнялось условие1 e A2 1 ,22найдем A | ln |.VII.5. ВычислениефункцийинтеграловотбыстроосциллирующихНачнем с простого примера. Пусть требуется вычислить1etsin kt dt0при большом значении k, например, k = 100.
Интегралы типа193VII. ЧИСЛЕННОЕ ИНТЕГРИРОВАНИЕ1 f (t ) sin kt dt ,0где f (t) – гладкая функция, часто приходится вычислять в некоторыхразделах физики. Сложность задачи состоит в том, что подынтегральнаяфункция совершает большое число колебаний. Вычисление интеграла постандартной формуле Симпсона, конечно, возможно, но требует сетки сочень малым шагом: каждая волна должна быть описана некоторым числом узлов сетки, а волн много.Дело осложняется еще и тем, что вычисление должно проводиться свысокой точностью, так как результат есть сумма большого числа близких величин с противоположными знаками (интегралов от отдельныхволн подынтегральной функции), происходит сильное сокращение знаков и для обеспечения точности остатка (результата), отдельные слагаемые должны вычисляться с существенно более высокой точностью. Длявычисления подобных интегралов используется следующий прием:гладкая функция f (t) аппроксимируется некоторой другой гладкойфункцией f *(t), такой, чтобы интеграл от f *(t) sin kt вычислялся аналитически.Итак, дело сводится к тождественному преобразованию000f (t )sin kt dt f * (t )sin kt dt f (t ) f * (t ) sin kt dt.Второе слагаемое является малым и отбрасывается.
Правда, если оценить отбрасываемую величину, опираясь только на оценку типа |f (t) –f *(t)| ≤ ε, т.е. в данном случае величиной πε, ничего хорошего (даже еслиε – точная оценка погрешности аппроксимации) не получится, так каквеличина πε может оказаться значительно большей интересующего насинтеграла. На самом деле погрешность существенно меньше. Это ведьинтеграл от гладкой функции, не превосходящей ε, умноженной набыстроосциллирующую функцию.
Естественно ожидать, что погрешность будет во столько раз меньше результата, во сколько раз |f – f *|меньше f. При f (t) = e – t интеграл вычисляется аналитически.VII.6. Задачи на доказательствоVII.6.1. Доказать, что вычисление интеграла от строго выпуклой функции f '' > 0 по формулам прямоугольников со средней точкой дает заниженное значение интеграла.194VII.7. Задачи с решениямиVII.6.2. Доказать, что вычисление интеграла от строго выпуклой функции f '' > 0 методом трапеций дает завышенное значение интеграла.VII.6.3. Доказать, что для погрешности квадратурной формулы трапецийсправедливо представлениеbb1ba()()()( x a)( x b) f ( x) dx.fxdxfafba22 aVII.6.4.
Доказать, что при применении правила Рунге к формуле трапеций получается формула Симпсона. Насколько при этом повышаетсяпорядок точности метода?aVII.6.5. Пусть для вычисления интегралаf ( x) p( x)dx с четной весовойaфункцией p(x) используется симметричный относительно нуля наборузлов xn 1 k xk , k 1,..., n .
Доказать, что веса, соответствующиесимметричным узлам, будут равны: cn 1 k ck , k 1,..., n .VII.6.6. Показать, что операция построения формулы (2.8)I p 1 I p (h / 2) ( I p (h / 2) I p (h)) / (2 p 1)является экстраполяцией, т.е. при Ip(h/2) ≠ Ip(h) величина Ip+1 всегда лежит вне отрезка с концами Ip(h/2) и Ip(h).bVII.6.7. Пусть Cq f (q) ( x) dx , q 4 . Получить оценку погрешaности формулы Симпсона R4 ( f ) k Cq hq , где k – абсолютная постоянная, h – шаг интегрированная.VII.7. Задачи с решениями1VII.7.1.
Оценить погрешность вычисления интеграла I 10 1x2dx какинтеграла от сеточной функции, вычисляемого методом трапеций последующим данным:x 0. 0.1250.250.3750.5 0.6250.75f(x) 1.0 0.984615 0.941176 0.876712 0.8 0.719101 0.641950.87510.566372 0.5VII. ЧИСЛЕННОЕ ИНТЕГРИРОВАНИЕРешение. Воспользуемся правилом Рунге практической оценки погрешности. Вычисление интеграла методом трапеций приводит к результату Ih = 0.784747, вычисление с двойным шагом I2h = 0.782794.
Выбранметод второго порядка точности, поэтому оценка погрешностиI Ih I h I 2h22 10.784747 0.782794 0.000651 .3Эту величину можно использовать для уточнения результата экстраполяцией Ричардсона IR = 0.785398. Заметим, что точное значение интеграла с шестью значащими цифрами I = π/4 = 0.785398, экстраполяцияРичардсона дает точный ответ при заданной точности входных данных.1VII.7.2. Вычислить интеграл I 10 1 sinx2dx от заданной сеточнойфункции методом трапеций, сделать уточнение экстраполяцией Ричардсона, сравнить результат с вычислениями методом Симпсона.
Предложить второй вариант уточнения результата метода трапеций на заданнойсетке.0.250.3750.50.6250.750.8751x 0. 0.125f(x) 1.0 0.984616 0.941213 0.877068 0.801665 0.724235 0.652187 0.590671 0.543044Решение. Вычисления методом трапеций приводят к значению интеграла Ih = 0.792891, вычисления с двойным шагом I2h = 0.791647, тогдаэкстраполяция Ричардсона даетI Ih I h I 2h22 1 0.792891 0.792891 0.791647 0.793306 .3Этот результат совпадает с вычислением методом Симпсона.Второй способ уточнения вычисления интеграла методом трапеций отизвестной функции связан с применением формулы Эйлера–Маклоренадля вычисления интеграла на каждом отдельном отрезке интегрирования:I xi , xi 1 xi 1xif ( x ) dx hh2 f ( xi ) f ( xi 1 ) f ( xi 1 ) f ( xi ) O(h5 ).212При фиксированном шаге сетки суммирование интегралов по каждомуотрезку приведет к составной формуле трапеций с поправкой Эйлера–196VII.7.
Задачи с решениямиМаклорена:hh2 f ( xi ) f ( xi 1 ) f (b) f (a) O(h5 ),12i 0 2bNI f ( x) dx axi a (b a) i / N .Вычислим производную подынтегральной функции:12 x cos x 2,f ( x) 2 (1 sin x 2 )2 1 sin x тогда f (0) = 0, f (1) = – 0.318667, и уточненное значение интегралаI Ih h2(0.125)2(0.318667 0) 0.793306 , f (1) f (0) 0.792891 1212таким образом, оба способа уточнения в данном случае приводят к одному и тому же результату.При задании значений функций с шестью значащими цифрамибольшую точность вычисления интеграла мы достичь не можем.VII.7.3.
Пример применения алгоритма Рунге–Ромберга. Для вычисления интеграла11 x2 1.1 x dxиспользуется таблица значений подынте-0гральной функции:x 0.0.1250.250.3750.50.6250.750.8751f(x) 0.909091 1.007266 1.139113 1.278655 1.443376 1.643421 1.889822 2.151657 0.00Специально взята функция, имеющая максимум ближе к одному из концов интегрирования (максимум в точке x = 1/1.1) и обращающаяся в нульна этом же конце. Поведение функции на правом краю интервала интегрирования сеточной функцией никак не прописано. Этот интеграл можно вычислить точно, и его значение равно0.42 0.11 I ex 1.1 1 arctgarctg 1.485130 .0.21 0.210.21 Применение метода трапеций дает приближенное значение интеграла Ih = 1.375982, вычисление с двойным шагом I2h = 1.231714, уточнение экстраполяцией Ричардсона совпадает с применением метода Симп197VII.
ЧИСЛЕННОЕ ИНТЕГРИРОВАНИЕсона IS = IR = 1.424071. Оценка точности любого из методов по производным в данном случае будет некорректна, т.к. максимальное значениепроизводных на правом конце равно бесконечности, по этой же причиненельзя использовать поправку Эйлера–Маклорена. Как было показано вдвух предыдущих задачах, для функций с меньшими значениями производной уточнение экстраполяцией Ричардсона дает очень хороший результат.Построим решение по алгоритму Рунге–Ромберга на основе экстраполяции вычисленных интегралов в нулевой шаг интегрирования.Ошибка численного интегрирования в силу симметричности используемых формул разлагается в ряд по h2, поэтому как для метода трапецийвычисления определенного интеграла, так и метода Симпсона экстраполяцию в нулевой шаг следует делать не по величине шага h, а по величине h2.Метод трапеций для различных шагов интегрирования позволяетвычислить соответствующие значения интеграла и составить таблицуразделенных разностей для экстраполяции:h2(0.5) = 0.25I0.948961(0.25)2 = 0.06251.231714(0.125)2 = 0.0156251.3759822Δ1Δ2–1.508016–3.07771736.697392Экстраполяция в нулевой шаг интегрирования дает I = 1.430612,ошибка по сравнению с точным решением довольно велика (3,67%).Составим таблицу разделенных разностей вычисления интеграламетодом Симпсона для различных шагов интегрирования:h2(0.5) = 0.25I1.113765833(0.25)2 = 0.06251.32596525(0.125)2 = 0.0156251.4240722Δ1–1.131730222–2.092944Δ24.1011788Заметим, что экстраполяция Ричардсона по двум последним значениям (метод шестого порядка аппроксимации) дает I = 1.430612, что198VII.7.
Задачи с решениямисовпадает с результатом экстраполяции в нулевой шаг интегрированияпо трем сеткам метода трапеций, что естественно.Экстраполяция в нулевой шаг интегрирования по трем значениямметода Симпсона (метод восьмого порядка) дает I = 1.460779, ошибкавдвое меньше (1,64%), но тоже велика.Измельчим шаг интегрирования еще вдвое, при этом впервые появляется единственная точка на убывающем участке функции:h2(0.25)2 = 0.0625IΔ1Δ2Δ31.32596525–2.092944–79.845679(0.125)2 = 0.0156251.42407223.7507015–3.48458(6)(0.0625)2 = 0.00390625 1.464907Экстраполяция в нулевой шаг интегрирования по трем последнимсеткам дает I = 1.479968, ошибка 0,35%.
Экстраполяция в нулевой шагинтегрирования по четырем последним шагам дает I = 1.480273, ошибкапримерно такая же –– 0,33%. Заметим, что вычисление методом Симпсона даже для самой мелкой сетки дает значительно большую ошибку –1.36%.Интересно сравнить этот результат с интегрированием по Симпсону, когда поведение функции на правом конце учтено еще лучше:h2(0.1) = 0.01I1.442913(0.05)2 = 0.00251.4709212Δ1–3.7344000–5.02613333Δ2137.78488888(0.025)2 = 0.000625 1.480345Экстраполяция в нулевой шаг интегрирования для последних трехмелких сеток дает результат I = 1.483702 с погрешностью (0,096%).1VII.7.4.