Численное дифференцирование
Раздел № 5
Задача приближенного вычисления производной может возникнуть в тех случаях, когда неизвестно аналитическое выражение для исследуемой функции. Функция может быть задана таблично, или известен только график функции, полученный, например, в результате показаний датчиков параметров технологического процесса.
Иногда, при решении некоторых задач на компьютере, из-за громоздкости выкладок может оказаться более удобным вычисление производных численным методом, чем аналитическим. При этом, разумеется, необходимо обосновать применяемый численный метод, т. е. убедиться в том, что погрешность численного метода находится в приемлемых границах.
Одним из эффективных методов решения дифференциальных уравнений является разностный метод, когда вместо искомой функции рассматривается таблица ее значений в определенных точках, при этом производные приближенно заменяются разностными формулами.
Пусть известен график функции у = f (х ) на отрезке [а ,b ].Можно построить график производной функции, вспомнив ее геометрический смысл. Воспользуемся тем фактом, что производная функции в точке х равна тангенсу угла наклона к оси абсцисс касательной к ее графику в этой точке.
Если х = х 0 ,найдем у 0 = f (x 0)с помощью графика и затем проведем касательную АВ к графику функции в точке (х 0 , y 0) (рис. 5.1). Проведем прямую, параллельную касательной АВ, через точку (-1, 0) и найдем точку у 1 ее пересечения с осью ординат. Тогда значение у 1 равно тангенсу угла наклона касательной к оси абсцисс, т. е. производной функции f (x )в точке х 0:
у 1 = = tgα = f ¢ (x 0), и точка М 0 (х 0 , у 1) принадлежит графику производной.
Чтобы построить график производной, необходимо разбить отрезок [а , b ]на несколько частей точками х i , затем для каждой точки графически построить значение производной и соединить полученные точки плавной кривой с помощью лекал.
На рис. 5.2 показано построение пяти точек М 1, М 2 ,... , М 5 и графика производной.
Алгоритм построения графика производной:
1. Строим касательную к графику функции у = f (x )в точке (х 1 , f (x 1));из точки (-1, 0) параллельно касательной в точке (х 1 , f (x 1)) проведем прямую до пересечения с осью ординат; эта точка пересечения дает значение производной f ¢ (х 1).Строим точку М 1 (х 1 , f ¢ (х 1)).
2. Аналогично построим остальные точки М 2 , М 3 , М 4 и М 5 .
3. Соединяем точки М 1 , М 2 , М 3 , М 4 , М 5 плавной кривой.
|
Полученная кривая является графиком производной.
Точность графического способа определения производной невысока. Мы приводим описание этого способа только в учебных целях.
Замечание . Если в алгоритме построения графика производной вместо точки (-1, 0) взять точку (-l ,0), где l > 0, то график будет построен в другом масштабе по оси ординат.
5 . 2 .Разностные формулы
а) Разностные формулы для обыкновенных производных
Разностные формулы для приближенного вычисления производной подсказаны самим определением производной. Пусть значения функции в точках x i обозначены через y i :
y i = f (x i ), x i = a+ ih , i = 0, 1, ... , n; h =
Мы рассматриваем случай равномерного распределения точек на отрезке [a , b ]. Для приближенного вычисления производных в точках x i можно использовать следующие разностные формулы , или разностные производные .
Так как предел отношения (5.1) при h ® 0 равен правой производной в точке х i , то это отношение иногда называют правой разностной производной в точке x i .По аналогичной причине отношение (5.2) называют левой разностной производной в точке x i .Отношение (5.3) называют центральной разностной производной в точке x i .
Оценим погрешность разностных формул (5.1)–(5.3), предполагая, что функция f (x ) разлагается в ряд Тейлора в окрестности точки x i :
f (x ) = f (x i )+ . (5.4)
Полагая в (5.4) х = x i + h или х = х i - h , получим
Непосредственной подстановкой разложений (5.5) и (5.6) в формулу (5.10) можно получить зависимость между второй производной функции и разностной формулой для производной второго порядка .
Графическое дифференцирование начинают с построения графика функции по заданным значениям. При экспериментальном исследовании такой график получают с помощью самопишущих приборов. Далее проводят касательные к кривой в фиксированных положениях и вычисляют значения производной по тангенсу угла, образованного касательной с осью абсцисс.
На рис. 5.8, а изображена кривая полученная экспериментально на установке (рис. 5.6). Определение углового ускорения (искомой функции) проводят графическим дифференцированием по соотношению:
(5.19)
Тангенс угла наклона касательной к кривой в некоторой точке i представляют в виде отношения отрезков , где К – выбранный отрезок интегрирования (рис. 5.8, б )
После подстановки этого соотношения в соотношение (5.19) получают
где - ордината искового графика углового ускорения;
Масштаб искомого графика ; единицы СИ: = мм; = мм/(рад с -2).
График функции строят по найденным значениям ординат для ряда позиций. Точки на кривой соединяют от руки плавной линией, а затем обводят с помощью лекала.
Графическое дифференцирование рассмотренным методом касательных имеет относительно низкую точность. Более высокую точность получают при графическом дифференцировании методом хорд (рис. 5.8, в и г ).
 |
На заданной кривой отмечают ряд точек 1 ", 2 ", 3" , которые соединяют хордами, т.е. заменяют заданную кривую ломаной линией. Принимают следующее, допущение: угол наклона касательных в точках, расположенных посередине каждого участка кривой, равен углу наклона соответствующей хорды. Это допущение вносит некоторую погрешность, но она относится только к данной точке. Эти погрешности не суммируются, что обеспечивает приемлемую точность метода.
Остальные построения аналогичны ранее описанным при графическом дифференцировании методом касательных. Выбирают отрезок (мм); проводят лучи, наклоненные под углами 
до пересечения с осью ординат в точках 1
", 2
", 3
" ... , которые переносят на ординаты, проведенные в середине каждого из интервалов. Полученные точки 1
*, 2
*, 3
* являются точками искомой функции 
.
Масштабы по осям координат при этом методе построения связаны таким же соотношением (5.21), которое было выведено для случая графического дифференцирования методом касательных.
Дифференцирование функции f(x) , заданной (либо вычисленной) в виде массива чисел, выполняют методом численного дифференцирования с применением ЭВМ.
Чем меньше шаг в массиве чисел, тем точнее можно вычислить значение производной функции в этом интервале
Для решения многих инженерных задач часто требуется вычисление производных. Когда есть формула, описывающая процесс, сложностей никаких нет: берем формулу и вычисляем производную, как учили еще в школе, находим значения производной в разных точках, и всё. Сложность, наверное, только в этом и состоит, чтобы вспомнить, как вычислять производные. А как быть, если у нас есть только несколько сотен или тысяч строк с данными, а никакой формулы нет? Чаще всего именно так на практике и бывает. Предлагаю два способа.
Первый заключается в том, что мы наш набор точек аппроксимируем стандартной функцией Excel, то есть подбираем функцию, которая лучше всего ложится на наши точки (в Excel это линейная функция, логарифмическая, экспоненциальная, полиномиальная и степенная). Второй способ – численное дифференцирование, для которого нам нужно будет только умение вводить формулы.
Вспомним, что такое производная вообще:
Производной функции f (x) в точке x называется предел отношения приращения Δf функции в точке x к приращению Δx аргумента, когда последнее стремится к нулю:
Вот и воспользуемся этим знанием: будем просто брать для расчета производной очень маленькие значения приращения аргумента, т.е. Δx.
Для того, чтобы найти приближённое значение производной в нужных нам точках (а у нас точки – это различные значения степени деформации ε) можно поступить вот как. Посмотрим еще раз на определение производной и видим, что при использовании малых значений приращения аргумента Δε (то есть малых приращений степени деформации, которые регистрируются при испытаниях) можно заменить значение реальной производной в точке x 0 (f’(x 0)=dy/dx (x 0)) на отношение Δy/Δx=(f (x 0 + Δx) – f (x 0))/Δx.
То есть вот что получается:
f’(x 0) ≈(f (x 0 + Δx) – f (x 0))/Δx (1)
Для вычисления этой производной в каждой точке мы производим вычисления с использованием двух соседних точек: первая с координатой ε 0 по горизонтальной оси, а вторая с координатой x 0 + Δx, т.е. одна – производную в которой вычисляем и та, что поправее. Вычисленная таким образом производная называется разностной производной вправо (вперед) с шагом Δ x .
Можем поступить наоборот, взяв уже другие две соседние точки: x 0 — Δx и x 0 , т.е интересующую нас и ту, что левее. Получаем формулу для вычисления разностной производной влево (назад) с шагом — Δ x .
f’(x 0) ≈(f (x 0) – f (x 0 — Δx))/Δx (2)
Предыдущие формулы были «левые» и «правые», а есть еще одна формула, которая позволяет вычислять центральную разностною производную с шагом 2 Δx, и которая чаще других используется для численного дифференцирования:
f’(x 0) ≈(f (x 0 + Δx) – f (x 0 — Δx))/2Δx (3)
Для проверки формулы рассмотрим простой пример с известной функцией y=x 3 . Построим таблицу в Excel с двумя с столбцами: x и y, а затем построим график по имеющимся точкам.
Производная функции y=x 3 это y=3x 2 , график которой, т.е. параболу, мы и должны получить с использованием наших формул.
Попробуем вычислить значения центральной разностной производной в точках х. Для этого. В ячейке второй строки нашей таблицы забиваем нашу формулу (3), т.е. следующую формулу в Excel:
Теперь строим график с использованием уже имеющихся значений х и полученных значений центральной разностной производной:
А вот и наша красненькая парабола! Значит, формула работает!
Ну а теперь можем перейти к конкретной инженерной задаче, про которую говорили в начале статьи – к нахождению изменения dσ/dε с увеличением деформации. Первая производная кривой «напряжение-деформация» σ=f (ε) в зарубежной литературе называется «скорость упрочнения» (strain hardening rate),а в нашей – «коэффициент упрочнения». Итак, в результате испытаний мы имеем массив данных, которой состоит из двух столбцов: один — со значениями деформаций ε и другой – со значениями напряжений σ в МПа. Возьмем холодную деформацию стали 1035 или наша 40Г (см. таблицу аналогов сталей) при 20°С.
| C | Mn | P | S | Si | N |
| 0.36 | 0.69 | 0.025 | 0.032 | 0.27 | 0.004 |
Вот наша кривая в координатах «истинное напряжение — истинная деформация» σ-ε:
Действуем так же, как и в предыдущем примере и получаем вот такую кривую:
Это и есть изменение скорости упрочнения по ходу деформации. Что с ней делать, это уже отдельный вопрос.