Вычисление функций

Материал из MachineLearning.

(Различия между версиями)
Перейти к: навигация, поиск
м
Текущая версия (07:17, 18 октября 2011) (править) (отменить)
(→Интерполяция степенным рядом)
 
(6 промежуточных версий не показаны.)
Строка 1: Строка 1:
== Введение ==
== Введение ==
-
 
+
Приближенное вычисление функций является важной практической задачей. Для большинства случаев практических вычислений бывает достаточно уже реализованных в программных системах функций, которые, кстати, также вычисляются приближенно. Однако, если необходимо вычислять эти функции с точностью, отличной от предлагаемой, или реализовать нестандартную функцию, то непременно возникает задача приближенного вычисления. Ниже описаны несколько подходов к решению данной задачи.
== Изложение метода. Анализ метода. Анализ ошибок ==
== Изложение метода. Анализ метода. Анализ ошибок ==
Строка 13: Строка 13:
К счастью, часто используемые в вычислениях элементарные функции разложимы в ряд Тейлора, хотя, как будет показано ниже, и для таких функций метод интерполяции ограниченным степенным рядом не дает желаемых результатов.
К счастью, часто используемые в вычислениях элементарные функции разложимы в ряд Тейлора, хотя, как будет показано ниже, и для таких функций метод интерполяции ограниченным степенным рядом не дает желаемых результатов.
-
Для разложения функции ''f(x)'' в ряд Тэйлора в окрестности точки ''x<sub>0</sub>''
+
Для разложения функции ''f(x)'' в ряд Тейлора в окрестности точки ''x<sub>0</sub>''
::<tex>f(x)=\sum_{n=0}^\infty f(x_0)^{(n)}\frac{(x-x_0)^n}{n!}</tex>
::<tex>f(x)=\sum_{n=0}^\infty f(x_0)^{(n)}\frac{(x-x_0)^n}{n!}</tex>
-
cпараведлива следующая оценка остаточного члена:
+
cправедлива следующая оценка остаточного члена:
::<tex>f(x)= f(x^*)^{(n+1)}\frac{(x-x_0)^{n+1}}{(n+1)!}, x\in (x_0,x)</tex>,
::<tex>f(x)= f(x^*)^{(n+1)}\frac{(x-x_0)^{n+1}}{(n+1)!}, x\in (x_0,x)</tex>,
Строка 22: Строка 22:
откуда можно найти число членов ряда, доставляющее необходимую точность приближения значения функции в точке ''х''.
откуда можно найти число членов ряда, доставляющее необходимую точность приближения значения функции в точке ''х''.
-
К недостаткам данного метода можно отнести наличие радиуса сходимости степенного ряда для разложения в окрестности точки ''x<sub>0</sub>''. Обычно радиус сходимости степенного ряда известен заранее, и с этой проблемой справляются путем замены переменной, внося ошибки округления. Но, гораздо более коварная ситуация возникает с всюду сходящимися рядами, так как зачастую они сходятся настолько медленно, что становятся непригодными для вычислений. Кроме того, абсолютная ошибка вычислений по указанной формуле распределена неравномерно и нарастает по мере удаления от точки ''x<sub>0</sub>''.
+
К недостаткам данного метода можно отнести наличие радиуса сходимости степенного ряда для разложения в окрестности точки ''x<sub>0</sub>''. Обычно радиус сходимости степенного ряда известен заранее, и с этой проблемой справляются путем замены переменной, внося ошибки округления. Но, гораздо более коварная ситуация возникает со всюду сходящимися рядами, так как зачастую они сходятся настолько медленно, что становятся непригодными для вычислений. Кроме того, абсолютная ошибка вычислений по указанной формуле распределена неравномерно и нарастает по мере удаления от точки ''x<sub>0</sub>''.
=== Полиномы Чебышева ===
=== Полиномы Чебышева ===
Строка 117: Строка 117:
== Числовой пример ==
== Числовой пример ==
-
=== Вычислние синуса на интервале [-1,1] ===
+
=== Вычисление синуса на интервале [-1,1] ===
Вычисление синуса произвольного угла можно свести к вычислению синуса угла в интервале <tex>[- \frac{\pi}2 , \frac{\pi}2 ]</tex> по формуле
Вычисление синуса произвольного угла можно свести к вычислению синуса угла в интервале <tex>[- \frac{\pi}2 , \frac{\pi}2 ]</tex> по формуле
Строка 143: Строка 143:
-
[[Изображение:Сhebyshev_error.png|Абсолютная ошибка вычислений при приближении рядом Чебышева ]]
+
[[Изображение:Teb_error.png‎|Абсолютная ошибка вычислений при приближении рядом Чебышева ]]
Метод экономизации можно применять не только для распространения максимальной ошибки ограничения по всему интервалу за счет уменьшения слагаемых в сумме, но и для уменьшения абсолютной величины ошибки, при неизменном количестве слагаемых.
Метод экономизации можно применять не только для распространения максимальной ошибки ограничения по всему интервалу за счет уменьшения слагаемых в сумме, но и для уменьшения абсолютной величины ошибки, при неизменном количестве слагаемых.
-
После корректировки коэффициенттов ряда Тэйлора с помощью метода экономизации получим:
+
После корректировки коэффициентов ряда Тэйлора с помощью метода экономизации получим:
::<tex>\sin\frac{\pi x}2 \approx 1.5707962x-0.64596332x^3+0.07968829x^5-0.004671873x^7+0.00015054436x^9</tex>
::<tex>\sin\frac{\pi x}2 \approx 1.5707962x-0.64596332x^3+0.07968829x^5-0.004671873x^7+0.00015054436x^9</tex>
Строка 154: Строка 154:
[[Изображение:Economization_error.png|Абсолютная ошибка вычислений при экономизайии ряда Тэйлора ]]
[[Изображение:Economization_error.png|Абсолютная ошибка вычислений при экономизайии ряда Тэйлора ]]
-
Economization_error.png
+
На следующем изображении представлены графики всех трех приближений и график функции <tex>\sin (\frac{\pi x}2)</tex>, значения приближенний увеличены в каждой точке пропорционально абсолютной ошибке. Имеет смысл сравнивать график для конкретного приближения с графиком синуса, но не между собой, так как порядок ошибок в каждом случае разный.
 +
 
 +
[[Изображение:Different_series.png|Графики приближений и график синуса]]
 +
 
 +
== См. также ==
 +
* [[Метод Монте-Карло]]
== Список Литературы ==
== Список Литературы ==
 +
*{{книга
 +
|автор = Д. Мак-Кракен, У. Дорн
 +
|заглавие = ЧИСЛЕННЫЕ МЕТОДЫ и ПРОГРАММИРОВАНИЕ на ФОРТРАНе
 +
|издание = 2-е изд
 +
|место = М.
 +
|издательство = Мир
 +
|год = 1977
 +
}}
 +
 +
*{{книга
 +
|автор = William H. Press, Saul A. Teukolsky, William T. Vetterling, Brian P. Flannery
 +
|заглавие = NUMERICAL RECIPES. The Art of Scientific Computing
 +
|издание = Third Edition
 +
|место = New York
 +
|издательство = Cambridge University Press
 +
|год = 2007
 +
|страниц = 1235
 +
|ссылка = http://www.nr.com
 +
|isbn = 0-521-88068-8
 +
}}
 +
 +
*{{книга
 +
|автор = F. B. HILDERBRAND
 +
|заглавие = Introduction to Numerical Analysis
 +
|издание = Second Edition
 +
|место = New York
 +
|издательство = Dover Publications, Inc.
 +
|год = 1974
 +
|страниц = 1235
 +
|ссылка = http://www.doverpublications.com
 +
|isbn = 0-486-65363-3
 +
}}
 +
 +
*{{книга
 +
|автор = Weisstein, Eric W.
 +
|часть = Chebyshev Polynomial of the First Kind
 +
|заглавие = MathWorld--A Wolfram Web Resource
 +
|ссылка = http://mathworld.wolfram.com/
 +
}}
 +
 +
*{{книга
 +
|автор = Scott A. Sarra
 +
|часть = Chebyshev Interpolation: An Interactive Tour
 +
|заглавие = Journal of Online Mathematics and Its Applications
 +
|год = 2006
 +
|том = 6
 +
|страницы = ID 1297
 +
|ссылка = http://mathdl.maa.org/
 +
}}
 +
[[Категория:Учебные задачи]]
[[Категория:Учебные задачи]]
 +
[[Категория:Численные методы]]

Текущая версия

Содержание

Введение

Приближенное вычисление функций является важной практической задачей. Для большинства случаев практических вычислений бывает достаточно уже реализованных в программных системах функций, которые, кстати, также вычисляются приближенно. Однако, если необходимо вычислять эти функции с точностью, отличной от предлагаемой, или реализовать нестандартную функцию, то непременно возникает задача приближенного вычисления. Ниже описаны несколько подходов к решению данной задачи.

Изложение метода. Анализ метода. Анализ ошибок

Интерполяция степенным рядом

Известно, что любую аналитическую функцию f(x) можно разложить в степенной ряд в окрестности некоторой точки x0.

f(x) = \sum_{n=0}^\infty a_n (x-x_0)^n

По известному разложению легко вычислить значение полинома, составленного из первых k членов степенного ряда в заданной точке. В общем случае вопрос о выборе числа k для достижения необходимой точности вычислений совсем не прост. Стандартная практика завершать вычисление суммы, когда модуль последнего добавленного слагаемого в определенное число раз меньше модуля накопленной суммы. К счастью, часто используемые в вычислениях элементарные функции разложимы в ряд Тейлора, хотя, как будет показано ниже, и для таких функций метод интерполяции ограниченным степенным рядом не дает желаемых результатов.

Для разложения функции f(x) в ряд Тейлора в окрестности точки x0

f(x)=\sum_{n=0}^\infty f(x_0)^{(n)}\frac{(x-x_0)^n}{n!}

cправедлива следующая оценка остаточного члена:

f(x)= f(x^*)^{(n+1)}\frac{(x-x_0)^{n+1}}{(n+1)!}, x\in (x_0,x),

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

К недостаткам данного метода можно отнести наличие радиуса сходимости степенного ряда для разложения в окрестности точки x0. Обычно радиус сходимости степенного ряда известен заранее, и с этой проблемой справляются путем замены переменной, внося ошибки округления. Но, гораздо более коварная ситуация возникает со всюду сходящимися рядами, так как зачастую они сходятся настолько медленно, что становятся непригодными для вычислений. Кроме того, абсолютная ошибка вычислений по указанной формуле распределена неравномерно и нарастает по мере удаления от точки x0.

Полиномы Чебышева

Полиномы Чебышева первого рода представляют собой ортогональную систему функций и определяются следующим образом:

T_n(x)= \cos{(n\,\arccos{\,x})}.

Например,

T_0(x)= 1,
T_1(x)= x,
T_2(x)= 2x^2-1,
T_3(x)= 4x^3-3x,
T_4(x)= 8x^4-8x^2+1,
T_5(x)= 16x^5-20x^3+5x.

Первые шесть полиномов Чебышева

Используя тригономитрические соотношения для косинусов суммы и разности, можно вывести рекуррентное соотношение для нахождения полиномов Чебышева:

T_{n+1}(x)= 2xT_n(x)-T_{n-1}(x).

Полином Tn(x) имеет на отрезке [-1,1] ровно n корней, расположенных в точках

x_k=\cos{\left(\frac{\pi(k+\frac12)}n \right)} \, k=0,1,\dots,n-1 .

Любую функцию f(x), определённую на отрезке [-1,1] можно приблизить следующей формулой:

f(x)\approx \sum_{j=0}^{N-1} c_j T_j(x) - \frac12 c_0 , где
c_j=\frac2N \sum_{k=0}^{N-1} f(x_k) T_j(x_k)=\frac2N \sum_{k=0}^{N-1} f \left(\cos{\left(\frac{\pi(k+\frac12)}N \right)} \right) \cos{\left(\frac{\pi j(k+\frac12)}N \right)} \, j=0,1,\dots,N-1 ,

которая в точности верна для всех x,являющихся корнями многочлена TN(x).

Преимущество разложения функции по полиномам Чебышева состоит в том, что при этом абсолютная ошибка вычислений знакопеременна и распределена более или менее равномерно по всему интервалу [-1,1]. Наилучшее приближение функции степенным рядом в том смысле, что максимальная ошибка при этом приближении минимальна, называется чебышевским приближением. Это приближение не совпадает с разложением по Чебышевским многочленам, но обычно его поиски, которые требуют больших вычислительных затрат, не оправдывают уменьшение ошибки. Т.е. приближение Чебышевскими многочленами практически совпадает с чебышевским приближением и гораздо более привлекательно в вычислительном плане.

Экономизация рядов

Метод экономизации рядов заключается в корректировке коэффициентов частичной суммы степенного ряда функции f(x). Он состоит из следующей последовательности действий:

  1. Вычислить необходимое число коэффициентов степенного ряда для приближения функции f(x) с требуемой точностью на отрезке [a,b];
  2. Сделать замену переменных для отображения интервала [a,b] в интервал [-1,1];
  3. Найти коэффициенты разложения полученного полинома по полиномам Чебышева;
  4. В полученном разложении взять первые k членов так, чтобы коэффициент при Tk+1 по абсолютной величине был меньше необходимой точности вычислений;
  5. Представить полученную сумму многочленом стандартного вида;
  6. Сделать обратную замену переменных.

Метод экономизации рядов позволяет распространить ошибку ограничения по всему интервалу, при этом уменьшив количество необходимых для вычисления слагаемых.

Вычисление рядов

Независимо от того, каким способом было получено разложение функции в степенной ряд — по формуле Тейлора, Чебышева или методом экономизации, — всегда бывает необходимо вычислить значение полинома вида

p(x)=a_0+a_1x+\dots+a_{n-1}x^{n-1}+a_nx^n .

Метод вычисления полинома оказывает большое влияние на распространение ошибок вычисления (ошибки исходной информации и ошибки округления).

Число умножений при вычислении по данной формуле составляет \frac{n(n+1)}2. Этот полином можно переписать в несколько ином виде, после чего его можно вычислить не только быстрее, но и во многих случаях с большей точностью.

Правило Горнера

Правило Горнера, широко известный метод вычисления значений полинома, записывается в виде:

p(x)=a_0+x(a_1+x(a_2+\dots+x(a_{n-1}+x(a_n))\,\dots)).

Вычисления по данной формуле производятся в соответствии с заданным порядком действий. Для вычисления полинома по правилу Горнера требуется n сложений и n умножений. Во многих практических расчетах применение правила Горнера не только экономит машинное время, но и повышает машинное время за счет уменьшения верхнего предела ошибки округления.

Параллельное вычисление многочлена

Полином степени n может быть вычислен за O(\log_2n) параллельных шагов. Процесс вычисления состоит в следующем:

  1. На вход алгоритму подается вектор коэффициентов исходного многочлена, при необходимости дополненный нулем, для того, чтобы вектор имел четную длину; k=0.
  2. Соседние элементы вектора попарно суммируются, причем второй член суммы умножается на x2k, где k - это номер итерации. Таким образом, на выходе мы получим вектор в два раза меньшей длины. Полученный вектор также дополняется нулем, в случае необходимости; k=k+1.
  3. Шаг 2 повторяется до тех пор, пока длина вектора на выходе не станет единичной.

Данный метод, в общем случае, имеет меньшую ошибку округления, чем метод последовательного вычисления.

Суммирование выражений с членами, связанными рекуррентными соотношениями

Данный метод предложен Clenshaw для вычисления сумм вида

f(x)=\sum_{k=0}^N c_k F_k(x),

где функции Fk(x) связаны рекуррентным соотношением

F_{n+1}(x)=\alpha (n,x) F_n(x) + \beta (n,x) F_{n-1}(x),

для некоторых функций α(n,x) и β(n,x).

Определим вспомогательные переменные yk следующими соотношениями:

y_{N+2}=y_{N+1}=0
y_k=\alpha (k,x) y_{k+1} + \beta (k+1,x)y_{k+2} + c_k \, (k=N,N-1,\dots,1).

Если разрешить эту систему относительно ck, и подставить полученные выражения в исходную формулу для f(x), то после сокращения получим:

f(x)=\beta (1,x) F_0(x)y_2 + F_1(x)y_1+F_0(x)c_0


Числовой пример

Вычисление синуса на интервале [-1,1]

Вычисление синуса произвольного угла можно свести к вычислению синуса угла в интервале [- \frac{\pi}2 , \frac{\pi}2 ] по формуле

sin(n\pi+y)=(-1)^n \sin y,

делая замену переменных y=\frac{\pi x}2 мы сведем диапазон значний аргумента к интервалу [-1,1].

При совершении данных операций возникает ошибка округления, но она оказывается существенно меньше ошибки ограничения. Проиллюстрируем это на примере.

Возьмем в качестве приближения функции \sin(\frac{\pi x}2) первые пять членов ряда Тэйлора:

\sin\frac{\pi x}2 \approx 1.5707963x-0.64596410x^3+0.079692626x^5-0.00468117541x^7+0.00016044118x^9

Ошибка ограничения меньше первого отброшенного члена:

|e_T| \le \frac 1{11!}( \frac{\pi x}2)^11 \le 3.6*10^{-6}

Абсолютная ошибка вычислений при ограничении ряда Тейлора

На рисунке изображен график абсолютной ошибки при вычислении по указанной формуле. Ошибка составляется из ошибки ограничения и ошибки округления, но в данном случае ошибка округления мала, по сравнению с ошибкой ограничения. Стоит отметить, что хотя ошибка практичски отсутствует вблизи нуля, она становится довольно большой при подходе к границам интервала.

При приближении той же функции полиномами Чебышева до 9 порядка получаем следующую формулу:

\sin\frac{\pi x}2 \approx 1.5707963x-0.64596336x^3+0.079688475x^5-0.0046722203x^7+0.00015081716x^9


Абсолютная ошибка вычислений при приближении рядом Чебышева


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

После корректировки коэффициентов ряда Тэйлора с помощью метода экономизации получим:

\sin\frac{\pi x}2 \approx 1.5707962x-0.64596332x^3+0.07968829x^5-0.004671873x^7+0.00015054436x^9

Абсолютная ошибка вычислений при экономизайии ряда Тэйлора

На следующем изображении представлены графики всех трех приближений и график функции \sin (\frac{\pi x}2), значения приближенний увеличены в каждой точке пропорционально абсолютной ошибке. Имеет смысл сравнивать график для конкретного приближения с графиком синуса, но не между собой, так как порядок ошибок в каждом случае разный.

Графики приближений и график синуса

См. также

Список Литературы

  • Д. Мак-Кракен, У. Дорн ЧИСЛЕННЫЕ МЕТОДЫ и ПРОГРАММИРОВАНИЕ на ФОРТРАНе. — 2-е изд. — М.: Мир, 1977.
Источник — «http://machinelearning.ru/wiki/index.php?title=%D0%92%D1%8B%D1%87%D0%B8%D1%81%D0%BB%D0%B5%D0%BD%D0%B8%D0%B5_%D1%84%D1%83%D0%BD%D0%BA%D1%86%D0%B8%D0%B9»
Личные инструменты