Very slow std :: pow () for databases very close to 1

I have a numerical code that solves an equation f(x) = 0in which I have to raise xto power p. I solve it using a bunch of things, but in the end I have Newton's method. The solution turns out to be equal x = 1, and therefore is the cause of my problems. When the iterative solution approaches 1, say x = 1 + 1e-13, the time required for the calculation std::pow(x, p)grows significantly, easily 100 times, which makes my code unusable.

The machine that works with this thing is AMD64 (Opteron 6172) on CentOS, and the team is simple y = std::pow(x, p);. This behavior appears on all my machines, all x64. As documented here , this is not only my problem (i.e., someone else got angry), appears only on x64 and only for xclose to 1.0. A similar thing happens with exp.

Solving this problem is vital to me. Does anyone know if there is a way around this slowness?

EDIT: John pointed out that this is due to denunciations. The question is how to fix this? C ++ code compiled with g++for use in GNU Octave. It seems that although I set it CXXFLAGSto turn on -mtune=nativeand -ffast-math, it does not help, and the code is just as slow.

PSEUDO-SOLUTION FOR NOW: For everyone who cares about this problem, the solutions below do not help me personally. I really need the usual speed std::pow(), but without the sluggishness around x = 1. The solution for me personally consists in the following hack:

inline double mpow(double x, double p) __attribute__ ((const));

inline double mpow(double x, double p)
{
    double y(x - 1.0);
    return (std::abs(y) > 1e-4) ? (std::pow(x, p)) : (1.0 + p * y * (1.0 + (p - 1.0) * y * (0.5 + (1.0 / 6.0) * (p - 2.0) * y)));
}

The score can be changed, but for -40 <p <40 the error is less than about 1e-11, which is good enough. The overhead is minimal from what I found, so this solves the problem.

+5
source share
3 answers

, a ** b == exp(log(a) * b) . , . : , .

, , ; exp(-2.4980018054066093e-15) , -2.4980018054066093e-15, , .

, :

sqrt(pow(a, b * 2))
pow(a * 2, b) / pow(2, b)
...

glibc: http://sourceware.org/bugzilla/show_bug.cgi?id=13932 - , , .

+9

64- Linux?

pow() FreeBSD.

Linux C (glibc) .

: http://entropymine.com/imageworsener/slowpow/

+1

. , - BFGS .

. , .

0

All Articles