Am Sun, 11 Sep 2016 15:00:12 +1000 schrieb Manu via Digitalmars-d <[email protected]>:
> On 9 September 2016 at 21:50, Stefan Koch via Digitalmars-d > <[email protected]> wrote: > > Hi, > > > > In short 80bit real are a real pain to support cross-platform. > > emulating them in software is prohibitively slow, and more importantly hard > > to get right. > > 64bit floating-point numbers are supported on more architectures and are > > much better supported. > > They are also trivial to use at ctfe. > > I vote for killing the 80bit handling at constant folding. > > > > Destroy! > > I just want CTFE '^^'. Please, can we have that? > It's impossible to CTFE any non-linear function. It's ridiculous, I > constantly want to generate a curve at compile time! I have experimented with a few iterative algorithms from around the web that are now in my module for "random stuff not in Phobos": /******************************************************************************* * * Computes the arcus tangens at compile time. * **************************************/ enum ctfeAtan(real x) in { assert(x == x && abs(x) != real.infinity); } body { if (abs(x) == 0) return x; // Reduce x to <0.5 for effective convergence of the Taylor series. x /= 1 + sqrt(1 + x * x); x /= 1 + sqrt(1 + x * x); // Sum up Taylor series to compute atan(). immutable xSqr = -x * x; real mul = x; real div = 1; real x_old; do { x_old = x; mul *= xSqr; div += 2; x += mul / div; } while (x !is x_old); // Compensate for the initial reduction by multiplying by 4. return 4 * x; } /******************************************************************************* * * Computes the arcs sinus at compile time. * **************************************/ enum ctfeAsin(real x) in { assert(x.isWithin(-1,+1)); } body { if (abs(x) == 0) return x; immutable div = 1 - x * x; return x / abs(x) * (div == 0 ? PI / 2 : ctfeAtan(sqrt(x * x / div))); } /******************************************************************************* * * Computes `x` to the power of `y` at compile-time. * * Params: * x = The base value. * y = The power. * * Source: * http://stackoverflow.com/a/7710097/4038614 * **************************************/ @safe @nogc pure nothrow real ctfePow(real x, real y) { if (y >= 1) { real temp = ctfePow( x, y / 2 ); return temp * temp; } real low = 0, high = 1; real sqr = sqrt( x ); real acc = sqr; real mid = high / 2; while (mid != y) { sqr = sqrt( sqr ); if (mid <= y) { low = mid; acc *= sqr; } else { high = mid; acc *= 1 / sqr; } mid = (low + high) / 2; } return acc; } /******************************************************************************* * * Computes the natural logarithm of `x`at compile time. * **************************************/ @safe @nogc pure nothrow FloatReg ctfeLog(FloatReg x) { if (x != x || x <= 0) return -FloatReg.nan; else if (x == 0) return -FloatReg.infinity; else if (x == FloatReg.infinity) return +FloatReg.infinity; uint m = 0; while (x <= ulong.max) { x *= 2; m++; } @safe @nogc pure nothrow static FloatReg agm(FloatReg x, FloatReg y) in { assert(x >= y); } body { real a, g; do { a = x; g = y; x = 0.5 * (a+g); y = sqrt(a*g); } while(x != a || y != g); return x; } return PI_2 / agm(1, 4/x) - m * LN2; } Especially the `log` function seemed like a good compromise between execution speed and accuracy. FloatReg is "the native float register type", but the way CTFE works it should be `real` anyways. The code is mostly not by me, but from StackOVerflow comments and Wikipedia. Most of it is common knowledge to every mathematician. -- Marco
