МОСКОВСКИЙ ГОСУДАРСТВЕННЫЙ СТРОИТЕЛЬНЫЙ УНИВЕРСИТЕТ
КУРСОВАЯ РАБОТА
«Программа приближенного вычисления определенного интеграла с помощью ф – лы Симпсона на компьютере»
Выполнил:
студент ф – та ЭОУС – 1 – 12
Валюгин А. С.
Принял:
Зоткин С. П.
Москва 2001
Введение
Определенный интеграл от функции, имеющей неэлементарную первообразную, можно вычислить с помощью той или иной приближенной формулы. Для решения этой задачи на компьютере, среди прочих, можно воспользоваться формулами прямоугольников, трапеций или формулой Симпсона. В данной работе рассматривается именно последняя.
Рассмотрим функцию 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), в итоге получаем малую фор – лу Симпсона
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, листинг - исходный код работающей программы с комментариями, а ручной счет предоставляет возможность проанализировать результаты выполнения программы.
Блок – схема программы
Ввод a, b, e, f(x)
n = 4, h = (b – a) / n
s_ab = f(a) + f(b)
s_even = 0, s_res = 0
i = 2, n – 1, 2
s_even = s_even + f(a + i * h)
s_odd = 0, s_pres = s_res
i = 1, n – 1, 2
s_odd = s_odd + f(a + i * h)
s_res = h / 3 * (s_ab + 2 * s_even + 4 * s_odd)
s_even = s_even + s_odd, n = n / 2, h = h / 2
| (s_pres – s_res) / s_res | > e
ДА
НЕТ
Вывод s_res
Спецификации
Имя переменной | Тип | Назначение |
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 | Предыдущий результат интегрирования |
Листинг программы
#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);
return 0;
}
/* Реализация фун – ции, задающей интегрируемую фун – цию */
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); /* Возвращаем результат */
}
Ручной счет
Таблица константных значений для n = 8
Имя переменной | Значение |
a | 0 |
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.78125 | 172.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.044 | 185.203 |
Подсчет s_res
o f(x) dx | s_res = h / 3 * (s_ab + 2 * s_even + 4 * s_odd) | Абсолютная ошибка |
324 | 325.266 | 1.266 |
1. Пример вычисления интерграла по формуле симпсона
2. Пример задачи по формуле симпсона
3. Вычислить по формуле симпсона интеграл
4. Вычисление интегралов примеры
6. Вычисление площади с помощью двойного интеграла
7. Вычисление объемов с помощью двойного интеграла
9. Коэффициент профессора ф в коньшина к примеры задач
10. Примеры вычисления тройных интегралов
11. Примеры сетевых графиков и вычисление
12. Примеры решения задач на эластичность спроса
13. Эластичность спроса примеры решения задач
14. Картель примеры
15. Вычисление формулой симпсона