Реферат: Метод Симпсона на компьютере
МОСКОВСКИЙГОСУДАРСТВЕННЫЙ СТРОИТЕЛЬНЫЙ УНИВЕРСИТЕТ
КУРСОВАЯ РАБОТА«Программа приближенноговычисления определенного интеграла с помощью ф – лы Симпсона на компьютере»
Выполнил:
студент ф – та ЭОУС – 1 – 12
Валюгин А. С.
Принял:
Зоткин С. П.
Москва2001
1. Введение
Определенный интеграл от функции,имеющей неэлементарную первообразную, можно вычислить с помощью той или иной приближенной формулы. Для решения этой задачи на компьютере, средипрочих, можно воспользоваться формуламипрямоугольников, трапеций или формулой Симпсона. В данной работе рассматриваетсяименно последняя.
Рассмотрим функцию y = f(x). Будем считать, что на отрезке [a, b] она положительна и непрерывна. Найдем площадькриволинейной трапеции aABb(рис. 1).
/>
рис. 1
Для этого разделим отрезок [a, b]точкой c = (a + b) / 2 пополам и в точке C(c, f(c))проведем касательную к линии y = f(x). После этого разделим [a, b] точками p и q на3 равные части и проведем через них прямые x = p и x = q. Пусть P и Q – точки пересечения этих прямых скасательной. Соединив A с P и B с Q,получим 3 прямолинейные трапеции aAPp, pPQq, qQBb. Тогда площадь трапеции aABb можно приближенно посчитать по следующей формуле
I » (aA + pP) / 2 * h +(pP + qQ) / 2 * h + (qQ + bB) / 2 * h, где h = (b – a) / 3.
Откуда получаем
I » (b – a) / 6 * (aA + 2 * (pP + qQ) + bB)
заметим, что aA = f(a), bB = f(b), а pP + qQ = 2 * f(c),в итоге получаем малую фор – лу Симпсона
/> <td/>I » (b – a) / 6 * (f(a) + 4 * f(c) + f(b)) (1)
Малая формула Симпсона дает интегралс хорошей точностью, когда график подинтегральной функции мало изогнут, вслучаях же, когда дана более сложная функция малая формула Симпсона непригодна.Тогда, чтобы посчитать интеграл заданной функции нужно разбить отрезок [a, b] на nчастей и к каждому из отрезков применить формулу (1). После указанных вышедействий получится “большая” формулаСимпсона, которая имеетвид,
I »h / 3 * (Yкр+ 2 * Yнеч+ 4 * Yчет) (2)
где Yкр = y1 + yn, Yнеч = y3 + y5 + … +yn – 1, Yчет = y2 + y4 + … + yn – 2, а h = (b – a) / n.
Задача. Пусть нужно проинтегрироватьфункцию f(x) = x³(x - 5)² на отрезке [0, 6] (рис.2). На этом отрезке функция непрерывна и принимает только неотрицательныезначения, т. е. знакопостоянна.
/>
рис. 2
Для выполнения поставленной задачисоставлена нижеописанная программа, приближенно вычисляющая определенныйинтеграл с помощью формулы Симпсона. Программа состоит из трех функций main, f и integral. Функция main вызываетфункцию integral для вычисления интеграла ираспечатывает на экране результат. Функция f принимает аргумент x типа floatи возвращает значение интегрируемой функции в этой точке. Integral – основная функция программы: онавыполняет все вычисления, связанные с нахождением определенного интеграла. Integral принимает четыре параметра: пределыинтегрирования типа float,допустимую относительную ошибку типа float и указатель на интегрируемую функцию. Вычисления выполняются до техпор, пока относительная ошибка, вычисляемая по формуле
| (In/2 – In) / In |,
где In интеграл причисле разбиений n, не будетменьше требуемой. Например, допустимая относительная ошибка e = 0.02 это значит, что максимальнаяпогрешность в вычислениях будет не больше, чем In * e =0.02 * In. Функция реализована с экономиейвычислений, т. е. учитывается, что Yкр постоянная,а Yнеч = Yнеч + Yчет,поэтому эти значения вычисляются единожды. Высокая точность и скоростьвычисления делают использование программы на основе формулы Симпсона болеежелательным при приближенном вычислении интегралов, чем использование программна основе формулы трапеции или метода прямоугольников.
Ниже предлагается блок – схема, спецификации, листинг и ручной счет программы на примере поставленнойвыше задачи. Блок – схема позволяет отследить и понятьособенности алгоритма программы, спецификации дают представление о назначениикаждой переменной в основной функции integral, листинг - исходный код работающейпрограммы с комментариями, а ручной счет предоставляет возможностьпроанализировать результаты выполнения программы.
2. Блок – схема программы
/>
/> ДА
/> /> /> /> /> /> /> /> /> />НЕТ
/>
3. Спецификации
Имя переменной
Тип
Назначение
n int Число разбиений отрезка [a, b] i int Счетчик циклов a float Нижний предел интегрирования b float Верхний предел интегрирования h float Шаг разбиения отрезка e float Допустимая относительная ошибка f float (*) Указатель на интегрируемую фун — цию s_ab float Сумма значений фун – ции в точках a и b s_even float Сумма значений фун – ции в нечетных точках s_odd float Сумма значений фун – ции в четных точках s_res float Текущий результат интегрирования s_pres float Предыдущий результат интегрирования4. Листинг программы
#include <stdio.h>
#include <math.h>
/*Прототип фун – ции, вычисляющей интеграл */
float integral(float, float, float, float (*)(float));
/* Прототип фун – ции, задающей интегрируемую фун – цию */
float f(float);
main()
{
float result;
result= integral(0, 6, .1, f);
printf("%f", result);
return0;
}
/*Реализация фун – ции, задающей интегрируемую фун – цию */
float f(float x)
{
/* Функция f(x) = x³(x - 5)² */
return pow(x, 3) * pow(x — 5,2);
}
/*Реализация фун – ции, вычисляющей интеграл */
float integral(float a, float b, float e, float (*f)(float))
{
int n = 4, i; /*Начальное число разбиений 4 */
float s_ab = f(a) + f(b); /*Сумма значений фун – ции в a и b */
float h = (b – a) / n; /* Вычисляем шаг*/
float s_even = 0, s_odd;
float s_res = 0, s_pres;
/*Сумма значений фун – ции в нечетных точках */
for (i = 2; i < n; i += 2) {
s_even += f(a + i * h);
}
do {
s_odd = 0;
s_pres = s_res;
/* Суммазначений фун – ции в четных точках */
for(i = 1; i < n; i += 2) {
s_odd += f(a + i * h);
}
/*Подсчет результата */
s_res = h / 3 * (s_ab + 2 * s_even + 4 *s_odd);
/* Избегаем деления на ноль */
if (s_res == 0) s_res = e;
s_even += s_odd;
n *= 2;
h /= 2;
} while (fabs((s_pres — s_res) / s_res) > e);/*Выполнять до тех пор, пока результат не будет удовлетворять допустимой ошибке*/
return fabs(s_res); /*Возвращаем результат */
}
5. Ручной счет
Таблицаконстантных значений для n= 8
Имя переменной
Значение
a b 6 e .1 s_ab 216 h .75Подсчетs_even
i
a + i * h
f(a + i * h)
s_even
2 1.5 41.34375 41.34375 4 3 108 149.34375 6 4.5 22.78125172.125
Подсчетs_odd
i
a + i * h
f(a + i * h)
s_odd
1 .75 7.62012 7.62012 3 2.25 86.14158 93.7617 5 3.75 82.3973 176.159 7 5.25 9.044185.203
Подсчетs_res
òf(x) dx
s_res = h / 3 * (s_ab + 2 * s_even + 4 * s_odd)
Абсолютнаяошибка
324 325.266 1.266