I'm not advocating that we use the opcodes themselves. That is pretty much tantamount to numerical precision suicide as most of them have pretty arcane restrictions on their valid input ranges that vary by platform and expect you to play games with scaling or newton-raphson as well. You need a fair bit of wrapping around the available floating point opcodes. e.g. in glibc on x86 expm1 looks like: .text ENTRY(__expm1) fldl 4(%esp) // x fxam // Is NaN or +-Inf? fstsw %ax movb $0x45, %ch andb %ah, %ch cmpb $0x40, %ch je 3f // If +-0, jump. #ifdef PIC LOAD_PIC_REG (dx) #endif cmpb $0x05, %ch je 2f // If +-Inf, jump. fldt MO(l2e) // log2(e) : x fmulp // log2(e)*x fld %st // log2(e)*x : log2(e)*x frndint // int(log2(e)*x) : log2(e)*x fsubr %st, %st(1) // int(log2(e)*x) : fract(log2(e)*x) fxch // fract(log2(e)*x) : int(log2(e)*x) f2xm1 // 2^fract(log2(e)*x)-1 : int(log2(e)*x) fscale // 2^(log2(e)*x)-2^int(log2(e)*x) : int(log2(e)*x) fxch // int(log2(e)*x) : 2^(log2(e)*x)-2^int(log2(e)*x) fldl MO(one) // 1 : int(log2(e)*x) : 2^(log2(e)*x)-2^int(log2(e)*x) fscale // 2^int(log2(e)*x) : int(log2(e)*x) : 2^(log2(e)*x)-2^int(log2(e)*x) fsubrl MO(one) // 1-2^int(log2(e)*x) : int(log2(e)*x) : 2^(log2(e)*x)-2^int(log2(e)*x) fstp %st(1) // 1-2^int(log2(e)*x) : 2^(log2(e)*x)-2^int(log2(e)*x) fsubrp %st, %st(1) // 2^(log2(e)*x) ret 2: testl $0x200, %eax // Test sign. jz 3f // If positive, jump. fstp %st fldl MO(minus1) // Set result to -1.0. 3: ret I don't think we want to get into replicating that logic. =) -Edward On Thu, Apr 17, 2014 at 2:48 PM, Henning Thielemann < schlepptop@henning-thielemann.de> wrote:
Am 17.04.2014 19:15, schrieb Edward Kmett:
log1p and expm1 are C standard library functions that are important for
work with exponentials and logarithms.
I propose adding them to the Floating class where it is defined in GHC.Float.
No question, these functions are useful. But I think there should be two proposals:
1) Add log1pFloat, log1pDouble, expm1Float, expm1Double to GHC.Float 2) Extend Floating class with log1p and expm1 methods.
I think the first item is unproblematic since it is a simple addition. Since FPUs sometimes directly implement log1p and expm1 functions, I wonder whether GHC also should support the according machine instructions. E.g. x86 has F2XM1 and FYL2XP1 and good old MC68882 had FETOXM1 and FLOGNP1.
The second item means to alter the Floating class which affects all custom Floating instances. I think one should add default implementations. They don't have an numerical advantage but they save programmers from code breakage.
We do not have to export these from Prelude. My knee-jerk reaction is
that we probably shouldn't.
Not exporting them from Prelude still means to export them from GHC.Float, right? I mean, users must be able to implement these methods for custom types like extended precision floating point numbers as provided by libqd. But there should also be a non-GHC module that exports the full Floating class.