Численные методы. Номер 1.. 1. Приближенное вычисление определенных интегралов приближенное вычисление интеграла по формулам трапеций и Симпсона, оценка погрешностей вычисления Формула трапеций имеет вид
 Скачать 4.47 Mb.
|
|
2 ТЕМА №1. ПРИБЛИЖЕННОЕ ВЫЧИСЛЕНИЕ ОПРЕДЕЛЕННЫХ ИНТЕГРАЛОВ Приближенное вычисление интеграла по формулам трапеций и Симпсона, оценка погрешностей вычисления Формула трапеций имеет вид f x dx h f a f b f a ih a b i n ( ) ( ) ( ) ( ) 2 1 1 (1.1), где h b a n (1.2). Формула Симпсона имеет вид 1 1 1 ) ( 1 3 ) ( ) ( 3 ) ( n i k b a ih a f b f a f h dx x f (1.3) Определение приближенного значения интеграла свелось к вычислению суммы. Все это просто реализовать в табличном процессоре Microsoft Excel, используя приемы «автозаполнение», копирование формул и «автосуммирование». Задание 1: Вычислить приближенное значение интеграла 2 1 2 1 1 dx x с числом разбиений интервала интегрирования на 10 равных частей по формулам трапеций и Симпсона в табличном процессоре Microsoft Excel ив пакете MathCAD. Сравнить полученные значения интегралов сточным, определенным по формуле Ньютона-Лейбница (используя первообразную функцию. Решение вычисляем приближенное значение интеграла по формуле трапеций (1.1) при n=10 в Microsoft Excel. В первый столбец заносим значения k от нуля до числа разбиений интервала интегрирования n, те. десять в рассматриваемом случае. В соседнем столбце вычисляем значения аргумента (x) в точках разбиения интервала интегрирования, те. числа от единицы – нижнего предела интегрирования до двух – верхнего предела интегрирования. В третьем столбце вычисляем значения подынтегральной функции в точках разбиения интервала интегрирования. Приближенное значение интеграла получается суммированием значений подынтегральной функции за вычетом полусуммы значений функции на нижнем и верхнем пределах интегрирования с последующим умножением результата на величину h. Выполняя вычисления в каждом столбце, рационально использовать прием копирования формул (рис. 1.1). 1 2 3 4 5 6 7 8 9 10 11 12 13 A B C D x f(x) h 0,1 1 0,5 1,1 0,452489 1,2 0,409836 1,3 0,371747 1,4 0,337838 1,5 0,307692 1,6 0,280899 1,7 0,257069 1,8 0,235849 1,9 0,21692 2 0,2 Рис. 1.1. Таблица Microsoft Excel c вычислением интеграла по формуле трапеций Табличка с решением в Microsoft Excel в режиме отображения формул приведена на рис. 1.2. 3 1 2 3 4 5 6 7 8 9 10 11 12 13 A B C D x f(x) h =(2-1)/10 1 =1/(1+A2^2) 1,1 =1/(1+A3^2) 1,2 =1/(1+A4^2) 1,3 =1/(1+A5^2) 1,4 =1/(1+A6^2) 1,5 =1/(1+A7^2) 1,6 =1/(1+A8^2) 1,7 =1/(1+A9^2) 1,8 =1/(1+A10^2) 1,9 =1/(1+A11^2) 2 =1/(1+A12^2) =(СУММ(B2:B12)-B2/2-B12/2)*D1 Рис. 1.2. Табличка с решением задания 1 в режиме отображения формул Вычисление по формуле Симпсона в Microsoft Excel выполняется аналогично вычислению по формуле трапеций. Решение приведено на рис. 1.3 и 1.4. 1 2 3 4 5 6 7 8 9 10 11 12 13 A B C D E k x f(x) h 0,1 0 1 0,5 0,5 1 1,1 0,452488688 1,809954751 2 1,2 0,409836066 0,819672131 3 1,3 0,371747212 1,486988848 4 1,4 0,337837838 0,675675676 5 1,5 0,307692308 1,230769231 6 1,6 0,280898876 0,561797753 7 1,7 0,257069409 1,028277635 8 1,8 0,235849057 0,471698113 9 1,9 0,21691974 0,867678959 10 2 0,2 0,2 11 0,322033919 Рис. 1.3. Таблица Microsoft Excel c вычислением интеграла по формуле Симпсона 1 2 3 4 5 6 7 8 9 10 11 12 13 A B C D E k x f(x) h =(B12-B2)/A12 0 1 =1/(1+B2^2) =C2 1 1,1 =1/(1+B3^2) =C3*(3+(-1)^(A3+1)) 2 1,2 =1/(1+B4^2) =C4*(3+(-1)^(A4+1)) 3 1,3 =1/(1+B5^2) =C5*(3+(-1)^(A5+1)) 4 1,4 =1/(1+B6^2) =C6*(3+(-1)^(A6+1)) 5 1,5 =1/(1+B7^2) =C7*(3+(-1)^(A7+1)) 6 1,6 =1/(1+B8^2) =C8*(3+(-1)^(A8+1)) 7 1,7 =1/(1+B9^2) =C9*(3+(-1)^(A9+1)) 8 1,8 =1/(1+B10^2) =C10*(3+(-1)^(A10+1)) 9 1,9 =1/(1+B11^2) =C11*(3+(-1)^(A11+1)) 10 2 =1/(1+B12^2) =C12 11 =(СУММ(C2:C12)-C2/2-C12/2)*E1 =СУММ(E2:E12)*E1/3 Рис. 1.4. Табличка с решением в режиме отображения формул по методу Симпсона Замечаем, что полученные значения несколько разнятся. В задании 1 взят интеграл, который можно вычислить точно, так что точность полученных результатов можно проверить, используя формулу Ньютона-Лейбница. ) 1 ( ) 2 ( 1 1 2 1 2 arctg arctg dx x 4 На рис. 1.5 - 1.6 приведены таблички с вычислением точного значения интеграла с помощью первообразной. точное значение интеграла Рис. 1.5. Таблица Microsoft Excel c вычислением точного значения интеграла точное значение интеграла =ATAN(2)-ATAN(1) Рис. 1.6. Табличка с вычислением точного значения интеграла в режиме отображения формул Как видим, более точно интеграл вычисляется по формуле Симпсона. Вычислим интеграл в пакете MathCAD. Вычисление интегралов в пакете MathCAD выполняется с использованием расширенных операторов пакета (палитра Calculus) и показано на рис. 1.7. Рис. 1.7. Вычисление интеграла в пакете MathCAD Отмечаем, что в пакете MathCAD интеграл вычислен также точно. Обычно результаты вычислений в MathCAD выводятся стремя значащими цифрами. Для вычислений в рамках данного примера этого недостаточно. Настройка на отображение нужного количества знаков после запятой в MathCAD находится в пункте меню «Format». Вычисление интеграла с заданной точностью Задание 2. Вычислить интеграл 6 2 2 1 2 6 0 x dx по формулам трапеций и Симпсона с точностью Решение в пакете Microsoft Excel. Для вычисления приближенного значения интеграла b a dx x f ) ( необходимо определить необходимое для достижения заданной точности значение количества точек, входящих в формулы трапеций 5 f x dx h f a f b f a ih a b i n ( ) ( ) ( ) ( ) 2 1 1 , где h b и Симпсона ) ( 4 2 2 4 2 4 ) ( 3 ) ( 1 2 4 3 2 1 b f y y y y y y a f h dx x f n n b a 1. Точность приближенного значения интеграла по формуле трапеций задается формулой 001 0 ) ( max 12 2 3 x f n a b . Для определения необходимой величины n воспользуемся этой формулой. В рассматриваемом случае a=1,2; b=2,6; 2 6 0 1 x x f . Находим первую и вторую производные подынтегральной функции 3 2 3 2 6 0 2 6 0 1 2 1 x x x x x f , 5 2 2 3 2 5 2 3 2 6 0 3 6 0 1 2 6 0 1 2 3 6 0 1 x x x x x x x x f 5 2 2 5 2 2 2 6 0 2 6 0 6 0 3 6 0 x x x x x . Максимальное значение второй производной на промежутке интегрирования можно оценить как максимальное значение числителя дробите. при x=2.6 и минимальным значением знаменателя дробите. при x=1.6 5 2 2 2 1 6 0 6 2 2 6 0 ) ( max x f . Подставляя все величины в формулу для определения погрешности вычисления, имеем 001 0 2 1 6 0 6 2 2 6 0 12 3 1 6 2 5 2 2 2 3 n . Из этого неравенства определим значения n 2 5 2 3 3 2 1 6 0 6 2 2 6 0 001 0 12 4 1 n , те. 0381375 497 2 n , n>22,294. Для расчета принимаем n=24. Интеграл вычисляем, как в первом задании (рис. 1.8). 6 Рис. 1.8. Вычисление интеграла с заданной точностью по формуле трапеций 2. Точность приближенного значения интеграла по формуле Симпсона задается формулой x f n a b IV max 180 4 5 . Для определения n, требуемого для достижения заданной точности, нужно оценить величину четвертой производной от подынтегральной функции. Выше вычислялась вторая производная, так что третья производная будет производной от второй 7 2 2 5 2 5 2 2 6 0 2 2 5 2 6 0 6 0 2 2 6 0 2 6 0 ) ( x x x x x x x x f . Упрощаем полученное выражение 2 2 7 2 2 6 0 10 6 0 4 6 0 1 ) ( x x x x x x f 7 2 2 6 0 6 4 5 x x x . Затем вычислим четвертую производную, как производную от третьей производной x x x x x x x f IV 2 2 7 6 0 6 4 5 6 0 18 4 5 9 2 3 7 2 2 . Приводя подобные члены, получаем 9 2 4 2 6 0 2 , 40 48 , 33 24 3 x x x x f IV Оценим максимальное значение полученной величины на промежутке изменения аргументах от 1.2 до 2.6 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 A B h= =(2,6-1,2)/24 x y 1,2 =1/(0,6+A3^2)^0,5 =A3+$B$1 =1/(0,6+A4^2)^0,5 =A4+$B$1 =1/(0,6+A5^2)^0,5 =A5+$B$1 =1/(0,6+A6^2)^0,5 =A6+$B$1 =1/(0,6+A7^2)^0,5 =A7+$B$1 =1/(0,6+A8^2)^0,5 =A8+$B$1 =1/(0,6+A9^2)^0,5 =A9+$B$1 =1/(0,6+A10^2)^0,5 =A10+$B$1 =1/(0,6+A11^2)^0,5 =A11+$B$1 =1/(0,6+A12^2)^0,5 =A12+$B$1 =1/(0,6+A13^2)^0,5 =A13+$B$1 =1/(0,6+A14^2)^0,5 =A14+$B$1 =1/(0,6+A15^2)^0,5 =A15+$B$1 =1/(0,6+A16^2)^0,5 =A16+$B$1 =1/(0,6+A17^2)^0,5 =A17+$B$1 =1/(0,6+A18^2)^0,5 =A18+$B$1 =1/(0,6+A19^2)^0,5 =A19+$B$1 =1/(0,6+A20^2)^0,5 =A20+$B$1 =1/(0,6+A21^2)^0,5 =A21+$B$1 =1/(0,6+A22^2)^0,5 =A22+$B$1 =1/(0,6+A23^2)^0,5 =A23+$B$1 =1/(0,6+A24^2)^0,5 =A24+$B$1 =1/(0,6+A25^2)^0,5 =A25+$B$1 =1/(0,6+A26^2)^0,5 =A26+$B$1 =1/(0,6+A27^2)^0,5 сумма =СУММ(B3:B27) интеграл =B1*(B28-B3/2-B27/2) 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 A B h= 0,058333 x y 1,2 0,70014 1,258333 0,676757 1,316667 0,654615 1,375 0,633645 1,433333 0,613781 1,491667 0,594957 1,55 0,57711 1,608333 0,560179 1,666667 0,544107 1,725 0,52884 1,783333 0,514326 1,841667 0,500517 1,9 0,48737 1,958333 0,474843 2,016667 0,462896 2,075 0,451495 2,133333 0,440605 2,191667 0,430196 2,25 0,420239 2,308333 0,410706 2,366667 0,401574 2,425 0,392818 2,483333 0,384418 2,541667 0,376353 2,6 0,368605 сумма 12,60109 интеграл 0,703892 7 1130 2 1 6 0 6 2 2 , 40 6 2 48 33 24 3 max 9 2 4 2 x f IV . Получим необходимое значение n для достижения требуемой точности приближенного значения интеграла 1130 180 2 , 1 6 , 2 001 , 0 4 5 n . Решаем неравенство относительно n 623 , 7 1130 001 , 0 180 2 , 1 6 , 2 5 4 n . Выполняем вычисления с n=8 (рис. 1.9, 1.10). Рис. 1.9. Вычисление интеграла с заданной точностью по формуле Симпсона (режим отображения чисел) Рис. 1.10. Вычисление интеграла с заданной точностью по формуле Симпсона (режим отображения формул) Как видим с точностью три знака после запятой значения интеграла совпадают Вычислить центр тяжести плоской фигуры Задание 3. Вычислить центр тяжести плоской фигуры (х с , ус, ограниченной кривой 4 1 3 ) ( 2 x x f и осями координат x>0 и y>0 по формулам 2 0 ) ( 1 dx x f x S x c , 2 0 2 ) ( 2 1 dx x f S y c , где 2 0 ) ( dx x f S . Решение. Для вычисления центра тяжести нужно вычислить значения трех интегралов. Интегралы вычисляем приближенно методом трапеций, как в задании 1. В первый столбец таблицы заносим значения аргументах, во втором вычисляем значения 1 2 3 4 5 6 7 8 9 10 11 12 13 C D E F h= 0,175 i x y 0 1,2 0,70014 0,70014 1 1,375 0,633645 2,534579 2 1,55 0,57711 1,15422 3 1,725 0,52884 2,115359 4 1,9 0,48737 0,97474 5 2,075 0,451495 1,80598 6 2,25 0,420239 0,840477 7 2,425 0,392818 1,571273 8 2,6 0,368605 0,368605 12,06537 0,703813 1 2 3 4 5 6 7 8 9 10 11 12 13 C D E F h= =(2,6-1,2)/8 i x y 0 1,2 =1/(0,6+D3^2)^0,5 =E3 СУММ 8 подынтегральной функции по формуле 4 1 3 ) ( 2 x x f . В третьем столбце вычисляем значения x f(x), в четвертом - f (x) 2 . Решение в Microsoft Excel приведено на рис. 1.11, 1.12. На рис. 1.11 в силу громоздкости формул вычисления интеграла при полной их идентичности не отображена формула вычисления интеграла 2 0 2 ) ( dx x f 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 B C D E x f(x) xf(x) f(x)f(x) 0 =3*(1-B20^2/4)^0,5 =B20*C20 =C20^2 0,2 =3*(1-B21^2/4)^0,5 =B21*C21 =C21^2 0,4 =3*(1-B22^2/4)^0,5 =B22*C22 =C22^2 0,6 =3*(1-B23^2/4)^0,5 =B23*C23 =C23^2 0,8 =3*(1-B24^2/4)^0,5 =B24*C24 =C24^2 1 =3*(1-B25^2/4)^0,5 =B25*C25 =C25^2 1,2 =3*(1-B26^2/4)^0,5 =B26*C26 =C26^2 1,4 =3*(1-B27^2/4)^0,5 =B27*C27 =C27^2 1,6 =3*(1-B28^2/4)^0,5 =B28*C28 =C28^2 1,8 =3*(1-B29^2/4)^0,5 =B29*C29 =C29^2 2 =3*(1-B30^2/4)^0,5 =B30*C30 =C30^2 =(СУММ(C20:C30)-C20/2-C30/2)*(B30-B20)/10 =(СУММ(D20:D30)-D30/2-D20/2)*(B30-B20)/10 =(СУММ(E20:E30)-E30/2-E20/2)*(B30-B20)/10 xc =D31/C31 yc =E31/(2*C31) Центр тяжести Рис. 1.13. Решение задания 3 в режиме отображения формул Рис. 1.14. Решение задания 3 в режиме отображения чисел Для проверки правдоподобности расчета результат вычислений можно проиллюстрировать графически. На основании значений аргумента и функции, расположенных в столбцах В и С, строим диаграмму типа точечная, состоящую из отдельных точек (рис. 1.15). Результат построения приведен на рис. 1.16. Далее нужно нанести положение вычисленного центра тяжести. Программа Microsoft Excel позволяет на график наносить дополнительные точки. Для этого делаем активной диаграмму, вызываем вспомогательное меню, из которого выбираем команду Исходные данные рис. 1.17). Выбираем пункт Добавить ив правые окна вписываем адреса ячеек, содержащих вычисленные координаты точки – центра тяжести (рис. 1.18). 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 Центр тяжести четверти эллипса x f(x) xf(x) f(x)f(x) 0 3 0 9 0,2 2,984962311 0,596992462 8,91 0,4 2,939387691 1,175755077 8,64 0,6 2,861817604 1,717090563 8,19 0,8 2,749545417 2,199636334 7,56 1 2,598076211 2,598076211 6,75 1,2 2,4 2,88 5,76 1,4 2,142428529 2,99939994 4,59 1,6 1,8 2,88 3,24 1,8 1,307669683 2,35380543 1,71 2 0 0 0 4,656777489 3,880151203 11,97 xc 0,833226671 yc 1,285223529 9 Рис. 1.15. Построение вида фигуры, центр тяжести которой определяется Рис. 1.16. Фигура, центр тяжести которой определяется Рис. 1.17. Вкладка Ряд команды Исходные данные 10 Рис. 1.18. Заполнение окна с адресом отображения точки Точка отображается пакетом Microsoft Excel на графике (рис. 1.19). Рис. 1.19. Результат построения графика с центром тяжести Решение задачи в пакете MathCAD выполняется с использованием расширенного оператора вычисления определенного интеграла (палитра Calculus) (рис. 1.20). 11 Рис. 1.20. Решение задания 3 в пакете MathCAD ВАРИАНТЫ ЗАДАНИЙ ПО ТЕМЕ 1 Задание 1. Задан интеграл Вычислить 1. приближенное значение по формуле трапеций (при числе разбиений промежутка интегрирования не менее 20) в Microsoft Excel; 2. приближенное значение по формуле Симпсона (при числе разбиений промежутка интегрирования не менее 20) в Microsoft Excel; 3. значение интеграла в пакете MathCAD; 4. точное значение интеграла через первообразную функцию. В отчете по выполнению задания привести формулы трапеций и Симпсона, по которым производился расчет график подынтегральной функции таблички Microsoft Excel с вычислениями в режимах отображения чисел и формул с сеткой и заголовками строки столбцов фрагмент листа MathCAD с решением. вид первообразной подынтегральной функции и точное значение интеграла погрешность приближенного значения интеграла. Вариант Задан интеграл x dx 2 0 4 1 Вариант Задан интеграл 1 9 2 Вариант 3 |