Материал: Larin_Anton_8383_cm_21_7

Рисунок 2 – График подынтегральной функции к отмеченными пределами интегрирования

Головная функция и функции для вычисления интеграла методом прямоугольников, трапеций и парабол предоставлены в приложении А. На вход каждая функция принимает следующие параметры: a,b – предел интегрирования

[ , ], eps – требуемая точность вычисления корня, f(x) – интегрируемая функция. Функции возвращают вычисленное значение определенного интеграла и количество отрезков, на которые понадобилось разбить интервал для достижения заданной точности. Вычисление количество отрезков интегрирования функции производится при помощи правила Рунге. Для этого в функциях производятся вычисления интегралов и /2, вычисляется

погрешность по формуле | /2− |, где для формул прямоугольников и трапеций

2 −1

= 2, для формулы Симпсона = 4, после чего данная величина сравнивается с требуемым значением , и если превосходит его, то количество отрезков

6

разбиения удваивается и происходит повторная оценка. В табл. 1 приведены расчеты значения интеграла и количества отрезков , необходимых для его расчета с варьируемой точностью eps, методом ‘*’, где ‘*’ принимает значения:

R – метод прямоугольников, T – метод трапеций, S – метод Симпсона(парабол).

Таблица 1 – Вычисление значения интеграла различными методами с варьированием значения eps

Выводы.

Проанализировав результаты расчетов, мы можем сделать вывод, что число итераций, необходимых для вычисления интеграла по составным формулам прямоугольника, трапеции и Симпсона возрастает с ростом требуемой точности результата.

По результатам оценки необходимого количества отрезков разбиения по правилу Рунге можем сделать вывод, что для данной функции наименьшее количество итераций для получения необходимой точности результата требуется при вычислении интеграла по формуле Симпсона(парабол), далее в порядке увеличения числа итераций идет формула прямоугольников, а затем формула трапеций. Следовательно, в данном случае метод парабол дает наилучшие результаты, так как дает наименьшую погрешность из исследованных методов.

7

ПРИЛОЖЕНИЕ А

ИСХОДНЫЙ КОД ПРОГРАММЫ

#include <stdio.h> #include <stdlib.h> #include <math.h>

#define CSQR(CX) ((CX)*(CX))

#ifndef M_PI

#define M_PI 3.14159265358979323846 #endif // !M_PI

#define A (M_PI/2) #define B M_PI

#define E 0.000000000000001

double F(double x) {

if (x < A)printf("Boo<%lf",x); if(x>B)printf("Foo>%lf",x); return sqrt(x)*exp(-CSQR(x));

}

double _RECT(double h, int n, double x0, double(*f)(double)) { double x = x0;

double s = 0;

for (int i = 0; i < n; i++) { s += f(x + h / 2);

x += h;

}

s *= h; return s;

}

double RECT(double a, double b, double e, double(*f)(double), int &n) { #define H ((b - a) / n)

double h= NAN, h2=NAN, Ih=NAN, Ih2= NAN; int k = 2;

double e0 = e * ((1 << k) - 1); double d = e0+1;

n = 1;

Ih2 = _RECT(H, n, a, f); while(d>=e0){

8

}

return (Ih2 + (Ih2-Ih)/3); #undef H

}

double _TRAP(double h, int n, double x0, double(*f)(double)) { #define X(CI) (x0+(CI)*h)

//double x = x0; double s = 0;

for (int i = 0; i < n; i++) { s += f(X(i))+f(X(i+1)); //x += h;

}

s *= h/2; return s;

}

double TRAP(double a, double b, double e, double(*f)(double), int &n) { #define H ((b - a) / n)

double h = NAN, h2 = NAN, Ih = NAN, Ih2 = NAN; int k = 2;

double e0 = e * ((1 << k) - 1); double d = e0 + 1;

n = 1;

Ih2 = _TRAP(H, n, a, f); while (d >= e0) {

n *= 2; Ih = Ih2;

Ih2 = _TRAP(H, n, a, f); d = fabsl(Ih2 - Ih);

}

return (Ih2 - (Ih2 - Ih) / 3); #undef H

}

double _SIMPS(double h, int n, double x0, double(*f)(double)) { //#define X(CI) (x0+(CI)*h)

int m = n / 2; double s = 0;

for (int i = 0; i < m; i++) {

9

s += f(X(2 * i)) + 4 * f(X(2 * i + 1)) + f(X(2 * i + 2));

}

s *= h / 3; return s;

}

double SIMPS(double a, double b, double e, double(*f)(double), int &n, bool echo = false) {

#define H ((b - a) / n)

double h = NAN, h2 = NAN, Ih = NAN, Ih2 = NAN; int k = 4;

double e0 = e * ((1 << k) - 1); double d = e0 + 1;

n = 1;

Ih2 = _SIMPS(H, n, a, f); while (d >= e0) {

n *= 2; Ih = Ih2;

Ih2 = _SIMPS(H, n, a, f); d = fabsl(Ih2 - Ih);

}

return (Ih2 - (Ih2 - Ih) / 15); #undef H

}

int main(){

double s_R, s_T, s_S; int k_R, k_T, k_S;

printf("Eps\t\ts_R\t\ts_T\t\ts_S\t\tk_R\tk_T\tk_S\n"); for (double eps = 0.1; eps > 0.000001; eps /= 10) {

s_R = RECT(A, B, eps, F, k_R); //s_T = TRAP_(eps, k_T);

s_T = TRAP(A, B, eps, F, k_T); s_S = SIMPS(A, B, eps, F, k_S);

printf("%lf\t%lf\t%lf\t%lf\t%d\t%d\t%d\n", eps, s_R, s_T, s_S, k_R, k_T, k_S);

}

system("pause"); return 0;

}

10