Численное интегрирование средствами MATLAB
В MATLAB реализованы множество современных методов численного интегрирования. Мы рассмотрим простейшие из них.
А. Метод трапеций. Метод трапеции реализован в MATLAB несколькими функциями:
- cumtrapz(у)
- cumtrapz(x,y)
- trapz(y)
- trapz(x,y)
Наиболее интересна последняя из них. Данная функция вычисляет интеграл от функции у(х) по х методом трапеций. Аргумент и функция задаются в виде векторов или х — в виде вектора, а у — в виде матрицы. Если у(х) — матрица, то функция возвращает вектор значений интеграла для каждого столбца матрицы.
Пример.Пусть подынтегральная функция имеет вид
у(х) = х*еx + ln(x) + 1
Необходимо вычислить определенный интеграл в диапазоне от 1 до 10 с шагом 0.5.
Решение:
x=1:0.5:10;
y=x.*exp(x)+log(x)+1;
trapz(x,y)
Решение:
ans =
2.032841320958599e+005
Б. Метод Симпсона. Метод Симпсона реализован в MATLAB также несколькими функциями. Мы рассмотрим простейшую из них:
quad('fun', a, b)
quad('fun', a, b, tol)
Обозначения:
- 'fun'— подынтегральная функция, взятая в одинарные кавычки;
- а, b — пределы интегрирования;
- tol — относительная погрешность, задаваемая пользователем; по умолчанию tol=10-3.
Пример. Вычислить интеграл от функции x+ex на отрезке [1, 2] с точностью 10-5.
quad('x+exp(x)', 1, 2, 1e-5)
Результат:
ans = 6.170774280645000e+000
Правило Рунге оценки погрешности интегрирования
В формулах для оценки погрешности квадратурных формул R используются значения производных подинтегральной функции, что требует дополнительного анализа и вычислений. В связи с этим получило распространение практическое правило Рунге оценки погрешности.
Пусть
- I – точное значение интеграла,
- I(n) – значение интеграла вычисленное при n узлах интегрирования h = (b-a)/n,
- I(2n) – значение интеграла вычисленное при 2*n узлах интегрирования, h = (b-a)/2n.
Необходимо определить, с какой точностью вычислен итеграл I(2n), т.е. найти абсолютную погрешность
Для непосредственно определения данной погрешности необходимо найти максимум модуля соответствующей производжной от интегрируемой функции на отрезке [a, b]. Часто это достаточно трудоемкий или вообще невозможный процесс. Напрмер если интегрируемая функция задана таблично. В таких случаях оценку погрешности величины I(2n) можно провести следующим образом:
Здесь m = 3 для методов средних прямоугольников и трапеций, m = 15 для метода Симпсона.
Примечание. Если решается задача численного вычисления интеграла с заданной точностью, процесс удвоения числа узлов интегрирования продолжается до тех пор, пока величена не станет меньше заданной погрешности.