Численное дифференцирование

Численное дифференцирование — совокупность методов вычисления значения производной дискретно заданной функции.

В основе численного дифференцирования лежит аппроксимация функции, от которой берется производная, интерполяционным многочленом.  Все основные формулы численного дифференцирования могут быть получены при помощи первого интерполяционного многочлена Ньютона (формулы Ньютона для начала таблицы).

Основными задачами являются вычисление производной на краях таблицы и в её середине. Для равномерной сетки формулы численного дифференцирования «в начале таблицы» можно представить в общем виде следующим образом:fi′=1bh∑jajfi+j+Δ(f),{f'_i} = \frac{1}{{bh}}\sum\limits_j {{a_j}} {f_{i + j}} + \Delta (f),

где Δ(f)\Delta (f) — погрешность формулы. Здесь коэффициенты aj{a_j} и bb зависят от степени nnn использовавшегося интерполяционного многочлена, то есть от необходимой точности (скорости сходимости к точному значению при уменьшении шага сетки) формулы. Коэффициенты представлены в таблице

na0a1a2a3a4a5b1−11000012−34−100023−1118−920064−2548−3616−30125−137300−300200−751260\begin{array}{cccccccccccccccccccc} n&{{a_0}}&{{a_1}}&{{a_2}}&{{a_3}}&{{a_4}}&{{a_5}}&b \\ 1&{ - 1}&1&0&0&0&0&1 \\ 2&{ - 3}&4&{ - 1}&0&0&0&2 \\ 3&{ - 11}&{18}&{ - 9}&2&0&0&6 \\ 4&{ - 25}&{48}&{ - 36}&{16}&{ - 3}&0&{12} \\ 5&{ - 137}&{300}&{ - 300}&{200}&{ - 75}&{12}&{60} \end{array}

Формулы численного дифференцирования

Один из универсальных способов построения формул численного дифференцирования состоит в том, что по значениям функции f(x)f\left( x \right) в некоторых узлах x0,x1,…,xn−1,xn{x_0},{x_1},…,{x_{n - 1}},{x_n} строят интерполяционный полином PN(x){P_N}\left( x \right)(в форме Лагранжа или в форме Ньютона) и приближенно полагают f(r)(x)≈PN(r)(x),0⩽r⩽N{f^{(r)}}(x) \approx P_N^{(r)}(x),0 \leqslant r \leqslant N

В ряде случаев, наряду с приближенным равенством удается (например, используя формулу Тейлора) получить точное равенство, содержащее остаточный член RR (погрешность численного дифференцирования)f(r)(x)=PN(r)(x)+R,0⩽r⩽N{f^{(r)}}(x) = P_N^{(r)}(x) + R,0 \leqslant r \leqslant N

Такие формулы называются формулами численного дифференцирования с остаточными членами.

Степень, с которой входит величина h=max⁡hi (hi=xi−xi−1)h = \max {h_i}{\mkern 1mu} ({h_i} = {x_i} - {x_{i - 1}}) в остаточный член, называется порядком погрешности формулы численного дифференцирования. Формулы с отброшенными остаточными членами называются просто формулами численного дифференцирования.

Ниже приводятся несколько распространенных формул численного дифференцирования с остаточными членами для первой (r=1)(r=1) и второй (r=2)(r=2) производных в узлах, расположенных с постоянным шагом  h>0h>0:

r=1,N=1r = 1,N = 1 (два узла):
f′(xo)=(f1−f0)/h−hf»(ξ)/2f′(x1)=(f1−f0)/h+hf»(ξ)/2\begin{aligned}&f'({x_o}) = ({f_1} - {f_0})/h - hf»\left( \xi \right)/2 \\ &f'({x_1}) = ({f_1} - {f_0})/h + hf»(\xi )/2\end{aligned}
r=1,N=2r = 1,N = 2 (три узла):
f′(xo)=(−3f0+4f1−f2)/2h+h2f»′(ξ)/3f′(x1)=(f2−f0)/2h−h2f»′(ξ)/6f′(x2)=(f0−4f1+3f2)/2h+h2f»′(ξ)/3\begin{aligned}&f'({x_o}) = ( - 3{f_0} + 4{f_1} - {f_2})/2h + {h^2}f»'(\xi )/3 \\ &f'({x_1}) = ({f_2} - {f_0})/2h - {h^2}f»'(\xi )/6 \\ &f'({x_2}) = ({f_0} - 4{f_1} + 3{f_2})/2h + {h^2}f»'(\xi )/3\end{aligned}
r=2,N=2r = 2,N = 2 (три узла):
f»(xo)=(f0−2f1+f2)/h2−hf»′(ξ)f»(x1)=(f0−2f1+f2)/h2−h2f(4)(ξ)/12f»(x2)=(f0−2f1+f2)/h2+hf»′(ξ)\begin{aligned}&f»({x_o}) = ({f_0} - 2{f_1} + {f_2})/{h^2} - hf»'(\xi ) \\ &f»({x_1}) = ({f_0} - 2{f_1} + {f_2})/{h^2} - {h^2}{f^{(4)}}(\xi )/12 \\ &f»({x_2}) = ({f_0} - 2{f_1} + {f_2})/{h^2} + hf»'(\xi )\end{aligned}
r=2,N=3r = 2,N = 3 (четыре узла):
f»(x0)=(2f0−5f1+4f2−f3)/h2+11h2f(4)(ξ)/12f»(x1)=(f0−2f1+f2)/h2−h2f(4)(ξ)/12f»(x2)=(f1−2f2+f3)/h2−h2f(4)(ξ)/12f»(x3)=(−f0+4f1−5f2+2f3)/h2+11h2f(4)(ξ)/12\begin{aligned}&f»({x_0}) = (2{f_0} - 5{f_1} + 4{f_2} - {f_3})/{h^2} + 11{h^2}{f^{(4)}}(\xi )/12 \\ &f»({x_1}) = ({f_0} - 2{f_1} + {f_2})/{h^2} - {h^2}{f^{(4)}}(\xi )/12 \\ &f»({x_2}) = ({f_1} - 2{f_2} + {f_3})/{h^2} - {h^2}{f^{(4)}}(\xi )/12 \\ &f»({x_3}) = ( - {f_0} + 4{f_1} - 5{f_2} + 2{f_3})/{h^2} + 11{h^2}{f^{(4)}}(\xi )/12\end{aligned}

где hh — шаг сетки, а точка ξ\xi — некоторая промежуточная точка.

В формулах численного дифференцирования с постоянным шагом hh значения функции f(x)f\left( x \right) делятся на hr{h^r}, где rr-порядок вычисляемой производной. Поэтому при малом hh неустранимые погрешности в значениях функции f(x)f\left( x \right) оказывают сильное влияние на результат численного дифференцирования. Таким образом, возникает задача выбора оптимального шага hh, так как погрешность собственно метода стремится к нулю при h→0h \to 0, а неустранимая погрешность растет. В результате общая погрешность, которая возникает при численном дифференцировании, может неограниченно возрастать при h→0h \to 0. Поэтому операцию численного дифференцирования считают некорректной.

Если функция задана таблицей значений в узлах с постоянным шагом hh, например по результатам измерений, производную вычисляют приближённо по соседним значениям f0,f1,f2f_0, f_1, f_2 в узлах x0,x1,x2x_0, x_1, x_2. В среднем узле x1x_1 правая разность равна f2−f1h\frac{f_2 - f_1}{h}, левая — f1−f0h\frac{f_1 - f_0}{h}, а центральная — f2−f02h\frac{f_2 - f_0}{2h}.

У первых двух формул порядок погрешности hh, а у центральной h2h^2: если уменьшить шаг вдвое, погрешность центральной формулы упадёт примерно в четыре раза, а у первых двух — только в два. Поэтому центральная формула, как правило, точнее. Вторую производную находят по формуле f′′(x1)≈f0−2f1+f2h2f''(x_1) \approx \frac{f_0 - 2f_1 + f_2}{h^2}. Слишком малый шаг вреден, если данные измерены с погрешностью: разности делятся на малое число (hh или h2h^2), и ошибка измерения сильно возрастает.

12345670,511,522,53xy
2
Правая разность
0,2247
Центральная разность
0,2588
Точно f′(4)
0,25
Кривая — f(x)=xf(x) = \sqrt x, чёрная прямая — касательная в точке x=4x = 4 с наклоном f′(4)=0,25f'(4) = 0{,}25. Наклон фиолетовой секущей равен правой разности, а зелёной — центральной. Уменьшайте hh: обе секущие приближаются к касательной, но зелёная гораздо быстрее, потому что её погрешность порядка h2h^2.

Именно так по последовательным положениям GPS оценивают скорость, а программы для моделирования и обработки измерений заменяют производные этими формулами.

Скорость бегуна

Часы с приёмником GPS записывают пройденный бегуном путь каждые 2 секунды. Измерено: s(10)=38,0s(10) = 38{,}0 м, s(12)=46,2s(12) = 46{,}2 м и s(14)=54,8s(14) = 54{,}8 м. Оцените скорость бегуна в момент t=12t = 12 с по правой, левой и центральной разностям, а также его ускорение.

Решение

  1. Узлы t0=10t_0 = 10, t1=12t_1 = 12, t2=14t_2 = 14, шаг h=2h = 2 с. Правая разность: 54,8−46,22=4,3\frac{54{,}8 - 46{,}2}{2} = 4{,}3 м/с. Левая: 46,2−38,02=4,1\frac{46{,}2 - 38{,}0}{2} = 4{,}1 м/с.
  2. Центральная разность: 54,8−38,02⋅2=16,84=4,2\frac{54{,}8 - 38{,}0}{2 \cdot 2} = \frac{16{,}8}{4} = 4{,}2 м/с. Это среднее арифметическое двух предыдущих значений и самая надёжная оценка.
  3. Ускорение — вторая производная: 38,0−2⋅46,2+54,822=0,44=0,1\frac{38{,}0 - 2 \cdot 46{,}2 + 54{,}8}{2^2} = \frac{0{,}4}{4} = 0{,}1 м/с².

Ответ: На 12-й секунде скорость бегуна около 4,2 м/с, т. е. около 15 км/ч, и он слегка ускоряется, с ускорением около 0,1 м/с².

Задача для тренировки

Вычислите правую и центральную разности для f(x)=xf(x) = \sqrt x в точке x=4x = 4 с шагом h=0,1h = 0{,}1 и сравните их с точным значением f′(4)f'(4).

Показать решение

Правая: 4,1−20,1≈0,2485\frac{\sqrt{4{,}1} - 2}{0{,}1} \approx 0{,}2485, центральная: 4,1−3,90,2≈0,25002\frac{\sqrt{4{,}1} - \sqrt{3{,}9}}{0{,}2} \approx 0{,}25002. Точно f′(4)=124=0,25f'(4) = \frac{1}{2\sqrt 4} = 0{,}25, поэтому центральная формула ошибается примерно на 0,000020{,}00002, а правая — примерно на 0,00150{,}0015.