В MATLAB нет функции, выполняющей численное дифференцирование. Исключение составляют функции вычисления производной от алгебраического полинома, заданного в символьной форме – diff и заданного в форме вектора – polyder. Функция polyder возвращает вектор коэффициентов производной от многочлена. В первой строке примера 37 был фактически определен многочлен 2x3 + 3x2 + 5x + 7 = 0. При аналитическом дифференцировании этого выражения степень многочлена понижается и число коэффициентов уменьшается на единицу (производная от последнего числа в левой части равна нулю).
Если в качестве аргумента функции diff передать вектор, то в результате будет вычислен вектор разностей. Причина, по которой численное дифференцирование в MATLAB реализовано только одной функцией, заключается в том, что вычисление производных основано на интерполяции дифференцируемой функции алгебраическими полиномами не очень высоких степеней. Простейший способ вычисления производной с использованием функции diff заключается в простом вычислении отношения приращения функции к приращению аргумента. Понятно, что в определении производной приращение аргумента обязано быть бесконечно малым, и, используя такое приближение для производной, точных результатов ожидать не стоит. Для получения более точных результатов следует использовать специальные формулы (например, Ньютона или Стирлинга) для численного дифференцирования. В примере 38 построен график квадратичной функции и её производной. Как легко убедиться, выполнив этот сценарий, производная будет линейной функцией. Аналитическое нахождение производной приведет к выражению 2*xd, и разность этого значения с приближенным значением производной dxdy будет равна 0.1
Численное интегрирование – вычисление площади под графиком функции (в одномерном случае), также реализуется при помощи интерполяции функции различными способами. Для эффективного вычисления интеграла необходимо представлять особенности подынтегральной функции. Например, в случае разрывов (устранимых) функции в некоторой точке следует разбить интеграл на два так, чтобы точка разрыва попала на границу. Для обеспечения нужной точности и скорости сходимости следует выбирать алгоритм интегрирования в соответствии с характером роста интегрируемой функции.
| Таблица 10. Функции MATLAB для вычисления определенных интегралов | |
| Название | Значение |
| trapz(y,x) | Вычисляет определенный интеграл методом трапеций. Аргументы: y - подынтегральная функция, x - вектор аргумент. Обеспечивает достаточную точность для линейных функций |
| quad(y(x), xmin, xmax, eps, trace) | Вычисляет определенный интеграл по формуле Симпсона, обеспечивает достаточную точность для функций, четвертая производная от которых равна нулю. Обязательными аргументами являются имя функции и границы интегрирования. В качестве дополнительного аргумента можно указать точность вычислений - eps, по умолчанию равна 10-6. Специальная опция trace позволяет включить вывод промежуточных значений на каждой итерации, для этого нужно задать любое ненулевое значение |
| quadl(y(x), xmin, xmax, eps, tarce) | Вычисляет определенный интеграл методом Лобатто. Смысл аргументов, такой же как у quad |
| quadgk(y(x), xmin, xmax, eps, tarce) | Вычисляет определенный интеграл методом Гаусса-Кронрода. Смысл аргументов такой же как у quad |
| dblquad(z(x,y), xmin, xmax, ymin, ymax, eps, method) | Вычисляет двойной определенный интеграл над прямоугольной областью. Необходимо указать границы изменения каждой переменной. Дополнительный параметр method может принимать два значения - quad и quadl |
| quad2d(z(x,y), xmin, xmax, ymin, ymax) | Вычисляет двойной определенный интеграл над плоской областью. Необходимо указать границы изменения каждой переменной |
| triplequad(z(x,y), xmin, xmax, ymin, ymax, eps, method) | Вычисляет тройной определенный интеграл. Необходимо указать границы изменения каждой переменной. Дополнительный параметр method может принимать два значения - quad и quadl |
Вызов функций quad и quadl с двумя выходными парметрами приводит к выводу в командное окно числа обращений к вычислению подынтегральной функции. Вызов quadgk с двумя параметрами приводит к выводу приближенного значения погрешности интегрирования
1 >> f=@(x)exp(−x).∗cos(x)
2 f =
3 @(x)exp(−x).∗cos(x)
4 >> format longEng
5 >> s=quad(f,0,1)
6 s =
7 555.396877498501e−003
8 >> s=quad(f,0,1,1e−8)
9 s =
10 555.396882649120e−003
11 >> s=quad(f,0,1,1e−12)
12 s =
13 555.396882653349e−003
14 >> [s,iter]=quad(f,0,1,1e−15)
15 s =
16 555.396882653350e−003
17 iter =
18 781.000000000000e+000
19 >> [s2,iter2]=quadl(f,0,1,1e−15)
20 s2 =
21 555.396882653350e−003
22 iter2 =
23 48.0000000000000e+000
24 >> [s3,err]=quadgk(f,0,1)
25 s3 =
26 555.396882653350e−003
27 err =
28 21.1636264069171e−018
29 >> [s2,iter2]=quadl(f,0,1,1e−15,2)
30 18 0.0000000000 5.00000000e−01 0.5553968827
31 23 0.0000000000 4.58758548e−02 0.0875482138
32 28 0.0917517095 9.23207464e−02 0.1510834176
33 33 0.2763932023 1.11803399e−01 0.1406211476
34 38 0.5000000000 1.11803399e−01 0.0995695619
35 43 0.7236067977 9.23207464e−02 0.0561174607
36 48 0.9082482905 4.58758548e−02 0.0204570810
37 s2 =
38 555.396882653350e−003
39 iter2 =
40 48.0000000000000e+000
41 >> f=@(x)exp(−x)
42 f =
43 @(x)exp(−x)
44 >> [s,err]=quadgk(f,2,Inf)
45 s =
46 135.335283236613e−003
47 err =
48 2.22626006030936e−012
Как видно из примера, функция quadgk обеспечивает вычисление интеграла от заданной функции с погрешностью, не превышающей 21.1636264069171 ⋅ 10−018, а функции quad и quadl обеспечивают точность 10−015 за 781 и 48 итераций соответственно. В фирменной документации утверждается, что quadl обеспечивает более высокую точность, чем quad; а quadgk обеспечивает высокую точность для периодических функций и, кроме этого, поддерживает бесконечные пределы интегрирования.