Быстрое вычисление exp(-x^2)

gala05

в задаче требуется считать очень много таких вот экспонент
нужна приближенная формула
диапазон изменения x - от -3 до 3

vovatroff

в задаче требуется считать очень много таких вот экспонент
нужна приближенная формула
диапазон изменения x - от -3 до 3
Юз фортран

zuzaka

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

spiritmc

Ряд Тейлора не устраивает?
Он же хорошо обрывается.
---
...Я работаю антинаучным аферистом...

gala05

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

spiritmc

У тебя какие требования по точности?
Равномерно приблизить с абс./отн. погрешностью <число>?
Поточнее в нуле?
Разложи в нуле, единице ... и восьмёрке.
Если придумаешь быстрый извращённый способ --- и в девятке.
И не забывай сравнивать с вычислением экспоненты.
---
...Я работаю антинаучным аферистом...

spiritmc

Знаю быстрый извращённый способ задействовать разложение в девятке для x87.
Недостаток --- губит вмешательство внешнего кода.
---
...Я работаю антинаучным аферистом...

gala05

равномерно на отрезке -3..3
p.s. я просил совета, а не непонятно по какому поводу язвительных слов.

afony

Так какая точность нужна?

afony

Если точности примерно в 0,002 хватит, то возьми рациональную функцию
(1+1/9*x^4-1/2*x^2-1/72*x^6+1/1008*x^8-1/30240*x^10)/(1+1/9*x^4+1/2*x^2+1/72*x^6+1/1008*x^8+1/30240*x^10)

zuzaka

ты думаешь, это будет заметно быстрее, чем e^(x^2)?

vovatroff

фортран не пойдет
Можно тогда уточнить постановку задачи: вы это вручную собираетесь
делать, на калькуляторе, или как? Почему не воспользоваться встроенной
экспонентой языка программирования, не Фортрана, так какого-нибудь еще?
Встроенные функции обычно быстрее все вычисляют, чем любые их эмуляции
доморощенными средствами. Или я не прав?

afony

А как ты собираешься считать e^{-x^2}? Если через ряд Тейлора, то быстрее.

kachokslava


сиреневая - Exp(-x^2 синяя - твоя дробь. жёлтая - разница между ними.

#include <stdio.h>
#include <math.h>
#include <stdlib.h>
double C1=0.5, C2=1/9.0, C3=1/72.0,C4=1/1008.0,C5=1/30240.0;
double f(double x)
{
double p1,p2,x2;
x2=x*x;
p1-C5*x2)+C4)*x2-C3*x2)+C2)*x2-C1)*x2+1;
p2+C5*x2)+C4)*x2+C3*x2)+C2)*x2+C1)*x2+1;
return p1/p2;
// (1-1/2*x^2+1/9*x^4-1/72*x^6+1/1008*x^8-1/30240*x^10)/
// (1+1/2*x^2+1/9*x^4+1/72*x^6+1/1008*x^8+1/30240*x^10)
}
double ff(double x)
{
return exp(-x*x);
}
int main
{
double x;
for(x=-3;x<3;x+=0.1)
printf("%lf\t%lf\t%lf\t%le\n",x, f(x ff(xf(x)-ff(x;
return 0;
}

т.е. она реально годится только на отрезке -1.5..1.5
ща ещё время померим

kachokslava

точность можно ещё такими ухищрениями повысить:

double f(double x)
{
double p1,p2,x2;
int f=0;
x2=x*x;
f=x2;
x2=x2-floor(x2);

p1-C5*x2)+C4)*x2-C3*x2)+C2)*x2-C1)*x2+1;
p2+C5*x2)+C4)*x2+C3*x2)+C2)*x2+C1)*x2+1;
p1/=p2;
for(;f>0;f--) p1/=E;
return p1;
}

на отрезке -3..3 погрешность достигает своего максимума в 0.9 и примерно равна 10^-5

kachokslava

exp(..) считается в сопроцессоре.
временные замеры показали следующую ситуацию:

1250 1.772441e+006
531 1.772415e+006


err=0;
t=clock;
for(i=0;i<10000;i++)
for(x=-3;x<3;x+=0.01)
err+=f(x);
t=clock-t;
printf("%d %le\n",t,err);
err=0;
t=clock;
for(i=0;i<10000;i++)
for(x=-3;x<3;x+=0.01)
err+=ff(x);
t=clock-t;
printf("%d %le\n",t,err);

f(x) - через дробь, ff(x) - через exp.
результат неутешительный: через Exp быстрее больше, чем в 2 раза.
AthlonXP 2600+ 2.08Ghz

afony

Неутешительный для кого? Я не ставил себе задачу оптимизации вычислений. Предложенная функция является дробью Паде - рациональной дробью наилучшего локального в окрестности 0 приближения. Естественно, что чем ближе к нулю, тем приближение лучше. Для отрезка [-3,3] она дает точность примерно в 0.002, а уж годится она или нет - зависит от требуемой степени точности. Можно вообще заранее посчитать e^{-x^2} во всех нужных точках и составить таблицу значений - тогда больше ничего считать не нужно. Все зависит от конкретной задачи.

gala05

спасибо за формулу
но задача была конкретной
нужно БЫСТРО вычислять exp(-x^2) в диапазоне [-3..3]

zuzaka

не, задача не была конкретной. Ты не указал точность. Полагаю, для точности в 0.1% табулирование еще будет довольно эффективно работать.

Sanych

Есть ещё такой вариант:
а)вычисляем -x^2
б)раскладываем в сумму табличных значений с нужной точностью
в)перемножнаем известные по таблице величины e^(слагаемые)
В зависимости от необходимой точности и прочих параметров, размер таблицы, число слагаемых и скорость могут оказаться разными. Можно ли это написать, чтобы было быстрее встроенной команды сопроцессора, я сомневаюсь .

afony

За формулу пожалуйста.
Но тогда тебе правильно ответил . Используй фортран - для счетных задач самый оптимальный вариант.

choconasty

Интерполяционный многочлен Лагранжа Ln(x) для f(x) := exp(-x) на отрезке [0, 9] степени (n — 1) по n узлам Чебышева.
m = 1, . . , n;
xm
:= (9/2)*( 1 + cos{pi*(2m — 1)/(2n)} )
Ln(x)
:= Sumk=1,..,n f(xk)Пs=1,..n,s≠k{(x — xs)/(xk — xs)}
Равномерная оценка модуля погрешности En:
En := 2 pow(9/4; n) / n!
В частности,
E11 = pow(3.,18.)/(pow(2.,29.)*25.*7.*11.) = +0.000374871
E10 = pow(3.,16.)/(pow(2.,27.)*25.*7.) = +0.00183270
E9 = pow(3.,14.)/(pow(2.,24.)*35.) = +0.00814535
E8 = pow(3.,14.)/(pow(2.,22.)*35.) = +0.0325814
Реальная погрешность заметно меньше.

vovatroff

Да человек, похоже, сам не знает, чего хочет. Такому помочь сложнее всего.

gala05

это просто ты думать или читать не хочешь
персонально для тебя развернутая формулировка: в компиляторах, в стандартных библиотеках ко всем средам программирования есть функция exp. мне нужно найти другой способ расчета экспоненты. чтобы миллион грубых вычислений экспоненты работало в несколько раз быстрее, чем встроенная функция exp. про точность и рабочий диапазон я писал
предложенные способы меня не устраивают (они все были внимательно прочитаны и перепробованы время выигрыша - процентов 5-10. это хорошо, но нужно больше

kachokslava

"встроенная" функция вычисляется на сопроцессоре. сомневаюсь, что есть методы посчитать быстрее.
может, не экспоненту надо ускорять в твоей программе?
кстати, табличный метод работает 60% времени вычисления exp(x при условии что приближение делается кусочно-постоянной функцией, а не кусочно-линейной

spiritmc

Если ты посмотришь в исходники библиотечной функции,
то увидишь, что там есть, где выиграть время.
---
"Vyroba umelych lidi, slecno, je tovarni tajemstvi."

vovatroff

предложенные способы меня не устраивают (они все были внимательно прочитаны и перепробованы время выигрыша - процентов 5-10. это хорошо, но нужно больше
Боюсь, что ты и вправду хочешь невозможного. Смирись.
Оставить комментарий
Имя или ник:
Комментарий: