Блюмин А.Г., Федотов А.А., Храпов П.В. - Численные методы (1078410), страница 3
Текст из файла (страница 3)
Может случиться, что c = 0. Тогдаглавная часть погрешности формулы средних прямоугольников является величиной порядка h4 . Но обычно c 6= 0.Если f ∈ C 4 [a, b], то можно получить также соотношениеI = Ihmp + c1 · h2 + O(h4 ),(1.2.15)где Ihmp — приближенное значение интеграла I, найденное по составнойформуле трапеций с шагом h; c1 — постоянная, не зависящая от h, c1 =Rb 001− 12 f (x)dx.aЕсли f ∈ C 6 [a, b], то аналогично выражениям (1.2.14) и (1.2.15) можнополучить следующее соотношениеI = IhC + c · h4 + O(h6 ),(1.2.16)где IhC — приближенное значение интеграла I, найденное по составной формуле Симпсона; c — некоторая не зависящая от h постоянная.Правило Рунге.
Пусть Ih — приближенное значение интеграла I, найденное по одной из трех рассмотренных составных формул (по формуламсредних прямоугольников, трапеций и Симпсона). Объединим соотношения(1.2.14), (1.2.15) и (1.2.16) в одно:I = Ih + c · hk + O(hk+2 ),19(1.2.17)Глава 1. Численные методы вычисления определенного интегралагде c не зависит от h, k — порядок точности квадратурной формулы (k = 2для составных формул средних прямоугольников и трапеций, k = 4 для составной формулы Симпсона).
Предполагается, что f ∈ C k+2 [a, b].На основании формулы (1.2.17) можем записать, что khI = Ih/2 + c ·+ O(hk+2 ).2(1.2.18)Вычитая равенство (1.2.18) из (1.2.17), находим khIh/2 − Ih = c ·(2k − 1) + O(hk+2 ).2Отсюда kIh/2 − Ihhc·= k+ O(hk+2 )22 −1и, следовательно, согласно формуле (1.2.18), с точностью до O(hk+2 ) имеемI − Ih/2 ≈Ih/2 − Ih.2k − 1(1.2.19)Вычисление приближенной оценки погрешности по формуле (1.2.19) привыполнении условия (1.2.17), т.е. при возможности представления значенияинтеграла I в виде (1.2.17), называется правилом Рунге.Вычитая из умноженного на 2k равенства (1.2.18) равенство (1.2.17), получаемI · (2k − 1) = 2k · Ih/2 − Ih + O(hk+2 ).Отсюда I = Ih∗ + O(hk+2 ), гдеIh∗2k · Ih/2 − Ih=.2k − 1Число Ih∗ называется уточненным по Ричардсону приближенным значением интеграла I.20Глава 1.
Численные методы вычисления определенного интеграла1.3. ЗаданиеДля предложенного варианта лабораторной работы интегралZbI=f (x)dxaвычислите:1) аналитически,2) численно с точностью до ε = 0.0001:• по формуле средних прямоугольников,• по формуле трапеций,• по формуле Симпсона.Точность вычислений определяется с помощью правила Рунге. Точность ε, с которой необходимо найти приближенное значение интеграла,считается достигнутой, когда в процессе вычислений будет выполненонеравенство| Ih/2 − Ih |< ε.2k − 1Алгоритм вычислений с использованием правила Рунге.
Приближенное вычисление интеграла с заданной точностью ε проводим методом итераций. На l-той итерации вычисляем значение Il = Ih интегралаI по одной из трех требуемых составных формул приближенного вычисления интегралов с шагом hl , затем находим значение Il+1 = Ih/2 потой же составной формуле, но с шагом hl+1 =hl /2. Если для найденныхзначений Il и Il+1 выполняется неравенство| Il+1 − Il |< ε,2k − 121(1.3.1)Глава 1. Численные методы вычисления определенного интегралато точность считается достигнутой.
В противном случае проводим следующую итерацию: Il присваиваем значение Il+1 , увеличиваем в двараза число разбиений n, находим новое значение Il+1 и опять проверяем выполнение условия (1.3.1).При вычислении начального приближения I0 (для l = 0) в качестве√шага h0 можно взять значение h0 ≈ k ε. Однако, при этом, соответствующее значению h0 первоначальное число разбиений n0 , если егоопределять по формуле n0 = (b − a)/h0 , скорее всего окажется не целым числом. Число разбиений n по своему смыслу на каждой итерацииl должно быть целым, поэтому вначале надо задавать число разбиений,а затем вычислять шаг, соответствующий данному числу разбиений.Это можно сделать следующим образом:b−ab−an0 = √+ 1, h0 =n0εдля формул средних прямоугольников и трапеции;b−ab−a√n0 =+1,h=02n024ε(1.3.2)(1.3.3)для формулы Симпсона.В этих формулах квадратные скобки [ ] обозначают целую частьзаключенного в них числа.3) дайте оценку сверху погрешности вычислений, используя формулы, выражающие Rn через соответствующие производные подынтегральнойфункции;4) оцените погрешность как разность между точным значением интегралаи значением, полученным численным методом;5) сравните между собой погрешности, полученные в п.п.
3 и 4;22Глава 1. Численные методы вычисления определенного интеграла6) оформите отчет по лабораторной работе. Отчет должен содержать описание использованного метода, результаты и текст программы.Варианты лабораторной работы и ответы представлены в таблице 1.2.Таблица 1.1. Координаты узловых точек и весовые коэффициенты квадратурной формулы ГауссаЧисло Номерузлов n точки i12345Координататочки xiКоэффициентqi1021x1 = −x2120.577350269211x1 = −x3q1 = q 3200.888888888930.77459666920.55555555561x1 = −x4q1 = q 42x2 = −x3q2 = q 330.33998104360.652145154940.86113631160.34785484511x1 = −x5q1 = q 52x2 = −x4q2 = q 3300.568888888940.53846931010.478628670550.90617984590.236926885123Глава 1. Численные методы вычисления определенного интегралаТаблица 1.2.
Варианты лабораторной работы№вар.ab12345101ex + 1e2012x + 1/ ln 22/ ln 23013x + 1/ ln 33/ ln 340.10.1 · eln(10 · x)0.150.20.2 · eln(5 · x)0.2612ex + 1/xe(e − 1) + ln 2701x · ex181ex2 + 16/x(e3 − 1)/3 + 169012x − e−x1/e10122x + 1/x3 + ln 211123x2 + 1/x7 + ln 212014x3 − e−x1/e13012x + exe14011/(1 + x2 )π/415011 − 2xe−xФункция f (x)224Ответ1/eГлава 1. Численные методы вычисления определенного интегралаОкончание таблицы 1.2.1234516012xexe−117011 − xe−x2/e181eln2 x/x1/31901x/(1 + x4 )2012e1/x /x2π/8√e− e21ln 22 ln 21/(ex − 1)ln(3/2)220π/2cos3 x · sin(2x)2/5230π/2(x + sin x)/(1 + cos x)π/224121/(x + x2 )ln(4/3)250π/2ex · cos x(eπ/2 − 1)/22601ex+e270.50.5 · eln(2x)1/228014x1/ ln 429015x + 1/ ln 55/ ln 5300110x + 1/ ln 1010/ ln 10225xee − eГлава 2Приближенное вычисление двойного интегралаЦель работы — изучение и применение численных методов для приближенного вычисления двойного интеграла.Продолжительность работы — 2–4 час.2.1.
Постановка задачиТребуется вычислить двойной интегралZZI=f (x, y) dx dy,(2.1.1)Dгде f (x, y) — непрерывная в области D функция двух переменных x и y.На практике редко удается выразить интеграл через элементарные функции и найти его точное значение. Поэтому обычно для вычисления интегралов применяются методы численного интегрирования. Они основаны назамене подынтегральной функции f (x, y) аппроксимирующими ее функциями, интегралы от которых легко вычисляются в элементарных функциях. Вкачестве аппроксимирующих функций, например, можно использовать многочлены.2.2. Численные методы вычисления двойного интегралаРассмотрим два способа численного интегрирования: метод ячеек и последовательное интегрирование.26Глава 2.
Приближенное вычисление двойного интеграла2.2.1. Метод ячеекПусть сначала область интегрирования является прямоугольником D == {a 6 x 6 b, c 6 x 6 d}. Среднее значение f (x, y) непрерывной в областиD функции f (x, y) по теореме о среднем представляется выражениемZZ1f (x, y) dx dy, S = (b − a) · (d − c).(2.2.1)f (x, y) =SDСчитая, что среднее значение приближенно равняется значению функциив центре прямоугольника: f (x, y) ≈ f (x, y), где x = (a + b)/2, y = (c + d)/2;из соотношения (2.1.1) получаем простейшую формулу для приближенноговычисления двойного интегралаZZf (x, y) dx dy ≈ S · f (x, y).(2.2.2)DНайдем погрешность формулы (2.2.2). Функцию f (x, y) будем считать достаточно гладкой, т.е. будем полагать, что она имеет все необходимые по ходурассуждения непрерывные производные. Разложим функцию f (x, y) по формуле Тейлора, выбирая центр прямоугольника (точку (x, y)) за точку разложенияf (x, y) = f (x, y) + (x − x) · f 0 x (x, y) + (y − y) · f 0 y (x, y)+1+ (x − x)2 · f 00 xx (x, y) + (x − x)(y − y) · f 00 xy (x, y)+21+ (y − y)2 · f 00 yy (x, y) + ...2(2.2.3).Погрешность есть разность точного и приближенного значений интеграла.27Глава 2.
Приближенное вычисление двойного интегралаПодставляя в (2.2.2) формулу (2.2.3), получим главный член погрешностиZb Zdf (x, y) dx dy − S · f (x, y) ≈R=a(2.2.4)c1· S · [(b − a)2 · f 00 xx (x, y) + (d − c)2 · f 00 yy (x, y)],24где члены, отброшенные при замене точного равенства приближенным, со≈держат производные старших порядков и более высокие степени длин сторонпрямоугольника D. Заметим, что все члены разложения, являющиеся нечетными функциями относительно центра прямоугольника, не вносят вклад впогрешность, поскольку интегралы от этих членов оказываются равными нулю.В общем случае длины сторон прямоугольника (b − a) и (d − c) не малы,поэтому главный член погрешности (2.2.4) может быть велик. Для повышения точности вычислений в области D (рис. 2.1) вводится сетка xi = a + ih1 ,yj = a + jh2 , i = 0, 1, ..., m, j = 0, 1, ..., n; h1 = (b − a)/m, h2 = (c − d)/n сдостаточно мелкими ячейками ∆Dij = {xi−1 6 x 6 xi , yj−1 6 y 6 yj , i =1, ..., m; j = 1, ..., n}.Вычисляя интеграл по каждой ячейке по формуле (2.2.2) и суммируя найденные значения по всем ячейкам, получаем формулу метода ячеекm XnX(2.2.5)I ≈ Ih =Sij · f (xi , y j ),i=1 j=1где Sij = h1 · h2 — площадь ячейки, xi = (xi−1 + xi )/2, y j = (yj−1 + yj )/2— координаты центра ячейки.
Здесь и далее пусть Ih будет приближеннымзначением интеграла (2.1.1), вычисленное по формуле (2.2.5) с шагами h1 иh2 .Справа в выражении (2.2.5) стоит интегральная сумма, поэтому для любой непрерывной функции f (x, y) эта сумма сходится к значению интеграла,28Глава 2. Приближенное вычисление двойного интегралаy( xi , yj )h1yn=dyjh2yj -1y0=cx0=a0xi -1xixm=bxРис. 2.1.когда периметры всех ячеек стремятся к нулю.Погрешность интегрирования (2.2.4) для одной ячейки ∆Dij представляется в видеRij ≈1Sij h21 · f 00 xx (xi , y j ) + h22 · f 00 yy (xi , y j ) .24(2.2.6)Суммируя выражения (2.2.6) по всем ячейкам, получаем погрешность метода ячеекR ≈ c1 · h21 + c2 · h22 ,(2.2.7)где1c1 =24ZZ1f 00 xx (x, y) dx dy, c2 =24DZZf 00 yy (x, y) dx dyDилиR = O(h21 + h22 ),(2.2.8)т.е.