Министерство Высшего Образования РФ.
Московский Институт Электронной Техники
(Технический Университет)
Лицей №1557
КУРСОВАЯ РАБОТА
“Вычисление интеграла методом
Ньютона-Котеса”
Написал: Коноплев А.А.
Проверил: доцент Колдаев В.Д.
Москва, 2001г.
1.
Введение..................................................................
................... 3 2. Теоретическая часть...................................................................4 3. Алгоритм работы....................................................................
....8 4. Код программы.................................................................
........17 . Модуль
K_graph............................................................17 . Модуль
Graphic.............................................................34 . Модуль
K_unit...............................................................38 . Основная программа....................................................40 5. Тестовые испытания.................................................................
42 6. Полезные советы по работе с программой.............................42 7. Окна ввода и вывода программы............................................. 8.
Вывод.....................................................................
.....................43 9. Список литературы................................................................
...44
Математика - одна из самых древних наук. Труды многих ученых вошли в мировой фонд и стали основой современных алгебры и геометрии. В конце XVII в., когда развитие науки шло быстрыми темпами, появились понятия дифференцирование, а вслед за ним и интегрирование. Многие правила нахождения неопределенного интеграла в то время не были известны, поэтому ученые пытались найти другие, обходные пути поиска значений. Первым методом явился метод Ньютона – поиск интеграла через график функции, т.е. нахождение площади под графиком, методом прямоугольников, в последствии усовершенствованный в метод трапеций. Позже был придуман параболический метод или метод Симпсона. Однако часть ученых терзал вопрос: А можно ли объединить все эти методы в один?? Ответ на него был дан одновременно двумя математиками Ньютоном и Котесом. Они вывели общую формулу, названную в их честь. Однако их метод был частично забыт. В этой работе будут изложены основные положения теории, рассмотрены различные примеры, приведены таблицы, полученные при различных погрешностях, и конечно описана работа и код программы, рассчитывающей интеграл методом Ньютона-Котеса.
Пусть некоторая функция f(x) задана в уздах интерполяции: [pic] (i=1,2,3…,n) на отрезке [а,b] таблицей значений:
|X0=a |X1 |X2 |… |XN=b | |Y0=f(x0) |Y1=f(x1) |Y2=f(x2) |… |YN=f(xN) |
Требуется найти значение интеграла [pic] . Для начала составим интерполяционный многочлен Лагранджа:
[pic]
Для равноотстоящих узлов интерполяционный многочлен имеет вид:
где q=(x-x0)/h – шаг интерполяции, заменим подынтегральную функцию f(x) интерполяционным многочленом Лагранжа:
Поменяем знак суммирования и интеграл и вынесем за знак интеграла постоянные элементы:
Так как dp=dx/h, то, заменив пределы интегрирования, имеем:
Для равноотстоящих узлов интерполяции на отрезке [a,b] величина шаг определяется как h=(a-b)/n. Представив это выражение для h в формулу (4) и вынося (b-a) за знак суммы, получим:
Положим, что
[pic] где i=0,1,2…,n; Числа Hi называют коэффициентами Ньютона-Котеса. Эти коэффиценты не зависят от вида f(x), а являются функцией только по n. Поэтому их можно вычислить заранее. Окончательная формула выглядит так:
[pic] Теперь рассмотрим несколько примеров.
Пример 1.
Вычислить с помощью метода Ньютона-Котаса: [pic] , при n=7.
Вычисление.
1) Определим шаг: h=(7-0)/7=1.
2)Найдем значения y:
|x0=0 |y0=1 | |x1=1 |y1=0.5 | |x2=2 |y2=0.2 | |x3=3 |y3=0.1 | |x4=4 |y4=0.0588 | |x5=5 |y5=0.0384 | |x6=6 |y6=0.0270 | |x7=7 |y7=0.02 |
3) Находим коэффициенты Ньютона-Котеса: H1=H7=0.0435, H1=H6=0.2040, H2=H5=0.0760 ,H3=H4=0.1730 Подставим значения в формулу и получим:
При подсчете с помощью формулы Ньютона-Лейбница получим:
Пример 2. Вычислить при помощи метода Ньютона-Котеса [pic] , взяв n=5; Вычисление: 1) Определим шаг h=(8-4)/5=0.8 2) Найдем значения y:
|x0=0 |y0=-2.61 | |x1=4.8 |y1=0.42 | |x2=5.6 |y2=4.34 | |x3=6.4 |y3=6.35 | |x4=7.2 |y4=4.38 | |x5=8 |y5=-0.16 |
3) Находим коэффициенты Ньютона –Котеса: H0=H5=0.065972 ;H1=H4=0.260417 ;H2=H3=0.173611 ; 4)Подставим значения в формулу и получим: [pic]
Рассмотрим частные случаи формулы Ньйтона-Котеса. Пусть n=1 тогда H0=H1=0.5 и конечная формула примет вид: [pic] Тем самым в качестве частного случая нашей формулы мы получили формулу трапеций. Взяв n=3, мы получим [pic] . Частный случай формулы Ньютона –Котеса – формула Симпсона
Теперь произведем анализ алгоритма и рассмотрим основной принцип работы программы.
Для вычисления интеграла сначала находятся коэффициенты Ньютона- Котеса. Их нахождение осуществляется в процедуре hkoef. Основной проблемой вычисления коэффициентов является интеграл от произведения множителей. Для его расчета необходимо:
А) посчитать коэффициенты при раскрытии скобок при q (процедура mnogoclen) Б) домножить их на 1/n , где n –степень при q (процедура koef) В) подставить вместо q значение n (функция integral)
Далее вычисляем факториалы (функция faktorial) и перемножаем полученные выражения (функция mainint). Для увеличения быстроты работы вводится вычисление половины от количества узлов интерполяции и последующей подстановкой их вместо неподсчитанных.
Процедура koef(w: массив;n:целый;var e:массив);
Процедура hkoef(n:целый;var h:массив);
Процедура mnogochlen(n,i:целые;var c:массив );
Процедура funktia(n:целая;a,b:вещест.;var y:массив;c:вещест.;f:строка);
Функция facktorial(n:целый):двойной;
Функция integral(w:массив;n:целый):двойной;
Функция mainint(n:целый;a,b:вещест.;y:массив):двойной;
Основная программа
Программа состоит из 8 файлов: . K_main.exe – файл загрузки основной программы . K_unit.tpu – модуль вычислительных процедур и функций . K_graph.tpu – модуль графических процедур . Graphic.tpu – модуль процедур для построения графика . Egavga.bgi – файл графической инициализации . Sans.chr, litt.chr – файлы шрифтов . Keyrus.com (не обязательно) – файл установки русского языка. Для работы программы с русским интерфайсом желательно запускать ее в режиме DOS.
================================================
==========МОДУЛЬ GRAPH==========
================================================ {$N+} unit k_graph; interface uses crt,graph,k_unit,graphic; procedure winwin1; procedure proline(ea:word); procedure winwwodab(ea:word); procedure error1(ea:word); procedure helpwin(ea:word); procedure error(ea:word); procedure newsctext(ea:word); procedure newsc(ea:word); procedure win1(ea:word); procedure win2(ea:word;var k:word); procedure wwodn(ea:word;var n:integer); procedure wwodab(ea:word;var a,b:real); procedure wwod1(ea:word;var y:array of double;var n:integer;var a,b:real); procedure wwod2(ea:word;var ea1:word;var n:integer;var a,b:real;var st:string); procedure win3(ea:word;n:integer;a,b:real;int:double;f:string;h:array of double;var k:word); implementation procedure proline(ea:word); {Проседура полосы процесса} var i:integer; f:string; c:char; begin newsc(ea); setcolor(15); setfillstyle(1,7); bar(160,150,460,260); rectangle(165,155,455,255); rectangle(167,157,453,253); case (ea mod 2) of
0: outtextxy(180,170,' Идет работа .Ждите..');
1: outtextxy(180,170,' Working.Please wait..'); end; setfillstyle(1,12); setcolor(0); rectangle(200,199,401,221); for i:=1 to 9 do line(200+i*20,200,200+i*20,220); delay(20000); for i:=1 to 100 do begin if ((i-1) mod 10)=0 then line(200+((i-1) div 10)*20,200,200+((i-1) div 10)*20,220); bar(round(200+2*(i-0.5)),200,200+2*i,220); delay(1100); setcolor(15); setfillstyle(1,7); bar(280,230,323,250); str(i,f); f:=f+'%'; outtextxy(290,235,f); if (i mod 25) =0 then bar(170,180,452,198); if (ea mod 2)=0 then case (i div 25) of
0: outtextxy(170,190,'Подготовка ');
1: outtextxy(170,190,'Расчет коеффициентов в многочлене');
2: outtextxy(170,190,'Расчет коеффициентов Ньютона-Котеса');
3: outtextxy(170,190,'Расчет интеграла'); end else case (i div 25) of
0: outtextxy(170,190,'Prepearing');
1: outtextxy(170,190,'Calculation of mnogochlen coeff.');
2: outtextxy(170,190,'Calculation of Newton-Cotes coeff. ');
3: outtextxy(170,190,'Calculation of integral'); end; setfillstyle(1,12); setcolor(0); end; end; procedure winwwodn(ea:word); {Окно ввода числа узлов интерполяции} var c:char; f:string; begin helpwin(ea); if (ea mod 2) =0 then begin outtextxy(360,140,' В этом окне необходимо '); outtextxy(360,155,' ввести количество узлов '); outtextxy(360,170,' интерполяции, от которого '); outtextxy(360,185,' будет зависить точность '); outtextxy(360,200,' вычисления интеграл и '); outtextxy(360,215,' количество зн чений функции.'); outtextxy(360,240,' ВНИМАНИЕ : НАСТОЯТЕЛЬНО '); outtextxy(360,250,' РЕКОМЕНДУЕТСЯ НЕ ВВОДИТЬ '); outtextxy(360,260,' ЗНАЧЕНИЕ N БОЛЬШЕ 12 !! '); end else begin outtextxy(360,140,' In this window you have to '); outtextxy(360,155,' put into the number. '); outtextxy(360,170,' The accuracy of calculation '); outtextxy(360,185,' and the number of function '); outtextxy(360,200,' parameters will depend on '); outtextxy(360,215,' this number. '); outtextxy(360,240,' WARNING: IT IS HARDLY '); outtextxy(360,250,' RECOMENDED NOT TO PUT IN '); outtextxy(360,260,' NUMBER MORE THEN 12 !! '); end; setcolor(2); setfillstyle(1,14); bar(70,200,340,300); rectangle(75,205,335,295); rectangle(77,207,333,293); if (ea mod 2) =0 then begin outtextxy(90,227,'Введите количество узлов(n):'); outtextxy(80,270,'ВНИМАНИЕ: При больших n возможна'); outtextxy(80,280,'некорректная работа компьютера!!'); end else begin outtextxy(80,217,'Put in number of'); outtextxy(80,227,' interpolation units:'); outtextxy(80,270,'WARNING:if you use big number '); outtextxy(80,280,'of units,PC wont work properly!'); end; setfillstyle(1,0); bar(190,240,230,255); end; procedure wwodn(ea:word;var n:integer); {Процедура ввода узлов n} var ec,p:integer; k,f:string; x:integer; c:char; begin newsc(ea); winwwodn(ea); repeat repeat winwwodn(ea); gotoxy(25,16); read(k); val(k,p,ec); if ec0 then begin error1(ea); readln; end; until ec=0; n:=p; if n>12 then begin if keypressed then c:=readkey; c:='r'; setcolor(15); setfillstyle(1,12); bar(140,210,490,300); rectangle(145,215,485,295); rectangle(147,217,483,293); if (ea mod 2) =0 then begin outtextxy(150,227,' Предупреждение!'); outtextxy(150,237,' Вы дейcтвительно хотите использовать'); outtextxy(150,250,' большое значение N ???'); end else begin outtextxy(150,227,' Warning!! '); outtextxy(150,237,' Do you realy want to use a big '); outtextxy(150,250,' number interpolation units(N)??? '); end; sound(600); delay(4000); nosound; setfillstyle(1,2); bar(320,260,350,280); setfillstyle(1,12); bar(250,260,280,280); repeat if keypressed then begin c:=readkey; if (c=#80) or (c=#72) or (c=#77) or (c=#75) then x:=x+1; setfillstyle(1,2); if (x mod 2)=0 then begin bar(250,260,280,280); setfillstyle(1,12); bar(320,260,350,280); end else begin bar(320,260,350,280); setfillstyle(1,12); bar(250,260,280,280);
Страницы: 1, 2, 3
