Главная Рефераты по рекламе Рефераты по физике Рефераты по философии Рефераты по финансам Рефераты по химии Рефераты по хозяйственному праву Рефераты по цифровым устройствам Рефераты по экологическому праву Рефераты по экономико-математическому моделированию Рефераты по экономической географии Рефераты по экономической теории Рефераты по этике Рефераты по юриспруденции Рефераты по языковедению Рефераты по юридическим наукам Рефераты по истории Рефераты по компьютерным наукам Рефераты по медицинским наукам Рефераты по финансовым наукам Рефераты по управленческим наукам Психология и педагогика Промышленность производство Биология и химия Языкознание филология Издательское дело и полиграфия Рефераты по краеведению и этнографии Рефераты по религии и мифологии Рефераты по медицине Рефераты по сексологии Рефераты по информатике программированию Краткое содержание произведений |
Реферат: Вычисление интеграла методом Ньютона-Котеса (теория и программа на Паскале)Реферат: Вычисление интеграла методом Ньютона-Котеса (теория и программа на Паскале)Министерство Высшего Образования РФ.Московский Институт Электронной Техники (Технический Университет) Лицей №1557КУРСОВАЯ РАБОТА“Вычисление интеграла методом Ньютона-Котеса” Написал: Коноплев А.А.Проверил: доцент Колдаев В.Д.Москва, 2001г.
Математика - одна из самых древних наук. Труды многих ученых вошли в мировой фонд и стали основой современных алгебры и геометрии. В конце XVII в., когда развитие науки шло быстрыми темпами, появились понятия дифференцирование, а вслед за ним и интегрирование. Многие правила нахождения неопределенного интеграла в то время не были известны, поэтому ученые пытались найти другие, обходные пути поиска значений. Первым методом явился метод Ньютона – поиск интеграла через график функции, т.е. нахождение площади под графиком, методом прямоугольников, в последствии усовершенствованный в метод трапеций. Позже был придуман параболический метод или метод Симпсона. Однако часть ученых терзал вопрос: А можно ли объединить все эти методы в один?? Ответ на него был дан одновременно двумя математиками Ньютоном и Котесом. Они вывели общую формулу, названную в их честь. Однако их метод был частично забыт. В этой работе будут изложены основные положения теории, рассмотрены различные примеры, приведены таблицы, полученные при различных погрешностях, и конечно описана работа и код программы, рассчитывающей интеграл методом Ньютона-Котеса.
Пусть некоторая функция f(x) задана в уздах интерполяции: (i=1,2,3…,n) на отрезке [а,b] таблицей значений:
Требуется найти значение интеграла . Для начала составим интерполяционный многочлен Лагранджа: Для равноотстоящих узлов интерполяционный многочлен имеет вид: где q=(x-x0)/h – шаг интерполяции, заменим подынтегральную функцию f(x) интерполяционным многочленом Лагранжа: Поменяем знак суммирования и интеграл и вынесем за знак интеграла постоянные элементы: Так как dp=dx/h, то, заменив пределы интегрирования, имеем:Для равноотстоящих узлов интерполяции на отрезке [a,b] величина шаг определяется как h=(a-b)/n. Представив это выражение для h в формулу (4) и вынося (b-a) за знак суммы, получим: Положим, чтогде i=0,1,2…,n; Числа Hi называют коэффициентами Ньютона-Котеса. Эти коэффиценты не зависят от вида f(x), а являются функцией только по n. Поэтому их можно вычислить заранее. Окончательная формула выглядит так: Теперь рассмотрим несколько примеров. Пример 1. Вычислить с помощью метода Ньютона-Котаса: , при n=7. Вычисление. 1) Определим шаг: h=(7-0)/7=1. 2)Найдем значения y:
3) Находим коэффициенты Ньютона-Котеса: H1=H7=0.0435, H1=H6=0.2040, H2=H5=0.0760 ,H3=H4=0.1730 Подставим значения в формулу и получим: П Пример 2. Вычислить при помощи метода Ньютона-Котеса , взяв n=5; Вычисление:
H0=H5=0.065972 ;H1=H4=0.260417 ;H2=H3=0.173611 ; 4)Подставим значения в формулу и получим: Рассмотрим частные случаи формулы Ньйтона-Котеса. Пусть n=1 тогда H0=H1=0.5 и конечная формула примет вид: Тем самым в качестве частного случая нашей формулы мы получили формулу трапеций. Взяв n=3, мы получим . Частный случай формулы Ньютона –Котеса – формула Симпсона
Теперь произведем анализ алгоритма и рассмотрим основной принцип работы программы. Для вычисления интеграла сначала находятся коэффициенты Ньютона-Котеса. Их нахождение осуществляется в процедуре 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 файлов:
Для работы программы с русским интерфайсом желательно запускать ее в режиме 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 ec<>0 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); END; end; if (ea mod 2) =0 then begin outtextxy(255,267,'ДА'); outtextxy(325,267,'НЕТ'); end else begin outtextxy(255,267,'YES'); outtextxy(325,267,'NO'); end; until c=#13; if abs(x mod 2)=1 then begin n:=0; setcolor(15); setfillstyle(1,2); bar(160,200,460,280); rectangle(165,205,455,275); rectangle(167,207,453,273); if (ea mod 2)=0 then begin outtextxy(180,227,'Для работы программы необходимо'); outtextxy(180,237,' заново ввести N.'); outtextxy(180,247,' Нажмите ENTER для продолжения.'); end else begin outtextxy(180,227,' To continue you have to '); outtextxy(180,237,' again put in N. '); outtextxy(180,247,' Press ENTER to continue.'); end; readln; readln; end; end; until n>0; end; procedure winwwodab(ea:word); {Окно ввода приделов интегрирования} var f:string; begin helpwin(ea); if (ea mod 2)=0 then begin outtextxy(360,140,' В этом окне необходимо'); outtextxy(360,155,' ввести сначала нижнее'); outtextxy(360,170,' значение интеграл и нажать'); outtextxy(360,185,' ENTER, а затем ввести'); outtextxy(360,200,' верхнее значение интеграла'); outtextxy(360,215,' и снова нажать ENTER.'); end else begin outtextxy(360,140,' In this window you have to:'); outtextxy(360,155,'firstly, put in lower value '); outtextxy(360,170,'of integral and press ENTER,'); outtextxy(360,185,'then put in higher value'); outtextxy(360,200,'of integral and press ENTER'); end; setcolor(2); setfillstyle(1,5); bar(10,210,335,320); rectangle(15,215,330,315); rectangle(17,217,328,313); settextstyle(0,0,0); if (ea mod 2)=0 then begin outtextxy(20,230,' Введите нижнее значение'); outtextxy(20,244,' интеграл :'); outtextxy(20,262,' Введите верхнее значение'); outtextxy(20,272,'интеграл :'); end else begin outtextxy(20,230,' Put in lower value of'); outtextxy(20,244,' integral:'); outtextxy(20,262,' Put in higher value of'); outtextxy(20,272,'integral:'); end; end; procedure wwodab(ea:word;var a,b:real); {Процедура ввода приделов интегрирования} var f:string; k:string; ec:integer; begin newsc(ea); winwwodab(ea); readln; repeat winwwodab(ea); gotoxy(16,16); read(k); val(k,a,ec); if ec<>0 then error1(ea); until ec=0; readln; repeat winwwodab(ea); str(a:4:2,f); outtextxy(120,244,f); gotoxy(16,18); read(k); val(k,b,ec); if ec<>0 then error1(ea); until ec=0; end; procedure helpwin(ea:word); {основа окна помощи} begin setfillstyle(1,3); bar(350,100,590,380); setcolor(0); rectangle(353,103,587,377); rectangle(355,105,585,375); setcolor(14); if (ea mod 2)=0 then outtextxy(360,115,' ОКНО ПОМОЩИ') else outtextxy(360,115,' HELP WINDOW'); end; procedure error1(ea:word); begin setcolor(15); setfillstyle(1,12); bar(140,210,490,280); rectangle(145,215,485,275); rectangle(147,217,483,273); if (ea mod 2)=0 then begin outtextxy(150,227,' Ошибка! '); outtextxy(150,237,' Вводимые параметр не число!! '); outtextxy(150,250,' Проверьте значение и заново введите его.'); end else begin outtextxy(150,227,' Error! '); outtextxy(150,237,' The value you entered isn`t a quantity!!'); outtextxy(150,250,' Check it and put it in again. '); end; sound(600); delay(4000); nosound; readln; readln; end; procedure error(ea:word); {Процедура ошибки} begin setcolor(15); setfillstyle(1,12); bar(140,210,490,260); rectangle(145,215,485,255); rectangle(147,217,483,253); if (ea mod 2)=0 then begin outtextxy(150,227,' Ошибка!'); outtextxy(150,237,' Недостаток вводимых параметров!!'); end else begin outtextxy(150,227,' Error!'); outtextxy(150,237,' Not all parameters are set!'); end; sound(600); delay(4000); nosound; readln; end; procedure newsctext(ea:word); {Текст для процедуры newsc} begin if ea mod 2 =0 then begin settextstyle(0,0,1); setcolor(15); outtextxy(400,440,'Язык - Русский. '); outtextxy(400,450,'Версия 1.0 Последнее издание'); outtextxy(400,460,'й Все права защищены.'); end else begin settextstyle(0,0,1); setcolor(15); outtextxy(400,440,'Language - English.'); outtextxy(400,450,'Version 1.0 Final release.'); outtextxy(400,460,'й All rights reserved.'); end; end; procedure newsc(ea:word); {Процедура обновления экрана} begin cleardevice; setfillstyle(10,8); floodfill(1,1,15); setcolor(0); setfillstyle(1,7); bar(80,10,580,80); rectangle(82,12,578,78); rectangle(85,15,575,75); settextstyle(0,0,2); setcolor(10); if ea mod 2 =0 then begin settextstyle(0,0,2); outtextxy(90,20,' Вычисление интеграл '); outtextxy(90,50,' методом Ньютона-Котеса.'); newsctext(ea); end else begin settextstyle(3,0,2); outtextxy(90,20,' Calculeting of integral'); outtextxy(90,47,' using the Newton-Cotes method.'); newsctext(ea); end; settextstyle(0,0,1); end; procedure winwin1; {Окно процедуры win1} begin setfillstyle(1,7); bar(160,110,460,380); setcolor(0); rectangle(162,113,457,377); rectangle(165,115,455,375); end; procedure win1(ea:word); {Вводное окно} begin settextstyle(0,0,1); setcolor(10); if (ea mod 2)=0 then begin outtextxy(168,135,'Министерство Высшего образования РФ); outtextxy(168,150,'Московский Государственный Институт'); outtextxy(168,160,' Электронной Техники '); outtextxy(168,170,' (Технический лниверситет) '); outtextxy(168,180,' Лицей №1557 '); outtextxy(168,210,' КУРСОВАЯ РАБО'А '); outtextxy(168,230,' «Вычисление интеграла '); outtextxy(168,245,' метедом Ньютона-Котеса» '); outtextxy(158,270,' Написал: Коноплев А.А. '); outtextxy(158,285,' Руководитель: доцент Колдаев В.Д.'); end else begin outtextxy(168,135,' Department of High Education '); outtextxy(168,150,' Moscow State Institute of '); outtextxy(168,160,' Electronic Technics '); outtextxy(168,170,' (Technics University) '); outtextxy(168,180,' Lyceum №1557 '); outtextxy(168,210,' COURSE WORK '); outtextxy(168,230,' «Calculation of integral '); outtextxy(168,245,' by Newton-Cotes method» '); outtextxy(158,270,' Author: Konoplev A.A. '); outtextxy(158,285,' Supervisor:senior lecturer '); outtextxy(158,300,' Koldaev V.D. '); end; end; procedure win2(ea:word;var k:word); {Окно выбора способа подсчета } var c:char; x:integer; f:string; begin setcolor(2); setfillstyle(1,5); bar(70,200,340,330); rectangle(75,205,335,325); rectangle(77,207,333,323); settextstyle(0,0,0); setfillstyle(1,15); bar(80,250,330,270); setfillstyle(1,5); bar(80,285,330,305); if ea mod 2 =0 then begin outtextxy(77,220,'Выбирете способ задания значений'); outtextxy(75,230,' функции. '); outtextxy(70,255,' По таблице(в ручную)'); outtextxy(70,295,' По расчетам(автом т.)'); end else begin outtextxy(77,220,' Choose a method of putting in'); outtextxy(75,230,' the values of function. '); outtextxy(70,255,' By the table(by hand)'); outtextxy(70,295,' By calculations(automat.)'); end; helpwin(ea); if ea mod 2 =0 then begin outtextxy(360,140,'В этом способе необходимо'); outtextxy(360,155,'самостоятельно вводить'); outtextxy(360,170,'значения функции.'); end else begin outtextxy(360,140,'In this method you have'); outtextxy(360,155,'to put in values of '); outtextxy(360,170,'function by yourself.'); end; x:=0; repeat if keypressed then begin c:=readkey; if (c=#80) or (c=#72) then x:=x+1; setfillstyle(1,15); if (x mod 2)=0 then begin bar(80,250,330,270); setfillstyle(1,5); bar(80,285,330,305); helpwin(ea); if ea mod 2 =0 then begin outtextxy(360,140,'В этом способе необходимо'); outtextxy(360,155,'самостоятельно вводить'); outtextxy(360,170,'значения функции.'); end else begin outtextxy(360,140,'In this method you have'); outtextxy(360,155,'to put in values of '); outtextxy(360,170,'function by yourself.'); end; end else begin bar(80,285,330,305); setfillstyle(1,5); bar(80,250,330,270); helpwin(ea); if ea mod 2 =0 then begin outtextxy(360,140,'В этом способе компьютер'); outtextxy(360,155,'сам вычесляет значения'); outtextxy(360,170,'функции по вводимой функции.'); end else begin outtextxy(360,140,'In this method PC will'); outtextxy(360,155,'automaticly count the value'); outtextxy(360,170,'of function by the function'); outtextxy(360,185,'you enter '); end; end; setcolor(2); if ea mod 2 =0 then begin outtextxy(70,255,' По таблице(в ручную)'); outtextxy(70,295,' По расчетам(автом т.)'); end else begin outtextxy(70,255,' By the table(by hand)'); outtextxy(70,295,' By calculations(automat.)'); end; end; until c=#13; k:=x mod 2; end; procedure wwod1(ea:word;var y:array of double;var n:integer;var a,b:real); {Окно ручного ввода функции} var i,p:integer; s,f:string; p1:real; c:char; begin wwodn(ea,n); if n=0 then wwodn(ea,n); newsc(ea); wwodab(ea,a,b); 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,'нажмите ENTER.'); end else begin outtextxy(360,140,'In this window you have'); outtextxy(360,155,'to gradually enter the'); outtextxy(360,170,'values of functions.'); outtextxy(360,185,'After each enter press'); outtextxy(360,200,'ENTER key.'); end; setfillstyle(1,9); bar(40,200,330,300); rectangle(45,205,325,295); rectangle(47,207,323,293); if ea mod 2 =0 then outtextxy(56,227,'Введите 0 -е значение финкции:') else outtextxy(56,227,' Enter 0 -th value of function:'); for i:=0 to n do begin setfillstyle(1,0); bar(137,250,180,273); gotoxy(19,17); setfillstyle(1,9); read(p1); y[i]:=p1; bar(120,227,134,240); str(i+1,s); outtextxy(120,227,s); bar(310,220,320,250); end; end; procedure wwod2(ea:word;var ea1:word;var n:integer;var a,b:real;var st:string); {Окно 2 меню автомат. подсчета} var i:integer; c,k:char; x:longint; f:string; begin repeat x:=-600000; if keypressed then c:=readkey; c:='t'; newsc(ea); setfillstyle(1,15); bar(70,120,342,330); setcolor(12); rectangle(75,125,337,325); rectangle(77,127,335,323); settextstyle(0,0,0); setfillstyle(1,11); bar(80,170,330,190); if ea mod 2 =0 then begin outtextxy(80,130,'Меню ввода параметров нахождения'); outtextxy(80,140,' интеграла'); outtextxy(80,180,' Ввести количество узлов(n)'); outtextxy(80,210,' Ввести приделы интегрирования'); outtextxy(80,240,' Ввести функцию'); outtextxy(80,270,' Считать интеграл'); outtextxy(80,300,' Выход '); end else begin outtextxy(80,130,'Menu of entering the parameters'); outtextxy(80,140,' of integral'); outtextxy(80,180,' Put in the number of units '); outtextxy(80,210,' Enter the bounds of integral'); outtextxy(80,240,' Enter function'); outtextxy(80,270,' Count integral'); outtextxy(80,300,' Exit '); end; helpwin(ea); if ea mod 2 =0 then begin outtextxy(360,140,' Нажмите Enter для'); outtextxy(360,155,' ввода количества узлов'); end else begin outtextxy(360,140,' Press Enter to put'); outtextxy(360,155,' in the number of units'); end; repeat if keypressed then begin c:=readkey; case c of #80: x:=x-1; #72: x:=x+1; end; setfillstyle(1,11); case (abs(x) mod 5) of 0: begin bar(80,170,330,190); setfillstyle(1,15); bar(80,200,330,220); bar(80,290,330,310); helpwin(ea); if ea mod 2 =0 then begin outtextxy(360,140,' Нажмите Enter для'); outtextxy(360,155,' ввода количества узлов'); end else begin outtextxy(360,140,' Press Enter to put'); outtextxy(360,155,'in the number of units.'); end; end; 1: begin bar(80,200,330,220); setfillstyle(1,15); bar(80,170,330,190); bar(80,230,330,250); helpwin(ea); if ea mod 2 =0 then begin outtextxy(360,140,' Нажмите ENTER для ввода'); outtextxy(360,155,'приделов интегрирования.'); end else begin outtextxy(360,140,' Press ENTER to put in'); outtextxy(360,155,'the bounds of integral.'); end; end; 2: begin bar(80,230,330,250); setfillstyle(1,15); bar(80,200,330,220); bar(80,260,330,280); helpwin(ea); if ea mod 2 =0 then begin outtextxy(360,140,' Нажмите ENTER для ввода'); outtextxy(360,155,'функции.'); end else begin outtextxy(360,140,' Press ENTER to enter'); outtextxy(360,155,'function.'); end; end; 3: begin bar(80,260,330,280); setfillstyle(1,15); bar(80,230,330,250); bar(80,290,330,310); helpwin(ea); if ea mod 2 =0 then begin outtextxy(360,140,' Нажмите ENTER для начала'); outtextxy(360,155,'подсчета самого интеграла.'); end else begin outtextxy(360,140,' Press ENTER to begin'); outtextxy(360,155,'integral calculations.'); end; end; 4: begin bar(80,290,330,310); setfillstyle(1,15); bar(80,260,330,280); bar(80,170,330,190); helpwin(ea); end; end; setcolor(12); if ea mod 2 =0 then begin outtextxy(80,130,'Меню ввода параметров нахождения'); outtextxy(80,140,' интеграла'); outtextxy(80,180,' Ввести количество узлов(n)'); outtextxy(80,210,' Ввести приделы интегрирования'); outtextxy(80,240,' Ввести функцию'); outtextxy(80,270,' Считать интеграл'); outtextxy(80,300,' Выход '); end else begin outtextxy(80,130,'Menu of entering the parameters'); outtextxy(80,140,' of integral'); outtextxy(80,180,' Put in the number of units '); outtextxy(80,210,' Enter the bounds of integral'); outtextxy(80,240,' Enter function'); outtextxy(80,270,' Count integral'); outtextxy(80,300,' Exit '); end; end; until c=#13; c:='t'; case (abs(x) mod 5) of 0: begin wwodn(ea,n); end; 1: wwodab(ea,a,b); 2: begin helpwin(ea); setcolor(15); setfillstyle(1,9); bar(70,200,340,300); rectangle(75,205,335,295); rectangle(77,207,333,293); if ea mod 2 =0 then begin outtextxy(86,227,'Введите функцию f(x):'); setcolor(14); outtextxy(360,140,' В этом окне необходимо'); outtextxy(360,155,' ввести саму функцию.'); outtextxy(360,200,'Примечание: 1.данная программа '); outtextxy(360,215,'распознает только '); outtextxy(360,230,'элементарные функции.'); outtextxy(360,245,'(x,cos(x) и др.)'); outtextxy(360,260,’2.При неправильном вводе’); outtextxy(360,275,’по умолчанию f(x)=x;’); outtextxy(360,275,’3.Если после нажатия ENTER’); outtextxy(360,275,’ничего не произошло, то outtextxy(360,275,’занововведите функцию.’); end else begin outtextxy(86,227,'Enter function f(x):'); setcolor(14); outtextxy(360,140,' In this window you have'); outtextxy(360,155,' to enter the function.'); outtextxy(360,200,'Note: This version of '); outtextxy(360,215,'programm can indentify only '); outtextxy(360,230,'simple functions, as'); outtextxy(360,245,'x,cos(x) and other.'); end; setfillstyle(1,0); bar(86,255,330,275); readln; gotoxy(13,17); read(st); writeln(st); readln; end; 3:if (n error(ea); 4: halt; end; until (n>0)and(a<>b)and(st<>'')and((abs(x) mod 5)=3); end; procedure win3(ea:word;n:integer;a,b:real;int:double;f:string;h:array of double;var k:word); {Последнее окно просмотра результатов} var i:integer; c:char; x:longint; p1,p:string; y:array[0..16] of double; begin funktia(n,a,b,y,1,f); f:='('+f+')'+'dx ='; repeat x:=-600000; newsc(ea); setfillstyle(1,2); bar(170,120,490,360); setcolor(14); rectangle(175,125,485,355); rectangle(177,127,483,353); settextstyle(0,0,0); setfillstyle(1,1); bar(180,170,480,190); if ea mod 2 =0 then begin outtextxy(180,135,Функция распознана.Интеграл подсчитан.'); outtextxy(180,180,' Посмотреть значение интеграла'); outtextxy(180,210,'Посмотреть коэффициенты Ньютона-Котеса'); outtextxy(180,240,' Посмотреть значения функции'); outtextxy(180,270,' Посмотреть график' ); outtextxy(180,300,' Считать снова'); outtextxy(180,330,' Выход '); end else begin outtextxy(180,135,'Function Indentified.Integral counted.'); outtextxy(180,180,' View value of integral'); outtextxy(180,210,' View Newton-Cotes coefficients'); outtextxy(180,240,' Veiw values of function'); outtextxy(180,270,' View graphik ' ); outtextxy(180,300,' Count again'); outtextxy(180,330,' Exit '); end; repeat if keypressed then begin c:=readkey; case c of #80: x:=x-1; #72: x:=x+1; end; setfillstyle(1,1); case (abs(x) mod 6) of 0: begin bar(180,170,480,190); setfillstyle(1,2); bar(180,200,480,220); bar(180,320,480,340); end; 1: begin bar(180,200,480,220); setfillstyle(1,2); bar(180,170,480,190); bar(180,230,480,250); end; 2: begin bar(180,230,480,250); setfillstyle(1,2); bar(180,200,480,220); bar(180,260,480,280); end; 3: begin bar(180,260,480,280); setfillstyle(1,2); bar(180,230,480,250); bar(180,290,480,310); end; 4: begin bar(180,290,480,310); setfillstyle(1,2); bar(180,260,480,280); bar(180,320,480,340); end; 5: begin bar(180,320,480,340); setfillstyle(1,2); bar(180,290,480,310); bar(180,170,480,190); end; end; if ea mod 2 =0 then begin outtextxy(180,135,'Функция распознана.Интеграл подсчитан.'); outtextxy(180,180,' Посмотреть значение интеграла'); outtextxy(180,210,'Посмотреть коэффициенты Ньютона-Котеса'); outtextxy(180,240,' Посмотреть значения функции'); outtextxy(180,270,' Посмотреть график ' ); outtextxy(180,300,' Считать снова'); outtextxy(180,330,' Выход '); end else begin outtextxy(180,135,'Function Indentified.Integral counted.'); outtextxy(180,180,' View value of integral'); outtextxy(180,210,' View Newton-Cotes coefficients'); outtextxy(180,240,' Veiw values of function'); outtextxy(180,270,' View graphik ' ); outtextxy(180,300,' Count again'); outtextxy(180,330,' Exit '); end; end; until c=#13; c:='t'; case (abs(x) mod 6) of 0:begin setcolor(15); setfillstyle(1,12); bar(140,200,490,280); rectangle(145,205,485,275); rectangle(147,207,483,273); settextstyle(2,0,1); setusercharsize(1,1,5,1); outtextxy(170,210,'S'); settextstyle(2,0,4); str(a:3:3,p); outtextxy(160,257,p); str(b:3:3,p); outtextxy(160,212,p); settextstyle(3,0,2); outtextxy(180,224,f); p:=''; str(abs(int):7:3,p); outtextxy(190+length(f)*12,224,p); readln; end; 1: begin newsc(ea); setfillstyle(1,2); bar(170,120,490,180+n*15); setcolor(14); rectangle(175,125,485,175+n*15); rectangle(177,127,483,173+n*15); if ea mod 2 =0 then begin outtextxy(180,130,'Коэффициенты Ньютона-Котеса:'); outtextxy(180,140+(n+1)*15,'Нажмите ENTER для продолжения'); end else begin outtextxy(180,130,'Newton-Cotes coefficients:'); outtextxy(180,140+(n+1)*15,'Press ENTER to continue'); end; hkoef(n,h); for i:=0 to n do begin str(i,p);str(h[i]:2:4,p1); p:='H'+p+' = '+p1; outtextxy(180,140+i*15,p); end; readln; end; 2:begin newsc(ea); setfillstyle(1,2); bar(170,120,490,180+n*15); setcolor(14); rectangle(175,125,485,175+n*15); rectangle(177,127,483,173+n*15); if ea mod 2 =0 then begin outtextxy(180,130,'Значения функции:'); outtextxy(180,140+(n+1)*15,'Нажмите ENTER для продолжения'); end else begin outtextxy(180,130,'Values of function:'); outtextxy(180,140+(n+1)*15,'Press ENTER to continue'); end; for i:=0 to n do begin str(i,p);str(y[i]:2:4,p1); p:='Y'+p+' = '+p1; p1:=''; outtextxy(180,140+i*15,p); str((a+i*(b-a)/n):2:4,p1); str(i,p); if ea mod 2 = 0 then p:=',При '+'X'+p+' = '+p1 else p:=',When '+'X'+p+' = '+p1; outtextxy(285,140+i*15,p); end; readln; end; 3: graphik(ea,a,b,f); 5: begin closegraph; halt; end; end; until (abs(x) mod 6)=4; k:=abs(x) mod 6; end; end. ================================================ ========МОДУЛЬ GRAPHIC======== ================================================ unit graphic; interface uses k_unit,crt,graph; procedure hwg(ea:word); procedure graphik(ea:word;a,b:real;f1:string); implementation procedure hwg(ea:word); {Процедура окна помощи при графике} var f:string; begin settextstyle(0,0,0); setfillstyle(1,3); bar(150,100,390,380); setcolor(0); rectangle(153,103,387,377); rectangle(155,105,385,375); setcolor(14); if ea mod 2 =0 then begin outtextxy(160,115,' ОКНО ПОМОЩИ'); outtextxy(160,140,' Для работы с графиком'); outtextxy(160,155,' используйте клавиши:'); outtextxy(160,180,' PAGE UP-первоначальный'); outtextxy(160,195,' вид графика;'); outtextxy(160,210,' HOME-начальный масштаб;'); outtextxy(160,225,' INSERT-включить/выключеть'); outtextxy(160,240,' заливку области;'); outtextxy(160,255,' DELETE-включить/выключеть'); outtextxy(160,270,' сетку;'); outtextxy(160,285,' END-показать/убрать цифры'); outtextxy(160,300,' F1- Помощь;'); outtextxy(160,315,' Стрелки ВВЕРХ/ВНИЗ- '); outtextxy(160,330,' увеличение/уменьшение'); outtextxy(160,345,' масштаб .'); outtextxy(160,360,'Для возрата нажмите ENTER.'); end else begin outtextxy(160,115,' HELP WINDOW'); outtextxy(160,140,' For the work with graphic'); outtextxy(160,155,' use this keys:'); outtextxy(160,180,' PAGE UP-Primery form of'); outtextxy(160,195,' graphik;'); outtextxy(160,210,' HOME-Primery scale;'); outtextxy(160,225,' INSERT-Turn on/off inking'); outtextxy(160,240,' the field;'); outtextxy(160,255,' DELETE-Turn on/off the'); outtextxy(160,270,' net;'); outtextxy(160,285,' END-View/delete the figures'); outtextxy(160,300,' F1- Help;'); outtextxy(160,315,' Arrows UP/DOWN-Increase/ '); outtextxy(160,330,' lower the scale;'); outtextxy(160,360,'Press ENTER to continue.'); end; readln; setcolor(15); end; procedure graphik(ea:word;a,b:real;f1:string); {процедура построения графиков} var f,f2:string; d:char; i,v,r:integer; x1,x2,n,p,x:integer; c,k,k1:longint; y:array[0..1] of double; begin x1:=-240; x2:=240; c:=24; setcolor(15); n:=0;v:=0;r:=0; repeat cleardevice; settextstyle(0,0,0); if ea mod 2 =0 then begin outtextxy(10,1,'Нажмите F1 для помощи'); str(c/24:2:2,f); f:='Масштаб '+f+':1'; end else begin outtextxy(10,1,'Press F1 for help'); str(c/24:2:2,f); f:='Scale '+f+':1'; end; outtextxy(200,1,f); settextstyle(3,0,1); outtextxy(307,10,'y'); outtextxy(574,235,'x'); outtextxy(310,240,'0'); setlinestyle(1,7,100); line(70,240,580,240); line(320,20,320,460); line(320,20,315,25); line(321,20,326,25); line(580,239,575,244); line(580,240,575,235); line(70,239,580,239); line(321,20,321,460); for i:=-9 to 10 do begin if ((320+i*24)71) then line(320+i*24,240,320+i*24,242); if ((240+i*24)19) then line(320,240+i*24,322,240+i*24); end; setcolor(15); for x:= -240+round((240+x1)/10) to 240+round((240+x1)/10) do begin funktia(1,x-1,x,y,c,f1); k:=round(240-(y[0])*c); k1:=round(240-(y[1])*c); if ((k0)or(k10)) then line(319-round((240+x1)/10)+x,k,320-round((240+x1)/10)+x,k1); end; if (v mod 2)=0 then begin funktia(1,a,b,y,1,f1); k:=round(240-(y[0])*c); k1:=round(240-(y[1])*c); line(320-round((240+x1)/10)+round(a*c),k,320-round((240+x1)/10)+round(a*c),240); line(320-round((240+x1)/10)+round(b*c),k1,320-round((240+x1)/10)+round(b*c),240); if 320-round((240+x1)/10)+a*c begin funktia(1,-240/c,240/c,y,1,f1); k:=round(240-(y[0])*c); line(80,k,80,240); end; if 320-round((240+x1)/10)+b*c>560 then begin funktia(1,(-240-round((240+x1)/10))/c,(240-round((240+x1)/10))/c,y,1,f1); k1:=round(240-(y[1])*c); line(560,k1,560,240); end; for x:= -240 to 240 do begin funktia(1,x-1,x,y,c,f1); k1:=round(240-(y[1])*c); if ((x/c)>a) and ((x/c)0) then begin if (abs(240-k1)>2) then begin if k1 k1:=k1+1 else k1:=k1-1; if c>7 then setfillstyle(6,3) else setfillstyle(1,3); floodfill(320-round((240+x1)/10)+x,k1,15); end; end; end; end; str(x1,f2); outtextxy(1,450,f2); if (n mod 2)=0 then for i:=-9 to 10 do begin settextstyle(2,0,2); setcolor(14); if ((320+i*24)71)and(i<>0) then begin str((i*24+round((240+x1)/10))/c:2:2,f); p:=247; outtextxy(310+i*24,p,f); str(-i*24/c:2:2,f); outtextxy(330,240+i*24,f); end; end; for i:=-9 to 10 do begin setcolor(15); if ((r mod 2)=1) and (i<>0) then begin if ((320+i*24)71) then line(320+i*24,20,320+i*24,460); if ((240+i*24)19) then line(80,240+i*24,560,240+i*24); end; end; setcolor(15); d:=readkey; case d of #75: begin x1:=x1-30; x2:=x2-30; end; #77: begin x1:=x1+30; x2:=x2+30; end; #80: if c>1 then c:=c-1; #72: c:=c+1; #71: c:=24; #79: n:=n+1; #83: r:=r+1; #82: v:=v+1; #73: begin c:=24; n:=0;r:=0;v:=0;x1:=-240;x2:=240; end; #59: hwg(ea); end; until d=#13; end; end. ================================================ ==========МОДУЛЬ UNIT========== ================================================ {$N+} Unit k_unit; {Модуль нахождения интеграл от многочлена q(q-1)..(q-i+1)(q-i-1)..(q-n),} {где n-точность интеграла ,i-номер коофициента. } interface procedure rasposn(f:string;x:real;var ec:word;var t:real); procedure hkoef(n:integer;var h:array of double); procedure funktia(n:integer;a,b:real;var y:array of double;c:real;f:string); procedure koef(w:array of double;n:integer;var e:array of double); procedure mnogochlen(n,i:integer;var c:array of double); function facktorial(n:integer):double; function integral(w:array of double;n:integer):double; function mainint(n:integer;a,b:real;y:array of double):double; implementation procedure rasposn(f:string;x:real;var ec:word;var t:real); {Процедура распознования функции} var k:word; begin k:=pos('x',f); if k<>0 then begin {Распознавание функции} ec:=1; {Код ошибки} t:=x; k:=pos('abs(x)',f); if k<>0 then t:=abs(x); k:=pos('sin(x)',f); if k<>0 then t:=sin(x); k:=pos('cos(x)',f); if k<>0 then t:=cos(x); k:=pos('arctg(x)',f); if k<>0 then t:=arctan(x); k:=pos('sqr(x)',f); if k<>0 then t:=x*x; k:=pos('exp(x)',f); if k<>0 then t:=exp(x); k:=pos('cos(x)*x',f); if k<>0 then t:=cos(x)*x; k:=pos('ln(x)',f); if k<>0 then begin if x>0 then t:=ln(x) else t:=0; end; k:=pos('sqrt(x)',f); if k<>0 then if x>=0 then t:=sqrt(x) else t:=0; k:=pos('arcctg(x)',f); if k<>0 then t:=pi/2-arctan(x); k:=pos('sin(x)/x',f); if k<>0 then if x<>0 then t:=sin(x)/x; end else ec:=0; end; procedure funktia(n:integer;a,b:real;var y:array of double;c:real;f:string); {Процедур подсчет Y-ков и распознавания функции} var t,h,x:real; k,i:integer; es:word; begin h:=(b-a)/n; for i:=0 to n do begin x:=(a+h*i)/c; rasposn(f,x,es,t); y[i]:=t; end; end; procedure koef(w:array of double;n:integer;var e:array of double); {Изменение коофициентов для интеграла} var t:integer; begin for t:=1 to n do e[t]:=w[t]/(n-t+2); end; procedure mnogochlen(n,i:integer;var c:array of double); {процедура нахождения коофициентов при Q^n(q в степени n )} var k,j:integer; d:array[1..100] of double; begin d[1]:=1; for j:=1 to n do begin {Вычисление коэффициентов при раскрытии q*(q-1)*(q-2)*..*(q-n)} d[j+1]:=d[j]*j*(-1); if j>1 then for k:=j downto 2 do d[k]:=d[k]+d[k-1]*j*(-1); end; c[1]:=d[1]; {Деление многочлена на (q-i) по схеме Горнера} for j:=1 to n+1 do c[j]:=i*c[j-1]+d[j]; koef(c,n,c); {Изменение коэффициентов при интегрировании} end; function facktorial(n:integer):double; {функция нахождения факториала } var t:integer; s:double; begin s:=1; if n=0 then s:=1 else for t:=1 to n do s:=s*t; facktorial:=s; end; function integral(w:array of double;n:integer):double; {функция подсчета самого интеграла} var t,p:integer; s,c:double; begin s:=0;p:=n; for t:=0 to p+1 do s:=s+w[t]*exp((p-t+2)*ln(p)); {Подсчет интеграла} integral:=s; end; procedure hkoef(n:integer;var h:array of double); {Процедура подсчета коэф. Ньютона-Котеса} var p,j,d,c,i:integer; kq:array[0..20] of double; s:array[0..20] of double; begin p:=n; if (p mod 2)=1 then {Вычисление половины от всех вычислений коэффициентов} d:=round((p-1)*0.5) else d:=round(0.5*p); for i:=0 to n do begin mnogochlen(p,i,kq); s[i]:=integral(kq,p); {Формирование массива из интегралов} end; for i:=0 to d do begin if ((p-i) mod 2) = 0 then c:=1 else c:=(-1); h[i]:=(c*s[i])/(facktorial(i)*facktorial(p-i)*p); h[p-i]:=h[i]; end; end; function mainint(n:integer;a,b:real;y:array of double):double; {функция подсчета основного интеграла} var sum:double; p,i:integer; kq,h:array[0..20] of double; begin p:=n; hkoef(n,h); sum:=0; for i:=0 to p do sum:=sum+h[i]*y[i]; {Сумма произведений y-ков на коэффициенты} mainint:=sum*(b-a); end; end. ================================================ =======ОСНОВНАЯ ПРОГРАММА======= ================================================ {$N+} program Newton_Cotes_metod;{Программа нахождения определенного интеграла} uses {методом Ньютона-Котеса } k_unit,k_graph,graph,crt; const t=15; var c:char; a1,b1,a,b:real; n1,v,r,n:integer; h,y:array[0..t] of double; ea,k:word; int:double; f:string; begin ea:=10; v:=detect; initgraph(v,r,''); cleardevice; newsc(ea); winwin1; setcolor(15); outtextxy(380,430,'Нажмите F2 для смены языка.'); repeat win1(ea); settextstyle(3,0,1); outtextxy(178,340,'Press Enter...'); delay(13000); bar(178,340,350,365); delay(13000); if keypressed then {Смена языка} begin c:=readkey; if c=#60 then begin ea:=ea+1; newsc(ea); winwin1; setcolor(15); if ea mod 2 =0 then outtextxy(380,430,'Нажмите F2 для смены языка.') else outtextxy(380,430,'Press F2 key to change language.'); end; end; until c=#13; repeat newsc(ea); win2(ea,k); {Ввод способа задания функции} case k of 0: wwod1(ea,y,n,a,b); 1: begin wwod2(ea,ea,n1,a1,b1,f); n:=n1;a:=a1;b:=b1; k:=4; end; end; if k=4 then funktia(n,a,b,y,1,f); int:=mainint(n,a,b,y); {Вычисление интеграла} hkoef(n,h); proline(ea); win3(ea,n,a,b,int,f,h,k); {Последнее меню вывода результатов} until k<>4; closegraph; end.
Рассмотрим результаты тестовых испытаний для функций sin(x) на интервале [-5;3] и exp(x) на интервале [2;8]
Видно, что при увеличении числа узлов интерполяции точность растет, однако при больших n (n>15) наблюдался обратный эффект. Рекомендуемый диапозон n: от 7 до 13.
k:=pos(‘Формула f(x)’,f); if k<>0 then t:= ‘Формула f(x)’; где ‘Формула f(x)’ – желаемая формула для распознования.
Для нахождения интеграла существует много методов, однако, метод Ньютона-Котеса один из самых быстрых: достаточно знать значения коэффициентов для n=4, чтобы с точностью до сотых мгновенно посчитать интеграл. Быстрота и простота –главные части этого метода.
В.И. Грызлов «Турбо Паскаль 7.0» Москва: ДМК 2000г. Данилина «Численные методы» Москва: Высшая школа 1978г. |
|||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
|